zhyang-dev/openfoam
1
0
Fork 0
mirror of https://develop.openfoam.com/Development/openfoam.git synced 2025-11-28 03:28:01 +00:00
Code Issues Actions Packages Projects Releases Wiki Activity

Feature grad update

Browse Source
This commit is contained in:
Mattijs Janssens
2025-10-22 15:46:53 +00:00
committed by Andrew Heather
parent f186b0098b
commit e100e4d084
16 changed files with 543 additions and 105 deletions
Download Patch File Download Diff File Expand all files Collapse all files

86
src/finiteVolume/finiteVolume/gradSchemes/LeastSquaresGrad/LeastSquaresGrad.C
View File

@ -1,4 +1,4 @@
/*---------------------------------------------------------------------------*\
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
@ -36,23 +36,18 @@ License
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
template<class Type, class Stencil>
Foam::tmp
<
Foam::GeometricField
<
typename Foam::outerProduct<Foam::vector, Type>::type,
Foam::fvPatchField,
Foam::volMesh
>
>
Foam::fv::LeastSquaresGrad<Type, Stencil>::calcGrad
void Foam::fv::LeastSquaresGrad<Type, Stencil>::calcGrad
(
const GeometricField<Type, fvPatchField, volMesh>& vtf,
const word& name
GeometricField
<
typename outerProduct<vector, Type>::type,
fvPatchField,
volMesh
>& lsGrad,
const GeometricField<Type, fvPatchField, volMesh>& vtf
) const
{
typedef typename outerProduct<vector, Type>::type GradType;
typedef GeometricField<GradType, fvPatchField, volMesh> GradFieldType;
const fvMesh& mesh = vtf.mesh();
@ -62,24 +57,7 @@ Foam::fv::LeastSquaresGrad<Type, Stencil>::calcGrad
mesh
);
tmp<GradFieldType> tlsGrad
(
new GradFieldType
(
IOobject
(
name,
vtf.instance(),
mesh,
IOobject::NO_READ,
IOobject::NO_WRITE
),
mesh,
dimensioned<GradType>(vtf.dimensions()/dimLength, Zero),
fvPatchFieldBase::extrapolatedCalculatedType()
)
);
GradFieldType& lsGrad = tlsGrad.ref();
lsGrad = dimensioned<GradType>(vtf.dimensions()/dimLength, Zero);
Field<GradType>& lsGradIf = lsGrad;
const extendedCentredCellToCellStencil& stencil = lsv.stencil();
@ -131,6 +109,50 @@ Foam::fv::LeastSquaresGrad<Type, Stencil>::calcGrad
// Correct the boundary conditions
lsGrad.correctBoundaryConditions();
gaussGrad<Type>::correctBoundaryConditions(vtf, lsGrad);
}
template<class Type, class Stencil>
Foam::tmp
<
Foam::GeometricField
<
typename Foam::outerProduct<Foam::vector, Type>::type,
Foam::fvPatchField,
Foam::volMesh
>
>
Foam::fv::LeastSquaresGrad<Type, Stencil>::calcGrad
(
const GeometricField<Type, fvPatchField, volMesh>& vtf,
const word& name
) const
{
typedef typename outerProduct<vector, Type>::type GradType;
typedef GeometricField<GradType, fvPatchField, volMesh> GradFieldType;
const fvMesh& mesh = vtf.mesh();
tmp<GradFieldType> tlsGrad
(
new GradFieldType
(
IOobject
(
name,
vtf.instance(),
mesh,
IOobject::NO_READ,
IOobject::NO_WRITE
),
mesh,
vtf.dimensions()/dimLength,
fvPatchFieldBase::extrapolatedCalculatedType()
)
);
GradFieldType& lsGrad = tlsGrad.ref();
calcGrad(lsGrad, vtf);
return tlsGrad;
}

8
src/finiteVolume/finiteVolume/gradSchemes/LeastSquaresGrad/LeastSquaresGrad.H
View File

@ -132,6 +132,14 @@ public:
const GeometricField<Type, fvPatchField, volMesh>& vsf,
const word& name
) const;
//- Calculate the grad of the given field into supplied field
virtual void calcGrad
(
GeometricField
<typename outerProduct<vector, Type>::type, fvPatchField, volMesh>& res,
const GeometricField<Type, fvPatchField, volMesh>&
) const;
};

147
src/finiteVolume/finiteVolume/gradSchemes/fourthGrad/fourthGrad.C
View File

@ -38,19 +38,21 @@ License
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
template<class Type>
Foam::tmp
<
Foam::GeometricField
<
typename Foam::outerProduct<Foam::vector, Type>::type,
Foam::fvPatchField,
Foam::volMesh
>
>
Foam::fv::fourthGrad<Type>::calcGrad
void Foam::fv::fourthGrad<Type>::calcGrad
(
GeometricField
<
typename outerProduct<vector, Type>::type,
fvPatchField,
volMesh
>& fGrad,
const GeometricField<Type, fvPatchField, volMesh>& vsf,
const word& name
const GeometricField
<
typename outerProduct<vector, Type>::type,
fvPatchField,
volMesh
>& secondfGrad
) const
{
// The fourth-order gradient is calculated in two passes. First,
@ -59,37 +61,11 @@ Foam::fv::fourthGrad<Type>::calcGrad
// gradient to complete the accuracy.
typedef typename outerProduct<vector, Type>::type GradType;
typedef GeometricField<GradType, fvPatchField, volMesh> GradFieldType;
const fvMesh& mesh = vsf.mesh();
// Assemble the second-order least-square gradient
// Calculate the second-order least-square gradient
tmp<GradFieldType> tsecondfGrad
= leastSquaresGrad<Type>(mesh).grad
(
vsf,
"leastSquaresGrad(" + vsf.name() + ")"
);
const GradFieldType& secondfGrad =
tsecondfGrad();
tmp<GradFieldType> tfGrad
(
new GradFieldType
(
IOobject
(
name,
vsf.instance(),
mesh,
IOobject::NO_READ,
IOobject::NO_WRITE
),
secondfGrad
)
);
GradFieldType& fGrad = tfGrad.ref();
// Initialise to gradient
fGrad = secondfGrad;
const vectorField& C = mesh.C();
@ -158,9 +134,102 @@ Foam::fv::fourthGrad<Type>::calcGrad
fGrad.correctBoundaryConditions();
gaussGrad<Type>::correctBoundaryConditions(vsf, fGrad);
}
template<class Type>
Foam::tmp
<
Foam::GeometricField
<
typename Foam::outerProduct<Foam::vector, Type>::type,
Foam::fvPatchField,
Foam::volMesh
>
>
Foam::fv::fourthGrad<Type>::calcGrad
(
const GeometricField<Type, fvPatchField, volMesh>& vsf,
const word& name
) const
{
// The fourth-order gradient is calculated in two passes. First,
// the standard least-square gradient is assembled. Then, the
// fourth-order correction is added to the second-order accurate
// gradient to complete the accuracy.
typedef typename outerProduct<vector, Type>::type GradType;
typedef GeometricField<GradType, fvPatchField, volMesh> GradFieldType;
const fvMesh& mesh = vsf.mesh();
// Assemble the second-order least-square gradient
// Calculate the second-order least-square gradient
tmp<GradFieldType> tsecondfGrad
= leastSquaresGrad<Type>(mesh).grad
(
vsf,
"leastSquaresGrad(" + vsf.name() + ")"
);
const GradFieldType& secondfGrad = tsecondfGrad();
tmp<GradFieldType> tfGrad
(
new GradFieldType
(
IOobject
(
name,
vsf.instance(),
mesh,
IOobject::NO_READ,
IOobject::NO_WRITE
),
mesh,
vsf.dimensions()/dimLength,
fvPatchFieldBase::extrapolatedCalculatedType()
)
);
calcGrad(tfGrad.ref(), vsf, secondfGrad);
return tfGrad;
}
template<class Type>
void Foam::fv::fourthGrad<Type>::calcGrad
(
GeometricField
<
typename outerProduct<vector, Type>::type,
fvPatchField,
volMesh
>& fGrad,
const GeometricField<Type, fvPatchField, volMesh>& vsf
) const
{
// The fourth-order gradient is calculated in two passes. First,
// the standard least-square gradient is assembled. Then, the
// fourth-order correction is added to the second-order accurate
// gradient to complete the accuracy.
typedef typename outerProduct<vector, Type>::type GradType;
typedef GeometricField<GradType, fvPatchField, volMesh> GradFieldType;
const fvMesh& mesh = vsf.mesh();
// Assemble the second-order least-square gradient
// Calculate the second-order least-square gradient
tmp<GradFieldType> tsecondfGrad
= leastSquaresGrad<Type>(mesh).grad
(
vsf,
"leastSquaresGrad(" + vsf.name() + ")"
);
calcGrad(fGrad, vsf, tsecondfGrad());
}
// ************************************************************************* //

26
src/finiteVolume/finiteVolume/gradSchemes/fourthGrad/fourthGrad.H
View File

@ -99,6 +99,32 @@ public:
const GeometricField<Type, fvPatchField, volMesh>& vsf,
const word& name
) const;
//- Calculate the grad of the given field into supplied field
virtual void calcGrad
(
GeometricField
<typename outerProduct<vector, Type>::type, fvPatchField, volMesh>& res,
const GeometricField<Type, fvPatchField, volMesh>&
) const;
//- Helper : calculate the grad of the given field into supplied field
void calcGrad
(
GeometricField
<
typename outerProduct<vector, Type>::type,
fvPatchField,
volMesh
>& fGrad,
const GeometricField<Type, fvPatchField, volMesh>& vsf,
const GeometricField
<
typename outerProduct<vector, Type>::type,
fvPatchField,
volMesh
>& secondfGrad
) const;
};

61
src/finiteVolume/finiteVolume/gradSchemes/gaussGrad/gaussGrad.C
View File

@ -140,6 +140,67 @@ Foam::fv::gaussGrad<Type>::calcGrad
}
template<class Type>
void Foam::fv::gaussGrad<Type>::calcGrad
(
GeometricField
<
typename outerProduct<vector, Type>::type,
fvPatchField,
volMesh
>& gGrad,
const GeometricField<Type, fvPatchField, volMesh>& vsf
) const
{
typedef typename outerProduct<vector, Type>::type GradType;
DebugPout<< "gaussGrad<Type>::calcGrad on " << vsf.name()
<< " into " << gGrad.name() << endl;
const fvMesh& mesh = vsf.mesh();
const labelUList& owner = mesh.owner();
const labelUList& neighbour = mesh.neighbour();
const vectorField& Sf = mesh.Sf();
const auto tssf(tinterpScheme_().interpolate(vsf));
const auto& ssf = tssf();
gGrad = dimensioned<GradType>(vsf.dimensions()/dimLength, Zero);
Field<GradType>& igGrad = gGrad;
const Field<Type>& issf = ssf;
forAll(owner, facei)
{
const GradType Sfssf = Sf[facei]*issf[facei];
igGrad[owner[facei]] += Sfssf;
igGrad[neighbour[facei]] -= Sfssf;
}
forAll(mesh.boundary(), patchi)
{
const labelUList& pFaceCells =
mesh.boundary()[patchi].faceCells();
const vectorField& pSf = mesh.Sf().boundaryField()[patchi];
const fvsPatchField<Type>& pssf = ssf.boundaryField()[patchi];
forAll(mesh.boundary()[patchi], facei)
{
igGrad[pFaceCells[facei]] += pSf[facei]*pssf[facei];
}
}
igGrad /= mesh.V();
gGrad.correctBoundaryConditions();
correctBoundaryConditions(vsf, gGrad);
}
template<class Type>
void Foam::fv::gaussGrad<Type>::correctBoundaryConditions
(

8
src/finiteVolume/finiteVolume/gradSchemes/gaussGrad/gaussGrad.H
View File

@ -145,6 +145,14 @@ public:
const word& name
) const;
//- Calculate the grad of the given field into supplied field
virtual void calcGrad
(
GeometricField
<typename outerProduct<vector, Type>::type, fvPatchField, volMesh>& res,
const GeometricField<Type, fvPatchField, volMesh>&
) const;
//- Correct the boundary values of the gradient using the patchField
//- snGrad functions
static void correctBoundaryConditions

28
src/finiteVolume/finiteVolume/gradSchemes/gradScheme/gradScheme.C
View File

@ -96,7 +96,7 @@ Foam::fv::gradScheme<Type>::grad
GradFieldType* pgGrad =
mesh().objectRegistry::template getObjectPtr<GradFieldType>(name);
if (!this->mesh().cache(name) || this->mesh().changing())
if (!this->mesh().cache(name) || this->mesh().topoChanging())
{
// Delete any old occurrences to avoid double registration
if (pgGrad && pgGrad->ownedByRegistry())
@ -119,17 +119,14 @@ Foam::fv::gradScheme<Type>::grad
}
else
{
if (pgGrad->upToDate(vsf))
if (pgGrad->upToDate(vsf) && this->mesh().upToDatePoints(*pgGrad))
{
solution::cachePrintMessage("Reusing", name, vsf);
}
else
{
solution::cachePrintMessage("Updating", name, vsf);
delete pgGrad;
pgGrad = calcGrad(vsf, name).ptr();
regIOobject::store(pgGrad);
calcGrad(*pgGrad, vsf);
}
}
@ -180,4 +177,23 @@ Foam::fv::gradScheme<Type>::grad
}
template<class Type>
void Foam::fv::gradScheme<Type>::calcGrad
(
GeometricField
<
typename outerProduct<vector, Type>::type,
fvPatchField,
volMesh
>& result,
const GeometricField<Type, fvPatchField, volMesh>& vsf
) const
{
DebugPout<< "gradScheme<Type>::calcGrad on " << vsf.name()
<< " into " << result.name() << endl;
result = calcGrad(vsf, result.name());
}
// ************************************************************************* //

8
src/finiteVolume/finiteVolume/gradSchemes/gradScheme/gradScheme.H
View File

@ -142,6 +142,14 @@ public:
const word& name
) const = 0;
//- Calculate the grad of the given field into supplied field
virtual void calcGrad
(
GeometricField
<typename outerProduct<vector, Type>::type, fvPatchField, volMesh>& res,
const GeometricField<Type, fvPatchField, volMesh>&
) const;
//- Calculate and return the grad of the given field
//- which may have been cached
tmp

84
src/finiteVolume/finiteVolume/gradSchemes/limitedGradSchemes/cellLimitedGrad/cellLimitedGrad.C
View File

@ -62,40 +62,25 @@ void Foam::fv::cellLimitedGrad<Type, Limiter>::limitGradient
template<class Type, class Limiter>
Foam::tmp
<
Foam::GeometricField
<
typename Foam::outerProduct<Foam::vector, Type>::type,
Foam::fvPatchField,
Foam::volMesh
>
>
Foam::fv::cellLimitedGrad<Type, Limiter>::calcGrad
void Foam::fv::cellLimitedGrad<Type, Limiter>::calcGrad
(
const GeometricField<Type, fvPatchField, volMesh>& vsf,
const word& name
) const
{
const fvMesh& mesh = vsf.mesh();
tmp
<
GeometricField
<typename outerProduct<vector, Type>::type, fvPatchField, volMesh>
> tGrad = basicGradScheme_().calcGrad(vsf, name);
if (k_ < SMALL)
{
return tGrad;
}
GeometricField
<
typename outerProduct<vector, Type>::type,
fvPatchField,
volMesh
>& g = tGrad.ref();
>& g,
const GeometricField<Type, fvPatchField, volMesh>& vsf
) const
{
const fvMesh& mesh = vsf.mesh();
basicGradScheme_().calcGrad(g, vsf);
if (k_ < SMALL)
{
return;
}
const labelUList& owner = mesh.owner();
const labelUList& neighbour = mesh.neighbour();
@ -227,6 +212,49 @@ Foam::fv::cellLimitedGrad<Type, Limiter>::calcGrad
limitGradient(limiter, g);
g.correctBoundaryConditions();
gaussGrad<Type>::correctBoundaryConditions(vsf, g);
}
template<class Type, class Limiter>
Foam::tmp
<
Foam::GeometricField
<
typename Foam::outerProduct<Foam::vector, Type>::type,
Foam::fvPatchField,
Foam::volMesh
>
>
Foam::fv::cellLimitedGrad<Type, Limiter>::calcGrad
(
const GeometricField<Type, fvPatchField, volMesh>& vsf,
const word& name
) const
{
typedef typename outerProduct<vector, Type>::type GradType;
typedef GeometricField<GradType, fvPatchField, volMesh> GradFieldType;
const fvMesh& mesh = vsf.mesh();
tmp<GradFieldType> tGrad
(
new GradFieldType
(
IOobject
(
name,
vsf.instance(),
mesh,
IOobject::NO_READ,
IOobject::NO_WRITE
),
mesh,
vsf.dimensions()/dimLength,
fvPatchFieldBase::extrapolatedCalculatedType()
)
);
calcGrad(tGrad.ref(), vsf);
return tGrad;
}

8
src/finiteVolume/finiteVolume/gradSchemes/limitedGradSchemes/cellLimitedGrad/cellLimitedGrad.H
View File

@ -150,6 +150,14 @@ public:
const GeometricField<Type, fvPatchField, volMesh>& vsf,
const word& name
) const;
//- Calculate the grad of the given field into supplied field
virtual void calcGrad
(
GeometricField
<typename outerProduct<vector, Type>::type, fvPatchField, volMesh>& res,
const GeometricField<Type, fvPatchField, volMesh>&
) const;
};

88
src/fused/finiteVolume/fusedGaussGrad.C
View File

@ -226,6 +226,94 @@ Foam::fv::fusedGaussGrad<Type>::calcGrad
}
template<class Type>
void Foam::fv::fusedGaussGrad<Type>::calcGrad
(
GeometricField
<
typename outerProduct<vector, Type>::type,
fvPatchField,
volMesh
>& gGrad,
const GeometricField<Type, fvPatchField, volMesh>& vf
) const
{
typedef typename outerProduct<vector, Type>::type GradType;
const fvMesh& mesh = vf.mesh();
DebugPout<< "fusedGaussGrad<Type>::calcGrad on " << vf.name()
<< " into " << gGrad.name() << endl;
gGrad = dimensioned<GradType>(vf.dimensions()/dimLength, Zero);
if (this->tinterpScheme_().corrected())
{
const auto tfaceCorr(this->tinterpScheme_().correction(vf));
auto& faceCorr = tfaceCorr();
DebugPout<< "fusedGaussGrad<Type>::calcGrad corrected interpScheme "
<< this->tinterpScheme_().type() << endl;
const auto interpolate = [&]
(
const vector& area,
const scalar lambda,
const Type& ownVal,
const Type& neiVal,
const Type& correction
) -> GradType
{
return area*((lambda*(ownVal - neiVal) + neiVal) + correction);
};
fvc::surfaceSum
(
this->tinterpScheme_().weights(vf),
vf,
faceCorr,
interpolate,
gGrad,
false
);
}
else
{
DebugPout<< "fusedGaussGrad<Type>::calcGrad uncorrected interpScheme "
<< this->tinterpScheme_().type() << endl;
const auto interpolate = [&]
(
const vector& area,
const scalar lambda,
const Type& ownVal,
const Type& neiVal
) -> GradType
{
return area*(lambda*(ownVal - neiVal) + neiVal);
};
fvc::surfaceSum
(
tinterpScheme_().weights(vf),
vf,
interpolate,
gGrad,
false
);
}
gGrad.primitiveFieldRef() /= mesh.V();
gGrad.correctBoundaryConditions();
correctBoundaryConditions(vf, gGrad);
}
template<class Type>
template<class GradType>
void Foam::fv::fusedGaussGrad<Type>::correctBoundaryConditions

8
src/fused/finiteVolume/fusedGaussGrad.H
View File

@ -145,6 +145,14 @@ public:
const word& name
) const;
//- Calculate the grad of the given field into supplied field
virtual void calcGrad
(
GeometricField
<typename outerProduct<vector, Type>::type, fvPatchField, volMesh>& res,
const GeometricField<Type, fvPatchField, volMesh>&
) const;
//- Correct the boundary values of the gradient using the patchField
//- snGrad functions
template<class GradType>

8
tutorials/incompressible/simpleFoam/pitzDaily_fused/Allrun
View File

@ -19,4 +19,12 @@ sed "s/GAUSS/fusedGauss/g" system/fvSchemes.template > system/fvSchemes
runParallel -s fusedGauss $(getApplication)
echo "Updating fvSolution to cache gradients"
sed 's/GRADIENT/"grad(.*)"/g' system/fvSolution.template > system/fvSolution
runParallel -s fusedGaussCached $(getApplication)
# Restore
cp system/fvSolution.template system/fvSolution
#------------------------------------------------------------------------------

9
tutorials/incompressible/simpleFoam/pitzDaily_fused/system/controlDict
View File

@ -16,6 +16,15 @@ FoamFile
libs (fusedFiniteVolume);
/*
DebugSwitches
{
fusedGauss 1;
Gauss 1;
solution 1;
}
*/
application simpleFoam;
startFrom startTime;

6
tutorials/incompressible/simpleFoam/pitzDaily_fused/system/fvSolution
View File

@ -55,5 +55,11 @@ relaxationFactors
}
}
cache
{
// Cache all grads
GRADIENT;
}
// ************************************************************************* //

65
tutorials/incompressible/simpleFoam/pitzDaily_fused/system/fvSolution.template Normal file
View File

@ -0,0 +1,65 @@
/*--------------------------------*- C++ -*----------------------------------*\
| ========= | |
| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \\ / O peration | Version: v2506 |
| \\ / A nd | Website: www.openfoam.com |
| \\/ M anipulation | |
\*---------------------------------------------------------------------------*/
FoamFile
{
version 2.0;
format ascii;
class dictionary;
object fvSolution;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
solvers
{
p
{
solver GAMG;
tolerance 1e-06;
relTol 0.1;
smoother GaussSeidel;
}
"(U|k|epsilon|omega|f|v2)"
{
solver smoothSolver;
smoother symGaussSeidel;
tolerance 1e-05;
relTol 0.1;
}
}
SIMPLE
{
nNonOrthogonalCorrectors 0;
consistent yes;
residualControl
{
p 1e-2;
U 1e-3;
"(k|epsilon|omega|f|v2)" 1e-3;
}
}
relaxationFactors
{
equations
{
U 0.9; // 0.9 is more stable but 0.95 more convergent
".*" 0.9; // 0.9 is more stable but 0.95 more convergent
}
}
cache
{
// Cache all grads
GRADIENT;
}
// ************************************************************************* //