field.py 22.4 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
#
Martin Reinecke's avatar
Martin Reinecke committed
14
# Copyright(C) 2013-2020 Max-Planck-Society
Theo Steininger's avatar
Theo Steininger committed
15
#
16
# NIFTy is being developed at the Max-Planck-Institut fuer Astrophysik.
17

Martin Reinecke's avatar
Martin Reinecke committed
18
from functools import reduce
csongor's avatar
csongor committed
19
import numpy as np
20

Martin Reinecke's avatar
stage 1    
Martin Reinecke committed
21
from . import utilities
Martin Reinecke's avatar
Martin Reinecke committed
22
from .domain_tuple import DomainTuple
Martin Reinecke's avatar
misc    
Martin Reinecke committed
23
from .operators.operator import Operator
24

25

Martin Reinecke's avatar
misc    
Martin Reinecke committed
26
class Field(Operator):
27
    """The discrete representation of a continuous field over multiple spaces.
Theo Steininger's avatar
Theo Steininger committed
28

29
    Stores data arrays and carries all the needed meta-information (i.e. the
30
    domain) for operators to be able to operate on them.
Theo Steininger's avatar
Theo Steininger committed
31

32
33
    Parameters
    ----------
34
    domain : DomainTuple
Philipp Arras's avatar
Docs    
Philipp Arras committed
35
        The domain of the new Field.
Martin Reinecke's avatar
Martin Reinecke committed
36
37
    val : numpy.ndarray
        This object's shape must match the domain shape
38
        After construction, the object will no longer be writeable!
Martin Reinecke's avatar
Martin Reinecke committed
39

Martin Reinecke's avatar
Martin Reinecke committed
40
41
    Notes
    -----
Martin Reinecke's avatar
Martin Reinecke committed
42
    If possible, do not invoke the constructor directly, but use one of the
43
    many convenience functions for instantiation!
44
    """
45

46
47
    _scalar_dom = DomainTuple.scalar_domain()

48
49
50
    def __init__(self, domain, val):
        if not isinstance(domain, DomainTuple):
            raise TypeError("domain must be of type DomainTuple")
Martin Reinecke's avatar
stage 1    
Martin Reinecke committed
51
        if not isinstance(val, np.ndarray):
Martin Reinecke's avatar
Martin Reinecke committed
52
            if np.isscalar(val):
Martin Reinecke's avatar
stage 1    
Martin Reinecke committed
53
                val = np.full(domain.shape, val)
Martin Reinecke's avatar
Martin Reinecke committed
54
            else:
Martin Reinecke's avatar
stage 1    
Martin Reinecke committed
55
                raise TypeError("val must be of type numpy.ndarray")
56
        if domain.shape != val.shape:
Martin Reinecke's avatar
Martin Reinecke committed
57
            raise ValueError("shape mismatch between val and domain")
58
59
        self._domain = domain
        self._val = val
Martin Reinecke's avatar
stage 1    
Martin Reinecke committed
60
        self._val.flags.writeable = False
Martin Reinecke's avatar
Martin Reinecke committed
61

Martin Reinecke's avatar
Martin Reinecke committed
62
63
64
65
    @staticmethod
    def scalar(val):
        return Field(Field._scalar_dom, val)

66
67
68
69
70
71
72
    # prevent implicit conversion to bool
    def __nonzero__(self):
        raise TypeError("Field does not support implicit conversion to bool")

    def __bool__(self):
        raise TypeError("Field does not support implicit conversion to bool")

73
    @staticmethod
74
    def full(domain, val):
Martin Reinecke's avatar
Martin Reinecke committed
75
76
77
78
79
        """Creates a Field with a given domain, filled with a constant value.

        Parameters
        ----------
        domain : Domain, tuple of Domain, or DomainTuple
Philipp Arras's avatar
Docs    
Philipp Arras committed
80
            Domain of the new Field.
Martin Reinecke's avatar
Martin Reinecke committed
81
        val : float/complex/int scalar
Philipp Arras's avatar
Docs    
Philipp Arras committed
82
            Fill value. Data type of the field is inferred from val.
Martin Reinecke's avatar
Martin Reinecke committed
83
84
85
86

        Returns
        -------
        Field
Philipp Arras's avatar
Docs    
Philipp Arras committed
87
            The newly created Field.
Martin Reinecke's avatar
Martin Reinecke committed
88
        """
89
90
        if not np.isscalar(val):
            raise TypeError("val must be a scalar")
91
92
93
        if not (np.isreal(val) or np.iscomplex(val)):
            raise TypeError("need arithmetic scalar")
        domain = DomainTuple.make(domain)
Martin Reinecke's avatar
Martin Reinecke committed
94
        return Field(domain, val)
95

96
    @staticmethod
Martin Reinecke's avatar
Martin Reinecke committed
97
    def from_raw(domain, arr):
Martin Reinecke's avatar
Martin Reinecke committed
98
99
100
101
102
        """Returns a Field constructed from `domain` and `arr`.

        Parameters
        ----------
        domain : DomainTuple, tuple of Domain, or Domain
Philipp Arras's avatar
Docs    
Philipp Arras committed
103
            The domain of the new Field.
Martin Reinecke's avatar
Martin Reinecke committed
104
105
106
107
        arr : numpy.ndarray
            The data content to be used for the new Field.
            Its shape must match the shape of `domain`.
        """
Martin Reinecke's avatar
stage 1    
Martin Reinecke committed
108
        return Field(DomainTuple.make(domain), arr)
Martin Reinecke's avatar
Martin Reinecke committed
109

Martin Reinecke's avatar
Martin Reinecke committed
110
    def cast_domain(self, new_domain):
Martin Reinecke's avatar
Martin Reinecke committed
111
112
113
114
115
116
117
118
119
120
121
        """Returns a field with the same data, but a different domain

        Parameters
        ----------
        new_domain : Domain, tuple of Domain, or DomainTuple
            The domain for the returned field. Must be shape-compatible to
            `self`.

        Returns
        -------
        Field
Philipp Arras's avatar
Philipp Arras committed
122
            Field defined on `new_domain`, but with the same data as `self`.
Martin Reinecke's avatar
Martin Reinecke committed
123
        """
124
        return Field(DomainTuple.make(new_domain), self._val)
125

Martin Reinecke's avatar
Martin Reinecke committed
126
127
    @staticmethod
    def from_random(random_type, domain, dtype=np.float64, **kwargs):
128
        """Draws a random field with the given parameters.
129
130
131

        Parameters
        ----------
Martin Reinecke's avatar
Martin Reinecke committed
132
133
        random_type : 'pm1', 'normal', or 'uniform'
            The random distribution to use.
Martin Reinecke's avatar
Martin Reinecke committed
134
        domain : DomainTuple
Philipp Arras's avatar
Docs    
Philipp Arras committed
135
            The domain of the output random Field.
136
        dtype : type
Philipp Arras's avatar
Docs    
Philipp Arras committed
137
            The datatype of the output random Field.
Theo Steininger's avatar
Theo Steininger committed
138

139
140
        Returns
        -------
Martin Reinecke's avatar
Martin Reinecke committed
141
        Field
Martin Reinecke's avatar
Martin Reinecke committed
142
            The newly created Field.
143
        """
Martin Reinecke's avatar
cleanup    
Martin Reinecke committed
144
        from .random import Random
Martin Reinecke's avatar
Martin Reinecke committed
145
        domain = DomainTuple.make(domain)
Martin Reinecke's avatar
cleanup    
Martin Reinecke committed
146
147
148
        generator_function = getattr(Random, random_type)
        arr = generator_function(dtype=dtype, shape=domain.shape, **kwargs)
        return Field(domain, arr)
149

theos's avatar
theos committed
150
151
    @property
    def val(self):
Martin Reinecke's avatar
Martin Reinecke committed
152
        """numpy.ndarray : the array storing the field's entries.
Martin Reinecke's avatar
Martin Reinecke committed
153

Martin Reinecke's avatar
Martin Reinecke committed
154
155
        Notes
        -----
Martin Reinecke's avatar
Martin Reinecke committed
156
        The returned array is read-only.
157
        """
Martin Reinecke's avatar
Martin Reinecke committed
158
        return self._val
csongor's avatar
csongor committed
159

Martin Reinecke's avatar
Martin Reinecke committed
160
161
162
163
164
    def val_rw(self):
        """numpy.ndarray : a copy of the array storing the field's entries.
        """
        return self._val.copy()

Martin Reinecke's avatar
Martin Reinecke committed
165
166
    @property
    def dtype(self):
Martin Reinecke's avatar
Martin Reinecke committed
167
        """type : the data type of the field's entries"""
Martin Reinecke's avatar
Martin Reinecke committed
168
169
        return self._val.dtype

Martin Reinecke's avatar
Martin Reinecke committed
170
171
    @property
    def domain(self):
Martin Reinecke's avatar
Martin Reinecke committed
172
        """DomainTuple : the field's domain"""
Martin Reinecke's avatar
Martin Reinecke committed
173
174
        return self._domain

175
176
    @property
    def shape(self):
Martin Reinecke's avatar
Martin Reinecke committed
177
        """tuple of int : the concatenated shapes of all sub-domains"""
Martin Reinecke's avatar
Martin Reinecke committed
178
        return self._domain.shape
csongor's avatar
csongor committed
179

180
    @property
Martin Reinecke's avatar
Martin Reinecke committed
181
    def size(self):
Martin Reinecke's avatar
Martin Reinecke committed
182
        """int : total number of pixels in the field"""
Martin Reinecke's avatar
Martin Reinecke committed
183
        return self._domain.size
csongor's avatar
csongor committed
184

Theo Steininger's avatar
Theo Steininger committed
185
186
    @property
    def real(self):
Martin Reinecke's avatar
Martin Reinecke committed
187
        """Field : The real part of the field"""
Martin Reinecke's avatar
Martin Reinecke committed
188
189
190
        if utilities.iscomplextype(self.dtype):
            return Field(self._domain, self._val.real)
        return self
Theo Steininger's avatar
Theo Steininger committed
191
192
193

    @property
    def imag(self):
Martin Reinecke's avatar
Martin Reinecke committed
194
        """Field : The imaginary part of the field"""
Martin Reinecke's avatar
Martin Reinecke committed
195
        if not utilities.iscomplextype(self.dtype):
196
            raise ValueError(".imag called on a non-complex Field")
Martin Reinecke's avatar
Martin Reinecke committed
197
        return Field(self._domain, self._val.imag)
Theo Steininger's avatar
Theo Steininger committed
198

199
    def scalar_weight(self, spaces=None):
Martin Reinecke's avatar
Martin Reinecke committed
200
201
202
203
204
        """Returns the uniform volume element for a sub-domain of `self`.

        Parameters
        ----------
        spaces : int, tuple of int or None
Philipp Arras's avatar
Docs    
Philipp Arras committed
205
            Indices of the sub-domains of the field's domain to be considered.
Martin Reinecke's avatar
Martin Reinecke committed
206
207
208
209
210
            If `None`, the entire domain is used.

        Returns
        -------
        float or None
Philipp Arras's avatar
Docs    
Philipp Arras committed
211
            If the requested sub-domain has a uniform volume element, it is
Martin Reinecke's avatar
Martin Reinecke committed
212
213
            returned. Otherwise, `None` is returned.
        """
214
        return self._domain.scalar_weight(spaces)
215

Martin Reinecke's avatar
Martin Reinecke committed
216
    def total_volume(self, spaces=None):
217
        """Returns the total volume of the field's domain or of a subspace of it.
Martin Reinecke's avatar
Martin Reinecke committed
218
219
220
221

        Parameters
        ----------
        spaces : int, tuple of int or None
Philipp Arras's avatar
Docs    
Philipp Arras committed
222
            Indices of the sub-domains of the field's domain to be considered.
223
            If `None`, the total volume of the whole domain is returned.
Martin Reinecke's avatar
Martin Reinecke committed
224
225
226
227

        Returns
        -------
        float
228
            the total volume of the requested (sub-)domain.
Martin Reinecke's avatar
Martin Reinecke committed
229
        """
230
        return self._domain.total_volume(spaces)
Martin Reinecke's avatar
Martin Reinecke committed
231

232
    def weight(self, power=1, spaces=None):
233
        """Weights the pixels of `self` with their invidual pixel volumes.
234
235
236
237

        Parameters
        ----------
        power : number
238
            The pixel values get multiplied with their volume-factor**power.
Theo Steininger's avatar
Theo Steininger committed
239

Martin Reinecke's avatar
Martin Reinecke committed
240
241
242
        spaces : None, int or tuple of int
            Determines on which sub-domain the operation takes place.
            If None, the entire domain is used.
Theo Steininger's avatar
Theo Steininger committed
243

244
245
        Returns
        -------
Martin Reinecke's avatar
Martin Reinecke committed
246
        Field
Theo Steininger's avatar
Theo Steininger committed
247
            The weighted field.
248
        """
Martin Reinecke's avatar
Martin Reinecke committed
249
        aout = self.val_rw()
csongor's avatar
csongor committed
250

Martin Reinecke's avatar
Martin Reinecke committed
251
        spaces = utilities.parse_spaces(spaces, len(self._domain))
csongor's avatar
csongor committed
252

253
254
        fct = 1.
        for ind in spaces:
Martin Reinecke's avatar
Martin Reinecke committed
255
            wgt = self._domain[ind].dvol
256
257
258
            if np.isscalar(wgt):
                fct *= wgt
            else:
Martin Reinecke's avatar
Martin Reinecke committed
259
                new_shape = np.ones(len(self.shape), dtype=np.int)
Martin Reinecke's avatar
Martin Reinecke committed
260
261
                new_shape[self._domain.axes[ind][0]:
                          self._domain.axes[ind][-1]+1] = wgt.shape
262
                wgt = wgt.reshape(new_shape)
263
                aout *= wgt**power
264
        fct = fct**power
Martin Reinecke's avatar
Martin Reinecke committed
265
        if fct != 1.:
266
            aout *= fct
267

Martin Reinecke's avatar
stage 1    
Martin Reinecke committed
268
        return Field(self._domain, aout)
csongor's avatar
csongor committed
269

270
    def outer(self, x):
271
        """Computes the outer product of 'self' with x.
272
273
274
275
276
277

        Parameters
        ----------
        x : Field

        Returns
Philipp Arras's avatar
Docs    
Philipp Arras committed
278
279
280
        -------
        Field
            Defined on the product space of self.domain and x.domain.
281
282
283
        """
        if not isinstance(x, Field):
            raise TypeError("The multiplier must be an instance of " +
284
                            "the Field class")
285
286
287
        from .operators.outer_product_operator import OuterProduct
        return OuterProduct(self, x.domain)(x)

Martin Reinecke's avatar
Martin Reinecke committed
288
    def vdot(self, x, spaces=None):
289
        """Computes the dot product of 'self' with x.
Theo Steininger's avatar
Theo Steininger committed
290

291
292
293
        Parameters
        ----------
        x : Field
Philipp Arras's avatar
Philipp Arras committed
294
            x must be defined on the same domain as `self`.
Theo Steininger's avatar
Theo Steininger committed
295

Philipp Arras's avatar
Philipp Arras committed
296
        spaces : None, int or tuple of int
297
298
            The dot product is only carried out over the sub-domains in this
            tuple. If None, it is carried out over all sub-domains.
Philipp Arras's avatar
Philipp Arras committed
299
            Default: None.
Theo Steininger's avatar
Theo Steininger committed
300

301
302
        Returns
        -------
Philipp Arras's avatar
Philipp Arras committed
303
        float, complex, either scalar (for full dot products) or Field (for partial dot products).
304
        """
305
        if not isinstance(x, Field):
306
            raise TypeError("The dot-partner must be an instance of " +
307
                            "the Field class")
theos's avatar
theos committed
308

Martin Reinecke's avatar
Martin Reinecke committed
309
        if x._domain != self._domain:
310
            raise ValueError("Domain mismatch")
theos's avatar
theos committed
311

Martin Reinecke's avatar
Martin Reinecke committed
312
        ndom = len(self._domain)
313
314
315
        spaces = utilities.parse_spaces(spaces, ndom)

        if len(spaces) == ndom:
Martin Reinecke's avatar
Martin Reinecke committed
316
            return Field.scalar(np.array(np.vdot(self._val, x._val)))
317
318
        # If we arrive here, we have to do a partial dot product.
        # For the moment, do this the explicit, non-optimized way
Martin Reinecke's avatar
Martin Reinecke committed
319
        return (self.conjugate()*x).sum(spaces=spaces)
theos's avatar
theos committed
320

Martin Reinecke's avatar
Martin Reinecke committed
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
    def s_vdot(self, x):
        """Computes the dot product of 'self' with x.

        Parameters
        ----------
        x : Field
            x must be defined on the same domain as `self`.

        Returns
        -------
        float or complex
            The dot product
        """
        if not isinstance(x, Field):
            raise TypeError("The dot-partner must be an instance of " +
                            "the Field class")

        if x._domain != self._domain:
            raise ValueError("Domain mismatch")

        return np.vdot(self._val, x._val)

Martin Reinecke's avatar
Martin Reinecke committed
343
    def norm(self, ord=2):
344
        """Computes the L2-norm of the field values.
csongor's avatar
csongor committed
345

Martin Reinecke's avatar
Martin Reinecke committed
346
347
        Parameters
        ----------
Philipp Arras's avatar
Philipp Arras committed
348
349
        ord : int
            Accepted values: 1, 2, ..., np.inf. Default: 2.
350
351
352
353

        Returns
        -------
        float
Martin Reinecke's avatar
Martin Reinecke committed
354
            The L2-norm of the field values.
355
        """
Martin Reinecke's avatar
stage 1    
Martin Reinecke committed
356
        return np.linalg.norm(self._val.reshape(-1), ord=ord)
357

Martin Reinecke's avatar
tweaks    
Martin Reinecke committed
358
    def conjugate(self):
359
        """Returns the complex conjugate of the field.
Theo Steininger's avatar
Theo Steininger committed
360

361
362
        Returns
        -------
Martin Reinecke's avatar
Martin Reinecke committed
363
364
        Field
            The complex conjugated field.
csongor's avatar
csongor committed
365
        """
Martin Reinecke's avatar
Martin Reinecke committed
366
        if utilities.iscomplextype(self._val.dtype):
Martin Reinecke's avatar
Martin Reinecke committed
367
368
            return Field(self._domain, self._val.conjugate())
        return self
csongor's avatar
csongor committed
369

theos's avatar
theos committed
370
    # ---General unary/contraction methods---
371

theos's avatar
theos committed
372
    def __pos__(self):
373
        return self
374

theos's avatar
theos committed
375
    def __neg__(self):
Martin Reinecke's avatar
Martin Reinecke committed
376
        return Field(self._domain, -self._val)
csongor's avatar
csongor committed
377

theos's avatar
theos committed
378
    def __abs__(self):
Martin Reinecke's avatar
Martin Reinecke committed
379
        return Field(self._domain, abs(self._val))
csongor's avatar
csongor committed
380

381
    def _contraction_helper(self, op, spaces):
theos's avatar
theos committed
382
        if spaces is None:
Martin Reinecke's avatar
Martin Reinecke committed
383
            return Field.scalar(getattr(self._val, op)())
384

Martin Reinecke's avatar
Martin Reinecke committed
385
        spaces = utilities.parse_spaces(spaces, len(self._domain))
csongor's avatar
csongor committed
386

Martin Reinecke's avatar
Martin Reinecke committed
387
        axes_list = tuple(self._domain.axes[sp_index] for sp_index in spaces)
388

Martin Reinecke's avatar
Martin Reinecke committed
389
        if len(axes_list) > 0:
theos's avatar
theos committed
390
            axes_list = reduce(lambda x, y: x+y, axes_list)
csongor's avatar
csongor committed
391

Martin Reinecke's avatar
stage1    
Martin Reinecke committed
392
        # perform the contraction on the data
Martin Reinecke's avatar
Martin Reinecke committed
393
        data = getattr(self._val, op)(axis=axes_list)
csongor's avatar
csongor committed
394

theos's avatar
theos committed
395
396
        # check if the result is scalar or if a result_field must be constr.
        if np.isscalar(data):
Martin Reinecke's avatar
Martin Reinecke committed
397
            return Field.scalar(data)
csongor's avatar
csongor committed
398
        else:
Martin Reinecke's avatar
Martin Reinecke committed
399
            return_domain = tuple(dom
Martin Reinecke's avatar
Martin Reinecke committed
400
                                  for i, dom in enumerate(self._domain)
theos's avatar
theos committed
401
                                  if i not in spaces)
402

403
            return Field(DomainTuple.make(return_domain), data)
csongor's avatar
csongor committed
404

405
    def sum(self, spaces=None):
Martin Reinecke's avatar
Martin Reinecke committed
406
407
408
409
        """Sums up over the sub-domains given by `spaces`.

        Parameters
        ----------
Philipp Arras's avatar
Philipp Arras committed
410
        spaces : None, int or tuple of int
Martin Reinecke's avatar
Martin Reinecke committed
411
            The summation is only carried out over the sub-domains in this
Martin Reinecke's avatar
cleanup    
Martin Reinecke committed
412
            tuple. If None, it is carried out over all sub-domains.
Martin Reinecke's avatar
Martin Reinecke committed
413
414
415

        Returns
        -------
Martin Reinecke's avatar
Martin Reinecke committed
416
417
        Field
            The result of the summation.
Martin Reinecke's avatar
Martin Reinecke committed
418
        """
419
        return self._contraction_helper('sum', spaces)
csongor's avatar
csongor committed
420

Martin Reinecke's avatar
Martin Reinecke committed
421
422
423
424
425
426
427
428
429
430
    def s_sum(self):
        """Returns the sum over all entries

        Returns
        -------
        scalar
            The result of the summation.
        """
        return self._val.sum()

431
    def integrate(self, spaces=None):
Martin Reinecke's avatar
Martin Reinecke committed
432
433
434
435
436
437
438
        """Integrates over the sub-domains given by `spaces`.

        Integration is performed by summing over `self` multiplied by its
        volume factors.

        Parameters
        ----------
Philipp Arras's avatar
Philipp Arras committed
439
        spaces : None, int or tuple of int
Martin Reinecke's avatar
Martin Reinecke committed
440
441
442
443
444
            The summation is only carried out over the sub-domains in this
            tuple. If None, it is carried out over all sub-domains.

        Returns
        -------
Martin Reinecke's avatar
Martin Reinecke committed
445
446
        Field
            The result of the integration.
Martin Reinecke's avatar
Martin Reinecke committed
447
        """
Martin Reinecke's avatar
Martin Reinecke committed
448
449
450
        swgt = self.scalar_weight(spaces)
        if swgt is not None:
            res = self.sum(spaces)
Martin Reinecke's avatar
fixes    
Martin Reinecke committed
451
            res = res*swgt
Martin Reinecke's avatar
Martin Reinecke committed
452
            return res
453
454
455
        tmp = self.weight(1, spaces=spaces)
        return tmp.sum(spaces)

Martin Reinecke's avatar
Martin Reinecke committed
456
457
458
459
460
461
462
463
464
465
466
467
468
469
470
471
472
    def s_integrate(self):
        """Integrates over the Field.

        Integration is performed by summing over `self` multiplied by its
        volume factors.

        Returns
        -------
        Scalar
            The result of the integration.
        """
        swgt = self.scalar_weight()
        if swgt is not None:
            return self.s_sum()*swgt
        tmp = self.weight(1)
        return tmp.s_sum()

473
    def prod(self, spaces=None):
Martin Reinecke's avatar
Martin Reinecke committed
474
475
476
477
        """Computes the product over the sub-domains given by `spaces`.

        Parameters
        ----------
Philipp Arras's avatar
Philipp Arras committed
478
        spaces : None, int or tuple of int
Martin Reinecke's avatar
Martin Reinecke committed
479
480
            The operation is only carried out over the sub-domains in this
            tuple. If None, it is carried out over all sub-domains.
Philipp Arras's avatar
Philipp Arras committed
481
            Default: None.
Martin Reinecke's avatar
Martin Reinecke committed
482
483
484

        Returns
        -------
Martin Reinecke's avatar
Martin Reinecke committed
485
        Field
Philipp Arras's avatar
Philipp Arras committed
486
            The result of the product.
Martin Reinecke's avatar
Martin Reinecke committed
487
        """
488
        return self._contraction_helper('prod', spaces)
csongor's avatar
csongor committed
489

Martin Reinecke's avatar
Martin Reinecke committed
490
491
492
    def s_prod(self):
        return self._val.prod()

493
494
    def all(self, spaces=None):
        return self._contraction_helper('all', spaces)
csongor's avatar
csongor committed
495

Martin Reinecke's avatar
Martin Reinecke committed
496
497
498
    def s_all(self):
        return self._val.all()

499
500
    def any(self, spaces=None):
        return self._contraction_helper('any', spaces)
csongor's avatar
csongor committed
501

Martin Reinecke's avatar
Martin Reinecke committed
502
503
504
    def s_any(self):
        return self._val.any()

505
506
507
508
509
510
511
512
513
514
515
#     def min(self, spaces=None):
#         """Determines the minimum over the sub-domains given by `spaces`.
#
#         Parameters
#         ----------
#         spaces : None, int or tuple of int (default: None)
#             The operation is only carried out over the sub-domains in this
#             tuple. If None, it is carried out over all sub-domains.
#
#         Returns
#         -------
Martin Reinecke's avatar
Martin Reinecke committed
516
517
#         Field
#             The result of the operation.
518
519
520
521
522
523
524
525
526
527
528
529
530
531
#         """
#         return self._contraction_helper('min', spaces)
#
#     def max(self, spaces=None):
#         """Determines the maximum over the sub-domains given by `spaces`.
#
#         Parameters
#         ----------
#         spaces : None, int or tuple of int (default: None)
#             The operation is only carried out over the sub-domains in this
#             tuple. If None, it is carried out over all sub-domains.
#
#         Returns
#         -------
Martin Reinecke's avatar
Martin Reinecke committed
532
533
#         Field
#             The result of the operation.
534
535
#         """
#         return self._contraction_helper('max', spaces)
csongor's avatar
csongor committed
536

537
    def mean(self, spaces=None):
Martin Reinecke's avatar
Martin Reinecke committed
538
539
540
541
542
543
544
        """Determines the mean over the sub-domains given by `spaces`.

        ``x.mean(spaces)`` is equivalent to
        ``x.integrate(spaces)/x.total_volume(spaces)``.

        Parameters
        ----------
Philipp Arras's avatar
Philipp Arras committed
545
        spaces : None, int or tuple of int
Martin Reinecke's avatar
Martin Reinecke committed
546
            The operation is only carried out over the sub-domains in this
Martin Reinecke's avatar
cleanup    
Martin Reinecke committed
547
            tuple. If None, it is carried out over all sub-domains.
Martin Reinecke's avatar
Martin Reinecke committed
548
549
550

        Returns
        -------
Martin Reinecke's avatar
Martin Reinecke committed
551
552
        Field
            The result of the operation.
Martin Reinecke's avatar
Martin Reinecke committed
553
        """
554
555
        if self.scalar_weight(spaces) is not None:
            return self._contraction_helper('mean', spaces)
Martin Reinecke's avatar
Martin Reinecke committed
556
        # MR FIXME: not very efficient
557
558
        # MR FIXME: do we need "spaces" here?
        tmp = self.weight(1, spaces)
Martin Reinecke's avatar
Martin Reinecke committed
559
        return tmp.sum(spaces)*(1./tmp.total_volume(spaces))
csongor's avatar
csongor committed
560

Martin Reinecke's avatar
Martin Reinecke committed
561
562
563
564
565
566
567
568
569
570
571
572
573
    def s_mean(self):
        """Determines the field mean

        ``x.s_mean()`` is equivalent to
        ``x.s_integrate()/x.total_volume()``.

        Returns
        -------
        scalar
            The result of the operation.
        """
        return self.s_integrate()/self.total_volume()

574
    def var(self, spaces=None):
Martin Reinecke's avatar
Martin Reinecke committed
575
576
577
578
        """Determines the variance over the sub-domains given by `spaces`.

        Parameters
        ----------
Philipp Arras's avatar
Philipp Arras committed
579
        spaces : None, int or tuple of int
Martin Reinecke's avatar
Martin Reinecke committed
580
581
            The operation is only carried out over the sub-domains in this
            tuple. If None, it is carried out over all sub-domains.
Philipp Arras's avatar
Philipp Arras committed
582
            Default: None.
Martin Reinecke's avatar
Martin Reinecke committed
583
584
585

        Returns
        -------
Martin Reinecke's avatar
Martin Reinecke committed
586
587
        Field
            The result of the operation.
Martin Reinecke's avatar
Martin Reinecke committed
588
        """
589
590
        if self.scalar_weight(spaces) is not None:
            return self._contraction_helper('var', spaces)
Martin Reinecke's avatar
Martin Reinecke committed
591
592
        # MR FIXME: not very efficient or accurate
        m1 = self.mean(spaces)
Martin Reinecke's avatar
fixes    
Martin Reinecke committed
593
594
595
        from .operators.contraction_operator import ContractionOperator
        op = ContractionOperator(self._domain, spaces)
        m1 = op.adjoint_times(m1)
Martin Reinecke's avatar
Martin Reinecke committed
596
        if utilities.iscomplextype(self.dtype):
597
            sq = abs(self-m1)**2
Martin Reinecke's avatar
Martin Reinecke committed
598
        else:
599
600
            sq = (self-m1)**2
        return sq.mean(spaces)
csongor's avatar
csongor committed
601

Martin Reinecke's avatar
Martin Reinecke committed
602
603
604
605
606
607
608
609
610
611
612
613
614
615
616
617
618
619
    def s_var(self):
        """Determines the field variance

        Returns
        -------
        scalar
            The result of the operation.
        """
        if self.scalar_weight() is not None:
            return self._val.var()
        # MR FIXME: not very efficient or accurate
        m1 = self.s_mean()
        if utilities.iscomplextype(self.dtype):
            sq = abs(self-m1)**2
        else:
            sq = (self-m1)**2
        return sq.s_mean()

620
    def std(self, spaces=None):
Martin Reinecke's avatar
Martin Reinecke committed
621
622
623
624
625
626
627
        """Determines the standard deviation over the sub-domains given by
        `spaces`.

        ``x.std(spaces)`` is equivalent to ``sqrt(x.var(spaces))``.

        Parameters
        ----------
Philipp Arras's avatar
Philipp Arras committed
628
        spaces : None, int or tuple of int
Martin Reinecke's avatar
Martin Reinecke committed
629
630
            The operation is only carried out over the sub-domains in this
            tuple. If None, it is carried out over all sub-domains.
Philipp Arras's avatar
Philipp Arras committed
631
            Default: None.
Martin Reinecke's avatar
Martin Reinecke committed
632
633
634

        Returns
        -------
Martin Reinecke's avatar
Martin Reinecke committed
635
636
        Field
            The result of the operation.
Martin Reinecke's avatar
Martin Reinecke committed
637
        """
638
639
        if self.scalar_weight(spaces) is not None:
            return self._contraction_helper('std', spaces)
Martin Reinecke's avatar
Martin Reinecke committed
640
        return self.var(spaces).ptw("sqrt")
csongor's avatar
csongor committed
641

Martin Reinecke's avatar
Martin Reinecke committed
642
643
644
645
646
647
648
649
650
651
652
653
654
655
    def s_std(self):
        """Determines the standard deviation of the Field.

        ``x.s_std()`` is equivalent to ``sqrt(x.s_var())``.

        Returns
        -------
        scalar
            The result of the operation.
        """
        if self.scalar_weight() is not None:
            return self._val.std()
        return np.sqrt(self.s_var())

theos's avatar
theos committed
656
    def __repr__(self):
Martin Reinecke's avatar
5->6    
Martin Reinecke committed
657
        return "<nifty6.Field>"
theos's avatar
theos committed
658
659

    def __str__(self):
Martin Reinecke's avatar
5->6    
Martin Reinecke committed
660
        return "nifty6.Field instance\n- domain      = " + \
661
               self._domain.__str__() + \
Martin Reinecke's avatar
Martin Reinecke committed
662
               "\n- val         = " + repr(self._val)
Martin Reinecke's avatar
cleanup    
Martin Reinecke committed
663

Martin Reinecke's avatar
more    
Martin Reinecke committed
664
    def extract(self, dom):
Martin Reinecke's avatar
Martin Reinecke committed
665
        if dom != self._domain:
Martin Reinecke's avatar
more    
Martin Reinecke committed
666
667
668
            raise ValueError("domain mismatch")
        return self

669
670
671
    def extract_part(self, dom):
        if dom != self._domain:
            raise ValueError("domain mismatch")
Martin Reinecke's avatar
more    
Martin Reinecke committed
672
673
674
        return self

    def unite(self, other):
Martin Reinecke's avatar
Martin Reinecke committed
675
676
677
678
        return self+other

    def flexible_addsub(self, other, neg):
        return self-other if neg else self+other
679

Martin Reinecke's avatar
Martin Reinecke committed
680
681
682
683
    def _binary_op(self, other, op):
        # if other is a field, make sure that the domains match
        f = getattr(self._val, op)
        if isinstance(other, Field):
Martin Reinecke's avatar
Martin Reinecke committed
684
            if other._domain != self._domain:
Martin Reinecke's avatar
Martin Reinecke committed
685
686
687
688
689
                raise ValueError("domains are incompatible.")
            return Field(self._domain, f(other._val))
        if np.isscalar(other):
            return Field(self._domain, f(other))
        return NotImplemented
Martin Reinecke's avatar
Martin Reinecke committed
690

Martin Reinecke's avatar
Martin Reinecke committed
691
692
693
694
    def _prep_args(self, args, kwargs):
        for arg in args + tuple(kwargs.values()):
            if not (arg is None or np.isscalar(arg) or arg.jac is None):
                raise TypeError("bad argument")
Martin Reinecke's avatar
misc    
Martin Reinecke committed
695
696
697
698
        argstmp = tuple(arg if arg is None or np.isscalar(arg) else arg._val
                        for arg in args)
        kwargstmp = {key: val if val is None or np.isscalar(val) else val._val
                     for key, val in kwargs.items()}
Martin Reinecke's avatar
Martin Reinecke committed
699
700
701
702
703
        return argstmp, kwargstmp

    def ptw(self, op, *args, **kwargs):
        from .pointwise import ptw_dict
        argstmp, kwargstmp = self._prep_args(args, kwargs)
Martin Reinecke's avatar
misc    
Martin Reinecke committed
704
705
706
707
        return Field(self._domain, ptw_dict[op][0](self._val, *argstmp, **kwargstmp))

    def ptw_with_deriv(self, op, *args, **kwargs):
        from .pointwise import ptw_dict
Martin Reinecke's avatar
Martin Reinecke committed
708
        argstmp, kwargstmp = self._prep_args(args, kwargs)
Martin Reinecke's avatar
misc    
Martin Reinecke committed
709
710
        tmp = ptw_dict[op][1](self._val, *argstmp, **kwargstmp)
        return (Field(self._domain, tmp[0]), Field(self._domain, tmp[1]))
Martin Reinecke's avatar
Martin Reinecke committed
711

Martin Reinecke's avatar
Martin Reinecke committed
712
713
714
for op in ["__add__", "__radd__",
           "__sub__", "__rsub__",
           "__mul__", "__rmul__",
715
716
717
           "__truediv__", "__rtruediv__",
           "__floordiv__", "__rfloordiv__",
           "__pow__", "__rpow__",
718
719
720
           "__lt__", "__le__", "__gt__", "__ge__", "__eq__", "__ne__"]:
    def func(op):
        def func2(self, other):
Martin Reinecke's avatar
Martin Reinecke committed
721
            return self._binary_op(other, op)
722
723
        return func2
    setattr(Field, op, func(op))
724
725
726
727
728
729
730
731
732

for op in ["__iadd__", "__isub__", "__imul__", "__idiv__",
           "__itruediv__", "__ifloordiv__", "__ipow__"]:
    def func(op):
        def func2(self, other):
            raise TypeError(
                "In-place operations are deliberately not supported")
        return func2
    setattr(Field, op, func(op))