Go to the documentation of this file.
29 #include "surfaceInterpolate.H"
46 tmp<volScalarField> CoEulerDdtScheme<Type>::CorDeltaT()
const
50 tmp<volScalarField> tcorDeltaT
57 cofrDeltaT.instance(),
62 extrapolatedCalculatedFvPatchScalarField::typeName
73 corDeltaT[owner[facei]] =
74 max(corDeltaT[owner[facei]], cofrDeltaT[facei]);
76 corDeltaT[neighbour[facei]] =
77 max(corDeltaT[neighbour[facei]], cofrDeltaT[facei]);
80 const surfaceScalarField::Boundary& cofrDeltaTbf =
81 cofrDeltaT.boundaryField();
83 forAll(cofrDeltaTbf, patchi)
86 const fvPatch&
p = pcofrDeltaT.patch();
87 const labelUList& faceCells =
p.patch().faceCells();
89 forAll(pcofrDeltaT, patchFacei)
91 corDeltaT[faceCells[patchFacei]] =
max
93 corDeltaT[faceCells[patchFacei]],
94 pcofrDeltaT[patchFacei]
99 corDeltaT.correctBoundaryConditions();
106 tmp<surfaceScalarField> CoEulerDdtScheme<Type>::CofrDeltaT()
const
111 static_cast<const objectRegistry&
>(
mesh())
112 .lookupObject<surfaceScalarField>(phiName_);
114 if (
phi.dimensions() == dimensionSet(0, 3, -1, 0, 0))
123 return max(Co/maxCo_, scalar(1))/deltaT;
125 else if (
phi.dimensions() == dimensionSet(1, 0, -1, 0, 0))
128 static_cast<const objectRegistry&
>(
mesh())
129 .lookupObject<volScalarField>(rhoName_).oldTime();
138 return max(Co/maxCo_, scalar(1))/deltaT;
142 <<
"Incorrect dimensions of phi: " <<
phi.dimensions()
150 tmp<GeometricField<Type, fvPatchField, volMesh>>
160 "ddt("+dt.
name()+
')',
161 mesh().time().timeName(),
177 tdtdt.
ref().primitiveFieldRef() =
179 *(1.0 -
mesh().Vsc0()/
mesh().Vsc());
210 "ddt("+vf.name()+
')',
211 mesh().time().timeName(),
223 rDeltaT.dimensions()*vf.dimensions(),
262 "ddt("+
rho.name()+
','+vf.name()+
')',
263 mesh().time().timeName(),
275 rDeltaT.dimensions()*
rho.dimensions()*vf.dimensions(),
314 "ddt("+
rho.name()+
','+vf.name()+
')',
315 mesh().time().timeName(),
327 rDeltaT.dimensions()*
rho.dimensions()*vf.dimensions(),
331 -
rho.oldTime().primitiveField()
337 -
rho.oldTime().boundaryField()
371 mesh().time().timeName(),
387 alpha.primitiveField()
388 *
rho.primitiveField()
391 -
alpha.oldTime().primitiveField()
392 *
rho.oldTime().primitiveField()
397 alpha.boundaryField()
401 -
alpha.oldTime().boundaryField()
402 *
rho.oldTime().boundaryField()
444 scalarField rDeltaT(CorDeltaT()().primitiveField());
446 fvm.diag() = rDeltaT*
mesh().Vsc();
479 scalarField rDeltaT(CorDeltaT()().primitiveField());
481 fvm.diag() = rDeltaT*
rho.value()*
mesh().Vsc();
516 scalarField rDeltaT(CorDeltaT()().primitiveField());
518 fvm.diag() = rDeltaT*
rho.primitiveField()*
mesh().Vsc();
523 *
rho.oldTime().primitiveField()
529 *
rho.oldTime().primitiveField()
556 scalarField rDeltaT(CorDeltaT()().primitiveField());
559 rDeltaT*
alpha.primitiveField()*
rho.primitiveField()*
mesh().Vsc();
564 *
alpha.oldTime().primitiveField()
565 *
rho.oldTime().primitiveField()
571 *
alpha.oldTime().primitiveField()
572 *
rho.oldTime().primitiveField()
602 "ddtCorr(" +
U.name() +
',' +
Uf.name() +
')',
603 mesh().time().timeName(),
606 this->fvcDdtPhiCoeff(
U.oldTime(), phiUf0, phiCorr)
634 "ddtCorr(" +
U.name() +
',' +
phi.name() +
')',
635 mesh().time().timeName(),
638 this->fvcDdtPhiCoeff(
U.oldTime(),
phi.oldTime(), phiCorr)
664 rho.oldTime()*
U.oldTime()
677 +
rho.name() +
',' +
U.name() +
',' +
Uf.name() +
')',
678 mesh().time().timeName(),
681 this->fvcDdtPhiCoeff(rhoU0, phiUf0, phiCorr,
rho.oldTime())
705 +
rho.name() +
',' +
U.name() +
',' +
Uf.name() +
')',
706 mesh().time().timeName(),
722 <<
"dimensions of Uf are not correct"
725 return fluxFieldType::null();
749 rho.oldTime()*
U.oldTime()
764 +
rho.name() +
',' +
U.name() +
',' +
phi.name() +
')',
765 mesh().time().timeName(),
796 +
rho.name() +
',' +
U.name() +
',' +
phi.name() +
')',
797 mesh().time().timeName(),
813 <<
"dimensions of phi are not correct"
816 return fluxFieldType::null();
845 tmeshPhi.
ref().setOriented();
Defines the attributes of an object for which implicit objectRegistry management is supported,...
static tmp< GeometricField< typename innerProduct< vector, Type >::type, fvsPatchField, surfaceMesh > > dotInterpolate(const surfaceVectorField &Sf, const GeometricField< Type, fvPatchField, volMesh > &tvf)
Interpolate field onto faces.
A class for managing temporary objects.
static constexpr const zero Zero
Global zero (0)
const dimensionSet dimVelocity
const dimensionSet dimDensity
tmp< GeometricField< Type, fvPatchField, volMesh > > fvcDdt(const dimensioned< Type > &)
fvsPatchField< scalar > fvsPatchScalarField
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
const Internal::FieldType & primitiveField() const
Return a const-reference to the internal field.
autoPtr< surfaceVectorField > Uf
const GeometricField< Type, PatchField, GeoMesh > & oldTime() const
Return old time field.
Calculate the divergence of the given field.
const word & name() const
Return const reference to name.
const Type & value() const
Return const reference to value.
tmp< fvMatrix< Type > > fvmDdt(const GeometricField< Type, fvPatchField, volMesh > &)
const surfaceScalarField & deltaCoeffs() const
Return reference to cell-centre difference coefficients.
#define forAll(list, i)
Loop across all elements in list.
A special matrix type and solver, designed for finite volume solutions of scalar equations.
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
const dimensionSet dimArea(sqr(dimLength))
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
Generic dimensioned Type class.
This boundary condition is not designed to be evaluated; it is assmued that the value is assigned via...
errorManip< error > abort(error &err)
static tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > interpolate(const GeometricField< Type, fvPatchField, volMesh > &tvf, const surfaceScalarField &faceFlux, Istream &schemeData)
Interpolate field onto faces using scheme given by Istream.
tmp< fluxFieldType > fvcDdtUfCorr(const GeometricField< Type, fvPatchField, volMesh > &U, const GeometricField< Type, fvsPatchField, surfaceMesh > &Uf)
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
Internal & ref(const bool updateAccessTime=true)
Return a reference to the dimensioned internal field.
tmp< fluxFieldType > fvcDdtPhiCorr(const GeometricField< Type, fvPatchField, volMesh > &U, const fluxFieldType &phi)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
A special matrix type and solver, designed for finite volume solutions of scalar equations....
const dimensionSet dimVol(dimVolume)
Older spelling for dimVolume.
UList< label > labelUList
A UList of labels.
tmp< surfaceScalarField > meshPhi(const GeometricField< Type, fvPatchField, volMesh > &)
const dimensionSet & dimensions() const
Return const reference to dimensions.
const dimensionSet dimVolume(pow3(dimLength))
const Boundary & boundaryField() const
Return const-reference to the boundary field.