add CANTERA_THERMO to control method to calculate h cp and update T. CANTERA_THERMO=1(current default) means using Cantera. CANTERA_THERMO=0 means using deepflame. (#555)

* change energy from ha to hs

* Update README.md

* Update README.md

* Update CanteraMixture.C

* change ha to h in fvSolution

* fix

* fix bug

* rever commit: change ha to hs

* fix EEqn

* Increase the test error of the dfLowMachFoam_2DaachenBomb case to 30.

---------

Co-authored-by: maorz1998 <[email protected]>
Co-authored-by: Zhi Chen <[email protected]>
Co-authored-by: pkuLmq <[email protected]>
Co-authored-by: OpenFOAMFans <[email protected]>

Author: 5 people

README.md

@@ -29,7 +29,7 @@ Detailed guide for installation and tutorials is available on [our documentation
 ## Features
 New in v1.5 (2025/1/8):
 - Provide a new boundary condition, totalFlowRateAdvectiveDiffusion (adopted from OpenFOAM v7), which accounts for diffusion effects at the boundary
-- Develop a new solver, `dfBuoyancyFoam`(adopted from fireFoam), a transient, compressible solver designed to model turbulent reacting flows while incorporating buoyancy effects
+- Develop a new solver, `dfBuoyancyFoam`(adopted from fireFoam in OpenFOAM v7), a transient, compressible solver designed to model turbulent reacting flows incorporating buoyancy effects
 - Add radiation models, which are integrated into the current `dfBuoyancyFoam` solver, and can be incorporated into other solvers if needed
 - Add new combustion models (including infinitelyFastChemistry and eddyDissipationModel) to the `dfBuoyancyFoam` solver
 - Imply DeepFlame on Arm architecture platforms

applications/solvers/df0DFoam/Make/options

@@ -5,6 +5,7 @@ EXE_INC = -std=c++14 \
     -Wno-unused-but-set-variable \
     -Wno-old-style-cast \
     $(PFLAGS) $(PINC) \
+	$(if $(CANTERA_THERMO),-DCANTERA_THERMO,) \
     $(if $(LIBTORCH_ROOT),-DUSE_LIBTORCH,) \
     $(if $(PYTHON_INC_DIR),-DUSE_PYTORCH,) \
     -I$(CANTERA_ROOT)/include \

applications/solvers/dfHighSpeedFoam/Make/options

@@ -5,6 +5,7 @@ EXE_INC = -std=c++14 \
     -Wno-unused-but-set-variable \
     -Wno-old-style-cast \
     $(PFLAGS) $(PINC) \
+	$(if $(CANTERA_THERMO),-DCANTERA_THERMO,) \
     $(if $(LIBTORCH_ROOT),-DUSE_LIBTORCH,) \
     $(if $(PYTHON_INC_DIR),-DUSE_PYTORCH,) \
     -I$(FOAM_APP)/solvers/compressible/rhoCentralFoam/BCs/lnInclude \

applications/solvers/dfHighSpeedFoam/rhoEEqn.H

@@ -39,7 +39,7 @@ else
 ea = rhoE/rho - 0.5*magSqr(U);
 ea.correctBoundaryConditions();
 
-ha = ea + p/rho;
+// ha = ea + p/rho;
 chemistry->correctThermo(); // before this, we must update ha = e + p/rho
 
 rhoE.boundaryFieldRef() == rho.boundaryField()*(ea.boundaryField() + 0.5*magSqr(U.boundaryField()));
@@ -49,15 +49,15 @@ if (!inviscid)
     if((ddtSchemes == "RK2SSP") || (ddtSchemes == "RK3SSP"))
     {
         rhoE_rhs =   fvc::laplacian(turbulence->alphaEff()*thermo.gamma(), ea)
-                  - diffAlphaD 
+                  - diffAlphaD
                   + fvc::div(hDiffCorrFlux);
-        
+
         rhoE +=  rkcoe3[nrk]*rhoE_rhs*runTime.deltaT();
 
         ea = rhoE/rho - 0.5*magSqr(U);
         ea.correctBoundaryConditions();
 
-        ha = ea + p/rho;
+        // ha = ea + p/rho;
         chemistry->correctThermo(); // before this, we must update ha = e + p/rho
     }
     else
@@ -84,12 +84,12 @@ if (!inviscid)
 
         eEqn.solve("ea");
 
-        ha = ea + p/rho;
+        // ha = ea + p/rho;
         chemistry->correctThermo();
         rhoE = rho*(ea + 0.5*magSqr(U));
         rhoE.boundaryFieldRef() == rho.boundaryField()*(ea.boundaryField() + 0.5*magSqr(U.boundaryField()));
     }
-    
+
 }
 
 Info<< "min/max(T) = "

applications/solvers/dfLowMachFoam/EEqn.H

@@ -1,5 +1,7 @@
 {
     volScalarField& he = thermo.he();
+
+
     tmp<fvScalarMatrix> thcSource(new fvScalarMatrix(he, dimEnergy/dimTime));
     fvScalarMatrix& hcSource = thcSource.ref();
     forAll(Y, i)
@@ -8,6 +10,7 @@
     }
     //hSource.source() -= parcels.hsTrans()/runTime.deltaT();// equivalent to parcels.Sh(he)
 
+
     start1 = std::clock();
     fvScalarMatrix EEqn
     (
@@ -16,33 +19,43 @@
     +   fvc::ddt(rho, K) + fvc::div(phi, K)
     -   dpdt
     ==
+    (
+        turbName == "laminar"
+        ?
         (
-            turbName == "laminar"
-            ?
-            (
-                fvm::laplacian(turbulence->alpha(), he)
-            -   diffAlphaD
-            +   fvc::div(hDiffCorrFlux)
-            // +   rho*(U&g)
-            +   parcels.Sh(he)
-            +   hcSource
-            )
-            :
-            (
-                fvm::laplacian(turbulence->alphaEff(), he)
-            // +   rho*(U&g)
-            +   parcels.Sh(he)
-            +   hcSource
-            )
+            fvm::laplacian(turbulence->alpha(), he)
+        -   diffAlphaD
+        +   fvc::div(hDiffCorrFlux)
         )
+        :
+        fvm::laplacian(turbulence->alphaEff(), he)
+    )
+    +   parcels.Sh(he)
+    // +   rho*(U&g)
     );
+
+    if(CanteraMixture::getEnergyName() == word("hs"))
+    {
+        EEqn -= combustion->Qdot();
+    }
+    else if(CanteraMixture::getEnergyName() == word("ha"))
+    {
+        EEqn -= hcSource;
+    }
+    else
+    {
+        FatalErrorInFunction
+            << "Wrong CanteraMixture enregy name: " << CanteraMixture::getEnergyName()
+            << abort(FatalError);
+    }
+
     end1 = std::clock();
     time_monitor_EEqn += double(end1 - start1) / double(CLOCKS_PER_SEC);
     time_monitor_EEqn_mtxAssembly += double(end1 - start1) / double(CLOCKS_PER_SEC);
 
     // EEqn.relax();
     start1 = std::clock();
-    EEqn.solve("ha");
+    EEqn.solve("h");
     end1 = std::clock();
     time_monitor_EEqn += double(end1 - start1) / double(CLOCKS_PER_SEC);
     time_monitor_EEqn_solve += double(end1 - start1) / double(CLOCKS_PER_SEC);

applications/solvers/dfLowMachFoam/EEqn_GPU.H

@@ -37,7 +37,7 @@
         )
     );
     // EEqn.relax();
-    EEqn.solve("ha");
+    EEqn.solve("h");
     // checkResult
     // TODO: for temp, now we compare ldu, finally we compare csr
     std::vector<double> h_internal_coeffs(dfDataBase.num_boundary_surfaces);

applications/solvers/dfLowMachFoam/Make/options

@@ -8,6 +8,7 @@ EXE_INC = -std=c++14 \
     -Wno-old-style-cast \
     -I. \
     $(PFLAGS) $(PINC) \
+	$(if $(CANTERA_THERMO),-DCANTERA_THERMO,) \
     $(if $(LIBTORCH_ROOT),-DUSE_LIBTORCH,) \
     $(if $(PYTHON_INC_DIR),-DUSE_PYTORCH,) \
     -I$(FOAM_APP)/solvers/lagrangian/reactingParcelFoam \

applications/solvers/dfSprayFoam/EEqn.H

@@ -9,8 +9,8 @@
         hcSource.source() -= parcels.rhoTrans(i)*chemistry->mixture().Hc(i)/runTime.deltaT();
     }
     //hSource.source() -= parcels.hsTrans()/runTime.deltaT();// equivalent to parcels.Sh(he)
-    
-    
+
+
     fvScalarMatrix EEqn
     (
         turbName == "laminar"
@@ -22,9 +22,8 @@
     -   fvm::laplacian(turbulence->alphaEff(), he)
     +   diffAlphaD
     ==
-        rho*(U&g)
+       rho*(U&g)
     +   parcels.Sh(he)
-    +   hcSource
     +   fvc::div(hDiffCorrFlux)
     //+ radiation->Sh(thermo, he)
     //+ fvOptions(rho, he)
@@ -38,20 +37,34 @@
     ==
         rho*(U&g)
     +   parcels.Sh(he)
-    +   hcSource
     //+ radiation->Sh(thermo, he)
     //+ fvOptions(rho, he)
         )
     );
 
+    if(CanteraMixture::getEnergyName() == word("hs"))
+    {
+        EEqn -= combustion->Qdot();
+    }
+    else if(CanteraMixture::getEnergyName() == word("ha"))
+    {
+        EEqn -= hcSource;
+    }
+    else
+    {
+        FatalErrorInFunction
+            << "Wrong CanteraMixture enregy name: " << CanteraMixture::getEnergyName()
+            << abort(FatalError);
+    }
+
     EEqn.relax();
 
     //fvOptions.constrain(EEqn);
 
-    EEqn.solve("ha");
+    EEqn.solve("h");
 
     //fvOptions.correct(he);
 
     //radiation->correct();
-    
+
 }

applications/solvers/dfSprayFoam/Make/options

@@ -35,6 +35,7 @@ EXE_INC = -std=c++14 \
     -I$(CANTERA_ROOT)/include \
     $(if $(LIBTORCH_ROOT),-I$(LIBTORCH_ROOT)/include,) \
     $(if $(LIBTORCH_ROOT),-I$(LIBTORCH_ROOT)/include/torch/csrc/api/include,) \
+	$(if $(CANTERA_THERMO),-DCANTERA_THERMO,) \
     $(if $(LIBTORCH_ROOT),-DUSE_LIBTORCH,) \
     $(if $(PYTHON_INC_DIR),-DUSE_PYTORCH,) \
     $(PYTHON_INC_DIR)

bashrc.in

@@ -5,6 +5,10 @@ export LIBTORCH_ROOT=LIBTORCH_DIR
 export PYTHON_INC_DIR="PYTORCH_INC"
 export PYTHON_LIB_DIR="PYTORCH_LIB"
 export CANTERA_ROOT=LIBCANTERA_DIR
+export CANTERA_MECH=1
+export CANTERA_REACTION=1
+export CANTERA_THERMO=1
+export CANTERA_TRANSPORT=1
 export CANTERA_DATA=$CANTERA_ROOT/share/cantera/data
 export LD_LIBRARY_PATH=$LIBTORCH_ROOT/lib:$LD_LIBRARY_PATH
 export LD_LIBRARY_PATH=$CANTERA_ROOT/lib:$LD_LIBRARY_PATH