openmc.deplete.cram.IPFCramSolver

class openmc.deplete.cram.IPFCramSolver(alpha, theta, alpha0)[source]

CRAM depletion solver that uses incomplete partial factorization

Provides a __call__() that utilizes an incomplete partial factorization (IPF) for the Chebyshev Rational Approximation Method (CRAM), as described in the following paper: M. Pusa, “Higher-Order Chebyshev Rational Approximation Method and Application to Burnup Equations,” Nucl. Sci. Eng., 182:3, 297-318.

Parameters
  • alpha (numpy.ndarray) – Complex residues of poles used in the factorization. Must be a vector with even number of items.

  • theta (numpy.ndarray) – Complex poles. Must have an equal size as alpha.

  • alpha0 (float) – Limit of the approximation at infinity

Variables
  • alpha (numpy.ndarray) – Complex residues of poles theta in the incomplete partial factorization. Denoted as \(\tilde{\alpha}\)

  • theta (numpy.ndarray) – Complex poles \(\theta\) of the rational approximation

  • alpha0 (float) – Limit of the approximation at infinity

__call__(A, n0, dt)[source]

Solve depletion equations using IPF CRAM

Parameters
  • A (scipy.sparse.csr_matrix) – Sparse transmutation matrix A[j, i] desribing rates at which isotope i transmutes to isotope j

  • n0 (numpy.ndarray) – Initial compositions, typically given in number of atoms in some material or an atom density

  • dt (float) – Time [s] of the specific interval to be solved

Returns

Final compositions after dt

Return type

numpy.ndarray