37 typename Foam::turbulence::IntegralScaleBox<Type>::kernelType
41 { kernelType::GAUSSIAN,
"Gaussian" },
42 { kernelType::EXPONENTIAL ,
"exponential" }
56 const dictionary&
dict
59 if (
dict.found(coordinateSystem::typeName_()))
61 return coordinateSystem::New
63 p_.patch().boundaryMesh().mesh(),
65 coordinateSystem::typeName_()
81 scalar minMag =
mag(nf[minCmpt]);
82 for (direction cmpt = 1; cmpt < pTraits<vector>::nComponents; ++cmpt)
84 const scalar
s =
mag(nf[cmpt]);
102 new coordSystem::cartesian
119 csysPtr_->localPosition
124 p_.patch().meshPoints()
130 const bool globalReduction =
true;
131 const boundBox bb(localPos, globalReduction);
133 return Vector2D<vector>(bb.span(), bb.min());
141 return Vector2D<scalar>
143 boundingBoxSpan_[1]/n_.x(),
144 boundingBoxSpan_[2]/n_.y()
152 if (!Pstream::master())
157 labelList spans(pTraits<TypeL>::nComponents, label(1));
158 const Vector<label> slice(label(1), n_.x(), n_.y());
159 const TypeL
L(convert(L_));
164 j = pTraits<Type>::nComponents;
167 for (label i = j; i < pTraits<TypeL>::nComponents; ++i)
169 const label slicei = label(i/pTraits<Type>::nComponents);
171 const label
n = ceil(
L[i]);
172 const label twiceN = 4*
n;
174 spans[i] = slice[slicei] + twiceN;
185 if (!Pstream::master())
192 pTraits<TypeL>::nComponents,
196 const TypeL
L(convert(L_));
201 i = pTraits<Type>::nComponents;
204 for (direction dir = i; dir < pTraits<TypeL>::nComponents; ++dir)
208 const label
n = ceil(
L[dir]);
212 const label twiceN = 4*
n;
219 const scalar initElem = -2*scalar(
n);
234 const scalar
C = -0.5*constant::mathematical::pi;
238 const scalar nSqr =
n*
n;
261 if (!Pstream::master())
269 for (direction dir = 0; dir < pTraits<Type>::nComponents; ++dir)
276 *spans_[dir+pTraits<TypeL>::nComponents/3]
277 *spans_[dir+2*(pTraits<TypeL>::nComponents/3)]
280 if (randomSet.size() > 1e8)
283 <<
"Size of random-number set is relatively high:" <<
nl
284 <<
" size = " << randomSet.size() <<
nl
285 <<
" Please consider to use the forward-stepwise method."
295 [&]{ return rndGen_.GaussNormal<scalar>(); }
307 if (!Pstream::master())
313 const label nx = n_.x();
314 const label ny = n_.y();
315 const label
nPoints = (nx + 1)*(ny + 1);
319 for (label j = 0; j <= ny; ++j)
321 for (label i = 0; i <= nx; ++i)
326 boundingBoxMin_[1] + i*delta_.x(),
327 boundingBoxMin_[2] + j*delta_.y()
343 if (!Pstream::master())
349 const label nx = n_.x();
350 const label ny = n_.y();
351 const label
nFaces = nx*ny;
355 for (label j = 0; j < ny; ++j)
357 for (label i = 0; i < nx; ++i)
359 const label
k = j*(nx+1) + i;
360 faces[m] = face({
k,
k+(nx+1),
k+(nx+2),
k+1});
372 if (debug && Pstream::master())
374 const auto& tm = p_.
patch().boundaryMesh().mesh().time();
375 OBJstream
os(tm.path()/
"patch.obj");
376 os.write(patchFaces_, patchPoints_,
false);
385 SubList<face>(patchFaces_, patchFaces_.size()),
394typename Foam::turbulence::IntegralScaleBox<Type>::TypeL
397 const typename Foam::turbulence::IntegralScaleBox<Type>::TypeL&
L
402 const scalar deltaT =
403 p_.patch().boundaryMesh().mesh().time().deltaTValue();
405 for (direction dir = 0; dir < pTraits<Type>::nComponents; ++dir)
411 Ls[dir + Ls.size()/3] /= delta_.x();
412 Ls[dir + 2*(Ls.size()/3)] /= delta_.y();
425 constexpr scalar
c1 = -0.25*constant::mathematical::pi;
436 constexpr scalar
c1 = -0.25*constant::mathematical::pi;
454 constexpr scalar
c2 = -0.5*constant::mathematical::pi;
465 constexpr scalar
c2 = -0.5*constant::mathematical::pi;
480 if (p_.patch().boundaryMesh().mesh().time().isAdjustTimeStep())
482 C1_ = calcC1(convert(L_));
483 C2_ = calcC2(convert(L_));
499 kernelType_(kernelType::GAUSSIAN),
503 boundingBoxSpan_(Zero),
504 boundingBoxMin_(Zero),
527 csysPtr_(
b.csysPtr_.clone()),
528 kernelType_(
b.kernelType_),
532 boundingBoxSpan_(
b.boundingBoxSpan_),
533 boundingBoxMin_(
b.boundingBoxMin_),
538 patchPoints_(
b.patchPoints_),
539 patchFaces_(
b.patchFaces_),
556 csysPtr_(calcCoordinateSystem(
dict)),
559 kernelTypeNames.getOrDefault
569 boundingBoxSpan_(
Zero),
570 boundingBoxMin_(
Zero),
571 L_(
dict.get<TypeL>(
"L")),
577 fsm_(
dict.getOrDefault(
"fsm", false)),
585 <<
"Integral scale set contains a very small input" <<
nl
590 if (
min(n_.
x(), n_.
y()) <= 0)
593 <<
"Number of faces on box inlet plane has non-positive input"
608 csysPtr_(
b.csysPtr_.clone()),
609 kernelType_(
b.kernelType_),
613 boundingBoxSpan_(
b.boundingBoxSpan_),
614 boundingBoxMin_(
b.boundingBoxMin_),
619 patchPoints_(
b.patchPoints_),
620 patchFaces_(
b.patchFaces_),
635 calcCoordinateSystem();
638 if (debug && csysPtr_)
640 Info<<
"Local coordinate system:" <<
nl
641 <<
" - origin = " << csysPtr_->origin() <<
nl
642 <<
" - e1-axis = " << csysPtr_->e1() <<
nl
643 <<
" - e2-axis = " << csysPtr_->e2() <<
nl
644 <<
" - e3-axis = " << csysPtr_->e3() <<
nl <<
endl;
650 boundingBoxSpan_ = bb.
x();
652 boundingBoxMin_ = bb.
y();
655 delta_ = calcDelta();
657 spans_ = calcSpans();
659 kernel_ = calcKernel();
663 patchPoints_ = calcPatchPoints();
665 patchFaces_ = calcPatchFaces();
671 C1_ = calcC1(convert(L_));
673 C2_ = calcC2(convert(L_));
683 for (
direction dir = 0; dir < pTraits<Type>::nComponents; ++dir)
689 const label sliceSpan =
702 for (
direction dir = 0; dir < pTraits<Type>::nComponents; ++dir)
706 const label sliceSpan =
711 for (label i = 0; i < sliceSpan; ++i)
713 slice[i] = rndGen_.GaussNormal<scalar>();
725 for (
direction dir = 0; dir < pTraits<Type>::nComponents; ++dir)
737 const label szkernel1 = kernel1.
size();
738 const label szkernel2 = kernel2.
size();
739 const label szkernel3 = kernel3.
size();
741 const label sz1 = spans_[dir];
744 const label sz23 = sz2*sz3;
745 const label sz123 = sz1*sz23;
747 const label validSlice2 = sz2 - (szkernel2 - 1);
748 const label validSlice3 = sz3 - (szkernel3 - 1);
753 const label filterCentre = label(szkernel2/label(2));
754 const label endIndex = sz2 - filterCentre;
758 for (label i = 0; i < sz1; ++i)
760 for (label j = 0; j < sz3; ++j)
764 for (label
k = filterCentre;
k < endIndex; ++
k)
768 for (label
p = szkernel2 - 1;
p >= 0; --
p, ++q)
770 tmp[i1] += in[i0 + q]*kernel2[
p];
775 i0 += 2*filterCentre;
784 const label filterCentre = label(szkernel3/label(2));
785 const label endIndex = sz3 - filterCentre;
789 for (label i = 0; i < sz1; ++i)
791 const label sl = i*sz23;
793 for (label j = 0; j < sz2; ++j)
798 for (label
k = 0;
k < endIndex - filterCentre; ++
k)
802 for (label
p = szkernel3 - 1, q = 0;
p >= 0; --
p, ++q)
804 tmp[i1] += tmp2[i2 + q*sz2]*kernel3[
p];
809 i1 += (sz2 + filterCentre);
810 i2 += (sz2 + filterCentre);
817 label i1 = (szkernel2 - label(1))/label(2);
818 label i2 = (szkernel2 - label(1))/label(2);
821 for (label i = 0; i < validSlice3; ++i)
823 for (label j = 0; j < validSlice2; ++j)
828 for (label
k = szkernel1 - 1;
k >= 0; --
k)
871 for (
direction dir = 0; dir < pTraits<vector>::nComponents; ++dir)
876 slice_.component(dir)*C1_[dir] +
fld.component(dir)*C2_[dir]
891 os.writeEntryIfDifferent<
bool>(
"fsm",
false, fsm_);
892 os.writeEntry(
"n", n_);
893 os.writeEntry(
"L", L_);
894 os.writeEntry(
"kernelType", kernelTypeNames[kernelType_]);
897 csysPtr_->writeEntry(
os);
static const Foam::dimensionedScalar C("", Foam::dimTemperature, 234.5)
Info<< nl<< "Wrote faMesh in vtk format: "<< writer.output().name()<< nl;}{ vtk::lineWriter writer(aMesh.points(), aMesh.edges(), fileName(aMesh.mesh().time().globalPath()/"finiteArea-edges"));writer.writeGeometry();writer.beginCellData(4);writer.writeProcIDs();{ Field< scalar > fld(faMeshTools::flattenEdgeField(aMesh.magLe(), true))
Enum is a wrapper around a list of names/values that represent particular enumeration (or int) values...
void replace(const direction, const UList< cmptType > &)
Replace a component field of the field.
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
void size(const label n)
Older name for setAddressableSize.
Templated 2D Vector derived from VectorSpace adding construction from 2 components,...
const Cmpt & y() const
Access to the vector y component.
const Cmpt & x() const
Access to the vector x component.
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
virtual bool write()
Write the output fields.
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
A traits class, which is primarily used for primitives.
Tensor of scalars, i.e. Tensor<scalar>.
A class for managing temporary objects.
Class to describe the integral-scale container being used in the turbulentDigitalFilterInletFvPatchFi...
void initialise()
Initialise integral-scale box properties.
void refill()
Add a new integral-scale box slice to the rear of the box.
const primitivePatch & patch()
Return const reference to integral-scale box inlet patch.
void correlate(scalarField &fld)
Apply forward-stepwise correlation for scalar fields.
Field< Type > convolve() const
Embed two-point correlations, i.e. L.
A Vector of values with scalar precision, where scalar is float/double depending on the compilation f...
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
OBJstream os(runTime.globalPath()/outputName)
const labelList nFaces(UPstream::listGatherValues< label >(aMesh.nFaces()))
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))
#define WarningInFunction
Report a warning using Foam::Warning.
const dimensionedScalar c2
Second radiation constant: default SI units: [m.K].
const dimensionedScalar c1
First radiation constant: default SI units: [W/m2].
dimensionedScalar exp(const dimensionedScalar &ds)
List< label > labelList
A List of labels.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &df)
void cmptMin(FieldField< Field, typename FieldField< Field, Type >::cmptType > &cf, const FieldField< Field, Type > &f)
messageStream Info
Information stream (stdout output on master, null elsewhere)
void inplaceRotateList(ListType< DataType > &list, label n)
Inplace reversal of a list using the Reversal Block Swapping algorithm.
vectorField pointField
pointField is a vectorField.
List< scalar > scalarList
A List of scalars.
vector point
Point is a vector.
Ostream & endl(Ostream &os)
Add newline and flush stream.
dimensionedScalar sqrt(const dimensionedScalar &ds)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
Type gAverage(const FieldField< Field, Type > &f)
static constexpr const zero Zero
Global zero (0)
List< face > faceList
A List of faces.
errorManipArg< error, int > exit(error &err, const int errNo=1)
List< scalarList > scalarListList
A List of scalarList.
constexpr char nl
The newline '\n' character (0x0a)
constexpr auto end(C &c) -> decltype(c.end())
Return iterator to the end of the container c.
constexpr auto begin(C &c) -> decltype(c.begin())
Return iterator to the beginning of the container c.
#define forAll(list, i)
Loop across all elements in list.
const vector L(dict.get< vector >("L"))