From ce20082e31ff12b2981b9f52e26bc4d5ec5f4d7a Mon Sep 17 00:00:00 2001 From: Michael Zingale Date: Thu, 30 May 2024 14:35:43 -0400 Subject: [PATCH] allow for an eccentricity --- Exec/science/circular_det/_prob_params | 9 ++++++--- Exec/science/circular_det/inputs | 7 ++++--- .../circular_det/problem_initialize_state_data.H | 12 ++++++++---- 3 files changed, 18 insertions(+), 10 deletions(-) diff --git a/Exec/science/circular_det/_prob_params b/Exec/science/circular_det/_prob_params index e781e95b79..44e0db3a3a 100644 --- a/Exec/science/circular_det/_prob_params +++ b/Exec/science/circular_det/_prob_params @@ -10,13 +10,16 @@ nfrac real 0.0_rt y ofrac real 0.0_rt y +# dimensionless width of transition between hot and cold w_T real 5.e-4_rt y -center_T real 0.3_rt y +# dimensionless semi-major axis +a_T real 0.3_rt y -smallx real 1.e-12_rt y +# eccentricity +ecc_T real 0.0 y -idir integer 1 y +smallx real 1.e-12_rt y ihe4 integer -1 diff --git a/Exec/science/circular_det/inputs b/Exec/science/circular_det/inputs index cfb75abb68..e86483d1d9 100644 --- a/Exec/science/circular_det/inputs +++ b/Exec/science/circular_det/inputs @@ -40,7 +40,7 @@ castro.use_flattening = 1 castro.riemann_solver = 1 # TIME STEP CONTROL -castro.cfl = 0.5 # cfl number for hyperbolic system +castro.cfl = 0.4 # cfl number for hyperbolic system castro.init_shrink = 0.1 # scale back initial timestep castro.change_max = 1.05 # scale back initial timestep @@ -82,8 +82,9 @@ problem.smallx = 1.e-10 problem.idir = 1 -problem.w_T = 1.e-3 -problem.center_T = 0.1 +problem.w_T = 1.e-2 +problem.a_T = 0.15 +problem.ecc_T = 0.7 # refinement diff --git a/Exec/science/circular_det/problem_initialize_state_data.H b/Exec/science/circular_det/problem_initialize_state_data.H index 9c820e2845..579e3346c8 100644 --- a/Exec/science/circular_det/problem_initialize_state_data.H +++ b/Exec/science/circular_det/problem_initialize_state_data.H @@ -15,18 +15,22 @@ void problem_initialize_state_data (int i, int j, int k, const Real* problo = geomdata.ProbLo(); const Real* probhi = geomdata.ProbHi(); - Real width = problem::w_T * (probhi[0] - problo[0]); - Real c_T = problem::center_T * (probhi[0] - problo[0]); + // semi-major axis + Real a = problem::a_T * (probhi[0] - problo[0]); + + // semi-minor axis + Real b = a * std::sqrt(1.0 - problem::ecc_T * problem::ecc_T); // compute distance from the center Real xx = problo[0] + dx[0] * (static_cast(i) + 0.5_rt) - problem::center[0]; Real yy = problo[1] + dx[1] * (static_cast(j) + 0.5_rt) - problem::center[1]; - Real r = std::sqrt(xx * xx + yy * yy); + Real r = std::sqrt((xx / a) * (xx / a) + + (yy / b) * (yy / b)); state(i,j,k,URHO) = problem::dens; - Real sigma = 1.0_rt / (1.0_rt + std::exp(-(c_T - r)/ width)); + Real sigma = 1.0_rt / (1.0_rt + std::exp(-(1.0_rt - r)/ problem::w_T)); state(i,j,k,UTEMP) = problem::T_l + (problem::T_r - problem::T_l) * (1.0_rt - sigma);