Merge pull request #191 from piotrbartman/thrust_storage

rollback indexed_storage.py changes

PySDM/backends/numba/impl/_algorithmic_methods.py

@@ -30,7 +30,7 @@ def calculate_displacement(dim, scheme, displacement, courant, cell_origin, posi
     @staticmethod
     @numba.njit(int64(int64[:], float64[:], int64[:], int64, float64[:, :], float64[:, :], float64[:], int64[:], numba.boolean, int64, float64[:]),
                 **{**conf.JIT_FLAGS, **{'parallel': False}})
-    # TODO: waits for https://github.com/numba/numba/issues/5279
+    # TODO: reopen https://github.com/numba/numba/issues/5279 with minimal rep. ex.
     def coalescence_body(n, volume, idx, length, intensive, extensive, gamma, healthy, adaptive, subs, adaptive_memory):
         result = 1
         for i in prange(length - 1):
@@ -86,10 +86,6 @@ def compute_gamma_body(prob, rand):
         """
         for i in prange(len(prob)):
             prob[i] = np.ceil(prob[i] - rand[i // 2])
-            # TODO: same in Thrust?
-#             prob[i] *= -1.
-#             prob[i] += rand[i // 2]
-#             prob[i] = -np.floor(prob[i])
 
     @staticmethod
     def compute_gamma(prob, rand):

PySDM/backends/numba/impl/_algorithmic_step_methods.py

@@ -25,7 +25,7 @@ def amin(row, idx, length):
         return result
 
     @staticmethod
-    # @numba.njit(**conf.JIT_FLAGS)  # TODO: "np.dot() only supported on float and complex arrays"
+    # @numba.njit(**conf.JIT_FLAGS)  # Note: in Numba 0.51 "np.dot() only supported on float and complex arrays"
     def cell_id(cell_id, cell_origin, strides):
         cell_id.data[:] = np.dot(strides.data, cell_origin.data)
 

PySDM/backends/numba/impl/_maths_methods.py

@@ -35,7 +35,7 @@ def floor(output):
     @staticmethod
     @numba.njit(void(int64[:, :], float64[:, :]), **conf.JIT_FLAGS)
     def floor_out_of_place(output, input_data):
-        output[:] = np.floor(input_data)  # TODO: Try input_data//1 instead of np.floor(input_data)
+        output[:] = np.floor(input_data)
 
     @staticmethod
     @numba.njit(**{**conf.JIT_FLAGS, **{'parallel': False}})

PySDM/backends/numba/storage/indexed_storage.py

@@ -64,7 +64,8 @@ def to_ndarray(self):
         return self.data[:self.length].copy()
 
     def read_row(self, i):
-        result = IndexedStorage(self.idx, self.data[i, :], (1, *self.shape[1:]), self.dtype)
+        # TODO: shape like in ThrustRTC
+        result = IndexedStorage(self.idx, self.data[i, :], *self.shape[1:], self.dtype)
         return result
 
     def remove_zeros(self):

PySDM/backends/numba/storage/storage.py

@@ -166,6 +166,6 @@ def to_ndarray(self):
     def upload(self, data):
         np.copyto(self.data, data, casting='safe')
 
-    # TODO: optimize
+    # TODO: remove
     def write_row(self, i, row):
         self.data[i, :] = row.data

PySDM/backends/thrustRTC/impl/_maths_methods.py

@@ -2,9 +2,7 @@
 Created at 10.12.2019
 """
 
-import numpy as np
 import ThrustRTC as trtc
-import CURandRTC as rndrtc
 from ._storage_methods import StorageMethods
 from PySDM.backends.thrustRTC.nice_thrust import nice_thrust
 from PySDM.backends.thrustRTC.conf import NICE_THRUST_FLAGS
@@ -119,31 +117,3 @@ def power(output, exponent):
     @nice_thrust(**NICE_THRUST_FLAGS)
     def subtract(output, subtrahend):
         MathsMethods.__subtract_body.launch_n(output.size(), [output, subtrahend])
-        # trtc.Transform_Binary(output, subtrahend, output, trtc.Minus())
-
-    __urand_init_rng_state_body = trtc.For(['rng', 'states', 'seed'], 'i', '''
-        rng.state_init(1234, i, 0, states[i]);
-        ''')
-
-    __urand_body = trtc.For(['states', 'vec_rnd'], 'i', '''
-        vec_rnd[i]=states[i].rand01();
-        ''')
-
-    __rng = rndrtc.DVRNG()
-    states = trtc.device_vector('RNGState', 2**19)
-    __urand_init_rng_state_body.launch_n(states.size(), [__rng, states, trtc.DVInt64(12)])
-
-    @staticmethod
-    @nice_thrust(**NICE_THRUST_FLAGS)
-    def urand(data, seed=None):
-        # TODO: print("Numpy import!: ThrustRTC.urand(...)")
-
-        seed = seed or np.random.randint(2**16)
-        dseed = trtc.DVInt64(seed)
-        # MathsMethods.__urand_init_rng_state_body.launch_n(MathsMethods.states.size(), [MathsMethods.__rng, MathsMethods.states, dseed])
-        MathsMethods.__urand_body.launch_n(data.size(), [MathsMethods.states, data])
-        # hdata = data.to_host()
-        # print(np.mean(hdata))
-        # np.random.seed(seed)
-        # output = np.random.uniform(0, 1, data.shape)
-        # StorageMethods.upload(output, data)

PySDM/backends/thrustRTC/impl/_storage_methods.py

@@ -9,7 +9,6 @@
 
 
 class StorageMethods:
-    # TODO check static For
     storage = trtc.DVVector.DVVector
     integer = np.int64
     double = np.float64
@@ -27,7 +26,6 @@ def array(shape, dtype):
             raise NotImplementedError
 
         data = trtc.device_vector(elem_cls, int(np.prod(shape)))
-        # TODO: trtc.Fill(data, trtc.DVConstant(np.nan))
 
         StorageMethods.__equip(data, shape, elem_dtype)
         return data
@@ -109,19 +107,10 @@ def shuffle_global(idx, length, u01):
     @nice_thrust(**NICE_THRUST_FLAGS)
     def shuffle_local(idx, u01, cell_start):
         StorageMethods.__shuffle_local_body.launch_n(cell_start.size() - 1, [cell_start, u01, idx])
-        # TODO: print("Numba import!: ThrustRTC.shuffle_local(...)")
-        # from PySDM.backends.numba.numba import Numba
-        # host_idx = StorageMethods.to_ndarray(idx)
-        # host_u01 = StorageMethods.to_ndarray(u01)
-        # host_cell_start = StorageMethods.to_ndarray(cell_start)
-        # Numba.shuffle_local(host_idx, host_u01, host_cell_start)
-        # device_idx = StorageMethods.from_ndarray(host_idx)
-        # trtc.Copy(device_idx, idx)
 
     @staticmethod
     @nice_thrust(**NICE_THRUST_FLAGS)
     def to_ndarray(data):
-        # TODO: move to __equip??
         if isinstance(data, StorageMethods.storage):
             pass
         elif isinstance(data, trtc.DVVector.DVRange):

PySDM/dynamics/condensation/condensation.py

@@ -42,7 +42,7 @@ def register(self, builder):
         self.max_substeps = int(self.core.dt)
         self.max_substeps = int(self.core.dt)
         self.substeps = self.core.Storage.empty(self.core.mesh.n_cell, dtype=int)
-        self.substeps[:] = np.maximum(1, int(self.core.dt))  # TODO: reset substeps
+        self.substeps[:] = np.maximum(1, int(self.core.dt))  # TODO: min substep length
         self.ripening_flags = self.core.Storage.empty(self.core.mesh.n_cell, dtype=int)
         self.ripening_flags[:] = 0
 

PySDM/mesh.py

@@ -37,7 +37,7 @@ def mesh_0d(dv=None):
 
     @staticmethod
     def __strides(grid):
-        domain = np.empty(tuple(grid))  # TODO optimize
+        domain = np.empty(tuple(grid))
         strides = np.array(domain.strides).reshape(1, -1) // domain.itemsize
         return strides
 

PySDM/state/products/aerosol_specific_concentration.py

@@ -22,7 +22,7 @@ def __init__(self, radius_threshold):
     def get(self):
         self.download_moment_to_buffer('volume', rank=0,
                                        filter_range=[0, phys.volume(self.radius_threshold)])
-        result = self.buffer.copy()  # TODO !!!
+        result = self.buffer.copy()  # TODO
         self.download_to_buffer(self.core.environment['rhod'])
         result[:] /= self.core.mesh.dv
         result[:] /= self.buffer

PySDM/state/products/total_particle_specific_concentration.py

@@ -19,7 +19,7 @@ def __init__(self):
 
     def get(self):
         self.download_moment_to_buffer('volume', rank=0)
-        result = self.buffer.copy()  # TODO !!!
+        result = self.buffer.copy()  # TODO
         self.download_to_buffer(self.core.environment['rhod'])
         result[:] /= self.core.mesh.dv
         result[:] /= self.buffer

PySDM_examples/ICMW_2012_case_1/demo_setup.py

@@ -81,8 +81,6 @@ def adaptive(self):
     def condensation_coord(self):
         return self.ui_condensation_coord.value
 
-    # TODO ui_ept = Checkbox(value=Setup.enable_particle_temperatures, description="enable particle temperatures")
-
     ui_processes = [Checkbox(value=Setup.processes[key], description=key) for key in Setup.processes.keys()]
 
     @property
@@ -144,17 +142,13 @@ def mpdata_iters(self):
     def box(self):
         layout = Accordion(children=[
             VBox([self.ui_th_std0, self.ui_qv0, self.ui_p0, self.ui_kappa, self.ui_amplitude]),
-            VBox([*self.ui_processes
-                  #   , self.ui_ept  # TODO
-                  ]),
+            VBox([*self.ui_processes]),
             VBox([self.ui_nx, self.ui_nz, self.ui_sdpg, self.ui_dt, self.ui_n_steps,
                   self.ui_condensation_rtol_x, self.ui_condensation_rtol_thd,
                   self.ui_adaptive, self.ui_condensation_coord,
                   *self.ui_mpdata_options]),
-            # VBox([])  # TODO
         ])
         layout.set_title(0, 'environment parameters')
         layout.set_title(1, 'processes')
         layout.set_title(2, 'discretisation')
-        # layout.set_title(3, 'parallelisation')  # TODO
-        return layout
+        return layout

PySDM_examples/ICMW_2012_case_1/setup.py

@@ -96,8 +96,6 @@ def rhod(self, zZ):
         arg = np.power(self.p0/const.p1000, kappa) - z * kappa * const.g / self.th_std0 / phys.R(self.qv0)
         p = const.p1000 * np.power(arg, 1/kappa)
 
-        # np.testing.assert_array_less(p, Setup.p0)  # TODO: less or equal
-
         # density using "dry" potential temp.
         pd = p * (1 - self.qv0 / (self.qv0 + const.eps))
         rhod = pd / (np.power(p / const.p1000, kappa) * const.Rd * self.th_std0)

PySDM_examples/ICMW_2012_case_1/simulation.py

@@ -93,14 +93,11 @@ def reinit(self):
             builder.add_dynamic(displacement)
         if self.setup.processes["coalescence"]:
             builder.add_dynamic(Coalescence(kernel=self.setup.kernel))
-        # TODO
-        # if self.setup.processes["relaxation"]:
-        #     raise NotImplementedError()
 
         attributes = {}
         moist_environment_init(attributes, builder.core.environment,
                                spatial_discretisation=spatial_sampling.pseudorandom,
-                               spectral_discretisation=spectral_sampling.constant_multiplicity,  # TODO: random
+                               spectral_discretisation=spectral_sampling.constant_multiplicity,
                                spectrum_per_mass_of_dry_air=self.setup.spectrum_per_mass_of_dry_air,
                                r_range=(self.setup.r_min, self.setup.r_max),
                                kappa=self.setup.kappa)

PySDM_examples/Yang_et_al_2018_Fig_2/simulation.py

@@ -22,7 +22,7 @@ class Simulation:
     def __init__(self, setup):
 
         dt_output = setup.total_time / setup.n_steps  # TODO: overwritten in jupyter example
-        self.n_substeps = 1  # TODO:
+        self.n_substeps = 1  # TODO
         while (dt_output / self.n_substeps >= setup.dt_max):
             self.n_substeps += 1
         self.bins_edges = phys.volume(setup.r_bins_edges)
@@ -61,8 +61,7 @@ def __init__(self, setup):
     def save(self, output):
         cell_id = 0
         output["r_bins_values"].append(self.particles.products["Particles Wet Size Spectrum"].get())
-        volume = self.particles.state['volume']
-        volume = volume.to_ndarray()  # TODO
+        volume = self.particles.state['volume'].to_ndarray()
         output["r"].append(phys.radius(volume=volume))
         output["S"].append(self.particles.environment["RH"][cell_id] - 1)
         output["qv"].append(self.particles.environment["qv"][cell_id])

PySDM_tests/smoke_tests/Arabas_and_Shima_2017_Fig_5/test_conservation.py

@@ -45,7 +45,7 @@ def test_water_mass_conservation(setup_idx, mass_of_dry_air, scheme):
 
     # Assert
     qt = simulation.core.environment["qv"].to_ndarray() + ql(simulation)
-    np.testing.assert_approx_equal(qt, qt0, 14)  # TODO: was 15 at some point...
+    np.testing.assert_approx_equal(qt, qt0, 14)
 
 
 @pytest.mark.parametrize("setup_idx", range(len(w_avgs)))

PySDM_tests/smoke_tests/ICMW_2012_case_1/test_initialisation.py

@@ -37,26 +37,26 @@ def test_initialisation(plot=False):
     histogram_dry = np.empty((len(r_bins) - 1, n_levels))
     histogram_wet = np.empty_like(histogram_dry)
 
-    moment_0 = setup.backend.array(n_cell, dtype=int)
-    moments = setup.backend.array((n_moments, n_cell), dtype=float)
+    moment_0 = setup.backend.Storage.empty(n_cell, dtype=int)
+    moments = setup.backend.Storage.empty((n_moments, n_cell), dtype=float)
     tmp = np.empty(n_cell)
     simulation.reinit()
 
     # Act (moments)
     simulation.run()
     particles = simulation.core
     environment = simulation.core.environment
-    rhod = setup.backend.to_ndarray(environment["rhod"]).reshape(setup.grid).mean(axis=0)
+    rhod = environment["rhod"].to_ndarray().reshape(setup.grid).mean(axis=0)
 
     for i in range(len(histogram_dry)):
         particles.state.moments(
             moment_0, moments, specs={}, attr_name='dry volume', attr_range=(v_bins[i], v_bins[i + 1]))
-        particles.backend.download(moment_0, tmp)
+        moment_0.download(tmp)
         histogram_dry[i, :] = tmp.reshape(setup.grid).sum(axis=0) / (particles.mesh.dv * setup.grid[0])
 
         particles.state.moments(
             moment_0, moments, specs={}, attr_name='volume', attr_range=(v_bins[i], v_bins[i + 1]))
-        particles.backend.download(moment_0, tmp)
+        moment_0.download(tmp)
         histogram_wet[i, :] = tmp.reshape(setup.grid).sum(axis=0) / (particles.mesh.dv * setup.grid[0])
 
     # Plot

PySDM_tests/unit_tests/backends/__parametrisation__.py

@@ -3,14 +3,13 @@
 """
 
 import pytest
-import os
+from PySDM.backends import ThrustRTC
 
 from PySDM.backends.default import Default
 
 backend = Default()
-backends = []  # TODO: add Pythran
-if os.environ.get('TRAVIS') != 'true':
-    from PySDM.backends import ThrustRTC
+backends = []
+if ThrustRTC.ENABLE:
     backends.append(ThrustRTC())
 
 

PySDM_tests/unit_tests/backends/test_algorithmic_methods.py

@@ -16,6 +16,8 @@
 from .__parametrisation__ import backend
 
 
+# TODO: not implemented
+@pytest.mark.skip()
 @pytest.mark.parametrize('sut', backends)
 class TestAlgorithmicMethods:
 

PySDM_tests/unit_tests/backends/test_algorithmic_step_methods.py

@@ -14,6 +14,8 @@
 from .utils import universal_test
 
 
+# TODO: not implemented
+@pytest.mark.skip()
 @pytest.mark.parametrize('sut', backends)
 class TestAlgorithmicStepMethods:
 

PySDM_tests/unit_tests/backends/test_maths_methods.py

@@ -13,6 +13,8 @@
 from .utils import universal_test, generate_data, generate_idx
 
 
+# TODO: not implemented
+@pytest.mark.skip()
 @pytest.mark.parametrize('sut', backends)
 class TestMathsMethods:
 

PySDM_tests/unit_tests/backends/test_storage_methods.py

@@ -13,6 +13,8 @@
 from .utils import universal_test, generate_data
 
 
+# TODO: not implemented
+@pytest.mark.skip()
 @pytest.mark.parametrize('sut', backends)
 class TestStorageMethods:
 

PySDM_tests/unit_tests/dynamics/coalescence/test_sdm_single_cell.py

@@ -110,7 +110,6 @@ def test_multi_droplet(self, v, n, p):
         assert np.amin(particles.state['n'].to_ndarray()) >= 0
         assert np.sum(particles.state['n'].to_ndarray() * particles.state['volume'].to_ndarray()) == np.sum(n * v)
 
-    # TODO integration test?
     def test_multi_step(self):
         # Arrange
         n_sd = 256
@@ -137,7 +136,6 @@ def test_multi_step(self):
         desired = np.sum(n * v)
         np.testing.assert_almost_equal(actual=actual, desired=desired)
 
-    # TODO: move to backend tests
     @staticmethod
     def test_compute_gamma():
         # Arrange
@@ -160,8 +158,6 @@ def test_compute_gamma():
                 # Assert
                 assert expected(p, r) == prob_arr.to_ndarray()[0]
 
-    # TODO test_compute_probability
-
     @staticmethod
     @pytest.mark.parametrize("optimized_random", (True, False))
     def test_rnd_reuse(optimized_random):

PySDM_tests/unit_tests/dynamics/displacement/test_advection.py

@@ -19,7 +19,7 @@ def test_single_cell(self):
         sut()
 
         # Assert
-        # TODO
+        pass
 
     def test_advection(self):
         # Arrange

PySDM_tests/unit_tests/initialisation/test_spectral_discretisation.py

@@ -33,12 +33,4 @@ def test_spectral_discretisation(discretisation):
 	actual = np.sum(n)
 	desired = spectrum.cumulative(m_range[1]) - spectrum.cumulative(m_range[0])
 	quotient = actual / desired
-	# TODO relative error
 	np.testing.assert_almost_equal(actual=quotient, desired=1.0, decimal=2)
-
-
-# TODO test_linear()
-
-# TODO test_logarithmic()
-
-# TODO test_constant_multiplicity()

PySDM_tests/unit_tests/state/test_state.py

@@ -33,8 +33,7 @@ def test_housekeeping(self, volume, n):
         attributes = {'n': n, 'volume': volume}
         particles.build(attributes)
         sut = particles.state
-        # TODO
-        sut.healthy = TestState.storage([0])
+        sut.healthy = False
 
         # Act
         n_sd = sut.SD_num