1 /*---------------------------------------------------------------------------*\
3 \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
5 \\ / A nd | Copyright (C) 2011 OpenFOAM Foundation
7 -------------------------------------------------------------------------------
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
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 * * * * * * * * * * * //
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()),
43 autoPtr<interpolation<Type> > interpolator
45 interpolation<Type>::New(interpolationScheme, field)
48 forAll(samplers, setI)
50 Field<Type>& values = this->operator[](setI);
51 const sampledSet& samples = samplers[setI];
53 values.setSize(samples.size());
54 forAll(samples, sampleI)
56 const point& samplePt = samples[sampleI];
57 label cellI = samples.cells()[sampleI];
58 label faceI = samples.faces()[sampleI];
60 if (cellI == -1 && faceI == -1)
62 // Special condition for illegal sampling points
63 values[sampleI] = pTraits<Type>::max;
67 values[sampleI] = interpolator().interpolate
80 Foam::sampledSets::volFieldSampler<Type>::volFieldSampler
82 const GeometricField<Type, fvPatchField, volMesh>& field,
83 const PtrList<sampledSet>& samplers
86 List<Field<Type> >(samplers.size()),
89 forAll(samplers, setI)
91 Field<Type>& values = this->operator[](setI);
92 const sampledSet& samples = samplers[setI];
94 values.setSize(samples.size());
95 forAll(samples, sampleI)
97 label cellI = samples.cells()[sampleI];
101 values[sampleI] = pTraits<Type>::max;
105 values[sampleI] = field[cellI];
112 template <class Type>
113 Foam::sampledSets::volFieldSampler<Type>::volFieldSampler
115 const List<Field<Type> >& values,
119 List<Field<Type> >(values),
125 void Foam::sampledSets::writeSampleFile
127 const coordSet& masterSampleSet,
128 const PtrList<volFieldSampler<Type> >& masterFields,
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)
139 valueSetNames[fieldi] = masterFields[fieldi].name();
140 valueSets[fieldi] = &masterFields[fieldi][setI];
145 timeDir/formatter.getFileName(masterSampleSet, valueSetNames)
163 "void Foam::sampledSets::writeSampleFile"
166 "const PtrList<volFieldSampler<Type> >&, "
169 "const writer<Type>&"
171 ) << "File " << ofs.name() << " could not be opened. "
172 << "No data will be written" << endl;
178 void Foam::sampledSets::combineSampledValues
180 const PtrList<volFieldSampler<T> >& sampledFields,
181 const labelListList& indexSets,
182 PtrList<volFieldSampler<T> >& masterFields
185 forAll(sampledFields, fieldi)
187 List<Field<T> > masterValues(indexSets.size());
189 forAll(indexSets, setI)
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())
200 ListListOps::combine<Field<T> >
203 Foam::accessOp<Field<T> >()
207 masterValues[setI] = UIndirectList<T>
218 new volFieldSampler<T>
221 sampledFields[fieldi].name()
229 void Foam::sampledSets::sampleAndWrite
231 fieldGroup<Type>& fields
236 bool interpolate = interpolationScheme_ != "cell";
238 // Create or use existing writer
239 if (fields.formatter.empty())
241 fields = writeFormat_;
244 // Storage for interpolated values
245 PtrList<volFieldSampler<Type> > sampledFields(fields.size());
247 forAll(fields, fieldi)
249 if (Pstream::master() && verbose_)
251 Pout<< "sampledSets::sampleAndWrite: "
252 << fields[fieldi] << endl;
257 GeometricField<Type, fvPatchField, volMesh> vf
262 mesh_.time().timeName(),
276 new volFieldSampler<Type>
278 interpolationScheme_,
289 new volFieldSampler<Type>(vf, *this)
300 new volFieldSampler<Type>
302 interpolationScheme_,
304 <GeometricField<Type, fvPatchField, volMesh> >
315 new volFieldSampler<Type>
318 <GeometricField<Type, fvPatchField, volMesh> >
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())
335 forAll(masterSampledSets_, setI)
339 masterSampledSets_[setI],
342 outputPath_/mesh_.time().timeName(),
351 // ************************************************************************* //