API
Model Interface
ClimaCalibrate.set_up_forward_model
— Functionset_up_forward_model(member, iteration, experiment_dir::AbstractString)
+set_up_forward_model(member, iteration, experiment_config::ExperimentConfig)
Set up and configure a single member's forward model. Used in conjunction with run_forward_model
.
This function must be overriden by a component's model interface and should set things like the parameter path and other member-specific settings.
ClimaCalibrate.run_forward_model
— Functionrun_forward_model(model_config)
Execute the forward model simulation with the given configuration.
This function should be overridden with model-specific implementation details. config
should be obtained from set_up_forward_model
: run_forward_model(set_up_forward_model(member, iter, experiment_dir))
ClimaCalibrate.observation_map
— Functionobservation_map(iteration)
Runs the observation map for the specified iteration. This function must be implemented for each calibration experiment.
Backend Interface
ClimaCalibrate.get_backend
— Functionget_backend()
Get ideal backend for deploying forward model runs. Each backend is found via gethostname()
. Defaults to JuliaBackend if none is found.
ClimaCalibrate.calibrate
— Functioncalibrate(::Type{AbstractBackend}, config::ExperimentConfig; kwargs...)
+calibrate(::Type{AbstractBackend}, experiment_dir; kwargs...)
Run a full calibration, scheduling the forward model runs on Caltech's HPC cluster.
Takes either an ExperimentConfig or an experiment folder.
Available Backends: CaltechHPCBackend, ClimaGPUBackend, DerechoBackend, JuliaBackend
Keyword Arguments
- `experiment_dir: Directory containing experiment configurations.
- `model_interface: Path to the model interface file.
hpc_kwargs
: Dictionary of resource arguments, passed to the job scheduler.reruns
: Number of times to retry a failed ensemble member.verbose::Bool
: Enable verbose logging.- Any keyword arguments for the EnsembleKalmanProcess constructor, such as
scheduler
Usage
Open julia: julia --project=experiments/surface_fluxes_perfect_model
using ClimaCalibrate
+
+experiment_dir = joinpath(pkgdir(ClimaCalibrate), "experiments", "surface_fluxes_perfect_model")
+model_interface = joinpath(experiment_dir, "model_interface.jl")
+
+# Generate observational data and load interface
+include(joinpath(experiment_dir, "generate_data.jl"))
+include(joinpath(experiment_dir, "observation_map.jl"))
+include(model_interface)
+
+hpc_kwargs = kwargs(time = 3)
+backend = get_backend()
+eki = calibrate(backend, experiment_dir; model_interface, hpc_kwargs);
ClimaCalibrate.model_run
— Functionmodel_run(backend, iter, member, output_dir, experiment_dir; model_interface, verbose, hpc_kwargs)
Construct and execute a command to run a single forward model on a given job scheduler.
Dispatches on backend
to run slurm_model_run
or pbs_model_run
.
Arguments:
- iter: Iteration number
- member: Member number
- output_dir: Calibration experiment output directory
- experiment_dir: Directory containing the experiment's Project.toml
- model_interface: File containing the model interface
- moduleloadstr: Commands which load the necessary modules
- hpc_kwargs: Dictionary containing the resources for the job. Easily generated using
kwargs
.
ClimaCalibrate.module_load_string
— Functionmodule_load_string(backend)
Return a string that loads the correct modules for a given backend when executed via bash.
Job Scheduler
ClimaCalibrate.wait_for_jobs
— Functionwait_for_jobs(jobids, output_dir, iter, experiment_dir, model_interface, module_load_str, model_run_func; verbose, hpc_kwargs, reruns=1)
Wait for a set of jobs to complete. If a job fails, it will be rerun up to reruns
times.
This function monitors the status of multiple jobs and handles failures by rerunning the failed jobs up to the specified number of reruns
. It logs errors and job completion status, ensuring all jobs are completed before proceeding.
Arguments:
jobids
: Vector of job IDs.output_dir
: Directory for output.iter
: Iteration number.experiment_dir
: Directory for the experiment.model_interface
: Interface to the model.module_load_str
: Commands to load necessary modules.model_run_func
: Function to run the model.verbose
: Print detailed logs if true.hpc_kwargs
: HPC job parameters.reruns
: Number of times to rerun failed jobs.
ClimaCalibrate.log_member_error
— Functionlog_member_error(output_dir, iteration, member, verbose=false)
Log a warning message when an error occurs. If verbose, includes the ensemble member's output.
ClimaCalibrate.kill_job
— Functionkill_job(jobid::SlurmJobID)
+kill_job(jobid::PBSJobID)
End a running job, catching errors in case the job can not be ended.
ClimaCalibrate.job_status
— Functionjob_status(jobid)
Parse the slurm jobid's state and return one of three status symbols: :COMPLETED, :FAILED, or :RUNNING.
ClimaCalibrate.kwargs
— Functionkwargs(; kwargs...)
Create a dictionary from keyword arguments.
ClimaCalibrate.slurm_model_run
— Functionslurm_model_run(iter, member, output_dir, experiment_dir, model_interface, module_load_str; hpc_kwargs)
Construct and execute a command to run a single forward model on Slurm. Helper function for model_run
.
ClimaCalibrate.generate_sbatch_script
— Functiongenerate_sbatch_script(iter, member, output_dir, experiment_dir, model_interface; module_load_str, hpc_kwargs)
Generate a string containing an sbatch script to run the forward model.
hpc_kwargs
is turned into a series of sbatch directives using generate_sbatch_directives
. module_load_str
is used to load the necessary modules and can be obtained via module_load_string
.
ClimaCalibrate.generate_sbatch_directives
— Functiongenerate_sbatch_directives(hpc_kwargs)
Generate Slurm sbatch directives from HPC kwargs.
ClimaCalibrate.submit_slurm_job
— Functionsubmit_slurm_job(sbatch_filepath; env=deepcopy(ENV))
Submit a job to the Slurm scheduler using sbatch, removing unwanted environment variables.
Unset variables: "SLURMMEMPERCPU", "SLURMMEMPERGPU", "SLURMMEMPER_NODE"
ClimaCalibrate.pbs_model_run
— Functionpbs_model_run(iter, member, output_dir, experiment_dir, model_interface, module_load_str; hpc_kwargs)
Construct and execute a command to run a single forward model on PBS Pro. Helper function for model_run
.
ClimaCalibrate.generate_pbs_script
— Functiongeneratepbsscript( iter, member, outputdir, experimentdir, modelinterface; moduleloadstr, hpckwargs, )
Generate a string containing a PBS script to run the forward model.
Returns:
qsub_contents::Function
: A function generating the content of the PBS script based on the provided arguments. This will run the contents of thejulia_script
, which have to be run from a file due to Derecho'sset_gpu_rank
.julia_script::String
: The Julia script string to be executed by the PBS job.
Helper function for pbs_model_run
.
ClimaCalibrate.submit_pbs_job
— Functionsubmit_pbs_job(sbatch_filepath; env=deepcopy(ENV))
Submit a job to the PBS Pro scheduler using qsub, removing unwanted environment variables.
Unset variables: "PBSMEMPERCPU", "PBSMEMPERGPU", "PBSMEMPER_NODE"
EnsembleKalmanProcesses Interface
ClimaCalibrate.initialize
— Functioninitialize(ensemble_size, observations, noise, prior, output_dir; kwargs...)
+initialize(ensemble_size, observations, prior, output_dir; kwargs...)
+initialize(eki::EnsembleKalmanProcess, prior, output_dir)
+initialize(config::ExperimentConfig; kwargs...)
+initialize(filepath::AbstractString; kwargs...)
Initialize the EnsembleKalmanProcess object and parameter files.
Can take in an existing EnsembleKalmanProcess which will be used to generate the initial parameter ensemble.
Noise is optional when the observation is an EKP.ObservationSeries.
Additional kwargs may be passed through to the EnsembleKalmanProcess constructor.
ClimaCalibrate.save_G_ensemble
— Functionsave_G_ensemble(config::ExperimentConfig, iteration, G_ensemble)
+save_G_ensemble(output_dir::AbstractString, iteration, G_ensemble)
Saves the ensemble's observation map output to the correct directory based on the provided configuration. Takes an output directory, either extracted from an ExperimentConfig or passed directly.
ClimaCalibrate.update_ensemble
— Functionupdate_ensemble(output_dir::AbstractString, iteration, prior)
+update_ensemble(config::ExperimentConfig, iteration)
+update_ensemble(config_file::AbstractString, iteration)
Updates the EnsembleKalmanProcess object and saves the parameters for the next iteration.
ClimaCalibrate.ExperimentConfig
— TypeExperimentConfig(
+ n_iterations::Integer,
+ ensemble_size::Integer,
+ observations,
+ noise,
+ prior::ParameterDistribution,
+ output_dir,
+)
+ExperimentConfig(filepath::AbstractString; kwargs...)
Construct an ExperimentConfig from a given YAML file or directory containing 'experiment_config.yml'.
ExperimentConfig holds the configuration for a calibration experiment. This can be constructed from a YAML configuration file or directly using individual parameters.
ClimaCalibrate.get_prior
— Functionget_prior(param_dict::AbstractDict; names = nothing)
+get_prior(prior_path::AbstractString; names = nothing)
Constructs the combined prior distribution from a param_dict
or a TOML configuration file specified by prior_path
. If names
is provided, only those parameters are used.
ClimaCalibrate.get_param_dict
— Functionget_param_dict(distribution; names)
Generates a dictionary for parameters based on the specified distribution, assumed to be of floating-point type. If names
is not provided, the distribution's names will be used.
ClimaCalibrate.path_to_iteration
— Functionpath_to_iteration(output_dir, iteration)
Creates the path to the directory for a specific iteration within the specified output directory.
ClimaCalibrate.path_to_ensemble_member
— Functionpath_to_ensemble_member(output_dir, iteration, member)
Constructs the path to an ensemble member's directory for a given iteration and member number.
ClimaCalibrate.path_to_model_log
— Functionpath_to_model_log(output_dir, iteration, member)
Constructs the path to an ensemble member's forward model log for a given iteration and member number.