PySDM.particulator

The very class exposing PySDM.particulator.Particulator.run() method for launching simulations

  1"""
  2The very class exposing `PySDM.particulator.Particulator.run()` method for launching simulations
  3"""
  4
  5import numpy as np
  6
  7from PySDM.backends.impl_common.backend_methods import BackendMethods
  8from PySDM.backends.impl_common.freezing_attributes import (
  9    SingularAttributes,
 10    TimeDependentAttributes,
 11)
 12from PySDM.backends.impl_common.index import make_Index
 13from PySDM.backends.impl_common.indexed_storage import make_IndexedStorage
 14from PySDM.backends.impl_common.pair_indicator import make_PairIndicator
 15from PySDM.backends.impl_common.pairwise_storage import make_PairwiseStorage
 16from PySDM.impl.particle_attributes import ParticleAttributes
 17
 18
 19class Particulator:  # pylint: disable=too-many-public-methods,too-many-instance-attributes
 20    def __init__(self, n_sd, backend):
 21        assert isinstance(backend, BackendMethods)
 22        self.__n_sd = n_sd
 23
 24        self.backend = backend
 25        self.formulae = backend.formulae
 26        self.environment = None
 27        self.attributes: (ParticleAttributes, None) = None
 28        self.dynamics = {}
 29        self.products = {}
 30        self.observers = []
 31
 32        self.n_steps = 0
 33
 34        self.sorting_scheme = "default"
 35        self.condensation_solver = None
 36
 37        self.Index = make_Index(backend)  # pylint: disable=invalid-name
 38        self.PairIndicator = make_PairIndicator(backend)  # pylint: disable=invalid-name
 39        self.PairwiseStorage = make_PairwiseStorage(  # pylint: disable=invalid-name
 40            backend
 41        )
 42        self.IndexedStorage = make_IndexedStorage(  # pylint: disable=invalid-name
 43            backend
 44        )
 45
 46        self.timers = {}
 47        self.null = self.Storage.empty(0, dtype=float)
 48
 49    def run(self, steps):
 50        for _ in range(steps):
 51            for key, dynamic in self.dynamics.items():
 52                with self.timers[key]:
 53                    dynamic()
 54            self.n_steps += 1
 55            self._notify_observers()
 56
 57    def _notify_observers(self):
 58        reversed_order_so_that_environment_is_last = reversed(self.observers)
 59        for observer in reversed_order_so_that_environment_is_last:
 60            observer.notify()
 61
 62    @property
 63    def Storage(self):
 64        return self.backend.Storage
 65
 66    @property
 67    def Random(self):
 68        return self.backend.Random
 69
 70    @property
 71    def n_sd(self) -> int:
 72        return self.__n_sd
 73
 74    @property
 75    def dt(self) -> float:
 76        if self.environment is not None:
 77            return self.environment.dt
 78        return None
 79
 80    @property
 81    def mesh(self):
 82        if self.environment is not None:
 83            return self.environment.mesh
 84        return None
 85
 86    def normalize(self, prob, norm_factor):
 87        self.backend.normalize(
 88            prob=prob,
 89            cell_id=self.attributes["cell id"],
 90            cell_idx=self.attributes.cell_idx,
 91            cell_start=self.attributes.cell_start,
 92            norm_factor=norm_factor,
 93            timestep=self.dt,
 94            dv=self.mesh.dv,
 95        )
 96
 97    def update_TpRH(self):
 98        self.backend.temperature_pressure_rh(
 99            # input
100            rhod=self.environment.get_predicted("rhod"),
101            thd=self.environment.get_predicted("thd"),
102            water_vapour_mixing_ratio=self.environment.get_predicted(
103                "water_vapour_mixing_ratio"
104            ),
105            # output
106            T=self.environment.get_predicted("T"),
107            p=self.environment.get_predicted("p"),
108            RH=self.environment.get_predicted("RH"),
109        )
110
111    def condensation(self, *, rtol_x, rtol_thd, counters, RH_max, success, cell_order):
112        """Updates droplet volumes by simulating condensation driven by prior changes
113          in environment thermodynamic state, updates the environment state.
114        In the case of parcel environment, condensation is driven solely by changes in
115          the dry-air density (theta and water_vapour_mixing_ratio should not be changed
116          by other dynamics).
117        In the case of prescribed-flow/kinematic environments, the dry-air density is
118          constant in time throughout the simulation.
119        This function should only change environment's predicted `thd` and
120          `water_vapour_mixing_ratio` (and not `rhod`).
121        """
122        self.backend.condensation(
123            solver=self.condensation_solver,
124            n_cell=self.mesh.n_cell,
125            cell_start_arg=self.attributes.cell_start,
126            water_mass=self.attributes["signed water mass"],
127            multiplicity=self.attributes["multiplicity"],
128            vdry=self.attributes["dry volume"],
129            idx=self.attributes._ParticleAttributes__idx,
130            rhod=self.environment["rhod"],
131            thd=self.environment["thd"],
132            water_vapour_mixing_ratio=self.environment["water_vapour_mixing_ratio"],
133            dv=self.environment.dv,
134            prhod=self.environment.get_predicted("rhod"),
135            pthd=self.environment.get_predicted("thd"),
136            predicted_water_vapour_mixing_ratio=self.environment.get_predicted(
137                "water_vapour_mixing_ratio"
138            ),
139            kappa=self.attributes["kappa"],
140            f_org=self.attributes["dry volume organic fraction"],
141            rtol_x=rtol_x,
142            rtol_thd=rtol_thd,
143            v_cr=self.attributes["critical volume"],
144            timestep=self.dt,
145            counters=counters,
146            cell_order=cell_order,
147            RH_max=RH_max,
148            success=success,
149            cell_id=self.attributes["cell id"],
150            reynolds_number=self.attributes["Reynolds number"],
151            air_density=self.environment["air density"],
152            air_dynamic_viscosity=self.environment["air dynamic viscosity"],
153        )
154        self.attributes.mark_updated("signed water mass")
155
156    def collision_coalescence_breakup(
157        self,
158        *,
159        enable_breakup,
160        gamma,
161        rand,
162        Ec,
163        Eb,
164        fragment_mass,
165        coalescence_rate,
166        breakup_rate,
167        breakup_rate_deficit,
168        is_first_in_pair,
169        warn_overflows,
170        max_multiplicity,
171    ):
172        # pylint: disable=too-many-locals
173        idx = self.attributes._ParticleAttributes__idx
174        healthy = self.attributes._ParticleAttributes__healthy_memory
175        cell_id = self.attributes["cell id"]
176        multiplicity = self.attributes["multiplicity"]
177        attributes = self.attributes.get_extensive_attribute_storage()
178        if enable_breakup:
179            self.backend.collision_coalescence_breakup(
180                multiplicity=multiplicity,
181                idx=idx,
182                attributes=attributes,
183                gamma=gamma,
184                rand=rand,
185                Ec=Ec,
186                Eb=Eb,
187                fragment_mass=fragment_mass,
188                healthy=healthy,
189                cell_id=cell_id,
190                coalescence_rate=coalescence_rate,
191                breakup_rate=breakup_rate,
192                breakup_rate_deficit=breakup_rate_deficit,
193                is_first_in_pair=is_first_in_pair,
194                warn_overflows=warn_overflows,
195                particle_mass=self.attributes["water mass"],
196                max_multiplicity=max_multiplicity,
197            )
198        else:
199            self.backend.collision_coalescence(
200                multiplicity=multiplicity,
201                idx=idx,
202                attributes=attributes,
203                gamma=gamma,
204                healthy=healthy,
205                cell_id=cell_id,
206                coalescence_rate=coalescence_rate,
207                is_first_in_pair=is_first_in_pair,
208            )
209        self.attributes.sanitize()
210        self.attributes.mark_updated("multiplicity")
211        for key in self.attributes.get_extensive_attribute_keys():
212            self.attributes.mark_updated(key)
213
214    def oxidation(
215        self,
216        *,
217        kinetic_consts,
218        timestep,
219        equilibrium_consts,
220        dissociation_factors,
221        do_chemistry_flag,
222    ):
223        self.backend.oxidation(
224            n_sd=self.n_sd,
225            cell_ids=self.attributes["cell id"],
226            do_chemistry_flag=do_chemistry_flag,
227            k0=kinetic_consts["k0"],
228            k1=kinetic_consts["k1"],
229            k2=kinetic_consts["k2"],
230            k3=kinetic_consts["k3"],
231            K_SO2=equilibrium_consts["K_SO2"],
232            K_HSO3=equilibrium_consts["K_HSO3"],
233            dissociation_factor_SO2=dissociation_factors["SO2"],
234            timestep=timestep,
235            # input
236            droplet_volume=self.attributes["volume"],
237            pH=self.attributes["pH"],
238            # output
239            moles_O3=self.attributes["moles_O3"],
240            moles_H2O2=self.attributes["moles_H2O2"],
241            moles_S_IV=self.attributes["moles_S_IV"],
242            moles_S_VI=self.attributes["moles_S_VI"],
243        )
244        for attr in ("moles_S_IV", "moles_S_VI", "moles_H2O2", "moles_O3"):
245            self.attributes.mark_updated(attr)
246
247    def dissolution(
248        self,
249        *,
250        gaseous_compounds,
251        system_type,
252        dissociation_factors,
253        timestep,
254        environment_mixing_ratios,
255        do_chemistry_flag,
256    ):
257        self.backend.dissolution(
258            n_cell=self.mesh.n_cell,
259            n_threads=1,
260            cell_order=np.arange(self.mesh.n_cell),
261            cell_start_arg=self.attributes.cell_start,
262            idx=self.attributes._ParticleAttributes__idx,
263            do_chemistry_flag=do_chemistry_flag,
264            mole_amounts={
265                key: self.attributes["moles_" + key] for key in gaseous_compounds.keys()
266            },
267            env_mixing_ratio=environment_mixing_ratios,
268            # note: assuming condensation was called
269            env_p=self.environment.get_predicted("p"),
270            env_T=self.environment.get_predicted("T"),
271            env_rho_d=self.environment.get_predicted("rhod"),
272            timestep=timestep,
273            dv=self.mesh.dv,
274            droplet_volume=self.attributes["volume"],
275            multiplicity=self.attributes["multiplicity"],
276            system_type=system_type,
277            dissociation_factors=dissociation_factors,
278        )
279        for key in gaseous_compounds.keys():
280            self.attributes.mark_updated(f"moles_{key}")
281
282    def chem_recalculate_cell_data(self, *, equilibrium_consts, kinetic_consts):
283        self.backend.chem_recalculate_cell_data(
284            equilibrium_consts=equilibrium_consts,
285            kinetic_consts=kinetic_consts,
286            temperature=self.environment.get_predicted("T"),
287        )
288
289    def chem_recalculate_drop_data(self, *, dissociation_factors, equilibrium_consts):
290        self.backend.chem_recalculate_drop_data(
291            dissociation_factors=dissociation_factors,
292            equilibrium_consts=equilibrium_consts,
293            pH=self.attributes["pH"],
294            cell_id=self.attributes["cell id"],
295        )
296
297    def recalculate_cell_id(self):
298        if not self.attributes.has_attribute("cell origin"):
299            return
300        self.backend.cell_id(
301            self.attributes["cell id"],
302            self.attributes["cell origin"],
303            self.backend.Storage.from_ndarray(self.environment.mesh.strides),
304        )
305        self.attributes._ParticleAttributes__sorted = False
306
307    def sort_within_pair_by_attr(self, is_first_in_pair, attr_name):
308        self.backend.sort_within_pair_by_attr(
309            self.attributes._ParticleAttributes__idx,
310            is_first_in_pair,
311            self.attributes[attr_name],
312        )
313
314    def moments(
315        self,
316        *,
317        moment_0,
318        moments,
319        specs: dict,
320        attr_name="signed water mass",
321        attr_range=(-np.inf, np.inf),
322        weighting_attribute="water mass",
323        weighting_rank=0,
324        skip_division_by_m0=False,
325    ):
326        """
327        Writes to `moment_0` and `moment` the zero-th and the k-th statistical moments
328        of particle attributes computed filtering by value of the attribute `attr_name`
329        to fall within `attr_range`. The moment ranks are defined by `specs`.
330
331        Parameters:
332            specs: e.g., `specs={'volume': (1,2,3), 'kappa': (1)}` computes three moments
333                of volume and one moment of kappa
334            skip_division_by_m0: if set to `True`, the values written to `moments` are
335                multiplied by the 0-th moment (e.g., total volume instead of mean volume)
336        """
337        if len(specs) == 0:
338            raise ValueError("empty specs passed")
339        attr_data, ranks = [], []
340        for attr in specs:
341            for rank in specs[attr]:
342                attr_data.append(self.attributes[attr])
343                ranks.append(rank)
344        assert len(set(attr_data)) <= 1
345        if len(attr_data) == 0:
346            attr_data = self.backend.Storage.empty((0,), dtype=float)
347        else:
348            attr_data = attr_data[0]
349
350        ranks = self.backend.Storage.from_ndarray(np.array(ranks, dtype=float))
351
352        self.backend.moments(
353            moment_0=moment_0,
354            moments=moments,
355            multiplicity=self.attributes["multiplicity"],
356            attr_data=attr_data,
357            cell_id=self.attributes["cell id"],
358            idx=self.attributes._ParticleAttributes__idx,
359            length=self.attributes.super_droplet_count,
360            ranks=ranks,
361            min_x=attr_range[0],
362            max_x=attr_range[1],
363            x_attr=self.attributes[attr_name],
364            weighting_attribute=self.attributes[weighting_attribute],
365            weighting_rank=weighting_rank,
366            skip_division_by_m0=skip_division_by_m0,
367        )
368
369    def spectrum_moments(
370        self,
371        *,
372        moment_0,
373        moments,
374        attr,
375        rank,
376        attr_bins,
377        attr_name="water mass",
378        weighting_attribute="water mass",
379        weighting_rank=0,
380    ):
381        attr_data = self.attributes[attr]
382        self.backend.spectrum_moments(
383            moment_0=moment_0,
384            moments=moments,
385            multiplicity=self.attributes["multiplicity"],
386            attr_data=attr_data,
387            cell_id=self.attributes["cell id"],
388            idx=self.attributes._ParticleAttributes__idx,
389            length=self.attributes.super_droplet_count,
390            rank=rank,
391            x_bins=attr_bins,
392            x_attr=self.attributes[attr_name],
393            weighting_attribute=self.attributes[weighting_attribute],
394            weighting_rank=weighting_rank,
395        )
396
397    def adaptive_sdm_end(self, dt_left):
398        return self.backend.adaptive_sdm_end(dt_left, self.attributes.cell_start)
399
400    def remove_precipitated(
401        self, *, displacement, precipitation_counting_level_index
402    ) -> float:
403        rainfall_mass = self.backend.flag_precipitated(
404            cell_origin=self.attributes["cell origin"],
405            position_in_cell=self.attributes["position in cell"],
406            water_mass=self.attributes["water mass"],
407            multiplicity=self.attributes["multiplicity"],
408            idx=self.attributes._ParticleAttributes__idx,
409            length=self.attributes.super_droplet_count,
410            healthy=self.attributes._ParticleAttributes__healthy_memory,
411            precipitation_counting_level_index=precipitation_counting_level_index,
412            displacement=displacement,
413        )
414        self.attributes.sanitize()
415        return rainfall_mass
416
417    def flag_out_of_column(self):
418        self.backend.flag_out_of_column(
419            cell_origin=self.attributes["cell origin"],
420            position_in_cell=self.attributes["position in cell"],
421            idx=self.attributes._ParticleAttributes__idx,
422            length=self.attributes.super_droplet_count,
423            healthy=self.attributes._ParticleAttributes__healthy_memory,
424            domain_top_level_index=self.mesh.grid[-1],
425        )
426        self.attributes.sanitize()
427
428    def calculate_displacement(
429        self, *, displacement, courant, cell_origin, position_in_cell, n_substeps
430    ):
431        for dim in range(len(self.environment.mesh.grid)):
432            self.backend.calculate_displacement(
433                dim=dim,
434                displacement=displacement,
435                courant=courant[dim],
436                cell_origin=cell_origin,
437                position_in_cell=position_in_cell,
438                n_substeps=n_substeps,
439            )
440
441    def isotopic_fractionation(self, heavy_isotopes: tuple):
442        self.backend.isotopic_fractionation()
443        for isotope in heavy_isotopes:
444            self.attributes.mark_updated(f"moles_{isotope}")
445
446    def seeding(
447        self,
448        *,
449        seeded_particle_index,
450        seeded_particle_multiplicity,
451        seeded_particle_extensive_attributes,
452        number_of_super_particles_to_inject,
453    ):
454        n_null = self.n_sd - self.attributes.super_droplet_count
455        if n_null == 0:
456            raise ValueError(
457                "No available seeds to inject. Please provide particles with nan filled attributes."
458            )
459
460        if number_of_super_particles_to_inject > n_null:
461            raise ValueError(
462                "Trying to inject more super particles than space available."
463            )
464
465        if number_of_super_particles_to_inject > len(seeded_particle_multiplicity):
466            raise ValueError(
467                "Trying to inject multiple super particles with the same attributes. \
468                Instead increase multiplicity of injected particles."
469            )
470
471        self.backend.seeding(
472            idx=self.attributes._ParticleAttributes__idx,
473            multiplicity=self.attributes["multiplicity"],
474            extensive_attributes=self.attributes.get_extensive_attribute_storage(),
475            seeded_particle_index=seeded_particle_index,
476            seeded_particle_multiplicity=seeded_particle_multiplicity,
477            seeded_particle_extensive_attributes=seeded_particle_extensive_attributes,
478            number_of_super_particles_to_inject=number_of_super_particles_to_inject,
479        )
480        self.attributes.reset_idx()
481        self.attributes.sanitize()
482
483        self.attributes.mark_updated("multiplicity")
484        for key in self.attributes.get_extensive_attribute_keys():
485            self.attributes.mark_updated(key)
486
487    def deposition(self):
488        self.backend.deposition(
489            multiplicity=self.attributes["multiplicity"],
490            signed_water_mass=self.attributes["signed water mass"],
491            current_temperature=self.environment["T"],
492            current_total_pressure=self.environment["p"],
493            current_relative_humidity=self.environment["RH"],
494            current_water_activity=self.environment["a_w_ice"],
495            current_vapour_mixing_ratio=self.environment["water_vapour_mixing_ratio"],
496            current_dry_air_density=self.environment["rhod"],
497            current_dry_potential_temperature=self.environment["thd"],
498            cell_volume=self.environment.mesh.dv,
499            time_step=self.dt,
500            cell_id=self.attributes["cell id"],
501            reynolds_number=self.attributes["Reynolds number"],
502            schmidt_number=self.environment["Schmidt number"],
503            predicted_vapour_mixing_ratio=self.environment.get_predicted(
504                "water_vapour_mixing_ratio"
505            ),
506            predicted_dry_potential_temperature=self.environment.get_predicted("thd"),
507        )
508        self.attributes.mark_updated("signed water mass")
509        # TODO #1524 - should we update here?
510        # self.update_TpRH(only_if_not_last='VapourDepositionOnIce')
511
512    def immersion_freezing_time_dependent(self, *, thaw: bool, rand: Storage):
513        self.backend.freeze_time_dependent(
514            rand=rand,
515            attributes=TimeDependentAttributes(
516                immersed_surface_area=self.attributes["immersed surface area"],
517                signed_water_mass=self.attributes["signed water mass"],
518            ),
519            timestep=self.dt,
520            cell=self.attributes["cell id"],
521            a_w_ice=self.environment["a_w_ice"],
522            temperature=self.environment["T"],
523            relative_humidity=self.environment["RH"],
524            thaw=thaw,
525        )
526        self.attributes.mark_updated("signed water mass")
527
528    def immersion_freezing_singular(self, *, thaw: bool):
529        self.backend.freeze_singular(
530            attributes=SingularAttributes(
531                freezing_temperature=self.attributes["freezing temperature"],
532                signed_water_mass=self.attributes["signed water mass"],
533            ),
534            temperature=self.environment["T"],
535            relative_humidity=self.environment["RH"],
536            cell=self.attributes["cell id"],
537            thaw=thaw,
538        )
539        self.attributes.mark_updated("signed water mass")
class Particulator:
 20class Particulator:  # pylint: disable=too-many-public-methods,too-many-instance-attributes
 21    def __init__(self, n_sd, backend):
 22        assert isinstance(backend, BackendMethods)
 23        self.__n_sd = n_sd
 24
 25        self.backend = backend
 26        self.formulae = backend.formulae
 27        self.environment = None
 28        self.attributes: (ParticleAttributes, None) = None
 29        self.dynamics = {}
 30        self.products = {}
 31        self.observers = []
 32
 33        self.n_steps = 0
 34
 35        self.sorting_scheme = "default"
 36        self.condensation_solver = None
 37
 38        self.Index = make_Index(backend)  # pylint: disable=invalid-name
 39        self.PairIndicator = make_PairIndicator(backend)  # pylint: disable=invalid-name
 40        self.PairwiseStorage = make_PairwiseStorage(  # pylint: disable=invalid-name
 41            backend
 42        )
 43        self.IndexedStorage = make_IndexedStorage(  # pylint: disable=invalid-name
 44            backend
 45        )
 46
 47        self.timers = {}
 48        self.null = self.Storage.empty(0, dtype=float)
 49
 50    def run(self, steps):
 51        for _ in range(steps):
 52            for key, dynamic in self.dynamics.items():
 53                with self.timers[key]:
 54                    dynamic()
 55            self.n_steps += 1
 56            self._notify_observers()
 57
 58    def _notify_observers(self):
 59        reversed_order_so_that_environment_is_last = reversed(self.observers)
 60        for observer in reversed_order_so_that_environment_is_last:
 61            observer.notify()
 62
 63    @property
 64    def Storage(self):
 65        return self.backend.Storage
 66
 67    @property
 68    def Random(self):
 69        return self.backend.Random
 70
 71    @property
 72    def n_sd(self) -> int:
 73        return self.__n_sd
 74
 75    @property
 76    def dt(self) -> float:
 77        if self.environment is not None:
 78            return self.environment.dt
 79        return None
 80
 81    @property
 82    def mesh(self):
 83        if self.environment is not None:
 84            return self.environment.mesh
 85        return None
 86
 87    def normalize(self, prob, norm_factor):
 88        self.backend.normalize(
 89            prob=prob,
 90            cell_id=self.attributes["cell id"],
 91            cell_idx=self.attributes.cell_idx,
 92            cell_start=self.attributes.cell_start,
 93            norm_factor=norm_factor,
 94            timestep=self.dt,
 95            dv=self.mesh.dv,
 96        )
 97
 98    def update_TpRH(self):
 99        self.backend.temperature_pressure_rh(
100            # input
101            rhod=self.environment.get_predicted("rhod"),
102            thd=self.environment.get_predicted("thd"),
103            water_vapour_mixing_ratio=self.environment.get_predicted(
104                "water_vapour_mixing_ratio"
105            ),
106            # output
107            T=self.environment.get_predicted("T"),
108            p=self.environment.get_predicted("p"),
109            RH=self.environment.get_predicted("RH"),
110        )
111
112    def condensation(self, *, rtol_x, rtol_thd, counters, RH_max, success, cell_order):
113        """Updates droplet volumes by simulating condensation driven by prior changes
114          in environment thermodynamic state, updates the environment state.
115        In the case of parcel environment, condensation is driven solely by changes in
116          the dry-air density (theta and water_vapour_mixing_ratio should not be changed
117          by other dynamics).
118        In the case of prescribed-flow/kinematic environments, the dry-air density is
119          constant in time throughout the simulation.
120        This function should only change environment's predicted `thd` and
121          `water_vapour_mixing_ratio` (and not `rhod`).
122        """
123        self.backend.condensation(
124            solver=self.condensation_solver,
125            n_cell=self.mesh.n_cell,
126            cell_start_arg=self.attributes.cell_start,
127            water_mass=self.attributes["signed water mass"],
128            multiplicity=self.attributes["multiplicity"],
129            vdry=self.attributes["dry volume"],
130            idx=self.attributes._ParticleAttributes__idx,
131            rhod=self.environment["rhod"],
132            thd=self.environment["thd"],
133            water_vapour_mixing_ratio=self.environment["water_vapour_mixing_ratio"],
134            dv=self.environment.dv,
135            prhod=self.environment.get_predicted("rhod"),
136            pthd=self.environment.get_predicted("thd"),
137            predicted_water_vapour_mixing_ratio=self.environment.get_predicted(
138                "water_vapour_mixing_ratio"
139            ),
140            kappa=self.attributes["kappa"],
141            f_org=self.attributes["dry volume organic fraction"],
142            rtol_x=rtol_x,
143            rtol_thd=rtol_thd,
144            v_cr=self.attributes["critical volume"],
145            timestep=self.dt,
146            counters=counters,
147            cell_order=cell_order,
148            RH_max=RH_max,
149            success=success,
150            cell_id=self.attributes["cell id"],
151            reynolds_number=self.attributes["Reynolds number"],
152            air_density=self.environment["air density"],
153            air_dynamic_viscosity=self.environment["air dynamic viscosity"],
154        )
155        self.attributes.mark_updated("signed water mass")
156
157    def collision_coalescence_breakup(
158        self,
159        *,
160        enable_breakup,
161        gamma,
162        rand,
163        Ec,
164        Eb,
165        fragment_mass,
166        coalescence_rate,
167        breakup_rate,
168        breakup_rate_deficit,
169        is_first_in_pair,
170        warn_overflows,
171        max_multiplicity,
172    ):
173        # pylint: disable=too-many-locals
174        idx = self.attributes._ParticleAttributes__idx
175        healthy = self.attributes._ParticleAttributes__healthy_memory
176        cell_id = self.attributes["cell id"]
177        multiplicity = self.attributes["multiplicity"]
178        attributes = self.attributes.get_extensive_attribute_storage()
179        if enable_breakup:
180            self.backend.collision_coalescence_breakup(
181                multiplicity=multiplicity,
182                idx=idx,
183                attributes=attributes,
184                gamma=gamma,
185                rand=rand,
186                Ec=Ec,
187                Eb=Eb,
188                fragment_mass=fragment_mass,
189                healthy=healthy,
190                cell_id=cell_id,
191                coalescence_rate=coalescence_rate,
192                breakup_rate=breakup_rate,
193                breakup_rate_deficit=breakup_rate_deficit,
194                is_first_in_pair=is_first_in_pair,
195                warn_overflows=warn_overflows,
196                particle_mass=self.attributes["water mass"],
197                max_multiplicity=max_multiplicity,
198            )
199        else:
200            self.backend.collision_coalescence(
201                multiplicity=multiplicity,
202                idx=idx,
203                attributes=attributes,
204                gamma=gamma,
205                healthy=healthy,
206                cell_id=cell_id,
207                coalescence_rate=coalescence_rate,
208                is_first_in_pair=is_first_in_pair,
209            )
210        self.attributes.sanitize()
211        self.attributes.mark_updated("multiplicity")
212        for key in self.attributes.get_extensive_attribute_keys():
213            self.attributes.mark_updated(key)
214
215    def oxidation(
216        self,
217        *,
218        kinetic_consts,
219        timestep,
220        equilibrium_consts,
221        dissociation_factors,
222        do_chemistry_flag,
223    ):
224        self.backend.oxidation(
225            n_sd=self.n_sd,
226            cell_ids=self.attributes["cell id"],
227            do_chemistry_flag=do_chemistry_flag,
228            k0=kinetic_consts["k0"],
229            k1=kinetic_consts["k1"],
230            k2=kinetic_consts["k2"],
231            k3=kinetic_consts["k3"],
232            K_SO2=equilibrium_consts["K_SO2"],
233            K_HSO3=equilibrium_consts["K_HSO3"],
234            dissociation_factor_SO2=dissociation_factors["SO2"],
235            timestep=timestep,
236            # input
237            droplet_volume=self.attributes["volume"],
238            pH=self.attributes["pH"],
239            # output
240            moles_O3=self.attributes["moles_O3"],
241            moles_H2O2=self.attributes["moles_H2O2"],
242            moles_S_IV=self.attributes["moles_S_IV"],
243            moles_S_VI=self.attributes["moles_S_VI"],
244        )
245        for attr in ("moles_S_IV", "moles_S_VI", "moles_H2O2", "moles_O3"):
246            self.attributes.mark_updated(attr)
247
248    def dissolution(
249        self,
250        *,
251        gaseous_compounds,
252        system_type,
253        dissociation_factors,
254        timestep,
255        environment_mixing_ratios,
256        do_chemistry_flag,
257    ):
258        self.backend.dissolution(
259            n_cell=self.mesh.n_cell,
260            n_threads=1,
261            cell_order=np.arange(self.mesh.n_cell),
262            cell_start_arg=self.attributes.cell_start,
263            idx=self.attributes._ParticleAttributes__idx,
264            do_chemistry_flag=do_chemistry_flag,
265            mole_amounts={
266                key: self.attributes["moles_" + key] for key in gaseous_compounds.keys()
267            },
268            env_mixing_ratio=environment_mixing_ratios,
269            # note: assuming condensation was called
270            env_p=self.environment.get_predicted("p"),
271            env_T=self.environment.get_predicted("T"),
272            env_rho_d=self.environment.get_predicted("rhod"),
273            timestep=timestep,
274            dv=self.mesh.dv,
275            droplet_volume=self.attributes["volume"],
276            multiplicity=self.attributes["multiplicity"],
277            system_type=system_type,
278            dissociation_factors=dissociation_factors,
279        )
280        for key in gaseous_compounds.keys():
281            self.attributes.mark_updated(f"moles_{key}")
282
283    def chem_recalculate_cell_data(self, *, equilibrium_consts, kinetic_consts):
284        self.backend.chem_recalculate_cell_data(
285            equilibrium_consts=equilibrium_consts,
286            kinetic_consts=kinetic_consts,
287            temperature=self.environment.get_predicted("T"),
288        )
289
290    def chem_recalculate_drop_data(self, *, dissociation_factors, equilibrium_consts):
291        self.backend.chem_recalculate_drop_data(
292            dissociation_factors=dissociation_factors,
293            equilibrium_consts=equilibrium_consts,
294            pH=self.attributes["pH"],
295            cell_id=self.attributes["cell id"],
296        )
297
298    def recalculate_cell_id(self):
299        if not self.attributes.has_attribute("cell origin"):
300            return
301        self.backend.cell_id(
302            self.attributes["cell id"],
303            self.attributes["cell origin"],
304            self.backend.Storage.from_ndarray(self.environment.mesh.strides),
305        )
306        self.attributes._ParticleAttributes__sorted = False
307
308    def sort_within_pair_by_attr(self, is_first_in_pair, attr_name):
309        self.backend.sort_within_pair_by_attr(
310            self.attributes._ParticleAttributes__idx,
311            is_first_in_pair,
312            self.attributes[attr_name],
313        )
314
315    def moments(
316        self,
317        *,
318        moment_0,
319        moments,
320        specs: dict,
321        attr_name="signed water mass",
322        attr_range=(-np.inf, np.inf),
323        weighting_attribute="water mass",
324        weighting_rank=0,
325        skip_division_by_m0=False,
326    ):
327        """
328        Writes to `moment_0` and `moment` the zero-th and the k-th statistical moments
329        of particle attributes computed filtering by value of the attribute `attr_name`
330        to fall within `attr_range`. The moment ranks are defined by `specs`.
331
332        Parameters:
333            specs: e.g., `specs={'volume': (1,2,3), 'kappa': (1)}` computes three moments
334                of volume and one moment of kappa
335            skip_division_by_m0: if set to `True`, the values written to `moments` are
336                multiplied by the 0-th moment (e.g., total volume instead of mean volume)
337        """
338        if len(specs) == 0:
339            raise ValueError("empty specs passed")
340        attr_data, ranks = [], []
341        for attr in specs:
342            for rank in specs[attr]:
343                attr_data.append(self.attributes[attr])
344                ranks.append(rank)
345        assert len(set(attr_data)) <= 1
346        if len(attr_data) == 0:
347            attr_data = self.backend.Storage.empty((0,), dtype=float)
348        else:
349            attr_data = attr_data[0]
350
351        ranks = self.backend.Storage.from_ndarray(np.array(ranks, dtype=float))
352
353        self.backend.moments(
354            moment_0=moment_0,
355            moments=moments,
356            multiplicity=self.attributes["multiplicity"],
357            attr_data=attr_data,
358            cell_id=self.attributes["cell id"],
359            idx=self.attributes._ParticleAttributes__idx,
360            length=self.attributes.super_droplet_count,
361            ranks=ranks,
362            min_x=attr_range[0],
363            max_x=attr_range[1],
364            x_attr=self.attributes[attr_name],
365            weighting_attribute=self.attributes[weighting_attribute],
366            weighting_rank=weighting_rank,
367            skip_division_by_m0=skip_division_by_m0,
368        )
369
370    def spectrum_moments(
371        self,
372        *,
373        moment_0,
374        moments,
375        attr,
376        rank,
377        attr_bins,
378        attr_name="water mass",
379        weighting_attribute="water mass",
380        weighting_rank=0,
381    ):
382        attr_data = self.attributes[attr]
383        self.backend.spectrum_moments(
384            moment_0=moment_0,
385            moments=moments,
386            multiplicity=self.attributes["multiplicity"],
387            attr_data=attr_data,
388            cell_id=self.attributes["cell id"],
389            idx=self.attributes._ParticleAttributes__idx,
390            length=self.attributes.super_droplet_count,
391            rank=rank,
392            x_bins=attr_bins,
393            x_attr=self.attributes[attr_name],
394            weighting_attribute=self.attributes[weighting_attribute],
395            weighting_rank=weighting_rank,
396        )
397
398    def adaptive_sdm_end(self, dt_left):
399        return self.backend.adaptive_sdm_end(dt_left, self.attributes.cell_start)
400
401    def remove_precipitated(
402        self, *, displacement, precipitation_counting_level_index
403    ) -> float:
404        rainfall_mass = self.backend.flag_precipitated(
405            cell_origin=self.attributes["cell origin"],
406            position_in_cell=self.attributes["position in cell"],
407            water_mass=self.attributes["water mass"],
408            multiplicity=self.attributes["multiplicity"],
409            idx=self.attributes._ParticleAttributes__idx,
410            length=self.attributes.super_droplet_count,
411            healthy=self.attributes._ParticleAttributes__healthy_memory,
412            precipitation_counting_level_index=precipitation_counting_level_index,
413            displacement=displacement,
414        )
415        self.attributes.sanitize()
416        return rainfall_mass
417
418    def flag_out_of_column(self):
419        self.backend.flag_out_of_column(
420            cell_origin=self.attributes["cell origin"],
421            position_in_cell=self.attributes["position in cell"],
422            idx=self.attributes._ParticleAttributes__idx,
423            length=self.attributes.super_droplet_count,
424            healthy=self.attributes._ParticleAttributes__healthy_memory,
425            domain_top_level_index=self.mesh.grid[-1],
426        )
427        self.attributes.sanitize()
428
429    def calculate_displacement(
430        self, *, displacement, courant, cell_origin, position_in_cell, n_substeps
431    ):
432        for dim in range(len(self.environment.mesh.grid)):
433            self.backend.calculate_displacement(
434                dim=dim,
435                displacement=displacement,
436                courant=courant[dim],
437                cell_origin=cell_origin,
438                position_in_cell=position_in_cell,
439                n_substeps=n_substeps,
440            )
441
442    def isotopic_fractionation(self, heavy_isotopes: tuple):
443        self.backend.isotopic_fractionation()
444        for isotope in heavy_isotopes:
445            self.attributes.mark_updated(f"moles_{isotope}")
446
447    def seeding(
448        self,
449        *,
450        seeded_particle_index,
451        seeded_particle_multiplicity,
452        seeded_particle_extensive_attributes,
453        number_of_super_particles_to_inject,
454    ):
455        n_null = self.n_sd - self.attributes.super_droplet_count
456        if n_null == 0:
457            raise ValueError(
458                "No available seeds to inject. Please provide particles with nan filled attributes."
459            )
460
461        if number_of_super_particles_to_inject > n_null:
462            raise ValueError(
463                "Trying to inject more super particles than space available."
464            )
465
466        if number_of_super_particles_to_inject > len(seeded_particle_multiplicity):
467            raise ValueError(
468                "Trying to inject multiple super particles with the same attributes. \
469                Instead increase multiplicity of injected particles."
470            )
471
472        self.backend.seeding(
473            idx=self.attributes._ParticleAttributes__idx,
474            multiplicity=self.attributes["multiplicity"],
475            extensive_attributes=self.attributes.get_extensive_attribute_storage(),
476            seeded_particle_index=seeded_particle_index,
477            seeded_particle_multiplicity=seeded_particle_multiplicity,
478            seeded_particle_extensive_attributes=seeded_particle_extensive_attributes,
479            number_of_super_particles_to_inject=number_of_super_particles_to_inject,
480        )
481        self.attributes.reset_idx()
482        self.attributes.sanitize()
483
484        self.attributes.mark_updated("multiplicity")
485        for key in self.attributes.get_extensive_attribute_keys():
486            self.attributes.mark_updated(key)
487
488    def deposition(self):
489        self.backend.deposition(
490            multiplicity=self.attributes["multiplicity"],
491            signed_water_mass=self.attributes["signed water mass"],
492            current_temperature=self.environment["T"],
493            current_total_pressure=self.environment["p"],
494            current_relative_humidity=self.environment["RH"],
495            current_water_activity=self.environment["a_w_ice"],
496            current_vapour_mixing_ratio=self.environment["water_vapour_mixing_ratio"],
497            current_dry_air_density=self.environment["rhod"],
498            current_dry_potential_temperature=self.environment["thd"],
499            cell_volume=self.environment.mesh.dv,
500            time_step=self.dt,
501            cell_id=self.attributes["cell id"],
502            reynolds_number=self.attributes["Reynolds number"],
503            schmidt_number=self.environment["Schmidt number"],
504            predicted_vapour_mixing_ratio=self.environment.get_predicted(
505                "water_vapour_mixing_ratio"
506            ),
507            predicted_dry_potential_temperature=self.environment.get_predicted("thd"),
508        )
509        self.attributes.mark_updated("signed water mass")
510        # TODO #1524 - should we update here?
511        # self.update_TpRH(only_if_not_last='VapourDepositionOnIce')
512
513    def immersion_freezing_time_dependent(self, *, thaw: bool, rand: Storage):
514        self.backend.freeze_time_dependent(
515            rand=rand,
516            attributes=TimeDependentAttributes(
517                immersed_surface_area=self.attributes["immersed surface area"],
518                signed_water_mass=self.attributes["signed water mass"],
519            ),
520            timestep=self.dt,
521            cell=self.attributes["cell id"],
522            a_w_ice=self.environment["a_w_ice"],
523            temperature=self.environment["T"],
524            relative_humidity=self.environment["RH"],
525            thaw=thaw,
526        )
527        self.attributes.mark_updated("signed water mass")
528
529    def immersion_freezing_singular(self, *, thaw: bool):
530        self.backend.freeze_singular(
531            attributes=SingularAttributes(
532                freezing_temperature=self.attributes["freezing temperature"],
533                signed_water_mass=self.attributes["signed water mass"],
534            ),
535            temperature=self.environment["T"],
536            relative_humidity=self.environment["RH"],
537            cell=self.attributes["cell id"],
538            thaw=thaw,
539        )
540        self.attributes.mark_updated("signed water mass")
Particulator(n_sd, backend)
21    def __init__(self, n_sd, backend):
22        assert isinstance(backend, BackendMethods)
23        self.__n_sd = n_sd
24
25        self.backend = backend
26        self.formulae = backend.formulae
27        self.environment = None
28        self.attributes: (ParticleAttributes, None) = None
29        self.dynamics = {}
30        self.products = {}
31        self.observers = []
32
33        self.n_steps = 0
34
35        self.sorting_scheme = "default"
36        self.condensation_solver = None
37
38        self.Index = make_Index(backend)  # pylint: disable=invalid-name
39        self.PairIndicator = make_PairIndicator(backend)  # pylint: disable=invalid-name
40        self.PairwiseStorage = make_PairwiseStorage(  # pylint: disable=invalid-name
41            backend
42        )
43        self.IndexedStorage = make_IndexedStorage(  # pylint: disable=invalid-name
44            backend
45        )
46
47        self.timers = {}
48        self.null = self.Storage.empty(0, dtype=float)
backend
formulae
environment
dynamics
products
observers
n_steps
sorting_scheme
condensation_solver
Index
PairIndicator
PairwiseStorage
IndexedStorage
timers
null
def run(self, steps):
50    def run(self, steps):
51        for _ in range(steps):
52            for key, dynamic in self.dynamics.items():
53                with self.timers[key]:
54                    dynamic()
55            self.n_steps += 1
56            self._notify_observers()
Storage
63    @property
64    def Storage(self):
65        return self.backend.Storage
Random
67    @property
68    def Random(self):
69        return self.backend.Random
n_sd: int
71    @property
72    def n_sd(self) -> int:
73        return self.__n_sd
dt: float
75    @property
76    def dt(self) -> float:
77        if self.environment is not None:
78            return self.environment.dt
79        return None
mesh
81    @property
82    def mesh(self):
83        if self.environment is not None:
84            return self.environment.mesh
85        return None
def normalize(self, prob, norm_factor):
87    def normalize(self, prob, norm_factor):
88        self.backend.normalize(
89            prob=prob,
90            cell_id=self.attributes["cell id"],
91            cell_idx=self.attributes.cell_idx,
92            cell_start=self.attributes.cell_start,
93            norm_factor=norm_factor,
94            timestep=self.dt,
95            dv=self.mesh.dv,
96        )
def update_TpRH(self):
 98    def update_TpRH(self):
 99        self.backend.temperature_pressure_rh(
100            # input
101            rhod=self.environment.get_predicted("rhod"),
102            thd=self.environment.get_predicted("thd"),
103            water_vapour_mixing_ratio=self.environment.get_predicted(
104                "water_vapour_mixing_ratio"
105            ),
106            # output
107            T=self.environment.get_predicted("T"),
108            p=self.environment.get_predicted("p"),
109            RH=self.environment.get_predicted("RH"),
110        )
def condensation(self, *, rtol_x, rtol_thd, counters, RH_max, success, cell_order):
112    def condensation(self, *, rtol_x, rtol_thd, counters, RH_max, success, cell_order):
113        """Updates droplet volumes by simulating condensation driven by prior changes
114          in environment thermodynamic state, updates the environment state.
115        In the case of parcel environment, condensation is driven solely by changes in
116          the dry-air density (theta and water_vapour_mixing_ratio should not be changed
117          by other dynamics).
118        In the case of prescribed-flow/kinematic environments, the dry-air density is
119          constant in time throughout the simulation.
120        This function should only change environment's predicted `thd` and
121          `water_vapour_mixing_ratio` (and not `rhod`).
122        """
123        self.backend.condensation(
124            solver=self.condensation_solver,
125            n_cell=self.mesh.n_cell,
126            cell_start_arg=self.attributes.cell_start,
127            water_mass=self.attributes["signed water mass"],
128            multiplicity=self.attributes["multiplicity"],
129            vdry=self.attributes["dry volume"],
130            idx=self.attributes._ParticleAttributes__idx,
131            rhod=self.environment["rhod"],
132            thd=self.environment["thd"],
133            water_vapour_mixing_ratio=self.environment["water_vapour_mixing_ratio"],
134            dv=self.environment.dv,
135            prhod=self.environment.get_predicted("rhod"),
136            pthd=self.environment.get_predicted("thd"),
137            predicted_water_vapour_mixing_ratio=self.environment.get_predicted(
138                "water_vapour_mixing_ratio"
139            ),
140            kappa=self.attributes["kappa"],
141            f_org=self.attributes["dry volume organic fraction"],
142            rtol_x=rtol_x,
143            rtol_thd=rtol_thd,
144            v_cr=self.attributes["critical volume"],
145            timestep=self.dt,
146            counters=counters,
147            cell_order=cell_order,
148            RH_max=RH_max,
149            success=success,
150            cell_id=self.attributes["cell id"],
151            reynolds_number=self.attributes["Reynolds number"],
152            air_density=self.environment["air density"],
153            air_dynamic_viscosity=self.environment["air dynamic viscosity"],
154        )
155        self.attributes.mark_updated("signed water mass")

Updates droplet volumes by simulating condensation driven by prior changes in environment thermodynamic state, updates the environment state. In the case of parcel environment, condensation is driven solely by changes in the dry-air density (theta and water_vapour_mixing_ratio should not be changed by other dynamics). In the case of prescribed-flow/kinematic environments, the dry-air density is constant in time throughout the simulation. This function should only change environment's predicted thd and water_vapour_mixing_ratio (and not rhod).

def collision_coalescence_breakup( self, *, enable_breakup, gamma, rand, Ec, Eb, fragment_mass, coalescence_rate, breakup_rate, breakup_rate_deficit, is_first_in_pair, warn_overflows, max_multiplicity):
157    def collision_coalescence_breakup(
158        self,
159        *,
160        enable_breakup,
161        gamma,
162        rand,
163        Ec,
164        Eb,
165        fragment_mass,
166        coalescence_rate,
167        breakup_rate,
168        breakup_rate_deficit,
169        is_first_in_pair,
170        warn_overflows,
171        max_multiplicity,
172    ):
173        # pylint: disable=too-many-locals
174        idx = self.attributes._ParticleAttributes__idx
175        healthy = self.attributes._ParticleAttributes__healthy_memory
176        cell_id = self.attributes["cell id"]
177        multiplicity = self.attributes["multiplicity"]
178        attributes = self.attributes.get_extensive_attribute_storage()
179        if enable_breakup:
180            self.backend.collision_coalescence_breakup(
181                multiplicity=multiplicity,
182                idx=idx,
183                attributes=attributes,
184                gamma=gamma,
185                rand=rand,
186                Ec=Ec,
187                Eb=Eb,
188                fragment_mass=fragment_mass,
189                healthy=healthy,
190                cell_id=cell_id,
191                coalescence_rate=coalescence_rate,
192                breakup_rate=breakup_rate,
193                breakup_rate_deficit=breakup_rate_deficit,
194                is_first_in_pair=is_first_in_pair,
195                warn_overflows=warn_overflows,
196                particle_mass=self.attributes["water mass"],
197                max_multiplicity=max_multiplicity,
198            )
199        else:
200            self.backend.collision_coalescence(
201                multiplicity=multiplicity,
202                idx=idx,
203                attributes=attributes,
204                gamma=gamma,
205                healthy=healthy,
206                cell_id=cell_id,
207                coalescence_rate=coalescence_rate,
208                is_first_in_pair=is_first_in_pair,
209            )
210        self.attributes.sanitize()
211        self.attributes.mark_updated("multiplicity")
212        for key in self.attributes.get_extensive_attribute_keys():
213            self.attributes.mark_updated(key)
def oxidation( self, *, kinetic_consts, timestep, equilibrium_consts, dissociation_factors, do_chemistry_flag):
215    def oxidation(
216        self,
217        *,
218        kinetic_consts,
219        timestep,
220        equilibrium_consts,
221        dissociation_factors,
222        do_chemistry_flag,
223    ):
224        self.backend.oxidation(
225            n_sd=self.n_sd,
226            cell_ids=self.attributes["cell id"],
227            do_chemistry_flag=do_chemistry_flag,
228            k0=kinetic_consts["k0"],
229            k1=kinetic_consts["k1"],
230            k2=kinetic_consts["k2"],
231            k3=kinetic_consts["k3"],
232            K_SO2=equilibrium_consts["K_SO2"],
233            K_HSO3=equilibrium_consts["K_HSO3"],
234            dissociation_factor_SO2=dissociation_factors["SO2"],
235            timestep=timestep,
236            # input
237            droplet_volume=self.attributes["volume"],
238            pH=self.attributes["pH"],
239            # output
240            moles_O3=self.attributes["moles_O3"],
241            moles_H2O2=self.attributes["moles_H2O2"],
242            moles_S_IV=self.attributes["moles_S_IV"],
243            moles_S_VI=self.attributes["moles_S_VI"],
244        )
245        for attr in ("moles_S_IV", "moles_S_VI", "moles_H2O2", "moles_O3"):
246            self.attributes.mark_updated(attr)
def dissolution( self, *, gaseous_compounds, system_type, dissociation_factors, timestep, environment_mixing_ratios, do_chemistry_flag):
248    def dissolution(
249        self,
250        *,
251        gaseous_compounds,
252        system_type,
253        dissociation_factors,
254        timestep,
255        environment_mixing_ratios,
256        do_chemistry_flag,
257    ):
258        self.backend.dissolution(
259            n_cell=self.mesh.n_cell,
260            n_threads=1,
261            cell_order=np.arange(self.mesh.n_cell),
262            cell_start_arg=self.attributes.cell_start,
263            idx=self.attributes._ParticleAttributes__idx,
264            do_chemistry_flag=do_chemistry_flag,
265            mole_amounts={
266                key: self.attributes["moles_" + key] for key in gaseous_compounds.keys()
267            },
268            env_mixing_ratio=environment_mixing_ratios,
269            # note: assuming condensation was called
270            env_p=self.environment.get_predicted("p"),
271            env_T=self.environment.get_predicted("T"),
272            env_rho_d=self.environment.get_predicted("rhod"),
273            timestep=timestep,
274            dv=self.mesh.dv,
275            droplet_volume=self.attributes["volume"],
276            multiplicity=self.attributes["multiplicity"],
277            system_type=system_type,
278            dissociation_factors=dissociation_factors,
279        )
280        for key in gaseous_compounds.keys():
281            self.attributes.mark_updated(f"moles_{key}")
def chem_recalculate_cell_data(self, *, equilibrium_consts, kinetic_consts):
283    def chem_recalculate_cell_data(self, *, equilibrium_consts, kinetic_consts):
284        self.backend.chem_recalculate_cell_data(
285            equilibrium_consts=equilibrium_consts,
286            kinetic_consts=kinetic_consts,
287            temperature=self.environment.get_predicted("T"),
288        )
def chem_recalculate_drop_data(self, *, dissociation_factors, equilibrium_consts):
290    def chem_recalculate_drop_data(self, *, dissociation_factors, equilibrium_consts):
291        self.backend.chem_recalculate_drop_data(
292            dissociation_factors=dissociation_factors,
293            equilibrium_consts=equilibrium_consts,
294            pH=self.attributes["pH"],
295            cell_id=self.attributes["cell id"],
296        )
def recalculate_cell_id(self):
298    def recalculate_cell_id(self):
299        if not self.attributes.has_attribute("cell origin"):
300            return
301        self.backend.cell_id(
302            self.attributes["cell id"],
303            self.attributes["cell origin"],
304            self.backend.Storage.from_ndarray(self.environment.mesh.strides),
305        )
306        self.attributes._ParticleAttributes__sorted = False
def sort_within_pair_by_attr(self, is_first_in_pair, attr_name):
308    def sort_within_pair_by_attr(self, is_first_in_pair, attr_name):
309        self.backend.sort_within_pair_by_attr(
310            self.attributes._ParticleAttributes__idx,
311            is_first_in_pair,
312            self.attributes[attr_name],
313        )
def moments( self, *, moment_0, moments, specs: dict, attr_name='signed water mass', attr_range=(-inf, inf), weighting_attribute='water mass', weighting_rank=0, skip_division_by_m0=False):
315    def moments(
316        self,
317        *,
318        moment_0,
319        moments,
320        specs: dict,
321        attr_name="signed water mass",
322        attr_range=(-np.inf, np.inf),
323        weighting_attribute="water mass",
324        weighting_rank=0,
325        skip_division_by_m0=False,
326    ):
327        """
328        Writes to `moment_0` and `moment` the zero-th and the k-th statistical moments
329        of particle attributes computed filtering by value of the attribute `attr_name`
330        to fall within `attr_range`. The moment ranks are defined by `specs`.
331
332        Parameters:
333            specs: e.g., `specs={'volume': (1,2,3), 'kappa': (1)}` computes three moments
334                of volume and one moment of kappa
335            skip_division_by_m0: if set to `True`, the values written to `moments` are
336                multiplied by the 0-th moment (e.g., total volume instead of mean volume)
337        """
338        if len(specs) == 0:
339            raise ValueError("empty specs passed")
340        attr_data, ranks = [], []
341        for attr in specs:
342            for rank in specs[attr]:
343                attr_data.append(self.attributes[attr])
344                ranks.append(rank)
345        assert len(set(attr_data)) <= 1
346        if len(attr_data) == 0:
347            attr_data = self.backend.Storage.empty((0,), dtype=float)
348        else:
349            attr_data = attr_data[0]
350
351        ranks = self.backend.Storage.from_ndarray(np.array(ranks, dtype=float))
352
353        self.backend.moments(
354            moment_0=moment_0,
355            moments=moments,
356            multiplicity=self.attributes["multiplicity"],
357            attr_data=attr_data,
358            cell_id=self.attributes["cell id"],
359            idx=self.attributes._ParticleAttributes__idx,
360            length=self.attributes.super_droplet_count,
361            ranks=ranks,
362            min_x=attr_range[0],
363            max_x=attr_range[1],
364            x_attr=self.attributes[attr_name],
365            weighting_attribute=self.attributes[weighting_attribute],
366            weighting_rank=weighting_rank,
367            skip_division_by_m0=skip_division_by_m0,
368        )

Writes to moment_0 and moment the zero-th and the k-th statistical moments of particle attributes computed filtering by value of the attribute attr_name to fall within attr_range. The moment ranks are defined by specs.

Parameters: specs: e.g., specs={'volume': (1,2,3), 'kappa': (1)} computes three moments of volume and one moment of kappa skip_division_by_m0: if set to True, the values written to moments are multiplied by the 0-th moment (e.g., total volume instead of mean volume)

def spectrum_moments( self, *, moment_0, moments, attr, rank, attr_bins, attr_name='water mass', weighting_attribute='water mass', weighting_rank=0):
370    def spectrum_moments(
371        self,
372        *,
373        moment_0,
374        moments,
375        attr,
376        rank,
377        attr_bins,
378        attr_name="water mass",
379        weighting_attribute="water mass",
380        weighting_rank=0,
381    ):
382        attr_data = self.attributes[attr]
383        self.backend.spectrum_moments(
384            moment_0=moment_0,
385            moments=moments,
386            multiplicity=self.attributes["multiplicity"],
387            attr_data=attr_data,
388            cell_id=self.attributes["cell id"],
389            idx=self.attributes._ParticleAttributes__idx,
390            length=self.attributes.super_droplet_count,
391            rank=rank,
392            x_bins=attr_bins,
393            x_attr=self.attributes[attr_name],
394            weighting_attribute=self.attributes[weighting_attribute],
395            weighting_rank=weighting_rank,
396        )
def adaptive_sdm_end(self, dt_left):
398    def adaptive_sdm_end(self, dt_left):
399        return self.backend.adaptive_sdm_end(dt_left, self.attributes.cell_start)
def remove_precipitated(self, *, displacement, precipitation_counting_level_index) -> float:
401    def remove_precipitated(
402        self, *, displacement, precipitation_counting_level_index
403    ) -> float:
404        rainfall_mass = self.backend.flag_precipitated(
405            cell_origin=self.attributes["cell origin"],
406            position_in_cell=self.attributes["position in cell"],
407            water_mass=self.attributes["water mass"],
408            multiplicity=self.attributes["multiplicity"],
409            idx=self.attributes._ParticleAttributes__idx,
410            length=self.attributes.super_droplet_count,
411            healthy=self.attributes._ParticleAttributes__healthy_memory,
412            precipitation_counting_level_index=precipitation_counting_level_index,
413            displacement=displacement,
414        )
415        self.attributes.sanitize()
416        return rainfall_mass
def flag_out_of_column(self):
418    def flag_out_of_column(self):
419        self.backend.flag_out_of_column(
420            cell_origin=self.attributes["cell origin"],
421            position_in_cell=self.attributes["position in cell"],
422            idx=self.attributes._ParticleAttributes__idx,
423            length=self.attributes.super_droplet_count,
424            healthy=self.attributes._ParticleAttributes__healthy_memory,
425            domain_top_level_index=self.mesh.grid[-1],
426        )
427        self.attributes.sanitize()
def calculate_displacement( self, *, displacement, courant, cell_origin, position_in_cell, n_substeps):
429    def calculate_displacement(
430        self, *, displacement, courant, cell_origin, position_in_cell, n_substeps
431    ):
432        for dim in range(len(self.environment.mesh.grid)):
433            self.backend.calculate_displacement(
434                dim=dim,
435                displacement=displacement,
436                courant=courant[dim],
437                cell_origin=cell_origin,
438                position_in_cell=position_in_cell,
439                n_substeps=n_substeps,
440            )
def isotopic_fractionation(self, heavy_isotopes: tuple):
442    def isotopic_fractionation(self, heavy_isotopes: tuple):
443        self.backend.isotopic_fractionation()
444        for isotope in heavy_isotopes:
445            self.attributes.mark_updated(f"moles_{isotope}")
def seeding( self, *, seeded_particle_index, seeded_particle_multiplicity, seeded_particle_extensive_attributes, number_of_super_particles_to_inject):
447    def seeding(
448        self,
449        *,
450        seeded_particle_index,
451        seeded_particle_multiplicity,
452        seeded_particle_extensive_attributes,
453        number_of_super_particles_to_inject,
454    ):
455        n_null = self.n_sd - self.attributes.super_droplet_count
456        if n_null == 0:
457            raise ValueError(
458                "No available seeds to inject. Please provide particles with nan filled attributes."
459            )
460
461        if number_of_super_particles_to_inject > n_null:
462            raise ValueError(
463                "Trying to inject more super particles than space available."
464            )
465
466        if number_of_super_particles_to_inject > len(seeded_particle_multiplicity):
467            raise ValueError(
468                "Trying to inject multiple super particles with the same attributes. \
469                Instead increase multiplicity of injected particles."
470            )
471
472        self.backend.seeding(
473            idx=self.attributes._ParticleAttributes__idx,
474            multiplicity=self.attributes["multiplicity"],
475            extensive_attributes=self.attributes.get_extensive_attribute_storage(),
476            seeded_particle_index=seeded_particle_index,
477            seeded_particle_multiplicity=seeded_particle_multiplicity,
478            seeded_particle_extensive_attributes=seeded_particle_extensive_attributes,
479            number_of_super_particles_to_inject=number_of_super_particles_to_inject,
480        )
481        self.attributes.reset_idx()
482        self.attributes.sanitize()
483
484        self.attributes.mark_updated("multiplicity")
485        for key in self.attributes.get_extensive_attribute_keys():
486            self.attributes.mark_updated(key)
def deposition(self):
488    def deposition(self):
489        self.backend.deposition(
490            multiplicity=self.attributes["multiplicity"],
491            signed_water_mass=self.attributes["signed water mass"],
492            current_temperature=self.environment["T"],
493            current_total_pressure=self.environment["p"],
494            current_relative_humidity=self.environment["RH"],
495            current_water_activity=self.environment["a_w_ice"],
496            current_vapour_mixing_ratio=self.environment["water_vapour_mixing_ratio"],
497            current_dry_air_density=self.environment["rhod"],
498            current_dry_potential_temperature=self.environment["thd"],
499            cell_volume=self.environment.mesh.dv,
500            time_step=self.dt,
501            cell_id=self.attributes["cell id"],
502            reynolds_number=self.attributes["Reynolds number"],
503            schmidt_number=self.environment["Schmidt number"],
504            predicted_vapour_mixing_ratio=self.environment.get_predicted(
505                "water_vapour_mixing_ratio"
506            ),
507            predicted_dry_potential_temperature=self.environment.get_predicted("thd"),
508        )
509        self.attributes.mark_updated("signed water mass")
510        # TODO #1524 - should we update here?
511        # self.update_TpRH(only_if_not_last='VapourDepositionOnIce')
def immersion_freezing_time_dependent(self, *, thaw: bool, rand: <property object>):
513    def immersion_freezing_time_dependent(self, *, thaw: bool, rand: Storage):
514        self.backend.freeze_time_dependent(
515            rand=rand,
516            attributes=TimeDependentAttributes(
517                immersed_surface_area=self.attributes["immersed surface area"],
518                signed_water_mass=self.attributes["signed water mass"],
519            ),
520            timestep=self.dt,
521            cell=self.attributes["cell id"],
522            a_w_ice=self.environment["a_w_ice"],
523            temperature=self.environment["T"],
524            relative_humidity=self.environment["RH"],
525            thaw=thaw,
526        )
527        self.attributes.mark_updated("signed water mass")
def immersion_freezing_singular(self, *, thaw: bool):
529    def immersion_freezing_singular(self, *, thaw: bool):
530        self.backend.freeze_singular(
531            attributes=SingularAttributes(
532                freezing_temperature=self.attributes["freezing temperature"],
533                signed_water_mass=self.attributes["signed water mass"],
534            ),
535            temperature=self.environment["T"],
536            relative_humidity=self.environment["RH"],
537            cell=self.attributes["cell id"],
538            thaw=thaw,
539        )
540        self.attributes.mark_updated("signed water mass")