Go to the documentation of this file.
31 #include "surfaceInterpolate.H"
45 void Foam::interfaceProperties::correctContactAngle
47 surfaceVectorField::Boundary& nHatb,
48 const surfaceVectorField::Boundary& gradAlphaf
51 const fvMesh&
mesh = alpha1_.mesh();
52 const volScalarField::Boundary& abf = alpha1_.
boundaryField();
58 if (isA<alphaContactAngleTwoPhaseFvPatchScalarField>(abf[patchi]))
60 alphaContactAngleTwoPhaseFvPatchScalarField& acap =
61 const_cast<alphaContactAngleTwoPhaseFvPatchScalarField&
>
63 refCast<const alphaContactAngleTwoPhaseFvPatchScalarField>
88 b2[facei] =
cos(
acos(a12[facei]) - theta[facei]);
96 nHatp = a*nf +
b*nHatp;
97 nHatp /= (
mag(nHatp) + deltaN_.
value());
99 acap.gradient() = (nf & nHatp)*
mag(gradAlphaf[patchi]);
106 void Foam::interfaceProperties::calculateK()
108 const fvMesh&
mesh = alpha1_.mesh();
128 correctContactAngle(nHatfv.boundaryFieldRef(), gradAlphaf.boundaryField());
131 nHatf_ = nHatfv & Sf;
152 Foam::interfaceProperties::interfaceProperties
156 const IOdictionary&
dict
159 transportPropertiesDict_(
dict),
192 "interfaceProperties:K",
209 return sigmaPtr_->sigma()*K_;
223 return pos0(alpha1_ - 0.01)*
pos0(0.99 - alpha1_);
235 alpha1_.mesh().solverDict(alpha1_.name()).readEntry(
"cAlpha", cAlpha_);
236 sigmaPtr_->readDict(transportPropertiesDict_);
tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > snGrad(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
A class for managing temporary objects.
static constexpr const zero Zero
Global zero (0)
Calculate the snGrad of the given volField.
tmp< volScalarField > sigmaK() const
Calculate the divergence of the given field.
Unit conversion functions.
tmp< GeometricField< Type, fvPatchField, volMesh > > div(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
const Type & value() const
Return const reference to value.
dimensionedScalar pos0(const dimensionedScalar &ds)
const volScalarField & alpha1
#define forAll(list, i)
Loop across all elements in list.
Field< vector > vectorField
Specialisation of Field<T> for vector.
tmp< volScalarField > nearInterface() const
Indicator of the proximity of the interface.
const dimensionSet dimArea(sqr(dimLength))
const dimensionedScalar b
Wien displacement law constant: default SI units: [m.K].
fvsPatchField< vector > fvsPatchVectorField
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
tmp< surfaceScalarField > surfaceTensionForce() const
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.
constexpr scalar degToRad(const scalar deg) noexcept
Conversion from degrees to radians.
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
Global function forwards to reuseTmpDimensionedField::New.
GeometricField< vector, fvPatchField, volMesh > volVectorField
dimensionedScalar acos(const dimensionedScalar &ds)
Calculate the gradient of the given field.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
const dimensionedScalar e
Elementary charge.
List< ReturnType > get(const UPtrList< T > &list, const AccessOp &aop)
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
dimensionedScalar det(const dimensionedSphericalTensor &dt)
GeometricField< vector, fvsPatchField, surfaceMesh > surfaceVectorField
dimensionedScalar cbrt(const dimensionedScalar &ds)
const Boundary & boundaryField() const
Return const-reference to the boundary field.
dimensioned< Type > average(const DimensionedField< Type, GeoMesh > &df)
const dimensionSet dimless
Dimensionless.
bool read()
Read transportProperties dictionary.
dimensionedScalar cos(const dimensionedScalar &ds)
const surfaceVectorField & Sf() const
Return cell face area vectors.