37template<
class RhoFieldType>
40 const RhoFieldType&
rho,
44 const bool divideVolume,
53 scalar AOAmin = GREAT;
54 scalar AOAmax = -GREAT;
61 if (
area_[i] > ROOTVSMALL)
63 const label celli =
cells_[i];
65 const scalar radius =
x_[i].x();
96 scalar alphaGeom = thetag[i] + twist;
128 scalar Cd = invDr*(Cd2 - Cd1) + Cd1;
129 scalar Cl = invDr*(Cl2 - Cl1) + Cl1;
148 force[celli] =
transform(Rcyl, localForce);
152 force[celli] /=
V[celli];
165 <<
" min/max(AOA) = " <<
radToDeg(AOAmin) <<
", "
167 <<
" Effective drag = " << dragEff <<
nl
168 <<
" Effective lift = " << liftEff <<
endl;
183 if (mesh_.time().writeTime() || writeNow)
190 mesh_.time().timeName(),
199 auto&
field = tfield.ref().primitiveFieldRef();
201 if (cells_.size() != values.size())
209 const label celli = cells_[i];
210 field[celli] = values[i];
Generic GeometricField class.
Defines the attributes of an object for which implicit objectRegistry management is supported,...
A 1D array of objects of type <T>, where the size of the vector is known and used for subscript bound...
static autoPtr< Time > New()
Construct (dummy) Time - no functionObjects or libraries.
const Cmpt & z() const
Access to the vector z component.
const Cmpt & y() const
Access to the vector y component.
const Cmpt & x() const
Access to the vector x component.
const List< label > & profileID() const
Return const access to the profile ID list.
virtual void interpolate(const scalar radius, scalar &twist, scalar &chord, label &i1, label &i2, scalar &invDr) const
Return the twist and chord for a given radius.
virtual tensor R(const point &global) const
Generic dimensioned Type class.
const volVectorField & C() const
Return cell centres as volVectorField.
const DimensionedField< scalar, volMesh > & V() const
Return cell volumes.
scalar V() const noexcept
Return const access to the total cell volume.
labelList cells_
Set of cells to apply source to.
const fvMesh & mesh_
Reference to the mesh database.
scalar omega_
Rotational speed [rad/s].
void calculate(const RhoFieldType &rho, const vectorField &U, const scalarField &thetag, vectorField &force, const bool divideVolume=true, const bool output=true) const
Calculate forces.
scalar rhoRef_
Reference density for incompressible case.
List< scalar > area_
Area [m2].
scalar tipEffect_
Tip effect [0-1].
void writeField(const word &name, const List< Type > &values, const bool writeNow=false) const
Helper function to write rotor values.
coordSystem::cylindrical coordSys_
Rotor local cylindrical coordinate system (r-theta-z)
scalar rMax_
Maximum radius.
List< point > x_
Cell centre positions in local rotor frame.
bladeModel blade_
Blade data.
List< tensor > Rcone_
Rotation tensor for flap angle.
label nBlades_
Number of blades.
autoPtr< tensorField > Rcyl_
Cached rotation tensors for cylindrical coordinates.
profileModelList profiles_
Profile data.
A class for handling words, derived from Foam::string.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
volScalarField alphaEff("alphaEff", turbulence->nu()/Pr+alphat)
constexpr scalar pi(M_PI)
constexpr scalar twoPi(2 *M_PI)
Different types of constants.
dimensionSet invTransform(const dimensionSet &ds)
Return the argument; transformations do not change the dimensions.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
const dimensionSet dimless
Dimensionless.
dimensionSet transform(const dimensionSet &ds)
Return the argument; transformations do not change the dimensions.
messageStream Info
Information stream (stdout output on master, null elsewhere)
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
dimensionedScalar atan2(const dimensionedScalar &x, const dimensionedScalar &y)
Ostream & endl(Ostream &os)
Add newline and flush stream.
constexpr scalar radToDeg() noexcept
Multiplication factor for radians to degrees conversion.
void reduce(const List< UPstream::commsStruct > &comms, T &value, const BinaryOp &bop, const int tag, const label comm)
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
errorManip< error > abort(error &err)
static Ostream & output(Ostream &os, const IntRange< T > &range)
dimensionedScalar neg(const dimensionedScalar &ds)
static constexpr const zero Zero
Global zero (0)
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
dimensioned< typename typeOfMag< Type >::type > magSqr(const dimensioned< Type > &dt)
constexpr char nl
The newline '\n' character (0x0a)
volScalarField pDyn(IOobject("pDyn", runTime.timeName(), mesh, IOobject::NO_READ, IOobject::AUTO_WRITE), mesh, dimensionedScalar(dimPressure, Zero))
#define forAll(list, i)
Loop across all elements in list.
Unit conversion functions.