36 namespace solidThermophysicalTransportModels
51 void Foam::solidThermophysicalTransportModels::anisotropic::
52 setZonesPatchFaces()
const
54 if (!zoneCoordinateSystems_.size())
return;
59 const fvBoundaryMesh&
patches = mesh.boundary();
62 zonesPatchFaces_.
setSize(zoneCoordinateSystems_.size());
67 PtrDictionary<coordinateSystem>,
68 zoneCoordinateSystems_,
74 const labelList& zoneCells = mesh.cellZones()[iter().name()];
77 boolList cellInZone(mesh.nCells(),
false);
81 cellInZone[zoneCells[i]] =
true;
90 label nZonePatchFaces = 0;
94 const label facei = pp.start() + patchFacei;
96 if (cellInZone[own[facei]])
98 zonesPatchFaces_[zonei][
patchi][nZonePatchFaces++] =
123 coeffDict().lookupOrDefault<
Switch>(
"boundaryAligned", false)
131 Info<<
" Reading coordinate system for zones:" <<
endl;
137 const word& zoneName = iter().keyword();
142 zoneCoordinateSystems_.
insert
151 setZonesPatchFaces();
171 const scalar alignmentFactor =
175 *
mag(nKappa -
n*(nKappa &
n))
178 if (alignmentFactor > 1
e-3)
183 Info<<
" Kappa is not aligned with patch "
185 <<
" by an alignment factor of " << alignmentFactor
186 <<
" (0=aligned, 1=unaligned)" <<
nl
187 <<
" heat-flux correction will be applied." <<
endl;
192 if (!aligned && boundaryAligned_)
198 " boundaryAligned is set true, "
199 "boundary alignment of kappa will be enforced." <<
endl;
227 Foam::solidThermophysicalTransportModels::anisotropic::Kappa()
const
246 coordinateSystem_.R(mesh.
C()).transformDiagTensor(materialKappa);
260 zoneCoordinateSystems_,
271 const label celli = zoneCells[i];
292 const label patchFacei =
293 zonesPatchFaces_[zonei][
patchi][zonePatchFacei];
295 KappaPf[patchFacei] =
299 materialKappaPf[patchFacei]
312 Foam::solidThermophysicalTransportModels::anisotropic::Kappa
322 tmp<symmTensorField> tKappa
324 coordinateSystem_.R(CPf).transformDiagTensor(materialKappaPf)
332 PtrDictionary<coordinateSystem>,
333 zoneCoordinateSystems_,
337 const coordinateSystem& cs = iter();
341 const label patchFacei =
342 zonesPatchFaces_[zonei][
patchi][zonePatchFacei];
344 KappaPf[patchFacei] =
345 cs.R().transformDiagTensor
348 materialKappaPf[patchFacei]
405 return -(nKappa -
n*(nKappa &
n)) & gradT().boundaryField()[
patchi];
443 if (!zonesPatchFaces_.size())
445 setZonesPatchFaces();
462 zonesPatchFaces_.clear();
472 zonesPatchFaces_.clear();
482 zonesPatchFaces_.clear();
#define forAll(list, i)
Loop across all elements in list.
#define forAllConstIter(Container, container, iter)
Iterate across all elements in the container object of type.
Macros for easy insertion into run-time selection tables.
void insert(const word &, T *)
Add at head of dictionary.
const dimensionSet & dimensions() const
Return dimensions.
Generic GeometricField class.
Internal & ref()
Return a reference to the dimensioned internal field.
Internal::FieldType & primitiveFieldRef()
Return a reference to the internal field.
static tmp< GeometricField< Type, PatchField, GeoMesh > > New(const word &name, const Internal &, const PtrList< PatchField< Type >> &)
Return a temporary field constructed from name,.
const Boundary & boundaryField() const
Return const-reference to the boundary field.
Boundary & boundaryFieldRef()
Return a reference to the boundary field.
const word & name() const
Return name.
void setSize(const label)
Reset size of List.
A simple wrapper around bool so that it can be read as a word: true/false, on/off,...
virtual const fvMesh & mesh() const =0
Return const access to the mesh.
virtual const volScalarField & T() const =0
Temperature [K].
virtual symmTensor transformDiagTensor(const vector &p, const vector &v) const =0
Transform diagTensor masquerading as a vector using transformation.
Base class for other coordinate system specifications.
static autoPtr< coordinateSystem > New(const objectRegistry &obr, const dictionary &dict)
Select constructed from dictionary and objectRegistry.
const coordinateRotation & R() const
Return const reference to co-ordinate rotation.
A list of keyword definitions, which are a keyword followed by any number of values (e....
bool isDict(const word &) const
Check if entry is a sub-dictionary.
const dictionary & subDict(const word &) const
Find and return a sub-dictionary.
bool found(const word &, bool recursive=false, bool patternMatch=true) const
Search dictionary for given keyword.
Mesh data needed to do the Finite Volume discretisation.
const volVectorField & C() const
Return cell centres.
const fvBoundaryMesh & boundary() const
Return reference to boundary mesh.
const fvSchemes & schemes() const
Return the fvSchemes.
const fvSolution & solution() const
Return the fvSchemes.
const surfaceScalarField & magSf() const
Return cell face area magnitudes.
void setFluxRequired(const word &name) const
Class containing mesh-to-mesh mapping information after a mesh distribution where we send parts of me...
Class containing mesh-to-mesh mapping information.
const meshCellZones & cellZones() const
Return cell zones.
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
Base-class for solid thermodynamic properties.
Abstract base class for solid thermophysical transport models.
virtual void printCoeffs(const word &type)
Print model coefficients.
virtual const dictionary & coeffDict() const
Const access to the coefficients dictionary.
virtual void predict()=0
Predict the thermophysical transport coefficients if possible.
virtual const solidThermo & thermo() const
Access function to solid thermophysical properties.
Solid thermophysical transport model for anisotropic thermal conductivity.
virtual bool movePoints()
Update for mesh motion.
virtual tmp< scalarField > qCorr(const label patchi) const
Return the patch heat flux correction [W/m^2].
virtual void predict()
Correct the anisotropic viscosity.
virtual void mapMesh(const polyMeshMap &map)
Update from another mesh using the given map.
anisotropic(const solidThermo &thermo)
Construct from solid thermophysical properties.
virtual void distribute(const polyDistributionMap &map)
Redistribute or update using the given distribution map.
virtual tmp< volScalarField > kappa() const
Thermal conductivity [W/m/K].
virtual tmp< surfaceScalarField > q() const
Return the heat flux [W/m^2].
virtual void topoChange(const polyTopoChangeMap &map)
Update topology using the given map.
virtual bool read()
Read thermophysicalTransport dictionary.
virtual tmp< fvScalarMatrix > divq(volScalarField &he) const
Return the source term for the energy equation.
void enableCache(const word &name) const
Enable caching of the given field.
A class for managing temporary objects.
A class for handling words, derived from string.
#define NotImplemented
Issue a FatalErrorIn for a function not currently implemented.
Calculate the laplacian of the given field.
Calculate the snGrad of the given volField.
Calculate the matrix for the laplacian of the field.
const polyBoundaryMesh & bMesh
const fvPatchList & patches
static tmp< SurfaceField< Type > > interpolate(const VolField< Type > &tvf, const surfaceScalarField &faceFlux, Istream &schemeData)
Interpolate field onto faces using scheme given by Istream.
tmp< VolField< Type > > laplacian(const VolField< Type > &vf, const word &name)
tmp< VolField< typename outerProduct< vector, Type >::type > > grad(const SurfaceField< Type > &ssf)
tmp< fvMatrix< Type > > laplacian(const VolField< Type > &vf, const word &name)
tmp< fvMatrix< Type > > laplacianCorrection(const VolField< scalar > &gamma, const VolField< Type > &vf)
defineTypeNameAndDebug(anisotropic, 0)
addToRunTimeSelectionTable(solidThermophysicalTransportModel, anisotropic, dictionary)
Type gSum(const FieldField< Field, Type > &f)
List< label > labelList
A List of labels.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
Ostream & endl(Ostream &os)
Add newline and flush stream.
List< bool > boolList
Bool container classes.
T returnReduce(const T &Value, const BinaryOp &bop, const int tag=Pstream::msgType(), const label comm=UPstream::worldComm)
dimensioned< scalar > mag(const dimensioned< Type > &)
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
dimensioned< symmTensor > dimensionedSymmTensor
Dimensioned tensor obtained from generic dimensioned type.
Field< vector > vectorField
Specialisation of Field<T> for vector.
Field< symmTensor > symmTensorField
Specialisation of Field<T> for symmTensor.
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
faceListList boundary(nPatches)
fluidMulticomponentThermo & thermo