hehe2

Author: Dimitrios Adam

examples/1D_reactive_shocktube/case.py

@@ -57,9 +57,9 @@
     "p": 0,
     "dt": float(dt),
     "t_step_start": 0,
-    "t_step_stop": NT,
-    "t_step_save": NS,
-    "t_step_print": NS,
+    "t_step_stop": 1,
+    "t_step_save": 1,
+    "t_step_print": 1,
     "parallel_io": "F" if args.mfc.get("mpi", True) else "F",
     # Simulation Algorithm Parameters
     "model_eqns": 2,

examples/2D_riemann_test/case.py

@@ -28,8 +28,8 @@
             "n": 499,
             "p": 0,
             "dt": 4.0e-08,
-            "t_step_start": 0,
-            "t_step_stop": 10000,
+            "t_step_start": 10000,
+            "t_step_stop": 17000,
             "t_step_save": 1000,
             # Simulation Algorithm Parameters
             "num_patches": 1,

src/common/m_chemistry.fpp

@@ -146,6 +146,7 @@ contains
 
                     rho = q_cons_qp(contxe)%sf(x, y, z)
                     T = q_T_sf%sf(x, y, z)
+                    print *, T
 
                     call get_net_production_rates(rho, T, Ys, omega)
 
@@ -167,11 +168,12 @@ contains
 
     end subroutine s_compute_chemistry_reaction_flux
 
-    subroutine s_compute_chemistry_diffusion_flux(idir, q_prim_qp, flux_src_vf, irx, iry, irz)
+    subroutine s_compute_chemistry_diffusion_flux(idir, q_prim_qp, flux_src_vf, irx, iry, irz, q_T_sf)
 
         type(scalar_field), dimension(sys_size), intent(in) :: q_prim_qp
         type(scalar_field), dimension(sys_size), intent(inout) :: flux_src_vf
         type(int_bounds_info), intent(in) :: irx, iry, irz
+        type(scalar_field), intent(in) :: q_T_sf
 
         integer, intent(in) :: idir
         #:if not MFC_CASE_OPTIMIZATION and USING_AMD
@@ -248,15 +250,23 @@ contains
                             rho_L = q_prim_qp(1)%sf(x, y, z)
                             rho_R = q_prim_qp(1)%sf(x + offsets(1), y + offsets(2), z + offsets(3))
 
-                            T_L = P_L/rho_L/Rgas_L      
-                            T_R = P_R/rho_R/Rgas_R
+                            if ( x .eq. -1 .or. x .eq. 499 .or. y .eq. -1 .or. y .eq. 499) then
+                              T_L = P_L/rho_L/Rgas_L     
+                              T_R =  P_R/rho_R/Rgas_R
+
+                            else 
+                             T_L =q_T_sf%sf(x, y, z)! P_L/rho_L/Rgas_L      
+                              T_R = q_T_sf%sf(x + offsets(1), y + offsets(2), z + offsets(3))! P_R/rho_R/Rgas_R
+                            end if
 
+          
                             rho_cell = 0.5_wp*(rho_L + rho_R)
                             dT_dxi = (T_R - T_L)/grid_spacing
 
-                            if (x .eq. 200 .and. idir .eq. 2) then
-                            print *, dT_dxi, y , T_L, T_R
-                            end if
+                           ! if (x .eq. 200 .and. idir == 2) then
+                           ! print *, T_L, T_R, y_cc(y), y
+                           ! end if
+                         
 
                             ! Get transport properties
                             call get_species_mass_diffusivities_mixavg(P_L, T_L, Ys_L, mass_diffusivities_mixavg1)

src/common/m_variables_conversion.fpp

@@ -814,19 +814,19 @@ contains
                                             qv_K, rhoYks, pres, T, pres_mag=pres_mag)
 
 
-                       if (.not. (t_step == 0 .and. stage == 1)) then
+                 !     if (.not. (t_step == 0 .and. stage == 1)) then
                             ! Safety check: Ensure indices are within the allocated domain (0:m, 0:n, 0:p)
                             ! to prevent segfaults if the loop iterates over buffer zones.
-                            if (j >= 0 .and. j <= m .and. k >= 0 .and. k <= n .and. l >= 0 .and. l <= p) then
-                                if (ghost_points_index%sf(j, k, l) == 1) then
-                                 pres = pressure_ghost_point%sf(j, k, l)
-                                               if (j .eq. 200) then
-                               !print *, pres, k, y_cc(k), T
-                               end if
-                              end if
+                 !           if (j >= 0 .and. j <= m .and. k >= 0 .and. k <= n .and. l >= 0 .and. l <= p) then
+                 !               if (ghost_points_index%sf(j, k, l) == 1) then
+                 !                pres = pressure_ghost_point%sf(j, k, l)
+                  !                             if (j .eq. 200) then
+                  !             !print *, pres, k, y_cc(k), T
+                  !             end if
+                  !            end if
 
-                           end if
-                       end if
+                !           end if
+                !       end if
 
                         qK_prim_vf(E_idx)%sf(j, k, l) = pres
 

src/simulation/m_rhs.fpp

@@ -954,7 +954,7 @@ contains
                 ! RHS for diffusion
                 if (chemistry .and. chem_params%diffusion) then
                     call nvtxStartRange("RHS-CHEM-DIFFUSION")
-                    call s_compute_chemistry_diffusion_flux(id, q_prim_qp%vf, flux_src_n(id)%vf, irx, iry, irz)
+                    call s_compute_chemistry_diffusion_flux(id, q_prim_qp%vf, flux_src_n(id)%vf, irx, iry, irz,q_T_sf)
                     call nvtxEndRange
                 end if