temperatureCoupledBase.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2011-2022 OpenFOAM Foundation
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
11  OpenFOAM is free software: you can redistribute it and/or modify it
12  under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19  for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23 
24 \*---------------------------------------------------------------------------*/
25 
26 #include "temperatureCoupledBase.H"
27 #include "fluidThermo.H"
28 #include "solidThermo.H"
30 
31 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
32 
34 (
35  const fvPatch& patch
36 )
37 :
38  patch_(patch)
39 {}
40 
41 
43 (
44  const fvPatch& patch,
45  const dictionary& dict
46 )
47 :
48  patch_(patch)
49 {}
50 
51 
53 (
54  const fvPatch& patch,
55  const temperatureCoupledBase& base
56 )
57 :
58  patch_(patch)
59 {}
60 
61 
62 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
63 
65 (
66  const fvPatchScalarField& Tp
67 ) const
68 {
69  const fvMesh& mesh = patch_.boundaryMesh().mesh();
70  const label patchi = patch_.index();
71 
72  const word& phase(Tp.internalField().group());
73 
74  const word fluidThermoName
75  (
76  IOobject::groupName(physicalProperties::typeName, phase)
77  );
78 
79  if (mesh.foundObject<fluidThermo>(fluidThermoName))
80  {
81  static word ttmName
82  (
83  IOobject::groupName
84  (
85  thermophysicalTransportModel::typeName,
86  phase
87  )
88  );
89 
90  if (mesh.foundObject<thermophysicalTransportModel>(ttmName))
91  {
92  const thermophysicalTransportModel& ttm =
94 
95  return ttm.kappaEff(patchi);
96  }
97  else
98  {
99  const fluidThermo& thermo =
100  mesh.lookupObject<fluidThermo>(fluidThermoName);
101 
102  return thermo.kappa().boundaryField()[patchi];
103  }
104  }
105  else if (mesh.foundObject<solidThermo>(physicalProperties::typeName))
106  {
107  const solidThermo& thermo =
108  mesh.lookupObject<solidThermo>(physicalProperties::typeName);
109 
110  if (!thermo.isotropic())
111  {
112  const symmTensorField kappa(thermo.KappaLocal(patchi));
113  const vectorField n(patch_.nf());
114 
115  return n & kappa & n;
116  }
117  else
118  {
119  return thermo.kappa().boundaryField()[patchi];
120  }
121  }
122  else
123  {
125  << "Cannot find a fluidThermo or solidThermo instance"
126  << exit(FatalError);
127 
128  return scalarField::null();
129  }
130 }
131 
132 
134 {}
135 
136 
137 // ************************************************************************* //
virtual tmp< symmTensorField > KappaLocal(const label patchi) const =0
Anisotropic thermal conductivity for patch.
temperatureCoupledBase(const fvPatch &patch)
Construct from patch.
fluidReactionThermo & thermo
Definition: createFields.H:28
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:156
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:306
const Boundary & boundaryField() const
Return const-reference to the boundary field.
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:63
bool foundObject(const word &name) const
Is the named Type found?
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:66
fvMesh & mesh
const Type & lookupObject(const word &name) const
Lookup and return the object of the given Type.
const dimensionedScalar kappa
Coulomb constant: default SI units: [N.m2/C2].
virtual bool isotropic() const =0
Return true if thermal conductivity is isotropic.
A class for handling words, derived from string.
Definition: word.H:59
Base-class for fluid thermodynamic properties.
Definition: fluidThermo.H:53
virtual const volScalarField & kappa() const =0
Thermal conductivity of mixture [W/m/K].
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:54
Base-class for solid thermodynamic properties.
Definition: solidThermo.H:53
Abstract base class for thermophysical transport models (RAS, LES and laminar).
Common functions used in temperature coupled boundaries.
label patchi
virtual tmp< volScalarField > kappaEff() const =0
Effective thermal turbulent diffusivity for temperature.
tmp< scalarField > kappa(const fvPatchScalarField &Tp) const
Given patch temperature calculate corresponding K field.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:95
Single incompressible phase derived from the phase-fraction. Used as part of the multiPhaseMixture fo...
Definition: phase.H:52
label n
A class for managing temporary objects.
Definition: PtrList.H:53
void write(Ostream &) const
Write.
const DimensionedField< Type, volMesh > & internalField() const
Return dimensioned internal field reference.
Definition: fvPatchField.H:357