zhyang-dev/OpenFOAM-6
1
0
Fork 0
mirror of https://github.com/OpenFOAM/OpenFOAM-6.git synced 2025-12-08 06:57:46 +00:00
Code Issues Actions Packages Projects Releases Wiki Activity

integrationSchemes: Corrections to coupled/non-coupled force splitting

Browse Source 
The integration splitting implemented in commit a5806207 has been shown
to be incorrect in some cases. A new procedure has been implemented
which can correctly split the implicit-explicit integral into a number
of pieces, in order to calculate the contribution of each. This is
intended for integrating coupled and non-coupled particle momentum and
heat transfers.

However, currently there is only ever one implicit coefficient used in
these transfers (there is no implicit non-coupled contribution). The
evaluation has therefore been short-cutted to only do the integration
with respect to the coupled contributions. The splitting functionality
has been retained in case additional separate implicit coefficients are
required in the future.

This change was made with help from Timo Niemi, VTT
This resolves bug report https://bugs.openfoam.org/view.php?id=2666
This commit is contained in:
Will Bainbridge
2017-11-28 10:23:52 +00:00
parent 9a3789f45a
commit 9ee3bbc943
8 changed files with 125 additions and 31 deletions
Download Patch File Download Diff File Expand all files Collapse all files

10
src/lagrangian/intermediate/integrationScheme/Euler/Euler.C
View File

@ -61,4 +61,14 @@ Foam::scalar Foam::integrationSchemes::Euler::dtEff
}
Foam::scalar Foam::integrationSchemes::Euler::sumDtEff
(
const scalar dt,
const scalar Beta
) const
{
return sqr(dt)/(1 + Beta*dt);
}
// ************************************************************************* //

4
src/lagrangian/intermediate/integrationScheme/Euler/Euler.H
View File

@ -79,6 +79,10 @@ public:
//- Return the integration effective time step
virtual scalar dtEff(const scalar dt, const scalar Beta) const;
//- Return the integral of the effective time step (using an Euler
// integration method)
virtual scalar sumDtEff(const scalar dt, const scalar Beta) const;
};

18
src/lagrangian/intermediate/integrationScheme/analytical/analytical.C
View File

@ -57,7 +57,23 @@ Foam::scalar Foam::integrationSchemes::analytical::dtEff
const scalar Beta
) const
{
return mag(Beta*dt) > SMALL ? (1 - exp(- Beta*dt))/Beta : dt;
return
mag(Beta*dt) > SMALL
? (1 - exp(- Beta*dt))/Beta
: dt;
}
Foam::scalar Foam::integrationSchemes::analytical::sumDtEff
(
const scalar dt,
const scalar Beta
) const
{
return
mag(Beta*dt) > SMALL
? dt/Beta - (1 - exp(- Beta*dt))/sqr(Beta)
: sqr(dt);
}

3
src/lagrangian/intermediate/integrationScheme/analytical/analytical.H
View File

@ -79,6 +79,9 @@ public:
//- Return the integration effective time step
virtual scalar dtEff(const scalar dt, const scalar Beta) const;
//- Return the integral of the effective time step
virtual scalar sumDtEff(const scalar dt, const scalar Beta) const;
};

48
src/lagrangian/intermediate/integrationScheme/integrationScheme/integrationScheme.H
View File

@ -34,7 +34,8 @@ Description
The methods are defined in terms of the effective time-step \f$\Delta
t_e\f$ by which the explicit rate is multipled. The effective time-step is
a function of the actual time step and the implicit coefficient.
a function of the actual time step and the implicit coefficient, which must
be implemented in each derived scheme.
\f[
\Delta t_e = f(\Delta t, B)
@ -49,14 +50,20 @@ Description
integration can be split up to detemine the effect of each contribution.
\f[
\frac{d \phi}{d t} = \left( \sum \alpha_i \right) -
\left( \sum \beta_i \right) \phi
\frac{d \phi_i}{d t} = \alpha_i - \beta_i \phi
\f]
\f[
\Delta \phi_i = (\alpha_i - \beta_i \phi^n) \Delta t_e
\Delta \phi_i = \alpha_i \Delta t -
\beta_i \int_0^{\Delta t} \phi d t
\f]
\f[
\Delta \phi_i = (\alpha_i - \beta_i \phi^n) \Delta t -
(A - B \phi^n) \int_0^{\Delta t} t_e dt
\f]
Note that the effective time step \f$\Delta t_e\f$ is not split up.
These partial calculations are defined in terms of the integral of the
effective time-step, \f$\int_0^{\Delta t} t_e dt\f$, which is also
implemented in every derivation.
SourceFiles
integrationScheme.C
@ -124,9 +131,19 @@ public:
// Member Functions
//- Perform the full integration
//- Perform the integration explicitly
template<class Type>
Type Delta
static Type explicitDelta
(
const Type& phi,
const scalar dtEff,
const Type& Alpha,
const scalar Beta
);
//- Perform the integration
template<class Type>
Type delta
(
const Type& phi,
const scalar dt,
@ -134,18 +151,23 @@ public:
const scalar Beta
) const;
//- Perform a stage of the integration
//- Perform a part of the integration
template<class Type>
static Type delta
Type partialDelta
(
const Type& phi,
const scalar dtEff,
const Type& alpha,
const scalar beta
);
const scalar dt,
const Type& Alpha,
const scalar Beta,
const Type& alphai,
const scalar betai
) const;
//- Return the integration effective time step
virtual scalar dtEff(const scalar dt, const scalar Beta) const = 0;
//- Return the integral of the effective time step
virtual scalar sumDtEff(const scalar dt, const scalar Beta) const = 0;
};

35
src/lagrangian/intermediate/integrationScheme/integrationScheme/integrationSchemeTemplates.C
View File

@ -28,15 +28,15 @@ License
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
template<class Type>
inline Type Foam::integrationScheme::Delta
inline Type Foam::integrationScheme::explicitDelta
(
const Type& phi,
const scalar dt,
const scalar dtEff,
const Type& Alpha,
const scalar Beta
) const
)
{
return delta(phi, dtEff(dt, Beta), Alpha, Beta);
return (Alpha - Beta*phi)*dtEff;
}
@ -44,12 +44,29 @@ template<class Type>
inline Type Foam::integrationScheme::delta
(
const Type& phi,
const scalar dtEff,
const Type& alpha,
const scalar beta
)
const scalar dt,
const Type& Alpha,
const scalar Beta
) const
{
return (alpha - beta*phi)*dtEff;
return explicitDelta(phi, dtEff(dt, Beta), Alpha, Beta);
}
template<class Type>
inline Type Foam::integrationScheme::partialDelta
(
const Type& phi,
const scalar dt,
const Type& Alpha,
const scalar Beta,
const Type& alphai,
const scalar betai
) const
{
return
explicitDelta(phi, dt, alphai, betai)
- explicitDelta(phi, sumDtEff(dt, Beta), Alpha, Beta)*betai;
}

30
src/lagrangian/intermediate/parcels/Templates/KinematicParcel/KinematicParcel.C
View File

@ -179,6 +179,8 @@ const Foam::vector Foam::KinematicParcel<ParcelType>::calcVelocity
const forceSuSp Fncp = forces.calcNonCoupled(p, ttd, dt, mass, Re, mu);
const scalar massEff = forces.massEff(p, ttd, mass);
/*
// Proper splitting ...
// Calculate the integration coefficients
const vector acp = (Fcp.Sp()*td.Uc() + Fcp.Su())/massEff;
const vector ancp = (Fncp.Sp()*td.Uc() + Fncp.Su() + Su)/massEff;
@ -186,12 +188,32 @@ const Foam::vector Foam::KinematicParcel<ParcelType>::calcVelocity
const scalar bncp = Fncp.Sp()/massEff;
// Integrate to find the new parcel velocity
const scalar dtEff = cloud.UIntegrator().dtEff(dt, bcp + bncp);
const vector deltaUcp = integrationScheme::delta(U_, dtEff, acp, bcp);
const vector deltaUncp = integrationScheme::delta(U_, dtEff, ancp, bncp);
const vector deltaUcp =
cloud.UIntegrator().partialDelta
(
U_, dt, acp + ancp, bcp + bncp, acp, bcp
);
const vector deltaUncp =
cloud.UIntegrator().partialDelta
(
U_, dt, acp + ancp, bcp + bncp, ancp, bncp
);
const vector deltaT = deltaUcp + deltaUncp;
*/
// Shortcut splitting assuming no implicit non-coupled force ...
// Calculate the integration coefficients
const vector acp = (Fcp.Sp()*td.Uc() + Fcp.Su())/massEff;
const vector ancp = (Fncp.Su() + Su)/massEff;
const scalar bcp = Fcp.Sp()/massEff;
// Integrate to find the new parcel velocity
const vector deltaU = cloud.UIntegrator().delta(U_, dt, acp + ancp, bcp);
const vector deltaUncp = ancp*dt;
const vector deltaUcp = deltaU - deltaUncp;
// Calculate the new velocity and the momentum transfer terms
vector Unew = U_ + deltaUcp + deltaUncp;
vector Unew = U_ + deltaU;
dUTrans -= massEff*deltaUcp;

8
src/lagrangian/intermediate/parcels/Templates/ThermoParcel/ThermoParcel.C
View File

@ -286,12 +286,12 @@ Foam::scalar Foam::ThermoParcel<ParcelType>::calcHeatTransfer
ancp /= m*Cp_;
// Integrate to find the new parcel temperature
const scalar dtEff = cloud.TIntegrator().dtEff(dt, bcp);
const scalar deltaTcp = integrationScheme::delta(T_, dtEff, acp, bcp);
const scalar deltaTncp = integrationScheme::delta(T_, dtEff, ancp, 0);
const scalar deltaT = cloud.TIntegrator().delta(T_, dt, acp + ancp, bcp);
const scalar deltaTncp = ancp*dt;
const scalar deltaTcp = deltaT - deltaTncp;
// Calculate the new temperature and the enthalpy transfer terms
scalar Tnew = T_ + deltaTcp + deltaTncp;
scalar Tnew = T_ + deltaT;
Tnew = min(max(Tnew, cloud.constProps().TMin()), cloud.constProps().TMax());
dhsTrans -= m*Cp_*deltaTcp;