Go to the documentation of this file.
41 namespace functionObjects
68 const scalar kappaNorm
83 Log <<
"Computing wave numbers" <<
endl;
88 kappa[kappai] = kappai*kappaNorm;
91 Log <<
"Computing energy spectrum" <<
endl;
97 label i = round((Cc.
x() - 0.5*deltaC.
x())/(deltaC.
x())*(
N.
x() - 1));
98 label j = round((Cc.
y() - 0.5*deltaC.
y())/(deltaC.
y())*(
N.
y() - 1));
99 label
k = round((Cc.
z() - 0.5*deltaC.
z())/(deltaC.
z())*(
N.
z() - 1));
100 label kappai = round(
Foam::sqrt(scalar(i*i + j*j +
k*
k)));
102 E[kappai] += Ec[celli];
107 Log <<
"Writing spectrum" <<
endl;
114 os <<
kappa[kappai] <<
tab << E[kappai] <<
nl;
130 cellAddr_(mesh_.nCells()),
148 dict.readIfPresent(
"U", UName_);
153 const cell&
c = mesh_.cells()[0];
157 const face&
f = mesh_.faces()[
c[facei]];
158 cellBb.
add(mesh_.points(),
f);
173 vector expectedMax(N_.x()*cellDx.x(), N_.y()*cellDx.y(), N_.z()*cellDx.z());
177 if (
mag(relativeSize[i] - 1) > 1
e-3)
180 <<
name() <<
" function object is only appropriate for "
181 <<
"isotropic structured IJK meshes. Mesh extents: " <<
L
182 <<
", computed IJK mesh extents: " << expectedMax
186 Log <<
"Mesh extents (deltax,deltay,deltaz): " <<
L <<
endl;
187 Log <<
"Number of cells (Nx,Ny,Nz): " << N_ <<
endl;
193 deltaC_ = cMax - c0_;
197 label i = round((
C[celli].
x() - c0_.x())/(deltaC_.x())*(N_.x() - 1));
198 label j = round((
C[celli].
y() - c0_.y())/(deltaC_.y())*(N_.y() - 1));
199 label
k = round((
C[celli].z() - c0_.z())/(deltaC_.z())*(N_.z() - 1));
201 cellAddr_[celli] =
k + j*N_.y() + i*N_.y()*N_.z();
228 toProc << Uc <<
C << cellAddr_;
244 fromProc >> Up >>
Cp >> cellAddrp;
249 calcAndWriteSpectrum(Uijk, Cijk, c0_, deltaC_, N_, kappaNorm_);
258 calcAndWriteSpectrum(Uijk, Cijk, c0_, deltaC_, N_, kappaNorm_);
const vector L(dict.get< vector >("L"))
const Cmpt & x() const
Access to the vector x component.
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
A class for handling words, derived from Foam::string.
void calcAndWriteSpectrum(const vectorField &U, const vectorField &C, const vector &c0, const vector &deltaC, const Vector< int > &N, const scalar kappaNorm)
Calculate and write the spectrum.
Output inter-processor communications stream operating on external buffer.
T returnReduce(const T &Value, const BinaryOp &bop, const int tag=Pstream::msgType(), const label comm=UPstream::worldComm)
static constexpr const zero Zero
Global zero (0)
static const boundBox invertedBox
A large inverted boundBox: min/max == +/- ROOTVGREAT.
static label nProcs(const label communicator=0)
Number of processes in parallel run.
virtual bool read(const dictionary &)
Read the field min/max data.
bool read(const char *buf, int32_t &val)
Same as readInt32.
Buffers for inter-processor communications streams (UOPstream, UIPstream).
static bool & parRun()
Is this a parallel run?
autoPtr< surfaceVectorField > Uf
Cmpt cmptProduct(const VectorSpace< Form, Cmpt, Ncmpts > &vs)
Ostream & endl(Ostream &os)
Add newline and flush stream.
const point & max() const
Maximum describing the bounding box.
List< treeBoundBox > meshBb(1, treeBoundBox(boundBox(coarseMesh.points(), false)).extend(rndGen, 1e-3))
const dimensionedScalar kappa
Coulomb constant: default SI units: [N.m2/C2].
Specialization of Foam::functionObject for an Foam::fvMesh, providing a reference to the Foam::fvMesh...
const Cmpt & z() const
Access to the vector z component.
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
const point & min() const
Minimum describing the bounding box.
#define forAll(list, i)
Loop across all elements in list.
dimensioned< typename typeOfMag< Type >::type > magSqr(const dimensioned< Type > &dt)
virtual void writeHeader(Ostream &os, const string &str) const
Write a commented header to stream.
complexField ReComplexField(const UList< scalar > &re)
Create complex field from a list of real (using imag == 0)
void cmptMax(FieldField< Field, typename FieldField< Field, Type >::cmptType > &cf, const FieldField< Field, Type > &f)
virtual bool read(const dictionary &dict)
Read.
virtual bool execute()
Execute, currently does nothing.
virtual void writeFileHeader(Ostream &os)
Output file header information.
word name(const complex &c)
Return string representation of complex.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
void finishedSends(const bool block=true)
Mark all sends as having been done. This will start receives.
constexpr scalar twoPi(2 *M_PI)
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
virtual bool read(const dictionary &dict)
Read optional controls.
Macros for easy insertion into run-time selection tables.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Output to file stream, using an OSstream.
static bool master(const label communicator=0)
Am I the master process.
energySpectrum(const energySpectrum &)=delete
No copy construct.
virtual bool write()
Write the energySpectrum.
dimensioned< Type > cmptDivide(const dimensioned< Type > &, const dimensioned< Type > &)
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
const Cmpt & y() const
Access to the vector y component.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
virtual void writeCommented(Ostream &os, const string &str) const
Write a commented string to stream.
T & ref()
Return reference to the managed object without nullptr checking.
addToRunTimeSelectionTable(functionObject, ObukhovLength, dictionary)
dimensionedScalar sqrt(const dimensionedScalar &ds)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
label k
Boltzmann constant.
const dimensionedScalar e
Elementary charge.
static tmp< complexField > forwardTransform(const tmp< complexField > &field, const UList< int > &nn)
const volScalarField & Cp
defineTypeNameAndDebug(ObukhovLength, 0)
A bounding box defined in terms of min/max extrema points.
const dimensionedScalar c
Speed of light in a vacuum.
A List with indirect addressing.
Base class for writing single files from the function objects.
A face is a list of labels corresponding to mesh vertices.
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
const Vector< label > N(dict.get< Vector< label >>("N"))
Input inter-processor communications stream operating on external buffer.
Graphite solid properties.
A cell is defined as a list of faces with extra functionality.
void add(const boundBox &bb)
Extend to include the second box.