PySDM_examples.seeding.simulation

  1import numpy as np
  2
  3from PySDM_examples.seeding.settings import Settings
  4
  5from PySDM import Builder
  6from PySDM.backends import CPU
  7from PySDM.environments import Parcel
  8from PySDM.dynamics import Condensation, AmbientThermodynamics, Coalescence, Seeding
  9from PySDM.dynamics.collisions.collision_kernels import Geometric
 10from PySDM.initialisation.sampling.spectral_sampling import ConstantMultiplicity
 11from PySDM import products
 12from PySDM.physics import si
 13
 14
 15class Simulation:
 16    def __init__(self, settings: Settings):
 17        builder = Builder(
 18            n_sd=settings.n_sd_seeding + settings.n_sd_initial,
 19            backend=CPU(
 20                formulae=settings.formulae, override_jit_flags={"parallel": False}
 21            ),
 22            environment=Parcel(
 23                dt=settings.timestep,
 24                mass_of_dry_air=settings.mass_of_dry_air,
 25                w=settings.updraft,
 26                initial_water_vapour_mixing_ratio=settings.initial_water_vapour_mixing_ratio,
 27                p0=settings.initial_total_pressure,
 28                T0=settings.initial_temperature,
 29            ),
 30        )
 31        builder.add_dynamic(AmbientThermodynamics())
 32        builder.add_dynamic(Condensation())
 33        if settings.enable_collisions:
 34            builder.add_dynamic(Coalescence(collision_kernel=Geometric()))
 35        builder.add_dynamic(
 36            Seeding(
 37                super_droplet_injection_rate=settings.super_droplet_injection_rate,
 38                seeded_particle_multiplicity=settings.seeded_particle_multiplicity,
 39                seeded_particle_extensive_attributes=settings.seeded_particle_extensive_attributes,
 40            )
 41        )
 42
 43        r_dry, n_in_dv = ConstantMultiplicity(
 44            settings.initial_aerosol_dry_radii
 45        ).sample(n_sd=settings.n_sd_initial, backend=builder.particulator.backend)
 46        attributes = builder.particulator.environment.init_attributes(
 47            n_in_dv=n_in_dv, kappa=settings.initial_aerosol_kappa, r_dry=r_dry
 48        )
 49        self.particulator = builder.build(
 50            attributes={
 51                k: np.pad(
 52                    array=v,
 53                    pad_width=(0, settings.n_sd_seeding),
 54                    mode="constant",
 55                    constant_values=np.nan if k == "multiplicity" else 0,
 56                )
 57                for k, v in attributes.items()
 58            },
 59            products=(
 60                products.SuperDropletCountPerGridbox(name="sd_count"),
 61                products.Time(),
 62                products.WaterMixingRatio(
 63                    radius_range=(settings.rain_water_radius_threshold, np.inf),
 64                    name="rain water mixing ratio",
 65                ),
 66                products.EffectiveRadius(
 67                    name="r_eff",
 68                    unit="um",
 69                    radius_range=(0.5 * si.um, 25 * si.um),
 70                ),
 71                products.ParticleConcentration(
 72                    name="n_drop",
 73                    unit="cm^-3",
 74                    radius_range=(0.5 * si.um, 25 * si.um),
 75                ),
 76            ),
 77        )
 78        self.n_steps = int(settings.t_max // settings.timestep)
 79
 80    def run(self):
 81        output = {
 82            "attributes": {"water mass": []},
 83            "products": {key: [] for key in self.particulator.products},
 84        }
 85        for step in range(self.n_steps + 1):
 86            if step != 0:
 87                self.particulator.run(steps=1)
 88            for key, attr in output["attributes"].items():
 89                data = self.particulator.attributes[key].to_ndarray(raw=True)
 90                data[data == 0] = np.nan
 91                attr.append(data)
 92            for key, prod in output["products"].items():
 93                value = self.particulator.products[key].get()
 94                if not isinstance(value, float):
 95                    (value,) = value
 96                prod.append(float(value))
 97        for out in ("attributes", "products"):
 98            for key, val in output[out].items():
 99                output[out][key] = np.array(val)
100        return output
class Simulation:
 16class Simulation:
 17    def __init__(self, settings: Settings):
 18        builder = Builder(
 19            n_sd=settings.n_sd_seeding + settings.n_sd_initial,
 20            backend=CPU(
 21                formulae=settings.formulae, override_jit_flags={"parallel": False}
 22            ),
 23            environment=Parcel(
 24                dt=settings.timestep,
 25                mass_of_dry_air=settings.mass_of_dry_air,
 26                w=settings.updraft,
 27                initial_water_vapour_mixing_ratio=settings.initial_water_vapour_mixing_ratio,
 28                p0=settings.initial_total_pressure,
 29                T0=settings.initial_temperature,
 30            ),
 31        )
 32        builder.add_dynamic(AmbientThermodynamics())
 33        builder.add_dynamic(Condensation())
 34        if settings.enable_collisions:
 35            builder.add_dynamic(Coalescence(collision_kernel=Geometric()))
 36        builder.add_dynamic(
 37            Seeding(
 38                super_droplet_injection_rate=settings.super_droplet_injection_rate,
 39                seeded_particle_multiplicity=settings.seeded_particle_multiplicity,
 40                seeded_particle_extensive_attributes=settings.seeded_particle_extensive_attributes,
 41            )
 42        )
 43
 44        r_dry, n_in_dv = ConstantMultiplicity(
 45            settings.initial_aerosol_dry_radii
 46        ).sample(n_sd=settings.n_sd_initial, backend=builder.particulator.backend)
 47        attributes = builder.particulator.environment.init_attributes(
 48            n_in_dv=n_in_dv, kappa=settings.initial_aerosol_kappa, r_dry=r_dry
 49        )
 50        self.particulator = builder.build(
 51            attributes={
 52                k: np.pad(
 53                    array=v,
 54                    pad_width=(0, settings.n_sd_seeding),
 55                    mode="constant",
 56                    constant_values=np.nan if k == "multiplicity" else 0,
 57                )
 58                for k, v in attributes.items()
 59            },
 60            products=(
 61                products.SuperDropletCountPerGridbox(name="sd_count"),
 62                products.Time(),
 63                products.WaterMixingRatio(
 64                    radius_range=(settings.rain_water_radius_threshold, np.inf),
 65                    name="rain water mixing ratio",
 66                ),
 67                products.EffectiveRadius(
 68                    name="r_eff",
 69                    unit="um",
 70                    radius_range=(0.5 * si.um, 25 * si.um),
 71                ),
 72                products.ParticleConcentration(
 73                    name="n_drop",
 74                    unit="cm^-3",
 75                    radius_range=(0.5 * si.um, 25 * si.um),
 76                ),
 77            ),
 78        )
 79        self.n_steps = int(settings.t_max // settings.timestep)
 80
 81    def run(self):
 82        output = {
 83            "attributes": {"water mass": []},
 84            "products": {key: [] for key in self.particulator.products},
 85        }
 86        for step in range(self.n_steps + 1):
 87            if step != 0:
 88                self.particulator.run(steps=1)
 89            for key, attr in output["attributes"].items():
 90                data = self.particulator.attributes[key].to_ndarray(raw=True)
 91                data[data == 0] = np.nan
 92                attr.append(data)
 93            for key, prod in output["products"].items():
 94                value = self.particulator.products[key].get()
 95                if not isinstance(value, float):
 96                    (value,) = value
 97                prod.append(float(value))
 98        for out in ("attributes", "products"):
 99            for key, val in output[out].items():
100                output[out][key] = np.array(val)
101        return output
Simulation(settings: PySDM_examples.seeding.settings.Settings)
17    def __init__(self, settings: Settings):
18        builder = Builder(
19            n_sd=settings.n_sd_seeding + settings.n_sd_initial,
20            backend=CPU(
21                formulae=settings.formulae, override_jit_flags={"parallel": False}
22            ),
23            environment=Parcel(
24                dt=settings.timestep,
25                mass_of_dry_air=settings.mass_of_dry_air,
26                w=settings.updraft,
27                initial_water_vapour_mixing_ratio=settings.initial_water_vapour_mixing_ratio,
28                p0=settings.initial_total_pressure,
29                T0=settings.initial_temperature,
30            ),
31        )
32        builder.add_dynamic(AmbientThermodynamics())
33        builder.add_dynamic(Condensation())
34        if settings.enable_collisions:
35            builder.add_dynamic(Coalescence(collision_kernel=Geometric()))
36        builder.add_dynamic(
37            Seeding(
38                super_droplet_injection_rate=settings.super_droplet_injection_rate,
39                seeded_particle_multiplicity=settings.seeded_particle_multiplicity,
40                seeded_particle_extensive_attributes=settings.seeded_particle_extensive_attributes,
41            )
42        )
43
44        r_dry, n_in_dv = ConstantMultiplicity(
45            settings.initial_aerosol_dry_radii
46        ).sample(n_sd=settings.n_sd_initial, backend=builder.particulator.backend)
47        attributes = builder.particulator.environment.init_attributes(
48            n_in_dv=n_in_dv, kappa=settings.initial_aerosol_kappa, r_dry=r_dry
49        )
50        self.particulator = builder.build(
51            attributes={
52                k: np.pad(
53                    array=v,
54                    pad_width=(0, settings.n_sd_seeding),
55                    mode="constant",
56                    constant_values=np.nan if k == "multiplicity" else 0,
57                )
58                for k, v in attributes.items()
59            },
60            products=(
61                products.SuperDropletCountPerGridbox(name="sd_count"),
62                products.Time(),
63                products.WaterMixingRatio(
64                    radius_range=(settings.rain_water_radius_threshold, np.inf),
65                    name="rain water mixing ratio",
66                ),
67                products.EffectiveRadius(
68                    name="r_eff",
69                    unit="um",
70                    radius_range=(0.5 * si.um, 25 * si.um),
71                ),
72                products.ParticleConcentration(
73                    name="n_drop",
74                    unit="cm^-3",
75                    radius_range=(0.5 * si.um, 25 * si.um),
76                ),
77            ),
78        )
79        self.n_steps = int(settings.t_max // settings.timestep)
particulator
n_steps
def run(self):
 81    def run(self):
 82        output = {
 83            "attributes": {"water mass": []},
 84            "products": {key: [] for key in self.particulator.products},
 85        }
 86        for step in range(self.n_steps + 1):
 87            if step != 0:
 88                self.particulator.run(steps=1)
 89            for key, attr in output["attributes"].items():
 90                data = self.particulator.attributes[key].to_ndarray(raw=True)
 91                data[data == 0] = np.nan
 92                attr.append(data)
 93            for key, prod in output["products"].items():
 94                value = self.particulator.products[key].get()
 95                if not isinstance(value, float):
 96                    (value,) = value
 97                prod.append(float(value))
 98        for out in ("attributes", "products"):
 99            for key, val in output[out].items():
100                output[out][key] = np.array(val)
101        return output