Go to the documentation of this file.
48 void Foam::cyclicFaPatch::calcTransforms()
50 const label sizeby2 = this->
size()/2;
55 for (label edgei=0; edgei < sizeby2; ++edgei)
58 half1Ctrs[edgei] = this->
edgeCentres()[edgei + sizeby2];
66 scalar maxMatchError = 0;
69 for (label edgei = 0; edgei < sizeby2; ++edgei)
71 half0Normals[edgei] = eN[edgei];
72 half1Normals[edgei] = eN[edgei + sizeby2];
74 scalar magLe =
mag(half0Normals[edgei]);
75 scalar nbrMagLe =
mag(half1Normals[edgei]);
76 scalar avLe = 0.5*(magLe + nbrMagLe);
78 if (magLe < ROOTVSMALL && nbrMagLe < ROOTVSMALL)
83 half0Normals[edgei] =
point(1, 0, 0);
84 half1Normals[edgei] = half0Normals[edgei];
95 half0Normals[edgei] /= magLe;
96 half1Normals[edgei] /= nbrMagLe;
103 label nbrEdgei = errorEdge + sizeby2;
104 scalar magLe =
mag(half0Normals[errorEdge]);
105 scalar nbrMagLe =
mag(half1Normals[errorEdge]);
106 scalar avLe = 0.5*(magLe + nbrMagLe);
109 <<
"edge " << errorEdge
110 <<
" area does not match neighbour "
111 << nbrEdgei <<
" by "
112 << 100*
mag(magLe - nbrMagLe)/avLe
113 <<
"% -- possible edge ordering problem." <<
endl
114 <<
"patch:" <<
name()
115 <<
" my area:" << magLe
116 <<
" neighbour area:" << nbrMagLe
119 <<
"Mesh edge:" <<
start() + errorEdge
121 <<
"Neighbour edge:" <<
start() + nbrEdgei
123 <<
"Other errors also exist, only the largest is reported. "
124 <<
"Please rerun with cyclic debug flag set"
143 <<
"Transformation tensor is not constant for the cyclic "
144 <<
"patch. Please reconsider your setup and definition of "
145 <<
"cyclic boundaries." <<
endl;
157 const label sizeby2 = deltas.size()/2;
159 scalar maxMatchError = 0;
160 label errorEdge = -1;
162 for (label edgei = 0; edgei < sizeby2; ++edgei)
164 scalar avL = 0.5*(magL[edgei] + magL[edgei + sizeby2]);
168 mag(magL[edgei] - magL[edgei + sizeby2])/avL
177 mag(magL[edgei] - magL[edgei + sizeby2])/avL
183 scalar di = deltas[edgei];
184 scalar dni = deltas[edgei + sizeby2];
186 w[edgei] = dni/(di + dni);
187 w[edgei + sizeby2] = 1 - w[edgei];
191 if (maxMatchError > matchTol_)
193 scalar avL = 0.5*(magL[errorEdge] + magL[errorEdge + sizeby2]);
196 <<
"edge " << errorEdge <<
" and " << errorEdge + sizeby2
197 <<
" areas do not match by "
198 << 100*
mag(magL[errorEdge] - magL[errorEdge + sizeby2])/avL
199 <<
"% -- possible edge ordering problem." <<
nl
200 <<
"Cyclic area match tolerance = "
201 << matchTol_ <<
" patch: " <<
name()
210 const label sizeby2 = deltas.size()/2;
212 for (label edgei = 0; edgei < sizeby2; ++edgei)
214 scalar di = deltas[edgei];
215 scalar dni = deltas[edgei + sizeby2];
217 dc[edgei] = 1.0/(di + dni);
218 dc[edgei + sizeby2] = dc[edgei];
252 const label sizeby2 = patchD.size()/2;
255 auto& pdv = tpdv.ref();
260 for (label edgei = 0; edgei < sizeby2; ++edgei)
262 const vector& ddi = patchD[edgei];
263 vector dni = patchD[edgei + sizeby2];
265 pdv[edgei] = ddi - dni;
266 pdv[edgei + sizeby2] = -pdv[edgei];
271 for (label edgei = 0; edgei < sizeby2; ++edgei)
273 const vector& ddi = patchD[edgei];
274 vector dni = patchD[edgei + sizeby2];
276 pdv[edgei] = ddi -
transform(forwardT()[0], dni);
277 pdv[edgei + sizeby2] = -
transform(reverseT()[0], pdv[edgei]);
290 return patchInternalField(internalData);
300 return patchInternalField(internalData, edgeFaces);
311 auto& pnf = tpnf.ref();
313 const label sizeby2 = this->size()/2;
315 for (label edgei=0; edgei < sizeby2; ++edgei)
317 pnf[edgei] = interfaceData[edgei + sizeby2];
318 pnf[edgei + sizeby2] = interfaceData[edgei];
331 return internalFieldTransfer(commsType, iF, this->
faceCells());
343 auto& pnf = tpnf.ref();
345 const label sizeby2 = this->size()/2;
347 for (label edgei=0; edgei < sizeby2; ++edgei)
349 pnf[edgei] = iF[edgeCells[edgei + sizeby2]];
350 pnf[edgei + sizeby2] = iF[edgeCells[edgei]];
addToRunTimeSelectionTable(decompositionMethod, kahipDecomp, dictionary)
vectorField pointField
pointField is a vectorField.
virtual void initMovePoints(const pointField &)
Initialise the patches for moving points.
virtual void calcGeometry()
Calculate the patch geometry.
void makeWeights(scalarField &) const
Make patch weighting factors.
bool parallel() const
Are the cyclic planes parallel.
A class for managing temporary objects.
virtual void movePoints(const pointField &)
Correct patch after moving points.
void makeDeltaCoeffs(scalarField &) const
Make patch face - neighbour cell distances.
virtual void calcGeometry()
Calculate the patch geometry.
Ostream & endl(Ostream &os)
Add newline and flush stream.
dimensionSet transform(const dimensionSet &ds)
Return the argument; transformations do not change the dimensions.
virtual tmp< vectorField > delta() const
Return cell-centre to face-centre vector.
Field< vector > vectorField
Specialisation of Field<T> for vector.
virtual void initMovePoints(const pointField &)
Initialise the patches for moving points.
const vectorField & edgeCentres() const
Return edge centres.
#define SeriousErrorInFunction
Report an error message using Foam::SeriousError.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
virtual tmp< labelField > interfaceInternalField(const labelUList &internalData) const
static const scalar matchTol_
Relative tolerance (for geometric matching). Is factor of.
virtual void initGeometry()
Initialise the calculation of the patch geometry.
Macros for easy insertion into run-time selection tables.
virtual const tensorField & forwardT() const
Return face transformation tensor.
errorManip< error > abort(error &err)
virtual const tensorField & reverseT() const
Return neighbour-cell transformation tensor.
label start() const
Patch start in edge list.
errorManipArg< error, int > exit(error &err, const int errNo=1)
commsTypes
Types of communications.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
void calcTransformTensors(const vector &Cf, const vector &Cr, const vector &nf, const vector &nr) const
Calculate the uniform transformation tensors.
const scalarField & magEdgeLengths() const
Return edge length magnitudes.
virtual tmp< labelField > transfer(const Pstream::commsTypes commsType, const labelUList &interfaceData) const
Transfer and return neighbour field.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
virtual void movePoints(const pointField &)
Correct patches after moving points.
const dimensionedScalar e
Elementary charge.
static tmp< T > New(Args &&... args)
Construct tmp of T with forwarding arguments.
tmp< vectorField > edgeNormals() const
Return edge normals.
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
virtual tmp< vectorField > delta() const
Return delta (P to N) vectors across coupled patch.
const word & name() const noexcept
The patch name.
virtual void initGeometry()
Initialise the calculation of the patch geometry.
vector point
Point is a vector.
virtual tmp< labelField > internalFieldTransfer(const Pstream::commsTypes commsType, const labelUList &internalData) const
Return neighbour field.
defineTypeNameAndDebug(combustionModel, 0)
Smooth ATC in cells next to a set of patches supplied by type.
virtual label size() const
Patch size is the number of edge labels.