62 void Foam::cyclicFaPatch::calcTransforms()
64 const label sizeby2 = this->size()/2;
69 for (label edgei=0; edgei < sizeby2; ++edgei)
72 half1Ctrs[edgei] = this->
edgeCentres()[edgei + sizeby2];
78 const vectorField eN(edgeNormals()*magEdgeLengths());
80 scalar maxMatchError = 0;
83 for (label edgei = 0; edgei < sizeby2; ++edgei)
85 half0Normals[edgei] = eN[edgei];
86 half1Normals[edgei] = eN[edgei + sizeby2];
88 scalar magLe =
mag(half0Normals[edgei]);
89 scalar nbrMagLe =
mag(half1Normals[edgei]);
90 scalar avLe = 0.5*(magLe + nbrMagLe);
92 if (magLe < ROOTVSMALL && nbrMagLe < ROOTVSMALL)
97 half0Normals[edgei] =
point(1, 0, 0);
98 half1Normals[edgei] = half0Normals[edgei];
100 else if (
mag(magLe - nbrMagLe)/avLe > matchTol_)
109 half0Normals[edgei] /= magLe;
110 half1Normals[edgei] /= nbrMagLe;
115 if (maxMatchError > matchTol_)
117 label nbrEdgei = errorEdge + sizeby2;
118 scalar magLe =
mag(half0Normals[errorEdge]);
119 scalar nbrMagLe =
mag(half1Normals[errorEdge]);
120 scalar avLe = 0.5*(magLe + nbrMagLe);
123 <<
"edge " << errorEdge
124 <<
" area does not match neighbour " 125 << nbrEdgei <<
" by " 126 << 100*
mag(magLe - nbrMagLe)/avLe
127 <<
"% -- possible edge ordering problem." <<
endl 128 <<
"patch:" <<
name()
129 <<
" my area:" << magLe
130 <<
" neighbour area:" << nbrMagLe
131 <<
" matching tolerance:" << matchTol_
133 <<
"Mesh edge:" << start() + errorEdge
135 <<
"Neighbour edge:" << start() + nbrEdgei
137 <<
"Other errors also exist, only the largest is reported. " 138 <<
"Please rerun with cyclic debug flag set" 154 if (forwardT().size() > 1 || reverseT().size() > 1)
157 <<
"Transformation tensor is not constant for the cyclic " 158 <<
"patch. Please reconsider your setup and definition of " 159 <<
"cyclic boundaries." <<
endl;
171 const label sizeby2 = deltas.size()/2;
173 scalar maxMatchError = 0;
174 label errorEdge = -1;
176 for (label edgei = 0; edgei < sizeby2; ++edgei)
178 scalar avL = 0.5*(magL[edgei] + magL[edgei + sizeby2]);
182 mag(magL[edgei] - magL[edgei + sizeby2])/avL
191 mag(magL[edgei] - magL[edgei + sizeby2])/avL
197 scalar di = deltas[edgei];
198 scalar dni = deltas[edgei + sizeby2];
200 w[edgei] = dni/(di + dni);
201 w[edgei + sizeby2] = 1 - w[edgei];
205 if (maxMatchError > matchTol_)
207 scalar avL = 0.5*(magL[errorEdge] + magL[errorEdge + sizeby2]);
210 <<
"edge " << errorEdge <<
" and " << errorEdge + sizeby2
211 <<
" areas do not match by " 212 << 100*
mag(magL[errorEdge] - magL[errorEdge + sizeby2])/avL
213 <<
"% -- possible edge ordering problem." <<
nl 214 <<
"Cyclic area match tolerance = " 215 << matchTol_ <<
" patch: " <<
name()
224 const label sizeby2 = deltas.size()/2;
226 for (label edgei = 0; edgei < sizeby2; ++edgei)
228 scalar di = deltas[edgei];
229 scalar dni = deltas[edgei + sizeby2];
231 dc[edgei] = 1.0/(di + dni);
232 dc[edgei + sizeby2] = dc[edgei];
274 const label sizeby2 = patchD.size()/2;
277 auto& pdv = tpdv.ref();
282 for (label edgei = 0; edgei < sizeby2; ++edgei)
284 const vector& ddi = patchD[edgei];
285 const vector& dni = patchD[edgei + sizeby2];
287 pdv[edgei] = ddi - dni;
288 pdv[edgei + sizeby2] = -pdv[edgei];
293 for (label edgei = 0; edgei < sizeby2; ++edgei)
295 const vector& ddi = patchD[edgei];
296 const vector& dni = patchD[edgei + sizeby2];
298 pdv[edgei] = ddi -
transform(forwardT()[0], dni);
299 pdv[edgei + sizeby2] = -
transform(reverseT()[0], pdv[edgei]);
312 return patchInternalField(internalData);
323 patchInternalField(internalData, edgeFaces, tpfld.ref());
335 auto& pnf = tpnf.ref();
337 const label sizeby2 = this->size()/2;
339 for (label edgei=0; edgei < sizeby2; ++edgei)
341 pnf[edgei] = interfaceData[edgei + sizeby2];
342 pnf[edgei + sizeby2] = interfaceData[edgei];
355 return internalFieldTransfer(commsType, iF, this->
faceCells());
367 auto& pnf = tpnf.ref();
369 const label sizeby2 = this->size()/2;
371 for (label edgei=0; edgei < sizeby2; ++edgei)
373 pnf[edgei] = iF[edgeCells[edgei + sizeby2]];
374 pnf[edgei + sizeby2] = iF[edgeCells[edgei]];
virtual void initMovePoints(PstreamBuffers &, const pointField &)
Initialise the patches for moving points.
virtual tmp< labelField > interfaceInternalField(const labelUList &internalData) const
Return the values of the given internal data adjacent to the interface as a field.
static const scalar matchTol_
Relative tolerance (for geometric matching). Is factor of.
errorManipArg< error, int > exit(error &err, const int errNo=1)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
commsTypes
Communications types.
error FatalError
Error stream (stdout output on all processes), with additional 'FOAM FATAL ERROR' header text and sta...
A list of keyword definitions, which are a keyword followed by a number of values (eg...
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
virtual tmp< labelField > internalFieldTransfer(const Pstream::commsTypes commsType, const labelUList &internalData) const
Return neighbour field.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
virtual void movePoints(PstreamBuffers &, const pointField &)
Correct patches after moving points.
constexpr char nl
The newline '\n' character (0x0a)
void makeDeltaCoeffs(scalarField &) const
Make patch face - neighbour cell distances.
virtual void initMovePoints(PstreamBuffers &, const pointField &)
Initialise the patches for moving points.
Ostream & endl(Ostream &os)
Add newline and flush stream.
virtual void calcGeometry(PstreamBuffers &)
Calculate the patch geometry.
Smooth ATC in cells next to a set of patches supplied by type.
virtual void initGeometry(PstreamBuffers &)
Initialise the calculation of the patch geometry.
virtual tmp< vectorField > delta() const =0
Return delta (P to N) vectors across coupled patch.
Macros for easy insertion into run-time selection tables.
UList< label > labelUList
A UList of labels.
#define SeriousErrorInFunction
Report an error message using Foam::SeriousError.
const dimensionedScalar e
Elementary charge.
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for expressions::valueTypeCode::INVALID.
virtual void calcGeometry(PstreamBuffers &)
Calculate the patch geometry.
Info<< nl<< "Wrote faMesh in vtk format: "<< writer.output().name()<< nl;}{ const Field< vector > edgeCentres(faMeshTools::flattenEdgeField(aMesh.edgeCentres(), true))
A class for handling words, derived from Foam::string.
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
void makeWeights(scalarField &) const
Make patch weighting factors.
static tmp< T > New(Args &&... args)
Construct tmp with forwarding arguments.
errorManip< error > abort(error &err)
A Vector of values with scalar precision, where scalar is float/double depending on the compilation f...
defineTypeNameAndDebug(combustionModel, 0)
Buffers for inter-processor communications streams (UOPstream, UIPstream).
cyclicFaPatch(const word &name, const dictionary &dict, const label index, const faBoundaryMesh &bm, const word &patchType)
Construct from dictionary.
vector point
Point is a vector.
coupledFaPatch is an abstract base class for patches that couple regions of the computational domain ...
Finite area boundary mesh.
Field< vector > vectorField
Specialisation of Field<T> for vector.
virtual tmp< labelField > transfer(const Pstream::commsTypes commsType, const labelUList &interfaceData) const
Transfer and return neighbour field.
A class for managing temporary objects.
dimensionSet transform(const dimensionSet &ds)
Return the argument; transformations do not change the dimensions.
virtual void movePoints(PstreamBuffers &, const pointField &)
Correct patch after moving points.
virtual tmp< vectorField > delta() const
Return delta (P to N) vectors across coupled patch.
virtual void initGeometry(PstreamBuffers &)
Initialise the calculation of the patch geometry.
addToRunTimeSelectionTable(functionObject, pointHistory, dictionary)