Skip to content

Commit b3da1a3

Browse files
zhichen3zhichen3
authored
scale T to T9 manually in NSE solve and allow hybrj to do optimzations internally (#1884)
This seems to help with the second approach of solving (rho, e, Ye) mode. Previously I think it also does the scaling but via hybrj.mode=2 with the diag values I passed in. But perhaps its better to scale T manually and then allow hybrj to scale internally via hybrj.mode=1.
1 parent 7bb013a commit b3da1a3

4 files changed

Lines changed: 24 additions & 28 deletions

File tree

‎nse_solver/nse_eqns.H‎

Lines changed: 9 additions & 7 deletions
Original file line numberDiff line numberDiff line change
@@ -250,7 +250,7 @@ void nse_fcn(amrex::Array1D<amrex::Real, 1, 3>& x,
250250
auto nse_state = state_data.state;
251251
nse_state.mu_p = x(1);
252252
nse_state.mu_n = x(2);
253-
nse_state.T = x(3);
253+
nse_state.T = x(3) * 1.0e9_rt;
254254

255255
// Finish evaluating NSE X_k
256256

@@ -384,7 +384,7 @@ void nse_jcn(amrex::Array1D<amrex::Real, 1, 3>& x,
384384
auto nse_state = state_data.state;
385385
nse_state.mu_p = x(1);
386386
nse_state.mu_n = x(2);
387-
nse_state.T = x(3);
387+
nse_state.T = x(3) * 1.0e9_rt;
388388

389389
const amrex::Real T_in = nse_state.T;
390390
const amrex::Real e_in = state_data.state.e;
@@ -471,20 +471,20 @@ void nse_jcn(amrex::Array1D<amrex::Real, 1, 3>& x,
471471

472472
amrex::Real dX_dmu_p = nse_state.xn[n] * zion[n] * ikTMeV;
473473
amrex::Real dX_dmu_n = nse_state.xn[n] * (aion[n] - zion[n]) * ikTMeV;
474-
amrex::Real dX_dT = nse_state.xn[n] *
474+
amrex::Real dX_dT9 = 1.0e9_rt * nse_state.xn[n] *
475475
(dpf_dT / pf - du_c_dT + (1.5_rt - (zion[n] * nse_state.mu_p +
476476
(aion[n] - zion[n]) * nse_state.mu_n +
477477
network::bion(n+1)) * ikTMeV) / T_in);
478478

479479
fjac(1, 1) += dX_dmu_p;
480480
fjac(1, 2) += dX_dmu_n;
481-
fjac(1, 3) += dX_dT;
481+
fjac(1, 3) += dX_dT9;
482482
fjac(2, 1) += zion[n] * aion_inv[n] * dX_dmu_p;
483483
fjac(2, 2) += zion[n] * aion_inv[n] * dX_dmu_n;
484-
fjac(2, 3) += zion[n] * aion_inv[n] * dX_dT;
484+
fjac(2, 3) += zion[n] * aion_inv[n] * dX_dT9;
485485
fjac(3, 1) += aion_inv[n] * dX_dmu_p;
486486
fjac(3, 2) += aion_inv[n] * dX_dmu_n;
487-
fjac(3, 3) += aion_inv[n] * dX_dT;
487+
fjac(3, 3) += aion_inv[n] * dX_dT9;
488488
}
489489

490490
//
@@ -501,13 +501,15 @@ void nse_jcn(amrex::Array1D<amrex::Real, 1, 3>& x,
501501
eos_state.T = T_in;
502502
eos(eos_input_rt, eos_state);
503503

504+
amrex::Real de_dT9 = 1.0e9_rt * eos_state.dedT;
505+
504506
fjac(3, 1) = - eos_state.abar * eos_state.abar *
505507
(eos_state.dedA + Ye * eos_state.dedZ) * fjac(3, 1) / e_in;
506508

507509
fjac(3, 2) = - eos_state.abar * eos_state.abar *
508510
(eos_state.dedA + Ye * eos_state.dedZ) * fjac(3, 2) / e_in;
509511

510-
fjac(3, 3) = (eos_state.dedT - eos_state.abar * eos_state.abar *
512+
fjac(3, 3) = (de_dT9 - eos_state.abar * eos_state.abar *
511513
(eos_state.dedA + Ye * eos_state.dedZ) * fjac(3, 3)) / e_in;
512514
}
513515

‎nse_solver/nse_solver.H‎

Lines changed: 8 additions & 14 deletions
Original file line numberDiff line numberDiff line change
@@ -56,17 +56,9 @@ void nse_hybrid_solver(nse_solver_data<T>& state_data,
5656
// mode = 1: scales internally
5757
// mode = 2: scales based on user-supplied diag
5858

59-
if constexpr(neqs == 2){
60-
hj.mode = 1;
61-
} else {
62-
hj.mode = 2;
63-
64-
// Scale variable, this is to work with hj.mode = 2
65-
66-
hj.diag(1) = 1.0_rt;
67-
hj.diag(2) = 1.0_rt;
68-
hj.diag(3) = 1e9_rt;
69-
}
59+
// Choose mode 1 for all.
60+
// If T is unknown, we manually scale to work with T9
61+
hj.mode = 1;
7062

7163
// random flag number for evaluation in for loop;
7264

@@ -85,7 +77,8 @@ void nse_hybrid_solver(nse_solver_data<T>& state_data,
8577
outer_x(1) = state_data.state.mu_p;
8678
outer_x(2) = state_data.state.mu_n;
8779
if constexpr (neqs == 3) {
88-
outer_x(3) = state_data.state.T;
80+
// For temperature, solve for T9.
81+
outer_x(3) = state_data.state.T * 1.e-9_rt;
8982
}
9083

9184
// fine tuning initial guesses
@@ -115,7 +108,7 @@ void nse_hybrid_solver(nse_solver_data<T>& state_data,
115108
bool converged = true;
116109

117110
if constexpr (neqs == 3) {
118-
if (hj.x(3) < 1.0e9_rt || hj.x(3) >= 2e10_rt) {
111+
if (hj.x(3) < 1.0_rt || hj.x(3) >= 20.0_rt) {
119112
// If temperature is too low or even negative.
120113
// Or if the temperature is too high
121114
// Then the solved Temperature is not correct.
@@ -134,7 +127,8 @@ void nse_hybrid_solver(nse_solver_data<T>& state_data,
134127
state_data.state.mu_p = hj.x(1);
135128
state_data.state.mu_n = hj.x(2);
136129
if constexpr (neqs == 3) {
137-
state_data.state.T = hj.x(3);
130+
// Convert to T from T9
131+
state_data.state.T = hj.x(3) * 1.e9_rt;
138132
}
139133
return;
140134
}

‎unit_test/nse_net_cell/burn_cell.H‎

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -126,7 +126,7 @@ void burn_cell_c()
126126

127127
// See the result of constraint eq and jacobian using the solved input
128128

129-
amrex::Array1D<amrex::Real, 1, 3> x {nse_state.mu_p, nse_state.mu_n, nse_state.T};
129+
amrex::Array1D<amrex::Real, 1, 3> x {nse_state.mu_p, nse_state.mu_n, nse_state.T * 1.e-9};
130130
amrex::Array1D<amrex::Real, 1, 3> f;
131131
amrex::Array2D<amrex::Real, 1, 3, 1, 3> jac;
132132
int flag = 0;

‎unit_test/nse_net_cell/ci-benchmarks/nse_net_unit_test.out‎

Lines changed: 6 additions & 6 deletions
Original file line numberDiff line numberDiff line change
@@ -1,5 +1,5 @@
1-
Initializing AMReX (25.09-20-g6e2d7bbde65e)...
2-
AMReX (25.09-20-g6e2d7bbde65e) initialized
1+
Initializing AMReX (24.10-342-g82fa0de69508)...
2+
AMReX (24.10-342-g82fa0de69508) initialized
33
starting the single zone burn...
44
reading in network electron-capture / beta-decay tables...
55
chemical potential of proton is -3
@@ -33,13 +33,13 @@ change in T: 6000000000 6069144727
3333
change in abar: 7.026880148 6.26238884
3434
recovered energy: 2.784234255e+18 2.784234255e+18
3535
Constraint Eqns using solved input are: -4.135580767e-15 0 -2.220446049e-16
36-
Constraint Jacobian using solved input are: 19.70206532 19.80498652 3.069602949e-09 9.82702147 9.875043853 1.533907452e-09 0.2616834419 0.2616834419 7.450188201e-10
36+
Constraint Jacobian using solved input are: 19.70206532 19.80498652 3.069602949 9.82702147 9.875043853 1.533907452 0.2616834419 0.2616834419 0.7450188201
3737
-----------------------------------------------
3838
Using nse.solve_nse_e_mode == 2 :
3939
change in T: 6000000000 6069144727
4040
change in abar: 7.026880148 6.26238884
4141
recovered energy: 2.784234255e+18 2.784234255e+18
42-
Constraint Eqns using solved input are: 1.218192214e-12 4.196643033e-14 8.792966355e-14
43-
Constraint Jacobian using solved input are: 19.70206532 19.80498652 3.069602949e-09 9.82702147 9.875043853 1.533907452e-09 0.2616834419 0.2616834419 7.450188201e-10
42+
Constraint Eqns using solved input are: -6.831202271e-13 -3.014255512e-14 -3.907985047e-13
43+
Constraint Jacobian using solved input are: 19.70206532 19.80498652 3.069602949 9.82702147 9.875043853 1.533907452 0.2616834419 0.2616834419 0.7450188201
4444
-----------------------------------------------
45-
AMReX (25.09-20-g6e2d7bbde65e) finalized
45+
AMReX (24.10-342-g82fa0de69508) finalized

0 commit comments

Comments
 (0)