WarpX
Public Member Functions | Protected Member Functions | Protected Attributes | List of all members
ImplicitSolver Class Referenceabstract

#include <ImplicitSolver.H>

Inheritance diagram for ImplicitSolver:
SemiImplicitEM ThetaImplicitEM

Public Member Functions

 ImplicitSolver ()=default
 
virtual ~ImplicitSolver ()=default
 
 ImplicitSolver (const ImplicitSolver &)=delete
 
ImplicitSolveroperator= (const ImplicitSolver &)=delete
 
 ImplicitSolver (ImplicitSolver &&)=delete
 
ImplicitSolveroperator= (ImplicitSolver &&)=delete
 
virtual void Define (WarpX *a_WarpX)=0
 Read user-provided parameters that control the implicit solver. Allocate internal arrays for intermediate field values needed by the solver. More...
 
bool IsDefined () const
 
virtual void PrintParameters () const =0
 
void GetParticleSolverParams (int &a_max_particle_iter, amrex::ParticleReal &a_particle_tol) const
 
virtual void OneStep (amrex::Real a_time, amrex::Real a_dt, int a_step)=0
 Advance fields and particles by one time step using the specified implicit algorithm. More...
 
virtual void ComputeRHS (WarpXSolverVec &a_RHS, const WarpXSolverVec &a_E, amrex::Real a_time, amrex::Real a_dt, int a_nl_iter, bool a_from_jacobian)=0
 Computes the RHS of the equation corresponding to the specified implicit algorithm. The discrete equations corresponding to numerical integration of ODEs are often written in the form U = b + RHS(U), where U is the variable to be solved for (e.g., the solution at the next time step), b is a constant (i.e., the solution from the previous time step), and RHS(U) is the right-hand-side of the equation. Iterative solvers, such as Picard and Newton, and higher-order Runge-Kutta methods, need to compute RHS(U) multiple times for different values of U. Thus, a routine that returns this value is needed. e.g., Ebar - E^n = cvac^2*0.5*dt*(curl(Bbar) - mu0*Jbar(Ebar,Bbar)) Here, U = Ebar, b = E^n, and the expression on the right is RHS(U). More...
 

Protected Member Functions

void parseNonlinearSolverParams (const amrex::ParmParse &pp)
 parse nonlinear solver parameters (if one is used) More...
 

Protected Attributes

WarpXm_WarpX
 Pointer back to main WarpX class. More...
 
bool m_is_defined = false
 
NonlinearSolverType m_nlsolver_type
 Nonlinear solver type and object. More...
 
std::unique_ptr< NonlinearSolver< WarpXSolverVec, ImplicitSolver > > m_nlsolver
 
amrex::ParticleReal m_particle_tolerance = 1.0e-10
 tolerance used by the iterative method used to obtain a self-consistent update of the particle positions and velocities for given E and B on the grid More...
 
int m_max_particle_iterations = 21
 maximum iterations for the iterative method used to obtain a self-consistent update of the particle positions and velocities for given E and B on the grid More...
 

Constructor & Destructor Documentation

◆ ImplicitSolver() [1/3]

ImplicitSolver::ImplicitSolver ( )
default

◆ ~ImplicitSolver()

virtual ImplicitSolver::~ImplicitSolver ( )
virtualdefault

◆ ImplicitSolver() [2/3]

ImplicitSolver::ImplicitSolver ( const ImplicitSolver )
delete

◆ ImplicitSolver() [3/3]

ImplicitSolver::ImplicitSolver ( ImplicitSolver &&  )
delete

Member Function Documentation

◆ ComputeRHS()

virtual void ImplicitSolver::ComputeRHS ( WarpXSolverVec a_RHS,
const WarpXSolverVec a_E,
amrex::Real  a_time,
amrex::Real  a_dt,
int  a_nl_iter,
bool  a_from_jacobian 
)
pure virtual

Computes the RHS of the equation corresponding to the specified implicit algorithm. The discrete equations corresponding to numerical integration of ODEs are often written in the form U = b + RHS(U), where U is the variable to be solved for (e.g., the solution at the next time step), b is a constant (i.e., the solution from the previous time step), and RHS(U) is the right-hand-side of the equation. Iterative solvers, such as Picard and Newton, and higher-order Runge-Kutta methods, need to compute RHS(U) multiple times for different values of U. Thus, a routine that returns this value is needed. e.g., Ebar - E^n = cvac^2*0.5*dt*(curl(Bbar) - mu0*Jbar(Ebar,Bbar)) Here, U = Ebar, b = E^n, and the expression on the right is RHS(U).

Implemented in ThetaImplicitEM, and SemiImplicitEM.

◆ Define()

virtual void ImplicitSolver::Define ( WarpX a_WarpX)
pure virtual

Read user-provided parameters that control the implicit solver. Allocate internal arrays for intermediate field values needed by the solver.

Implemented in ThetaImplicitEM, and SemiImplicitEM.

◆ GetParticleSolverParams()

void ImplicitSolver::GetParticleSolverParams ( int a_max_particle_iter,
amrex::ParticleReal &  a_particle_tol 
) const
inline

◆ IsDefined()

bool ImplicitSolver::IsDefined ( ) const
inline

◆ OneStep()

virtual void ImplicitSolver::OneStep ( amrex::Real  a_time,
amrex::Real  a_dt,
int  a_step 
)
pure virtual

Advance fields and particles by one time step using the specified implicit algorithm.

Implemented in ThetaImplicitEM, and SemiImplicitEM.

◆ operator=() [1/2]

ImplicitSolver& ImplicitSolver::operator= ( const ImplicitSolver )
delete

◆ operator=() [2/2]

ImplicitSolver& ImplicitSolver::operator= ( ImplicitSolver &&  )
delete

◆ parseNonlinearSolverParams()

void ImplicitSolver::parseNonlinearSolverParams ( const amrex::ParmParse pp)
inlineprotected

parse nonlinear solver parameters (if one is used)

◆ PrintParameters()

virtual void ImplicitSolver::PrintParameters ( ) const
pure virtual

Implemented in ThetaImplicitEM, and SemiImplicitEM.

Member Data Documentation

◆ m_is_defined

bool ImplicitSolver::m_is_defined = false
protected

◆ m_max_particle_iterations

int ImplicitSolver::m_max_particle_iterations = 21
protected

maximum iterations for the iterative method used to obtain a self-consistent update of the particle positions and velocities for given E and B on the grid

◆ m_nlsolver

std::unique_ptr<NonlinearSolver<WarpXSolverVec,ImplicitSolver> > ImplicitSolver::m_nlsolver
protected

◆ m_nlsolver_type

NonlinearSolverType ImplicitSolver::m_nlsolver_type
protected

Nonlinear solver type and object.

◆ m_particle_tolerance

amrex::ParticleReal ImplicitSolver::m_particle_tolerance = 1.0e-10
protected

tolerance used by the iterative method used to obtain a self-consistent update of the particle positions and velocities for given E and B on the grid

◆ m_WarpX

WarpX* ImplicitSolver::m_WarpX
protected

Pointer back to main WarpX class.


The documentation for this class was generated from the following file: