ENH: autoLayerDriver: better layering information message
[OpenFOAM-2.0.x.git] / src / sampling / sampledSet / sampledSets / sampledSetsTemplates.C
blobab96643d7df0e49899c5a83699dd030f27bcbae9
1 /*---------------------------------------------------------------------------*\
2   =========                 |
3   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
4    \\    /   O peration     |
5     \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
6      \\/     M anipulation  |
7 -------------------------------------------------------------------------------
8 License
9     This file is part of OpenFOAM.
11     OpenFOAM is free software: you can redistribute it and/or modify it
12     under the terms of the GNU General Public License as published by
13     the Free Software Foundation, either version 3 of the License, or
14     (at your option) any later version.
16     OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17     ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18     FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
19     for more details.
21     You should have received a copy of the GNU General Public License
22     along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
24 \*---------------------------------------------------------------------------*/
26 #include "sampledSets.H"
27 #include "volFields.H"
28 #include "ListListOps.H"
30 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
32 template <class Type>
33 Foam::sampledSets::volFieldSampler<Type>::volFieldSampler
35     const word& interpolationScheme,
36     const GeometricField<Type, fvPatchField, volMesh>& field,
37     const PtrList<sampledSet>& samplers
40     List<Field<Type> >(samplers.size()),
41     name_(field.name())
43     autoPtr<interpolation<Type> > interpolator
44     (
45         interpolation<Type>::New(interpolationScheme, field)
46     );
48     forAll(samplers, setI)
49     {
50         Field<Type>& values = this->operator[](setI);
51         const sampledSet& samples = samplers[setI];
53         values.setSize(samples.size());
54         forAll(samples, sampleI)
55         {
56             const point& samplePt = samples[sampleI];
57             label cellI = samples.cells()[sampleI];
58             label faceI = samples.faces()[sampleI];
60             if (cellI == -1 && faceI == -1)
61             {
62                 // Special condition for illegal sampling points
63                 values[sampleI] = pTraits<Type>::max;
64             }
65             else
66             {
67                 values[sampleI] = interpolator().interpolate
68                 (
69                     samplePt,
70                     cellI,
71                     faceI
72                 );
73             }
74         }
75     }
79 template <class Type>
80 Foam::sampledSets::volFieldSampler<Type>::volFieldSampler
82     const GeometricField<Type, fvPatchField, volMesh>& field,
83     const PtrList<sampledSet>& samplers
86     List<Field<Type> >(samplers.size()),
87     name_(field.name())
89     forAll(samplers, setI)
90     {
91         Field<Type>& values = this->operator[](setI);
92         const sampledSet& samples = samplers[setI];
94         values.setSize(samples.size());
95         forAll(samples, sampleI)
96         {
97             label cellI = samples.cells()[sampleI];
99             if (cellI ==-1)
100             {
101                 values[sampleI] = pTraits<Type>::max;
102             }
103             else
104             {
105                 values[sampleI] = field[cellI];
106             }
107         }
108     }
112 template <class Type>
113 Foam::sampledSets::volFieldSampler<Type>::volFieldSampler
115     const List<Field<Type> >& values,
116     const word& name
119     List<Field<Type> >(values),
120     name_(name)
124 template<class Type>
125 void Foam::sampledSets::writeSampleFile
127     const coordSet& masterSampleSet,
128     const PtrList<volFieldSampler<Type> >& masterFields,
129     const label setI,
130     const fileName& timeDir,
131     const writer<Type>& formatter
134     wordList valueSetNames(masterFields.size());
135     List<const Field<Type>*> valueSets(masterFields.size());
137     forAll(masterFields, fieldi)
138     {
139         valueSetNames[fieldi] = masterFields[fieldi].name();
140         valueSets[fieldi] = &masterFields[fieldi][setI];
141     }
143     fileName fName
144     (
145         timeDir/formatter.getFileName(masterSampleSet, valueSetNames)
146     );
148     OFstream ofs(fName);
149     if (ofs.opened())
150     {
151         formatter.write
152         (
153             masterSampleSet,
154             valueSetNames,
155             valueSets,
156             ofs
157         );
158     }
159     else
160     {
161         WarningIn
162         (
163             "void Foam::sampledSets::writeSampleFile"
164             "("
165                 "const coordSet&, "
166                 "const PtrList<volFieldSampler<Type> >&, "
167                 "const label, "
168                 "const fileName&, "
169                 "const writer<Type>&"
170             ")"
171         )   << "File " << ofs.name() << " could not be opened. "
172             << "No data will be written" << endl;
173     }
177 template<class T>
178 void Foam::sampledSets::combineSampledValues
180     const PtrList<volFieldSampler<T> >& sampledFields,
181     const labelListList& indexSets,
182     PtrList<volFieldSampler<T> >& masterFields
185     forAll(sampledFields, fieldi)
186     {
187         List<Field<T> > masterValues(indexSets.size());
189         forAll(indexSets, setI)
190         {
191             // Collect data from all processors
192             List<Field<T> > gatheredData(Pstream::nProcs());
193             gatheredData[Pstream::myProcNo()] = sampledFields[fieldi][setI];
194             Pstream::gatherList(gatheredData);
196             if (Pstream::master())
197             {
198                 Field<T> allData
199                 (
200                     ListListOps::combine<Field<T> >
201                     (
202                         gatheredData,
203                         Foam::accessOp<Field<T> >()
204                     )
205                 );
207                 masterValues[setI] = UIndirectList<T>
208                 (
209                     allData,
210                     indexSets[setI]
211                 )();
212             }
213         }
215         masterFields.set
216         (
217             fieldi,
218             new volFieldSampler<T>
219             (
220                 masterValues,
221                 sampledFields[fieldi].name()
222             )
223         );
224     }
228 template<class Type>
229 void Foam::sampledSets::sampleAndWrite
231     fieldGroup<Type>& fields
234     if (fields.size())
235     {
236         bool interpolate = interpolationScheme_ != "cell";
238         // Create or use existing writer
239         if (fields.formatter.empty())
240         {
241             fields = writeFormat_;
242         }
244         // Storage for interpolated values
245         PtrList<volFieldSampler<Type> > sampledFields(fields.size());
247         forAll(fields, fieldi)
248         {
249             if (Pstream::master() && verbose_)
250             {
251                 Pout<< "sampledSets::sampleAndWrite: "
252                     << fields[fieldi] << endl;
253             }
255             if (loadFromFiles_)
256             {
257                 GeometricField<Type, fvPatchField, volMesh> vf
258                 (
259                     IOobject
260                     (
261                         fields[fieldi],
262                         mesh_.time().timeName(),
263                         mesh_,
264                         IOobject::MUST_READ,
265                         IOobject::NO_WRITE,
266                         false
267                     ),
268                     mesh_
269                 );
271                 if (interpolate)
272                 {
273                     sampledFields.set
274                     (
275                         fieldi,
276                         new volFieldSampler<Type>
277                         (
278                             interpolationScheme_,
279                             vf,
280                             *this
281                         )
282                     );
283                 }
284                 else
285                 {
286                     sampledFields.set
287                     (
288                         fieldi,
289                         new volFieldSampler<Type>(vf, *this)
290                     );
291                 }
292             }
293             else
294             {
295                 if (interpolate)
296                 {
297                     sampledFields.set
298                     (
299                         fieldi,
300                         new volFieldSampler<Type>
301                         (
302                             interpolationScheme_,
303                             mesh_.lookupObject
304                             <GeometricField<Type, fvPatchField, volMesh> >
305                             (fields[fieldi]),
306                             *this
307                         )
308                     );
309                 }
310                 else
311                 {
312                     sampledFields.set
313                     (
314                         fieldi,
315                         new volFieldSampler<Type>
316                         (
317                             mesh_.lookupObject
318                             <GeometricField<Type, fvPatchField, volMesh> >
319                             (fields[fieldi]),
320                             *this
321                         )
322                     );
323                 }
324             }
325         }
327         // Combine sampled fields from processors.
328         // Note: only master results are valid
330         PtrList<volFieldSampler<Type> > masterFields(sampledFields.size());
331         combineSampledValues(sampledFields, indexSets_, masterFields);
333         if (Pstream::master())
334         {
335             forAll(masterSampledSets_, setI)
336             {
337                 writeSampleFile
338                 (
339                     masterSampledSets_[setI],
340                     masterFields,
341                     setI,
342                     outputPath_/mesh_.time().timeName(),
343                     fields.formatter()
344                 );
345             }
346         }
347     }
351 // ************************************************************************* //