Go to the documentation of this file.
39 template<
class CloudType>
51 this->owner().
name() +
":alpha",
52 this->owner().db().time().
timeName(),
59 zeroGradientFvPatchScalarField::typeName
63 applyLimiting_(this->coeffDict().
lookup(
"applyLimiting")),
64 applyGravity_(this->coeffDict().
lookup(
"applyGravity")),
65 alphaMin_(this->coeffDict().getScalar(
"alphaMin")),
66 rhoMin_(this->coeffDict().getScalar(
"rhoMin"))
68 alpha_ = this->owner().theta();
73 template<
class CloudType>
81 phiCorrect_(cm.phiCorrect_()),
82 uCorrect_(cm.uCorrect_()),
83 applyLimiting_(cm.applyLimiting_),
84 applyGravity_(cm.applyGravity_),
85 alphaMin_(cm.alphaMin_),
94 template<
class CloudType>
101 template<
class CloudType>
131 mesh.setFluxRequired(alpha_.name());
137 alpha_ =
max(this->owner().theta(), alphaMin_);
138 alpha_.correctBoundaryConditions();
146 this->owner().db().time().
timeName(),
156 rho.correctBoundaryConditions();
167 this->owner().db().time().
timeName(),
178 this->particleStressModel_->dTaudTheta
180 alpha_.primitiveField(),
181 rho.primitiveField(),
225 alphaEqn +=
fvm::div(phiGByA(), alpha_);
252 this->owner().db().time().
timeName(),
262 U.correctBoundaryConditions();
272 phiCorrect_.ref() -= phiGByA();
275 forAll(phiCorrect_(), facei)
278 const scalar phiCurr =
phi[facei];
279 scalar& phiCorr = phiCorrect_.ref()[facei];
283 if (phiCurr*phiCorr < 0)
289 else if (phiCorr > 0)
291 phiCorr =
max(phiCorr - phiCurr, 0);
295 phiCorr =
min(phiCorr - phiCurr, 0);
301 phiCorrect_.ref() += phiGByA();
315 uCorrect_->correctBoundaryConditions();
332 template<
class CloudType>
342 const label celli =
p.cell();
343 const label facei =
p.tetFace();
346 const vector U = uCorrect_()[celli];
350 const scalar nMag =
mag(nHat);
355 const label patchi =
mesh.boundaryMesh().whichPatch(facei);
358 phi = phiCorrect_()[facei];
363 phiCorrect_().boundaryField()[patchi]
365 mesh.boundaryMesh()[patchi].whichFace(facei)
370 const scalar t =
p.coordinates()[0];
374 return U + (1.0 - t)*nHat*(
phi/nMag - (
U & nHat));
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> reconstruct(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
const dimensionSet dimPressure
Defines the attributes of an object for which implicit objectRegistry management is supported,...
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Dimensionless.
const word cloudName(propsDict.get< word >("cloud"))
A class for handling words, derived from Foam::string.
A class for managing temporary objects.
static constexpr const zero Zero
Global zero (0)
const dimensionSet dimVelocity
const dimensionSet dimDensity
virtual tmp< Field< Type > > primitiveField() const =0
Return an internal field of the average.
virtual ~Implicit()
Destructor.
virtual void cacheFields(const bool store)
Calculate the inter particles stresses.
tmp< GeometricField< Type, faPatchField, areaMesh > > div(const GeometricField< Type, faePatchField, edgeMesh > &ssf)
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
This boundary condition applies a zero-gradient condition from the patch internal field onto the patc...
This boundary condition supplies a fixed value constraint, and is the base class for a number of othe...
Calculate the matrix for the divergence of the given field and flux.
#define forAll(list, i)
Loop across all elements in list.
dimensioned< vector > dimensionedVector
Dimensioned vector obtained from generic dimensioned type.
virtual vector velocityCorrection(typename CloudType::parcelType &p, const scalar deltaT) const
Calculate the velocity correction.
const fvMesh & mesh() const
Return reference to the mesh.
word name(const complex &c)
Return string representation of complex.
const cellModel * lookup(const word &modelName)
Deprecated(2017-11) equivalent to cellModel::ptr static method.
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.
Templated base class for dsmc cloud.
Base class for packing models.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Mesh data needed to do the Finite Volume discretisation.
const uniformDimensionedVectorField & g
Implicit(const dictionary &dict, CloudType &owner)
Construct from components.
Internal::FieldType & primitiveFieldRef(const bool updateAccessTime=true)
Return a reference to the internal field.
void correctBoundaryConditions()
Correct boundary field.
tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > linearInterpolate(const GeometricField< Type, fvPatchField, volMesh > &vf)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
tmp< GeometricField< Type, faPatchField, areaMesh > > ddt(const dimensioned< Type > dt, const faMesh &mesh)
Implicit model for applying an inter-particle stress to the particles.
A special matrix type and solver, designed for finite volume solutions of scalar equations....
SolverPerformance< Type > solve(const dictionary &)
Solve returning the solution statistics.
Reconstruct volField from a face flux field.
ParcelType parcelType
Type of parcel the cloud was instantiated for.
tmp< GeometricField< Type, faPatchField, areaMesh > > laplacian(const GeometricField< Type, faPatchField, areaMesh > &vf, const word &name)
tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > flux() const
Return the face-flux field from the matrix.
static tmp< GeometricField< Type, faePatchField, edgeMesh > > interpolate(const GeometricField< Type, faPatchField, areaMesh > &tvf, const edgeScalarField &faceFlux, Istream &schemeData)
Interpolate field onto faces using scheme given by Istream.
Calculate the matrix for the first temporal derivative.