Commit e7bc9096 authored by Martin Reinecke's avatar Martin Reinecke

more

parent d6894d23
......@@ -825,6 +825,8 @@ class Field(Loggable, Versionable, object):
if dtype is None:
dtype = self.dtype
if x is not None:
if np.isscalar(x):
return np.full(self.shape,x, dtype=dtype)
return np.asarray(x, dtype=dtype).reshape(self.shape)
else:
return np.empty(self.shape, dtype=dtype)
......
......@@ -114,7 +114,7 @@ class LaplaceOperator(EndomorphicOperator):
dpos = self._dpos.reshape((1,)*axis + (nval-1,))
dposc = self._dposc.reshape((1,)*axis + (nval,))
deriv = (x.val[sl_r]-x.val[sl_l])/dpos # defined between points
ret = x.val.copy_empty()
ret = np.empty_like(x.val)
ret[sl_l] = deriv
ret[prefix + (-1,)] = 0.
ret[sl_r] -= deriv
......@@ -144,7 +144,7 @@ class LaplaceOperator(EndomorphicOperator):
y[prefix + (slice(None, 2),)] = 0.
y[prefix + (-1,)] = 0.
deriv = (y[sl_r]-y[sl_l])/dpos # defined between points
ret = x.val.copy_empty()
ret = np.empty_like(x.val)
ret[sl_l] = deriv
ret[prefix + (-1,)] = 0.
ret[sl_r] -= deriv
......
......@@ -4,8 +4,6 @@ from __future__ import division
from builtins import range
import numpy as np
from d2o import STRATEGIES
from .smoothing_operator import SmoothingOperator
......@@ -152,61 +150,16 @@ class DirectSmoothingOperator(SmoothingOperator):
affected_axis = affected_axes[0]
distance_array = x.domain[spaces[0]].get_distance_array(
distribution_strategy='not')
distance_array = distance_array.get_local_data(copy=False)
distance_array = x.domain[spaces[0]].get_distance_array()
#MR FIXME: this causes calls of log(0.) which should probably be avoided
if self.log_distances:
np.log(np.maximum(distance_array,1e-15), out=distance_array)
# collect the local data + ghost cells
local_data_Q = False
if x.distribution_strategy == 'not':
local_data_Q = True
elif x.distribution_strategy in STRATEGIES['slicing']:
# infer the local start/end based on the slicing information of
# x's d2o. Only gets non-trivial for axis==0.
if 0 != affected_axis:
local_data_Q = True
else:
start_index = x.val.distributor.local_start
start_distance = distance_array[start_index]
augmented_start_distance = \
(start_distance -
self.effective_smoothing_width*self.sigma)
augmented_start_index = \
np.searchsorted(distance_array, augmented_start_distance)
true_start = start_index - augmented_start_index
end_index = x.val.distributor.local_end
end_distance = distance_array[end_index-1]
augmented_end_distance = \
(end_distance + self.effective_smoothing_width*self.sigma)
augmented_end_index = \
np.searchsorted(distance_array, augmented_end_distance)
true_end = true_start + x.val.distributor.local_length
augmented_slice = slice(augmented_start_index,
augmented_end_index)
augmented_data = x.val.get_data(augmented_slice,
local_keys=True,
copy=False)
augmented_data = augmented_data.get_local_data(copy=False)
augmented_distance_array = distance_array[augmented_slice]
else:
raise ValueError("Direct smoothing not implemented for given"
"distribution strategy.")
if local_data_Q:
# if the needed data resides on the nodes already, the necessary
# are the same; no matter what the distribution strategy was.
augmented_data = x.val.get_local_data(copy=False)
augmented_distance_array = distance_array
true_start = 0
true_end = x.shape[affected_axis]
augmented_data = x.val
augmented_distance_array = distance_array
true_start = 0
true_end = x.shape[affected_axis]
# perform the convolution along the affected axes
# currently only one axis is supported
......@@ -226,6 +179,4 @@ class DirectSmoothingOperator(SmoothingOperator):
smooth_length=true_sigma,
smoothing_width=self.effective_smoothing_width)
result = x.copy_empty()
result.val.set_local_data(local_result, copy=False)
return result
return local_result
......@@ -34,26 +34,22 @@ class FFTSmoothingOperator(SmoothingOperator):
codomain = transformed_x.domain[spaces[0]]
coaxes = transformed_x.domain_axes[spaces[0]]
# create the kernel using the knowledge of codomain about itself
axes_local_distribution_strategy = \
transformed_x.val.get_axes_local_distribution_strategy(axes=coaxes)
kernel = codomain.get_distance_array(
distribution_strategy=axes_local_distribution_strategy)
kernel = codomain.get_distance_array()
#MR FIXME: this causes calls of log(0.) which should probably be avoided
if self.log_distances:
kernel.apply_scalar_function(np.log, inplace=True)
kernel=np.log(kernel)
kernel.apply_scalar_function(
codomain.get_fft_smoothing_kernel_function(self.sigma),
inplace=True)
#kernel.apply_scalar_function(
# codomain.get_fft_smoothing_kernel_function(self.sigma),
# inplace=True)
kernel = codomain.get_fft_smoothing_kernel_function(self.sigma)(kernel)
# now, apply the kernel to transformed_x
# this is done node-locally utilizing numpys reshaping in order to
# apply the kernel to the correct axes
local_transformed_x = transformed_x.val.get_local_data(copy=False)
local_kernel = kernel.get_local_data(copy=False)
local_transformed_x = transformed_x.val
local_kernel = kernel
reshaper = [local_transformed_x.shape[i] if i in coaxes else 1
for i in range(len(transformed_x.shape))]
......@@ -66,13 +62,13 @@ class FFTSmoothingOperator(SmoothingOperator):
else:
local_transformed_x *= local_kernel
transformed_x.val.set_local_data(local_transformed_x, copy=False)
transformed_x.val=local_transformed_x
smoothed_x = transformator.adjoint_times(transformed_x,
spaces=spaces)
result = x.copy_empty()
result.set_val(smoothed_x, copy=False)
result=smoothed_x
return result
......
......@@ -92,7 +92,6 @@ class SmoothingOperator_Tests(unittest.TestCase):
inp = Field.from_random(domain=sp, random_type='normal', std=1, mean=4,
dtype=tp)
out = smo(inp)
inp = inp.val.get_full_data()
assert_allclose(inp.sum(), out.sum(), rtol=tol, atol=tol)
@expand(product([10, 15], [7, 10], [1, 0.4], [2, 0.3], [0., 1., 3.7],
......@@ -104,7 +103,6 @@ class SmoothingOperator_Tests(unittest.TestCase):
inp = Field.from_random(domain=sp, random_type='normal', std=1, mean=4,
dtype=tp)
out = smo(inp)
inp = inp.val.get_full_data()
assert_allclose(inp.sum(), out.sum(), rtol=tol, atol=tol)
@expand(product([100, 200], [False, True], [0., 1., 3.7],
......@@ -117,7 +115,6 @@ class SmoothingOperator_Tests(unittest.TestCase):
inp = Field.from_random(domain=ps, random_type='normal', std=1, mean=4,
dtype=tp)
out = smo(inp)
inp = inp.val.get_full_data()
assert_allclose(inp.sum(), out.sum(), rtol=tol, atol=tol)
@expand(product([10, 15], [7, 10], [False, True], [0., 1., 3.7],
......@@ -130,5 +127,4 @@ class SmoothingOperator_Tests(unittest.TestCase):
inp = Field.from_random(domain=ps, random_type='normal', std=1, mean=4,
dtype=tp)
out = smo(inp)
inp = inp.val.get_full_data()
assert_allclose(inp.sum(), out.sum(), rtol=tol, atol=tol)
......@@ -127,4 +127,4 @@ class LMSpaceFunctionalityTests(unittest.TestCase):
@expand(get_distance_array_configs())
def test_distance_array(self, lmax, expected):
l = LMSpace(lmax)
assert_almost_equal(l.get_distance_array('not').data, expected)
assert_almost_equal(l.get_distance_array().data, expected)
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