diff --git a/.github/workflows/main.yml b/.github/workflows/main.yml index 0f576df1..8eb621e9 100644 --- a/.github/workflows/main.yml +++ b/.github/workflows/main.yml @@ -1,4 +1,4 @@ -name: Build Status +name: PySDM-examples defaults: run: @@ -27,11 +27,11 @@ jobs: python-version: ${{ matrix.python-version }} - run: | - pip install pytest nbconvert + pip install -r test-time-requirements.txt - run: | pip install -e . - # TODO #184 (https://github.com/numba/numba/issues/6350#issuecomment-728174860) + # https://github.com/numba/numba/issues/6350#issuecomment-728174860 - if: matrix.platform == 'ubuntu-latest' env: NUMBA_THREADING_LAYER: omp diff --git a/PySDM_examples/Arabas_and_Shima_2017/settings.py b/PySDM_examples/Arabas_and_Shima_2017/settings.py index b06abbb4..e71dbcc6 100644 --- a/PySDM_examples/Arabas_and_Shima_2017/settings.py +++ b/PySDM_examples/Arabas_and_Shima_2017/settings.py @@ -5,20 +5,24 @@ from PySDM.physics.constants import si from PySDM.physics import constants as const, formulae as phys from PySDM.dynamics import condensation +from PySDM.physics.formulae import Formulae import numpy as np from pystrict import strict @strict class Settings: - def __init__(self, w_avg: float, N_STP: float, r_dry: float, mass_of_dry_air: float): + def __init__(self, w_avg: float, N_STP: float, r_dry: float, mass_of_dry_air: float, coord: str = 'VolumeLogarithm'): + self.formulae = Formulae(saturation_vapour_pressure='AugustRocheMagnus', condensation_coordinate=coord) + self.p0 = 1000 * si.hectopascals self.RH0 = .98 self.kappa = .2 # TODO #441 self.T0 = 300 * si.kelvin self.z_half = 150 * si.metres - self.q0 = const.eps / (self.p0 / self.RH0 / phys.pvs(self.T0) - 1) + pvs = self.formulae.saturation_vapour_pressure.pvs_Celsius(self.T0 - const.T0) + self.q0 = const.eps / (self.p0 / self.RH0 / pvs - 1) self.w_avg = w_avg self.r_dry = r_dry self.N_STP = N_STP @@ -29,7 +33,7 @@ def __init__(self, w_avg: float, N_STP: float, r_dry: float, mass_of_dry_air: fl self.rtol_x = condensation.default_rtol_x self.rtol_thd = condensation.default_rtol_thd self.coord = 'volume logarithm' - self.dt_cond_range = (1e-10 * si.s, 1 * si.s) + self.dt_cond_range = condensation.default_cond_range @property def dt_max(self): diff --git a/PySDM_examples/Arabas_and_Shima_2017/simulation.py b/PySDM_examples/Arabas_and_Shima_2017/simulation.py index f2ccd229..86a5809c 100644 --- a/PySDM_examples/Arabas_and_Shima_2017/simulation.py +++ b/PySDM_examples/Arabas_and_Shima_2017/simulation.py @@ -24,7 +24,7 @@ def __init__(self, settings, backend=CPU): while dt_output / self.n_substeps >= settings.dt_max: # TODO #334 dt_max self.n_substeps += 1 - builder = Builder(backend=backend, n_sd=1) + builder = Builder(backend=backend, n_sd=1, formulae=settings.formulae) builder.set_environment(Parcel( dt=dt_output / self.n_substeps, mass_of_dry_air=settings.mass_of_dry_air, @@ -43,11 +43,11 @@ def __init__(self, settings, backend=CPU): )) attributes = {} r_dry = np.array([settings.r_dry]) - attributes['dry volume'] = phys.volume(radius=r_dry) + attributes['dry volume'] = settings.formulae.trivia.volume(radius=r_dry) attributes['n'] = np.array([settings.n_in_dv], dtype=np.int64) environment = builder.core.environment - r_wet = r_wet_init(r_dry, environment, np.zeros_like(attributes['n']), settings.kappa) - attributes['volume'] = phys.volume(radius=r_wet) + r_wet = r_wet_init(r_dry, environment, kappa=settings.kappa) + attributes['volume'] = settings.formulae.trivia.volume(radius=r_wet) products = [ PySDM_products.ParticleMeanRadius(), PySDM_products.CondensationTimestepMin(), @@ -56,7 +56,8 @@ def __init__(self, settings, backend=CPU): PySDM_products.Time(), PySDM_products.ActivatingRate(), PySDM_products.DeactivatingRate(), - PySDM_products.RipeningRate() + PySDM_products.RipeningRate(), + PySDM_products.PeakSupersaturation() ] self.core = builder.build(attributes, products) @@ -67,7 +68,7 @@ def save(self, output): cell_id = 0 output["r"].append(self.core.products['radius_m1'].get(unit=const.si.metre)[cell_id]) output["dt_cond_min"].append(self.core.products['dt_cond_min'].get()[cell_id]) - output["z"].append(self.core.products["z"].get()) + output["z"].append(self.core.products["z"].get()[cell_id]) output["S"].append(self.core.products["RH_env"].get()[cell_id]/100 - 1) output["t"].append(self.core.products["t"].get()) diff --git a/PySDM_examples/Arabas_et_al_2015/demo_plots.py b/PySDM_examples/Arabas_et_al_2015/demo_plots.py index c769a4bb..529a2f9b 100644 --- a/PySDM_examples/Arabas_et_al_2015/demo_plots.py +++ b/PySDM_examples/Arabas_et_al_2015/demo_plots.py @@ -5,7 +5,7 @@ import matplotlib import matplotlib.pyplot as plt import numpy as np - +from PySDM import products class _Plot: @@ -37,9 +37,36 @@ def __init__(self, fig, ax, grid, size, product, show=False, lines=False, cmap=' self.lines['X'][1] = plt.plot([-1] * 2, zlim, **self.line_args)[0] self.lines['Z'][1] = plt.plot(xlim, [-1] * 2, **self.line_args)[0] - data = np.full_like(self.nans, product.range[0]) + product_range, product_scale = { + products.WaterMixingRatio: ((0, 1), 'linear'), + products.TotalParticleSpecificConcentration: ((20, 50), 'linear'), + products.TotalParticleConcentration: ((20, 50), 'linear'), + products.SuperDropletCount: ((0, 10), 'linear'), + products.ParticlesVolumeSpectrum: ((20, 50), 'linear'), + products.AerosolConcentration: ((1e0, 1e2), 'linear'), + products.CloudDropletConcentration: ((1e0, 1e2), 'linear'), + products.DrizzleConcentration: ((1e-3, 1e1), 'log'), + products.ParticleMeanRadius: ((1, 25), 'linear'), + products.CloudDropletEffectiveRadius: ((0, 25), 'linear'), + products.AerosolSpecificConcentration: ((0, 3e2), 'linear'), + products.WaterVapourMixingRatio: ((5, 7.5), 'linear'), + products.Temperature: ((275, 305), 'linear'), + products.RelativeHumidity: ((75, 105), 'linear'), + products.Pressure: ((90000, 100000), 'linear'), + products.DryAirPotentialTemperature: ((275, 300), 'linear'), + products.DryAirDensity: ((0.95, 1.3), 'linear'), + products.PeakSupersaturation: ((-1, 1), 'linear'), + products.RipeningRate: ((1e-1, 1e1), 'log'), + products.ActivatingRate: ((1e-1, 1e1), 'log'), + products.DeactivatingRate: ((1e-1, 1e1), 'log'), + products.CollisionRateDeficit: ((0, 1e10), 'linear'), + products.CollisionRate: ((0, 1e10), 'linear'), + products.CondensationTimestepMin: ((.01, 10), 'log'), + products.CondensationTimestepMax: ((.01, 10), 'log'), + products.CoalescenceTimestepMin: ((.01, 10), 'log'), + }[product.__class__] + data = np.full_like(self.nans, product_range[0]) label = f"{product.description} [{product.unit}]" - scale = product.scale self.ax.set_xlabel('X [m]') self.ax.set_ylabel('Z [m]') @@ -48,11 +75,11 @@ def __init__(self, fig, ax, grid, size, product, show=False, lines=False, cmap=' origin='lower', extent=(*xlim, *zlim), cmap=cmap, - norm=matplotlib.colors.LogNorm() if scale == 'log' and np.isfinite( + norm=matplotlib.colors.LogNorm() if product_scale == 'log' and np.isfinite( data).all() else None # TODO #37 this is always None!!! ) plt.colorbar(self.im, ax=self.ax).set_label(label) - self.im.set_clim(vmin=product.range[0], vmax=product.range[1]) + self.im.set_clim(vmin=product_range[0], vmax=product_range[1]) if show: plt.show() @@ -117,5 +144,5 @@ def update(self, data): if data is not None: self.ydata[0:len(data)] = data[:] else: - self.ydata[0:len(data)] = np.nan + self.ydata[:] = np.nan self.timeseries.set_ydata(self.ydata) diff --git a/PySDM_examples/Arabas_et_al_2015/demo_settings.py b/PySDM_examples/Arabas_et_al_2015/demo_settings.py index 01f2f52a..df421336 100644 --- a/PySDM_examples/Arabas_et_al_2015/demo_settings.py +++ b/PySDM_examples/Arabas_et_al_2015/demo_settings.py @@ -4,9 +4,11 @@ from ..utils.widgets import IntSlider, FloatSlider, VBox, Checkbox, Accordion, Dropdown from PySDM_examples.Arabas_et_al_2015.settings import Settings +from PySDM import physics import numpy as np import numba import os +import inspect class DemoSettings: @@ -14,6 +16,7 @@ class DemoSettings: def __init__(self): settings = Settings() + self.ui_th_std0 = FloatSlider(description="$\\theta_0$ [K]", value=settings.th_std0, min=280, max=300) self.ui_qv0 = FloatSlider(description="q$_{v0}$ [g/kg]", value=settings.qv0 * 1000, min=5, max=10) self.ui_p0 = FloatSlider(description="p$_0$ [hPa]", value=settings.p0 / 100, min=900, max=1100) @@ -43,10 +46,6 @@ def __init__(self): self.ui_coalescence_adaptive = Checkbox( value=settings.condensation_adaptive, description='coalescence adaptive time-step') - self.ui_condensation_coord = Dropdown( - options=['volume', 'volume logarithm'], - value=settings.condensation_coord, - description='condensational variable coordinate') self.ui_processes = [Checkbox(value=settings.processes[key], description=key) for key in settings.processes.keys()] self.ui_sdpg = IntSlider(value=settings.n_sd_per_gridbox, description="n_sd/gridbox", min=1, max=1000) self.fct_description = "MPDATA: flux-corrected transport option" @@ -59,27 +58,65 @@ def __init__(self): Checkbox(value=settings.mpdata_iga, description=self.iga_description), IntSlider(value=settings.mpdata_iters, description=self.nit_description, min=1, max=5) ] + self.ui_formulae_options = [ + Dropdown( + description=k, + options=physics.formulae._choices(getattr(physics, k)).keys(), + value=v.default + ) + for k, v in inspect.signature(physics.Formulae.__init__).parameters.items() + if k not in ('self', 'fastmath', 'seed') + ] + self.ui_formulae_options.append( + Checkbox( + value=inspect.signature(physics.Formulae.__init__).parameters['fastmath'].default, + description='fastmath' + ) + ) + self.ui_output_options = { + 'interval': IntSlider(description='interval [s]', value=settings.output_interval, min=30, max=60*15), + 'aerosol_radius_threshold': FloatSlider( + description='aerosol/cloud threshold [um]', + value=settings.aerosol_radius_threshold / physics.si.um, + min=.1, max=1, step=.1), + 'drizzle_radius_threshold': IntSlider( + description='cloud/drizzle threshold [um]', + value=settings.drizzle_radius_threshold / physics.si.um, + min=10, max=100, step=5) + } # TODO #37 - self.v_bins = settings.v_bins - self.n_sd = settings.n_sd + self.r_bins_edges = settings.r_bins_edges self.size = settings.size - self.rhod = settings.rhod - self.field_values = settings.field_values - self.aerosol_radius_threshold = settings.aerosol_radius_threshold - self.drizzle_radius_threshold = settings.drizzle_radius_threshold - self.stream_function = settings.stream_function self.condensation_substeps = settings.condensation_substeps self.condensation_dt_cond_range = settings.condensation_dt_cond_range self.condensation_schedule = settings.condensation_schedule self.kernel = settings.kernel self.spectrum_per_mass_of_dry_air = settings.spectrum_per_mass_of_dry_air - self.n_spin_up = settings.n_spin_up self.coalescence_dt_coal_range = settings.coalescence_dt_coal_range self.coalescence_optimized_random = settings.coalescence_optimized_random self.coalescence_substeps = settings.coalescence_substeps - self.output_interval = settings.output_interval - self.versions = settings.versions + + for attr in ('n_sd', 'rhod', 'field_values', 'versions', 'n_spin_up', 'stream_function'): + setattr(self, attr, getattr(settings, attr)) + + @property + def aerosol_radius_threshold(self): + return self.ui_output_options['aerosol_radius_threshold'].value * physics.si.um + + @property + def drizzle_radius_threshold(self): + return self.ui_output_options['drizzle_radius_threshold'].value * physics.si.um + + @property + def output_interval(self): + return self.ui_output_options['interval'].value + + @property + def formulae(self) -> physics.Formulae: + return physics.Formulae( + **{widget.description: widget.value for widget in self.ui_formulae_options} + ) @property def steps_per_output_interval(self) -> int: @@ -137,10 +174,6 @@ def condensation_adaptive(self): def coalescence_adaptive(self): return self.ui_coalescence_adaptive.value - @property - def condensation_coord(self): - return self.ui_condensation_coord.value - @property def processes(self): result = {} @@ -186,10 +219,14 @@ def box(self): VBox([*self.ui_processes]), VBox([self.ui_nx, self.ui_nz, self.ui_sdpg, self.ui_dt, self.ui_simulation_time, self.ui_condensation_rtol_x, self.ui_condensation_rtol_thd, - self.ui_condensation_adaptive, self.ui_coalescence_adaptive, self.ui_condensation_coord, + self.ui_condensation_adaptive, self.ui_coalescence_adaptive, *self.ui_mpdata_options]), + VBox([*self.ui_formulae_options]), + VBox([*self.ui_output_options.values()]) ]) - layout.set_title(0, 'environment parameters') + layout.set_title(0, 'environment') layout.set_title(1, 'processes') layout.set_title(2, 'discretisation') + layout.set_title(3, 'physics') + layout.set_title(4, 'output') return layout diff --git a/PySDM_examples/Arabas_et_al_2015/demo_viewer.py b/PySDM_examples/Arabas_et_al_2015/demo_viewer.py index 019a2823..3706fab6 100644 --- a/PySDM_examples/Arabas_et_al_2015/demo_viewer.py +++ b/PySDM_examples/Arabas_et_al_2015/demo_viewer.py @@ -43,7 +43,7 @@ def reinit(self, products): if len(val.shape) == 2 ] - r_bins = phys.radius(volume=self.settings.v_bins) + r_bins = self.settings.r_bins_edges.copy() const.convert_to(r_bins, const.si.micrometres) self.spectrumOutput = Output() with self.spectrumOutput: diff --git a/PySDM_examples/Arabas_et_al_2015/mpdata.py b/PySDM_examples/Arabas_et_al_2015/mpdata.py index 3a448f5c..8c978dc5 100644 --- a/PySDM_examples/Arabas_et_al_2015/mpdata.py +++ b/PySDM_examples/Arabas_et_al_2015/mpdata.py @@ -63,8 +63,5 @@ def wait(self): self.thread.join() def step(self): - try: # TODO #417: move this to within PyMPDATA - for mpdata in self.mpdatas.values(): - mpdata.advance(1) - except NumbaExperimentalFeatureWarning: - pass + for mpdata in self.mpdatas.values(): + mpdata.advance(1) diff --git a/PySDM_examples/Arabas_et_al_2015/netcdf_exporter.py b/PySDM_examples/Arabas_et_al_2015/netcdf_exporter.py index 12452bd0..396de7e0 100644 --- a/PySDM_examples/Arabas_et_al_2015/netcdf_exporter.py +++ b/PySDM_examples/Arabas_et_al_2015/netcdf_exporter.py @@ -20,7 +20,8 @@ def _create_dimensions(self, ncdf): ncdf.createDimension("T", len(self.settings.output_steps)) for index, label in enumerate(self.XZ): ncdf.createDimension(label, self.settings.grid[index]) - ncdf.createDimension("ParticleVolume", len(self.settings.v_bins) - 1) + ncdf.createDimension("ParticleVolume", + len(self.settings.formulae.trivia.volume(self.settings.r_bins_edges)) - 1) def _create_variables(self, ncdf): self.vars["T"] = ncdf.createVariable("T", "f", ["T"]) @@ -28,7 +29,8 @@ def _create_variables(self, ncdf): for index, label in enumerate(self.XZ): self.vars[label] = ncdf.createVariable(label, "f", (label,)) - self.vars[label][:] = (self.settings.size[index] / self.settings.grid[index]) * (1 / 2 + np.arange(self.settings.grid[index])) + self.vars[label][:] = ((self.settings.size[index] / self.settings.grid[index]) + * (1 / 2 + np.arange(self.settings.grid[index]))) self.vars[label].units = "metres" # TODO #340 ParticleVolume var diff --git a/PySDM_examples/Arabas_et_al_2015/settings.py b/PySDM_examples/Arabas_et_al_2015/settings.py index 1fc34506..0bfef588 100644 --- a/PySDM_examples/Arabas_et_al_2015/settings.py +++ b/PySDM_examples/Arabas_et_al_2015/settings.py @@ -1,38 +1,31 @@ -""" -Created at 25.09.2019 -""" - from typing import Iterable - import numba -import numpy -import numpy as np +import numpy, numpy as np import scipy from pystrict import strict - -import PySDM -from PySDM.dynamics import condensation -from PySDM.dynamics.coalescence import coalescence -from PySDM.dynamics.coalescence.kernels import Geometric -from PySDM.initialisation.spectra import Lognormal -from PySDM.initialisation.spectra import Sum -from PySDM.physics import constants as const -from PySDM.physics import formulae as phys -from PySDM.physics.constants import si +import PySDM, PyMPDATA +from PySDM.dynamics import condensation, coalescence +from PySDM.physics.coalescence_kernels import Geometric +from PySDM.physics import spectra, si, Formulae, constants as const +from PySDM.backends.numba.conf import JIT_FLAGS -# from PyMPDATA import __version__ as TODO #339 - @strict class Settings: def __dir__(self) -> Iterable[str]: - return 'dt', 'grid', 'size', 'n_spin_up', 'versions', 'steps_per_output_interval' + return 'dt', 'grid', 'size', 'n_spin_up', 'versions', 'steps_per_output_interval', 'formulae' - def __init__(self): - key_packages = (PySDM, numba, numpy, scipy) - self.versions = str({pkg.__name__: pkg.__version__ for pkg in key_packages}) + def __init__(self, fastmath: bool = JIT_FLAGS['fastmath']): + key_packages = (PySDM, PyMPDATA, numba, numpy, scipy) + self.versions = {} + for pkg in key_packages: + try: + self.versions[pkg.__name__] = pkg.__version__ + except AttributeError: + pass + self.versions = str(self.versions) - self.condensation_coord = 'volume logarithm' + self.formulae = Formulae(condensation_coordinate='VolumeLogarithm', fastmath=fastmath) self.condensation_rtol_x = condensation.default_rtol_x self.condensation_rtol_thd = condensation.default_rtol_thd @@ -57,19 +50,21 @@ def __init__(self): self.dt = 5 * si.second self.spin_up_time = 1 * si.hour - self.v_bins = phys.volume(np.logspace(np.log10(0.001 * si.micrometre), np.log10(100 * si.micrometre), 101, endpoint=True)) + self.r_bins_edges = np.logspace(np.log10(0.001 * si.micrometre), + np.log10(100 * si.micrometre), + 101, endpoint=True) - self.mode_1 = Lognormal( + self.mode_1 = spectra.Lognormal( norm_factor=60 / si.centimetre ** 3 / const.rho_STP, m_mode=0.04 * si.micrometre, s_geom=1.4 ) - self.mode_2 = Lognormal( + self.mode_2 = spectra.Lognormal( norm_factor=40 / si.centimetre**3 / const.rho_STP, m_mode=0.15 * si.micrometre, s_geom=1.6 ) - self.spectrum_per_mass_of_dry_air = Sum((self.mode_1, self.mode_2)) + self.spectrum_per_mass_of_dry_air = spectra.Sum((self.mode_1, self.mode_2)) self.processes = { "particle advection": True, @@ -80,8 +75,6 @@ def __init__(self): # "relaxation": False # TODO #338 } - self.enable_particle_temperatures = False - self.mpdata_iters = 2 self.mpdata_iga = True self.mpdata_fct = True @@ -115,7 +108,7 @@ def output_steps(self) -> np.ndarray: @property def field_values(self): return { - 'th': phys.th_dry(self.th_std0, self.qv0), + 'th': self.formulae.state_variable_triplet.th_dry(self.th_std0, self.qv0), 'qv': self.qv0 } @@ -128,6 +121,6 @@ def stream_function(self, xX, zZ): return - self.rho_w_max * X / np.pi * np.sin(np.pi * zZ) * np.cos(2 * np.pi * xX) def rhod(self, zZ): - p = phys.Hydrostatic.p_of_z_assuming_const_th_and_qv(self.g, self.p0, self.th_std0, self.qv0, z=zZ * self.size[-1]) - rhod = phys.ThStd.rho_d(p, self.qv0, self.th_std0) + p = self.formulae.hydrostatics.p_of_z_assuming_const_th_and_qv(self.g, self.p0, self.th_std0, self.qv0, z=zZ * self.size[-1]) + rhod = self.formulae.state_variable_triplet.rho_d(p, self.qv0, self.th_std0) return rhod diff --git a/PySDM_examples/Arabas_et_al_2015/simulation.py b/PySDM_examples/Arabas_et_al_2015/simulation.py index 8d6f76ac..97420427 100644 --- a/PySDM_examples/Arabas_et_al_2015/simulation.py +++ b/PySDM_examples/Arabas_et_al_2015/simulation.py @@ -32,8 +32,7 @@ def products(self): return self.core.products def reinit(self, products=None): - - builder = Builder(n_sd=self.settings.n_sd, backend=self.backend) + builder = Builder(n_sd=self.settings.n_sd, backend=self.backend, formulae=self.settings.formulae) environment = Kinematic2D(dt=self.settings.dt, grid=self.settings.grid, size=self.settings.size, @@ -45,8 +44,9 @@ def reinit(self, products=None): if products is not None: products = list(products) products = products or [ - PySDM_products.ParticlesWetSizeSpectrum(v_bins=self.settings.v_bins, normalise_by_dv=True), - PySDM_products.ParticlesDrySizeSpectrum(v_bins=self.settings.v_bins, normalise_by_dv=True), # Note: better v_bins + # Note: consider better radius_bins_edges + PySDM_products.ParticlesWetSizeSpectrum(radius_bins_edges=self.settings.r_bins_edges, normalise_by_dv=True), + PySDM_products.ParticlesDrySizeSpectrum(radius_bins_edges=self.settings.r_bins_edges, normalise_by_dv=True), PySDM_products.TotalParticleConcentration(), PySDM_products.TotalParticleSpecificConcentration(), PySDM_products.AerosolConcentration(radius_threshold=self.settings.aerosol_radius_threshold), @@ -78,7 +78,6 @@ def reinit(self, products=None): kappa=self.settings.kappa, rtol_x=self.settings.condensation_rtol_x, rtol_thd=self.settings.condensation_rtol_thd, - coord=self.settings.condensation_coord, adaptive=self.settings.condensation_adaptive, substeps=self.settings.condensation_substeps, dt_cond_range=self.settings.condensation_dt_cond_range, @@ -99,7 +98,6 @@ def reinit(self, products=None): if self.settings.processes["particle advection"]: displacement = Displacement( courant_field=fields.courant_field, - scheme='FTBS', enable_sedimentation=self.settings.processes["sedimentation"]) builder.add_dynamic(displacement) products.append(PySDM_products.SurfacePrecipitation()) # TODO #37 ditto @@ -127,7 +125,7 @@ def reinit(self, products=None): if self.storage is not None: self.storage.init(self.settings) - def run(self, controller=DummyController()): + def run(self, controller=DummyController(), vtk_exporter=None): with controller: for step in self.settings.output_steps: if controller.panic: @@ -136,6 +134,9 @@ def run(self, controller=DummyController()): self.core.run(step - self.core.n_steps) self.store(step) + + if vtk_exporter is not None: + vtk_exporter.export_particles(self.core) controller.set_percent(step / self.settings.output_steps[-1]) diff --git a/PySDM_examples/Arabas_et_al_2015/spin_up.py b/PySDM_examples/Arabas_et_al_2015/spin_up.py index 5872563c..db1ff3d5 100644 --- a/PySDM_examples/Arabas_et_al_2015/spin_up.py +++ b/PySDM_examples/Arabas_et_al_2015/spin_up.py @@ -25,6 +25,5 @@ def set(self, dynamic, attr, value): key = dynamic.__name__ if key in self.particles.dynamics: setattr(self.particles.dynamics[key], attr, value) - # TODO #292 log it else: - warnings.warn(f"{key} not found!") \ No newline at end of file + warnings.warn(f"{key} not found!") diff --git a/PySDM_examples/Bartman_2020_MasterThesis/fig_5_BDF_VS_ADAPTIVE.py b/PySDM_examples/Bartman_2020_MasterThesis/fig_5_BDF_VS_ADAPTIVE.py index 72bef6a0..f1250569 100644 --- a/PySDM_examples/Bartman_2020_MasterThesis/fig_5_BDF_VS_ADAPTIVE.py +++ b/PySDM_examples/Bartman_2020_MasterThesis/fig_5_BDF_VS_ADAPTIVE.py @@ -4,8 +4,8 @@ from PySDM_examples.Arabas_and_Shima_2017.simulation import Simulation from PySDM_examples.Arabas_and_Shima_2017.settings import setups -from PySDM.backends.numba import bdf - +from PySDM.backends.numba.test_helpers import bdf +from PySDM.backends import CPU, GPU import numpy as np import matplotlib.pyplot as plt import matplotlib @@ -23,7 +23,7 @@ def data(n_output, rtols, schemes, setups_num): settings = setups[settings_idx] settings.n_output = n_output simulation = Simulation(settings) - bdf.patch_core(simulation.core, settings.coord) + bdf.patch_core(simulation.core) results = simulation.run() for rtol in rtols: resultant_data[scheme][rtol].append(results) @@ -35,7 +35,7 @@ def data(n_output, rtols, schemes, setups_num): settings.rtol_x = rtol settings.rtol_thd = rtol settings.n_output = n_output - simulation = Simulation(settings) + simulation = Simulation(settings, backend=CPU if scheme=='CPU' else GPU) results = simulation.run() resultant_data[scheme][rtol].append(results) return resultant_data @@ -98,7 +98,7 @@ def plot(data, rtols, schemes, setups_num, path=None): def main(save=None): rtols = [1e-7, 1e-11] - schemes = ['default', 'BDF'] + schemes = ['CPU', 'BDF'] setups_num = len(setups) input_data = data(80, rtols, schemes, setups_num) plot(input_data, rtols, schemes, setups_num, save) diff --git a/PySDM_examples/Bartman_et_al_2021/demo.ipynb b/PySDM_examples/Bartman_et_al_2021/demo.ipynb index 4bf338b7..1ce9bed0 100644 --- a/PySDM_examples/Bartman_et_al_2021/demo.ipynb +++ b/PySDM_examples/Bartman_et_al_2021/demo.ipynb @@ -14,14 +14,14 @@ }, { "cell_type": "code", - "execution_count": 3, + "execution_count": 2, "metadata": {}, "outputs": [ { "name": "stderr", "output_type": "stream", "text": [ - "/home/slayoo/devel/PySDM/PySDM_examples/Bartman_et_al_2021/../../PySDM/backends/__init__.py:29: UserWarning: CUDA library found but cuInit() failed (error code: 999; message: unknown error)\n", + "/home/slayoo/devel/PySDM/PySDM/backends/__init__.py:29: UserWarning: CUDA library found but cuInit() failed (error code: 999; message: unknown error)\n", " warnings.warn(\n" ] } @@ -49,7 +49,7 @@ }, { "cell_type": "code", - "execution_count": 4, + "execution_count": 3, "metadata": {}, "outputs": [], "source": [ @@ -66,7 +66,7 @@ "settings.n_sd_per_gridbox = 128 if 'CI' not in os.environ else 32\n", "settings.grid = (32, 32)\n", "settings.dt = 32 * si.second\n", - "settings.simulation_time = 1.75 * settings.spin_up_time\n", + "settings.simulation_time = .175 * settings.spin_up_time\n", "settings.output_interval = 1 * si.minute\n", "settings.condensation_rtol_x = 1e-6\n", "settings.condensation_rtol_thd = 5e-7\n", @@ -85,13 +85,13 @@ }, { "cell_type": "code", - "execution_count": 5, + "execution_count": null, "metadata": {}, "outputs": [ { "data": { "application/vnd.jupyter.widget-view+json": { - "model_id": "d3674df536034ff29822a815844bfc81", + "model_id": "9112d3fa87404d9998ba6aa91f4c0913", "version_major": 2, "version_minor": 0 }, @@ -109,24 +109,9 @@ }, { "cell_type": "code", - "execution_count": 6, + "execution_count": null, "metadata": {}, - "outputs": [ - { - "data": { - "application/vnd.jupyter.widget-view+json": { - "model_id": "101f71a13f394040955d93dde8904bad", - "version_major": 2, - "version_minor": 0 - }, - "text/plain": [ - "FloatProgress(value=0.0, max=1.0)" - ] - }, - "metadata": {}, - "output_type": "display_data" - } - ], + "outputs": [], "source": [ "temp_file = TemporaryFile('.nc')\n", "exporter = NetCDFExporter(storage, settings, simulation, temp_file.absolute_path)\n", @@ -135,19 +120,9 @@ }, { "cell_type": "code", - "execution_count": 7, + "execution_count": null, "metadata": {}, - "outputs": [ - { - "name": "stderr", - "output_type": "stream", - "text": [ - ":43: DeprecationWarning: `np.int` is a deprecated alias for the builtin `int`. To silence this warning, use `int` by itself. Doing this will not modify any behavior and is safe. When replacing `np.int`, you may wish to use e.g. `np.int64` or `np.int32` to specify the precision. If you wish to review your current use, check the release note link for additional information.\n", - "Deprecated in NumPy 1.20; for more details and guidance: https://numpy.org/devdocs/release/1.20.0-notes.html#deprecations\n", - " frame_list = np.arange(ncdf.variables['T'].shape[0], dtype=np.int)\n" - ] - } - ], + "outputs": [], "source": [ "default_figsize = rcParams[\"figure.figsize\"]\n", "figsize = (1.75 * default_figsize[0], 3.1* default_figsize[1])\n", @@ -191,7 +166,7 @@ "plots[-1].ax.set_ylim(0, .001)\n", "\n", "interval = 100 #ms\n", - "frame_list = np.arange(ncdf.variables['T'].shape[0], dtype=np.int)\n", + "frame_list = np.arange(ncdf.variables['T'].shape[0], dtype=int)\n", "\n", "def update(frame_num):\n", " step = frame_num*ncdf.steps_per_output_interval\n", @@ -213,54 +188,9 @@ }, { "cell_type": "code", - "execution_count": 8, + "execution_count": null, "metadata": {}, - "outputs": [ - { - "name": "stderr", - "output_type": "stream", - "text": [ - "/home/slayoo/devel/PySDM/PySDM_examples/Bartman_et_al_2021/../../PySDM_examples/Arabas_et_al_2015/demo_plots.py:68: RuntimeWarning: All-NaN slice encountered\n", - " self.ax.set_title(f\"min:{np.nanmin(data): .3g} max:{np.nanmax(data): .3g} t/dt:{step: >6}\")\n" - ] - }, - { - "data": { - "application/vnd.jupyter.widget-view+json": { - "model_id": "ac485dedc4b04dd4b18602b7342ad2d2", - "version_major": 2, - "version_minor": 0 - }, - "text/plain": [ - "HTML(value='