ift issueshttps://gitlab.mpcdf.mpg.de/groups/ift/-/issues2017-09-23T07:14:41Zhttps://gitlab.mpcdf.mpg.de/ift/nifty/-/issues/192Drawing gaussian random fields over multiple domains2017-09-23T07:14:41ZPumpe, Daniel (dpumpe)Drawing gaussian random fields over multiple domainsHi,
currently I am facing some inconsistencies for drawing gaussian random fields over multiple domains. Yet I did not find the root of these inconsistencies. However in the minimal test example below, `s_1` and `s_2` do not match, nei...Hi,
currently I am facing some inconsistencies for drawing gaussian random fields over multiple domains. Yet I did not find the root of these inconsistencies. However in the minimal test example below, `s_1` and `s_2` do not match, neither mean nor variance even though they should be drawn from the same gaussian just in two different ways. Further none of the drawn fields, `s_1` `s_2`, obey the desired spectra.
What am I doing wrong or missing?
```
from nifty import *
x_0 = RGSpace(50, zerocenter=False, distances=.5)
k_0 = RGRGTransformation.get_codomain(x_0)
p_0 = PowerSpace(harmonic_partner=k_0)
x_1 = RGSpace(75, zerocenter=False, distances=.25)
k_1 = RGRGTransformation.get_codomain(x_1)
p_1 = PowerSpace(harmonic_partner=k_1)
fft_0 = FFTOperator(domain=x_0, target=k_0)
fft_1 = FFTOperator(domain=x_1, target=k_1)
p_spec_0 = Field(p_0, val=lambda l: 1/(1+l)**3)
p_spec_1 = Field(p_1, val=lambda k: 10/(1+k)**5)
def random_field_method_1():
w = Field((x_0, x_1), val=0.)
for ii in xrange(100):
outer = Field((p_0, p_1), val=np.outer(p_spec_0.val, p_spec_1.val))
sk = outer.power_synthesize()
w += fft_0.inverse_times(fft_1.inverse_times(sk, spaces=1), spaces=0)
return w/100.
def random_field_method_2():
w = Field((x_0, x_1), val=0.)
S_0 = create_power_operator(k_0, sqrt(p_spec_0))
S_1 = create_power_operator(k_1, sqrt(p_spec_1))
for ii in xrange(100):
rand = Field.from_random('normal', domain=(x_0, x_1))
rand_k = fft_0.times(fft_1.times(rand, spaces=1), spaces=0)
sk = S_0.times(S_1.times(rand_k, spaces=1), spaces=0)
w += fft_0.inverse_times(fft_1.inverse_times(sk, spaces=1), spaces=0)
return w/100.
def analyze_power(signal):
sk = fft_0.times(fft_1.times(signal,spaces=1), spaces=0)
power = sk.power_analyze()
p_drawn_0 = power.sum(spaces=1)/p_spec_1.sum()
p_drawn_1 = power.sum(spaces=0)/p_spec_0.sum()
return p_drawn_0, p_drawn_1
s_1 = random_field_method_1()
s_2 = random_field_method_2()
method_1_power_0, method_1_power_1 = analyze_power(s_1)
methdo_2_power_0, method_2_power_2 = analyze_power(s_2)
```https://gitlab.mpcdf.mpg.de/ift/nifty/-/issues/191Field.power_synthesize(): some clarifications needed2017-12-09T14:10:49ZMartin ReineckeField.power_synthesize(): some clarifications neededI'm trying to understand the intricacies of `Field.power_synthesize()` and am encountering a few points that are unclear to me:
- which combinations of `real_signal` and `real_power` are allowed? Specifically, is it allowed/sensible to ...I'm trying to understand the intricacies of `Field.power_synthesize()` and am encountering a few points that are unclear to me:
- which combinations of `real_signal` and `real_power` are allowed? Specifically, is it allowed/sensible to have `real_signal==True` and `real_power==False`?
- if `self.dtype==float`, does it make sense to have `real_power==False`? In that case, `local_rescaler.imag` will become zero.
Once I understand all of this better, I'd volunteer to extend the docstring, and (if necessary) add a few sanity checks to the code.https://gitlab.mpcdf.mpg.de/ift/nifty/-/issues/190DiagonalOperator calls nonexistent attribute adjoint()2017-09-15T23:36:07ZMartin ReineckeDiagonalOperator calls nonexistent attribute adjoint()`DiagonalOperator` contains the following code:
```
def _adjoint_times(self, x, spaces):
return self._times_helper(x, spaces,
operation=lambda z: z.adjoint().__mul__)
def _adjoint_invers...`DiagonalOperator` contains the following code:
```
def _adjoint_times(self, x, spaces):
return self._times_helper(x, spaces,
operation=lambda z: z.adjoint().__mul__)
def _adjoint_inverse_times(self, x, spaces):
return self._times_helper(x, spaces,
operation=lambda z: z.adjoint().__rtruediv__)
```
The `z.adjoint()` expression does not work, since nothing in NIFTy has an .adjoint() attribute. So far we have not noticed this, apparently because all of our `DiagonalOperator`s were self-adjoint and these methods were not called.
I don't know what the correct version is. @theos?https://gitlab.mpcdf.mpg.de/ift/nifty/-/issues/189PowerSpace volume factors2017-12-01T20:18:08ZMartin ReineckePowerSpace volume factorsAs mentioned in the NIFTy paper, there are many different ways to define volume factors for the PowerSpace.
If I remember correctly, the agreement that we reached some time ago was that `PowerSpace.weight()` should use a weight of 1, i....As mentioned in the NIFTy paper, there are many different ways to define volume factors for the PowerSpace.
If I remember correctly, the agreement that we reached some time ago was that `PowerSpace.weight()` should use a weight of 1, i.e. it should not do anything. However, the current implementation uses `rho` as weights.
I suggest to change this to 1 as discussed. The advantage of using this weight over the others is that the user does not need to do any additional corrections when they have to use a different weight; they can just apply their own volume factors without needing to undo anything the space did by default.https://gitlab.mpcdf.mpg.de/ift/nifty/-/issues/188Field.real/imag copy arrays, contrary to documentation2018-01-11T14:25:17ZMartin ReineckeField.real/imag copy arrays, contrary to documentation```
import nifty as ift
a=ift.RGSpace(10)
f=ift.Field(a,val=1.+1j)
f.real.val+=2
f.imag.val+=2
print f.val
```
This prints '1+1j' as the field value, i.e. the manipulations of real and imaginary parts did not affect the original field....```
import nifty as ift
a=ift.RGSpace(10)
f=ift.Field(a,val=1.+1j)
f.real.val+=2
f.imag.val+=2
print f.val
```
This prints '1+1j' as the field value, i.e. the manipulations of real and imaginary parts did not affect the original field. This seems to contradict the documentation, which states that the data is not copied.https://gitlab.mpcdf.mpg.de/ift/nifty/-/issues/187Field.dim returns 0 for fields without domain2017-09-28T01:37:34ZMartin ReineckeField.dim returns 0 for fields without domain```
import nifty as ift
f=ift.Field((),val=42)
print f.val
print f.dim
```
As expected, the Field contains a single value, but its `dim` property is 0, which seems incorrect and is also inconsistent with `numpy.ndarray`'s behavior for z...```
import nifty as ift
f=ift.Field((),val=42)
print f.val
print f.dim
```
As expected, the Field contains a single value, but its `dim` property is 0, which seems incorrect and is also inconsistent with `numpy.ndarray`'s behavior for zero-dimensional arrays.https://gitlab.mpcdf.mpg.de/ift/D2O/-/issues/24Missing file in last commit?2018-03-18T10:55:45ZMartin ReineckeMissing file in last commit?Is it possible that you forgot to add a file `random.py` in your last commit?
As things are, D2O now imports Python's default `random` module, but this probably doesn't address the problems with MPI-parallel seeding you mentioned.Is it possible that you forgot to add a file `random.py` in your last commit?
As things are, D2O now imports Python's default `random` module, but this probably doesn't address the problems with MPI-parallel seeding you mentioned.https://gitlab.mpcdf.mpg.de/ift/nifty/-/issues/186TransformationCache re-uses FFTs in situations where it shouldn't2017-09-28T01:59:35ZMartin ReineckeTransformationCache re-uses FFTs in situations where it shouldn't(This is mostly a reminder ... the problem should go away automatically once the `byebye_zerocenter` branch is merged.)
Currently the TransformationCache for FFTs caches and re-uses transforms even when the harmonic base has changes in ...(This is mostly a reminder ... the problem should go away automatically once the `byebye_zerocenter` branch is merged.)
Currently the TransformationCache for FFTs caches and re-uses transforms even when the harmonic base has changes in the meantime, which leads to wrong results. This can be avoided by simply not changing the harmonic base during a run, but this approach seems impractical, since we violate this constraint already in test_fft_operator.py.
Apparently, the whole test environment is not torn down and set up again between individual tests. I verified this by printing the size of the transformation cache whenever its create() method is called (line 27 in transformation_cache.py) and then ran
`nosetests test/test_operators/test_fft_operator.py -s`
I would have expected the length of the cache to be 1 all the time, but it went up beyond 100.
This is possibly also the cause for the strange behaviour in `test_power_synthesize_analyze()`.
To verify this, I suggest running the tests on maste ran on `byebye_zerocenter` and compare the convergence behaviour.https://gitlab.mpcdf.mpg.de/ift/nifty/-/issues/185Binbounds calculation in PowerSpace needs fixing2017-09-15T23:26:29ZMartin ReineckeBinbounds calculation in PowerSpace needs fixingCurrently, the initialization of a `PowerSpace` has several problems:
- when specifying `logarithmic=False`, the code behaves differently from the code before the `PowerSpace` reorganization two months ago (this is because I'm now inter...Currently, the initialization of a `PowerSpace` has several problems:
- when specifying `logarithmic=False`, the code behaves differently from the code before the `PowerSpace` reorganization two months ago (this is because I'm now interpreting the parameter `logarithmic=False` to mean "please use linear binning", instead of "please use natural binning", which was its original meaning).
- the "magic" formulae for computing nbins and binbounds are buggy; this was not really obvious before the reorganization, because they were not used in most circumstances).
- the algorithm tries to derive binbounds from insufficient user information (e.g. if the user specifies the number of bins, the code simply guesses the position of the first and last binbounds, and the guess is not necessarily intuitive).
I would strongly favor an approach where the user has to specify a set of parameters that defines the binbounds absolutely unambiguously. E.g.:
- by providing a binbounds array
- by providing a type of spacing (linear or logarithmic), a number of bins, and the first and last binbound
- by providing absolutely nothing, which leads to natural binning
There are very few cases where it could make sense for the code to guess the binbounds; the only situations that I can imagine are 1D RGSpaces and LMSpaces with linear binning. Here, there is an obvious choice for the number of bins and the placing of the binbounds.
Would it be OK if I introduced the necessary additional constraints?https://gitlab.mpcdf.mpg.de/ift/nifty/-/issues/184PowerSpace probing calls sugar.create_composed_fft_operator and dies2017-08-29T08:26:05ZPhilipp Arrasparras@mpa-garching.mpg.dePowerSpace probing calls sugar.create_composed_fft_operator and diesI narrowed the problem down to commit a09c6ed by @theos. The attached file executes until this commit and is broken afterwards.
Apparently, the prober tries to fft something although the probing takes place in a PowerSpace which does no...I narrowed the problem down to commit a09c6ed by @theos. The attached file executes until this commit and is broken afterwards.
Apparently, the prober tries to fft something although the probing takes place in a PowerSpace which does not admit for a fft.
[Critical+Filter.py](/uploads/7b0e6c37d18c6cdf932457da0ba682d1/Critical+Filter.py)https://gitlab.mpcdf.mpg.de/ift/nifty/-/issues/183Harmonize ConjugateGradient with other minimizers?2017-09-28T01:17:59ZMartin ReineckeHarmonize ConjugateGradient with other minimizers?In principle, ConjugateGradient is not fundamentally different from our descent minimizers, but its interface (including the form of the callbacks) is very different.
Wouldn't it be helpful to adjust it in a way that it also takes an `E...In principle, ConjugateGradient is not fundamentally different from our descent minimizers, but its interface (including the form of the callbacks) is very different.
Wouldn't it be helpful to adjust it in a way that it also takes an `Energy` object for minimization and passes its current energy to the callback in the same way the descent minimizers do? Then we could, for example, replace `ConjugateGradient` objects with, say, `VL_BFGS` ones and experiment much more flexibly with different minimizers.
Of course, `ConjugateGradient` would need a special subclass of `Energy` to work properly: it would, for example, need to invoke the linear operator within this energy directly. But that should be no fundamental problem ... after all, we already have some energy classes that provide `curvature` and some that don't.
I would be very interested to hear feedback!https://gitlab.mpcdf.mpg.de/ift/nifty/-/issues/182Completely remove inverse smoothing2017-09-15T23:25:14ZMartin ReineckeCompletely remove inverse smoothingThe inverse_times() of both `SmoothingOerator`s is currently broken.
`DirectSmoothingOperator` considers the inverse operation of "smoothing with sigma" to be "smoothing with 1/sigma", which is simply wrong.
`FFTSmoothingOperator` does...The inverse_times() of both `SmoothingOerator`s is currently broken.
`DirectSmoothingOperator` considers the inverse operation of "smoothing with sigma" to be "smoothing with 1/sigma", which is simply wrong.
`FFTSmoothingOperator` does in principle perform the correct operations, but does not check for division by zero, which occurs for most sigma values and indicates that the smoothing is in fact non-invertible with finite precision arithmetics.
Unfortunately the test `test_inverse_adjoint_times` does not recognize this, because the results it compares are both NaN and therefore the same...
Since this operation appears to be unused, I'd suggest to remove it ... doing it properly may be possible (at least for FFT smoothing), but is only worth the effort if someone really needs it.https://gitlab.mpcdf.mpg.de/ift/nifty/-/issues/181Logarithmic smoothing on non-PowerSpaces?2017-09-15T23:25:48ZMartin ReineckeLogarithmic smoothing on non-PowerSpaces?Currently, NIFTy allows to construct a `SmoothingOperator` like this:
```
import nifty as ift
op=ift.SmoothingOperator(ift.RGSpace(10),sigma=1.,log_distances=True)
```
I do not think this makes sense, as the RGSace is periodic and ther...Currently, NIFTy allows to construct a `SmoothingOperator` like this:
```
import nifty as ift
op=ift.SmoothingOperator(ift.RGSpace(10),sigma=1.,log_distances=True)
```
I do not think this makes sense, as the RGSace is periodic and there are no absolute distances.
This combination of parameters also leads to overflows because the code in FFTSmoothingOperator tries to compute log(0.).
Therefore I propose to ensure that `log_distance==False` whenever an FFTSmoothingOperator is constructed and raise an exception otherwise.https://gitlab.mpcdf.mpg.de/ift/nifty/-/issues/180Use MPCDF's runners for continuous integration?2017-11-13T09:39:16ZMartin ReineckeUse MPCDF's runners for continuous integration?It seems that MPCDF now offers shared runners (http://www.mpcdf.mpg.de/about-mpcdf/publications/bits-n-bytes?BB-View=196&BB-Doc=187).
Should we try to use these instead of Theo's machine?It seems that MPCDF now offers shared runners (http://www.mpcdf.mpg.de/about-mpcdf/publications/bits-n-bytes?BB-View=196&BB-Doc=187).
Should we try to use these instead of Theo's machine?https://gitlab.mpcdf.mpg.de/ift/nifty/-/issues/179Field.from_random crashes mpi2018-02-15T14:31:18ZJakob KnollmuellerField.from_random crashes mpiHi,
I found another bug with mpi:
[minimal_example.py](/uploads/3863c7e25ac0661caab65a64bd964264/minimal_example.py)
the error I get is:
```
[MainThread][ERROR ] root Uncaught exception
Traceback (most recent call last):
F...Hi,
I found another bug with mpi:
[minimal_example.py](/uploads/3863c7e25ac0661caab65a64bd964264/minimal_example.py)
the error I get is:
```
[MainThread][ERROR ] root Uncaught exception
Traceback (most recent call last):
File "minimal_example-3.py", line 16, in <module>
std=some_diagonal.val)
File "/usr/local/lib/python2.7/dist-packages/ift_nifty-3.1.0-py2.7.egg/nifty/field.py", line 251, in from_random
lambda shape: generator_function(dtype=f.dtype,
File "/usr/local/lib/python2.7/dist-packages/d2o/distributed_data_object.py", line 469, in apply_generator
self.set_local_data(generator(self.distributor.local_shape), copy=copy)
File "/usr/local/lib/python2.7/dist-packages/ift_nifty-3.1.0-py2.7.egg/nifty/field.py", line 253, in <lambda>
**random_arguments))
File "/usr/local/lib/python2.7/dist-packages/ift_nifty-3.1.0-py2.7.egg/nifty/random.py", line 50, in normal
x *= dtype.type(std)
ValueError: operands could not be broadcast together with shapes (16,64,64) (64,64,64) (16,64,64)
-------------------------------------------------------
Primary job terminated normally, but 1 process returned
a non-zero exit code.. Per user-direction, the job has been aborted.
-------------------------------------------------------
--------------------------------------------------------------------------
mpirun detected that one or more processes exited with non-zero status, thus causing
the job to be terminated. The first process to do so was:
Process name: [[29272,1],1]
Exit code: 1
--------------------------------------------------------------------------
```
I only get it if the std in from_random is a array, for numbers it seems to work
Jakobhttps://gitlab.mpcdf.mpg.de/ift/nifty/-/issues/178Replace the current implementation of MPI FFTs2018-01-11T14:25:44ZMartin ReineckeReplace the current implementation of MPI FFTsThe current way of performing MPI-parallel FFTs has several problems:
- it relies on a fork of `pyfftw` that may or may not be merged with the official version in the future, causing potential confusion for users
- it does not work for ...The current way of performing MPI-parallel FFTs has several problems:
- it relies on a fork of `pyfftw` that may or may not be merged with the official version in the future, causing potential confusion for users
- it does not work for all array sizes due to FFTW limitations.
The second point already makes broad regression testing of many different FFT sizes quite tricky.
My suggestion to overcome both these drawbacks is based on the fact that MPI communication during an FFT is _only_ needed if the first field dimension needs to be transformed, and that a multi-D FFT can be separated in to FFTs along individual axes.
If an FFT along the first axis is required, NIFTy (or D2O) could just do the following:
- MPI-tanspose the field in such a way that the first dimension is no longer distributed across CPUs
- perform the FFT along this dimension (no MPI required)
- revert the transpose again
- (perform FFTs along the other requested axes)
Internally this is exactly how FFTW handles this problem as well.
The transposition algorithm is not trivial, but certainly implementable without too many difficulties.
Doing this would also get rid of the special "fftw" distribution strategy.
Additional bonus: MPI FFTs would then also work with numpy FFT. The dependence on MPI-enabled FFTW would vanish, making NIFTy configuration simpler.https://gitlab.mpcdf.mpg.de/ift/nifty/-/issues/177Documentation for FFTSmoothingOperator and DirectSmoothingOperator missing2018-02-21T21:18:52ZPhilipp Arrasparras@mpa-garching.mpg.deDocumentation for FFTSmoothingOperator and DirectSmoothingOperator missinghttps://gitlab.mpcdf.mpg.de/ift/nifty/-/issues/176mpi crashes with power_analyze2017-08-14T22:04:56ZJakob Knollmuellermpi crashes with power_analyzeCalling power_analyze of a fftw distributed field crashes the mpirun. I attached a minimal example:
[minimal_example.py](/uploads/3253d1e327fb18e6ec1c6326779cea6f/minimal_example.py)
I run it with:
mpirun -n 2 python minimal_example.py ...Calling power_analyze of a fftw distributed field crashes the mpirun. I attached a minimal example:
[minimal_example.py](/uploads/3253d1e327fb18e6ec1c6326779cea6f/minimal_example.py)
I run it with:
mpirun -n 2 python minimal_example.py
and I get:
```
[MainThread][ERROR ] root Uncaught exception
Traceback (most recent call last):
File "minimal_example.py", line 22, in <module>
pp=sh.power_analyze()
File "/usr/local/lib/python2.7/dist-packages/ift_nifty-3.0.4-py2.7.egg/nifty/field.py", line 377, in power_analyze
for part in parts]
File "/usr/local/lib/python2.7/dist-packages/ift_nifty-3.0.4-py2.7.egg/nifty/field.py", line 413, in _single_power_analyze
axes=work_field.domain_axes[space_index])
File "/usr/local/lib/python2.7/dist-packages/ift_nifty-3.0.4-py2.7.egg/nifty/field.py", line 439, in _calculate_power_spectrum
axes=axes)
File "/usr/local/lib/python2.7/dist-packages/ift_nifty-3.0.4-py2.7.egg/nifty/field.py", line 469, in _shape_up_pindex
semiscaled_local_data = local_data.reshape(semiscaled_shape)
ValueError: cannot reshape array of size 131072 into shape (64,64,64)
-------------------------------------------------------
Primary job terminated normally, but 1 process returned
a non-zero exit code.. Per user-direction, the job has been aborted.
-------------------------------------------------------
--------------------------------------------------------------------------
mpirun detected that one or more processes exited with non-zero status, thus causing
the job to be terminated. The first process to do so was:
Process name: [[12289,1],0]
Exit code: 1
--------------------------------------------------------------------------
```https://gitlab.mpcdf.mpg.de/ift/nifty/-/issues/175Indexing bug in Laplace operator2017-08-17T23:22:06ZMartin ReineckeIndexing bug in Laplace operatorI'm running the following test case on current master:
```
from nifty import *
s_space = RGSpace([256,256])
FFT = FFTOperator(s_space)
h_space = FFT.target[0]
p_space = PowerSpace(h_space, logarithmic=True, distribution_strategy='not')...I'm running the following test case on current master:
```
from nifty import *
s_space = RGSpace([256,256])
FFT = FFTOperator(s_space)
h_space = FFT.target[0]
p_space = PowerSpace(h_space, logarithmic=True, distribution_strategy='not')
t0 = Field(p_space, val=log(0.1/(1+p_space.kindex)**3))
sop = SmoothnessOperator(p_space)
#print t0.val.get_full_data()
print t0.vdot(sop(t0))
#print t0.val.get_full_data()
print t0.vdot(sop(t0))
```
This produces the output:
```
0.645594411887
139746.578501
```
which is definitely incorrect, since the same number should be reported twice.
Things become even weirder. If I remove the commenting "#" from the first print statement, I get:
```
[ -2.30258509 -4.38202663 -4.94670585 -5.75181152 -6.49001612
-7.17883085 -8.05624052 -8.99425111 -9.93125436 -10.87947582
-11.85503513 -12.85994974 -16.24599361]
2.93684941121
139746.578501
```
So the print statement changes the result of the next statement?!
Removing the next "#" as well produces:
```
[ -2.30258509 -4.38202663 -4.94670585 -5.75181152 -6.49001612
-7.17883085 -8.05624052 -8.99425111 -9.93125436 -10.87947582
-11.85503513 -12.85994974 -16.24599361]
2.93684941121
[ -2.30258509 -4.38202663 -4.94670585 -5.75181152 -6.49001612
-7.17883085 -8.05624052 -8.99425111 -9.93125436 -10.87947582
-11.85503513 -12.85994974 -16.24599361]
2.93684941121
```
First question: can anyone reproduce this, or is it a bug in my local Python installation?
If this is reproducible, it definitely needs to be fixed quickly, but I'm completely at my wits' end ...https://gitlab.mpcdf.mpg.de/ift/nifty/-/issues/174More zerocenter weirdness...2017-08-15T00:50:33ZMartin ReineckeMore zerocenter weirdness...```
from nifty import *
s=RGSpace(3,zerocenter=True)
f=Field(s,val=[0.,1.,2.])
op=FFTOperator(s)
print "op domain:", op.domain[0]
print "op target:", op.target[0]
print "op times real field", op.times(f).val
```
This prints
```
op doma...```
from nifty import *
s=RGSpace(3,zerocenter=True)
f=Field(s,val=[0.,1.,2.])
op=FFTOperator(s)
print "op domain:", op.domain[0]
print "op target:", op.target[0]
print "op times real field", op.times(f).val
```
This prints
```
op domain: RGSpace(shape=(3,), zerocenter=(True,), distances=(0.33333333333333331,), harmonic=False)
op target: RGSpace(shape=(3,), zerocenter=(True,), distances=(1.0,), harmonic=True)
op times real field [-0.33333333+0.j -0.16666667+0.8660254j 0.16666667+0.8660254j]
```
What I don't understand: if the target domain of the operator is zero-centered, why is the first value of the result real-only, and not the middle one?