targetCoeffTrim Class Reference

Trim model where the operating characteristics of rotor (e.g. blade pitch angle) can vary to reach a specified target thrust and torque. More...

Inheritance diagram for targetCoeffTrim:
[legend]
Collaboration diagram for targetCoeffTrim:
[legend]

Public Member Functions

 TypeName ("targetCoeffTrim")
 Run-time type information. More...
 
 targetCoeffTrim (const fv::rotorDiskSource &rotor, const dictionary &dict)
 Constructor from rotor and dictionary. More...
 
 targetCoeffTrim (const targetCoeffTrim &)=delete
 No copy construct. More...
 
void operator= (const targetCoeffTrim &)=delete
 No copy assignment. More...
 
virtual ~targetCoeffTrim ()=default
 Destructor. More...
 
void read (const dictionary &dict)
 Read. More...
 
virtual tmp< scalarFieldthetag () const
 Return the geometric angle of attack [rad]. More...
 
virtual void correct (const vectorField &U, vectorField &force)
 Correct the model. More...
 
virtual void correct (const volScalarField rho, const vectorField &U, vectorField &force)
 Correct the model for compressible flow. More...
 
template<class RhoFieldType >
Foam::vector calcCoeffs (const RhoFieldType &rho, const vectorField &U, const scalarField &thetag, vectorField &force) const
 
- Public Member Functions inherited from trimModel
 TypeName ("trimModel")
 Run-time type information. More...
 
 declareRunTimeSelectionTable (autoPtr, trimModel, dictionary,(const fv::rotorDiskSource &rotor, const dictionary &dict),(rotor, dict))
 
 trimModel (const fv::rotorDiskSource &rotor, const dictionary &dict, const word &name)
 Construct from components. More...
 
virtual ~trimModel ()=default
 Destructor. More...
 
virtual void read (const dictionary &dict)
 Read. More...
 
virtual tmp< scalarFieldthetag () const =0
 Return the geometric angle of attack [rad]. More...
 
virtual void correct (const vectorField &U, vectorField &force)=0
 Correct the model. More...
 
virtual void correct (const volScalarField rho, const vectorField &U, vectorField &force)=0
 Correct the model for compressible flow. More...
 

Protected Member Functions

template<class RhoFieldType >
vector calcCoeffs (const RhoFieldType &rho, const vectorField &U, const scalarField &alphag, vectorField &force) const
 Calculate the rotor force and moment coefficients vector. More...
 
template<class RhoFieldType >
void correctTrim (const RhoFieldType &rho, const vectorField &U, vectorField &force)
 Correct the model. More...
 

Protected Attributes

label calcFrequency_
 Number of iterations between calls to 'correct'. More...
 
bool useCoeffs_
 Flag to indicate whether to solve coeffs (true) or forces (false) More...
 
vector target_
 Target coefficient vector (thrust force, roll moment, pitch moment) More...
 
vector theta_
 Pitch angles (collective, roll, pitch) [rad]. More...
 
label nIter_
 Maximum number of iterations in trim routine. More...
 
scalar tol_
 Convergence tolerance. More...
 
scalar relax_
 Under-relaxation coefficient. More...
 
scalar dTheta_
 Perturbation angle used to determine jacobian. More...
 
scalar alpha_
 Coefficient to allow for conversion between US and EU definitions. More...
 
- Protected Attributes inherited from trimModel
const fv::rotorDiskSourcerotor_
 Reference to the rotor source model. More...
 
const word name_
 Name of model. More...
 
dictionary coeffs_
 Coefficients dictionary. More...
 

Additional Inherited Members

- Static Public Member Functions inherited from trimModel
static autoPtr< trimModelNew (const fv::rotorDiskSource &rotor, const dictionary &dict)
 Return a reference to the selected trim model. More...
 

Detailed Description

Trim model where the operating characteristics of rotor (e.g. blade pitch angle) can vary to reach a specified target thrust and torque.

Solves:

\[ c^{old} + J.d(\theta) = c^{target} \]

where

\( n \) = Time level
\( c \) = Coefficient vector (thrust force, roll moment, pitch moment)
\( \theta \) = Pitch angle vector (collective, roll, pitch)
\( J \) = Jacobian [3x3] matrix

The trimmed pitch angles are found via solving the above with a Newton-Raphson iterative method. The solver tolerance can be user-input, using the 'tol' entry.

If coefficients are requested (useCoeffs = true), the force and moments are normalised using:

\[ c = \frac{F}{\alpha \pi \rho \omega^2 R^4} \]

and

\[ c = \frac{M}{\alpha \pi \rho \omega^2 R^5} \]

where

\( F \) = Force
\( M \) = Moment
\( \alpha \) = User-input conversion coefficient (default = 1)
\( \rho \) = Fluid desity
\( \omega \) = Rotor angular velocity
\( \pi \) = Mathematical pi
\( R \) = Rotor radius
Usage
Minimal example by using constant/fvOptions: rotorDiskSource1 { Mandatory/Optional (inherited) entries ...

Mandatory entries (runtime modifiable) trimModel targetForce;

targetForceCoeffs { Conditional mandatory entries (runtime modifiable)

when trimModel=targetForce target { Mandatory entries (runtime modifiable) thrust 0.5; pitch 0.5; roll 0.5;

Optional entries (runtime modifiable) useCoeffs true; }

pitchAngles { theta0Ini 0.1; theta1cIni 0.1; theta1sIni 0.1; }

calcFrequency 1;

Optional entries (runtime modifiable) nIter 50; tol 1e-8; relax 1; dTheta 0.1; alpha 1.0; } }

where the entries mean:

Property Description Type Reqd Dflt
calcFrequency

Number of iterations between calls to 'correct'

label yes -
useCoeffs

Flag to indicate whether to solve coeffs

(true) or forces (false)

bool no true
thrust Target thrust (coefficient) scalar yes -
pitch Target pitch (coefficient) scalar yes -
roll Target roll moment (coefficient) scalar yes -
theta0Ini Initial pitch angle [deg] scalar yes -
theta1cIni

Initial lateral pitch angle (cos coeff) [deg]

scalar yes -
theta1sIni

Initial longitudinal pitch angle (sin coeff) [deg]

scalar yes -
nIter

Maximum number of iterations in trim routine

label no 50
tol Convergence tolerance in trim routine scalar no 1e-8
relax Relaxation factor in trim routine scalar no 1
dTheta

Perturbation angle used to determine Jacobian [deg]

scalar no 0.1
alpha

Coefficient to allow for conversion between US

and EU definitions

scalar no 1
See also
Source files

Definition at line 313 of file targetCoeffTrim.H.

Constructor & Destructor Documentation

◆ targetCoeffTrim() [1/2]

targetCoeffTrim ( const fv::rotorDiskSource rotor,
const dictionary dict 
)

Constructor from rotor and dictionary.

Definition at line 195 of file targetCoeffTrim.C.

References dict, and targetCoeffTrim::read().

Here is the call graph for this function:

◆ targetCoeffTrim() [2/2]

targetCoeffTrim ( const targetCoeffTrim )
delete

No copy construct.

◆ ~targetCoeffTrim()

virtual ~targetCoeffTrim ( )
virtualdefault

Destructor.

Member Function Documentation

◆ calcCoeffs() [1/2]

vector calcCoeffs ( const RhoFieldType &  rho,
const vectorField U,
const scalarField alphag,
vectorField force 
) const
protected

Calculate the rotor force and moment coefficients vector.

◆ correctTrim()

void correctTrim ( const RhoFieldType &  rho,
const vectorField U,
vectorField force 
)
protected

Correct the model.

Definition at line 101 of file targetCoeffTrim.C.

References Foam::endl(), Foam::Info, Foam::inv(), Foam::mag(), Foam::nl, Foam::radToDeg(), rho, Foam::type(), U, and Foam::Zero.

Here is the call graph for this function:

◆ TypeName()

TypeName ( "targetCoeffTrim"  )

Run-time type information.

◆ operator=()

void operator= ( const targetCoeffTrim )
delete

No copy assignment.

◆ read()

void read ( const dictionary dict)
virtual

Read.

Reimplemented from trimModel.

Definition at line 218 of file targetCoeffTrim.C.

References Foam::degToRad(), dict, dictionary::get(), dictionary::getOrDefault(), kEpsilonLopesdaCosta< BasicTurbulenceModel >::read(), and dictionary::readEntry().

Referenced by targetCoeffTrim::targetCoeffTrim().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ thetag()

Foam::tmp< Foam::scalarField > thetag ( ) const
virtual

Return the geometric angle of attack [rad].

Implements trimModel.

Definition at line 254 of file targetCoeffTrim.C.

References Foam::cos(), forAll, Time::New(), psi, Foam::sin(), and x.

Referenced by targetCoeffTrim::calcCoeffs().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ correct() [1/2]

void correct ( const vectorField U,
vectorField force 
)
virtual

Correct the model.

Implements trimModel.

Definition at line 271 of file targetCoeffTrim.C.

References U.

◆ correct() [2/2]

void correct ( const volScalarField  rho,
const vectorField U,
vectorField force 
)
virtual

Correct the model for compressible flow.

Implements trimModel.

Definition at line 281 of file targetCoeffTrim.C.

References rho, and U.

◆ calcCoeffs() [2/2]

Foam::vector calcCoeffs ( const RhoFieldType &  rho,
const vectorField U,
const scalarField thetag,
vectorField force 
) const

Member Data Documentation

◆ calcFrequency_

label calcFrequency_
protected

Number of iterations between calls to 'correct'.

Definition at line 322 of file targetCoeffTrim.H.

◆ useCoeffs_

bool useCoeffs_
protected

Flag to indicate whether to solve coeffs (true) or forces (false)

Definition at line 325 of file targetCoeffTrim.H.

Referenced by targetCoeffTrim::calcCoeffs().

◆ target_

vector target_
protected

Target coefficient vector (thrust force, roll moment, pitch moment)

Definition at line 328 of file targetCoeffTrim.H.

◆ theta_

vector theta_
protected

Pitch angles (collective, roll, pitch) [rad].

Definition at line 331 of file targetCoeffTrim.H.

◆ nIter_

label nIter_
protected

Maximum number of iterations in trim routine.

Definition at line 334 of file targetCoeffTrim.H.

◆ tol_

scalar tol_
protected

Convergence tolerance.

Definition at line 337 of file targetCoeffTrim.H.

◆ relax_

scalar relax_
protected

Under-relaxation coefficient.

Definition at line 340 of file targetCoeffTrim.H.

◆ dTheta_

scalar dTheta_
protected

Perturbation angle used to determine jacobian.

Definition at line 343 of file targetCoeffTrim.H.

◆ alpha_

scalar alpha_
protected

Coefficient to allow for conversion between US and EU definitions.

Definition at line 346 of file targetCoeffTrim.H.

Referenced by targetCoeffTrim::calcCoeffs().


The documentation for this class was generated from the following files: