linear_operator.py 8.97 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

Martin Reinecke's avatar
Martin Reinecke committed
19
from builtins import str
20
import abc
21
from ...nifty_meta import NiftyMeta
22

23
24
from ...field import Field
from ... import nifty_utilities as utilities
Martin Reinecke's avatar
Martin Reinecke committed
25
from future.utils import with_metaclass
26
27


28
class LinearOperator(with_metaclass(
Martin Reinecke's avatar
Martin Reinecke committed
29
        NiftyMeta, type('NewBase', (object,), {}))):
30
    """NIFTY base class for linear operators.
Theo Steininger's avatar
Theo Steininger committed
31

32
33
34
35
    The base NIFTY operator class is an abstract class from which
    other specific operator subclasses, including those preimplemented
    in NIFTY (e.g. the EndomorphicOperator, ProjectionOperator,
    DiagonalOperator, SmoothingOperator, ResponseOperator,
Theo Steininger's avatar
Theo Steininger committed
36
    PropagatorOperator, ComposedOperator) are derived.
37
38
39

    Parameters
    ----------
Theo Steininger's avatar
Theo Steininger committed
40
41
42
    default_spaces : tuple of ints *optional*
        Defines on which space(s) of a given field the Operator acts by
        default (default: None)
43
44
45

    Attributes
    ----------
Theo Steininger's avatar
Theo Steininger committed
46
47
48
49
    domain : tuple of DomainObjects, i.e. Spaces and FieldTypes
        The domain on which the Operator's input Field lives.
    target : tuple of DomainObjects, i.e. Spaces and FieldTypes
        The domain in which the Operators result lives.
50
    unitary : boolean
Theo Steininger's avatar
Theo Steininger committed
51
        Indicates whether the Operator is unitary or not.
52
53
54
55
56
57
58
59
60
61
62

    Raises
    ------
    NotImplementedError
        Raised if
            * domain is not defined
            * target is not defined
            * unitary is not set to (True/False)

    Notes
    -----
Theo Steininger's avatar
Theo Steininger committed
63
64
65
    All Operators wihtin NIFTy are linear and must therefore be a subclasses of
    the LinearOperator. A LinearOperator must have the attributes domain,
    target and unitary to be properly defined.
66
67
68

    """

69
    def __init__(self, default_spaces=None):
70
        self._default_spaces = default_spaces
71

72
73
    @staticmethod
    def _parse_domain(domain):
74
        return utilities.parse_domain(domain)
75

76
    @abc.abstractproperty
77
    def domain(self):
78
        """
Theo Steininger's avatar
Theo Steininger committed
79
80
81
        domain : tuple of DomainObjects, i.e. Spaces and FieldTypes
            The domain on which the Operator's input Field lives.
            Every Operator which inherits from the abstract LinearOperator
82
83
84
            base class must have this attribute.

        """
Theo Steininger's avatar
Theo Steininger committed
85

86
        raise NotImplementedError
87
88
89

    @abc.abstractproperty
    def target(self):
90
        """
Theo Steininger's avatar
Theo Steininger committed
91
92
93
        target : tuple of DomainObjects, i.e. Spaces and FieldTypes
            The domain on which the Operator's output Field lives.
            Every Operator which inherits from the abstract LinearOperator
94
95
96
            base class must have this attribute.

        """
Theo Steininger's avatar
Theo Steininger committed
97

98
99
        raise NotImplementedError

100
101
    @abc.abstractproperty
    def unitary(self):
102
103
104
105
106
107
108
        """
        unitary : boolean
            States whether the Operator is unitary or not.
            Every Operator which inherits from the abstract LinearOperator
            base class must have this attribute.

        """
Theo Steininger's avatar
Theo Steininger committed
109

110
111
        raise NotImplementedError

112
113
114
115
    @property
    def default_spaces(self):
        return self._default_spaces

Martin Reinecke's avatar
cleanup    
Martin Reinecke committed
116
117
    def __call__(self, x, spaces=None):
        return self.times(x, spaces)
118

Theo Steininger's avatar
Theo Steininger committed
119
    def times(self, x, spaces=None):
120
121
122
123
124
125
        """ Applies the Operator to a given Field.

        Operator and Field have to live over the same domain.

        Parameters
        ----------
Theo Steininger's avatar
Theo Steininger committed
126
127
128
129
        x : Field
            The input Field.
        spaces : tuple of ints
            Defines on which space(s) of the given Field the Operator acts.
130
131
132

        Returns
        -------
Theo Steininger's avatar
Theo Steininger committed
133
134
        out : Field
            The processed Field living on the target-domain.
135
136
137

        """

138
        spaces = self._check_input_compatibility(x, spaces)
Theo Steininger's avatar
Theo Steininger committed
139
        y = self._times(x, spaces)
140
141
        return y

Theo Steininger's avatar
Theo Steininger committed
142
    def inverse_times(self, x, spaces=None):
143
144
145
146
147
148
        """ Applies the inverse-Operator to a given Field.

        Operator and Field have to live over the same domain.

        Parameters
        ----------
Theo Steininger's avatar
Theo Steininger committed
149
150
151
152
        x : Field
            The input Field.
        spaces : tuple of ints
            Defines on which space(s) of the given Field the Operator acts.
153
154
155

        Returns
        -------
Theo Steininger's avatar
Theo Steininger committed
156
157
        out : Field
            The processed Field living on the target-domain.
158
159
160

        """

161
        spaces = self._check_input_compatibility(x, spaces, inverse=True)
162

163
        try:
Theo Steininger's avatar
Theo Steininger committed
164
            y = self._inverse_times(x, spaces)
165
166
        except(NotImplementedError):
            if (self.unitary):
Theo Steininger's avatar
Theo Steininger committed
167
                y = self._adjoint_times(x, spaces)
168
169
            else:
                raise
170
171
        return y

Theo Steininger's avatar
Theo Steininger committed
172
    def adjoint_times(self, x, spaces=None):
173
174
175
176
177
178
        """ Applies the adjoint-Operator to a given Field.

        Operator and Field have to live over the same domain.

        Parameters
        ----------
Theo Steininger's avatar
Theo Steininger committed
179
        x : Field
180
            applies the Operator to the given Field
Theo Steininger's avatar
Theo Steininger committed
181
        spaces : tuple of ints
182
183
184
185
            defines on which space of the given Field the Operator acts

        Returns
        -------
Theo Steininger's avatar
Theo Steininger committed
186
187
        out : Field
            The processed Field living on the target-domain.
188
189

        """
Theo Steininger's avatar
Theo Steininger committed
190

191
        if self.unitary:
192
            return self.inverse_times(x, spaces)
193

194
        spaces = self._check_input_compatibility(x, spaces, inverse=True)
195

196
        try:
Theo Steininger's avatar
Theo Steininger committed
197
            y = self._adjoint_times(x, spaces)
198
199
        except(NotImplementedError):
            if (self.unitary):
Theo Steininger's avatar
Theo Steininger committed
200
                y = self._inverse_times(x, spaces)
201
202
            else:
                raise
203
204
        return y

Theo Steininger's avatar
Theo Steininger committed
205
    def adjoint_inverse_times(self, x, spaces=None):
206
207
208
209
210
211
        """ Applies the adjoint-inverse Operator to a given Field.

        Operator and Field have to live over the same domain.

        Parameters
        ----------
Theo Steininger's avatar
Theo Steininger committed
212
        x : Field
213
            applies the Operator to the given Field
Theo Steininger's avatar
Theo Steininger committed
214
        spaces : tuple of ints
215
216
217
218
            defines on which space of the given Field the Operator acts

        Returns
        -------
Theo Steininger's avatar
Theo Steininger committed
219
220
        out : Field
            The processed Field living on the target-domain.
221

Theo Steininger's avatar
Theo Steininger committed
222
223
224
225
226
        Notes
        -----
        If the operator has an `inverse` then the inverse adjoint is identical
        to the adjoint inverse. We provide both names for convenience.

227
        """
228

229
        spaces = self._check_input_compatibility(x, spaces)
230

231
        try:
Theo Steininger's avatar
Theo Steininger committed
232
            y = self._adjoint_inverse_times(x, spaces)
233
234
        except(NotImplementedError):
            if self.unitary:
Theo Steininger's avatar
Theo Steininger committed
235
                y = self._times(x, spaces)
236
237
            else:
                raise
238
239
        return y

Theo Steininger's avatar
Theo Steininger committed
240
241
    def inverse_adjoint_times(self, x, spaces=None):
        return self.adjoint_inverse_times(x, spaces)
242

243
    def _times(self, x, spaces):
244
245
        raise NotImplementedError(
            "no generic instance method 'times'.")
246

247
    def _adjoint_times(self, x, spaces):
248
249
        raise NotImplementedError(
            "no generic instance method 'adjoint_times'.")
250

251
    def _inverse_times(self, x, spaces):
252
253
        raise NotImplementedError(
            "no generic instance method 'inverse_times'.")
254

255
    def _adjoint_inverse_times(self, x, spaces):
256
257
        raise NotImplementedError(
            "no generic instance method 'adjoint_inverse_times'.")
258

259
    def _check_input_compatibility(self, x, spaces, inverse=False):
260
        if not isinstance(x, Field):
Martin Reinecke's avatar
updates    
Martin Reinecke committed
261
            raise ValueError("supplied object is not a `Field`.")
262

263
264
265
266
267
        if spaces is None and self.default_spaces is not None:
            if not inverse:
                spaces = self.default_spaces
            else:
                spaces = self.default_spaces[::-1]
268

Martin Reinecke's avatar
Martin Reinecke committed
269
270
271
        # sanitize the `spaces` input
        if spaces is not None:
            spaces = utilities.cast_iseq_to_tuple(spaces)
272
273
274
275
276

        # if the operator's domain is set to something, there are two valid
        # cases:
        # 1. Case:
        #   The user specifies with `spaces` that the operators domain should
277
        #   be applied to certain spaces in the domain-tuple of x.
278
279
280
        # 2. Case:
        #   The domains of self and x match completely.

Martin Reinecke's avatar
updates    
Martin Reinecke committed
281
        self_domain = self.target if inverse else self.domain
282

283
        if spaces is None:
284
            if self_domain != x.domain:
Martin Reinecke's avatar
Martin Reinecke committed
285
286
                raise ValueError("The operator's and and field's domains "
                                 "don't match.")
287
        else:
288
289
            for i, space_index in enumerate(spaces):
                if x.domain[space_index] != self_domain[i]:
Martin Reinecke's avatar
cleanup    
Martin Reinecke committed
290
291
                    raise ValueError("The operator's and and field's domains "
                                     "don't match.")
292

293
        return spaces
294
295
296

    def __repr__(self):
        return str(self.__class__)