diff --git a/src/dfCombustionModels/FGM/baseFGM/baseFGM.C b/src/dfCombustionModels/FGM/baseFGM/baseFGM.C
index 1f8cbca14..850656a23 100644
--- a/src/dfCombustionModels/FGM/baseFGM/baseFGM.C
+++ b/src/dfCombustionModels/FGM/baseFGM/baseFGM.C
@@ -283,6 +283,7 @@ Foam::combustionModels::baseFGM<ReactionThermo>::baseFGM
     ZcvarMax_(0.25),
     ZcvarMin_(-0.25),
     rho_(const_cast<volScalarField&>(this->mesh().objectRegistry::lookupObject<volScalarField>("rho"))),
+    rho_inRhoThermo_(dynamic_cast<rhoThermo&>(thermo).rho()),
     p_(this->thermo().p()),
     T_(this->thermo().T()),
     U_(this->mesh().objectRegistry::lookupObject<volVectorField>("U")),
diff --git a/src/dfCombustionModels/FGM/baseFGM/baseFGM.H b/src/dfCombustionModels/FGM/baseFGM/baseFGM.H
index 19eb378c1..59863fbd9 100644
--- a/src/dfCombustionModels/FGM/baseFGM/baseFGM.H
+++ b/src/dfCombustionModels/FGM/baseFGM/baseFGM.H
@@ -185,6 +185,8 @@ class baseFGM
 
     //- fluid density
      volScalarField& rho_;
+    //- fluid density
+     volScalarField& rho_inRhoThermo_;
 
     //- fluid pressure
     const volScalarField& p_;
diff --git a/src/dfCombustionModels/FGM/flareFGM/flareFGM.C b/src/dfCombustionModels/FGM/flareFGM/flareFGM.C
index 62c1b7f06..cafdf3357 100644
--- a/src/dfCombustionModels/FGM/flareFGM/flareFGM.C
+++ b/src/dfCombustionModels/FGM/flareFGM/flareFGM.C
@@ -436,6 +436,9 @@ void Foam::combustionModels::flareFGM<ReactionThermo>::retrieval()
     this->T_.max(TMin);  
     this->T_.min(TMax);  
 
+    dimensionedScalar R_uniGas("R_uniGas",dimensionSet(1,2,-2,-1,-1,0,0),8.314e3);
+    this->psi_ = this->Wt_/(R_uniGas*this->T_);
+  
     if(this->mesh().time().timeIndex() > 0)   
     {
         dimensionedScalar p_operateDim("p_operateDim", dimensionSet(1,-1,-2,0,0,0,0),this->incompPref_);  
@@ -449,9 +452,8 @@ void Foam::combustionModels::flareFGM<ReactionThermo>::retrieval()
             this->rho_ = this->p_*this->psi_;
         }
     }
+    this->rho_inRhoThermo_ = this->rho_;
 
-    dimensionedScalar R_uniGas("R_uniGas",dimensionSet(1,2,-2,-1,-1,0,0),8.314e3);
-    this->psi_ = this->Wt_/(R_uniGas*this->T_);
 }
 
 // ************************************************************************* //
diff --git a/test/corrtest.cpp b/test/corrtest.cpp
index a3bf3d407..3bb325793 100644
--- a/test/corrtest.cpp
+++ b/test/corrtest.cpp
@@ -56,11 +56,11 @@ TEST(corrtest,df0DFoam_H2){
 //}
 
 TEST(corrtest,dfLowMachFoam_TGV){
-    EXPECT_FLOAT_EQ(TGV500,1533.35);   // compare the maximum temperature along y direction in 2D TGV after 500 time steps
-    EXPECT_FLOAT_EQ(TGV400,1298.49);   //  ..........400 time steps
-    EXPECT_FLOAT_EQ(TGV300,874.988);
-    EXPECT_FLOAT_EQ(TGV200,541.233);
-    EXPECT_FLOAT_EQ(TGV100,363.754);
+    EXPECT_FLOAT_EQ(TGV500,1532.92);   // compare the maximum temperature along y direction in 2D TGV after 500 time steps
+    EXPECT_FLOAT_EQ(TGV400,1297.64);   //  ..........400 time steps
+    EXPECT_FLOAT_EQ(TGV300,871.092);
+    EXPECT_FLOAT_EQ(TGV200,537.614);
+    EXPECT_FLOAT_EQ(TGV100,363.504);
 }
 
 TEST(corrtest,2DSandia){