Go to the documentation of this file.
50 scalar maxRatio = 1 + coeff;
60 const label own = owner[facei];
61 const label nbr = neighbour[facei];
66 changedFaces.append(facei);
71 changedFaces.append(facei);
85 label facei =
patch.start() + patchFacei;
86 label own =
mesh.faceOwner()[facei];
88 changedFaces.append(facei);
94 changedFaces.shrink();
95 changedFacesInfo.shrink();
102 cellData[celli] =
field[celli];
109 td.maxRatio = maxRatio;
119 mesh.globalData().nTotalCells(),
125 field[celli] = cellData[celli].value();
128 field.correctBoundaryConditions();
137 const scalar alphaDiff,
139 const scalar alphaMin
152 cellData[celli] =
field[celli];
163 const label own = owner[facei];
164 const label nbr = neighbour[facei];
168 changedFaces.append(facei);
169 changedFacesInfo.append
185 label facei =
patch.start() + patchFacei;
186 label own =
mesh.faceOwner()[facei];
190 alpha.boundaryField()[patchi].patchNeighbourField()
193 if (
mag(
alpha[own] - alphapn[patchFacei]) > alphaDiff)
195 changedFaces.append(facei);
202 changedFaces.shrink();
203 changedFacesInfo.shrink();
217 smoothData.setFaceInfo(changedFaces, changedFacesInfo);
223 field[celli] = cellData[celli].value();
226 field.correctBoundaryConditions();
235 const scalar alphaDiff
255 const label own = owner[facei];
256 const label nbr = neighbour[facei];
260 changedFaces.append(facei);
261 changedFacesInfo.append
277 label facei =
patch.start() + patchFacei;
278 label own =
mesh.faceOwner()[facei];
282 alpha.boundaryField()[patchi].patchNeighbourField()
285 if (
mag(
alpha[own] - alphapn[patchFacei]) > alphaDiff)
287 changedFaces.append(facei);
288 changedFacesInfo.append
297 changedFaces.shrink();
298 changedFacesInfo.shrink();
308 sweepData.setFaceInfo(changedFaces, changedFacesInfo);
314 if (cellData[celli].valid(
sweepData.data()))
320 field.correctBoundaryConditions();
341 mesh.time().timeName(),
372 if (
alpha1[celli] < cutoff)
374 intvDotVapor.
value() +=
375 alpha2[celli]*mDotSmear[celli]*Vol[celli];
377 else if (
alpha1[celli] > 1.0 - cutoff)
379 intvDotLiquid.
value() +=
380 alpha1[celli]*mDotSmear[celli]*Vol[celli];
393 if (intvDotVapor.
value() > VSMALL)
395 Nv = intmSource0/intvDotVapor;
397 if (intvDotLiquid.
value() > VSMALL)
399 Nl = intmSource0/intvDotLiquid;
405 if (
alpha1[celli] < cutoff)
407 mDotOut[celli] = Nv.
value()*(1 -
alpha1[celli])*mDotSmear[celli];
409 else if (
alpha1[celli] > 1.0 - cutoff)
412 mDotOut[celli] = -Nl.value()*
alpha1[celli]*mDotSmear[celli];
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Class used to pass additional data in.
Helper class used by fvc::sweep function.
static constexpr const zero Zero
Global zero (0)
const volScalarField & alpha2
dimensioned< Type > domainIntegrate(const GeometricField< Type, fvPatchField, volMesh > &vf)
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
const Type & value() const
Return const reference to value.
const volScalarField & alpha1
void spreadSource(volScalarField &mDotOut, const volScalarField &mDotIn, const volScalarField &alpha1, const volScalarField &alpha2, const dimensionedScalar &D, const scalar cutoff)
Helper class used by the fvc::smooth and fvc::spread functions.
This boundary condition applies a zero-gradient condition from the patch internal field onto the patc...
#define forAll(list, i)
Loop across all elements in list.
A special matrix type and solver, designed for finite volume solutions of scalar equations.
void spread(volScalarField &field, const volScalarField &alpha, const label nLayers, const scalar alphaDiff=0.2, const scalar alphaMax=0.99, const scalar alphaMin=0.01)
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
A patch is a list of labels that address the faces in the global face list.
Calculate the matrix for the laplacian of the field.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
dimensionedScalar alphaMax("alphaMax", dimless/dimTime, laminarTransport)
reduce(hasMovingMesh, orOp< bool >())
Mesh data needed to do the Finite Volume discretisation.
const dimensionSet dimMass(1, 0, 0, 0, 0, 0, 0)
scalar maxRatio
Cut off distance.
Calculate the matrix for implicit and explicit sources.
Wave propagation of information through grid. Every iteration information goes through one layer of c...
const std::string patch
OpenFOAM patch number as a std::string.
A 1D array of objects of type <T>, where the size of the vector is known and used for subscript bound...
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
Volume integrate volField creating a volField.
const dimensionedScalar & D
A special matrix type and solver, designed for finite volume solutions of scalar equations....
SolverPerformance< Type > solve(const dictionary &)
Solve returning the solution statistics.
tmp< GeometricField< Type, faPatchField, areaMesh > > laplacian(const GeometricField< Type, faPatchField, areaMesh > &vf, const word &name)
void sweep(volScalarField &field, const volScalarField &alpha, const label nLayers, const scalar alphaDiff=0.2)
Provides functions smooth spread and sweep which use the FaceCellWave algorithm to smooth and redistr...
const dimensionSet dimless
Dimensionless.
void smooth(volScalarField &field, const scalar coeff)