Skip to content
New issue

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

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

Already on GitHub? # to your account

added diagnostic_interval int for nonlinear solver diagnostic. #5748

Open
wants to merge 3 commits into
base: development
Choose a base branch
from
Open
Show file tree
Hide file tree
Changes from all commits
Commits
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
7 changes: 7 additions & 0 deletions Source/FieldSolver/ImplicitSolvers/ImplicitSolver.H
Original file line number Diff line number Diff line change
Expand Up @@ -96,6 +96,8 @@ public:
[[nodiscard]] amrex::Array<amrex::LinOpBCType,AMREX_SPACEDIM> GetLinOpBCLo () const;
[[nodiscard]] amrex::Array<amrex::LinOpBCType,AMREX_SPACEDIM> GetLinOpBCHi () const;

[[nodiscard]] amrex::Real GetTheta () const { return m_theta; }

protected:

/**
Expand All @@ -115,6 +117,11 @@ protected:
*/
mutable amrex::Real m_dt = 0.0;

/**
* \brief Time-biasing parameter for fields used on RHS to advance system.
*/
amrex::Real m_theta = 0.5;

/**
* \brief Nonlinear solver type and object
*/
Expand Down
2 changes: 1 addition & 1 deletion Source/FieldSolver/ImplicitSolvers/SemiImplicitEM.H
Original file line number Diff line number Diff line change
Expand Up @@ -63,7 +63,7 @@ public:

void ComputeRHS ( WarpXSolverVec& a_RHS,
const WarpXSolverVec& a_E,
amrex::Real half_time,
amrex::Real start_time,
int a_nl_iter,
bool a_from_jacobian ) override;

Expand Down
5 changes: 3 additions & 2 deletions Source/FieldSolver/ImplicitSolvers/SemiImplicitEM.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -80,7 +80,7 @@ void SemiImplicitEM::OneStep ( amrex::Real start_time,
// Solve nonlinear system for Eg at t_{n+1/2}
// Particles will be advanced to t_{n+1/2}
m_E.Copy(m_Eold); // initial guess for Eg^{n+1/2}
m_nlsolver->Solve( m_E, m_Eold, half_time, 0.5_rt*m_dt, a_step );
m_nlsolver->Solve( m_E, m_Eold, start_time, m_dt, a_step );

// Update WarpX owned Efield_fp to t_{n+1/2}
m_WarpX->SetElectricFieldAndApplyBCs( m_E, half_time );
Expand All @@ -103,12 +103,13 @@ void SemiImplicitEM::OneStep ( amrex::Real start_time,

void SemiImplicitEM::ComputeRHS ( WarpXSolverVec& a_RHS,
const WarpXSolverVec& a_E,
amrex::Real half_time,
amrex::Real start_time,
int a_nl_iter,
bool a_from_jacobian )
{
// Update WarpX-owned Efield_fp using current state of Eg from
// the nonlinear solver at time n+theta
const amrex::Real half_time = start_time + 0.5_rt*m_dt;
m_WarpX->SetElectricFieldAndApplyBCs( a_E, half_time );

// Update particle positions and velocities using the current state
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -64,7 +64,7 @@ public:

void ComputeRHS ( WarpXSolverVec& a_RHS,
const WarpXSolverVec& a_E,
amrex::Real half_time,
amrex::Real start_time,
int a_nl_iter,
bool a_from_jacobian ) override;

Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -81,7 +81,7 @@ void StrangImplicitSpectralEM::OneStep ( amrex::Real start_time,

// Solve nonlinear system for E at t_{n+1/2}
// Particles will be advanced to t_{n+1/2}
m_nlsolver->Solve( m_E, m_Eold, half_time, 0.5_rt*m_dt, a_step );
m_nlsolver->Solve( m_E, m_Eold, start_time, m_dt, a_step );

// Update WarpX owned Efield_fp and Bfield_fp to t_{n+1/2}
UpdateWarpXFields( m_E, half_time );
Expand All @@ -101,12 +101,13 @@ void StrangImplicitSpectralEM::OneStep ( amrex::Real start_time,

void StrangImplicitSpectralEM::ComputeRHS ( WarpXSolverVec& a_RHS,
WarpXSolverVec const & a_E,
amrex::Real half_time,
amrex::Real start_time,
int a_nl_iter,
bool a_from_jacobian )
{
// Update WarpX-owned Efield_fp and Bfield_fp using current state of
// E from the nonlinear solver at time n+1/2
const amrex::Real half_time = start_time + 0.5_rt*m_dt;
UpdateWarpXFields( a_E, half_time );

// Self consistently update particle positions and velocities using the
Expand Down
5 changes: 0 additions & 5 deletions Source/FieldSolver/ImplicitSolvers/ThetaImplicitEM.H
Original file line number Diff line number Diff line change
Expand Up @@ -86,11 +86,6 @@ public:

private:

/**
* \brief Time-biasing parameter for fields used on RHS to advance system
*/
amrex::Real m_theta = 0.5;

/**
* \brief Solver vectors to be used in the nonlinear solver to solve for the
* electric field E. The main logic for determining which variables should be
Expand Down
2 changes: 1 addition & 1 deletion Source/FieldSolver/ImplicitSolvers/ThetaImplicitEM.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -106,7 +106,7 @@ void ThetaImplicitEM::OneStep ( const amrex::Real start_time,
// Solve nonlinear system for Eg at t_{n+theta}
// Particles will be advanced to t_{n+1/2}
m_E.Copy(m_Eold); // initial guess for Eg^{n+theta}
m_nlsolver->Solve( m_E, m_Eold, start_time, m_theta*m_dt, a_step );
m_nlsolver->Solve( m_E, m_Eold, start_time, m_dt, a_step );

// Update WarpX owned Efield_fp and Bfield_fp to t_{n+theta}
UpdateWarpXFields( m_E, start_time );
Expand Down
6 changes: 3 additions & 3 deletions Source/NonlinearSolvers/CurlCurlMLMGPC.H
Original file line number Diff line number Diff line change
Expand Up @@ -272,8 +272,8 @@ void CurlCurlMLMGPC<T,Ops>::Update (const T& a_U)
amrex::ignore_unused(a_U);

// set the coefficients alpha and beta for curl-curl op
// (m_dt here is actually theta<=0.5 times simulation dt)
const RT alpha = (this->m_dt*PhysConst::c) * (this->m_dt*PhysConst::c);
const RT thetaDt = m_ops->GetTheta()*this->m_dt;
const RT alpha = (thetaDt*PhysConst::c) * (thetaDt*PhysConst::c);
const RT beta = RT(1.0);

// currently not implemented in 1D
Expand All @@ -283,7 +283,7 @@ void CurlCurlMLMGPC<T,Ops>::Update (const T& a_U)

if (m_verbose) {
amrex::Print() << "Updating " << amrex::getEnumNameString(PreconditionerType::pc_curl_curl_mlmg)
<< ": theta*dt = " << this->m_dt << ", "
<< ": theta*dt = " << thetaDt << ", "
<< " coefficients: "
<< "alpha = " << alpha << ", "
<< "beta = " << beta << "\n";
Expand Down
8 changes: 5 additions & 3 deletions Source/NonlinearSolvers/NewtonSolver.H
Original file line number Diff line number Diff line change
Expand Up @@ -245,6 +245,7 @@ void NewtonSolver<Vec,Ops>::ParseParameters ()
pp_newton.query("max_iterations", m_maxits);
pp_newton.query("require_convergence", m_require_convergence);
pp_newton.query("diagnostic_file", this->m_diagnostic_file);
pp_newton.query("diagnostic_interval", this->m_diagnostic_interval);

const amrex::ParmParse pp_gmres("gmres");
pp_gmres.query("verbose_int", m_gmres_verbose_int);
Expand Down Expand Up @@ -358,12 +359,13 @@ void NewtonSolver<Vec,Ops>::Solve ( Vec& a_U,
}
}

if (!this->m_diagnostic_file.empty() && amrex::ParallelDescriptor::IOProcessor()) {
if (!this->m_diagnostic_file.empty() && amrex::ParallelDescriptor::IOProcessor() &&
((a_step+1)%this->m_diagnostic_interval==0 || a_step==0)) {
std::ofstream diagnostic_file{this->m_diagnostic_file, std::ofstream::out | std::ofstream::app};
diagnostic_file << std::setprecision(14);
diagnostic_file << a_step;
diagnostic_file << a_step+1;
diagnostic_file << " ";;
diagnostic_file << a_time;
diagnostic_file << a_time + a_dt;
diagnostic_file << " ";;
diagnostic_file << iter;
diagnostic_file << " ";;
Expand Down
1 change: 1 addition & 0 deletions Source/NonlinearSolvers/NonlinearSolver.H
Original file line number Diff line number Diff line change
Expand Up @@ -83,6 +83,7 @@ protected:
bool m_is_defined = false;
mutable bool m_verbose = true;
std::string m_diagnostic_file;
int m_diagnostic_interval = 1;

};

Expand Down
9 changes: 5 additions & 4 deletions Source/NonlinearSolvers/PicardSolver.H
Original file line number Diff line number Diff line change
Expand Up @@ -156,6 +156,7 @@ void PicardSolver<Vec,Ops>::ParseParameters ()
pp_picard.query("max_iterations", m_maxits);
pp_picard.query("require_convergence", m_require_convergence);
pp_picard.query("diagnostic_file", this->m_diagnostic_file);
pp_picard.query("diagnostic_interval", this->m_diagnostic_interval);
}

template <class Vec, class Ops>
Expand All @@ -169,7 +170,6 @@ void PicardSolver<Vec,Ops>::Solve ( Vec& a_U,
WARPX_ALWAYS_ASSERT_WITH_MESSAGE(
this->m_is_defined,
"PicardSolver::Solve() called on undefined object");
amrex::ignore_unused(a_dt);
using namespace amrex::literals;

//
Expand Down Expand Up @@ -244,12 +244,13 @@ void PicardSolver<Vec,Ops>::Solve ( Vec& a_U,
}
}

if (!this->m_diagnostic_file.empty() && amrex::ParallelDescriptor::IOProcessor()) {
if (!this->m_diagnostic_file.empty() && amrex::ParallelDescriptor::IOProcessor() &&
(a_step%this->m_diagnostic_interval==0 || a_step==0)) {
std::ofstream diagnostic_file{this->m_diagnostic_file, std::ofstream::out | std::ofstream::app};
diagnostic_file << std::setprecision(14);
diagnostic_file << a_step;
diagnostic_file << a_step+1;
diagnostic_file << " ";
diagnostic_file << a_time;
diagnostic_file << a_time + a_dt;
diagnostic_file << " ";
diagnostic_file << iter;
diagnostic_file << " ";
Expand Down
Loading