From b7c7bc57a13e59d881a0e87d34441c809c18d782 Mon Sep 17 00:00:00 2001 From: pkuLmq <1900011062@pku.edu.cn> Date: Tue, 26 Jul 2022 11:01:49 +0800 Subject: [PATCH 1/4] correct fvc to fvm in rhoEEqn.H and add inviscid switch --- .../solvers/dfHighSpeedFoam/createFields.H | 21 ++++++++--- .../solvers/dfHighSpeedFoam/rhoEEqn.H | 2 +- .../constant/thermophysicalProperties | 2 +- .../det_1D_H2_01mm/system/fvSchemes | 26 +++++++------- .../det_1D_H2_01mm/system/fvSolution | 36 +++++++++---------- 5 files changed, 49 insertions(+), 38 deletions(-) diff --git a/applications/solvers/dfHighSpeedFoam/createFields.H b/applications/solvers/dfHighSpeedFoam/createFields.H index 5dded90b..55ea3a79 100644 --- a/applications/solvers/dfHighSpeedFoam/createFields.H +++ b/applications/solvers/dfHighSpeedFoam/createFields.H @@ -12,11 +12,22 @@ const volScalarField& T = thermo.T(); const volScalarField& psi = thermo.psi(); const volScalarField& mu = thermo.mu(); -bool inviscid(true); -if (max(mu.primitiveField()) > 0.0) -{ - inviscid = false; -} +dictionary thermoDict +( + IOdictionary + ( + IOobject + ( + "thermophysicalProperties", + runTime.constant(), + mesh, + IOobject::MUST_READ_IF_MODIFIED, + IOobject::NO_WRITE + ) + ) +); + +bool inviscid(thermoDict.lookup("inviscid")); Info<< "Reading field U\n" << endl; volVectorField U diff --git a/applications/solvers/dfHighSpeedFoam/rhoEEqn.H b/applications/solvers/dfHighSpeedFoam/rhoEEqn.H index c1ba18d6..b69e0095 100644 --- a/applications/solvers/dfHighSpeedFoam/rhoEEqn.H +++ b/applications/solvers/dfHighSpeedFoam/rhoEEqn.H @@ -29,7 +29,7 @@ if (!inviscid) ( fvm::ddt(rho, e) - fvc::ddt(rho, e) // alpha in deepflame is considered to calculate h by default (kappa/Cp), so multiply gamma to correct alpha - - fvc::laplacian(turbulence->alphaEff()*thermo.gamma(), e) + - fvm::laplacian(turbulence->alphaEff()*thermo.gamma(), e) + diffAlphaD == fvc::div(hDiffCorrFlux) diff --git a/examples/dfHighSpeedFoam/det_1D_H2_01mm/constant/thermophysicalProperties b/examples/dfHighSpeedFoam/det_1D_H2_01mm/constant/thermophysicalProperties index 7f207bb7..e983a193 100755 --- a/examples/dfHighSpeedFoam/det_1D_H2_01mm/constant/thermophysicalProperties +++ b/examples/dfHighSpeedFoam/det_1D_H2_01mm/constant/thermophysicalProperties @@ -14,5 +14,5 @@ FoamFile object thermophysicalProperties; } // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // - +inviscid true; // ************************************************************************* // diff --git a/examples/dfHighSpeedFoam/det_1D_H2_01mm/system/fvSchemes b/examples/dfHighSpeedFoam/det_1D_H2_01mm/system/fvSchemes index 3f654d84..b710a8e5 100755 --- a/examples/dfHighSpeedFoam/det_1D_H2_01mm/system/fvSchemes +++ b/examples/dfHighSpeedFoam/det_1D_H2_01mm/system/fvSchemes @@ -19,39 +19,39 @@ fluxScheme Kurganov; ddtSchemes { - default Euler; + default Euler; } gradSchemes { - default Gauss linear; + default Gauss linear; } divSchemes { - default none; - div(tauMC) Gauss linear; - div(hDiffCorrFlux) Gauss cubic; - div(phi,Yi_h) Gauss vanLeer; + default none; + div(tauMC) Gauss linear; + div(hDiffCorrFlux) Gauss cubic; + div(phi,Yi_h) Gauss vanLeer; } laplacianSchemes { - default Gauss linear uncorrected; + default Gauss linear uncorrected; } interpolationSchemes { - default linear; - reconstruct(rho) Minmod; - reconstruct(U) MinmodV; - reconstruct(T) Minmod; - reconstruct(Yi) Minmod; + default linear; + reconstruct(rho) Minmod; + reconstruct(U) MinmodV; + reconstruct(T) Minmod; + reconstruct(Yi) Minmod; } snGradSchemes { - default uncorrected; + default uncorrected; } diff --git a/examples/dfHighSpeedFoam/det_1D_H2_01mm/system/fvSolution b/examples/dfHighSpeedFoam/det_1D_H2_01mm/system/fvSolution index 17c77b5d..7ee0f7f0 100755 --- a/examples/dfHighSpeedFoam/det_1D_H2_01mm/system/fvSolution +++ b/examples/dfHighSpeedFoam/det_1D_H2_01mm/system/fvSolution @@ -19,41 +19,41 @@ solvers { U { - //solver smoothSolver; - //smoother GaussSeidel; - //nSweeps 2; - //tolerance 1e-17; - //relTol 0; + //solver smoothSolver; + //smoother GaussSeidel; + //nSweeps 2; + //tolerance 1e-17; + //relTol 0; - solver PBiCGStab; - preconditioner DIC; - tolerance 1e-11; - relTol 0; + solver PBiCGStab; + preconditioner DIC; + tolerance 1e-11; + relTol 0; } h { - $U; - tolerance 1e-11; - relTol 0; + $U; + tolerance 1e-11; + relTol 0; } e { - $U; - tolerance 1e-11; - relTol 0; + $U; + tolerance 1e-11; + relTol 0; } "rho.*" { - solver diagonal; + solver diagonal; } "(O2|N2|Yi)" { - solver PBiCGStab; - preconditioner DILU; + solver PBiCGStab; + preconditioner DILU; tolerance 1e-11; relTol 0; } From 0fe9cee9eb0add945c248ee9a407b6d674e166c2 Mon Sep 17 00:00:00 2001 From: pkuLmq <1900011062@pku.edu.cn> Date: Tue, 26 Jul 2022 11:06:08 +0800 Subject: [PATCH 2/4] correct fvc to fvm in rhoEEqn.H and add inviscid switch --- examples/dfHighSpeedFoam/det_1D_H2_01mm/system/fvSchemes | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/examples/dfHighSpeedFoam/det_1D_H2_01mm/system/fvSchemes b/examples/dfHighSpeedFoam/det_1D_H2_01mm/system/fvSchemes index b710a8e5..15ea6d13 100755 --- a/examples/dfHighSpeedFoam/det_1D_H2_01mm/system/fvSchemes +++ b/examples/dfHighSpeedFoam/det_1D_H2_01mm/system/fvSchemes @@ -46,7 +46,7 @@ interpolationSchemes reconstruct(rho) Minmod; reconstruct(U) MinmodV; reconstruct(T) Minmod; - reconstruct(Yi) Minmod; + reconstruct(Yi) Minmod; } snGradSchemes From e784571231d12adf601640fea7fb18bea5b5ecb5 Mon Sep 17 00:00:00 2001 From: pkuLmq <1900011062@pku.edu.cn> Date: Tue, 26 Jul 2022 17:33:08 +0800 Subject: [PATCH 3/4] fix inviscid in rhoYEqn.H and rhoEEqn.H --- .../solvers/dfHighSpeedFoam/rhoYEqn.H | 29 ++++++++++--------- 1 file changed, 15 insertions(+), 14 deletions(-) diff --git a/applications/solvers/dfHighSpeedFoam/rhoYEqn.H b/applications/solvers/dfHighSpeedFoam/rhoYEqn.H index f0a50aec..105042db 100644 --- a/applications/solvers/dfHighSpeedFoam/rhoYEqn.H +++ b/applications/solvers/dfHighSpeedFoam/rhoYEqn.H @@ -1,7 +1,3 @@ -hDiffCorrFlux = Zero; -diffAlphaD = Zero; -sumYDiffError = Zero; - tmp> mvConvection ( fv::convectionScheme::New @@ -13,12 +9,6 @@ tmp> mvConvection ) ); -forAll(Y, i) -{ - sumYDiffError += chemistry.rhoD(i)*fvc::grad(Y[i]); -} -const surfaceScalarField phiUc = linearInterpolate(sumYDiffError) & mesh.Sf(); - { start = std::clock(); chemistry.solve(mesh.time().deltaTValue()); @@ -34,12 +24,9 @@ const surfaceScalarField phiUc = linearInterpolate(sumYDiffError) & mesh.Sf(); forAll(Y, i) { volScalarField& Yi = Y[i]; - hDiffCorrFlux += chemistry.hai(i)*(chemistry.rhoD(i)*fvc::grad(Yi) - Yi*sumYDiffError); - diffAlphaD += fvc::laplacian(thermo.alpha()*chemistry.hai(i), Yi); if (i != inertIndex) - { - tmp DEff = chemistry.rhoD(i) + turbulence->mut()/Sct; + { fvScalarMatrix YiEqn ( fvm::ddt(rho, Yi) @@ -50,6 +37,20 @@ const surfaceScalarField phiUc = linearInterpolate(sumYDiffError) & mesh.Sf(); if (!inviscid) { + hDiffCorrFlux = Zero; + diffAlphaD = Zero; + sumYDiffError = Zero; + + forAll(Y, i) + { + sumYDiffError += chemistry.rhoD(i)*fvc::grad(Y[i]); + } + const surfaceScalarField phiUc = linearInterpolate(sumYDiffError) & mesh.Sf(); + + hDiffCorrFlux += chemistry.hai(i)*(chemistry.rhoD(i)*fvc::grad(Yi) - Yi*sumYDiffError); + diffAlphaD += fvc::laplacian(thermo.alpha()*chemistry.hai(i), Yi); + tmp DEff = chemistry.rhoD(i) + turbulence->mut()/Sct; + YiEqn -= fvm::laplacian(DEff(), Yi) - mvConvection->fvmDiv(phiUc, Yi); } From 89e6a09abe3060e677b5e8a01e921acd8490f0c1 Mon Sep 17 00:00:00 2001 From: pkuLmq <1900011062@pku.edu.cn> Date: Wed, 27 Jul 2022 17:22:16 +0800 Subject: [PATCH 4/4] correct rhoYEqn.H --- applications/solvers/dfHighSpeedFoam/rhoYEqn.H | 11 +++++++---- 1 file changed, 7 insertions(+), 4 deletions(-) diff --git a/applications/solvers/dfHighSpeedFoam/rhoYEqn.H b/applications/solvers/dfHighSpeedFoam/rhoYEqn.H index 105042db..2a19deec 100644 --- a/applications/solvers/dfHighSpeedFoam/rhoYEqn.H +++ b/applications/solvers/dfHighSpeedFoam/rhoYEqn.H @@ -1,3 +1,10 @@ +if (!inviscid) +{ + hDiffCorrFlux = Zero; + diffAlphaD = Zero; + sumYDiffError = Zero; +} + tmp> mvConvection ( fv::convectionScheme::New @@ -37,10 +44,6 @@ tmp> mvConvection if (!inviscid) { - hDiffCorrFlux = Zero; - diffAlphaD = Zero; - sumYDiffError = Zero; - forAll(Y, i) { sumYDiffError += chemistry.rhoD(i)*fvc::grad(Y[i]);