cyclicACMIFvPatchField.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | Copyright (C) 2013-2014 OpenFOAM Foundation
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
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.
15 
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.
20 
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/>.
23 
24 \*---------------------------------------------------------------------------*/
25 
26 #include "cyclicACMIFvPatchField.H"
27 #include "transformField.H"
28 
29 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
30 
31 template<class Type>
33 (
34  const fvPatch& p,
36 )
37 :
40  cyclicACMIPatch_(refCast<const cyclicACMIFvPatch>(p))
41 {}
42 
43 
44 template<class Type>
46 (
48  const fvPatch& p,
50  const fvPatchFieldMapper& mapper
51 )
52 :
54  coupledFvPatchField<Type>(ptf, p, iF, mapper),
55  cyclicACMIPatch_(refCast<const cyclicACMIFvPatch>(p))
56 {
57  if (!isA<cyclicACMIFvPatch>(this->patch()))
58  {
60  (
61  "cyclicACMIFvPatchField<Type>::cyclicACMIFvPatchField"
62  "("
63  "const cyclicACMIFvPatchField<Type>& ,"
64  "const fvPatch&, "
65  "const DimensionedField<Type, volMesh>&, "
66  "const fvPatchFieldMapper&"
67  ")"
68  ) << " patch type '" << p.type()
69  << "' not constraint type '" << typeName << "'"
70  << "\n for patch " << p.name()
71  << " of field " << this->dimensionedInternalField().name()
72  << " in file " << this->dimensionedInternalField().objectPath()
73  << exit(FatalIOError);
74  }
75 }
76 
77 
78 template<class Type>
80 (
81  const fvPatch& p,
83  const dictionary& dict
84 )
85 :
88  cyclicACMIPatch_(refCast<const cyclicACMIFvPatch>(p))
89 {
90  if (!isA<cyclicACMIFvPatch>(p))
91  {
93  (
94  "cyclicACMIFvPatchField<Type>::cyclicACMIFvPatchField"
95  "("
96  "const fvPatch&, "
97  "const DimensionedField<Type, volMesh>&, "
98  "const dictionary&"
99  ")",
100  dict
101  ) << " patch type '" << p.type()
102  << "' not constraint type '" << typeName << "'"
103  << "\n for patch " << p.name()
104  << " of field " << this->dimensionedInternalField().name()
105  << " in file " << this->dimensionedInternalField().objectPath()
106  << exit(FatalIOError);
107  }
108 
109  if (!dict.found("value") && this->coupled())
110  {
111  this->evaluate(Pstream::blocking);
112  }
113 }
114 
115 
116 template<class Type>
118 (
120 )
121 :
124  cyclicACMIPatch_(ptf.cyclicACMIPatch_)
125 {}
126 
127 
128 template<class Type>
130 (
131  const cyclicACMIFvPatchField<Type>& ptf,
133 )
134 :
136  coupledFvPatchField<Type>(ptf, iF),
137  cyclicACMIPatch_(ptf.cyclicACMIPatch_)
138 {}
139 
140 
141 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
142 
143 template<class Type>
145 {
146  return cyclicACMIPatch_.coupled();
147 }
148 
149 
150 template<class Type>
153 {
154  const Field<Type>& iField = this->internalField();
155  const labelUList& nbrFaceCellsCoupled =
156  cyclicACMIPatch_.cyclicACMIPatch().neighbPatch().faceCells();
157  const labelUList& faceCellsNonOverlap =
158  cyclicACMIPatch_.cyclicACMIPatch().nonOverlapPatch().faceCells();
159 
160  Field<Type> pnfCoupled(iField, nbrFaceCellsCoupled);
161  Field<Type> pfNonOverlap(iField, faceCellsNonOverlap);
162 
163  tmp<Field<Type> > tpnf
164  (
165  new Field<Type>
166  (
167  cyclicACMIPatch_.interpolate
168  (
169  pnfCoupled,
170  pfNonOverlap
171  )
172  )
173  );
174 
175  if (doTransform())
176  {
177  tpnf() = transform(forwardT(), tpnf());
178  }
179 
180  return tpnf;
181 }
182 
183 
184 template<class Type>
187 {
190  (
191  this->internalField()
192  );
193 
194  return refCast<const cyclicACMIFvPatchField<Type> >
195  (
196  fld.boundaryField()[cyclicACMIPatch_.neighbPatchID()]
197  );
198 }
199 
200 
201 template<class Type>
204 {
207  (
208  this->internalField()
209  );
210 
211  return fld.boundaryField()[cyclicACMIPatch_.nonOverlapPatchID()];
212 }
213 
214 
215 template<class Type>
217 (
218  scalarField& result,
219  const scalarField& psiInternal,
220  const scalarField& coeffs,
221  const direction cmpt,
222  const Pstream::commsTypes
223 ) const
224 {
225  // note: only applying coupled contribution
226 
227  const labelUList& nbrFaceCellsCoupled =
228  cyclicACMIPatch_.cyclicACMIPatch().neighbPatch().faceCells();
229 
230  scalarField pnf(psiInternal, nbrFaceCellsCoupled);
231 
232  // Transform according to the transformation tensors
233  transformCoupleField(pnf, cmpt);
234 
235  const labelUList& faceCells = cyclicACMIPatch_.faceCells();
236 
237  pnf = cyclicACMIPatch_.interpolate(pnf);
238 
239  forAll(faceCells, elemI)
240  {
241  result[faceCells[elemI]] -= coeffs[elemI]*pnf[elemI];
242  }
243 }
244 
245 
246 template<class Type>
248 (
249  Field<Type>& result,
250  const Field<Type>& psiInternal,
251  const scalarField& coeffs,
252  const Pstream::commsTypes
253 ) const
254 {
255  // note: only applying coupled contribution
256 
257  const labelUList& nbrFaceCellsCoupled =
258  cyclicACMIPatch_.cyclicACMIPatch().neighbPatch().faceCells();
259 
260  Field<Type> pnf(psiInternal, nbrFaceCellsCoupled);
261 
262  // Transform according to the transformation tensors
263  transformCoupleField(pnf);
264 
265  const labelUList& faceCells = cyclicACMIPatch_.faceCells();
266 
267  pnf = cyclicACMIPatch_.interpolate(pnf);
268 
269  forAll(faceCells, elemI)
270  {
271  result[faceCells[elemI]] -= coeffs[elemI]*pnf[elemI];
272  }
273 }
274 
275 
276 template<class Type>
278 (
279  const scalarField& deltaCoeffs
280 ) const
281 {
282  // note: only applying coupled contribution
283  return coupledFvPatchField<Type>::snGrad(deltaCoeffs);
284 }
285 
286 
287 template<class Type>
289 {
290  // update non-overlap patch - some will implement updateCoeffs, and
291  // others will implement evaluate
292 
293  // scale neighbour field by (1 - mask)
294 
295  const scalarField& mask = cyclicACMIPatch_.cyclicACMIPatch().mask();
296  const fvPatchField<Type>& npf = nonOverlapPatchField();
297  const_cast<fvPatchField<Type>&>(npf).updateCoeffs(1.0 - mask);
298 }
299 
300 
301 template<class Type>
303 (
304  const Pstream::commsTypes comms
305 )
306 {
307  // update non-overlap patch (if not already updated by updateCoeffs)
308 
309  // scale neighbour field by (1 - mask)
310 
311  fvPatchField<Type>& npf =
312  const_cast<fvPatchField<Type>&>(nonOverlapPatchField());
313 
314  if (!npf.updated())
315  {
316  const scalarField& mask = cyclicACMIPatch_.cyclicACMIPatch().mask();
317 
318  npf.evaluate(comms);
319 
320  npf *= 1.0 - mask;
321  }
322 }
323 
324 
325 template<class Type>
327 (
328  const Pstream::commsTypes comms
329 )
330 {
331  // blend contributions from the coupled and non-overlap patches
332 
333  // neighbour patch field is updated via updateCoeffs or initEvaluate
334  // and is already scaled by (1 - mask)
335  const fvPatchField<Type>& npf = nonOverlapPatchField();
336 
338  const Field<Type>& cpf = *this;
339 
340  const scalarField& mask = cyclicACMIPatch_.cyclicACMIPatch().mask();
341  Field<Type>::operator=(mask*cpf + npf);
342 
344 }
345 
346 
347 template<class Type>
350 (
351  const tmp<scalarField>& w
352 ) const
353 {
354  // note: do not blend based on mask field
355  // - when applied this is scaled by the areas which are already scaled
357 }
358 
359 
360 template<class Type>
363 (
364  const tmp<scalarField>& w
365 ) const
366 {
367  // note: do not blend based on mask field
368  // - when applied this is scaled by the areas which are already scaled
370 }
371 
372 
373 template<class Type>
376 (
377  const scalarField& deltaCoeffs
378 ) const
379 {
380  // note: do not blend based on mask field
381  // - when applied this is scaled by the areas which are already scaled
383 }
384 
385 
386 template<class Type>
389 {
390  // note: do not blend based on mask field
391  // - when applied this is scaled by the areas which are already scaled
393 }
394 
395 
396 template<class Type>
399 (
400  const scalarField& deltaCoeffs
401 ) const
402 {
403  // note: do not blend based on mask field
404  // - when applied this is scaled by the areas which are already scaled
406 }
407 
408 
409 template<class Type>
412 {
413  // note: do not blend based on mask field
414  // - when applied this is scaled by the areas which are already scaled
416 }
417 
418 
419 template<class Type>
421 (
422  fvMatrix<Type>& matrix
423 )
424 {
425  const scalarField& mask = cyclicACMIPatch_.cyclicACMIPatch().mask();
426 
427  // nothing to be done by the AMI, but re-direct to non-overlap patch
428  // with non-overlap patch weights
429  const fvPatchField<Type>& npf = nonOverlapPatchField();
430 
431  const_cast<fvPatchField<Type>&>(npf).manipulateMatrix(matrix, 1.0 - mask);
432 }
433 
434 
435 template<class Type>
437 {
439  this->writeEntry("value", os);
440 }
441 
442 
443 // ************************************************************************* //
virtual void write(Ostream &os) const
Write.
virtual tmp< Field< Type > > valueInternalCoeffs(const tmp< scalarField > &) const
Return the matrix diagonal coefficients corresponding to the.
bool found(const word &, bool recursive=false, bool patternMatch=true) const
Search dictionary for given keyword.
Definition: dictionary.C:306
virtual void updateInterfaceMatrix(scalarField &result, const scalarField &psiInternal, const scalarField &coeffs, const direction cmpt, const Pstream::commsTypes commsType) const
Update result field based on interface functionality.
unsigned char direction
Definition: direction.H:43
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
void updateCoeffs()
Update the coefficients associated with the patch field.
const fvPatchField< Type > & nonOverlapPatchField() const
Return reference to non-overlapping patchField.
GeometricBoundaryField & boundaryField()
Return reference to GeometricBoundaryField.
virtual void evaluate(const Pstream::commsTypes commsType)
Evaluate the patch field.
virtual void initEvaluate(const Pstream::commsTypes commsType)
Initialise the evaluation of the patch field.
virtual tmp< Field< Type > > snGrad() const
Return patch-normal gradient.
virtual tmp< Field< Type > > gradientBoundaryCoeffs() const
Return the matrix source coefficients corresponding to the.
This boundary condition enforces a cyclic condition between a pair of boundaries, whereby communicati...
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
Foam::fvPatchFieldMapper.
A special matrix type and solver, designed for finite volume solutions of scalar equations. Face addressing is used to make all matrix assembly and solution loops vectorise.
Definition: fvPatchField.H:68
Abstract base class for cyclic ACMI coupled interfaces.
virtual void manipulateMatrix(fvMatrix< Type > &matrix)
Manipulate matrix.
bool updated() const
Return true if the boundary condition has already been updated.
Definition: fvPatchField.H:348
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
runTime write()
virtual bool coupled() const
Return true if coupled. Note that the underlying patch.
dictionary dict
virtual tmp< Field< Type > > gradientInternalCoeffs() const
Return the matrix diagonal coefficients corresponding to the.
IOerror FatalIOError
Abstract base class for coupled patches.
volScalarField & p
Definition: createFields.H:51
dimensionSet transform(const dimensionSet &)
Definition: dimensionSet.C:465
virtual tmp< Field< Type > > valueBoundaryCoeffs(const tmp< scalarField > &) const
Return the matrix source coefficients corresponding to the.
commsTypes
Types of communications.
Definition: UPstream.H:64
#define forAll(list, i)
Definition: UList.H:421
Spatial transformation functions for primitive fields.
cyclicACMIFvPatchField(const fvPatch &, const DimensionedField< Type, volMesh > &)
Construct from patch and internal field.
tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > snGrad(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvcSnGrad.C:45
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:65
Pre-declare SubField and related Field type.
Definition: Field.H:57
virtual void evaluate(const Pstream::commsTypes commsType=Pstream::blocking)
Evaluate the patch field, sets Updated to false.
Definition: fvPatchField.C:323
gmvFile<< "tracers "<< particles.size()<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){gmvFile<< iter().position().x()<< ' ';}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){gmvFile<< iter().position().y()<< ' ';}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){gmvFile<< iter().position().z()<< ' ';}gmvFile<< nl;forAll(lagrangianScalarNames, i){const word &name=lagrangianScalarNames[i];IOField< scalar > fld(IOobject( name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE ))
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
Definition: error.H:314
rDeltaT dimensionedInternalField()
Generic GeometricField class.
const cyclicACMIFvPatchField< Type > & neighbourPatchField() const
Return reference to neighbour patchField.
A 1D vector of objects of type <T>, where the size of the vector is known and can be used for subscri...
Definition: HashTable.H:60
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:61
virtual tmp< Field< Type > > patchNeighbourField() const
Return neighbour coupled internal cell data.
#define FatalIOErrorIn(functionName, ios)
Report an error message using Foam::FatalIOError.
Definition: error.H:325
A class for managing temporary objects.
Definition: PtrList.H:118
conserve internalField()+