Unverified Commit b9efcf23 authored by Stephan Rave's avatar Stephan Rave Committed by GitHub

Merge pull request #960 from TiKeil/avoid_nested_parameter_functionals

Avoid nested parameter functionals and functions for sums and products
parents 817951ee 1f701744
Pipeline #64839 passed with stages
in 58 minutes and 33 seconds
......@@ -14,6 +14,7 @@
* Tim Keil, tim.keil@uni-muenster.de
* Energy product in elliptic discretizer
* rename estimate --> estimate_error and estimator -> error_estimator
* avoid nested Product and Lincomb Functionals and Functions
## pyMOR 2020.1
......
......@@ -52,7 +52,7 @@ class Function(ParametricObject):
"""Shorthand for :meth:`~Function.evaluate`."""
return self.evaluate(x, mu)
def __add__(self, other):
def _add_sub(self, other, sign):
if isinstance(other, Number) and other == 0:
return self
elif not isinstance(other, Function):
......@@ -61,22 +61,68 @@ class Function(ParametricObject):
if np.all(other == 0.):
return self
other = ConstantFunction(other, dim_domain=self.dim_domain)
return LincombFunction([self, other], [1., 1.])
__radd__ = __add__
if self.name != 'LincombFunction' or not isinstance(self, LincombFunction):
if other.name == 'LincombFunction' and isinstance(other, LincombFunction):
functions = (self,) + other.functions
coefficients = (1.,) + (other.coefficients if sign == 1. else tuple(-c for c in other.coefficients))
else:
functions, coefficients = (self, other), (1., sign)
elif other.name == 'LincombFunction' and isinstance(other, LincombFunction):
functions = self.functions + other.functions
coefficients = self.coefficients + (other.coefficients if sign == 1.
else tuple(-c for c in other.coefficients))
else:
functions, coefficients = self.functions + (other,), self.coefficients + (sign,)
def __sub__(self, other):
if isinstance(other, Function):
return LincombFunction([self, other], [1., -1.])
return LincombFunction(functions, coefficients)
def _radd_sub(self, other, sign):
assert not isinstance(other, Function) # handled by __add__/__sub__
if isinstance(other, Number) and other == 0:
return self
other = np.array(other)
assert other.shape == self.shape_range
if np.all(other == 0.):
return self
other = ConstantFunction(other, dim_domain=self.dim_domain)
if self.name != 'LincombFunction' or not isinstance(self, LincombFunction):
functions, coefficients = (other, self), (1., sign)
else:
return self + (- np.array(other))
functions = (other,) + self.functions
coefficients = (1.,) + (self.coefficients if sign == 1. else tuple(-c for c in self.coefficients))
return LincombFunction(functions, coefficients)
def __add__(self, other):
return self._add_sub(other, 1.)
def __sub__(self, other):
return self._add_sub(other, -1.)
def __radd__(self, other):
return self._radd_sub(other, 1.)
def __rsub__(self, other):
return self._radd_sub(other, -1.)
def __mul__(self, other):
if not isinstance(other, (Number, ParameterFunctional, Function)):
return NotImplemented
if isinstance(other, (Number, ParameterFunctional)):
return LincombFunction([self], [other])
if isinstance(other, Function):
return ProductFunction([self, other])
return NotImplemented
if self.name != 'ProductFunction' or not isinstance(self, ProductFunction):
if isinstance(other, ProductFunction) and other.name == 'ProductFunction':
return other.with_(functions=other.functions + [self])
else:
return ProductFunction([self, other])
elif isinstance(other, ProductFunction) and other.name == 'ProductFunction':
functions = self.functions + other.functions
return ProductFunction(functions)
else:
return self.with_(functions=self.functions + [other])
__rmul__ = __mul__
......@@ -247,6 +293,9 @@ class LincombFunction(Function):
assert all(isinstance(c, (ParameterFunctional, Number)) for c in coefficients)
assert all(f.dim_domain == functions[0].dim_domain for f in functions[1:])
assert all(f.shape_range == functions[0].shape_range for f in functions[1:])
functions = tuple(functions)
coefficients = tuple(coefficients)
self.__auto_init(locals())
self.dim_domain = functions[0].dim_domain
self.shape_range = functions[0].shape_range
......
......@@ -221,57 +221,6 @@ class LincombOperator(Operator):
def as_source_array(self, mu=None):
return self._as_array(True, mu)
def _add_sub(self, other, sign):
if not isinstance(other, Operator):
return NotImplemented
if self.name != 'LincombOperator':
if isinstance(other, LincombOperator) and other.name == 'LincombOperator':
operators = (self,) + other.operators
coefficients = (1.,) + (other.coefficients if sign == 1. else tuple(-c for c in other.coefficients))
else:
operators, coefficients = (self, other), (1., sign)
elif isinstance(other, LincombOperator) and other.name == 'LincombOperator':
operators = self.operators + other.operators
coefficients = self.coefficients + (other.coefficients if sign == 1.
else tuple(-c for c in other.coefficients))
else:
operators, coefficients = self.operators + (other,), self.coefficients + (sign,)
return LincombOperator(operators, coefficients, solver_options=self.solver_options)
def _radd_sub(self, other, sign):
if not isinstance(other, Operator):
return NotImplemented
# note that 'other' can never be a LincombOperator
if self.name != 'LincombOperator':
operators, coefficients = (other, self), (1., sign)
else:
operators = (other,) + self.operators
coefficients = (1.,) + (self.coefficients if sign == 1. else tuple(-c for c in self.coefficients))
return LincombOperator(operators, coefficients, solver_options=other.solver_options)
def __add__(self, other):
return self._add_sub(other, 1.)
def __sub__(self, other):
return self._add_sub(other, -1.)
def __radd__(self, other):
return self._radd_sub(other, 1.)
def __rsub__(self, other):
return self._radd_sub(other, -1.)
def __mul__(self, other):
assert isinstance(other, (Number, ParameterFunctional))
if self.name != 'LincombOperator':
return LincombOperator((self,), (other,))
else:
return self.with_(coefficients=tuple(c * other for c in self.coefficients))
class ConcatenationOperator(Operator):
"""|Operator| representing the concatenation of two |Operators|.
......
......@@ -520,35 +520,50 @@ class Operator(ParametricObject):
"""
raise NotImplementedError
def __add__(self, other):
"""Sum of two operators."""
if other == 0:
return self
def _add_sub(self, other, sign):
if not isinstance(other, Operator):
return NotImplemented
from pymor.operators.constructions import LincombOperator
if isinstance(other, LincombOperator):
return NotImplemented
if self.name != 'LincombOperator' or not isinstance(self, LincombOperator):
if other.name == 'LincombOperator' and isinstance(other, LincombOperator):
operators = (self,) + other.operators
coefficients = (1.,) + (other.coefficients if sign == 1. else tuple(-c for c in other.coefficients))
else:
operators, coefficients = (self, other), (1., sign)
elif other.name == 'LincombOperator' and isinstance(other, LincombOperator):
operators = self.operators + other.operators
coefficients = self.coefficients + (other.coefficients if sign == 1.
else tuple(-c for c in other.coefficients))
else:
return LincombOperator([self, other], [1., 1.])
operators, coefficients = self.operators + (other,), self.coefficients + (sign,)
__radd__ = __add__
return LincombOperator(operators, coefficients, solver_options=self.solver_options)
def __sub__(self, other):
def _radd_sub(self, other, sign):
if other == 0:
return self
if not isinstance(other, Operator):
return NotImplemented
from pymor.operators.constructions import LincombOperator
if isinstance(other, LincombOperator):
return NotImplemented
else:
return LincombOperator([self, other], [1., -1.])
def __add__(self, other):
return self._add_sub(other, 1.)
def __sub__(self, other):
return self._add_sub(other, -1.)
def __radd__(self, other):
return self._radd_sub(other, 1.)
def __rsub__(self, other):
return self._radd_sub(other, -1.)
def __mul__(self, other):
"""Product of operator by a scalar."""
if not isinstance(other, (Number, ParameterFunctional)):
return NotImplemented
assert isinstance(other, (Number, ParameterFunctional))
from pymor.operators.constructions import LincombOperator
return LincombOperator([self], [other])
if self.name != 'LincombOperator' or not isinstance(self, LincombOperator):
return LincombOperator((self,), (other,))
else:
return self.with_(coefficients=tuple(c * other for c in self.coefficients))
def __rmul__(self, other):
return self * other
......
......@@ -45,28 +45,72 @@ class ParameterFunctional(ParametricObject):
def __call__(self, mu=None):
return self.evaluate(mu)
def __add__(self, other):
def _add_sub(self, other, sign):
if not isinstance(other, (ParameterFunctional, Number)):
return NotImplemented
if isinstance(other, Number):
if other == 0:
return self
other = ConstantParameterFunctional(other)
if isinstance(other, ParameterFunctional):
return LincombParameterFunctional([self, other], [1., 1.])
if self.name != 'LincombParameterFunctional' or not isinstance(self, LincombParameterFunctional):
if other.name == 'LincombParameterFunctional' and isinstance(other, LincombParameterFunctional):
functionals = (self,) + other.functionals
coefficients = (1.,) + (other.coefficients if sign == 1. else tuple(-c for c in other.coefficients))
else:
functionals, coefficients = (self, other), (1., sign)
elif other.name == 'LincombParameterFunctional' and isinstance(other, LincombParameterFunctional):
functionals = self.functionals + other.functionals
coefficients = self.coefficients + (other.coefficients if sign == 1.
else tuple(-c for c in other.coefficients))
else:
return NotImplemented
functionals, coefficients = self.functionals + (other,), self.coefficients + (sign,)
__radd__ = __add__
return LincombParameterFunctional(functionals, coefficients)
def __sub__(self, other):
if isinstance(other, ParameterFunctional):
return LincombParameterFunctional([self, other], [1., -1.])
def _radd_sub(self, other, sign):
assert not isinstance(other, ParameterFunctional) # handled by __add__/__sub__
if not isinstance(other, Number):
return NotImplemented
if other == 0:
return self
other = ConstantParameterFunctional(other)
if self.name != 'LincombParameterFunctional' or not isinstance(self, LincombParameterFunctional):
functionals, coefficients = (other, self), (1., sign)
else:
return self + (- other)
functionals = (other,) + self.functionals
coefficients = (1.,) + (self.coefficients if sign == 1. else tuple(-c for c in self.coefficients))
return LincombParameterFunctional(functionals, coefficients)
def __add__(self, other):
return self._add_sub(other, 1.)
def __sub__(self, other):
return self._add_sub(other, -1.)
def __radd__(self, other):
return self._radd_sub(other, 1.)
def __rsub__(self, other):
return self._radd_sub(other, -1.)
def __mul__(self, other):
if not isinstance(other, (Number, ParameterFunctional)):
return NotImplemented
return ProductParameterFunctional([self, other])
if self.name != 'ProductParameterFunctional' or not isinstance(self, ProductParameterFunctional):
if isinstance(other, ProductParameterFunctional) and other.name == 'ProductParameterFunctional':
return other.with_(factors=other.factors + [self])
else:
return ProductParameterFunctional([self, other])
elif isinstance(other, ProductParameterFunctional) and other.name == 'ProductParameterFunctional':
factors = self.factors + other.factors
return ProductParameterFunctional(factors)
else:
return self.with_(factors=self.factors + [other])
__rmul__ = __mul__
......
......@@ -41,15 +41,19 @@ def test_lincomb_function():
GenericFunction(lambda X: np.zeros(X.shape[:-1]), dim_domain=steps)):
for one in (ConstantFunction(1.0, dim_domain=steps),
GenericFunction(lambda X: np.ones(X.shape[:-1]), dim_domain=steps), 1.0):
add = (zero + one) + 0
add = (zero + one) + 1 - 1
add_ = 1 - 1 + (zero + one)
sub = (zero - one) + np.zeros(())
neg = - zero
assert np.allclose(sub(x), [-1])
assert np.allclose(add(x), [1.0])
assert np.allclose(add_(x), [1.0])
assert np.allclose(neg(x), [0.0])
(repr(add), str(add), repr(one), str(one)) # just to cover the respective special funcs too
mul = neg * 1.
mul = neg * 1. * 1.
mul_ = 1. * 1. * neg
assert np.allclose(mul(x), [0.0])
assert np.allclose(mul_(x), [0.0])
with pytest.raises(AssertionError):
zero + ConstantFunction(dim_domain=steps + 1)
with pytest.raises(AssertionError):
......
......@@ -23,17 +23,38 @@ def test_LincombParameterFunctional():
zero = pf - pf
two_pf = pf + pf
two_pf_named = two_pf.with_(name='some_interesting_quantity')
three_pf = pf + 2*pf
three_pf_ = pf + two_pf
three_pf_named = pf + two_pf_named
pf_plus_one = pf + 1
sum_ = epf + pf
pf_squared = (pf + 2*epf) * (pf - 2*epf) + 4 * epf * epf
pf_plus_one_ = 1 + pf
sum_ = epf + pf + 1 - 3
pf_squared_ = pf * pf
pf_squared_named = pf_squared_.with_(name='some_interesting_quantity')
pf_times_pf_squared = pf * pf_squared_
pf_times_pf_squared_named = pf * pf_squared_named
pf_squared = 4 * epf * epf * 1 + (pf + 2*epf) * (pf - 2*epf)
assert zero(mu) == 0
assert two_pf(mu) == 2 * pf(mu)
assert three_pf(mu) == 3 * pf(mu)
assert three_pf_(mu) == 3 * pf(mu)
assert pf_plus_one(mu) == pf(mu) + 1
assert sum_(mu) == epf(mu) + pf(mu)
assert pf_plus_one_(mu) == pf(mu) + 1
assert sum_(mu) == epf(mu) + pf(mu) + 1 - 3
assert pf_squared_(mu) == pf(mu) * pf(mu)
assert pf_squared(mu) == pf(mu) * pf(mu)
assert pf_times_pf_squared(mu) == pf(mu) * pf(mu) * pf(mu)
# derivatives are linear
assert sum_.d_mu('mu', 0)(mu) == epf.d_mu('mu', 0)(mu) + pf.d_mu('mu', 0)(mu)
assert sum_.d_mu('mu', 1)(mu) == epf.d_mu('mu', 1)(mu) + pf.d_mu('mu', 1)(mu)
assert sum_.d_mu('nu', 0)(mu) == epf.d_mu('nu', 0)(mu) + pf.d_mu('nu', 0)(mu)
# sums and products are not nested
assert len(sum_.coefficients) == 4
assert len(pf_squared.functionals[0].factors) == 4
# named functions will not be merged and still be nested
assert three_pf_named(mu) == three_pf_(mu)
assert len(three_pf_named.coefficients) != len(three_pf_.coefficients)
assert pf_times_pf_squared_named(mu) == pf_times_pf_squared(mu)
assert len(pf_times_pf_squared_named.factors) != len(pf_times_pf_squared.factors)
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment