Skip to content

Commit

Permalink
add renormalization of the species after the update
Browse files Browse the repository at this point in the history
in some sense, this mirrors what we would do with Strang or
simplified-SDC, since we always renormalize after leaving the
burner.
  • Loading branch information
zingale committed Oct 10, 2023
1 parent b512a60 commit d860c8f
Show file tree
Hide file tree
Showing 2 changed files with 22 additions and 12 deletions.
26 changes: 16 additions & 10 deletions Source/sdc/Castro_sdc.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -133,6 +133,12 @@ Castro::do_sdc_update(int m_start, int m_end, Real dt)
const Box& bx = mfi.tilebox();
const Box& bx1 = mfi.growntilebox(1);

// this is the starting data
Array4<const Real> const& k_new_m_start_arr = (k_new[m_start])->array(mfi);

// this is where the update will be stored
Array4<Real> const& k_new_m_end_arr = (k_new[m_end])->array(mfi);

#ifdef REACTIONS
// advection + reactions
if (sdc_order == 2)
Expand Down Expand Up @@ -171,16 +177,16 @@ Castro::do_sdc_update(int m_start, int m_end, Real dt)

}

auto k_m = (*k_new[m_start]).array(mfi);
auto k_n = (*k_new[m_end]).array(mfi);
auto A_m = (*A_new[m_start]).array(mfi);
auto A_n = (*A_new[m_end]).array(mfi);
auto C_arr = C2.array();

amrex::ParallelFor(bx,
[=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept
{
sdc_update_o2(i, j, k, k_m, k_n, A_m, A_n, C_arr, dt_m, sdc_iteration, m_start);
sdc_update_o2(i, j, k,
k_new_m_start_arr, k_new_m_end_arr,
A_m, A_n, C_arr, dt_m, sdc_iteration, m_start);
});
}
else
Expand All @@ -189,10 +195,7 @@ Castro::do_sdc_update(int m_start, int m_end, Real dt)
// fourth order SDC reaction update -- we need to respect the
// difference between cell-centers and averages

Array4<const Real> const& k_new_m_start_arr=
(k_new[m_start])->array(mfi);
Array4<Real> const& k_new_m_end_arr=(k_new[m_end])->array(mfi);
Array4<const Real> const& C_source_arr=C_source.array(mfi);
Array4<const Real> const& C_source_arr = C_source.array(mfi);

// convert the starting U to cell-centered on a fab-by-fab basis
// -- including one ghost cell
Expand Down Expand Up @@ -264,9 +267,6 @@ Castro::do_sdc_update(int m_start, int m_end, Real dt)

}
#else
Array4<const Real> const& k_new_m_start_arr=
(k_new[m_start])->array(mfi);
Array4<Real> const& k_new_m_end_arr=(k_new[m_end])->array(mfi);
Array4<const Real> const& A_new_arr=(A_new[m_start])->array(mfi);
Array4<const Real> const& A_old_0_arr=(A_old[0])->array(mfi);
Array4<const Real> const& A_old_1_arr=(A_old[1])->array(mfi);
Expand Down Expand Up @@ -317,6 +317,12 @@ Castro::do_sdc_update(int m_start, int m_end, Real dt)
}
#endif

amrex::ParallelFor(bx,
[=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept
{
normalize_species_sdc(i, j, k, k_new_m_end_arr);
});

}
}

Expand Down
8 changes: 6 additions & 2 deletions Source/sdc/sdc_newton_solve.H
Original file line number Diff line number Diff line change
Expand Up @@ -109,7 +109,7 @@ sdc_newton_solve(const Real dt_m,

Array1D<Real, 1, NumSpec+1> f;

const int MAX_ITER = 100;
const int MAX_ITER = 200;

ierr = newton::NEWTON_SUCCESS;

Expand Down Expand Up @@ -251,7 +251,7 @@ sdc_newton_subdivide(const Real dt_m,
// converges or reaches our limit on the number of
// subintervals.

const int MAX_NSUB = 64;
const int MAX_NSUB = 128;
GpuArray<Real, NUM_STATE> U_begin;

// subdivide the timestep and do multiple Newtons. We come
Expand Down Expand Up @@ -308,6 +308,10 @@ sdc_newton_subdivide(const Real dt_m,
}
nsub *= 2;
}

if (ierr != newton::NEWTON_SUCCESS) {
std::cout << "Falled. " << ierr << std::endl;
}
}
#endif

Expand Down

0 comments on commit d860c8f

Please sign in to comment.