linear_operator.py 5.44 KB
Newer Older
1
2
3
4
5
6
7
8
9
10
11
12
# This program is free software: you can redistribute it and/or modify
# it under the terms of the GNU General Public License as published by
# the Free Software Foundation, either version 3 of the License, or
# (at your option) any later version.
#
# This program is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
# GNU General Public License for more details.
#
# You should have received a copy of the GNU General Public License
# along with this program.  If not, see <http://www.gnu.org/licenses/>.
Theo Steininger's avatar
Theo Steininger committed
13
14
15
16
17
#
# Copyright(C) 2013-2017 Max-Planck-Society
#
# NIFTy is being developed at the Max-Planck-Institut fuer Astrophysik
# and financially supported by the Studienstiftung des deutschen Volkes.
18

19
import abc
Martin Reinecke's avatar
Martin Reinecke committed
20
from ..utilities import NiftyMeta
Martin Reinecke's avatar
Martin Reinecke committed
21
from ..field import Field
Martin Reinecke's avatar
Martin Reinecke committed
22
from future.utils import with_metaclass
23
import numpy as np
24
25


26
class LinearOperator(with_metaclass(
Martin Reinecke's avatar
Martin Reinecke committed
27
        NiftyMeta, type('NewBase', (object,), {}))):
Theo Steininger's avatar
Theo Steininger committed
28

Martin Reinecke's avatar
Martin Reinecke committed
29
30
31
32
33
34
    _validMode = (False, True, True, False, True, False, False, False, True)
    _inverseMode = (0, 4, 8, 0, 1, 0, 0, 0, 2)
    _inverseCapability = (0, 4, 8, 12, 1, 5, 9, 13, 2, 6, 10, 14, 3, 7, 11, 15)
    _adjointMode = (0, 2, 1, 0, 8, 0, 0, 0, 4)
    _adjointCapability = (0, 2, 1, 3, 8, 10, 9, 11, 4, 6, 5, 7, 12, 14, 13, 15)
    _addInverse = (0, 5, 10, 15, 5, 5, 15, 15, 10, 15, 10, 15, 15, 15, 15, 15)
35

Martin Reinecke's avatar
Martin Reinecke committed
36
37
    def _dom(self, mode):
        return self.domain if (mode & 9) else self.target
38

Martin Reinecke's avatar
Martin Reinecke committed
39
40
    def _tgt(self, mode):
        return self.domain if (mode & 6) else self.target
41

42
43
    def __init__(self):
        pass
44

45
    @abc.abstractproperty
46
    def domain(self):
47
        """
Martin Reinecke's avatar
Martin Reinecke committed
48
        domain : DomainTuple
Theo Steininger's avatar
Theo Steininger committed
49
50
            The domain on which the Operator's input Field lives.
            Every Operator which inherits from the abstract LinearOperator
51
52
            base class must have this attribute.
        """
53
        raise NotImplementedError
54
55
56

    @abc.abstractproperty
    def target(self):
57
        """
Martin Reinecke's avatar
Martin Reinecke committed
58
        target : DomainTuple
Theo Steininger's avatar
Theo Steininger committed
59
60
            The domain on which the Operator's output Field lives.
            Every Operator which inherits from the abstract LinearOperator
61
62
            base class must have this attribute.
        """
63
64
        raise NotImplementedError

65
    @property
Martin Reinecke's avatar
Martin Reinecke committed
66
67
    def TIMES(self):
        return 1
68

Martin Reinecke's avatar
Martin Reinecke committed
69
70
71
    @property
    def ADJOINT_TIMES(self):
        return 2
72

Martin Reinecke's avatar
Martin Reinecke committed
73
74
75
    @property
    def INVERSE_TIMES(self):
        return 4
76

Martin Reinecke's avatar
Martin Reinecke committed
77
78
79
    @property
    def ADJOINT_INVERSE_TIMES(self):
        return 8
80

Martin Reinecke's avatar
Martin Reinecke committed
81
82
83
    @property
    def INVERSE_ADJOINT_TIMES(self):
        return 8
84

Martin Reinecke's avatar
Martin Reinecke committed
85
86
87
88
    @property
    def inverse(self):
        from .inverse_operator import InverseOperator
        return InverseOperator(self)
89

Martin Reinecke's avatar
Martin Reinecke committed
90
91
92
93
    @property
    def adjoint(self):
        from .adjoint_operator import AdjointOperator
        return AdjointOperator(self)
94

95
96
97
98
99
100
101
102
103
104
    @staticmethod
    def _toOperator(thing, dom):
        from .diagonal_operator import DiagonalOperator
        from .scaling_operator import ScalingOperator
        if isinstance(thing, LinearOperator):
            return thing
        if isinstance(thing, Field):
            return DiagonalOperator(thing)
        if np.isscalar(thing):
            return ScalingOperator(thing, dom)
Martin Reinecke's avatar
Martin Reinecke committed
105
        return NotImplemented
106

Martin Reinecke's avatar
Martin Reinecke committed
107
108
    def __mul__(self, other):
        from .chain_operator import ChainOperator
109
        other = self._toOperator(other, self.domain)
Martin Reinecke's avatar
Martin Reinecke committed
110
        return ChainOperator(self, other)
111

112
113
114
115
116
    def __rmul__(self, other):
        from .chain_operator import ChainOperator
        other = self._toOperator(other, self.target)
        return ChainOperator(other, self)

Martin Reinecke's avatar
Martin Reinecke committed
117
118
    def __add__(self, other):
        from .sum_operator import SumOperator
119
        other = self._toOperator(other, self.domain)
Martin Reinecke's avatar
Martin Reinecke committed
120
        return SumOperator(self, other)
121

122
123
124
    def __radd__(self, other):
        return self.__add__(other)

Martin Reinecke's avatar
Martin Reinecke committed
125
126
    def __sub__(self, other):
        from .sum_operator import SumOperator
127
        other = self._toOperator(other, self.domain)
Martin Reinecke's avatar
Martin Reinecke committed
128
        return SumOperator(self, other, neg=True)
129

130
131
132
133
134
    def __rsub__(self, other):
        from .sum_operator import SumOperator
        other = self._toOperator(other, self.domain)
        return SumOperator(other, self, neg=True)

Martin Reinecke's avatar
Martin Reinecke committed
135
136
    def supports(self, ops):
        return False
137

Martin Reinecke's avatar
Martin Reinecke committed
138
139
140
    @abc.abstractproperty
    def capability(self):
        raise NotImplementedError
141

Martin Reinecke's avatar
Martin Reinecke committed
142
143
144
145
146
147
    @abc.abstractmethod
    def apply(self, x, mode):
        raise NotImplementedError

    def __call__(self, x):
        return self.apply(x, self.TIMES)
148

Martin Reinecke's avatar
Martin Reinecke committed
149
150
151
152
153
154
155
156
    def times(self, x):
        return self.apply(x, self.TIMES)

    def inverse_times(self, x):
        return self.apply(x, self.INVERSE_TIMES)

    def adjoint_times(self, x):
        return self.apply(x, self.ADJOINT_TIMES)
157

Martin Reinecke's avatar
Martin Reinecke committed
158
159
    def adjoint_inverse_times(self, x):
        return self.apply(x, self.ADJOINT_INVERSE_TIMES)
160

Martin Reinecke's avatar
Martin Reinecke committed
161
162
    def inverse_adjoint_times(self, x):
        return self.apply(x, self.ADJOINT_INVERSE_TIMES)
163

Martin Reinecke's avatar
Martin Reinecke committed
164
165
166
167
168
    def _check_mode(self, mode):
        if not self._validMode[mode]:
            raise ValueError("invalid operator mode specified")
        if mode & self.capability == 0:
            raise ValueError("requested operator mode is not supported")
169

Martin Reinecke's avatar
Martin Reinecke committed
170
    def _check_input(self, x, mode):
171
        if not isinstance(x, Field):
Martin Reinecke's avatar
updates    
Martin Reinecke committed
172
            raise ValueError("supplied object is not a `Field`.")
173

Martin Reinecke's avatar
Martin Reinecke committed
174
175
176
177
        self._check_mode(mode)
        if x.domain != self._dom(mode):
                raise ValueError("The operator's and and field's domains "
                                 "don't match.")