PySDM_examples.Shima_et_al_2009.example

View Source

 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.builder import Builder
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(dv=settings.dv, dt=settings.dt)
18    builder = Builder(
19        n_sd=settings.n_sd,
20        backend=backend(formulae=settings.formulae),
21        environment=env,
22        dynamics=(
23            Coalescence(collision_kernel=settings.kernel, adaptive=settings.adaptive),
24        ),
25    )
26    attributes = {}
27    sampling = ConstantMultiplicity(settings.spectrum)
28    attributes["volume"], attributes["multiplicity"] = sampling.sample_deterministic(
29        settings.n_sd
30    )
31    products = (
32        ParticleVolumeVersusRadiusLogarithmSpectrum(
33            settings.radius_bins_edges, name="dv/dlnr"
34        ),
35        WallTime(),
36    )
37    particulator = builder.build(attributes, products)
38
39    for observer in observers:
40        particulator.observers.append(observer)
41
42    vals = {}
43    particulator.products["wall time"].reset()
44    for step in settings.output_steps:
45        particulator.run(step - particulator.n_steps)
46        vals[step] = particulator.products["dv/dlnr"].get()[0]
47        vals[step][:] *= settings.rho
48
49    exec_time = particulator.products["wall time"].get()
50    return vals, exec_time
51
52
53def main(plot: bool, save: Optional[str]):
54    with np.errstate(all="raise"):
55        settings = Settings()
56
57        settings.n_sd = 2**15
58
59        states, _ = run(settings)
60
61    with np.errstate(invalid="ignore"):
62        plotter = SpectrumPlotter(settings)
63        plotter.smooth = True
64        for step, vals in states.items():
65            _ = plotter.plot(vals, step * settings.dt)
66            # assert _ < 200  # TODO #327
67        if save is not None:
68            n_sd = settings.n_sd
69            plotter.save(save + "/" + f"{n_sd}_shima_fig_2" + "." + plotter.format)
70        if plot:
71            plotter.show()
72
73
74if __name__ == "__main__":
75    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=()): View Source

17def run(settings, backend=CPU, observers=()):
18    env = Box(dv=settings.dv, dt=settings.dt)
19    builder = Builder(
20        n_sd=settings.n_sd,
21        backend=backend(formulae=settings.formulae),
22        environment=env,
23        dynamics=(
24            Coalescence(collision_kernel=settings.kernel, adaptive=settings.adaptive),
25        ),
26    )
27    attributes = {}
28    sampling = ConstantMultiplicity(settings.spectrum)
29    attributes["volume"], attributes["multiplicity"] = sampling.sample_deterministic(
30        settings.n_sd
31    )
32    products = (
33        ParticleVolumeVersusRadiusLogarithmSpectrum(
34            settings.radius_bins_edges, name="dv/dlnr"
35        ),
36        WallTime(),
37    )
38    particulator = builder.build(attributes, products)
39
40    for observer in observers:
41        particulator.observers.append(observer)
42
43    vals = {}
44    particulator.products["wall time"].reset()
45    for step in settings.output_steps:
46        particulator.run(step - particulator.n_steps)
47        vals[step] = particulator.products["dv/dlnr"].get()[0]
48        vals[step][:] *= settings.rho
49
50    exec_time = particulator.products["wall time"].get()
51    return vals, exec_time

def main(plot: bool, save: Optional[str]): View Source

54def main(plot: bool, save: Optional[str]):
55    with np.errstate(all="raise"):
56        settings = Settings()
57
58        settings.n_sd = 2**15
59
60        states, _ = run(settings)
61
62    with np.errstate(invalid="ignore"):
63        plotter = SpectrumPlotter(settings)
64        plotter.smooth = True
65        for step, vals in states.items():
66            _ = plotter.plot(vals, step * settings.dt)
67            # assert _ < 200  # TODO #327
68        if save is not None:
69            n_sd = settings.n_sd
70            plotter.save(save + "/" + f"{n_sd}_shima_fig_2" + "." + plotter.format)
71        if plot:
72            plotter.show()