Skip to content

Commit

Permalink
Merge branch 'development' into 2d_spherical_ptheta
Browse files Browse the repository at this point in the history
  • Loading branch information
@zhichen3
zhichen3 committed Oct 1, 2024
2 parents 765969f + 6ea903d commit 810d9f9
Show file tree
Hide file tree
Showing 14 changed files with 321 additions and 153 deletions.
10 changes: 7 additions & 3 deletions Exec/hydro_tests/gresho_vortex/_prob_params
View file
Original file line number Diff line number Diff line change
@@ -1,11 +1,15 @@
p0 real 1.0_rt y
# initial Mach number
M0 real 0.1_rt y

# ambient density
rho0 real 1.0_rt y

# reference pressure -- this is set automatically
p0 real -1.0_rt

# reference timescale
t_r real 1.0_rt y

x_r real 0.0_rt

q_r real 0.0_rt

nsub integer 4 y
12 changes: 2 additions & 10 deletions Exec/hydro_tests/gresho_vortex/inputs.2d
View file
Original file line number Diff line number Diff line change
Expand Up @@ -79,16 +79,8 @@ amr.derive_plot_vars=ALL
# PROBLEM PARAMETERS
problem.rho0 = 1.0

# p0 = rho0 u_phi**2/(gamma * M**2) - 1/2 -- specify M to get p0

# M ~ 0.001
problem.p0 = 599999.5

# M ~ 0.01
# problem.p0 = 5999.5

# M ~ 0.1
# problem.p0 = 59.5
# Mach number of the flow
problem.M0 = 0.1

# EOS
eos.eos_assume_neutral = 1
11 changes: 11 additions & 0 deletions Exec/hydro_tests/gresho_vortex/problem_initialize.H
View file
Original file line number Diff line number Diff line change
Expand Up @@ -27,6 +27,17 @@ void problem_initialize ()
problem::x_r = probhi[0] - problo[0];
problem::q_r = 0.4_rt * M_PI * problem::x_r / problem::t_r;

// pressure peaks at r = 1/5 from the center
// where p = p0 + 25/2 r^2, so peak pressure is p = p0 + 1/2
//
// Mach number is |u_phi| / sqrt(gamma p / rho), so we get solve
// for p0. Note that the peak velocity is q,
//
// p0 = rho q^2 / (gamma M^2) - 1/2

problem::p0 = problem::rho0 * std::pow(problem::q_r, 2.0) /
(eos_rp::eos_gamma * std::pow(problem::M0, 2.0)) - 0.5_rt;

}

#endif
55 changes: 15 additions & 40 deletions Exec/hydro_tests/gresho_vortex/problem_initialize_state_data.H
View file
Original file line number Diff line number Diff line change
Expand Up @@ -22,71 +22,46 @@ void problem_initialize_state_data (int i, int j, int k,
const Real* dx = geomdata.CellSize();

Real x = problo[0] + dx[0] * (static_cast<Real>(i) + 0.5_rt);
Real xl = problo[0] + dx[0] * (static_cast<Real>(i));

Real y = 0.0;
Real yl = 0.0;

#if AMREX_SPACEDIM >= 2
y = problo[1] + dx[1] * (static_cast<Real>(j) + 0.5_rt);
yl = problo[1] + dx[1] * (static_cast<Real>(j));
#endif

Real z = 0.0;
Real zl = 0.0;

#if AMREX_SPACEDIM == 3
z = problo[2] + dx[2] * (static_cast<Real>(k) + 0.5_rt);
zl = problo[2] + dx[2] * (static_cast<Real>(k));
#endif

Real reint = 0.0;
Real u_tot = 0.0;

for (int kk = 0; kk < problem::nsub; kk++) {
Real zz = zl + dx[2] * (static_cast<Real>(kk) + 0.5_rt) / problem::nsub;

for (int jj = 0; jj < problem::nsub; jj++) {
Real yy = yl + dx[1] * (static_cast<Real>(jj) + 0.5_rt) / problem::nsub;

for (int ii = 0; ii < problem::nsub; ii++) {
Real xx = xl + dx[0] * (static_cast<Real>(ii) + 0.5_rt) / problem::nsub;

Real r = std::sqrt((xx - problem::center[0]) * (xx - problem::center[0]) +
(yy - problem::center[1]) * (yy - problem::center[1]));

Real u_phi;
Real p;
Real r = std::sqrt(amrex::Math::powi<2>(x - problem::center[0]) +
amrex::Math::powi<2>(y - problem::center[1]));

if (r < 0.2_rt) {
u_phi = 5.0_rt * r;
p = problem::p0 + 12.5_rt * r * r;
Real u_phi;
Real p;

} else if (r < 0.4_rt) {
u_phi = 2.0_rt - 5.0_rt * r;
p = problem::p0 + 12.5_rt * r * r + 4.0_rt *
(1.0_rt - 5.0_rt * r - std::log(0.2_rt) + std::log(r));
if (r < 0.2_rt) {
u_phi = 5.0_rt * r;
p = problem::p0 + 12.5_rt * r * r;

} else {
u_phi = 0.0_rt;
p = problem::p0 - 2.0_rt + 4.0_rt * std::log(2.0_rt);
}
} else if (r < 0.4_rt) {
u_phi = 2.0_rt - 5.0_rt * r;
p = problem::p0 + 12.5_rt * r * r + 4.0_rt *
(1.0_rt - 5.0_rt * r - std::log(0.2_rt) + std::log(r));

u_tot += u_phi;
reint += p/(eos_rp::eos_gamma - 1.0_rt);
}
}
} else {
u_phi = 0.0_rt;
p = problem::p0 - 2.0_rt + 4.0_rt * std::log(2.0_rt);
}

Real u_phi = u_tot / (problem::nsub * problem::nsub * problem::nsub);
reint = reint / (problem::nsub * problem::nsub * problem::nsub);
reint = p / (eos_rp::eos_gamma - 1.0_rt);

state(i,j,k,URHO) = problem::rho0;

// phi unit vector: \hat{\phi} = -sin(phi) \hat{x} + cos(phi) \hat{y}
// with cos(phi) = x/r; sin(phi) = y/r
Real r = std::sqrt((x - problem::center[0]) * (x - problem::center[0]) +
(y - problem::center[1]) * (y - problem::center[1]));

// -sin(phi) = y/r
state(i,j,k,UMX) = -problem::rho0 * problem::q_r * u_phi * ((y - problem::center[1]) / r);
Expand Down
4 changes: 2 additions & 2 deletions Exec/science/flame_wave/ci-benchmarks/grid_diag.out
View file
Original file line number Diff line number Diff line change
@@ -1,5 +1,5 @@
# COLUMN 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23
# TIMESTEP TIME MASS XMOM YMOM ZMOM ANG. MOM. X ANG. MOM. Y ANG. MOM. Z KIN. ENERGY INT. ENERGY GAS ENERGY GRAV. ENERGY TOTAL ENERGY CENTER OF MASS X-LOC CENTER OF MASS Y-LOC CENTER OF MASS Z-LOC CENTER OF MASS X-VEL CENTER OF MASS Y-VEL CENTER OF MASS Z-VEL MAXIMUM TEMPERATURE MAXIMUM DENSITY MAXIMUM T_S / T_E
0 0.0000000000000000e+00 1.4700685690736437e+20 0.0000000000000000e+00 0.0000000000000000e+00 0.0000000000000000e+00 0.0000000000000000e+00 0.0000000000000000e+00 0.0000000000000000e+00 0.0000000000000000e+00 2.9163021935541997e+37 2.9163021935541997e+37 0.0000000000000000e+00 2.9163021935541997e+37 2.7306641645125415e+04 6.8087525965685256e+02 0.0000000000000000e+00 0.0000000000000000e+00 0.0000000000000000e+00 0.0000000000000000e+00 1.3928678114307070e+09 3.0715027784567930e+07 0.0000000000000000e+00
1 5.9806047036494849e-08 1.4700700110495903e+20 3.3564221294760178e+23 7.2876252590129013e+23 1.2612378576902532e+20 3.1318392026636616e+23 -2.7141686393328650e+24 1.9066097490832445e+28 7.8605333893026802e+29 2.9163081863790971e+37 2.9163082649844082e+37 0.0000000000000000e+00 2.9163082649844082e+37 2.7306641717925890e+04 6.8087482283157362e+02 0.0000000000000000e+00 2.2831716205676648e+03 4.9573321027137572e+03 8.5794407627549896e-01 1.3929561430929687e+09 3.0715326643214259e+07 3.3873763068343888e-04
2 1.2559269877663918e-07 1.4700717750115456e+20 7.0484895627215180e+23 1.5304335700608418e+24 5.5619062084401319e+20 1.3810901029080806e+24 -1.1969287510263666e+25 4.0039912973182075e+28 3.0853296206485770e+30 2.9163145904134766e+37 2.9163148989464290e+37 0.0000000000000000e+00 2.9163148989464290e+37 2.7306641954244220e+04 6.8087461718998884e+02 0.0000000000000000e+00 4.7946567525018709e+03 1.0410604407725754e+04 3.7834249340624542e+00 1.3933944931271181e+09 3.0715562157661017e+07 3.3857586853705914e-04
1 5.9806047036494849e-08 1.4700700110495903e+20 3.3564221217849212e+23 7.2876252590226133e+23 1.2612378556916341e+20 3.1318391966955516e+23 -2.7141686347841742e+24 1.9066097493127676e+28 7.8605324800324878e+29 2.9163081863790706e+37 2.9163082649843738e+37 0.0000000000000000e+00 2.9163082649843738e+37 2.7306641717925955e+04 6.8087482283157317e+02 0.0000000000000000e+00 2.2831716153358752e+03 4.9573321027203638e+03 8.5794407491595881e-01 1.3929561431310942e+09 3.0715326643212341e+07 3.3873763041849706e-04
2 1.2559269877663918e-07 1.4700717750115462e+20 7.0484895445039834e+23 1.5304335700823335e+24 5.5619061962888426e+20 1.3810900989697918e+24 -1.1969287484132617e+25 4.0039912979259806e+28 3.0853275723280981e+30 2.9163145904135894e+37 2.9163148989463411e+37 0.0000000000000000e+00 2.9163148989463411e+37 2.7306641954244522e+04 6.8087461718998338e+02 0.0000000000000000e+00 4.7946567401095926e+03 1.0410604407871944e+04 3.7834249257966728e+00 1.3933944939675827e+09 3.0715562157563787e+07 3.3857586854903492e-04
4 changes: 2 additions & 2 deletions Exec/science/flame_wave/ci-benchmarks/species_diag.out
View file
Original file line number Diff line number Diff line change
@@ -1,5 +1,5 @@
# COLUMN 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15
# TIMESTEP TIME Mass He4 Mass C12 Mass O16 Mass Ne20 Mass Mg24 Mass Si28 Mass S32 Mass Ar36 Mass Ca40 Mass Ti44 Mass Cr48 Mass Fe52 Mass Ni56
0 0.0000000000000000e+00 2.8142378112400895e-15 7.3932235330395157e-24 7.3932235330395157e-24 7.3932235330395157e-24 7.3932235330395157e-24 7.3932235330395157e-24 7.3932235330395157e-24 7.3932235330395157e-24 7.3932235330395157e-24 7.3932235330395157e-24 7.3932235330395157e-24 7.3932235330395157e-24 7.1117997526547418e-14
1 5.9806047036494849e-08 2.8142323856670296e-15 5.4329394513687282e-21 7.4194476167127106e-24 7.3934515289195991e-24 7.3935541950036108e-24 7.3933935990475538e-24 7.3932321468746301e-24 7.3932312019180565e-24 7.3932308328287436e-24 7.3932307869393029e-24 7.3932307850741563e-24 7.3932307849933102e-24 7.1118070045957091e-14
2 1.2559269877663918e-07 2.8142264165380905e-15 1.1401993461614438e-20 7.4929577845086585e-24 7.3944696393124482e-24 7.3939305113509452e-24 7.3935838899709814e-24 7.3932425308155010e-24 7.3932405326302037e-24 7.3932397568503172e-24 7.3932396603620439e-24 7.3932396564405536e-24 7.3932396562705448e-24 7.1118158758588142e-14
1 5.9806047036494849e-08 2.8142323856670446e-15 5.4329394398418374e-21 7.4194443455948811e-24 7.3934514829832856e-24 7.3935541949306949e-24 7.3933935989353088e-24 7.3932321468740865e-24 7.3932312019180345e-24 7.3932308328287333e-24 7.3932307869392985e-24 7.3932307850741519e-24 7.3932307849933044e-24 7.1118070045957103e-14
2 1.2559269877663918e-07 2.8142264165381050e-15 1.1401993450709496e-20 7.4929528148974785e-24 7.3944694290388374e-24 7.3939305042586412e-24 7.3935838894001731e-24 7.3932425308122493e-24 7.3932405326301655e-24 7.3932397568503098e-24 7.3932396603620498e-24 7.3932396564405566e-24 7.3932396562705492e-24 7.1118158758588142e-14
46 changes: 23 additions & 23 deletions Exec/science/wdmerger/ci-benchmarks/wdmerger_collision_2D.out
View file
Original file line number Diff line number Diff line change
@@ -1,29 +1,29 @@
plotfile = plt00086
time = 1.25
variables minimum value maximum value
density 8.693611703e-05 19565534.299
xmom -5.4964100651e+14 1.3559128302e+14
ymom -2.5530096328e+15 2.5530122744e+15
density 8.6936468335e-05 19568007.309
xmom -5.4959005475e+14 1.3559838981e+14
ymom -2.5540206053e+15 2.5540214496e+15
zmom 0 0
rho_E 7.4982062146e+11 5.0669247219e+24
rho_e 7.1077581849e+11 5.0640768326e+24
Temp 242288.68588 1409652233.6
rho_He4 8.693611703e-17 3.599903302
rho_C12 3.4774446812e-05 7825956.6934
rho_O16 5.2161670217e-05 11739149.75
rho_Ne20 8.693611703e-17 181951.0571
rho_Mg24 8.693611703e-17 1192.7969729
rho_Si28 8.693611703e-17 6.6913702949
rho_S32 8.693611703e-17 0.00019493291655
rho_Ar36 8.693611703e-17 1.9565534609e-05
rho_Ca40 8.693611703e-17 1.9565534331e-05
rho_Ti44 8.693611703e-17 1.9565534308e-05
rho_Cr48 8.693611703e-17 1.9565534308e-05
rho_Fe52 8.693611703e-17 1.9565534308e-05
rho_Ni56 8.693611703e-17 1.9565534308e-05
phiGrav -5.8708033902e+17 -2.3375498341e+16
grav_x -684991644 -51428.243166
grav_y -739606241.84 739606820.44
rho_E 7.4987846833e+11 5.0676543192e+24
rho_e 7.1083104678e+11 5.0648015049e+24
Temp 242292.44287 1409581841.1
rho_He4 8.6936468335e-17 3.5973411848
rho_C12 3.4774587333e-05 7826946.398
rho_O16 5.2161881e-05 11740633.889
rho_Ne20 8.6936468335e-17 181819.49413
rho_Mg24 8.6936468335e-17 1190.7712386
rho_Si28 8.6936468335e-17 6.6817951193
rho_S32 8.6936468335e-17 0.00019451843176
rho_Ar36 8.6936468335e-17 1.9568007618e-05
rho_Ca40 8.6936468335e-17 1.9568007341e-05
rho_Ti44 8.6936468335e-17 1.9568007318e-05
rho_Cr48 8.6936468335e-17 1.9568007318e-05
rho_Fe52 8.6936468335e-17 1.9568007318e-05
rho_Ni56 8.6936468335e-17 1.9568007318e-05
phiGrav -5.8709462562e+17 -2.3375498549e+16
grav_x -685026429.13 -51428.265677
grav_y -739654246.49 739654206.24
grav_z 0 0
rho_enuc -2.7633982574e+12 7.6429034885e+23
rho_enuc -4.7815621457e+12 7.6360058391e+23

38 changes: 37 additions & 1 deletion Source/driver/timestep.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -35,6 +35,9 @@ Castro::estdt_cfl (int is_new)
// Courant-condition limited timestep

const auto dx = geom.CellSizeArray();
const auto problo = geom.ProbLoArray();
const auto coord = geom.Coord();
amrex::ignore_unused(problo, coord);

const MultiFab& stateMF = is_new ? get_new_data(State_Type) : get_old_data(State_Type);

Expand Down Expand Up @@ -82,6 +85,11 @@ Castro::estdt_cfl (int is_new)
Real dt2;
#if AMREX_SPACEDIM >= 2
dt2 = dx[1]/(c + std::abs(uy));
if (coord == 2) {
// dx[1] in Spherical2D is just dtheta, need rdtheta for physical length
// so just multiply by the smallest r
dt2 *= problo[0] + 0.5_rt * dx[0];
}
#else
dt2 = dt1;
#endif
Expand Down Expand Up @@ -127,6 +135,9 @@ Castro::estdt_mhd (int is_new)

// MHD timestep limiter
const auto dx = geom.CellSizeArray();
const auto problo = geom.ProbLoArray();
const auto coord = geom.Coord();
amrex::ignore_unused(problo, coord);

const MultiFab& U_state = is_new ? get_new_data(State_Type) : get_old_data(State_Type);

Expand Down Expand Up @@ -206,6 +217,9 @@ Castro::estdt_mhd (int is_new)
Real dt2;
#if AMREX_SPACEDIM >= 2
dt2 = dx[1]/(cy + std::abs(uy));
if (coord == 2) {
dt2 *= problo[0] + 0.5_rt * dx[0];
}
#else
dt2 = dt1;
#endif
Expand Down Expand Up @@ -238,6 +252,9 @@ Castro::estdt_temp_diffusion (int is_new)
// where D = k/(rho c_v), and k is the conductivity

const auto dx = geom.CellSizeArray();
const auto problo = geom.ProbLoArray();
const auto coord = geom.Coord();
amrex::ignore_unused(problo, coord);

const MultiFab& stateMF = is_new ? get_new_data(State_Type) : get_old_data(State_Type);

Expand Down Expand Up @@ -287,6 +304,10 @@ Castro::estdt_temp_diffusion (int is_new)
Real dt2;
#if AMREX_SPACEDIM >= 2
dt2 = 0.5_rt * dx[1]*dx[1] / D;
if (coord == 2) {
Real r = problo[0] + 0.5_rt * dx[0];
dt2 *= r * r;
}
#else
dt2 = dt1;
#endif
Expand Down Expand Up @@ -320,6 +341,9 @@ Castro::estdt_burning (int is_new)
}

const auto dx = geom.CellSizeArray();
const auto problo = geom.ProbLoArray();
const auto coord = geom.Coord();
amrex::ignore_unused(problo, coord);

MultiFab& stateMF = is_new ? get_new_data(State_Type) : get_old_data(State_Type);

Expand Down Expand Up @@ -368,7 +392,13 @@ Castro::estdt_burning (int is_new)
#if AMREX_SPACEDIM == 1
burn_state.dx = dx[0];
#else
burn_state.dx = amrex::min(AMREX_D_DECL(dx[0], dx[1], dx[2]));
Real dx1 = dx[1];
#if AMREX_SPACEDIM >= 2
if (coord == 2) {
dx1 *= problo[0] + 0.5_rt * dx[0];
}
#endif
burn_state.dx = amrex::min(AMREX_D_DECL(dx[0], dx1, dx[2]));
#endif

burn_state.rho = S(i,j,k,URHO);
Expand Down Expand Up @@ -464,6 +494,9 @@ Real
Castro::estdt_rad (int is_new)
{
auto dx = geom.CellSizeArray();
const auto problo = geom.ProbLoArray();
const auto coord = geom.Coord();
amrex::ignore_unused(problo, coord);

const MultiFab& stateMF = is_new ? get_new_data(State_Type) : get_old_data(State_Type);
const MultiFab& radMF = is_new ? get_new_data(Rad_Type) : get_old_data(Rad_Type);
Expand Down Expand Up @@ -523,6 +556,9 @@ Castro::estdt_rad (int is_new)
Real dt1 = dx[0] / (c + std::abs(ux));
#if AMREX_SPACEDIM >= 2
Real dt2 = dx[1] / (c + std::abs(uy));
if (coord == 2) {
dt2 *= problo[0] + 0.5_rt * dx[0];
}
#else
Real dt2 = std::numeric_limits<Real>::max();
#endif
Expand Down
8 changes: 8 additions & 0 deletions Source/hydro/Castro_hydro.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -237,12 +237,20 @@ Castro::check_for_cfl_violation(const MultiFab& State, const Real dt)
int cfl_violation = 0;

auto dx = geom.CellSizeArray();
const auto problo = geom.ProbLoArray();
const auto coord = geom.Coord();
amrex::ignore_unused(problo, coord);

Real dtdx = dt / dx[0];

Real dtdy = 0.0_rt;
if (AMREX_SPACEDIM >= 2) {
dtdy = dt / dx[1];
if (coord == 2) {
// dx[1] in Spherical2D is just rdtheta, need rdtheta for physical length
// Just choose to divide by the smallest r
dtdy /= problo[0] + 0.5_rt * dx[0];
}
}

Real dtdz = 0.0_rt;
Expand Down
Loading

0 comments on commit 810d9f9

Please sign in to comment.