Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

fix 4th order diffusion convergence + diffusion_test fixes #2592

Merged
merged 3 commits into from
Oct 5, 2023
Merged

fix 4th order diffusion convergence + diffusion_test fixes #2592

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
3 commits
Select commit Hold shift + click to select a range
aa22acf
Some diffusion_test fixes
zingale Oct 2, 2023
fb42e6e
fix convergence
zingale Oct 5, 2023
2d42a29
Merge branch 'development' into fix_4th_diffusion
zingale Oct 5, 2023
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
20 changes: 13 additions & 7 deletions Exec/unit_tests/diffusion_test/Prob.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -21,8 +21,8 @@ void Castro::problem_post_simulation(Vector<std::unique_ptr<AmrLevel> >& amr_lev
Castro& castro = dynamic_cast<Castro&>(*amr_level[n]);
Real time = castro.get_state_data(State_Type).curTime();

const int* domain_lo = castro.geom.Domain().loVect();
const int* domain_hi = castro.geom.Domain().hiVect();
auto domain_lo = castro.geom.Domain().loVect3d();
auto domain_hi = castro.geom.Domain().hiVect3d();

// the state data
MultiFab& S = castro.get_new_data(State_Type);
Expand All @@ -33,12 +33,18 @@ void Castro::problem_post_simulation(Vector<std::unique_ptr<AmrLevel> >& amr_lev
#ifdef TRUE_SDC
// if we are fourth-order, we need to convert to averages
if (sdc_order == 4) {
for (MFIter mfi(*analytic); mfi.isValid(); ++mfi) {
FArrayBox tmp;

const Box& gbx = mfi.growntilebox(1);
ca_make_fourth_in_place(AMREX_INT_ANYD(gbx.loVect()), AMREX_INT_ANYD(gbx.hiVect()),
BL_TO_FORTRAN_FAB((*analytic)[mfi]),
AMREX_INT_ANYD(domain_lo), AMREX_INT_ANYD(domain_hi));
for (MFIter mfi(*analytic); mfi.isValid(); ++mfi) {

const Box& bx = mfi.tilebox();
tmp.resize(bx, 1);
auto tmp_arr = tmp.array();

castro.make_fourth_in_place(bx,
analytic->array(mfi),
tmp_arr,
domain_lo, domain_hi);

}
}
Expand Down
8 changes: 4 additions & 4 deletions Exec/unit_tests/diffusion_test/README.md
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -69,10 +69,10 @@ analytic solution, giving:
A convergence test of the 4th-order SDC algorithm can be run as:

```
./Castro1d.gnu.ex inputs.1d castro.time_integration_method=2 castro.sdc_order=4 amr.n_cell=64
./Castro1d.gnu.ex inputs.1d castro.time_integration_method=2 castro.sdc_order=4 amr.n_cell=128
./Castro1d.gnu.ex inputs.1d castro.time_integration_method=2 castro.sdc_order=4 amr.n_cell=256
./Castro1d.gnu.ex inputs.1d castro.time_integration_method=2 castro.sdc_order=4 amr.n_cell=512
./Castro1d.gnu.TRUESDC.ex inputs.1d castro.time_integration_method=2 castro.sdc_order=4 amr.n_cell=64
./Castro1d.gnu.TRUESDC.ex inputs.1d castro.time_integration_method=2 castro.sdc_order=4 amr.n_cell=128
./Castro1d.gnu.TRUESDC.ex inputs.1d castro.time_integration_method=2 castro.sdc_order=4 amr.n_cell=256
./Castro1d.gnu.TRUESDC.ex inputs.1d castro.time_integration_method=2 castro.sdc_order=4 amr.n_cell=512
```

Note: this is Cartesian, not spherical (since we don't have spherical
Expand Down
30 changes: 15 additions & 15 deletions Source/hydro/fourth_order.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -883,40 +883,40 @@ Castro::fourth_add_diffusive_flux(const Box& bx,
if (idir == 0) {

if (is_avg) {
// we are working with the cell-center state
dTdx = (-q_arr(i+1,j,k,temp_comp) + 27*q_arr(i,j,k,temp_comp) -
27*q_arr(i-1,j,k,temp_comp) + q_arr(i-2,j,k,temp_comp)) / (24.0_rt * dx[0]);

} else {
// we are working with the cell-average state
dTdx = (-q_arr(i+1,j,k,temp_comp) + 15*q_arr(i,j,k,temp_comp) -
15*q_arr(i-1,j,k,temp_comp) + q_arr(i-2,j,k,temp_comp)) / (12.0_rt * dx[0]);

} else {
// we are working with the cell-center state
dTdx = (-q_arr(i+1,j,k,temp_comp) + 27*q_arr(i,j,k,temp_comp) -
27*q_arr(i-1,j,k,temp_comp) + q_arr(i-2,j,k,temp_comp)) / (24.0_rt * dx[0]);
}

} else if (idir == 1) {

if (is_avg) {
// we are working with the cell-center state
dTdx = (-q_arr(i,j+1,k,temp_comp) + 27*q_arr(i,j,k,temp_comp) -
27*q_arr(i,j-1,k,temp_comp) + q_arr(i,j-2,k,temp_comp)) / (24.0_rt * dx[1]);

} else {
// we are working with the cell-average state
dTdx = (-q_arr(i,j+1,k,temp_comp) + 15*q_arr(i,j,k,temp_comp) -
15*q_arr(i,j-1,k,temp_comp) + q_arr(i,j-2,k,temp_comp)) / (12.0_rt * dx[1]);

} else {
// we are working with the cell-center state
dTdx = (-q_arr(i,j+1,k,temp_comp) + 27*q_arr(i,j,k,temp_comp) -
27*q_arr(i,j-1,k,temp_comp) + q_arr(i,j-2,k,temp_comp)) / (24.0_rt * dx[1]);
}

} else {

if (is_avg) {
// we are working with the cell-center state
dTdx = (-q_arr(i,j,k+1,temp_comp) + 27*q_arr(i,j,k,temp_comp) -
27*q_arr(i,j,k-1,temp_comp) + q_arr(i,j,k-2,temp_comp)) / (24.0_rt * dx[2]);

} else {
// we are working with the cell-average state
dTdx = (-q_arr(i,j,k+1,temp_comp) + 15*q_arr(i,j,k,temp_comp) -
15*q_arr(i,j,k-1,temp_comp) + q_arr(i,j,k-2,temp_comp)) / (12.0_rt * dx[2]);

} else {
// we are working with the cell-center state
dTdx = (-q_arr(i,j,k+1,temp_comp) + 27*q_arr(i,j,k,temp_comp) -
27*q_arr(i,j,k-1,temp_comp) + q_arr(i,j,k-2,temp_comp)) / (24.0_rt * dx[2]);
}

}
Expand Down