From 91362ee1053f592f4b23cdcdb155fc39bc588828 Mon Sep 17 00:00:00 2001
From: Isuru Fernando <isuruf@gmail.com>
Date: Fri, 28 Jun 2019 18:19:16 -0400
Subject: [PATCH] Use fast spmv in p2m
---
sumpy/expansion/__init__.py | 29 ++++++++++++++++++++---------
sumpy/expansion/local.py | 6 +++---
sumpy/expansion/multipole.py | 6 +++---
sumpy/p2e.py | 6 +++---
4 files changed, 29 insertions(+), 18 deletions(-)
diff --git a/sumpy/expansion/__init__.py b/sumpy/expansion/__init__.py
index 2e741eec..baf09da4 100644
--- a/sumpy/expansion/__init__.py
+++ b/sumpy/expansion/__init__.py
@@ -115,12 +115,14 @@ class ExpansionBase(object):
"""
raise NotImplementedError
- def coefficients_from_source(self, avec, bvec, rscale):
+ def coefficients_from_source(self, avec, bvec, rscale, sac=None):
"""Form an expansion from a source point.
:arg avec: vector from source to center.
:arg bvec: vector from center to target. Not usually necessary,
except for line-Taylor expansion.
+ :arg sac: a symbolic assignment collection where temporary
+ expressions are stored.
:returns: a list of :mod:`sympy` expressions representing
the coefficients of the expansion.
@@ -174,11 +176,11 @@ class ExpansionTermsWrangler(object):
raise NotImplementedError
def get_full_kernel_derivatives_from_stored(self, stored_kernel_derivatives,
- rscale):
+ rscale, sac=None):
raise NotImplementedError
def get_stored_mpole_coefficients_from_full(self, full_mpole_coefficients,
- rscale):
+ rscale, sac=None):
raise NotImplementedError
@memoize_method
@@ -224,7 +226,7 @@ class FullExpansionTermsWrangler(ExpansionTermsWrangler):
ExpansionTermsWrangler.get_full_coefficient_identifiers)
def get_full_kernel_derivatives_from_stored(self, stored_kernel_derivatives,
- rscale):
+ rscale, sac=None):
return stored_kernel_derivatives
get_stored_mpole_coefficients_from_full = (
@@ -283,7 +285,7 @@ def _fast_spmv(reconstruct_matrix, vec, sac, transpose=False):
else:
res = []
expr_all = vec.copy()
- for row, deps in reversed(enumerate(reconstruct_matrix)):
+ for row, deps in reversed(list(enumerate(reconstruct_matrix))):
if len(deps) == 0:
res.append(expr_all[row])
continue
@@ -335,7 +337,7 @@ class LinearPDEBasedExpansionTermsWrangler(ExpansionTermsWrangler):
result[col] += coeff * val
return result
return _fast_spmv(reconstruct_matrix, full_mpole_coefficients, sac,
- tranpose=True)
+ transpose=True)
@property
def stored_identifiers(self):
@@ -367,7 +369,7 @@ class LinearPDEBasedExpansionTermsWrangler(ExpansionTermsWrangler):
full_coeffs = self.get_full_coefficient_identifiers()
matrix_rows = []
- count_nonzero_coeff = 0
+ count_nonzero_coeff = -len(stored_identifiers)
for irow, row in six.iteritems(coeff_matrix):
# For eg: (u_xxx / rscale**3) = (u_yy / rscale**2) * coeff1 +
# (u_xx / rscale**2) * coeff2
@@ -422,8 +424,17 @@ class LinearPDEBasedExpansionTermsWrangler(ExpansionTermsWrangler):
pde = self.get_single_pde()
assert len(pde.eqs) == 1
pde_dict = pde.eqs[0]
- max_mi_idx = max(coeff_ident_enumerate_dict[ident.mi] for ident in
- pde_dict.keys() if ident.mi in coeff_ident_enumerate_dict)
+ for ident in pde_dict.keys():
+ if ident.mi not in coeff_ident_enumerate_dict:
+ coeff_matrix = defaultdict(list)
+ reconstruct_matrix = []
+ for i in range(len(mis)):
+ coeff_matrix[i] = [(i, 1)]
+ reconstruct_matrix.append([])
+ return mis, coeff_matrix, reconstruct_matrix
+
+ max_mi_idx = max(coeff_ident_enumerate_dict[ident.mi] for
+ ident in pde_dict.keys())
max_mi = mis[max_mi_idx]
max_mi_coeff = pde_dict[CoeffIdentifier(max_mi, 0)]
max_mi_mult = -1/sym.sympify(max_mi_coeff)
diff --git a/sumpy/expansion/local.py b/sumpy/expansion/local.py
index e2150bfe..526da89d 100644
--- a/sumpy/expansion/local.py
+++ b/sumpy/expansion/local.py
@@ -57,7 +57,7 @@ class LineTaylorLocalExpansion(LocalExpansionBase):
def get_coefficient_identifiers(self):
return list(range(self.order+1))
- def coefficients_from_source(self, avec, bvec, rscale):
+ def coefficients_from_source(self, avec, bvec, rscale, sac=None):
# no point in heeding rscale here--just ignore it
if bvec is None:
raise RuntimeError("cannot use line-Taylor expansions in a setting "
@@ -115,7 +115,7 @@ class VolumeTaylorLocalExpansionBase(LocalExpansionBase):
Coefficients represent derivative values of the kernel.
"""
- def coefficients_from_source(self, avec, bvec, rscale):
+ def coefficients_from_source(self, avec, bvec, rscale, sac=None):
from sumpy.tools import MiDerivativeTaker
ppkernel = self.kernel.postprocess_at_source(
self.kernel.get_expression(avec), avec)
@@ -291,7 +291,7 @@ class _FourierBesselLocalExpansion(LocalExpansionBase):
def get_coefficient_identifiers(self):
return list(range(-self.order, self.order+1))
- def coefficients_from_source(self, avec, bvec, rscale):
+ def coefficients_from_source(self, avec, bvec, rscale, sac=None):
if not self.use_rscale:
rscale = 1
diff --git a/sumpy/expansion/multipole.py b/sumpy/expansion/multipole.py
index 6f31164b..c602e356 100644
--- a/sumpy/expansion/multipole.py
+++ b/sumpy/expansion/multipole.py
@@ -54,7 +54,7 @@ class VolumeTaylorMultipoleExpansionBase(MultipoleExpansionBase):
Coefficients represent the terms in front of the kernel derivatives.
"""
- def coefficients_from_source(self, avec, bvec, rscale):
+ def coefficients_from_source(self, avec, bvec, rscale, sac=None):
from sumpy.kernel import DirectionalSourceDerivative
kernel = self.kernel
@@ -95,7 +95,7 @@ class VolumeTaylorMultipoleExpansionBase(MultipoleExpansionBase):
for mi in self.get_full_coefficient_identifiers()]
return (
self.expansion_terms_wrangler.get_stored_mpole_coefficients_from_full(
- result, rscale))
+ result, rscale, sac=sac))
def get_scaled_multipole(self, expr, bvec, rscale, nderivatives,
nderivatives_for_scaling=None):
@@ -242,7 +242,7 @@ class _HankelBased2DMultipoleExpansion(MultipoleExpansionBase):
def get_coefficient_identifiers(self):
return list(range(-self.order, self.order+1))
- def coefficients_from_source(self, avec, bvec, rscale):
+ def coefficients_from_source(self, avec, bvec, rscale, sac=None):
if not self.use_rscale:
rscale = 1
diff --git a/sumpy/p2e.py b/sumpy/p2e.py
index daa7d93d..f5472fc4 100644
--- a/sumpy/p2e.py
+++ b/sumpy/p2e.py
@@ -86,9 +86,9 @@ class P2EBase(KernelCacheWrapper):
sac = SymbolicAssignmentCollection()
coeff_names = [
- sac.assign_unique("coeff%d" % i, coeff_i)
- for i, coeff_i in enumerate(
- self.expansion.coefficients_from_source(avec, None, rscale))]
+ sac.assign_unique("coeff%d" % i, coeff_i)
+ for i, coeff_i in enumerate(
+ self.expansion.coefficients_from_source(avec, None, rscale, sac))]
sac.run_global_cse()
--
GitLab