From 7e0370bde813b35727cebaba16525455c9e60706 Mon Sep 17 00:00:00 2001 From: Ethan Peterson Date: Wed, 8 Nov 2023 11:14:46 -0500 Subject: [PATCH] change from csr to csc one more change to CSC --- openmc/deplete/abc.py | 6 +++--- openmc/deplete/chain.py | 10 +++++----- openmc/deplete/cram.py | 4 ++-- 3 files changed, 10 insertions(+), 10 deletions(-) diff --git a/openmc/deplete/abc.py b/openmc/deplete/abc.py index d135b2dae5a..5d906efd840 100644 --- a/openmc/deplete/abc.py +++ b/openmc/deplete/abc.py @@ -543,7 +543,7 @@ class Integrator(ABC): User-supplied functions are expected to have the following signature: ``solver(A, n0, t) -> n1`` where - * ``A`` is a :class:`scipy.sparse.csr_matrix` making up the + * ``A`` is a :class:`scipy.sparse.csc_matrix` making up the depletion matrix * ``n0`` is a 1-D :class:`numpy.ndarray` of initial compositions for a given material in atoms/cm3 @@ -921,7 +921,7 @@ class SIIntegrator(Integrator): User-supplied functions are expected to have the following signature: ``solver(A, n0, t) -> n1`` where - * ``A`` is a :class:`scipy.sparse.csr_matrix` making up the + * ``A`` is a :class:`scipy.sparse.csc_matrix` making up the depletion matrix * ``n0`` is a 1-D :class:`numpy.ndarray` of initial compositions for a given material in atoms/cm3 @@ -1023,7 +1023,7 @@ def __call__(self, A, n0, dt): Parameters ---------- - A : scipy.sparse.csr_matrix + A : scipy.sparse.csc_matrix Sparse transmutation matrix ``A[j, i]`` describing rates at which isotope ``i`` transmutes to isotope ``j`` n0 : numpy.ndarray diff --git a/openmc/deplete/chain.py b/openmc/deplete/chain.py index dc01a1b34f3..e06cb572530 100644 --- a/openmc/deplete/chain.py +++ b/openmc/deplete/chain.py @@ -596,7 +596,7 @@ def form_matrix(self, rates, fission_yields=None): Returns ------- - scipy.sparse.csr_matrix + scipy.sparse.csc_matrix Sparse matrix representing depletion. See Also @@ -680,11 +680,11 @@ def form_matrix(self, rates, fission_yields=None): # Clear set of reactions reactions.clear() - # Use DOK matrix as intermediate representation, then convert to CSR and return + # Use DOK matrix as intermediate representation, then convert to CSC and return n = len(self) matrix_dok = sp.dok_matrix((n, n)) dict.update(matrix_dok, matrix) - return matrix_dok.tocsr() + return matrix_dok.tocsc() def form_rr_term(self, tr_rates, mats): """Function to form the transfer rate term matrices. @@ -713,7 +713,7 @@ def form_rr_term(self, tr_rates, mats): Returns ------- - scipy.sparse.csr_matrix + scipy.sparse.csc_matrix Sparse matrix representing transfer term. """ @@ -746,7 +746,7 @@ def form_rr_term(self, tr_rates, mats): n = len(self) matrix_dok = sp.dok_matrix((n, n)) dict.update(matrix_dok, matrix) - return matrix_dok.tocsr() + return matrix_dok.tocsc() def get_branch_ratios(self, reaction="(n,gamma)"): """Return a dictionary with reaction branching ratios diff --git a/openmc/deplete/cram.py b/openmc/deplete/cram.py index f5f80dfdc48..53de83bb68e 100644 --- a/openmc/deplete/cram.py +++ b/openmc/deplete/cram.py @@ -75,9 +75,9 @@ def __call__(self, A, n0, dt): Final compositions after ``dt`` """ - A = sp.csr_matrix(A * dt, dtype=np.float64) + A = dt * sp.csc_matrix(A, dtype=np.float64) y = n0.copy() - ident = sp.eye(A.shape[0]) + ident = sp.eye(A.shape[0], format='csc') for alpha, theta in zip(self.alpha, self.theta): y += 2*np.real(alpha*sla.spsolve(A - theta*ident, y)) return y * self.alpha0