/
HeatFlux.C
205 lines (171 loc) · 5.8 KB
/
HeatFlux.C
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
#include "HeatFlux.H"
#include "primitivePatchInterpolation.H"
#include "fvCFD.H"
using namespace Foam;
//----- preciceAdapter::CHT::HeatFlux -----------------------------------------
preciceAdapter::CHT::HeatFlux::HeatFlux(
const Foam::fvMesh& mesh,
const std::string nameT)
: T_(
const_cast<volScalarField*>(
&mesh.lookupObject<volScalarField>(nameT))),
mesh_(mesh)
{
dataType_ = scalar;
}
void preciceAdapter::CHT::HeatFlux::write(double* buffer, bool meshConnectivity, const unsigned int dim)
{
int bufferIndex = 0;
// For every boundary patch of the interface
for (uint j = 0; j < patchIDs_.size(); j++)
{
int patchID = patchIDs_.at(j);
const scalarField gradientPatch(
(T_->boundaryField()[patchID])
.snGrad());
// Extract the effective conductivity on the patch
extractKappaEff(patchID, meshConnectivity);
// If we use the mesh connectivity, we interpolate from the centres to the nodes
if (meshConnectivity)
{
//Setup Interpolation object
primitivePatchInterpolation patchInterpolator(mesh_.boundaryMesh()[patchID]);
scalarField gradientPoints;
//Interpolate
gradientPoints = patchInterpolator.faceToPointInterpolate(gradientPatch);
// For every cell of the patch
forAll(gradientPoints, i)
{
// Copy the heat flux into the buffer
// Q = - k * gradient(T)
//TODO: Interpolate kappa in case of a turbulent calculation
buffer[bufferIndex++] =
-getKappaEffAt(i) * gradientPoints[i];
}
}
else
{
// For every cell of the patch
forAll(gradientPatch, i)
{
// Copy the heat flux into the buffer
// Q = - k * gradient(T)
//TODO: Interpolate kappa in case of a turbulent calculation
buffer[bufferIndex++] =
-getKappaEffAt(i) * gradientPatch[i];
}
}
}
}
void preciceAdapter::CHT::HeatFlux::read(double* buffer, const unsigned int dim)
{
int bufferIndex = 0;
// For every boundary patch of the interface
for (uint j = 0; j < patchIDs_.size(); j++)
{
int patchID = patchIDs_.at(j);
// Extract the effective conductivity on the patch
// TODO: At the moment, reading with connectivity is not supported
extractKappaEff(patchID, /*meshConnectivity=*/false);
// Get the temperature gradient boundary patch
scalarField& gradientPatch(
refCast<fixedGradientFvPatchScalarField>(
T_->boundaryFieldRef()[patchID])
.gradient());
// For every cell of the patch
forAll(gradientPatch, i)
{
// Compute and assign the gradient from the buffer.
// The sign of the heat flux needs to be inversed,
// as the buffer contains the flux that enters the boundary:
// gradient(T) = -Q / -k
gradientPatch[i] =
buffer[bufferIndex++] / getKappaEffAt(i);
}
}
}
bool preciceAdapter::CHT::HeatFlux::isLocationTypeSupported(const bool meshConnectivity) const
{
// For cases with mesh connectivity, we support:
// - face nodes, only for writing
// - face centers, only for reading
// However, since we do not distinguish between reading and writing in the code, we
// always return true and offload the handling to the user.
if (meshConnectivity)
{
return true;
}
else
{
return (this->locationType_ == LocationType::faceCenters);
}
}
std::string preciceAdapter::CHT::HeatFlux::getDataName() const
{
return "HeatFlux";
}
//----- preciceAdapter::CHT::HeatFlux_Compressible ----------------------------
preciceAdapter::CHT::HeatFlux_Compressible::HeatFlux_Compressible(
const Foam::fvMesh& mesh,
const std::string nameT)
: HeatFlux(mesh, nameT),
Kappa_(new KappaEff_Compressible(mesh))
{
}
preciceAdapter::CHT::HeatFlux_Compressible::~HeatFlux_Compressible()
{
delete Kappa_;
}
void preciceAdapter::CHT::HeatFlux_Compressible::extractKappaEff(uint patchID, bool meshConnectivity)
{
Kappa_->extract(patchID, meshConnectivity);
}
scalar preciceAdapter::CHT::HeatFlux_Compressible::getKappaEffAt(int i)
{
return Kappa_->getAt(i);
}
//----- preciceAdapter::CHT::HeatFlux_Incompressible --------------------------
preciceAdapter::CHT::HeatFlux_Incompressible::HeatFlux_Incompressible(
const Foam::fvMesh& mesh,
const std::string nameT,
const std::string nameRho,
const std::string nameCp,
const std::string namePr,
const std::string nameAlphat)
: HeatFlux(mesh, nameT),
Kappa_(new KappaEff_Incompressible(mesh, nameRho, nameCp, namePr, nameAlphat))
{
}
preciceAdapter::CHT::HeatFlux_Incompressible::~HeatFlux_Incompressible()
{
delete Kappa_;
}
void preciceAdapter::CHT::HeatFlux_Incompressible::extractKappaEff(uint patchID, bool meshConnectivity)
{
Kappa_->extract(patchID, meshConnectivity);
}
scalar preciceAdapter::CHT::HeatFlux_Incompressible::getKappaEffAt(int i)
{
return Kappa_->getAt(i);
}
//----- preciceAdapter::CHT::HeatFlux_Basic -----------------------------------
preciceAdapter::CHT::HeatFlux_Basic::HeatFlux_Basic(
const Foam::fvMesh& mesh,
const std::string nameT,
const std::string nameKappa)
: HeatFlux(mesh, nameT),
Kappa_(new KappaEff_Basic(mesh, nameKappa))
{
}
preciceAdapter::CHT::HeatFlux_Basic::~HeatFlux_Basic()
{
delete Kappa_;
}
void preciceAdapter::CHT::HeatFlux_Basic::extractKappaEff(uint patchID, bool meshConnectivity)
{
Kappa_->extract(patchID, meshConnectivity);
}
scalar preciceAdapter::CHT::HeatFlux_Basic::getKappaEffAt(int i)
{
return Kappa_->getAt(i);
}