37 template<
class AlphaFieldType,
class RhoFieldType>
38 void Foam::fv::actuationDiskSource::calc
40 const AlphaFieldType&
alpha,
41 const RhoFieldType&
rho,
47 case forceMethodType::FROUDE:
53 case forceMethodType::VARIABLE_SCALING:
55 calcVariableScalingMethod(
alpha,
rho, eqn);
65 template<
class AlphaFieldType,
class RhoFieldType>
66 void Foam::fv::actuationDiskSource::calcFroudeMethod
68 const AlphaFieldType&
alpha,
69 const RhoFieldType&
rho,
80 label szMonitorCells = monitorCells_.size();
82 for (
const label celli : monitorCells_)
85 rhoRef = rhoRef +
rho[celli];
87 reduce(Uref, sumOp<vector>());
88 reduce(rhoRef, sumOp<scalar>());
89 reduce(szMonitorCells, sumOp<label>());
91 if (szMonitorCells == 0)
94 <<
"No cell is available for incoming velocity monitoring." 98 Uref /= szMonitorCells;
99 rhoRef /= szMonitorCells;
101 const scalar Ct = sink_*UvsCtPtr_->value(
mag(Uref));
102 const scalar
Cp = sink_*UvsCpPtr_->value(
mag(Uref));
104 if (
Cp <= VSMALL || Ct <= VSMALL)
107 <<
"Cp and Ct must be greater than zero." <<
nl 108 <<
"Cp = " <<
Cp <<
", Ct = " << Ct
113 const vector diskDir = this->diskDir();
114 const scalar a = 1.0 -
Cp/Ct;
115 const scalar
T = 2.0*rhoRef*diskArea_*
magSqr(Uref & diskDir)*a*(1 - a);
117 for (
const label celli : cells_)
119 Usource[celli] += ((cellsV[celli]/V())*
T)*diskDir;
124 mesh_.time().timeOutputValue() >= writeFileStart_
125 && mesh_.time().timeOutputValue() <= writeFileEnd_
128 Ostream&
os = file();
129 writeCurrentTime(
os);
137 template<
class AlphaFieldType,
class RhoFieldType>
138 void Foam::fv::actuationDiskSource::calcVariableScalingMethod
140 const AlphaFieldType&
alpha,
141 const RhoFieldType&
rho,
142 fvMatrix<vector>& eqn
152 label szMonitorCells = monitorCells_.size();
154 for (
const label celli : monitorCells_)
157 rhoRef = rhoRef +
rho[celli];
159 reduce(Uref, sumOp<vector>());
160 reduce(rhoRef, sumOp<scalar>());
161 reduce(szMonitorCells, sumOp<label>());
163 if (szMonitorCells == 0)
166 <<
"No cell is available for incoming velocity monitoring." 170 Uref /= szMonitorCells;
171 const scalar magUref =
mag(Uref);
172 rhoRef /= szMonitorCells;
176 scalar rhoDisk = 0.0;
179 for (
const label celli : cells_)
181 Udisk +=
U[celli]*cellsV[celli];
182 rhoDisk +=
rho[celli]*cellsV[celli];
183 totalV += cellsV[celli];
185 reduce(Udisk, sumOp<vector>());
186 reduce(rhoDisk, sumOp<scalar>());
187 reduce(totalV, sumOp<scalar>());
192 <<
"No cell in the actuator disk." 197 const scalar magUdisk =
mag(Udisk);
200 if (
mag(Udisk) < SMALL)
203 <<
"Velocity spatial-averaged on actuator disk is zero." <<
nl 204 <<
"Please check if the initial U field is zero." 209 const scalar Ct = sink_*UvsCtPtr_->value(magUref);
210 const scalar
Cp = sink_*UvsCpPtr_->value(magUref);
212 if (
Cp <= VSMALL || Ct <= VSMALL)
215 <<
"Cp and Ct must be greater than zero." <<
nl 216 <<
"Cp = " <<
Cp <<
", Ct = " << Ct
221 const scalar CtStar = Ct*
sqr(magUref/magUdisk);
222 const scalar CpStar =
Cp*
pow3(magUref/magUdisk);
225 const vector diskDir = this->diskDir();
226 const scalar
T = 0.5*rhoRef*diskArea_*
magSqr(Udisk & diskDir)*CtStar;
227 const scalar P = 0.5*rhoRef*diskArea_*
pow3(
mag(Udisk & diskDir))*CpStar;
229 for (
const label celli : cells_)
231 Usource[celli] += (cellsV[celli]/totalV*
T)*diskDir;
236 mesh_.time().timeOutputValue() >= writeFileStart_
237 && mesh_.time().timeOutputValue() <= writeFileEnd_
240 Ostream&
os = file();
241 writeCurrentTime(
os);
244 << Udisk <<
tab << CpStar <<
tab << CtStar <<
tab <<
T <<
tab << P
errorManipArg< error, int > exit(error &err, const int errNo=1)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
error FatalError
Error stream (stdout output on all processes), with additional 'FOAM FATAL ERROR' header text and sta...
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
constexpr char nl
The newline '\n' character (0x0a)
Ostream & endl(Ostream &os)
Add newline and flush stream.
constexpr char tab
The tab '\t' character(0x09)
enum forceMethodType forceMethod_
The type of the force computation method.
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
const volScalarField & Cp
OBJstream os(runTime.globalPath()/outputName)
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
dimensionedScalar pow3(const dimensionedScalar &ds)
Field< vector > vectorField
Specialisation of Field<T> for vector.
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
void reduce(T &value, const BinaryOp &bop, const int tag=UPstream::msgType(), const label comm=UPstream::worldComm)
Reduce inplace (cf. MPI Allreduce) using linear/tree communication schedule.
dimensioned< typename typeOfMag< Type >::type > magSqr(const dimensioned< Type > &dt)
static constexpr const zero Zero
Global zero (0)