diff --git a/src/finiteVolume/Make/files b/src/finiteVolume/Make/files index f904fecd73..a3ca09c0af 100644 --- a/src/finiteVolume/Make/files +++ b/src/finiteVolume/Make/files @@ -257,6 +257,7 @@ $(derivedFvPatchFields)/PrghPressure/prghPressureFvPatchScalarFields.C $(derivedFvPatchFields)/prghTotalHydrostaticPressure/prghTotalHydrostaticPressureFvPatchScalarField.C $(derivedFvPatchFields)/fixedValueInletOutlet/fixedValueInletOutletFvPatchFields.C $(derivedFvPatchFields)/transonicEntrainmentPressure/transonicEntrainmentPressureFvPatchScalarField.C +$(derivedFvPatchFields)/prghCyclicPressure/prghCyclicPressureFvPatchScalarField.C fvsPatchFields = fields/fvsPatchFields $(fvsPatchFields)/fvsPatchField/fvsPatchFields.C diff --git a/src/finiteVolume/fields/fvPatchFields/derived/prghCyclicPressure/prghCyclicPressureFvPatchScalarField.C b/src/finiteVolume/fields/fvPatchFields/derived/prghCyclicPressure/prghCyclicPressureFvPatchScalarField.C new file mode 100644 index 0000000000..a00078c5f5 --- /dev/null +++ b/src/finiteVolume/fields/fvPatchFields/derived/prghCyclicPressure/prghCyclicPressureFvPatchScalarField.C @@ -0,0 +1,207 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | Website: https://openfoam.org + \\ / A nd | Copyright (C) 2024 OpenFOAM Foundation + \\/ M anipulation | +------------------------------------------------------------------------------- +License + This file is part of OpenFOAM. + + OpenFOAM is free software: you can redistribute it and/or modify it + under the terms of the GNU General Public License as published by + the Free Software Foundation, either version 3 of the License, or + (at your option) any later version. + + OpenFOAM is distributed in the hope that it will be useful, but WITHOUT + ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + for more details. + + You should have received a copy of the GNU General Public License + along with OpenFOAM. If not, see . + +\*---------------------------------------------------------------------------*/ + +#include "prghCyclicPressureFvPatchScalarField.H" +#include "addToRunTimeSelectionTable.H" +#include "fieldMapper.H" +#include "volFields.H" +#include "surfaceFields.H" +#include "fvcSnGrad.H" + +// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // + +Foam::prghCyclicPressureFvPatchScalarField::prghCyclicPressureFvPatchScalarField +( + const fvPatch& p, + const DimensionedField& iF, + const dictionary& dict +) +: + jumpCyclicFvPatchScalarField(p, iF, dict), + rhoName_ + ( + cyclicPatch().owner() + ? dict.lookupOrDefault("rho", "rho") + : word::null + ), + rhoInf_(cyclicPatch().owner() ? dict.lookup("rhoInf") : NaN), + jump_(p.size(), Zero) +{ + if (dict.found("jump")) + { + jump_ = scalarField("jump", dict, p.size()); + } + + evaluateNoUpdateCoeffs(); +} + + +Foam::prghCyclicPressureFvPatchScalarField::prghCyclicPressureFvPatchScalarField +( + const prghCyclicPressureFvPatchScalarField& ptf, + const fvPatch& p, + const DimensionedField& iF, + const fieldMapper& mapper +) +: + jumpCyclicFvPatchScalarField(ptf, p, iF, mapper), + rhoName_(ptf.rhoName_), + rhoInf_(ptf.rhoInf_), + jump_(mapper(ptf.jump_)) +{} + + +Foam::prghCyclicPressureFvPatchScalarField::prghCyclicPressureFvPatchScalarField +( + const prghCyclicPressureFvPatchScalarField& ptf, + const DimensionedField& iF +) +: + jumpCyclicFvPatchScalarField(ptf, iF), + rhoName_(ptf.rhoName_), + rhoInf_(ptf.rhoInf_), + jump_(ptf.jump_) +{} + + +// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // + +Foam::tmp +Foam::prghCyclicPressureFvPatchScalarField::jump() const +{ + return jump_; +} + + +void Foam::prghCyclicPressureFvPatchScalarField::map +( + const fvPatchScalarField& ptf, + const fieldMapper& mapper +) +{ + jumpCyclicFvPatchScalarField::map(ptf, mapper); + + const prghCyclicPressureFvPatchScalarField& tiptf = + refCast(ptf); + + mapper(jump_, tiptf.jump_); +} + + +void Foam::prghCyclicPressureFvPatchScalarField::reset +( + const fvPatchScalarField& ptf +) +{ + jumpCyclicFvPatchScalarField::reset(ptf); + + const prghCyclicPressureFvPatchScalarField& tiptf = + refCast(ptf); + + jump_.reset(tiptf.jump_); +} + + +void Foam::prghCyclicPressureFvPatchScalarField::updateCoeffs() +{ + if (updated()) + { + return; + } + + const volScalarField& vf = + static_cast(internalField()); + + const prghCyclicPressureFvPatchScalarField& nbrPf = + refCast(nbrPatchField()); + + const label patchi = patch().index(), nbrPatchi = nbrPf.patch().index(); + + // Buoyancy fields + const volScalarField& rhoVf = + db().lookupObject + (cyclicPatch().owner() ? rhoName_ : nbrPf.rhoName_); + const volScalarField::Boundary& rhoBf = rhoVf.boundaryField(); + const surfaceScalarField::Boundary& ghfBf = + db().lookupObject("ghf").boundaryField(); + + // Pressure solution fields + const surfaceScalarField::Boundary& rAUfBf = + db().lookupObject("rAUf").boundaryField(); + const surfaceScalarField::Boundary& phiHbyABf = + db().lookupObject("phiHbyA").boundaryField(); + + // Delta coefficients for this field + const tmp deltaCoeffsSf = + fv::snGradScheme::New + ( + vf.mesh(), + vf.mesh().schemes().snGrad(internalField().name()) + )->deltaCoeffs(vf); + const scalarField& deltaCoeffsPf = + deltaCoeffsSf->boundaryField()[patch().index()]; + + // Calculate the jump + jump_ = + (rhoBf[patchi] - (cyclicPatch().owner() ? rhoInf_ : nbrPf.rhoInf_)) + *(ghfBf[nbrPatchi] - ghfBf[patchi]) + + ( + phiHbyABf[patchi]/rAUfBf[patchi]/patch().magSf() + + phiHbyABf[nbrPatchi]/rAUfBf[nbrPatchi]/nbrPf.patch().magSf() + )*patch().weights()/deltaCoeffsPf; + + jumpCyclicFvPatchScalarField::updateCoeffs(); +} + + +void Foam::prghCyclicPressureFvPatchScalarField::write +( + Ostream& os +) const +{ + fvPatchScalarField::write(os); + + if (cyclicPatch().owner()) + { + writeEntryIfDifferent(os, "rho", "rho", rhoName_); + writeEntry(os, "rhoInf", rhoInf_); + } + + writeEntry(os, "jump", jump_); +} + + +// * * * * * * * * * * * * * * Build Macro Function * * * * * * * * * * * * // + +namespace Foam +{ + makePatchTypeField + ( + fvPatchScalarField, + prghCyclicPressureFvPatchScalarField + ); +} + +// ************************************************************************* // diff --git a/src/finiteVolume/fields/fvPatchFields/derived/prghCyclicPressure/prghCyclicPressureFvPatchScalarField.H b/src/finiteVolume/fields/fvPatchFields/derived/prghCyclicPressure/prghCyclicPressureFvPatchScalarField.H new file mode 100644 index 0000000000..9c603e2e1c --- /dev/null +++ b/src/finiteVolume/fields/fvPatchFields/derived/prghCyclicPressure/prghCyclicPressureFvPatchScalarField.H @@ -0,0 +1,217 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | Website: https://openfoam.org + \\ / A nd | Copyright (C) 2024 OpenFOAM Foundation + \\/ M anipulation | +------------------------------------------------------------------------------- +License + This file is part of OpenFOAM. + + OpenFOAM is free software: you can redistribute it and/or modify it + under the terms of the GNU General Public License as published by + the Free Software Foundation, either version 3 of the License, or + (at your option) any later version. + + OpenFOAM is distributed in the hope that it will be useful, but WITHOUT + ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + for more details. + + You should have received a copy of the GNU General Public License + along with OpenFOAM. If not, see . + +Class + Foam::prghCyclicPressureFvPatchScalarField + +Description + This boundary condition provides a cyclic condition for p_rgh. It applies + corrections to the value and gradient on both sides of the cyclic to + account for the non-cylicity of the gravitational force. + + This condition is only needed when the cyclic patches have a transformation + and a normal component in the direction of gravity. If the cyclic patches + are orthogonal to the direction gravity, then a normal cyclic boundary + condition can be used instead. + + Care must be taken when using this boundary condition that the simulation + is actually cyclic. The following constraints apply: + + - Both cyclic patches must be oriented in the same way with respect to + gravity. In practice this means that applicability is limited to cyclics + with translational transformations. + + - The model cannot have any dependence on the absolute value of the + pressure field. The absolute value of the pressure, in reality, varies + between each repetition of the geometry; it is not actually formally + cyclic. Only the gradient of the pressure field can be truly cyclic. This + model is therefore only valid if the absolute value of the pressure is + arbitrary, and only the gradient has an effect on the solution. This is + the case for incompressible multiphase solutions or incompressible + Boussinesq-like models of density variation. It is not true if (for + example) a compressible thermodynamic model is being used. + +Usage + \table + Property | Description | Required | Default value + patchType | underlying patch type (should be \c cyclic) | yes | + rhoInf | far-field density | yes | + \endtable + + Example of the boundary condition specification: + \verbatim + + { + type prghCyclicPressure; + patchType cyclic; + rhoInf 1; + } + \endverbatim + +SourceFiles + prghCyclicPressureFvPatchScalarField.C + +\*---------------------------------------------------------------------------*/ + +#ifndef prghCyclicPressureFvPatchScalarField_H +#define prghCyclicPressureFvPatchScalarField_H + +#include "jumpCyclicFvPatchFields.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ + +/*---------------------------------------------------------------------------*\ + Class prghCyclicPressureFvPatchScalarField Declaration +\*---------------------------------------------------------------------------*/ + +class prghCyclicPressureFvPatchScalarField +: + public jumpCyclicFvPatchScalarField +{ + // Private Data + + //- Name of the density field + const word rhoName_; + + //- Far-field density + const scalar rhoInf_; + + //- Jump in value from the other patch to this one + scalarField jump_; + + + // Private Member Functions + + //- Get the far-field density from the owner patch + scalar rhoInf() const; + + //- Update result field based on interface functionality + template + void updateInterfaceMatrix + ( + scalarField& result, + const scalarField& psiInternal, + const scalarField& coeffs, + const Pstream::commsTypes commsType, + const Cmpt ... cmpt + ) const; + + +public: + + //- Runtime type information + TypeName("prghCyclicPressure"); + + + // Constructors + + //- Construct from patch, internal field and dictionary + prghCyclicPressureFvPatchScalarField + ( + const fvPatch&, + const DimensionedField&, + const dictionary& + ); + + //- Construct by mapping given fixedValueTypeFvPatchField + // onto a new patch + prghCyclicPressureFvPatchScalarField + ( + const prghCyclicPressureFvPatchScalarField&, + const fvPatch&, + const DimensionedField&, + const fieldMapper& + ); + + //- Disallow copy without setting internal field reference + prghCyclicPressureFvPatchScalarField + ( + const prghCyclicPressureFvPatchScalarField& + ) = delete; + + //- Copy constructor setting internal field reference + prghCyclicPressureFvPatchScalarField + ( + const prghCyclicPressureFvPatchScalarField&, + const DimensionedField& + ); + + //- Construct and return a clone setting internal field reference + virtual tmp clone + ( + const DimensionedField& iF + ) const + { + return tmp + ( + new prghCyclicPressureFvPatchScalarField + ( + *this, + iF + ) + ); + } + + + // Member Functions + + // Access + + //- Return the "jump" + virtual tmp jump() const; + + + // Mapping functions + + //- Map the given fvPatchField onto this fvPatchField + virtual void map(const fvPatchScalarField&, const fieldMapper&); + + //- Reset the fvPatchField to the given fvPatchField + // Used for mesh to mesh mapping + virtual void reset(const fvPatchScalarField&); + + + // Evaluation functions + + //- Update the patch pressure gradient field + virtual void updateCoeffs(); + + + //- Write + virtual void write(Ostream&) const; +}; + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace Foam + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + + +#endif + +// ************************************************************************* // diff --git a/tutorials/incompressibleVoF/trayedPipe/0/U b/tutorials/incompressibleVoF/trayedPipe/0/U new file mode 100644 index 0000000000..60183d17c9 --- /dev/null +++ b/tutorials/incompressibleVoF/trayedPipe/0/U @@ -0,0 +1,32 @@ +/*--------------------------------*- C++ -*----------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | Website: https://openfoam.org + \\ / A nd | Version: dev + \\/ M anipulation | +\*---------------------------------------------------------------------------*/ +FoamFile +{ + format ascii; + class volVectorField; + location "0"; + object U; +} +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +dimensions [0 1 -1 0 0 0 0]; + +internalField uniform (0 0 0); + +boundaryField +{ + #includeEtc "caseDicts/setConstraintTypes" + + wall + { + type noSlip; + } +} + + +// ************************************************************************* // diff --git a/tutorials/incompressibleVoF/trayedPipe/0/alpha.water.orig b/tutorials/incompressibleVoF/trayedPipe/0/alpha.water.orig new file mode 100644 index 0000000000..a450f9a9a4 --- /dev/null +++ b/tutorials/incompressibleVoF/trayedPipe/0/alpha.water.orig @@ -0,0 +1,32 @@ +/*--------------------------------*- C++ -*----------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | Website: https://openfoam.org + \\ / A nd | Version: dev + \\/ M anipulation | +\*---------------------------------------------------------------------------*/ +FoamFile +{ + format ascii; + class volScalarField; + location "0"; + object alpha.water; +} +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +dimensions []; + +internalField uniform 0; + +boundaryField +{ + #includeEtc "caseDicts/setConstraintTypes" + + wall + { + type zeroGradient; + } +} + + +// ************************************************************************* // diff --git a/tutorials/incompressibleVoF/trayedPipe/0/p_rgh b/tutorials/incompressibleVoF/trayedPipe/0/p_rgh new file mode 100644 index 0000000000..9a3313c123 --- /dev/null +++ b/tutorials/incompressibleVoF/trayedPipe/0/p_rgh @@ -0,0 +1,44 @@ +/*--------------------------------*- C++ -*----------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | Website: https://openfoam.org + \\ / A nd | Version: dev + \\/ M anipulation | +\*---------------------------------------------------------------------------*/ +FoamFile +{ + format ascii; + class volScalarField; + location "0"; + object p_rgh; +} +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +dimensions [1 -1 -2 0 0 0 0]; + +internalField uniform 0; + +boundaryField +{ + #includeEtc "caseDicts/setConstraintTypes" + + wall + { + type fixedFluxPressure; + value $internalField; + } + bottom + { + type prghCyclicPressure; + patchType cyclic; + rhoInf 1; + } + top + { + type prghCyclicPressure; + patchType cyclic; + } +} + + +// ************************************************************************* // diff --git a/tutorials/incompressibleVoF/trayedPipe/Allclean b/tutorials/incompressibleVoF/trayedPipe/Allclean new file mode 100755 index 0000000000..170a3cfe1e --- /dev/null +++ b/tutorials/incompressibleVoF/trayedPipe/Allclean @@ -0,0 +1,9 @@ +#!/bin/sh +cd ${0%/*} || exit 1 # Run from this directory + +# Source tutorial clean functions +. $WM_PROJECT_DIR/bin/tools/CleanFunctions + +cleanVoFCase + +#------------------------------------------------------------------------------ diff --git a/tutorials/incompressibleVoF/trayedPipe/Allrun b/tutorials/incompressibleVoF/trayedPipe/Allrun new file mode 100755 index 0000000000..d33b9817f0 --- /dev/null +++ b/tutorials/incompressibleVoF/trayedPipe/Allrun @@ -0,0 +1,13 @@ +#!/bin/sh +cd ${0%/*} || exit 1 # Run from this directory + +# Source tutorial run functions +. $WM_PROJECT_DIR/bin/tools/RunFunctions + +runApplication blockMesh +runApplication topoSet +runApplication createBaffles -overwrite +runApplication setFields +runApplication $(getApplication) + +#------------------------------------------------------------------------------ diff --git a/tutorials/incompressibleVoF/trayedPipe/constant/g b/tutorials/incompressibleVoF/trayedPipe/constant/g new file mode 100644 index 0000000000..770a56192e --- /dev/null +++ b/tutorials/incompressibleVoF/trayedPipe/constant/g @@ -0,0 +1,21 @@ +/*--------------------------------*- C++ -*----------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | Website: https://openfoam.org + \\ / A nd | Version: dev + \\/ M anipulation | +\*---------------------------------------------------------------------------*/ +FoamFile +{ + format ascii; + class uniformDimensionedVectorField; + location "constant"; + object g; +} +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +dimensions [0 1 -2 0 0 0 0]; +value (0 -9.81 0); + + +// ************************************************************************* // diff --git a/tutorials/incompressibleVoF/trayedPipe/constant/hRef b/tutorials/incompressibleVoF/trayedPipe/constant/hRef new file mode 100644 index 0000000000..2aeddd72d4 --- /dev/null +++ b/tutorials/incompressibleVoF/trayedPipe/constant/hRef @@ -0,0 +1,23 @@ +/*--------------------------------*- C++ -*----------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | Website: https://openfoam.org + \\ / A nd | Version: dev + \\/ M anipulation | +\*---------------------------------------------------------------------------*/ +FoamFile +{ + format ascii; + class uniformDimensionedScalarField; + location "constant"; + object hRef; +} +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#include "$FOAM_CASE/system/blockMeshDict" + +dimensions [0 1 0 0 0 0 0]; +value #calc "${height}/2"; + + +// ************************************************************************* // diff --git a/tutorials/incompressibleVoF/trayedPipe/constant/momentumTransport b/tutorials/incompressibleVoF/trayedPipe/constant/momentumTransport new file mode 100644 index 0000000000..8278c989ec --- /dev/null +++ b/tutorials/incompressibleVoF/trayedPipe/constant/momentumTransport @@ -0,0 +1,20 @@ +/*--------------------------------*- C++ -*----------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | Website: https://openfoam.org + \\ / A nd | Version: dev + \\/ M anipulation | +\*---------------------------------------------------------------------------*/ +FoamFile +{ + format ascii; + class dictionary; + location "constant"; + object momentumTransport; +} +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +simulationType laminar; + + +// ************************************************************************* // diff --git a/tutorials/incompressibleVoF/trayedPipe/constant/phaseProperties b/tutorials/incompressibleVoF/trayedPipe/constant/phaseProperties new file mode 100644 index 0000000000..234504e0d1 --- /dev/null +++ b/tutorials/incompressibleVoF/trayedPipe/constant/phaseProperties @@ -0,0 +1,22 @@ +/*--------------------------------*- C++ -*----------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | Website: https://openfoam.org + \\ / A nd | Version: dev + \\/ M anipulation | +\*---------------------------------------------------------------------------*/ +FoamFile +{ + format ascii; + class dictionary; + location "constant"; + object phaseProperties; +} +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +phases (water air); + +sigma 0.07; + + +// ************************************************************************* // diff --git a/tutorials/incompressibleVoF/trayedPipe/constant/physicalProperties.air b/tutorials/incompressibleVoF/trayedPipe/constant/physicalProperties.air new file mode 100644 index 0000000000..fab7bde58b --- /dev/null +++ b/tutorials/incompressibleVoF/trayedPipe/constant/physicalProperties.air @@ -0,0 +1,24 @@ +/*--------------------------------*- C++ -*----------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | Website: https://openfoam.org + \\ / A nd | Version: dev + \\/ M anipulation | +\*---------------------------------------------------------------------------*/ +FoamFile +{ + format ascii; + class dictionary; + location "constant"; + object physicalProperties.air; +} +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +viscosityModel constant; + +nu 1.48e-05; + +rho 1; + + +// ************************************************************************* // diff --git a/tutorials/incompressibleVoF/trayedPipe/constant/physicalProperties.water b/tutorials/incompressibleVoF/trayedPipe/constant/physicalProperties.water new file mode 100644 index 0000000000..e9b9d7dc4e --- /dev/null +++ b/tutorials/incompressibleVoF/trayedPipe/constant/physicalProperties.water @@ -0,0 +1,24 @@ +/*--------------------------------*- C++ -*----------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | Website: https://openfoam.org + \\ / A nd | Version: dev + \\/ M anipulation | +\*---------------------------------------------------------------------------*/ +FoamFile +{ + format ascii; + class dictionary; + location "constant"; + object physicalProperties.water; +} +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +viscosityModel constant; + +nu 1e-06; + +rho 1000; + + +// ************************************************************************* // diff --git a/tutorials/incompressibleVoF/trayedPipe/system/blockMeshDict b/tutorials/incompressibleVoF/trayedPipe/system/blockMeshDict new file mode 100644 index 0000000000..6f94b56c6d --- /dev/null +++ b/tutorials/incompressibleVoF/trayedPipe/system/blockMeshDict @@ -0,0 +1,89 @@ +/*--------------------------------*- C++ -*----------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | Website: https://openfoam.org + \\ / A nd | Version: dev + \\/ M anipulation | +\*---------------------------------------------------------------------------*/ +FoamFile +{ + format ascii; + class dictionary; + location "system"; + object blockMeshDict; +} +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +height 0.045; +radius 0.015; +nCellsPerMetre 2666.67; + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +halfThickness #calc "scalar(1)/$nCellsPerMetre"; + +nCellsHeight #calc "round($height*$nCellsPerMetre)"; +nCellsRadius #calc "round($radius*$nCellsPerMetre)"; + +convertToMeters 1; + +vertices +( + (0 0 0) (0 $height 0) + ($radius 0 #neg $halfThickness) ($radius $height #neg $halfThickness) + ($radius 0 $halfThickness) ($radius $height $halfThickness) +); + +blocks +( + hex (1 0 2 3 1 0 4 5) ($nCellsHeight $nCellsRadius 1) simpleGrading (1 1 1) +); + +boundary +( + walls + { + type wall; + faces + ( + (2 3 5 4) + ); + } + bottom + { + type cyclic; + faces + ( + (0 0 2 4) + ); + neighbourPatch top; + } + top + { + type cyclic; + faces + ( + (1 1 3 5) + ); + neighbourPatch bottom; + } + wedgeBack + { + type wedge; + faces + ( + (0 1 3 2) + ); + } + wedgeFront + { + type wedge; + faces + ( + (0 1 5 4) + ); + } +); + + +// ************************************************************************* // diff --git a/tutorials/incompressibleVoF/trayedPipe/system/controlDict b/tutorials/incompressibleVoF/trayedPipe/system/controlDict new file mode 100644 index 0000000000..7efdce8cab --- /dev/null +++ b/tutorials/incompressibleVoF/trayedPipe/system/controlDict @@ -0,0 +1,57 @@ +/*--------------------------------*- C++ -*----------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | Website: https://openfoam.org + \\ / A nd | Version: dev + \\/ M anipulation | +\*---------------------------------------------------------------------------*/ +FoamFile +{ + format ascii; + class dictionary; + location "system"; + object controlDict; +} +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +application foamRun; + +solver incompressibleVoF; + +startFrom startTime; + +startTime 0; + +stopAt endTime; + +endTime 2; + +deltaT 0.0001; + +writeControl adjustableRunTime; + +writeInterval 0.005; + +purgeWrite 0; + +writeFormat binary; + +writePrecision 6; + +writeCompression off; + +timeFormat general; + +timePrecision 6; + +runTimeModifiable yes; + +adjustTimeStep yes; + +maxCo 1; +maxAlphaCo 1; + +maxDeltaT 1; + + +// ************************************************************************* // diff --git a/tutorials/incompressibleVoF/trayedPipe/system/createBafflesDict b/tutorials/incompressibleVoF/trayedPipe/system/createBafflesDict new file mode 100644 index 0000000000..328139ee5f --- /dev/null +++ b/tutorials/incompressibleVoF/trayedPipe/system/createBafflesDict @@ -0,0 +1,39 @@ +/*--------------------------------*- C++ -*----------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | Website: https://openfoam.org + \\ / A nd | Version: dev + \\/ M anipulation | +\*---------------------------------------------------------------------------*/ +FoamFile +{ + format ascii; + class dictionary; + location "system"; + object createBafflesDict; +} +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +internalFacesOnly true; + +baffles +{ + wall + { + type faceZone; + zoneName baffles; + + owner + { + name baffles; + type wall; + } + neighbour + { + $owner; + } + } +} + + +// ************************************************************************* // diff --git a/tutorials/incompressibleVoF/trayedPipe/system/fvSchemes b/tutorials/incompressibleVoF/trayedPipe/system/fvSchemes new file mode 100644 index 0000000000..bc0a22f201 --- /dev/null +++ b/tutorials/incompressibleVoF/trayedPipe/system/fvSchemes @@ -0,0 +1,50 @@ +/*--------------------------------*- C++ -*----------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | Website: https://openfoam.org + \\ / A nd | Version: dev + \\/ M anipulation | +\*---------------------------------------------------------------------------*/ +FoamFile +{ + format ascii; + class dictionary; + location "system"; + object fvSchemes; +} +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +ddtSchemes +{ + default Euler; +} + +gradSchemes +{ + default Gauss linear; +} + +divSchemes +{ + div(phi,alpha) Gauss interfaceCompression vanLeer 1; + div(rhoPhi,U) Gauss linearUpwind grad(U); + div(((rho*nuEff)*dev2(T(grad(U))))) Gauss linear; +} + +laplacianSchemes +{ + default Gauss linear uncorrected; +} + +interpolationSchemes +{ + default linear; +} + +snGradSchemes +{ + default uncorrected; +} + + +// ************************************************************************* // diff --git a/tutorials/incompressibleVoF/trayedPipe/system/fvSolution b/tutorials/incompressibleVoF/trayedPipe/system/fvSolution new file mode 100644 index 0000000000..3c20ad7e52 --- /dev/null +++ b/tutorials/incompressibleVoF/trayedPipe/system/fvSolution @@ -0,0 +1,74 @@ +/*--------------------------------*- C++ -*----------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | Website: https://openfoam.org + \\ / A nd | Version: dev + \\/ M anipulation | +\*---------------------------------------------------------------------------*/ +FoamFile +{ + format ascii; + class dictionary; + location "system"; + object fvSolution; +} +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +solvers +{ + "alpha.water.*" + { + nAlphaCorr 2; + nAlphaSubCycles 1; + + MULESCorr yes; + nLimiterIter 3; + + solver smoothSolver; + smoother symGaussSeidel; + tolerance 1e-8; + relTol 0; + } + + "pcorr.*" + { + solver PCG; + preconditioner DIC; + tolerance 1e-5; + relTol 0; + } + + p_rgh + { + $pcorr; + tolerance 1e-7; + relTol 0.05; + } + + p_rghFinal + { + $p_rgh; + relTol 0; + } +} + +PIMPLE +{ + momentumPredictor no; + nOuterCorrectors 1; + nCorrectors 3; + nNonOrthogonalCorrectors 0; + pRefPoint (0 0.015 0); + pRefValue 0; +} + +relaxationFactors +{ + equations + { + ".*" 1; + } +} + + +// ************************************************************************* // diff --git a/tutorials/incompressibleVoF/trayedPipe/system/setFieldsDict b/tutorials/incompressibleVoF/trayedPipe/system/setFieldsDict new file mode 100644 index 0000000000..8883316095 --- /dev/null +++ b/tutorials/incompressibleVoF/trayedPipe/system/setFieldsDict @@ -0,0 +1,43 @@ +/*--------------------------------*- C++ -*----------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | Website: https://openfoam.org + \\ / A nd | Version: dev + \\/ M anipulation | +\*---------------------------------------------------------------------------*/ +FoamFile +{ + format ascii; + class dictionary; + location "system"; + object setFieldsDict; +} +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +alphaMean 0.25; + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#include "blockMeshDict" + +defaultFieldValues +( + volScalarFieldValue alpha.water 0 +); + +regions +( + cylinderToCell + { + point1 (0 0 0); + point2 (0 $height 0); + radius #calc "$radius*sqrt($alphaMean)"; + fieldValues + ( + volScalarFieldValue alpha.water 1 + ); + } +); + + +// ************************************************************************* // diff --git a/tutorials/incompressibleVoF/trayedPipe/system/topoSetDict b/tutorials/incompressibleVoF/trayedPipe/system/topoSetDict new file mode 100644 index 0000000000..bd5832daef --- /dev/null +++ b/tutorials/incompressibleVoF/trayedPipe/system/topoSetDict @@ -0,0 +1,128 @@ +/*--------------------------------*- C++ -*----------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | Website: https://openfoam.org + \\ / A nd | Version: dev + \\/ M anipulation | +\*---------------------------------------------------------------------------*/ +FoamFile +{ + format ascii; + class dictionary; + location "system"; + object topoSetDict; +} +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +nBaffles 6; +baffleOverlap 0.4; + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#include "blockMeshDict" + +innerBaffleOuterRadius #calc "${radius}/(2 - $baffleOverlap)"; +outerBaffleInnerRadius #calc "$radius - $innerBaffleOuterRadius"; + +actions +( + // Create baffles along the length of the pipe, alternating the zone + // between "innerBaffles" and "outerBaffles" + #codeStream + { + code + #{ + for (label i = 0; i < $