PySDM_examples.Shima_et_al_2009.example

 1import os
 2from typing import Optional
 3
 4import numpy as np
 5from PySDM_examples.Shima_et_al_2009.settings import Settings
 6from PySDM_examples.Shima_et_al_2009.spectrum_plotter import SpectrumPlotter
 7
 8from PySDM.backends import CPU
 9from PySDM import Particulator
10from PySDM.dynamics import Coalescence
11from PySDM.environments import Box
12from PySDM.initialisation.sampling.spectral_sampling import ConstantMultiplicity
13from PySDM.products import ParticleVolumeVersusRadiusLogarithmSpectrum, WallTime
14
15
16def run(settings, backend=CPU, observers=()):
17    env = Box(
18        dv=settings.dv,
19        dt=settings.dt,
20        backend=backend(formulae=settings.formulae),
21    )
22    attributes = {}
23    sampling = ConstantMultiplicity(settings.spectrum)
24    attributes["volume"], attributes["multiplicity"] = sampling.sample_deterministic(
25        settings.n_sd
26    )
27    products = (
28        ParticleVolumeVersusRadiusLogarithmSpectrum(
29            settings.radius_bins_edges, name="dv/dlnr"
30        ),
31        WallTime(),
32    )
33    particulator = Particulator(
34        n_sd=settings.n_sd,
35        environment=env,
36        dynamics=(
37            Coalescence(collision_kernel=settings.kernel, adaptive=settings.adaptive),
38        ),
39        attributes=attributes,
40        products=products,
41    )
42
43    for observer in observers:
44        particulator.observers.append(observer)
45
46    vals = {}
47    particulator.products["wall time"].reset()
48    for step in settings.output_steps:
49        particulator.run(step - particulator.n_steps)
50        vals[step] = particulator.products["dv/dlnr"].get()[0]
51        vals[step][:] *= settings.rho
52
53    exec_time = particulator.products["wall time"].get()
54    return vals, exec_time
55
56
57def main(plot: bool, save: Optional[str]):
58    with np.errstate(all="raise"):
59        settings = Settings()
60
61        settings.n_sd = 2**15
62
63        states, _ = run(settings)
64
65    with np.errstate(invalid="ignore"):
66        plotter = SpectrumPlotter(settings)
67        plotter.smooth = True
68        for step, vals in states.items():
69            _ = plotter.plot(vals, step * settings.dt)
70            # assert _ < 200  # TODO #327
71        if save is not None:
72            n_sd = settings.n_sd
73            plotter.save(save + "/" + f"{n_sd}_shima_fig_2" + "." + plotter.format)
74        if plot:
75            plotter.show()
76
77
78if __name__ == "__main__":
79    main(plot="CI" not in os.environ, save=None)
def run( settings, backend=functools.partial(<function _cached_backend>, backend_class=<class 'PySDM.backends.Numba'>), observers=()):
17def run(settings, backend=CPU, observers=()):
18    env = Box(
19        dv=settings.dv,
20        dt=settings.dt,
21        backend=backend(formulae=settings.formulae),
22    )
23    attributes = {}
24    sampling = ConstantMultiplicity(settings.spectrum)
25    attributes["volume"], attributes["multiplicity"] = sampling.sample_deterministic(
26        settings.n_sd
27    )
28    products = (
29        ParticleVolumeVersusRadiusLogarithmSpectrum(
30            settings.radius_bins_edges, name="dv/dlnr"
31        ),
32        WallTime(),
33    )
34    particulator = Particulator(
35        n_sd=settings.n_sd,
36        environment=env,
37        dynamics=(
38            Coalescence(collision_kernel=settings.kernel, adaptive=settings.adaptive),
39        ),
40        attributes=attributes,
41        products=products,
42    )
43
44    for observer in observers:
45        particulator.observers.append(observer)
46
47    vals = {}
48    particulator.products["wall time"].reset()
49    for step in settings.output_steps:
50        particulator.run(step - particulator.n_steps)
51        vals[step] = particulator.products["dv/dlnr"].get()[0]
52        vals[step][:] *= settings.rho
53
54    exec_time = particulator.products["wall time"].get()
55    return vals, exec_time
def main(plot: bool, save: Optional[str]):
58def main(plot: bool, save: Optional[str]):
59    with np.errstate(all="raise"):
60        settings = Settings()
61
62        settings.n_sd = 2**15
63
64        states, _ = run(settings)
65
66    with np.errstate(invalid="ignore"):
67        plotter = SpectrumPlotter(settings)
68        plotter.smooth = True
69        for step, vals in states.items():
70            _ = plotter.plot(vals, step * settings.dt)
71            # assert _ < 200  # TODO #327
72        if save is not None:
73            n_sd = settings.n_sd
74            plotter.save(save + "/" + f"{n_sd}_shima_fig_2" + "." + plotter.format)
75        if plot:
76            plotter.show()