Go to the documentation of this file.
44 cyclicACMIPatch_(refCast<const cyclicACMIFvPatch>(
p))
58 cyclicACMIPatch_(refCast<const cyclicACMIFvPatch>(
p,
dict))
60 if (!isA<cyclicACMIFvPatch>(
p))
63 <<
" patch type '" <<
p.type()
64 <<
"' not constraint type '" << typeName <<
"'"
65 <<
"\n for patch " <<
p.name()
66 <<
" of field " << this->internalField().name()
67 <<
" in file " << this->internalField().objectPath()
79 this->primitiveField()
81 if (!
fld.boundaryField().set(cyclicACMIPatch_.nonOverlapPatchID()))
84 <<
" patch " <<
p.name()
85 <<
" of field " << this->internalField().name()
86 <<
" refers to non-overlap patch "
87 << cyclicACMIPatch_.cyclicACMIPatch().nonOverlapPatchName()
88 <<
" which is not constructed yet." <<
nl
89 <<
" Either supply an initial value or change the ordering"
94 this->
evaluate(Pstream::commsTypes::blocking);
110 cyclicACMIPatch_(refCast<const cyclicACMIFvPatch>(
p))
112 if (!isA<cyclicACMIFvPatch>(this->
patch()))
115 <<
"' not constraint type '" << typeName <<
"'"
116 <<
"\n for patch " <<
p.name()
117 <<
" of field " << this->internalField().name()
118 <<
" in file " << this->internalField().objectPath()
133 cyclicACMIPatch_(ptf.cyclicACMIPatch_)
146 cyclicACMIPatch_(ptf.cyclicACMIPatch_)
155 return cyclicACMIPatch_.coupled();
163 const Field<Type>& iField = this->primitiveField();
173 cyclicACMIPatch_.interpolate
200 this->primitiveField()
203 return refCast<const cyclicACMIFvPatchField<Type>>
205 fld.boundaryField()[cyclicACMIPatch_.neighbPatchID()]
217 this->primitiveField()
221 return fld.boundaryField()[cyclicACMIPatch_.nonOverlapPatchID()];
247 lduAddr.
patchAddr(cyclicACMIPatch_.neighbPatchID());
252 transformCoupleField(pnf, cmpt);
254 pnf = cyclicACMIPatch_.interpolate(pnf);
258 this->addToInternalField(result, !
add,
faceCells, coeffs, pnf);
277 lduAddr.
patchAddr(cyclicACMIPatch_.neighbPatchID());
282 transformCoupleField(pnf);
284 pnf = cyclicACMIPatch_.interpolate(pnf);
288 this->addToInternalField(result, !
add,
faceCells, coeffs, pnf);
298 const scalarField& mask = cyclicACMIPatch_.cyclicACMIPatch().mask();
316 if (this->cyclicACMIPatch().owner())
318 label index = this->
patch().index();
320 const label globalPatchID =
338 const labelUList& u = matrix.lduAddr().upperAddr();
339 const labelUList& l = matrix.lduAddr().lowerAddr();
348 label globalFaceI =
faceMap[j];
350 const scalar boundCorr = -bndCoeffs[subFaceI];
351 const scalar intCorr = -intCoeffs[subFaceI];
353 matrix.upper()[globalFaceI] += boundCorr;
354 matrix.diag()[u[globalFaceI]] -= intCorr;
355 matrix.diag()[l[globalFaceI]] -= boundCorr;
357 if (matrix.asymmetric())
359 matrix.lower()[globalFaceI] += intCorr;
367 if (matrix.
psi(mat).mesh().fluxRequired(this->internalField().
name()))
378 const label nbrPathID =
379 cyclicACMIPatch_.cyclicACMIPatch().neighbPatchID();
381 const label nbrGlobalPatchID =
407 const label index(this->
patch().index());
409 const label nSubFaces
417 cyclicACMIPatch_.cyclicACMIPatch().AMI().srcWeights();
419 const scalarField& mask = cyclicACMIPatch_.cyclicACMIPatch().mask();
421 const scalar tol = cyclicACMIPolyPatch::tolerance();
426 for(label i=0; i<w.size(); i++)
428 if (mask[faceI] > tol)
430 const label localFaceId =
432 [mat][index][subFaceI];
433 mapCoeffs[subFaceI] = w[i]*coeffs[localFaceId];
439 return tmp<Field<scalar>>(
new Field<scalar>(mapCoeffs));
452 const scalarField& mask = cyclicACMIPatch_.cyclicACMIPatch().mask();
462 this->writeEntry(
"value",
os);
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
The class contains the addressing required by the lduMatrix: upper, lower and losort.
const lduPrimitiveMeshAssembly & lduMeshAssembly()
Return optional lduAdressing.
Pair< int > faceMap(const label facePi, const face &faceP, const label faceNi, const face &faceN)
const fvPatchField< Type > & nonOverlapPatchField() const
Return reference to non-overlapping patchField.
A class for managing temporary objects.
static constexpr const zero Zero
Global zero (0)
This boundary condition enforces a cyclic condition between a pair of boundaries, whereby communicati...
const FieldField< Field, Type > & boundaryCoeffs() const
dimensionSet transform(const dimensionSet &ds)
Return the argument; transformations do not change the dimensions.
Abstract base class for cyclic ACMI coupled interfaces.
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
#define forAll(list, i)
Loop across all elements in list.
const cyclicACMIFvPatchField< Type > & neighbourPatchField() const
Return reference to neighbour patchField.
Generic templated field type.
virtual bool coupled() const
Return true if coupled. Note that the underlying patch.
virtual void manipulateMatrix(fvMatrix< Type > &matrix)
Manipulate matrix.
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
virtual const cyclicACMIFvPatch & neighbPatch() const
Return neighbour fvPatch.
gmvFile<< "tracers "<< particles.size()<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().x()<< ' ';}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().y()<< ' ';}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().z()<< ' ';}gmvFile<< nl;for(const word &name :lagrangianScalarNames){ IOField< scalar > fld(IOobject(name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
const GeometricField< Type, fvPatchField, volMesh > & psi(const label i=0) const
Return psi.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
void add(FieldField< Field1, typename typeOfSum< Type1, Type2 >::type > &f, const FieldField< Field1, Type1 > &f1, const FieldField< Field2, Type2 > &f2)
OBJstream os(runTime.globalPath()/outputName)
virtual void updateInterfaceMatrix(solveScalarField &result, const bool add, const lduAddressing &lduAddr, const label patchId, const solveScalarField &psiInternal, const scalarField &coeffs, const direction cmpt, const Pstream::commsTypes commsType) const
Update result field based on interface functionality.
errorManipArg< error, int > exit(error &err, const int errNo=1)
virtual const labelUList & faceCells() const
Return faceCell addressing.
commsTypes
Types of communications.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
const std::string patch
OpenFOAM patch number as a std::string.
Abstract base class for coupled patches.
A traits class, which is primarily used for primitives.
virtual void write(Ostream &os) const
Write.
cyclicACMIFvPatchField(const fvPatch &, const DimensionedField< Type, volMesh > &)
Construct from patch and internal field.
virtual tmp< Field< Type > > patchNeighbourField() const
Return neighbour coupled internal cell data.
A special matrix type and solver, designed for finite volume solutions of scalar equations....
void write(vtk::formatter &fmt, const Type &val, const label n=1)
Component-wise write of a value (N times)
virtual const labelUList & patchAddr(const label patchNo) const =0
Return patch to internal addressing given patch number.
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
Foam::fvPatchFieldMapper.
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
bool coupled(solutionDict.getOrDefault("coupledEnergyField", false))
const FieldField< Field, Type > & internalCoeffs() const
Smooth ATC in cells next to a set of patches supplied by type.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
Cyclic patch for Arbitrarily Coupled Mesh Interface (ACMI)
string evaluate(label fieldWidth, const std::string &s, size_t pos=0, size_t len=std::string::npos)
String evaluation with specified (positive, non-zero) field width.