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/applications/solvers/dfHighSpeedFoam/rhoYEqn.H b/applications/solvers/dfHighSpeedFoam/rhoYEqn.H index f0a50aec..2a19deec 100644 --- a/applications/solvers/dfHighSpeedFoam/rhoYEqn.H +++ b/applications/solvers/dfHighSpeedFoam/rhoYEqn.H @@ -1,6 +1,9 @@ -hDiffCorrFlux = Zero; -diffAlphaD = Zero; -sumYDiffError = Zero; +if (!inviscid) +{ + hDiffCorrFlux = Zero; + diffAlphaD = Zero; + sumYDiffError = Zero; +} tmp> mvConvection ( @@ -13,12 +16,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 +31,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 +44,16 @@ const surfaceScalarField phiUc = linearInterpolate(sumYDiffError) & mesh.Sf(); if (!inviscid) { + 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); } 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..15ea6d13 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; }