37 typename Foam::turbulence::IntegralScaleBox<Type>::kernelType
41 { kernelType::GAUSSIAN,
"Gaussian" },
42 { kernelType::EXPONENTIAL ,
"exponential" }
56 const dictionary&
dict 59 return coordinateSystem::NewIfPresent(
dict);
71 scalar minMag =
mag(nf[minCmpt]);
74 const scalar
s =
mag(nf[cmpt]);
92 new coordSystem::cartesian
109 csysPtr_->localPosition
114 p_.patch().meshPoints()
120 const bool globalReduction =
true;
121 const boundBox bb(localPos, globalReduction);
123 return Vector2D<vector>(bb.span(), bb.min());
131 return Vector2D<scalar>
133 boundingBoxSpan_[1]/n_.x(),
134 boundingBoxSpan_[2]/n_.y()
142 if (!Pstream::master())
148 const Vector<label> slice(label(1), n_.x(), n_.y());
149 const TypeL
L(convert(L_));
161 const label
n = ceil(
L[i]);
162 const label twiceN = 4*
n;
164 spans[i] = slice[slicei] + twiceN;
175 if (!Pstream::master())
186 const TypeL
L(convert(L_));
198 const label
n = ceil(
L[dir]);
202 const label twiceN = 4*
n;
209 const scalar initElem = -2*scalar(
n);
228 const scalar nSqr =
n*
n;
251 if (!Pstream::master())
270 if (randomSet.size() > 1e8)
273 <<
"Size of random-number set is relatively high:" <<
nl 274 <<
" size = " << randomSet.size() <<
nl 275 <<
" Please consider to use the forward-stepwise method." 285 [&]{
return rndGen_.GaussNormal<scalar>(); }
297 if (!Pstream::master())
303 const label nx = n_.x();
304 const label ny = n_.y();
305 const label
nPoints = (nx + 1)*(ny + 1);
309 for (label j = 0; j <= ny; ++j)
311 for (label i = 0; i <= nx; ++i)
316 boundingBoxMin_[1] + i*delta_.x(),
317 boundingBoxMin_[2] + j*delta_.y()
333 if (!Pstream::master())
339 const label nx = n_.x();
340 const label ny = n_.y();
341 const label nFaces = nx*ny;
345 for (label j = 0; j < ny; ++j)
347 for (label i = 0; i < nx; ++i)
349 const label
k = j*(nx+1) + i;
350 faces[m] = face({
k,
k+(nx+1),
k+(nx+2),
k+1});
362 if (
debug && Pstream::master())
364 const auto& tm = p_.
patch().boundaryMesh().mesh().time();
365 OBJstream
os(tm.path()/
"patch.obj");
366 os.write(patchFaces_, patchPoints_,
false);
375 SubList<face>(patchFaces_, patchFaces_.size()),
384 typename Foam::turbulence::IntegralScaleBox<Type>::TypeL
387 const typename Foam::turbulence::IntegralScaleBox<Type>::TypeL&
L 392 const scalar deltaT =
393 p_.patch().boundaryMesh().mesh().time().deltaTValue();
401 Ls[dir + Ls.size()/3] /= delta_.x();
402 Ls[dir + 2*(Ls.size()/3)] /= delta_.y();
470 if (p_.patch().boundaryMesh().mesh().time().isAdjustTimeStep())
472 C1_ = calcC1(convert(L_));
473 C2_ = calcC2(convert(L_));
489 kernelType_(kernelType::GAUSSIAN),
493 boundingBoxSpan_(
Zero),
494 boundingBoxMin_(
Zero),
517 csysPtr_(
b.csysPtr_.clone()),
518 kernelType_(
b.kernelType_),
522 boundingBoxSpan_(
b.boundingBoxSpan_),
523 boundingBoxMin_(
b.boundingBoxMin_),
528 patchPoints_(
b.patchPoints_),
529 patchFaces_(
b.patchFaces_),
546 csysPtr_(calcCoordinateSystem(
dict)),
549 kernelTypeNames.getOrDefault
559 boundingBoxSpan_(
Zero),
560 boundingBoxMin_(
Zero),
567 fsm_(
dict.getOrDefault(
"fsm", false)),
575 <<
"Integral scale set contains a very small input" <<
nl 580 if (
min(n_.
x(), n_.
y()) <= 0)
583 <<
"Number of faces on box inlet plane has non-positive input" 593 const IntegralScaleBox&
b 598 csysPtr_(
b.csysPtr_.clone()),
599 kernelType_(
b.kernelType_),
603 boundingBoxSpan_(
b.boundingBoxSpan_),
604 boundingBoxMin_(
b.boundingBoxMin_),
609 patchPoints_(
b.patchPoints_),
610 patchFaces_(
b.patchFaces_),
625 calcCoordinateSystem();
628 if (
debug && csysPtr_)
630 Info<<
"Local coordinate system:" <<
nl 631 <<
" - origin = " << csysPtr_->origin() <<
nl 632 <<
" - e1-axis = " << csysPtr_->e1() <<
nl 633 <<
" - e2-axis = " << csysPtr_->e2() <<
nl 634 <<
" - e3-axis = " << csysPtr_->e3() <<
nl <<
endl;
638 const Vector2D<vector> bb(calcBoundBox());
640 boundingBoxSpan_ = bb.x();
642 boundingBoxMin_ = bb.y();
645 delta_ = calcDelta();
647 spans_ = calcSpans();
649 kernel_ = calcKernel();
653 patchPoints_ = calcPatchPoints();
655 patchFaces_ = calcPatchFaces();
661 C1_ = calcC1(convert(L_));
663 C2_ = calcC2(convert(L_));
679 const label sliceSpan =
696 const label sliceSpan =
701 for (label i = 0; i < sliceSpan; ++i)
703 slice[i] = rndGen_.GaussNormal<scalar>();
713 Field<Type> outFld(n_.x()*n_.y(),
Zero);
719 Field<scalar> out(n_.x()*n_.y(),
Zero);
727 const label szkernel1 = kernel1.
size();
728 const label szkernel2 = kernel2.size();
729 const label szkernel3 = kernel3.size();
731 const label sz1 = spans_[dir];
734 const label sz23 = sz2*sz3;
735 const label sz123 = sz1*sz23;
737 const label validSlice2 = sz2 - (szkernel2 - 1);
738 const label validSlice3 = sz3 - (szkernel3 - 1);
743 const label filterCentre = label(szkernel2/label(2));
744 const label endIndex = sz2 - filterCentre;
748 for (label i = 0; i < sz1; ++i)
750 for (label j = 0; j < sz3; ++j)
754 for (label
k = filterCentre;
k < endIndex; ++
k)
758 for (label
p = szkernel2 - 1;
p >= 0; --
p, ++q)
760 tmp[i1] += in[i0 + q]*kernel2[
p];
765 i0 += 2*filterCentre;
774 const label filterCentre = label(szkernel3/label(2));
775 const label endIndex = sz3 - filterCentre;
779 for (label i = 0; i < sz1; ++i)
781 const label sl = i*sz23;
783 for (label j = 0; j < sz2; ++j)
788 for (label
k = 0;
k < endIndex - filterCentre; ++
k)
792 for (label
p = szkernel3 - 1, q = 0;
p >= 0; --
p, ++q)
794 tmp[i1] += tmp2[i2 + q*sz2]*kernel3[
p];
799 i1 += (sz2 + filterCentre);
800 i2 += (sz2 + filterCentre);
807 label i1 = (szkernel2 - label(1))/label(2);
808 label i2 = (szkernel2 - label(1))/label(2);
811 for (label i = 0; i < validSlice3; ++i)
813 for (label j = 0; j < validSlice2; ++j)
818 for (label
k = szkernel1 - 1;
k >= 0; --
k)
820 sum += tmp[i1]*kernel1[
k];
830 outFld.replace(dir, out);
866 slice_.component(dir)*C1_[dir] +
fld.component(dir)*C2_[dir]
881 os.writeEntryIfDifferent<
bool>(
"fsm",
false, fsm_);
882 os.writeEntry(
"n", n_);
883 os.writeEntry(
"L", L_);
884 os.writeEntry(
"kernelType", kernelTypeNames[kernelType_]);
887 csysPtr_->writeEntry(
os);
List< ReturnType > get(const UPtrList< T > &list, const AccessOp &aop)
List of values generated by applying the access operation to each list item.
List< scalar > scalarList
List of scalar.
void size(const label n)
Older name for setAddressableSize.
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &f1)
Field< Type > convolve() const
Embed two-point correlations, i.e. L.
errorManipArg< error, int > exit(error &err, const int errNo=1)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
const Cmpt & y() const noexcept
Access to the vector y component.
A list of keyword definitions, which are a keyword followed by a number of values (eg...
const vector L(dict.get< vector >("L"))
dimensionedSymmTensor sqr(const dimensionedVector &dv)
constexpr char nl
The newline '\n' character (0x0a)
const dimensionedScalar c2
Second radiation constant: default SI units: [m.K].
dimensionedScalar sqrt(const dimensionedScalar &ds)
Ostream & endl(Ostream &os)
Add newline and flush stream.
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
A traits class, which is primarily used for primitives and vector-space.
::Foam::direction nComponents(const expressions::valueTypeCode) noexcept
The number of components associated with given valueTypeCode.
label k
Boltzmann constant.
Class to describe the integral-scale container being used in the turbulentDigitalFilterInletFvPatchFi...
#define forAll(list, i)
Loop across all elements in list.
List< face > faceList
List of faces.
vectorField pointField
pointField is a vectorField.
void refill()
Add a new integral-scale box slice to the rear of the box.
dimensionedScalar exp(const dimensionedScalar &ds)
const dimensionedScalar b
Wien displacement law constant: default SI units: [m.K].
List< scalarList > scalarListList
List of scalarList.
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
void shift()
Discard current time-step integral-scale box slice (the closest to the patch) by shifting from the ba...
void initialise()
Initialise integral-scale box properties.
constexpr scalar pi(M_PI)
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
void correlate(scalarField &fld)
Apply forward-stepwise correlation for scalar fields.
void write(Ostream &) const
Write integral-scale box settings.
An Ostream is an abstract base class for all output systems (streams, files, token lists...
A Vector of values with scalar precision, where scalar is float/double depending on the compilation f...
int debug
Static debugging option.
constexpr auto end(C &c) -> decltype(c.end())
Return iterator to the end of the container c.
OBJstream os(runTime.globalPath()/outputName)
const Cmpt & x() const noexcept
Access to the vector x component.
PrimitivePatch< SubList< face >, const pointField & > primitivePatch
A PrimitivePatch with a SubList addressing for the faces, const reference for the point field...
gmvFile<< "tracers "<< particles.size()<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().x()<< ' ';}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().y()<< ' ';}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().z()<< ' ';}gmvFile<< nl;for(const word &name :lagrangianScalarNames){ IOField< scalar > fld(IOobject(name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
vector point
Point is a vector.
void inplaceRotateList(ListType< DataType > &list, label n)
Inplace reversal of a list using the Reversal Block Swapping algorithm.
#define WarningInFunction
Report a warning using Foam::Warning.
Enum is a wrapper around a list of names/values that represent particular enumeration (or int) values...
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
Type gAverage(const FieldField< Field, Type > &f)
const dimensionedScalar c1
First radiation constant: default SI units: [W/m2].
messageStream Info
Information stream (stdout output on master, null elsewhere)
List< label > labelList
A List of labels.
void cmptMin(FieldField< Field, typename FieldField< Field, Type >::cmptType > &cf, const FieldField< Field, Type > &f)
IntegralScaleBox(const fvPatch &p)
Construct from patch.
gmvFile<< "tracers "<< particles.size()<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().x()<< " ";}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().y()<< " ";}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().z()<< " ";}gmvFile<< nl;forAll(lagrangianScalarNames, i){ word name=lagrangianScalarNames[i];IOField< scalar > s(IOobject(name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
Templated 2D Vector derived from VectorSpace adding construction from 2 components, element access using x() and y() member functions and the inner-product (dot-product).
constexpr auto begin(C &c) -> decltype(c.begin())
Return iterator to the beginning of the container c.
const primitivePatch & patch()
Return const reference to integral-scale box inlet patch.
IOerror FatalIOError
Error stream (stdout output on all processes), with additional 'FOAM FATAL IO ERROR' header text and ...
static constexpr const zero Zero
Global zero (0)