diff --git a/Source/Make.package b/Source/Make.package index a59ff13..dad2a9c 100644 --- a/Source/Make.package +++ b/Source/Make.package @@ -1,4 +1,4 @@ CEXE_sources += main.cpp mflo.cpp mflo_evolve.cpp mflo_plot.cpp -CEXE_sources += mflo_gridwork.cpp mflo_timestep.cpp +CEXE_sources += mflo_gridwork.cpp mflo_timestep.cpp mflo_splitchem.cpp CEXE_headers += mflo.H kernels_3d.H mflo_tagging.H mflo_bcfill.H globalDefines.H mflo_constants.H diff --git a/Source/Src_3d/compute_flux_3d.H b/Source/Src_3d/compute_flux_3d.H index 6aee2c9..db3e2a3 100644 --- a/Source/Src_3d/compute_flux_3d.H +++ b/Source/Src_3d/compute_flux_3d.H @@ -188,6 +188,13 @@ void compute_flux(int i, int j, int k,int sweepdir, facevel[XDIR] = half*(phi(left,VELX_INDX) + phi(right,VELX_INDX)); facevel[YDIR] = half*(phi(left,VELY_INDX) + phi(right,VELY_INDX)); facevel[ZDIR] = half*(phi(left,VELZ_INDX) + phi(right,VELZ_INDX)); + + Real facevfrac; + facevfrac=half*(phi(left,VFRAC_INDX)+phi(right,VFRAC_INDX)); + + //note: + //if facevfrac >=0.5, this face is a gas cell-fractional cell interface + //if facevfrac < 0.5, this face is solid cell-fractional cell interface for(int c=0;c= 0.5) { - amrex::Abort("unknown hyperbolics scheme"); - } + if(hyperbolics_order==1) + { + ausmp_up_flux(ul,ur,fhalf,normal); + } + else if(hyperbolics_order==2) + { + get_higherorder_states(ul2,ur2,ul,ur,ulm1,urp1); + ausmp_up_flux(ul2,ur2,fhalf,normal); + } + else if(hyperbolics_order==4) + { + central_with_diss_flux(ul,ur,ulm1,urp1,fhalf,normal,dissfactor); + } + else if(hyperbolics_order==5) + { + weno_flux(ulm2,ulm1,ul,ur,urp1,urp2,fhalf,normal); + } + else + { + amrex::Abort("unknown hyperbolics scheme"); + } - for(int c=0;c const& flxy, Array4 const& flxz), const GpuArray& dx, - int spec_in_solid) + int spec_in_solid, + int conj_ht) { // remember, we are solve dudt + del.F = 0 dsdt(i, j, k, n) = (flxx(i, j, k, n) - flxx(i + 1, j, k, n)) / dx[0] + @@ -33,15 +34,41 @@ void update_residual( int nvars_to_nullify=(spec_in_solid==1)?FLO_NVARS:FLO_NVARS+NUM_SPECIES; - - if(n < nvars_to_nullify) + //check for rhoe because we may have conjugate heat transfer + if(n!=RHOE_INDX) + { + if(n < nvars_to_nullify) + { + if(phi(i,j,k,VFRAC_INDX) < one) + { + dsdt(i, j, k, n)=zeroval; + } + } + } + //if n is rhoe + else { - if(phi(i,j,k,VFRAC_INDX) < one) - { - dsdt(i, j, k, n)=zeroval; + //if no conjugate heat transfer + //null the dsdt + if(!conj_ht) + { + if(phi(i,j,k,VFRAC_INDX) < one) + { + dsdt(i, j, k, n)=zeroval; + } + } + else + { + //when conjugate heat transfer is there + //null the dsdt in fractional cells + //as these are recalculated from cut cells + if(phi(i,j,k,VFRAC_INDX) > zeroval && phi(i,j,k,VFRAC_INDX) < one) + { + dsdt(i, j, k, n)=zeroval; + } } } - + if(dsdt(i,j,k,n)!=dsdt(i,j,k,n)) { amrex::Abort("NaN detected in dsdt\n"); diff --git a/Source/Src_3d/ib_utils.H b/Source/Src_3d/ib_utils.H index feb3b92..604d9fb 100644 --- a/Source/Src_3d/ib_utils.H +++ b/Source/Src_3d/ib_utils.H @@ -273,7 +273,8 @@ update_cutcells(int i, int j, int k, Array4 const& phi, //state variable const GpuArray& dx, const GpuArray& plo, - Real time, int spec_in_solid, int use_bg_gas) + Real time, int spec_in_solid, int use_bg_gas, + int chtflag) { Real lsgrad[3]={0.0}; @@ -334,8 +335,34 @@ update_cutcells(int i, int j, int k, phi(i,j,k,PRES_INDX)=do_interpolation_hom_neumann(i,j,k,STENCILGROW,xp,phi,PRES_INDX,dx); - //adiabatic wall for now - phi(i,j,k,TEMP_INDX)=do_interpolation_hom_neumann(i,j,k,STENCILGROW,xp,phi,TEMP_INDX,dx); + if(!chtflag) + { + //zero gradient + phi(i,j,k,TEMP_INDX)=do_interpolation_hom_neumann(i,j,k,STENCILGROW,xp,phi,TEMP_INDX,dx); + } + else + { + //3 by 3 stencil + Real sumtemp=0.0; + int cnt=0; + for(int n=-1;n<=1;n++) + { + for(int m=-1;m<=1;m++) + { + for(int l=-1;l<=1;l++) + { + if(!(l==0 && m==0 && n==0)) + { + sumtemp+=phi(i+l,j+m,k+n,TEMP_INDX); + cnt++; + } + } + } + } + //cnt should be 26 + phi(i,j,k,TEMP_INDX)=sumtemp/Real(cnt); + //amrex::Print()<<"sumtemp,cnt,temp:"< -void mflo::Evolve_split(Real t_ss,Real t_react,Real dt_react,Real dt_react_rk,Real dt_ss) -{ - Real cur_time = t_new[0]; - bool only_flow=true; - Real react_time=0.0; - Real flow_time; - Real tchem; - - amrex::Print()<<"running flow to steady-state\n"; - EvolveAMR(t_ss,only_flow); - flow_time=t_ss; - - while(react_time < t_react) - { - amrex::Print()<<"advance chemistry\n"; - //temporarily store inert gas concentration in dens index - if(using_bg_inertgas) - { - for (int l = 0; l <= finest_level; l++) - { - store_inertgas_conc(l); - } - - } - for (int l = 0; l <= finest_level; l++) - { - tchem=0.0; - while(tchem < dt_react) - { - Advance_chemistry(l, react_time+tchem, dt_react_rk); - tchem += dt_react_rk; - } - } - AverageDown(); - - for (int l = 0; l <= finest_level; l++) - { - //do update of flow variables - update_vars_after_chemsolve(l); - //update_primitive_vars(l); - } - - amrex::Print()<<"advance flow\n"; - EvolveAMR(flow_time+dt_ss,only_flow); - - react_time += dt_react; - flow_time += dt_ss; - } - -} - // advance solution to final time void mflo::EvolveAMR(Real final_time, bool only_flow) { @@ -74,6 +23,7 @@ void mflo::EvolveAMR(Real final_time, bool only_flow) { amrex::Print() << "\nCoarse STEP " << step + 1 << " starts ..." << std::endl; + ComputeDt(); int lev = 0; @@ -161,11 +111,11 @@ void mflo::compute_residual_norms() } } -void mflo::update_vars_after_chemsolve(int lev) +void mflo::update_primitive_vars(int lev) { MultiFab& S_new = phi_new[lev]; - int using_inert_gas=using_bg_inertgas; bool nsflag=(do_ns==1)?true:false; + int chtflag=conj_ht; #ifdef _OPENMP #pragma omp parallel if (Gpu::notInLaunchRegion()) @@ -178,28 +128,10 @@ void mflo::update_vars_after_chemsolve(int lev) amrex::ParallelFor(bx, [=] AMREX_GPU_DEVICE(int i, int j, int k) { - if(snew_arr(i,j,k,VFRAC_INDX) >= one) + //cant be greater than 1 + if(snew_arr(i,j,k,VFRAC_INDX) == one) { - Real u[NCVARS+NUM_SPECIES], p[NCVARS],spec[NUM_SPECIES]; - for (int c = 0; c < NUM_SPECIES; c++) - { - spec[c] = snew_arr(i, j, k, c + FLO_NVARS); - } - - //also pass background gas concentration - if(using_inert_gas) - { - snew_arr(i,j,k,RHO_INDX) = mflo_thermo::get_r_from_c(spec,snew_arr(i,j,k,DENS_INDX)); - snew_arr(i,j,k,DENS_INDX) = snew_arr(i,j,k,RHO_INDX); - } - else - { - snew_arr(i,j,k,RHO_INDX) = mflo_thermo::get_r_from_c(spec); - snew_arr(i,j,k,DENS_INDX) = snew_arr(i,j,k,RHO_INDX); - } - snew_arr(i,j,k,RHOU_INDX) = snew_arr(i,j,k,RHO_INDX)*snew_arr(i,j,k,VELX_INDX); - snew_arr(i,j,k,RHOV_INDX) = snew_arr(i,j,k,RHO_INDX)*snew_arr(i,j,k,VELY_INDX); - snew_arr(i,j,k,RHOW_INDX) = snew_arr(i,j,k,RHO_INDX)*snew_arr(i,j,k,VELZ_INDX); + Real u[NCVARS+NUM_SPECIES], p[NCVARS]; for (int c = 0; c < NCVARS; c++) { @@ -217,85 +149,13 @@ void mflo::update_vars_after_chemsolve(int lev) snew_arr(i, j, k, TEMP_INDX)=mflo_thermo::get_t_from_rpc(snew_arr(i, j, k, DENS_INDX), snew_arr(i, j, k, PRES_INDX),u+NCVARS); } - }); - - Real minpressure=S_new[mfi].min(PRES_INDX); - if(minpressure < 0 and nsflag) - { - Print() << "Minimum pressure in domain:" << minpressure << "\n"; - amrex::Abort("Pressure has gone negative"); - } - } - } -} - -void mflo::store_inertgas_conc(int lev) -{ - MultiFab& S_new = phi_new[lev]; - -#ifdef _OPENMP -#pragma omp parallel if (Gpu::notInLaunchRegion()) -#endif - { - for (MFIter mfi(S_new, TilingIfNotGPU()); mfi.isValid(); ++mfi) - { - const Box& bx = mfi.tilebox(); - Array4 snew_arr = S_new.array(mfi); - - amrex::ParallelFor(bx, [=] AMREX_GPU_DEVICE(int i, int j, int k) { - - if(snew_arr(i,j,k,VFRAC_INDX) >= one) + else { - Real spec[NUM_SPECIES]; - for (int c = 0; c < NUM_SPECIES; c++) + //there is conjugate heat transfer and I am in the solid.. + if(chtflag && snew_arr(i,j,k,VFRAC_INDX)==zeroval) { - spec[c] = snew_arr(i, j, k, c + FLO_NVARS); + snew_arr(i,j,k,TEMP_INDX)=mflo_thermo::get_solid_t_from_rhoe(snew_arr(i,j,k,RHOE_INDX)); } - snew_arr(i, j, k, DENS_INDX)= - mflo_thermo::get_bgasconc_from_rc(snew_arr(i,j,k,RHO_INDX),spec); - } - }); - } - } -} - - - -void mflo::update_primitive_vars(int lev) -{ - MultiFab& S_new = phi_new[lev]; - bool nsflag=(do_ns==1)?true:false; - -#ifdef _OPENMP -#pragma omp parallel if (Gpu::notInLaunchRegion()) -#endif - { - for (MFIter mfi(S_new, TilingIfNotGPU()); mfi.isValid(); ++mfi) - { - const Box& bx = mfi.tilebox(); - Array4 snew_arr = S_new.array(mfi); - - amrex::ParallelFor(bx, [=] AMREX_GPU_DEVICE(int i, int j, int k) { - - if(snew_arr(i,j,k,VFRAC_INDX) >= one) - { - Real u[NCVARS+NUM_SPECIES], p[NCVARS]; - - for (int c = 0; c < NCVARS; c++) - { - u[c + RHO_IND] = snew_arr(i, j, k, c + RHO_INDX); - } - for (int c = 0; c < NUM_SPECIES; c++) - { - u[c + NCVARS] = snew_arr(i, j, k, c + FLO_NVARS); - } - cons_to_prim(u, p); - for (int c = 0; c < NCVARS; c++) - { - snew_arr(i, j, k, c + DENS_INDX) = p[c + DENS_IND]; - } - snew_arr(i, j, k, TEMP_INDX)=mflo_thermo::get_t_from_rpc(snew_arr(i, j, k, DENS_INDX), - snew_arr(i, j, k, PRES_INDX),u+NCVARS); } }); @@ -307,6 +167,7 @@ void mflo::update_primitive_vars(int lev) } } } + } // advance a level by dt @@ -346,6 +207,7 @@ void mflo::timeStep(int lev, Real time, int iteration, bool only_flow) } } } + if (Verbose()) { @@ -590,33 +452,6 @@ void mflo::Advance_coupled(int lev, Real time, Real dt_lev, int iteration, int n MultiFab::Saxpy(S_new, dt_lev, dsdt_chemistry, 0, 0, S_new.nComp(), 0); }*/ } - - -void mflo::Advance_chemistry(int lev, Real time, Real dt_lev) -{ - constexpr int num_grow = 3; - std::swap(phi_old[lev], phi_new[lev]); // old becomes new and new becomes old - MultiFab& S_new = phi_new[lev]; // this is the old value, beware! - MultiFab& S_old = phi_old[lev]; // current value - - // source term - MultiFab dsdt_chemistry(grids[lev], dmap[lev], S_new.nComp(), 0); - - - // stage 1 - // compute dsdt for 1/2 timestep - compute_dsdt_chemistry(lev, num_grow, S_new, dsdt_chemistry, time); - // S_new=S_old+0.5*dt*dsdt //sold is the current value - MultiFab::LinComb(S_new, one, S_old, 0, half * dt_lev, dsdt_chemistry, 0, 0, S_new.nComp(), 0); - - // stage 2 - // dsdt for full time-step - compute_dsdt_chemistry(lev, num_grow, S_new, dsdt_chemistry, time + half * dt_lev); - // S_new=S_old+dt*dsdt - MultiFab::LinComb(S_new, one, S_old, 0, dt_lev, dsdt_chemistry, 0, 0, S_new.nComp(), 0); - -} - /* Advances the chemisty state from time -> time + dt_lev This is done using TimeIntegrator from AMReX, which can @@ -653,45 +488,8 @@ void mflo::Advance_chemistry_implicit(int lev, Real time, Real dt_lev) TimeIntegrator> integrator(state_old); integrator.set_rhs(rhs_function); // Advance from time to time + dt_lev - integrator.advance(state_old, state_new, time, dt_lev); //S_new/phi_new should have the new state -} - -void mflo::compute_dsdt_chemistry( - int lev, - const int num_grow, - MultiFab& S, - MultiFab& dsdt, - Real time) -{ - dsdt.setVal(zeroval); - -#ifdef _OPENMP -#pragma omp parallel if (Gpu::notInLaunchRegion()) -#endif - { - for (MFIter mfi(dsdt, TilingIfNotGPU()); mfi.isValid(); ++mfi) - { - const Box& bx = mfi.tilebox(); - - FArrayBox source_fab(bx,TOTAL_NVARS); //external sources - source_fab.setVal(zeroval); - Elixir source_fab_eli = source_fab.elixir(); - - Array4 s_arr = S.array(mfi); - Array4 dsdt_arr = dsdt.array(mfi); - Array4 source_arr = source_fab.array(); - - auto prob_lo = geom[lev].ProbLoArray(); - const auto dx = geom[lev].CellSizeArray(); - - amrex::ParallelFor(bx, [=] AMREX_GPU_DEVICE(int i, int j, int k) { - mflo_chem_reactions::compute_spec_source( - i, j, k, s_arr, source_arr, prob_lo, dx, time); - }); - - dsdt[mfi].plus(source_fab); - } - } + //S_new/phi_new should have the new state + integrator.advance(state_old, state_new, time, dt_lev); } void mflo::update_cutcell_data( @@ -706,6 +504,7 @@ void mflo::update_cutcell_data( const auto dx = geom[lev].CellSizeArray(); int using_inert_gas=using_bg_inertgas; int spec_in_solid=species_in_solid; + int chtflag=conj_ht; for (MFIter mfi(Sborder, TilingIfNotGPU()); mfi.isValid(); ++mfi) { @@ -714,7 +513,7 @@ void mflo::update_cutcell_data( amrex::ParallelFor(bx, [=] AMREX_GPU_DEVICE(int i, int j, int k) { update_cutcells(i,j,k,sborder_arr,dx,prob_lo,time, - spec_in_solid,using_inert_gas); + spec_in_solid,using_inert_gas,chtflag); }); } } @@ -755,6 +554,7 @@ void mflo::compute_dsdt_flow( const auto dx = geom[lev].CellSizeArray(); const auto prob_lo = geom[lev].ProbLoArray(); bool nsflag=(do_ns==1)?true:false; + int conjhtflag=conj_ht; int ncomp = Sborder.nComp(); int hyperbolics_order = order_hyp; @@ -847,7 +647,7 @@ void mflo::compute_dsdt_flow( update_residual( i, j, k, n, dsdt_arr, source_arr, sborder_arr, AMREX_D_DECL(flux_arr[0], flux_arr[1], flux_arr[2]), - dx,spec_in_solid); + dx,spec_in_solid,conjhtflag); }); } } @@ -879,5 +679,42 @@ void mflo::compute_dsdt_flow( } } } +} +void mflo::compute_dsdt_chemistry( + int lev, + const int num_grow, + MultiFab& S, + MultiFab& dsdt, + Real time) +{ + dsdt.setVal(zeroval); + +#ifdef _OPENMP +#pragma omp parallel if (Gpu::notInLaunchRegion()) +#endif + { + for (MFIter mfi(dsdt, TilingIfNotGPU()); mfi.isValid(); ++mfi) + { + const Box& bx = mfi.tilebox(); + + FArrayBox source_fab(bx,TOTAL_NVARS); //external sources + source_fab.setVal(zeroval); + Elixir source_fab_eli = source_fab.elixir(); + + Array4 s_arr = S.array(mfi); + Array4 dsdt_arr = dsdt.array(mfi); + Array4 source_arr = source_fab.array(); + + auto prob_lo = geom[lev].ProbLoArray(); + const auto dx = geom[lev].CellSizeArray(); + + amrex::ParallelFor(bx, [=] AMREX_GPU_DEVICE(int i, int j, int k) { + mflo_chem_reactions::compute_spec_source( + i, j, k, s_arr, source_arr, prob_lo, dx, time); + }); + + dsdt[mfi].plus(source_fab); + } + } } diff --git a/Source/mflo_splitchem.cpp b/Source/mflo_splitchem.cpp new file mode 100644 index 0000000..e22f82d --- /dev/null +++ b/Source/mflo_splitchem.cpp @@ -0,0 +1,195 @@ +#include +#include +#include +#include +#include +#include +#include +#include +#include + +#ifdef AMREX_MEM_PROFILING +#include +#endif + +#include + +//Note: This split evolve is used when chemistry is happening +//at a very slow timescale compared to flow. +//the strategy is to solve the flow to steady-state, advance the +//chemistry, and then solve the flow to steady state again. +//the technique is better described in +//Sitaraman et al., Chem Engg Res. and Design, 197, 2023 +//Sitaraman et al., Chemical Engineering Science 206 (2019): 348-360 +void mflo::Evolve_split(Real t_ss,Real t_react,Real dt_react,Real dt_react_rk,Real dt_ss) +{ + Real cur_time = t_new[0]; + bool only_flow=true; + Real react_time=0.0; + Real flow_time; + Real tchem; + + amrex::Print()<<"running flow to steady-state\n"; + EvolveAMR(t_ss,only_flow); + flow_time=t_ss; + + while(react_time < t_react) + { + amrex::Print()<<"advance chemistry\n"; + //temporarily store inert gas concentration in dens index + if(using_bg_inertgas) + { + for (int l = 0; l <= finest_level; l++) + { + store_inertgas_conc(l); + } + } + for (int l = 0; l <= finest_level; l++) + { + tchem=0.0; + while(tchem < dt_react) + { + Advance_chemistry(l, react_time+tchem, dt_react_rk); + tchem += dt_react_rk; + } + } + AverageDown(); + + for (int l = 0; l <= finest_level; l++) + { + //do update of flow variables + update_vars_after_chemsolve(l); + //update_primitive_vars(l); + } + + amrex::Print()<<"advance flow\n"; + EvolveAMR(flow_time+dt_ss,only_flow); + + react_time += dt_react; + flow_time += dt_ss; + } + +} + +void mflo::update_vars_after_chemsolve(int lev) +{ + MultiFab& S_new = phi_new[lev]; + int using_inert_gas=using_bg_inertgas; + bool nsflag=(do_ns==1)?true:false; + +#ifdef _OPENMP +#pragma omp parallel if (Gpu::notInLaunchRegion()) +#endif + { + for (MFIter mfi(S_new, TilingIfNotGPU()); mfi.isValid(); ++mfi) + { + const Box& bx = mfi.tilebox(); + Array4 snew_arr = S_new.array(mfi); + + amrex::ParallelFor(bx, [=] AMREX_GPU_DEVICE(int i, int j, int k) { + + if(snew_arr(i,j,k,VFRAC_INDX) >= one) + { + Real u[NCVARS+NUM_SPECIES], p[NCVARS],spec[NUM_SPECIES]; + for (int c = 0; c < NUM_SPECIES; c++) + { + spec[c] = snew_arr(i, j, k, c + FLO_NVARS); + } + + //also pass background gas concentration + if(using_inert_gas) + { + snew_arr(i,j,k,RHO_INDX) = mflo_thermo::get_r_from_c(spec,snew_arr(i,j,k,DENS_INDX)); + snew_arr(i,j,k,DENS_INDX) = snew_arr(i,j,k,RHO_INDX); + } + else + { + snew_arr(i,j,k,RHO_INDX) = mflo_thermo::get_r_from_c(spec); + snew_arr(i,j,k,DENS_INDX) = snew_arr(i,j,k,RHO_INDX); + } + snew_arr(i,j,k,RHOU_INDX) = snew_arr(i,j,k,RHO_INDX)*snew_arr(i,j,k,VELX_INDX); + snew_arr(i,j,k,RHOV_INDX) = snew_arr(i,j,k,RHO_INDX)*snew_arr(i,j,k,VELY_INDX); + snew_arr(i,j,k,RHOW_INDX) = snew_arr(i,j,k,RHO_INDX)*snew_arr(i,j,k,VELZ_INDX); + + for (int c = 0; c < NCVARS; c++) + { + u[c + RHO_IND] = snew_arr(i, j, k, c + RHO_INDX); + } + for (int c = 0; c < NUM_SPECIES; c++) + { + u[c + NCVARS] = snew_arr(i, j, k, c + FLO_NVARS); + } + cons_to_prim(u, p); + for (int c = 0; c < NCVARS; c++) + { + snew_arr(i, j, k, c + DENS_INDX) = p[c + DENS_IND]; + } + snew_arr(i, j, k, TEMP_INDX)=mflo_thermo::get_t_from_rpc(snew_arr(i, j, k, DENS_INDX), + snew_arr(i, j, k, PRES_INDX),u+NCVARS); + } + }); + + Real minpressure=S_new[mfi].min(PRES_INDX); + if(minpressure < 0 and nsflag) + { + Print() << "Minimum pressure in domain:" << minpressure << "\n"; + amrex::Abort("Pressure has gone negative"); + } + } + } +} + +void mflo::store_inertgas_conc(int lev) +{ + MultiFab& S_new = phi_new[lev]; + +#ifdef _OPENMP +#pragma omp parallel if (Gpu::notInLaunchRegion()) +#endif + { + for (MFIter mfi(S_new, TilingIfNotGPU()); mfi.isValid(); ++mfi) + { + const Box& bx = mfi.tilebox(); + Array4 snew_arr = S_new.array(mfi); + + amrex::ParallelFor(bx, [=] AMREX_GPU_DEVICE(int i, int j, int k) { + + if(snew_arr(i,j,k,VFRAC_INDX) >= one) + { + Real spec[NUM_SPECIES]; + for (int c = 0; c < NUM_SPECIES; c++) + { + spec[c] = snew_arr(i, j, k, c + FLO_NVARS); + } + snew_arr(i, j, k, DENS_INDX)= + mflo_thermo::get_bgasconc_from_rc(snew_arr(i,j,k,RHO_INDX),spec); + } + }); + } + } +} + +void mflo::Advance_chemistry(int lev, Real time, Real dt_lev) +{ + constexpr int num_grow = 3; + std::swap(phi_old[lev], phi_new[lev]); // old becomes new and new becomes old + MultiFab& S_new = phi_new[lev]; // this is the old value, beware! + MultiFab& S_old = phi_old[lev]; // current value + + // source term + MultiFab dsdt_chemistry(grids[lev], dmap[lev], S_new.nComp(), 0); + + + // stage 1 + // compute dsdt for 1/2 timestep + compute_dsdt_chemistry(lev, num_grow, S_new, dsdt_chemistry, time); + // S_new=S_old+0.5*dt*dsdt //sold is the current value + MultiFab::LinComb(S_new, one, S_old, 0, half * dt_lev, dsdt_chemistry, 0, 0, S_new.nComp(), 0); + + // stage 2 + // dsdt for full time-step + compute_dsdt_chemistry(lev, num_grow, S_new, dsdt_chemistry, time + half * dt_lev); + // S_new=S_old+dt*dsdt + MultiFab::LinComb(S_new, one, S_old, 0, dt_lev, dsdt_chemistry, 0, 0, S_new.nComp(), 0); + +} diff --git a/models/blockCatalyst1d/thermo.H b/models/blockCatalyst1d/thermo.H index 1cacb40..4bc27d6 100644 --- a/models/blockCatalyst1d/thermo.H +++ b/models/blockCatalyst1d/thermo.H @@ -15,6 +15,24 @@ using namespace amrex; namespace mflo_thermo { + AMREX_GPU_HOST_DEVICE AMREX_INLINE Real + get_solid_t_from_rhoe(Real rhoe) + { + //modify if conjugate heat transfer + //included + amrex::Real temp=300.0; + return(temp); + } + + AMREX_GPU_HOST_DEVICE AMREX_INLINE Real + get_solid_rhoe_from_t(Real temp) + { + //modify if conjugate + //heat transfer included + amrex::Real Cv_solid=1.0; + amrex::Real rho_solid=1.0; + return(rho_solid*Cv_solid*temp); + } AMREX_GPU_HOST_DEVICE AMREX_INLINE Real get_r_from_c(Real spec[NUM_SPECIES],Real bgasconc=zeroval) { diff --git a/models/catalystCylinder/thermo.H b/models/catalystCylinder/thermo.H index a23103b..9f87490 100644 --- a/models/catalystCylinder/thermo.H +++ b/models/catalystCylinder/thermo.H @@ -15,6 +15,24 @@ using namespace amrex; namespace mflo_thermo { + AMREX_GPU_HOST_DEVICE AMREX_INLINE Real + get_solid_t_from_rhoe(Real rhoe) + { + //modify if conjugate heat transfer + //included + amrex::Real temp=300.0; + return(temp); + } + + AMREX_GPU_HOST_DEVICE AMREX_INLINE Real + get_solid_rhoe_from_t(Real temp) + { + //modify if conjugate + //heat transfer included + amrex::Real Cv_solid=1.0; + amrex::Real rho_solid=1.0; + return(rho_solid*Cv_solid*temp); + } AMREX_GPU_HOST_DEVICE AMREX_INLINE Real get_r_from_c(Real spec[NUM_SPECIES],Real bgasconc=zeroval) { diff --git a/models/conjHtransferTest/GNUmakefile b/models/conjHtransferTest/GNUmakefile new file mode 100644 index 0000000..c55b6aa --- /dev/null +++ b/models/conjHtransferTest/GNUmakefile @@ -0,0 +1,37 @@ +PRECISION = DOUBLE +PROFILE = FALSE + +DEBUG = TRUE +DEBUG = FALSE + +DIM = 3 + +COMP = gnu + +USE_MPI = TRUE +USE_OMP = FALSE +USE_CUDA = FALSE + +Bpack := ./Make.package +Blocs := . + +include ../Make.mflo + +ifeq ($(USE_SUNDIALS),TRUE) +# NOTE: SUNDIALS_ROOT must point to the directory where sundials is installed +# A good check is to see if $(SUNDIALS_ROOT)/lib has a bunch of libsundials_ files +# To run with sundials, enabled, please compile with USE_SUNDIALS = TRUE +SUNDIALS_ROOT ?= $(TOP)/../../../../../sundials/instdir +SUNDIALS_LIB_DIR ?= $(SUNDIALS_ROOT)/lib + +USE_CVODE_LIBS ?= TRUE +USE_ARKODE_LIBS ?= TRUE + +DEFINES += -DAMREX_USE_SUNDIALS +INCLUDE_LOCATIONS += $(SUNDIALS_ROOT)/include +LIBRARY_LOCATIONS += $(SUNDIALS_LIB_DIR) + +LIBRARIES += -L$(SUNDIALS_LIB_DIR) -lsundials_cvode +LIBRARIES += -L$(SUNDIALS_LIB_DIR) -lsundials_arkode +LIBRARIES += -L$(SUNDIALS_LIB_DIR) -lsundials_nvecmanyvector +endif diff --git a/models/conjHtransferTest/Make.package b/models/conjHtransferTest/Make.package new file mode 100644 index 0000000..aaf4145 --- /dev/null +++ b/models/conjHtransferTest/Make.package @@ -0,0 +1,5 @@ +CEXE_headers += prob.H userfuncs.H externbc.H species.H transport.H chemistry.H +CEXE_sources += userfuncs.cpp species.cpp + +#f90EXE_sources += face_velocity_$(DIM)d.f90 +#f90EXE_source += Prob.f90 diff --git a/models/conjHtransferTest/Prob.H b/models/conjHtransferTest/Prob.H new file mode 100644 index 0000000..e9bb216 --- /dev/null +++ b/models/conjHtransferTest/Prob.H @@ -0,0 +1,124 @@ +#ifndef _PROB_H_ +#define _PROB_H_ + +#include +#include +#include +#include +#include +#include + +using namespace amrex; + +AMREX_GPU_HOST_DEVICE AMREX_INLINE +void initdata(Box const& bx, Array4 const& phi, GeometryData const& geomdata) +{ + const auto lo = lbound(bx); + const auto hi = ubound(bx); + + const Real* AMREX_RESTRICT prob_lo = geomdata.ProbLo(); + const Real* AMREX_RESTRICT prob_hi = geomdata.ProbHi(); + const Real* AMREX_RESTRICT dx = geomdata.CellSize(); + + Real xlen,ylen,zlen,maxlen; + Real spec[NUM_SPECIES]={zeroval}; + int dir=0; + + xlen=prob_hi[0]-prob_lo[0]; + ylen=prob_hi[1]-prob_lo[1]; + zlen=prob_hi[2]-prob_lo[2]; + + maxlen=std::max(xlen,std::max(ylen,zlen)); + + if(maxlen==xlen) dir=0; + if(maxlen==ylen) dir=1; + if(maxlen==zlen) dir=2; + + Real xyz[3]={zeroval}; + +#ifdef _OPENMP +#pragma omp parallel for collapse(2) if (GPU::notInLaunchRegion) +#endif + for (int k = lo.z; k <= hi.z; ++k) + { + for (int j = lo.y; j <= hi.y; ++j) + { + xyz[2] = prob_lo[2] + (0.5+k) * dx[2]; + xyz[1] = prob_lo[1] + (0.5+j) * dx[1]; + + AMREX_PRAGMA_SIMD + for (int i = lo.x; i <= hi.x; ++i) + { + xyz[0] = prob_lo[0] + (0.5+i) * dx[0]; + + + Real vfrac=0.0; + for(int kk=0;kk<2;kk++) + { + for(int jj=0;jj<2;jj++) + { + for(int ii=0;ii<2;ii++) + { + Real xx[3]={prob_lo[0]+(i+ii)*dx[0], + prob_lo[1]+(j+jj)*dx[1], + prob_lo[2]+(k+kk)*dx[2]}; + + if(xx[dir]zeroval) + { + phi(i,j,k,TEMP_INDX) = mflo_thermo::get_t_from_rpc(phi(i,j,k,DENS_INDX), + phi(i,j,k,PRES_INDX),spec); + } + else + { + phi(i,j,k,TEMP_INDX)=300.0; + } + + phi(i,j,k,RHO_INDX) = phi(i,j,k,DENS_INDX); + phi(i,j,k,RHOU_INDX) = phi(i,j,k,DENS_INDX)*phi(i,j,k,VELX_INDX); + phi(i,j,k,RHOV_INDX) = phi(i,j,k,DENS_INDX)*phi(i,j,k,VELY_INDX); + phi(i,j,k,RHOW_INDX) = phi(i,j,k,DENS_INDX)*phi(i,j,k,VELZ_INDX); + + if(vfrac>zeroval) + { + phi(i,j,k,RHOE_INDX) = phi(i,j,k,DENS_INDX)* + mflo_thermo::get_e_from_rpc(phi(i,j,k,DENS_INDX), + phi(i,j,k,PRES_INDX),spec); + phi(i,j,k,RHOE_INDX) += half*phi(i,j,k,DENS_INDX) * + ( phi(i,j,k,VELX_INDX)*phi(i,j,k,VELX_INDX) + + phi(i,j,k,VELY_INDX)*phi(i,j,k,VELY_INDX) + + phi(i,j,k,VELZ_INDX)*phi(i,j,k,VELZ_INDX) ); + } + else + { + phi(i,j,k,RHOE_INDX) = mflo_thermo::get_solid_rhoe_from_t(phi(i,j,k,TEMP_INDX)); + } + + for(int sp=0;sp +#include +#include +#include +#include + +using namespace amrex; +namespace mflo_chem_reactions +{ + AMREX_GPU_DEVICE AMREX_INLINE + void compute_spec_source(int i, int j, int k, + Array4 const& phi, + Array4 const& specsource, + GpuArray prob_lo, + GpuArray dx, + const Real time) + { + specsource(i,j,k,FLO_NVARS+AIR_ID) = zeroval; + specsource(i,j,k,FLO_NVARS+SOLSPEC_ID) = zeroval; + } +} +#endif diff --git a/models/conjHtransferTest/externbc.H b/models/conjHtransferTest/externbc.H new file mode 100644 index 0000000..f06e7b0 --- /dev/null +++ b/models/conjHtransferTest/externbc.H @@ -0,0 +1,54 @@ +#ifndef _EXTERNBC_H_ +#define _EXTERNBC_H_ + +#include +#include +#include +#include +#include +#include +#include + +using namespace amrex; + +AMREX_GPU_DEVICE +AMREX_FORCE_INLINE +void externalbc(const amrex::Real x[AMREX_SPACEDIM], + const amrex::Real s_int[], + amrex::Real s_ext[], + const int idir, + const int sgn, + const amrex::Real time, + amrex::GeometryData const& geomdata) +{ + //default adiabatic wall + s_ext[PRES_INDX] = s_int[PRES_INDX]; + s_ext[DENS_INDX] = s_int[DENS_INDX]; + s_ext[TEMP_INDX] = s_int[TEMP_INDX]; + s_ext[VELX_INDX] = -s_int[VELX_INDX]; + s_ext[VELY_INDX] = -s_int[VELY_INDX]; + s_ext[VELZ_INDX] = -s_int[VELZ_INDX]; + + s_ext[RHO_INDX] = s_int[RHO_INDX]; + s_ext[RHOU_INDX] = -s_int[RHOU_INDX]; + s_ext[RHOV_INDX] = -s_int[RHOV_INDX]; + s_ext[RHOW_INDX] = -s_int[RHOW_INDX]; + s_ext[RHOE_INDX] = s_int[RHOE_INDX]; + s_ext[VFRAC_INDX] = s_int[VFRAC_INDX]; + + for(int sp=0;sp 0) + { + s_ext[TEMP_INDX]=1100.0; + } + else + { + s_ext[TEMP_INDX]=300.0; + } + +} +#endif diff --git a/models/conjHtransferTest/inputs_x b/models/conjHtransferTest/inputs_x new file mode 100644 index 0000000..90580a2 --- /dev/null +++ b/models/conjHtransferTest/inputs_x @@ -0,0 +1,59 @@ +max_step = 300000 +stop_time = 1.35 + +# PROBLEM SIZE & GEOMETRY +geometry.is_periodic = 0 0 0 +geometry.coord_sys = 0 # 0 => cart +geometry.prob_lo = 0.0 0.0 0.0 +geometry.prob_hi = 0.1 0.0125 0.0125 +amr.n_cell = 64 4 4 + +# VERBOSITY +amr.v = 1 # verbosity in Amr + +# REFINEMENT +amr.max_level = 1 # maximum level number allowed +amr.blocking_factor = 4 # block factor in grid generation +amr.max_grid_size = 16 + +amr.regrid_int = 2 # how often to regrid + +# TIME STEP CONTROL +mflo.cfl = 0.2 # cfl number for hyperbolic system + # In this test problem, the velocity is + # time-dependent. We could use 0.9 in + # the 3D test, but need to use 0.7 in 2D + # to satisfy CFL condition. +mflo.order_hyp=1 +mflo.do_reflux = 1 +mflo.solve_navier_stokes=1 + +#periodic 0 +#extdir 3 +#foextrap 2 +#wall 6 + +mflo.lo_bc = 3 6 6 +mflo.hi_bc = 3 6 6 + +mflo.conjugate_heat_transfer=1 + +# Tagging +mflo.tagged_vars = volfrac +mflo.volfrac_refine = 1e20 +mflo.volfrac_refinegrad = 0.1 + +# PLOTFILES +amr.plot_file = plt # root name of plot file +amr.plot_int = 10000 # number of timesteps between plot files + +# CHECKPOINT +amr.chk_file = chk # root name of checkpoint file +amr.chk_int = -1 # number of timesteps between checkpoint files + +# INTEGRATION options (for chemistry) +## Uncomment below two lines to enable implicit chemistry +# integration.type = SUNDIALS +# integration.sundials.strategy = CVODE +## Note that this requires a working sundials installation and +## the code must be compiled with USE_SUNDIALS=TRUE, see GNUmakefile for additional info diff --git a/models/conjHtransferTest/species.H b/models/conjHtransferTest/species.H new file mode 100644 index 0000000..fba2a8d --- /dev/null +++ b/models/conjHtransferTest/species.H @@ -0,0 +1,26 @@ +#ifndef _SPECIES_H_ +#define _SPECIES_H_ + +#include +#include +#include +#include +#include +#include +#include + +#define NUM_SPECIES 2 +#define AIR_ID 0 +#define NUM_GAS_SPECIES 1 +#define SOLSPEC_ID 1 +namespace mflo_species +{ + extern amrex::Vector specnames; + extern AMREX_GPU_DEVICE_MANAGED amrex::Real molwts[NUM_SPECIES]; + extern AMREX_GPU_DEVICE_MANAGED amrex::Real advect_flags[NUM_SPECIES]; + void init(); + void close(); + int find_id(std::string specname); +} + +#endif diff --git a/models/conjHtransferTest/species.cpp b/models/conjHtransferTest/species.cpp new file mode 100644 index 0000000..744332a --- /dev/null +++ b/models/conjHtransferTest/species.cpp @@ -0,0 +1,40 @@ +#include +#include + +namespace mflo_species +{ + amrex::Vector specnames(NUM_SPECIES); + AMREX_GPU_DEVICE_MANAGED amrex::Real molwts[NUM_SPECIES]={zeroval}; + AMREX_GPU_DEVICE_MANAGED amrex::Real advect_flags[NUM_SPECIES]={one}; + + void init() + { + specnames[AIR_ID]="air"; + specnames[SOLSPEC_ID]="solidmat"; + + for(int sp=0;sp +#include +#include +#include +#include +#include +#include +#include +#include +#include + +using namespace amrex; +namespace mflo_thermo +{ + AMREX_GPU_HOST_DEVICE AMREX_INLINE Real + get_solid_t_from_rhoe(Real rhoe) + { + amrex::Real Cv_solid=100.0; + amrex::Real rho_solid=2000.0; + amrex::Real temp=rhoe/(rho_solid*Cv_solid); + //amrex::Print()<<"rhoe,temp:"< +#include +#include +#include + +using namespace amrex; +namespace mflo_chem_transport +{ + AMREX_GPU_DEVICE AMREX_INLINE + void compute_spec_dcoeff(int i, int j, int k, + Array4 const& phi, + Array4 const& dcoeff, + GpuArray prob_lo, + GpuArray dx, + const Real time) + { + dcoeff(i,j,k,AIR_ID)=1e-6; + dcoeff(i,j,k,SOLSPEC_ID)=zeroval; + } +} +#endif diff --git a/models/conjHtransferTest/userfuncs.H b/models/conjHtransferTest/userfuncs.H new file mode 100644 index 0000000..518731f --- /dev/null +++ b/models/conjHtransferTest/userfuncs.H @@ -0,0 +1,48 @@ +#ifndef _USERFUNCS_H_ +#define _USERFUNCS_H_ + +#include +#include +#include +#include +#include +#include +#include +#include +#include + +using namespace amrex; +namespace mflo_user_funcs +{ + extern AMREX_GPU_DEVICE_MANAGED Real iboffset; + void initialize_problem(); + + AMREX_GPU_DEVICE AMREX_INLINE + void compute_fluid_transport(int i, int j, int k, + Array4 const& phi, + Array4 const& transpcoeffs, + GpuArray prob_lo, + GpuArray dx, + const Real time) + { + Real volfrac=phi(i,j,k,VFRAC_INDX); + transpcoeffs(i,j,k,VISC_INDX ) = 0.01; //viscosity + transpcoeffs(i,j,k,THCOND_INDX) = 50.0*volfrac+500.0*(1.0-volfrac); //thermal conductivity + } + + AMREX_GPU_DEVICE AMREX_INLINE + void compute_fluid_source(int i, int j, int k, + Array4 const& phi, + Array4 const& source, + GpuArray prob_lo, + GpuArray dx, + const Real time) + { + for(int nc=0;nc +#include + +namespace mflo_user_funcs +{ + AMREX_GPU_DEVICE_MANAGED Real iboffset=0.7; + void initialize_problem() + { + ParmParse pp("user"); + pp.query("iboffset",iboffset); + + Print()<<"Initializing problem\n"; + + //This is a good place to have some global + //parameters defined for initializing solution + //vector + } +} diff --git a/models/porousParticle/thermo.H b/models/porousParticle/thermo.H index 2208328..37447c1 100644 --- a/models/porousParticle/thermo.H +++ b/models/porousParticle/thermo.H @@ -15,6 +15,24 @@ using namespace amrex; namespace mflo_thermo { + AMREX_GPU_HOST_DEVICE AMREX_INLINE Real + get_solid_t_from_rhoe(Real rhoe) + { + //modify if conjugate heat transfer + //included + amrex::Real temp=300.0; + return(temp); + } + + AMREX_GPU_HOST_DEVICE AMREX_INLINE Real + get_solid_rhoe_from_t(Real temp) + { + //modify if conjugate + //heat transfer included + amrex::Real Cv_solid=1.0; + amrex::Real rho_solid=1.0; + return(rho_solid*Cv_solid*temp); + } AMREX_GPU_HOST_DEVICE AMREX_INLINE Real get_r_from_c(Real spec[NUM_SPECIES],Real bgasconc=zeroval) { diff --git a/models/stagFlowCatalyst/thermo.H b/models/stagFlowCatalyst/thermo.H index 2b8988e..4360673 100644 --- a/models/stagFlowCatalyst/thermo.H +++ b/models/stagFlowCatalyst/thermo.H @@ -15,6 +15,24 @@ using namespace amrex; namespace mflo_thermo { + AMREX_GPU_HOST_DEVICE AMREX_INLINE Real + get_solid_t_from_rhoe(Real rhoe) + { + //modify if conjugate heat transfer + //included + amrex::Real temp=300.0; + return(temp); + } + + AMREX_GPU_HOST_DEVICE AMREX_INLINE Real + get_solid_rhoe_from_t(Real temp) + { + //modify if conjugate + //heat transfer included + amrex::Real Cv_solid=1.0; + amrex::Real rho_solid=1.0; + return(rho_solid*Cv_solid*temp); + } AMREX_GPU_HOST_DEVICE AMREX_INLINE Real get_r_from_c(Real spec[NUM_SPECIES],Real bgasconc=zeroval) { diff --git a/tests/Cartesian/Channel/thermo.H b/tests/Cartesian/Channel/thermo.H index 8b07c21..441094d 100644 --- a/tests/Cartesian/Channel/thermo.H +++ b/tests/Cartesian/Channel/thermo.H @@ -15,6 +15,24 @@ using namespace amrex; namespace mflo_thermo { + AMREX_GPU_HOST_DEVICE AMREX_INLINE Real + get_solid_t_from_rhoe(Real rhoe) + { + //modify if conjugate heat transfer + //included + amrex::Real temp=300.0; + return(temp); + } + + AMREX_GPU_HOST_DEVICE AMREX_INLINE Real + get_solid_rhoe_from_t(Real temp) + { + //modify if conjugate + //heat transfer included + amrex::Real Cv_solid=1.0; + amrex::Real rho_solid=1.0; + return(rho_solid*Cv_solid*temp); + } AMREX_GPU_HOST_DEVICE AMREX_INLINE Real get_r_from_c(Real spec[NUM_SPECIES],Real bgasconc=zeroval) { diff --git a/tests/Cartesian/ShockTube_N2/GNUmakefile b/tests/Cartesian/ShockTube_N2/GNUmakefile index 45afffc..26fc006 100644 --- a/tests/Cartesian/ShockTube_N2/GNUmakefile +++ b/tests/Cartesian/ShockTube_N2/GNUmakefile @@ -12,8 +12,6 @@ USE_MPI = TRUE USE_OMP = FALSE USE_CUDA = FALSE -USE_SEC_ORDER_FLUX = TRUE - Bpack := ./Make.package Blocs := . diff --git a/tests/Cartesian/ShockTube_N2/thermo.H b/tests/Cartesian/ShockTube_N2/thermo.H index 9e4e136..a46460c 100644 --- a/tests/Cartesian/ShockTube_N2/thermo.H +++ b/tests/Cartesian/ShockTube_N2/thermo.H @@ -16,6 +16,24 @@ using namespace amrex; namespace mflo_thermo { + AMREX_GPU_HOST_DEVICE AMREX_INLINE Real + get_solid_t_from_rhoe(Real rhoe) + { + //modify if conjugate heat transfer + //included + amrex::Real temp=300.0; + return(temp); + } + + AMREX_GPU_HOST_DEVICE AMREX_INLINE Real + get_solid_rhoe_from_t(Real temp) + { + //modify if conjugate + //heat transfer included + amrex::Real Cv_solid=1.0; + amrex::Real rho_solid=1.0; + return(rho_solid*Cv_solid*temp); + } extern AMREX_GPU_DEVICE_MANAGED Real R_univ; AMREX_GPU_HOST_DEVICE AMREX_INLINE Real get_r_from_c(Real spec[NUM_SPECIES],Real bgasconc=zeroval) diff --git a/tests/Cartesian/ShockTube_N2He/thermo.H b/tests/Cartesian/ShockTube_N2He/thermo.H index 06f20f8..9af3cec 100644 --- a/tests/Cartesian/ShockTube_N2He/thermo.H +++ b/tests/Cartesian/ShockTube_N2He/thermo.H @@ -15,6 +15,24 @@ using namespace amrex; namespace mflo_thermo { + AMREX_GPU_HOST_DEVICE AMREX_INLINE Real + get_solid_t_from_rhoe(Real rhoe) + { + //modify if conjugate heat transfer + //included + amrex::Real temp=300.0; + return(temp); + } + + AMREX_GPU_HOST_DEVICE AMREX_INLINE Real + get_solid_rhoe_from_t(Real temp) + { + //modify if conjugate + //heat transfer included + amrex::Real Cv_solid=1.0; + amrex::Real rho_solid=1.0; + return(rho_solid*Cv_solid*temp); + } extern AMREX_GPU_DEVICE_MANAGED Real R_univ; AMREX_GPU_HOST_DEVICE AMREX_INLINE Real get_r_from_c(Real spec[NUM_SPECIES],Real bgasconc=zeroval) diff --git a/tests/Cartesian/SpeciesEqTests/thermo.H b/tests/Cartesian/SpeciesEqTests/thermo.H index 5fdf99f..8696851 100644 --- a/tests/Cartesian/SpeciesEqTests/thermo.H +++ b/tests/Cartesian/SpeciesEqTests/thermo.H @@ -15,6 +15,24 @@ using namespace amrex; namespace mflo_thermo { + AMREX_GPU_HOST_DEVICE AMREX_INLINE Real + get_solid_t_from_rhoe(Real rhoe) + { + //modify if conjugate heat transfer + //included + amrex::Real temp=300.0; + return(temp); + } + + AMREX_GPU_HOST_DEVICE AMREX_INLINE Real + get_solid_rhoe_from_t(Real temp) + { + //modify if conjugate + //heat transfer included + amrex::Real Cv_solid=1.0; + amrex::Real rho_solid=1.0; + return(rho_solid*Cv_solid*temp); + } AMREX_GPU_HOST_DEVICE AMREX_INLINE Real get_r_from_c(Real spec[NUM_SPECIES],Real bgasconc=zeroval) { diff --git a/tests/Cartesian/TaylorGreen/thermo.H b/tests/Cartesian/TaylorGreen/thermo.H index f01d7e1..5483e20 100644 --- a/tests/Cartesian/TaylorGreen/thermo.H +++ b/tests/Cartesian/TaylorGreen/thermo.H @@ -15,6 +15,24 @@ using namespace amrex; namespace mflo_thermo { + AMREX_GPU_HOST_DEVICE AMREX_INLINE Real + get_solid_t_from_rhoe(Real rhoe) + { + //modify if conjugate heat transfer + //included + amrex::Real temp=300.0; + return(temp); + } + + AMREX_GPU_HOST_DEVICE AMREX_INLINE Real + get_solid_rhoe_from_t(Real temp) + { + //modify if conjugate + //heat transfer included + amrex::Real Cv_solid=1.0; + amrex::Real rho_solid=1.0; + return(rho_solid*Cv_solid*temp); + } AMREX_GPU_HOST_DEVICE AMREX_INLINE Real get_r_from_c(Real spec[NUM_SPECIES],Real bgasconc=zeroval) { diff --git a/tests/ImmBoundary/subsonicCylinder/thermo.H b/tests/ImmBoundary/subsonicCylinder/thermo.H index f01d7e1..5483e20 100644 --- a/tests/ImmBoundary/subsonicCylinder/thermo.H +++ b/tests/ImmBoundary/subsonicCylinder/thermo.H @@ -15,6 +15,24 @@ using namespace amrex; namespace mflo_thermo { + AMREX_GPU_HOST_DEVICE AMREX_INLINE Real + get_solid_t_from_rhoe(Real rhoe) + { + //modify if conjugate heat transfer + //included + amrex::Real temp=300.0; + return(temp); + } + + AMREX_GPU_HOST_DEVICE AMREX_INLINE Real + get_solid_rhoe_from_t(Real temp) + { + //modify if conjugate + //heat transfer included + amrex::Real Cv_solid=1.0; + amrex::Real rho_solid=1.0; + return(rho_solid*Cv_solid*temp); + } AMREX_GPU_HOST_DEVICE AMREX_INLINE Real get_r_from_c(Real spec[NUM_SPECIES],Real bgasconc=zeroval) { diff --git a/tests/ImmBoundary/supersonicSphere/thermo.H b/tests/ImmBoundary/supersonicSphere/thermo.H index f01d7e1..5483e20 100644 --- a/tests/ImmBoundary/supersonicSphere/thermo.H +++ b/tests/ImmBoundary/supersonicSphere/thermo.H @@ -15,6 +15,24 @@ using namespace amrex; namespace mflo_thermo { + AMREX_GPU_HOST_DEVICE AMREX_INLINE Real + get_solid_t_from_rhoe(Real rhoe) + { + //modify if conjugate heat transfer + //included + amrex::Real temp=300.0; + return(temp); + } + + AMREX_GPU_HOST_DEVICE AMREX_INLINE Real + get_solid_rhoe_from_t(Real temp) + { + //modify if conjugate + //heat transfer included + amrex::Real Cv_solid=1.0; + amrex::Real rho_solid=1.0; + return(rho_solid*Cv_solid*temp); + } AMREX_GPU_HOST_DEVICE AMREX_INLINE Real get_r_from_c(Real spec[NUM_SPECIES],Real bgasconc=zeroval) {