diff --git a/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/PatchInjection/PatchInjection.C b/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/PatchInjection/PatchInjection.C index 934f1343cf..6c1cc715ad 100644 --- a/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/PatchInjection/PatchInjection.C +++ b/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/PatchInjection/PatchInjection.C @@ -29,6 +29,17 @@ License #include "PatchInjection.H" #include "distributionModel.H" +// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * // + +template +const Foam::Enum::velocityType> +Foam::PatchInjection::velocityTypeNames_ +({ + { vtFixedValue, "fixedValue" }, + { vtPatchValue, "patchValue" }, + { vtZeroGradient, "zeroGradient" }, +}); + // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // template @@ -46,7 +57,21 @@ Foam::PatchInjection::PatchInjection ( this->coeffDict().getScalar("parcelsPerSecond") ), - U0_(this->coeffDict().lookup("U0")), + velocityType_ + ( + velocityTypeNames_.getOrDefault + ( + "velocityType", + this->coeffDict(), + vtFixedValue + ) + ), + U0_ + ( + velocityType_ == vtFixedValue + ? this->coeffDict().template get("U0") + : Zero + ), flowRateProfile_ ( Function1::New @@ -63,7 +88,9 @@ Foam::PatchInjection::PatchInjection this->coeffDict().subDict("sizeDistribution"), owner.rndGen() ) - ) + ), + currentParceli_(-1), + currentFacei_(-1) { // Convert from user time to reduce the number of time conversion calls const Time& time = owner.db().time(); @@ -87,16 +114,12 @@ Foam::PatchInjection::PatchInjection patchInjectionBase(im), duration_(im.duration_), parcelsPerSecond_(im.parcelsPerSecond_), + velocityType_(im.velocityType_), U0_(im.U0_), flowRateProfile_(im.flowRateProfile_.clone()), - sizeDistribution_(im.sizeDistribution_.clone()) -{} - - -// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * // - -template -Foam::PatchInjection::~PatchInjection() + sizeDistribution_(im.sizeDistribution_.clone()), + currentParceli_(im.currentParceli_), + currentFacei_(im.currentFacei_) {} @@ -167,16 +190,18 @@ Foam::scalar Foam::PatchInjection::volumeToInject template void Foam::PatchInjection::setPositionAndCell ( - const label, - const label, - const scalar, + const label parcelI, + const label nParcels, + const scalar time, vector& position, label& cellOwner, label& tetFacei, label& tetPti ) { - patchInjectionBase::setPositionAndCell + currentParceli_ = parcelI; + + currentFacei_ = patchInjectionBase::setPositionAndCell ( this->owner().mesh(), this->owner().rndGen(), @@ -191,16 +216,74 @@ void Foam::PatchInjection::setPositionAndCell template void Foam::PatchInjection::setProperties ( - const label, - const label, - const scalar, + const label parcelI, + const label nParcels, + const scalar time, typename CloudType::parcelType& parcel ) { - // set particle velocity - parcel.U() = U0_; + // Set particle velocity + switch (velocityType_) + { + case vtFixedValue: + { + parcel.U() = U0_; + break; + } + case vtPatchValue: + { + if (parcelI != currentParceli_) + { + WarningInFunction + << "Synchronisation problem: " + << "attempting to set injected parcel " << parcelI + << " properties using cached parcel " << currentParceli_ + << " properties" << endl; + } - // set particle diameter + const label patchFacei = currentFacei_; + + if (patchFacei < 0) + { + FatalErrorInFunction + << "Unable to set parcel velocity using patch value " + << "due to missing face index: patchFacei=" << patchFacei + << abort(FatalError); + } + + const volVectorField& U = this->owner().U(); + const label patchi = this->patchId_; + parcel.U() = U.boundaryField()[patchi][patchFacei]; + break; + } + case vtZeroGradient: + { + const label celli = parcel.cell(); + + if (celli < 0) + { + FatalErrorInFunction + << "Unable to set parcel velocity using zeroGradient " + << "due to missing cell index" + << abort(FatalError); + } + + const volVectorField& U = this->owner().U(); + parcel.U() = U[celli]; + break; + } + default: + { + FatalErrorInFunction + << "Unhandled velocityType " + << velocityTypeNames_[velocityType_] + << ". Available options are:" + << velocityTypeNames_.sortedToc() + << abort(FatalError); + } + } + + // Set particle diameter parcel.d() = sizeDistribution_->sample(); } diff --git a/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/PatchInjection/PatchInjection.H b/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/PatchInjection/PatchInjection.H index 25c1a007cf..517ebf7e02 100644 --- a/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/PatchInjection/PatchInjection.H +++ b/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/PatchInjection/PatchInjection.H @@ -72,6 +72,24 @@ class PatchInjection public InjectionModel, public patchInjectionBase { +public: + + // Public Data Types + + //- Velocity type enumeration + enum velocityType + { + vtFixedValue, //!< User supplied fixed value + vtPatchValue, //!< Patch face values + vtZeroGradient //!< Patch internal cell values + }; + + //- Velocity type names + static const Enum velocityTypeNames_; + + +private: + // Private data //- Injection duration [s] @@ -80,7 +98,11 @@ class PatchInjection //- Number of parcels to introduce per second [] const label parcelsPerSecond_; + //- Velocity type + const velocityType velocityType_; + //- Initial parcel velocity [m/s] + // Note: Only used for vtFixedValue const vector U0_; //- Flow rate profile relative to SOI [] @@ -89,6 +111,12 @@ class PatchInjection //- Parcel size distribution model const autoPtr sizeDistribution_; + //- Current parcel being processed + label currentParceli_; + + //- Current face being processed + label currentFacei_; + public: @@ -120,7 +148,7 @@ public: //- Destructor - virtual ~PatchInjection(); + virtual ~PatchInjection() = default; // Member Functions diff --git a/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/PatchInjection/patchInjectionBase.C b/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/PatchInjection/patchInjectionBase.C index 9832c478e3..e3c495ec93 100644 --- a/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/PatchInjection/patchInjectionBase.C +++ b/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/PatchInjection/patchInjectionBase.C @@ -144,7 +144,7 @@ void Foam::patchInjectionBase::updateMesh(const polyMesh& mesh) } -void Foam::patchInjectionBase::setPositionAndCell +Foam::label Foam::patchInjectionBase::setPositionAndCell ( const fvMesh& mesh, const scalar fraction01, @@ -155,6 +155,8 @@ void Foam::patchInjectionBase::setPositionAndCell label& tetPti ) { + label facei = -1; + if (cellOwners_.size() > 0) { // Determine which processor to inject from @@ -177,7 +179,7 @@ void Foam::patchInjectionBase::setPositionAndCell } // Set cellOwner - label facei = triToFace_[trii]; + facei = triToFace_[trii]; cellOwner = cellOwners_[facei]; // Find random point in triangle @@ -261,10 +263,12 @@ void Foam::patchInjectionBase::setPositionAndCell // Dummy position position = pTraits::max; } + + return facei; } -void Foam::patchInjectionBase::setPositionAndCell +Foam::label Foam::patchInjectionBase::setPositionAndCell ( const fvMesh& mesh, Random& rnd, @@ -276,7 +280,7 @@ void Foam::patchInjectionBase::setPositionAndCell { scalar fraction01 = rnd.globalSample01(); - setPositionAndCell + return setPositionAndCell ( mesh, fraction01, diff --git a/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/PatchInjection/patchInjectionBase.H b/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/PatchInjection/patchInjectionBase.H index aad5931ddf..6c2336670a 100644 --- a/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/PatchInjection/patchInjectionBase.H +++ b/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/PatchInjection/patchInjectionBase.H @@ -118,7 +118,8 @@ public: //- Set the injection position and owner cell, tetFace and tetPt // Supply the fraction used to determine the location on the patch - void setPositionAndCell + // Returns the seed patch face index + label setPositionAndCell ( const fvMesh& mesh, const scalar fraction01, @@ -130,7 +131,8 @@ public: ); //- Set the injection position and owner cell, tetFace and tetPt - virtual void setPositionAndCell + // Returns the seed patch face index + virtual label setPositionAndCell ( const fvMesh& mesh, Random& rnd, diff --git a/tutorials/lagrangian/kinematicParcelFoam/spinningDisk/0/U b/tutorials/lagrangian/kinematicParcelFoam/spinningDisk/0/U new file mode 100644 index 0000000000..f03c35985d --- /dev/null +++ b/tutorials/lagrangian/kinematicParcelFoam/spinningDisk/0/U @@ -0,0 +1,56 @@ +/*--------------------------------*- C++ -*----------------------------------*\ +| ========= | | +| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox | +| \\ / O peration | Version: v2106 | +| \\ / A nd | Website: www.openfoam.com | +| \\/ M anipulation | | +\*---------------------------------------------------------------------------*/ +FoamFile +{ + version 2.0; + format ascii; + class volVectorField; + object U; +} +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +dimensions [0 1 -1 0 0 0 0]; + +internalField uniform (0.1 0 0); + +boundaryField +{ + inlet + { + type fixedValue; + value $internalField; + } + + outlet + { + type inletOutlet; + inletValue uniform (0 0 0); + value uniform (0 0 0); + } + + cylinder + { + type rotatingWallVelocity; + origin (0 0 0); + axis (0 0 1); + omega 1; + } + + walls + { + type noSlip; + } + + frontAndBack + { + type empty; + } +} + + +// ************************************************************************* // diff --git a/tutorials/lagrangian/kinematicParcelFoam/spinningDisk/0/p b/tutorials/lagrangian/kinematicParcelFoam/spinningDisk/0/p new file mode 100644 index 0000000000..31a86b8b89 --- /dev/null +++ b/tutorials/lagrangian/kinematicParcelFoam/spinningDisk/0/p @@ -0,0 +1,46 @@ +/*--------------------------------*- C++ -*----------------------------------*\ +| ========= | | +| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox | +| \\ / O peration | Version: v2106 | +| \\ / A nd | Website: www.openfoam.com | +| \\/ M anipulation | | +\*---------------------------------------------------------------------------*/ +FoamFile +{ + version 2.0; + format ascii; + class volScalarField; + object p; +} +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +dimensions [0 2 -2 0 0 0 0]; + +internalField uniform 0; + +boundaryField +{ + inlet + { + type zeroGradient; + } + + outlet + { + type fixedValue; + value uniform 0; + } + + "(walls|cylinder)" + { + type zeroGradient; + } + + frontAndBack + { + type empty; + } +} + + +// ************************************************************************* // diff --git a/tutorials/lagrangian/kinematicParcelFoam/spinningDisk/constant/g b/tutorials/lagrangian/kinematicParcelFoam/spinningDisk/constant/g new file mode 100644 index 0000000000..e340b6cb92 --- /dev/null +++ b/tutorials/lagrangian/kinematicParcelFoam/spinningDisk/constant/g @@ -0,0 +1,21 @@ +/*--------------------------------*- C++ -*----------------------------------*\ +| ========= | | +| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox | +| \\ / O peration | Version: v2106 | +| \\ / A nd | Website: www.openfoam.com | +| \\/ M anipulation | | +\*---------------------------------------------------------------------------*/ +FoamFile +{ + version 2.0; + format ascii; + class uniformDimensionedVectorField; + object g; +} +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +dimensions [0 1 -2 0 0 0 0]; +value (0 0 -9.8); + + +// ************************************************************************* // diff --git a/tutorials/lagrangian/kinematicParcelFoam/spinningDisk/constant/kinematicCloudProperties b/tutorials/lagrangian/kinematicParcelFoam/spinningDisk/constant/kinematicCloudProperties new file mode 100644 index 0000000000..aeae28fbb1 --- /dev/null +++ b/tutorials/lagrangian/kinematicParcelFoam/spinningDisk/constant/kinematicCloudProperties @@ -0,0 +1,108 @@ +/*--------------------------------*- C++ -*----------------------------------*\ +| ========= | | +| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox | +| \\ / O peration | Version: v1912 | +| \\ / A nd | Website: www.openfoam.com | +| \\/ M anipulation | | +\*---------------------------------------------------------------------------*/ +FoamFile +{ + version 2.0; + format ascii; + class dictionary; + object kinematicCloudProperties; +} +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +solution +{ + active true; + coupled yes; + transient yes; + cellValueSourceCorrection no; + maxCo 0.3; + + sourceTerms + { + schemes + { + U semiImplicit 1; + } + } + + interpolationSchemes + { + rho cell; + U cellPoint; + muc cell; + p cell; + } + + integrationSchemes + { + U Euler; + } +} + +constantProperties +{ + rho0 1000; +} + +subModels +{ + particleForces + { + sphereDrag; + gravity; + } + + injectionModels + { + model1 + { + type patchInjection; + massTotal 1; + SOI 0; + parcelBasisType mass; + patch cylinder; + duration 10; + parcelsPerSecond 100; + velocityType patchValue; + //velocityType zeroGradient; + //U0 (-10 0 0); + flowRateProfile constant 1; + sizeDistribution + { + type normal; + normalDistribution + { + expectation 1e-3; + variance 1e-4; + minValue 1e-5; + maxValue 2e-3; + } + } + } + } + + dispersionModel none; + + patchInteractionModel standardWallInteraction; + + stochasticCollisionModel none; + + surfaceFilmModel none; + + standardWallInteractionCoeffs + { + type rebound; + } +} + + +cloudFunctions +{} + + +// ************************************************************************* // diff --git a/tutorials/lagrangian/kinematicParcelFoam/spinningDisk/constant/transportProperties b/tutorials/lagrangian/kinematicParcelFoam/spinningDisk/constant/transportProperties new file mode 100644 index 0000000000..2f4edd39aa --- /dev/null +++ b/tutorials/lagrangian/kinematicParcelFoam/spinningDisk/constant/transportProperties @@ -0,0 +1,24 @@ +/*--------------------------------*- C++ -*----------------------------------*\ +| ========= | | +| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox | +| \\ / O peration | Version: v2106 | +| \\ / A nd | Website: www.openfoam.com | +| \\/ M anipulation | | +\*---------------------------------------------------------------------------*/ +FoamFile +{ + version 2.0; + format ascii; + class dictionary; + object transportProperties; +} +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +transportModel Newtonian; + +nu 1e-05; + +rhoInf 1.2; + + +// ************************************************************************* // diff --git a/tutorials/lagrangian/kinematicParcelFoam/spinningDisk/constant/turbulenceProperties b/tutorials/lagrangian/kinematicParcelFoam/spinningDisk/constant/turbulenceProperties new file mode 100644 index 0000000000..15113f55f9 --- /dev/null +++ b/tutorials/lagrangian/kinematicParcelFoam/spinningDisk/constant/turbulenceProperties @@ -0,0 +1,20 @@ +/*--------------------------------*- C++ -*----------------------------------*\ +| ========= | | +| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox | +| \\ / O peration | Version: v2106 | +| \\ / A nd | Website: www.openfoam.com | +| \\/ M anipulation | | +\*---------------------------------------------------------------------------*/ +FoamFile +{ + version 2.0; + format ascii; + class dictionary; + object turbulenceProperties; +} +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +simulationType laminar; + + +// ************************************************************************* // diff --git a/tutorials/lagrangian/kinematicParcelFoam/spinningDisk/system/blockMeshDict b/tutorials/lagrangian/kinematicParcelFoam/spinningDisk/system/blockMeshDict new file mode 100644 index 0000000000..a38ed26339 --- /dev/null +++ b/tutorials/lagrangian/kinematicParcelFoam/spinningDisk/system/blockMeshDict @@ -0,0 +1,114 @@ +/*--------------------------------*- C++ -*----------------------------------*\ +| ========= | | +| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox | +| \\ / O peration | Version: v2106 | +| \\ / A nd | Website: www.openfoam.com | +| \\/ M anipulation | | +\*---------------------------------------------------------------------------*/ +FoamFile +{ + version 2.0; + format ascii; + class dictionary; + object blockMeshDict; +} +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +scale 1; + +r 1; +a #eval{1./sqrt(2)*$r}; +b #eval{3*$r}; + +n 16; + +vertices +( + (-$b -$b 0) // 0 + ( $b -$b 0) // 1 + (-$a -$a 0) // 2 + ( $a -$a 0) // 3 + (-$a $a 0) // 4 + ( $a $a 0) // 5 + (-$b $b 0) // 6 + ( $b $b 0) // 7 + + (-$b -$b 1) // 8 + ( $b -$b 1) // 9 + (-$a -$a 1) // 10 + ( $a -$a 1) // 11 + (-$a $a 1) // 12 + ( $a $a 1) // 13 + (-$b $b 1) // 14 + ( $b $b 1) // 15 +); + +blocks +( + hex (0 1 3 2 8 9 11 10) ($n $n 1) simpleGrading (1 1 1) + hex (1 7 5 3 9 15 13 11) ($n $n 1) simpleGrading (1 1 1) + hex (7 6 4 5 15 14 12 13) ($n $n 1) simpleGrading (1 1 1) + hex (6 0 2 4 14 8 10 12) ($n $n 1) simpleGrading (1 1 1) +); + +edges +( + arc 2 3 origin (0 0 0) + arc 3 5 origin (0 0 0) + arc 5 4 origin (0 0 0) + arc 4 2 origin (0 0 0) + + arc 10 11 origin (0 0 1) + arc 11 13 origin (0 0 1) + arc 13 12 origin (0 0 1) + arc 12 10 origin (0 0 1) +); + +boundary +( + walls + { + type wall; + faces + ( + (0 2) + (2 2) + ); + } + cylinder + { + type wall; + faces + ( + (0 3) + (1 3) + (2 3) + (3 3) + ); + } + outlet + { + type patch; + faces + ( + (1 2) + ); + } + inlet + { + type patch; + faces + ( + (3 2) + ); + } +); + +defaultPatch +{ + type empty; + name frontAndBack; +} + + +// ************************************************************************* // diff --git a/tutorials/lagrangian/kinematicParcelFoam/spinningDisk/system/controlDict b/tutorials/lagrangian/kinematicParcelFoam/spinningDisk/system/controlDict new file mode 100644 index 0000000000..dbace75dfe --- /dev/null +++ b/tutorials/lagrangian/kinematicParcelFoam/spinningDisk/system/controlDict @@ -0,0 +1,54 @@ +/*--------------------------------*- C++ -*----------------------------------*\ +| ========= | | +| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox | +| \\ / O peration | Version: v2106 | +| \\ / A nd | Website: www.openfoam.com | +| \\/ M anipulation | | +\*---------------------------------------------------------------------------*/ +FoamFile +{ + version 2.0; + format ascii; + class dictionary; + object controlDict; +} +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +application kinematicParcelFoam; + +startFrom startTime; + +startTime 0; + +stopAt endTime; + +endTime 10; + +deltaT 0.005; + +writeControl adjustable; + +writeInterval 0.25; + +purgeWrite 0; + +writeFormat binary; + +writePrecision 6; + +writeCompression off; + +timeFormat general; + +timePrecision 6; + +runTimeModifiable yes; + +adjustTimeStep yes; + +maxCo 0.9; + +maxDeltaT 0.1; + + +// ************************************************************************* // diff --git a/tutorials/lagrangian/kinematicParcelFoam/spinningDisk/system/decomposeParDict b/tutorials/lagrangian/kinematicParcelFoam/spinningDisk/system/decomposeParDict new file mode 100644 index 0000000000..586c037a8a --- /dev/null +++ b/tutorials/lagrangian/kinematicParcelFoam/spinningDisk/system/decomposeParDict @@ -0,0 +1,27 @@ +/*--------------------------------*- C++ -*----------------------------------*\ +| ========= | | +| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox | +| \\ / O peration | Version: v2106 | +| \\ / A nd | Website: www.openfoam.com | +| \\/ M anipulation | | +\*---------------------------------------------------------------------------*/ +FoamFile +{ + version 2.0; + format ascii; + class dictionary; + object decomposeParDict; +} +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +numberOfSubdomains 4; + +method hierarchical; + +coeffs +{ + n ( 2 2 1 ); +} + + +// ************************************************************************* // \ No newline at end of file diff --git a/tutorials/lagrangian/kinematicParcelFoam/spinningDisk/system/fvSchemes b/tutorials/lagrangian/kinematicParcelFoam/spinningDisk/system/fvSchemes new file mode 100644 index 0000000000..6333789c9a --- /dev/null +++ b/tutorials/lagrangian/kinematicParcelFoam/spinningDisk/system/fvSchemes @@ -0,0 +1,54 @@ +/*--------------------------------*- C++ -*----------------------------------*\ +| ========= | | +| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox | +| \\ / O peration | Version: v2106 | +| \\ / A nd | Website: www.openfoam.com | +| \\/ M anipulation | | +\*---------------------------------------------------------------------------*/ +FoamFile +{ + version 2.0; + format ascii; + class dictionary; + object fvSchemes; +} +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +ddtSchemes +{ + default Euler; +} + +gradSchemes +{ + default Gauss linear; + grad(p) Gauss linear; + grad(U) Gauss linear; +} + +divSchemes +{ + default none; + + div(phi,U) Gauss upwind; + + div((nuEff*dev2(T(grad(U))))) Gauss linear; +} + +laplacianSchemes +{ + default Gauss linear corrected; +} + +interpolationSchemes +{ + default linear; +} + +snGradSchemes +{ + default corrected; +} + + +// ************************************************************************* // diff --git a/tutorials/lagrangian/kinematicParcelFoam/spinningDisk/system/fvSolution b/tutorials/lagrangian/kinematicParcelFoam/spinningDisk/system/fvSolution new file mode 100644 index 0000000000..a4b3321817 --- /dev/null +++ b/tutorials/lagrangian/kinematicParcelFoam/spinningDisk/system/fvSolution @@ -0,0 +1,60 @@ +/*--------------------------------*- C++ -*----------------------------------*\ +| ========= | | +| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox | +| \\ / O peration | Version: v2106 | +| \\ / A nd | Website: www.openfoam.com | +| \\/ M anipulation | | +\*---------------------------------------------------------------------------*/ +FoamFile +{ + version 2.0; + format ascii; + class dictionary; + object fvSolution; +} +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +solvers +{ + p + { + solver PBiCGStab; + tolerance 0; + relTol 0.1; + smoother GaussSeidel; + preconditioner DIC; + } + + pFinal + { + $p; + tolerance 1e-06; + relTol 0; + } + + U + { + solver smoothSolver; + smoother GaussSeidel; + tolerance 1e-05; + relTol 0.01; + } + + UFinal + { + $U; + tolerance 1e-06; + relTol 0; + } +} + +PIMPLE +{ + momentumPredictor yes; + nOuterCorrectors 1; + nCorrectors 2; + nNonOrthogonalCorrectors 0; +} + + +// ************************************************************************* //