Commit 8e485f86 authored by Pumpe, Daniel (dpumpe)'s avatar Pumpe, Daniel (dpumpe)
Browse files

tests DiagonalOperator

parent 2c897522
Pipeline #12239 failed with stage
in 4 minutes and 45 seconds
......@@ -239,6 +239,8 @@ class LinearOperator(Loggable, object):
raise
return y
# If the operator supports inverse() then the inverse adjoint is identical
# to the adjoint inverse. We provide both names for convenience.
def adjoint_inverse_times(self, x, spaces=None, **kwargs):
""" Applies the adjoint-inverse Operator to a given Field.
......@@ -313,6 +315,13 @@ class LinearOperator(Loggable, object):
y = self._times(x, spaces, **kwargs)
else:
raise
try:
y = self._inverse_adjoint_times(x, spaces, **kwargs)
except(NotImplementedError):
if self.unitary:
y = self._times(x, spaces, **kwargs)
else:
raise
return y
def _times(self, x, spaces):
......
......@@ -74,15 +74,14 @@ class RGSpace(Space):
zerocenter : {bool, numpy.ndarray}, *optional*
Whether the Fourier zero-mode is located in the center of the
grid (or the center of each axis speparately) or not
MR FIXME: this also does something if the space is not harmonic!
(default: False).
distances : {float, numpy.ndarray}, *optional*
Distance between two grid points along each axis
(default: None).
If distances==None:
if harmonic==True, all distances will be set to 1
if harmonic==False, the distance along each axis will be
set to the inverse of the number of points along that axis
set to the inverse of the number of points along that
axis.
harmonic : bool, *optional*
Whether the space represents a Fourier or a position grid
(default: False).
......
import unittest
import numpy as np
from numpy.testing import assert_equal,\
assert_allclose,\
assert_approx_equal
from nifty import Field,\
RGSpace,\
LMSpace,\
HPSpace,\
GLSpace,\
DiagonalOperator
from itertools import product
from test.common import expand
def _get_spaces():
rg_space = RGSpace(20)
hp_space = HPSpace(nside=6)
lm_space = LMSpace(lmax=10)
gl_space = GLSpace(nlat=6, nlon=10)
return [rg_space, hp_space, lm_space, gl_space]
class DiagonalOperator_Tests(unittest.TestCase):
spaces = _get_spaces()
@expand(product(spaces, [True, False], [True, False]))
def test_property(self, space, bare, copy):
diag = Field.from_random('normal', domain=space)
D = DiagonalOperator(space, diagonal=diag)
if D.domain[0] != space:
raise TypeError
if D.unitary != True:
raise TypeError
if D.self_adjoint != True:
raise TypeError
@expand(product(spaces, [True, False], [True, False]))
def test_times_adjoint(self, space, bare, copy):
rand1 = Field.from_random('normal', domain=space)
rand2 = Field.from_random('normal', domain=space)
diag = Field.from_random('normal', domain=space)
D = DiagonalOperator(space, diagonal=diag, bare=bare, copy=copy)
tt1 = rand1.dot(D.times(rand2))
tt2 = rand2.dot(D.times(rand1))
assert_approx_equal(tt1, tt2)
@expand(product(spaces, [True, False], [True, False]))
def test_times_inverse(self, space, bare, copy):
rand1 = Field.from_random('normal', domain=space)
rand2 = Field.from_random('normal', domain=space)
diag = Field.from_random('normal', domain=space)
D = DiagonalOperator(space, diagonal=diag, bare=bare, copy=copy)
tt1 = D.times(D.inverse_times(rand1))
assert_allclose(rand1, tt1)
tt2 = D.inverse_times(D.times(rand2))
assert_allclose(rand2, tt2)
@expand(product(spaces, [True, False], [True, False]))
def test_times(self, space, bare, copy):
rand1 = Field.from_random('normal', domain=space)
diag = Field.from_random('normal', domain=space)
D = DiagonalOperator(space, diagonal=diag, bare=bare, copy=copy)
tt = D.times(rand1)
assert_equal(tt.domain, space)
@expand(product(spaces, [True, False], [True, False]))
def test_adjoint_times(self, space, bare, copy):
rand1 = Field.from_random('normal', domain=space)
diag = Field.from_random('normal', domain=space)
D = DiagonalOperator(space, diagonal=diag, bare=bare, copy=copy)
tt = D.adjoint_times(rand1)
assert_equal(tt.domain, space)
@expand(product(spaces, [True, False], [True, False]))
def test_inverse_times(self, space, bare, copy):
rand1 = Field.from_random('normal', domain=space)
diag = Field.from_random('normal', domain=space)
D = DiagonalOperator(space, diagonal=diag, bare=bare, copy=copy)
tt = D.inverse_times(rand1)
assert_equal(tt.domain, space)
@expand(product(spaces, [True, False], [True, False]))
def test_adjoint_inverse_times(self, space, bare, copy):
rand1 = Field.from_random('normal', domain=space)
diag = Field.from_random('normal', domain=space)
D = DiagonalOperator(space, diagonal=diag, bare=bare, copy=copy)
tt = D.adjoint_inverse_times(rand1)
assert_equal(tt.domain, space)
@expand(product(spaces, [True, False], [True, False]))
def test_diagonal(self, space, bare, copy):
diag = Field.from_random('normal', domain=space)
D = DiagonalOperator(space, diagonal=diag, bare=bare, copy=copy)
diag_op = D.diagonal(bare=bare, copy=copy)
assert_allclose(diag.val.get_full_data(), diag_op.val.get_full_data())
@expand(product(spaces, [True, False], [True, False]))
def test_inverse(self, space, bare, copy):
diag = Field.from_random('normal', domain=space)
D = DiagonalOperator(space, diagonal=diag, bare=bare, copy=copy)
diag_op = D.inverse_diagonal(bare=bare)
assert_allclose(1./diag.val.get_full_data(), diag_op.val.get_full_data())
@expand(product(spaces, [True, False], [True, False]))
def test_trace(self, space, bare, copy):
diag = Field.from_random('normal', domain=space)
D = DiagonalOperator(space, diagonal=diag, bare=bare, copy=copy)
trace_op = D.trace(bare=bare)
assert_allclose(trace_op, diag.sum())
@expand(product(spaces, [True, False], [True, False]))
def test_inverse_trace(self, space, bare, copy):
diag = Field.from_random('normal', domain=space)
D = DiagonalOperator(space, diagonal=diag, bare=bare, copy=copy)
trace_op = D.inverse_trace(bare=bare)
assert_allclose(trace_op, 1./diag.sum())
@expand(product(spaces, [True, False], [True, False]))
def test_trace_log(self, space, bare, copy):
diag = Field.from_random('normal', domain=space)
D = DiagonalOperator(space, diagonal=diag, bare=bare, copy=copy)
trace_log = D.trace_log()
assert_allclose(trace_log, diag.apply_scalar_function(np.log).sum())
@expand(product(spaces, [True, False], [True, False]))
def test_determinant(self, space, bare, copy):
diag = Field.from_random('normal', domain=space)
D = DiagonalOperator(space, diagonal=diag, bare=bare, copy=copy)
det = D.determinant()
assert_allclose(det, diag.val.get_full_data.prod())
@expand(product(spaces, [True, False], [True, False]))
def test_inverse_determinant(self, space, bare, copy):
diag = Field.from_random('normal', domain=space)
D = DiagonalOperator(space, diagonal=diag, bare=bare, copy=copy)
inv_det = D.inverse_determinant()
assert_allclose(inv_det, 1./D.determinant())
@expand(product(spaces, [True, False], [True, False]))
def test_log_determinant(self, space, bare, copy):
diag = Field.from_random('normal', domain=space)
D = DiagonalOperator(space, diagonal=diag, bare=bare, copy=copy)
log_det = D.log_determinant()
assert_allclose(log_det, np.log(D.determinant()))
\ No newline at end of file
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