diff --git a/src/lagrangian/intermediate/integrationScheme/Euler/Euler.C b/src/lagrangian/intermediate/integrationScheme/Euler/Euler.C index a0fa45f7b7..3380e1f361 100644 --- a/src/lagrangian/intermediate/integrationScheme/Euler/Euler.C +++ b/src/lagrangian/intermediate/integrationScheme/Euler/Euler.C @@ -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); +} + + // ************************************************************************* // diff --git a/src/lagrangian/intermediate/integrationScheme/Euler/Euler.H b/src/lagrangian/intermediate/integrationScheme/Euler/Euler.H index 21c417efa0..5ae87cbda9 100644 --- a/src/lagrangian/intermediate/integrationScheme/Euler/Euler.H +++ b/src/lagrangian/intermediate/integrationScheme/Euler/Euler.H @@ -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; }; diff --git a/src/lagrangian/intermediate/integrationScheme/analytical/analytical.C b/src/lagrangian/intermediate/integrationScheme/analytical/analytical.C index ed9362d670..4c9dc9bfd7 100644 --- a/src/lagrangian/intermediate/integrationScheme/analytical/analytical.C +++ b/src/lagrangian/intermediate/integrationScheme/analytical/analytical.C @@ -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); } diff --git a/src/lagrangian/intermediate/integrationScheme/analytical/analytical.H b/src/lagrangian/intermediate/integrationScheme/analytical/analytical.H index 5bfef1df4b..ea1176b4cb 100644 --- a/src/lagrangian/intermediate/integrationScheme/analytical/analytical.H +++ b/src/lagrangian/intermediate/integrationScheme/analytical/analytical.H @@ -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; }; diff --git a/src/lagrangian/intermediate/integrationScheme/integrationScheme/integrationScheme.H b/src/lagrangian/intermediate/integrationScheme/integrationScheme/integrationScheme.H index ead1eb98b0..01c41e6a7d 100644 --- a/src/lagrangian/intermediate/integrationScheme/integrationScheme/integrationScheme.H +++ b/src/lagrangian/intermediate/integrationScheme/integrationScheme/integrationScheme.H @@ -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 - Type Delta + static Type explicitDelta + ( + const Type& phi, + const scalar dtEff, + const Type& Alpha, + const scalar Beta + ); + + //- Perform the integration + template + 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 - 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; }; diff --git a/src/lagrangian/intermediate/integrationScheme/integrationScheme/integrationSchemeTemplates.C b/src/lagrangian/intermediate/integrationScheme/integrationScheme/integrationSchemeTemplates.C index 7202f2c126..0c7365ba06 100644 --- a/src/lagrangian/intermediate/integrationScheme/integrationScheme/integrationSchemeTemplates.C +++ b/src/lagrangian/intermediate/integrationScheme/integrationScheme/integrationSchemeTemplates.C @@ -28,15 +28,15 @@ License // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // template -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 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 +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; } diff --git a/src/lagrangian/intermediate/parcels/Templates/KinematicParcel/KinematicParcel.C b/src/lagrangian/intermediate/parcels/Templates/KinematicParcel/KinematicParcel.C index e654b56d3f..26206c4540 100644 --- a/src/lagrangian/intermediate/parcels/Templates/KinematicParcel/KinematicParcel.C +++ b/src/lagrangian/intermediate/parcels/Templates/KinematicParcel/KinematicParcel.C @@ -179,6 +179,8 @@ const Foam::vector Foam::KinematicParcel::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::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; diff --git a/src/lagrangian/intermediate/parcels/Templates/ThermoParcel/ThermoParcel.C b/src/lagrangian/intermediate/parcels/Templates/ThermoParcel/ThermoParcel.C index d44526e750..5d0dfec41e 100644 --- a/src/lagrangian/intermediate/parcels/Templates/ThermoParcel/ThermoParcel.C +++ b/src/lagrangian/intermediate/parcels/Templates/ThermoParcel/ThermoParcel.C @@ -286,12 +286,12 @@ Foam::scalar Foam::ThermoParcel::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;