1 #ifndef AMREX_PCG_SOLVER_H_
2 #define AMREX_PCG_SOLVER_H_
3 #include <AMReX_Config.H>
22 template <
int N,
typename T,
typename M,
typename P>
25 M const& mat,
P const& precond,
int maxiter, T rel_tol)
27 static_assert(std::is_floating_point_v<T>);
30 for (
int i = 0; i < N; ++i) {
33 if (rnorm0 == 0) {
return 0; }
38 for (iter = 1; iter <= maxiter; ++iter) {
42 for (
int i = 0; i < N; ++i) { rho +=
r[i]*z[i]; }
43 if (rho == 0) {
break; }
45 for (
int i = 0; i < N; ++i) { p[i] = z[i]; }
47 auto rr = rho * (T(1.0)/rho_prev);
48 for (
int i = 0; i < N; ++i) {
49 p[i] = z[i] + rr * p[i];
55 for (
int i = 0; i < N; ++i) { pq += p[i]*q[i]; }
56 if (pq == 0) {
break; }
57 T alpha = rho * (T(1.0)/pq);
59 for (
int i = 0; i < N; ++i) {
64 if (rnorm <= rnorm0*rel_tol) {
break; }
#define AMREX_FORCE_INLINE
Definition: AMReX_Extension.H:119
#define AMREX_RESTRICT
Definition: AMReX_Extension.H:37
#define AMREX_GPU_HOST_DEVICE
Definition: AMReX_GpuQualifiers.H:20
@ max
Definition: AMReX_ParallelReduce.H:17
static constexpr int M
Definition: AMReX_OpenBC.H:13
static constexpr int P
Definition: AMReX_OpenBC.H:14
Definition: AMReX_Amr.cpp:49
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE T abs(const GpuComplex< T > &a_z) noexcept
Return the absolute value of a complex number.
Definition: AMReX_GpuComplex.H:356
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE int pcg_solve(T *AMREX_RESTRICT x, T *AMREX_RESTRICT r, M const &mat, P const &precond, int maxiter, T rel_tol)
Preconditioned conjugate gradient solver.
Definition: AMReX_PCGSolver.H:24