ampsci
High-precision calculations for one- and two-valence atomic systems
AdamsMoulton Namespace Reference

Detailed Description

Contains classes and functions which use general N-step Adams Moulton method to solve systems of 2x2 ODEs, up to N=12.

Classes

struct  AM_Coefs
 Holds the K+1 Adams-Moulton coefficients for the K-step AM method. More...
 
struct  DerivativeMatrix
 Pure-virtual struct defining the derivative matrix for a 2x2 ODE system. More...
 
class  ODESolver2D
 Solves a 2x2 system of ODEs using a K-step Adams-Moulton method. More...
 

Functions

template<typename T , typename U , std::size_t N, std::size_t M>
constexpr T inner_product (const std::array< T, N > &a, const std::array< U, M > &b)
 Inner product of two std::arrays: sum_i a_i * b_i.
 

Variables

template<typename T >
constexpr bool is_complex_v = is_complex<T>::value
 Type trait: true if T is std::complex<U> for some U, false otherwise.
 

Function Documentation

◆ inner_product()

template<typename T , typename U , std::size_t N, std::size_t M>
constexpr T AdamsMoulton::inner_product ( const std::array< T, N > &  a,
const std::array< U, M > &  b 
)
constexpr

Inner product of two std::arrays: sum_i a_i * b_i.

\[ \text{inner\_product}(a,\, b) = \sum_{i=0}^{N-1} a_i \, b_i \]

where \( N = \min(\text{a.size()}, \text{b.size()}) \). The array types may differ (T and U), but U must be convertible to T. The return type is T (same as the first array).

Variable Documentation

◆ is_complex_v

template<typename T >
constexpr bool AdamsMoulton::is_complex_v = is_complex<T>::value
constexpr

Type trait: true if T is std::complex<U> for some U, false otherwise.

Examples:

static_assert(!is_complex_v<double>);
static_assert(!is_complex_v<float>);
static_assert(is_complex_v<std::complex<double>>);
static_assert(is_complex_v<std::complex<float>>);
static_assert(is_complex<std::complex<float>>::value);