Skip to content

Commit

Permalink
fix: adia0D example can run
Browse files Browse the repository at this point in the history
  • Loading branch information
yoctoyotta1024 committed Apr 17, 2024
1 parent 3c56d80 commit 32dfbe5
Show file tree
Hide file tree
Showing 11 changed files with 134 additions and 86 deletions.
30 changes: 14 additions & 16 deletions examples/adiabaticparcel/as2017.py
Original file line number Diff line number Diff line change
Expand Up @@ -8,7 +8,7 @@
Author: Clara Bayley (CB)
Additional Contributors:
-----
Last Modified: Friday 12th April 2024
Last Modified: Thursday 18th April 2024
Modified By: CB
-----
License: BSD 3-Clause "New" or "Revised" License
Expand Down Expand Up @@ -80,30 +80,28 @@

# parameters to edit in model configuration and plotting
params1 = {
"W_AVG": 1,
"T_HALF": 150,
"W_avg": 1,
"TAU_half": 150,
"T_END": 300,
"COUPLTSTEP": 1,
"OBSTSTEP": 2,
"lwdth": 2,
}
params2 = {
"W_AVG": 0.5,
"T_HALF": 300,
"W_avg": 0.5,
"TAU_half": 300,
"T_END": 600,
"COUPLTSTEP": 1,
"OBSTSTEP": 2,
"lwdth": 1,
}
params3 = {
"W_AVG": 0.002,
"T_HALF": 75000,
"W_avg": 0.002,
"TAU_half": 75000,
"T_END": 150000,
"COUPLTSTEP": 3,
"OBSTSTEP": 750,
"lwdth": 0.5,
}
paramslist = [params1, params2, params3]
lwdths = [2, 1, 0.5]

def displacement(time, w_avg, thalf):
'''displacement z given velocity, w, is sinusoidal
Expand Down Expand Up @@ -163,16 +161,16 @@ def displacement(time, w_avg, thalf):
plt.close()

fig, axs = plt.subplots(nrows=3, ncols=1, figsize=(5, 16))
for params in paramslist:
for params, lwdth in zip(paramslist, lwdths):

### edit relevant setup file parameters
params["zarrbasedir"] = binpath+"as2017_sol"+str(runnum)+".zarr"
params["setuptxt"] = binpath+"as2017_setup.txt"
params["setup_filename"] = binpath+"as2017_setup.txt"
editconfigfile.edit_config_params(configfile, params)

### delete any existing dataset
os.system("rm -rf "+params["zarrbasedir"])
os.system("rm "+params["setuptxt"])
os.system("rm "+params["setup_filename"])

### run model
os.chdir(path2build)
Expand All @@ -197,18 +195,18 @@ def displacement(time, w_avg, thalf):
supersat = thermo.supersaturation()
time = pyzarr.get_time(dataset).secs
sddata = pyzarr.get_supers(dataset, consts)
zprof = displacement(time, config["W_AVG"], config["T_HALF"])
zprof = displacement(time, config["W_avg"], config["TAU_half"])

attrs = ["radius", "xi", "msol"]
sd0 = sdtracing.attributes_for1superdroplet(sddata, 0, attrs)
numconc = np.sum(sddata["xi"][0])/gbxs["domainvol"]/1e6 # [/cm^3]

### plot results
wlab = "<w> = {:.1f}".format(config["W_AVG"]*100)+"cm s$^{-1}$"
wlab = "<w> = {:.1f}".format(config["W_avg"]*100)+"cm s$^{-1}$"
axs = as2017fig.condensation_validation_subplots(axs, time, sd0["radius"],
supersat[:, 0, 0, 0],
zprof,
lwdth=params["lwdth"],
lwdth=lwdth,
lab=wlab)

runnum += 1
Expand Down
6 changes: 3 additions & 3 deletions examples/adiabaticparcel/cuspbifurc.py
Original file line number Diff line number Diff line change
Expand Up @@ -6,7 +6,7 @@
Author: Clara Bayley (CB)
Additional Contributors:
-----
Last Modified: Friday 12th April 2024
Last Modified: Thursday 18th April 2024
Modified By: CB
-----
License: BSD 3-Clause "New" or "Revised" License
Expand Down Expand Up @@ -149,7 +149,7 @@ def displacement(time, w_avg, thalf):
supersat = thermo.supersaturation()
time = pyzarr.get_time(dataset).secs
sddata = pyzarr.get_supers(dataset, consts)
zprof = displacement(time, config["W_AVG"], config["T_HALF"])
zprof = displacement(time, config["W_avg"], config["TAU_half"])

### plot results
# sample drops to plot from whole range of SD ids
Expand All @@ -169,5 +169,5 @@ def displacement(time, w_avg, thalf):
thermo.temp[:, 0, 0, 0],
supersat[:, 0, 0, 0],
sddata.IONIC, sddata.MR_SOL,
config["W_AVG"], numconc,
config["W_avg"], numconc,
savename2)
61 changes: 30 additions & 31 deletions examples/adiabaticparcel/src/config/as2017_config.yaml
Original file line number Diff line number Diff line change
@@ -1,4 +1,3 @@
---
# ----- CLEO -----
# File: config.yaml
# Project: config
Expand All @@ -24,52 +23,52 @@

### Initialisation Parameters ###
inputfiles:
constants_filename : ../libs/cleoconstants.hpp # name of file for values of physical constants
grid_filename : ./share/as2017_dimlessGBxboundaries.dat # binary filename for initialisation of GBxs / GbxMaps
constants_filename: ../libs/cleoconstants.hpp # name of file for values of physical constants
grid_filename: ./share/as2017_dimlessGBxboundaries.dat # binary filename for initialisation of GBxs / GbxMaps

initsupers:
type : frombinary # type of initialisation of super-droplets
initsupers_filename : ./share/as2017_dimlessSDsinit.dat # binary filename for initialisation of SDs
totnsupers : 64 # initial total no. of SDs
type: frombinary # type of initialisation of super-droplets
initsupers_filename: ./share/as2017_dimlessSDsinit.dat # binary filename for initialisation of SDs
totnsupers: 64 # initial total no. of SDs

### Output Parameters ###
outputdata:
setup_filename : /home/m/m300950/CLEO/build_adia0D//bin/as2017_setup.txt # .txt filename to copy configuration to
stats_filename : /home/m/m300950/CLEO/build_adia0D//bin/as2017_stats.txt # .txt file to output runtime statistics to
zarrbasedir : /home/m/m300950/CLEO/build_adia0D//bin/as2017_sol8.zarr # zarr store base directory
maxchunk : 2500000 # maximum no. of elements in chunks of zarr store array
setup_filename: /home/m/m300950/CLEO/build_adia0D//bin/as2017_setup.txt # .txt filename to copy configuration to
stats_filename: /home/m/m300950/CLEO/build_adia0D//bin/as2017_stats.txt # .txt file to output runtime statistics to
zarrbasedir: /home/m/m300950/CLEO/build_adia0D//bin/as2017_sol2.zarr # zarr store base directory
maxchunk: 2500000 # maximum no. of elements in chunks of zarr store array

### SDM Runtime Parameters ###
domain:
nspacedims : 0 # no. of spatial dimensions to model
ngbxs : 1 # total number of Gbxs
nspacedims: 0 # no. of spatial dimensions to model
ngbxs: 1 # total number of Gbxs

timesteps:
CONDTSTEP : 1 # time between SD condensation [s]
COLLTSTEP : 1 # time between SD collision [s]
MOTIONTSTEP : 1 # time between SDM motion [s]
COUPLTSTEP : 3 # time between dynamic couplings [s]
OBSTSTEP : 750 # time between SDM observations [s]
T_END : 150000 # time span of integration from 0s to T_END [s]
CONDTSTEP: 1 # time between SD condensation [s]
COLLTSTEP: 1 # time between SD collision [s]
MOTIONTSTEP: 1 # time between SDM motion [s]
COUPLTSTEP: 3 # time between dynamic couplings [s]
OBSTSTEP: 750 # time between SDM observations [s]
T_END: 150000 # time span of integration from 0s to T_END [s]

### Microphysics Parameters ###
microphysics:
condensation:
do_alter_thermo : true # true = cond/evap alters the thermodynamic state
niters : 2 # no. iterations of Newton Raphson Method before testing for convergence
SUBTSTEP : 0.1 # smallest subtimestep in cases of substepping [s]
rtol : 0.0 # relative tolerance for implicit Euler integration
atol : 0.01 # abolute tolerance for implicit Euler integration
do_alter_thermo: true # true = cond/evap alters the thermodynamic state
niters: 2 # no. iterations of Newton Raphson Method before testing for convergence
SUBTSTEP: 0.1 # smallest subtimestep in cases of substepping [s]
rtol: 0.0 # relative tolerance for implicit Euler integration
atol: 0.01 # abolute tolerance for implicit Euler integration

### Coupled Dynamics Parameters ###
coupled_dynamics:
type : cvode # type of coupled dynamics to configure
type: cvode # type of coupled dynamics to configure
# initial (uniform) thermodynamic conditions #
P_init : 100000.0 # initial pressure [Pa]
TEMP_init : 273.15 # initial temperature [T]
relh_init : 98.0 # initial relative humidity (%)
P_init: 100000.0 # initial pressure [Pa]
TEMP_init: 273.15 # initial temperature [T]
relh_init: 98.0 # initial relative humidity (%)
# ODE solver parameters #
W_avg : 0.002 # average amplitude of sinusoidal w [m/s] (dP/dt ~ w*dP/dz)
TAU_half : 75000 # timescale for w sinusoid, tau_half = TAU_HALF/pi [s]
rtol : 1e-6 # relative tolerance for integration of [P, T, qv, qc] ODEs
atol : 1e-6 # relative tolerance for integration of [P, T, qv, qc] ODEs
W_avg: 0.002 # average amplitude of sinusoidal w [m/s] (dP/dt ~ w*dP/dz)
TAU_half: 75000 # timescale for w sinusoid, tau_half = TAU_half/pi [s]
rtol: 1e-6 # relative tolerance for integration of [P, T, qv, qc] ODEs
atol: 1e-6 # relative tolerance for integration of [P, T, qv, qc] ODEs
2 changes: 1 addition & 1 deletion examples/adiabaticparcel/src/config/cuspbifurc_config.yaml
Original file line number Diff line number Diff line change
Expand Up @@ -70,6 +70,6 @@ coupled_dynamics:
relh_init : 98.0 # initial relative humidity (%)
# ODE solver parameters #
W_avg : 0.002 # average amplitude of sinusoidal w [m/s] (dP/dt ~ w*dP/dz)
TAU_half : 75000 # timescale for w sinusoid, tau_half = TAU_HALF/pi [s]
TAU_half : 75000 # timescale for w sinusoid, tau_half = TAU_half/pi [s]
rtol : 1e-6 # relative tolerance for integration of [P, T, qv, qc] ODEs
atol : 1e-6 # relative tolerance for integration of [P, T, qv, qc] ODEs
6 changes: 3 additions & 3 deletions examples/exampleplotting/plotssrc/as2017fig.py
Original file line number Diff line number Diff line change
Expand Up @@ -6,7 +6,7 @@
Author: Clara Bayley (CB)
Additional Contributors:
-----
Last Modified: Friday 17th November 2023
Last Modified: Wednesday 17th April 2024
Modified By: CB
-----
License: BSD 3-Clause "New" or "Revised" License
Expand Down Expand Up @@ -100,7 +100,7 @@ def condensation_validation_subplots(axs, time, radius, supersat, zprof,
return axs

def arabas_shima_2017_fig(time, zprof, radius, msol, temp, supersat,
IONIC, MR_SOL, W_AVG, numconc, savename=""):
IONIC, MR_SOL, W_avg, numconc, savename=""):
''' plots the same plots as in Figure 5 of
"On the CCN (de)activation nonlinearities"
S. Arabas and S. Shima 2017 to check radius
Expand All @@ -118,7 +118,7 @@ def arabas_shima_2017_fig(time, zprof, radius, msol, temp, supersat,

textlab = "N = "+str(numconc)+"cm$^{-3}$\n" +\
"r$_{dry}$ = "+"{:.2g}\u03BCm\n".format(radius[0]) +\
"<w> = {:.1f}".format(W_AVG*100)+"cm s$^{-1}$"
"<w> = {:.1f}".format(W_avg*100)+"cm s$^{-1}$"
axs[0].text(0.03, 0.82, textlab, transform=axs[0].transAxes)

axs[0].set_xlim([-1, 1])
Expand Down
4 changes: 2 additions & 2 deletions libs/initialise/optional_config_params.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -9,7 +9,7 @@
* Author: Clara Bayley (CB)
* Additional Contributors:
* -----
* Last Modified: Wednesday 17th April 2024
* Last Modified: Thursday 18th April 2024
* Modified By: CB
* -----
* License: BSD 3-Clause "New" or "Revised" License
Expand Down Expand Up @@ -100,7 +100,7 @@ struct OptionalConfigParams {
double relh_init = NaNVals::dbl(); /**< initial relative humidity (%) */
/* ODE solver parameters */
double W_avg = NaNVals::dbl(); /**< average amplitude of sinusoidal w [m/s] (dP/dt ~ w*dP/dz) */
double TAU_half = NaNVals::dbl(); /**< timescale for w sinusoid, tau_half = TAU_HALF/pi [s] */
double TAU_half = NaNVals::dbl(); /**< timescale for w sinusoid, tau_half = TAU_half/pi [s] */
double rtol = NaNVals::dbl(); /**< relative tolerance for integration of [P, T, qv, qc] ODEs */
double atol = NaNVals::dbl(); /**< absolute tolerances for integration of [P, T, qv, qc] ODEs */
} cvodedynamics;
Expand Down
19 changes: 10 additions & 9 deletions pySD/editconfigfile.py
Original file line number Diff line number Diff line change
Expand Up @@ -9,7 +9,7 @@
Author: Clara Bayley (CB)
Additional Contributors:
-----
Last Modified: Wednesday 17th April 2024
Last Modified: Thursday 18th April 2024
Modified By: CB
-----
License: BSD 3-Clause "New" or "Revised" License
Expand All @@ -24,19 +24,21 @@ def update_param(node, param, new_value):
'''Function to recursively searches for 'param' key in YAML node
and updates it's value to with 'new_value' when found '''

is_success = False
if isinstance(node, dict):
if param in node:
node[param] = new_value # update value
is_success = True
return True
else:
for key, val in node.items():
update_param(val, param, new_value)
is_success = update_param(val, param, new_value)
if is_success:
return True
elif isinstance(node, list):
for item in node:
update_param(item, param, new_value)

return is_success
is_success = update_param(item, param, new_value)
if is_success:
return True
return False

def edit_config_params(filename, params2change):
''' rewrites config YAML file with key,value pairs listed in params2change updated to new values
Expand All @@ -51,9 +53,8 @@ def edit_config_params(filename, params2change):
# Update the parameters from the YAML file
for param, new_value in params2change.items():
is_success = update_param(data, param, new_value)

if not is_success:
errmsg = param+"could not be updated to new value: "+str(new_value)
errmsg = param+" could not be updated to new value: "+str(new_value)
raise ValueError(errmsg)

# Overwrite the YAML file
Expand Down
5 changes: 2 additions & 3 deletions pySD/gbxboundariesbinary_src/create_gbxboundaries.py
Original file line number Diff line number Diff line change
Expand Up @@ -6,7 +6,7 @@
Author: Clara Bayley (CB)
Additional Contributors:
-----
Last Modified: Tuesday 24th October 2023
Last Modified: Wednesday 17th April 2024
Modified By: CB
-----
License: BSD 3-Clause "New" or "Revised" License
Expand All @@ -23,8 +23,7 @@
from .. import cxx2py, writebinary

def get_COORD0_from_constsfile(constsfile, returnconsts=False):
''' create values from constants file & config file
required as inputs to create initial
''' create values from constants file required as inputs to create initial
superdroplet conditions '''

consts = cxx2py.read_cxxconsts_into_floats(constsfile)
Expand Down
24 changes: 21 additions & 3 deletions pySD/gbxboundariesbinary_src/read_gbxboundaries.py
Original file line number Diff line number Diff line change
@@ -1,3 +1,23 @@
'''
Copyright (c) 2024 MPI-M, Clara Bayley
----- CLEO -----
File: read_gbxboundaries.py
Project: gbxboundariesbinary_src
Created Date: Wednesday 17th January 2024
Author: Clara Bayley (CB)
Additional Contributors:
-----
Last Modified: Wednesday 17th April 2024
Modified By: CB
-----
License: BSD 3-Clause "New" or "Revised" License
https://opensource.org/licenses/BSD-3-Clause
-----
File Description:
'''

import numpy as np
import matplotlib.pyplot as plt

Expand Down Expand Up @@ -248,9 +268,7 @@ def grid_dimensions(gbxbounds):
return extents, spacings, griddims

def print_domain_info(constsfile, gridfile):
''' create values from constants file & config file
required as inputs to create initial
superdroplet conditions '''
''' prints information about domain reda from gridfile and constants file '''

isprint=True
COORD0 = get_COORD0_from_constsfile(constsfile)
Expand Down
4 changes: 2 additions & 2 deletions pySD/sdmout_src/ensembzarr.py
Original file line number Diff line number Diff line change
Expand Up @@ -6,7 +6,7 @@
Author: Clara Bayley (CB)
Additional Contributors:
-----
Last Modified: Friday 19th January 2024
Last Modified: Thursday 18th April 2024
Modified By: CB
-----
License: BSD 3-Clause "New" or "Revised" License
Expand Down Expand Up @@ -54,7 +54,7 @@ def write_ensemb_setupfile(ensembsetupfile, setupfile, datasets):
os.system('cp '+setupfile+" "+ensembsetupfile)
params = {
"initsupers_filename" : "[ensemble, see below]",
"setuptxt" : "[ensemble, see below]",
"setup_filename" : "[ensemble, see below]",
"zarrbasedir" : "[ensemble, see below]",
"stats_filename" : "[ensemble, see below]"
}
Expand Down
Loading

0 comments on commit 32dfbe5

Please sign in to comment.