PySDM_examples.seeding.simulation

View Source

  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            dynamics=[
 31                AmbientThermodynamics(),
 32                Condensation(),
 33            ]
 34            + (
 35                []
 36                if not settings.enable_collisions
 37                else [
 38                    Coalescence(collision_kernel=Geometric()),
 39                ]
 40            )
 41            + [
 42                Seeding(
 43                    **{
 44                        k: getattr(settings, k)
 45                        for k in (
 46                            "super_droplet_injection_rate",
 47                            "seeded_particle_multiplicity",
 48                            "seeded_particle_extensive_attributes",
 49                        )
 50                    }
 51                ),
 52            ],
 53        )
 54
 55        r_dry, n_in_dv = ConstantMultiplicity(
 56            settings.initial_aerosol_dry_radii
 57        ).sample_deterministic(
 58            n_sd=settings.n_sd_initial, backend=builder.particulator.backend
 59        )
 60        attributes = builder.particulator.environment.init_attributes(
 61            n_in_dv=n_in_dv, kappa=settings.initial_aerosol_kappa, r_dry=r_dry
 62        )
 63        self.particulator = builder.build(
 64            attributes={
 65                k: np.pad(
 66                    array=v,
 67                    pad_width=(0, settings.n_sd_seeding),
 68                    mode="constant",
 69                    constant_values=np.nan if k == "multiplicity" else 0,
 70                )
 71                for k, v in attributes.items()
 72            },
 73            products=(
 74                products.SuperDropletCountPerGridbox(name="sd_count"),
 75                products.Time(),
 76                products.WaterMixingRatio(
 77                    radius_range=(settings.rain_water_radius_threshold, np.inf),
 78                    name="rain water mixing ratio",
 79                ),
 80                products.EffectiveRadius(
 81                    name="r_eff",
 82                    unit="um",
 83                    radius_range=(0.5 * si.um, 25 * si.um),
 84                ),
 85                products.ParticleConcentration(
 86                    name="n_drop",
 87                    unit="cm^-3",
 88                    radius_range=(0.5 * si.um, 25 * si.um),
 89                ),
 90            ),
 91        )
 92        self.n_steps = int(settings.t_max // settings.timestep)
 93
 94    def run(self):
 95        output = {
 96            "attributes": {"water mass": []},
 97            "products": {key: [] for key in self.particulator.products},
 98        }
 99        for step in range(self.n_steps + 1):
100            if step != 0:
101                self.particulator.run(steps=1)
102            for key, attr in output["attributes"].items():
103                data = self.particulator.attributes[key].to_ndarray(raw=True)
104                data[data == 0] = np.nan
105                attr.append(data)
106            for key, prod in output["products"].items():
107                value = self.particulator.products[key].get()
108                if not isinstance(value, float):
109                    (value,) = value
110                prod.append(float(value))
111        for out in ("attributes", "products"):
112            for key, val in output[out].items():
113                output[out][key] = np.array(val)
114        return output

class Simulation: View Source

 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            dynamics=[
 32                AmbientThermodynamics(),
 33                Condensation(),
 34            ]
 35            + (
 36                []
 37                if not settings.enable_collisions
 38                else [
 39                    Coalescence(collision_kernel=Geometric()),
 40                ]
 41            )
 42            + [
 43                Seeding(
 44                    **{
 45                        k: getattr(settings, k)
 46                        for k in (
 47                            "super_droplet_injection_rate",
 48                            "seeded_particle_multiplicity",
 49                            "seeded_particle_extensive_attributes",
 50                        )
 51                    }
 52                ),
 53            ],
 54        )
 55
 56        r_dry, n_in_dv = ConstantMultiplicity(
 57            settings.initial_aerosol_dry_radii
 58        ).sample_deterministic(
 59            n_sd=settings.n_sd_initial, backend=builder.particulator.backend
 60        )
 61        attributes = builder.particulator.environment.init_attributes(
 62            n_in_dv=n_in_dv, kappa=settings.initial_aerosol_kappa, r_dry=r_dry
 63        )
 64        self.particulator = builder.build(
 65            attributes={
 66                k: np.pad(
 67                    array=v,
 68                    pad_width=(0, settings.n_sd_seeding),
 69                    mode="constant",
 70                    constant_values=np.nan if k == "multiplicity" else 0,
 71                )
 72                for k, v in attributes.items()
 73            },
 74            products=(
 75                products.SuperDropletCountPerGridbox(name="sd_count"),
 76                products.Time(),
 77                products.WaterMixingRatio(
 78                    radius_range=(settings.rain_water_radius_threshold, np.inf),
 79                    name="rain water mixing ratio",
 80                ),
 81                products.EffectiveRadius(
 82                    name="r_eff",
 83                    unit="um",
 84                    radius_range=(0.5 * si.um, 25 * si.um),
 85                ),
 86                products.ParticleConcentration(
 87                    name="n_drop",
 88                    unit="cm^-3",
 89                    radius_range=(0.5 * si.um, 25 * si.um),
 90                ),
 91            ),
 92        )
 93        self.n_steps = int(settings.t_max // settings.timestep)
 94
 95    def run(self):
 96        output = {
 97            "attributes": {"water mass": []},
 98            "products": {key: [] for key in self.particulator.products},
 99        }
100        for step in range(self.n_steps + 1):
101            if step != 0:
102                self.particulator.run(steps=1)
103            for key, attr in output["attributes"].items():
104                data = self.particulator.attributes[key].to_ndarray(raw=True)
105                data[data == 0] = np.nan
106                attr.append(data)
107            for key, prod in output["products"].items():
108                value = self.particulator.products[key].get()
109                if not isinstance(value, float):
110                    (value,) = value
111                prod.append(float(value))
112        for out in ("attributes", "products"):
113            for key, val in output[out].items():
114                output[out][key] = np.array(val)
115        return output

Simulation(settings: PySDM_examples.seeding.settings.Settings) View Source

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            dynamics=[
32                AmbientThermodynamics(),
33                Condensation(),
34            ]
35            + (
36                []
37                if not settings.enable_collisions
38                else [
39                    Coalescence(collision_kernel=Geometric()),
40                ]
41            )
42            + [
43                Seeding(
44                    **{
45                        k: getattr(settings, k)
46                        for k in (
47                            "super_droplet_injection_rate",
48                            "seeded_particle_multiplicity",
49                            "seeded_particle_extensive_attributes",
50                        )
51                    }
52                ),
53            ],
54        )
55
56        r_dry, n_in_dv = ConstantMultiplicity(
57            settings.initial_aerosol_dry_radii
58        ).sample_deterministic(
59            n_sd=settings.n_sd_initial, backend=builder.particulator.backend
60        )
61        attributes = builder.particulator.environment.init_attributes(
62            n_in_dv=n_in_dv, kappa=settings.initial_aerosol_kappa, r_dry=r_dry
63        )
64        self.particulator = builder.build(
65            attributes={
66                k: np.pad(
67                    array=v,
68                    pad_width=(0, settings.n_sd_seeding),
69                    mode="constant",
70                    constant_values=np.nan if k == "multiplicity" else 0,
71                )
72                for k, v in attributes.items()
73            },
74            products=(
75                products.SuperDropletCountPerGridbox(name="sd_count"),
76                products.Time(),
77                products.WaterMixingRatio(
78                    radius_range=(settings.rain_water_radius_threshold, np.inf),
79                    name="rain water mixing ratio",
80                ),
81                products.EffectiveRadius(
82                    name="r_eff",
83                    unit="um",
84                    radius_range=(0.5 * si.um, 25 * si.um),
85                ),
86                products.ParticleConcentration(
87                    name="n_drop",
88                    unit="cm^-3",
89                    radius_range=(0.5 * si.um, 25 * si.um),
90                ),
91            ),
92        )
93        self.n_steps = int(settings.t_max // settings.timestep)

particulator

n_steps

def run(self): View Source

 95    def run(self):
 96        output = {
 97            "attributes": {"water mass": []},
 98            "products": {key: [] for key in self.particulator.products},
 99        }
100        for step in range(self.n_steps + 1):
101            if step != 0:
102                self.particulator.run(steps=1)
103            for key, attr in output["attributes"].items():
104                data = self.particulator.attributes[key].to_ndarray(raw=True)
105                data[data == 0] = np.nan
106                attr.append(data)
107            for key, prod in output["products"].items():
108                value = self.particulator.products[key].get()
109                if not isinstance(value, float):
110                    (value,) = value
111                prod.append(float(value))
112        for out in ("attributes", "products"):
113            for key, val in output[out].items():
114                output[out][key] = np.array(val)
115        return output