A Python framework for modelling and simulating dynamical systems. A dynamical system is a system whose state evolves over time according to a fixed set of rules, often expressed as mathematical equations, that depend on its current state and, optionally, external inputs.
Models are structured using state transitions for encoding differential equations, or any other logic, as an example. Simulations are configured using methods such as parameter sweeps, Monte Carlo runs, and A/B testing. See cadCAD.education for the most comprehensive cadCAD beginner course.
Goals:
- simple API for ease of use
- performance driven (more speed = more experiments, larger parameter sweeps, in less time)
- cadCAD compatible (standard functions, data structures, and simulation results)
- maintainable, testable codebase
Have a question not answered in this README or the cadCAD Edu courses? Post a comment in the issues or check out the cadCAD Discord community.
- Open-source Models Using radCAD
- Example Models
- Features
- Installation
- Documentation
- Development
- Testing
- Jupyter Notebooks
- Benchmarking
- Acknowledgements
- Ethereum Economic Model by CADLabs: A modular dynamical-systems model of Ethereum's validator economics
- Beacon Runner by Ethereum RIG: An agent-based model of Ethereum's Proof-of-Stake consensus layer
- GEB Controller Simulations by Reflexer Protocol: A Proportional-Integral-Derivative (PID) controller based upon a reference document approach for the Maker DAI market that was never implemented
- Fei Protocol Model by CADLabs: A modular dynamical-systems model of Fei Protocol
Using Models as live in-the-loop digital twins, creating your own model pipelines, and streaming simulation results to update a visualization. That's what an iterable Model class enables.
Live radCAD demo model on Streamlit
A simple game where at each timestep, the following transitions occur:
- Any live cell with fewer than two live neighbours dies, as if by underpopulation.
- Any live cell with two or three live neighbours lives on to the next generation.
- Any live cell with more than three live neighbours dies, as if by overpopulation.
- Any dead cell with exactly three live neighbours becomes a live cell, as if by reproduction.
See examples/game_of_life/game-of-life.ipynb
A simple model that applies the two Lotka-Volterra differential equations, frequently used to describe the dynamics of biological systems in which two species interact:
Original models thanks to Danilo @danlessa!
- System dynamics model: examples/predator_prey_sd/predator-prey-sd.ipynb
- Agent based model: examples/predator_prey_abm/predator-prey-abm.ipynb
- Parameter sweeps
params = {
'a': [1, 2, 3],
'b': [1, 2],
'c': [1]
}
# Creates a parameter sweep of:
# [{'a': 1, 'b': 1, 'c': 1}, {'a': 2, 'b': 2, 'c': 1}, {'a': 3, 'b': 2, 'c': 1}]- Monte Carlo runs
RUNS = 100 # Set to the number of Monte Carlo Runs
Simulation(model=model, timesteps=TIMESTEPS, runs=RUNS)- A/B tests
model_a = Model(initial_state=states_a, state_update_blocks=state_update_blocks_a, params=params_a)
model_b = Model(initial_state=states_b, state_update_blocks=state_update_blocks_b, params=params_b)
simulation_1 = Simulation(model=model_a, timesteps=TIMESTEPS, runs=RUNS)
simulation_2 = Simulation(model=model_b, timesteps=TIMESTEPS, runs=RUNS)
# Simulate any number of models in parallel
experiment = Experiment([simulation_1, simulation_2])
result = experiment.run()- cadCAD compatibility and familiar data structure
a b simulation subset run substep timestep
0 1.000000 2.0 0 0 1 0 0
1 0.540302 2.0 0 0 1 1 1
2 0.540302 7.0 0 0 1 2 1
3 0.463338 7.0 0 0 1 1 2
4 0.463338 12.0 0 0 1 2 2
... ... ... ... ... ... ... ...
799999 0.003162 999982.0 1 1 1 2 99998
800000 0.003162 999982.0 1 1 1 1 99999
800001 0.003162 999992.0 1 1 1 2 99999
800002 0.003162 999992.0 1 1 1 1 100000
800003 0.003162 1000002.0 1 1 1 2 100000- Tested against Python 3.8 to 3.12 using Ubuntu, MacOS, and Windows
- Disable
deepcopyoption for improved performance (at cost of mutability) - Robust exception handling with partial results, and tracebacks
- Parallel processing with multiple backend options:
multiprocessing,pathos,ray - Distributed computing and remote execution in a cluster (AWS, GCP, Kubernetes, ...) using Ray - Fast and Simple Distributed Computing
- Hooks to easily extend the functionality - e.g. save results to HDF5 file format after completion
- Model classes are iterable, so you can iterate over them step-by-step from one state to the next (useful for gradient descent, live digital twins)
- Parameters can be configured using nested dataclasses! This enables typing and dot notation for accessing parameters, and the creation of parameter namespaces.
pip install radcadradCAD provides the following classes:
- A system is represented in some form as a
Model - A
Modelcan be simulated using aSimulation - An
Experimentconsists of one or moreSimulations - An
Experimentor aSimulationis run by theEngine
So, the hierarchy is as follows Model > Simulation > Experiment > Engine.
from radcad import Model, Simulation, Experiment
model = Model(initial_state=initial_state, state_update_blocks=state_update_blocks, params=params)
simulation = Simulation(model=model, timesteps=100_000, runs=1)
result = simulation.run()
# Or, multiple simulations:
# experiment = Experiment([simulation_1, simulation_2, simulation_3])
# result = experiment.run()
df = pd.DataFrame(result)A unique feature of radCAD is being able to use nested dataclasses to configure the model's parameters. This enables typing and dot notation for accessing parameters, and the creation of parameter namespaces, demonstrated below.
from dataclasses import dataclass
from radcad.utils import default
from radcad.types import StateVariables, PolicySignal
...
@dataclass
class LiquidityPoolParameters:
initial_liquidity: int = 100_000
@dataclass
class ProtocolParameters:
liquidity_pool: LiquidityPoolParameters = default(LiquidityPoolParameters())
...
@dataclass
class Parameters:
protocol: ProtocolParameters = default(ProtocolParameters())
agents: AgentParameters = default(AgentParameters())
...
models.params = Parameters()
...
def update_liquidity(
params: Parameters,
substep: int,
state_history: List[List[StateVariables]],
previous_state: StateVariables,
policy_input: PolicySignal
) -> Tuple[str, int]:
if not previous_state["liquidity"]:
updated_liquidity = params.protocol.liquidity_pool.initial_liquidity
else:
updated_liquidity = ...
return "liquidity", updated_liquidity# cadCAD configuration modules
from cadCAD.configuration.utils import config_sim
from cadCAD.configuration import Experiment
# cadCAD simulation engine modules
from cadCAD.engine import ExecutionMode, ExecutionContext
from cadCAD.engine import Executor
# cadCAD global simulation configuration list
from cadCAD import configs
# Clears any prior configs
del configs[:]
sim_config = config_sim({
'N': 1, # Number of Monte Carlo Runs
'T': range(100), # Number of timesteps
'M': system_params # System Parameters
})
experiment.append_configs(
# Model initial state
initial_state=initial_state,
# Model Partial State Update Blocks
partial_state_update_blocks=partial_state_update_blocks,
# Simulation configuration
sim_configs=sim_config
)
# ExecutionContext instance (used for more advanced cadCAD config)
exec_context = ExecutionContext()
# Creates a simulation Executor instance
simulation = Executor(
exec_context=exec_context,
# cadCAD configuration list
configs=configs
)
# Executes the simulation, and returns the raw results
result, _tensor_field, _sessions = simulation.execute()
df = pd.DataFrame(result)from radcad import Model, Simulation, Experiment
model = Model(
# Model initial state
initial_state=initial_state,
# Model Partial State Update Blocks
state_update_blocks=state_update_blocks,
# System Parameters
params=params
)
simulation = Simulation(
model=model,
timesteps=100_000, # Number of timesteps
runs=1 # Number of Monte Carlo Runs
)
# Executes the simulation, and returns the raw results
result = simulation.run()
df = pd.DataFrame(result)radCAD is already compatible with the cadCAD generalized dynamical systems model structure; existing state update blocks, policies, and state update functions should work as is. But to more easily refactor existing cadCAD models to use radCAD without changing the cadCAD API and configuration process, there is a compatibility mode. The compatibility mode doesn't guarrantee to handle all cadCAD options, but should work for most cadCAD models by translating the configuration and execution processes into radCAD behind the scenes.
To use the compatibility mode, install radCAD with the compat dependencies:
pip install -e .[compat]
# Or
poetry install -E compat
# Or
pdm install -G compatThen, update the cadCAD imports from cadCAD._ to radcad.compat.cadCAD._
from radcad.compat.cadCAD.configuration import Experiment
from radcad.compat.cadCAD.engine import Executor, ExecutionMode, ExecutionContext
from radcad.compat.cadCAD.configuration.utils import config_sim
from radcad.compat.cadCAD import configsNow run your existing cadCAD model using radCAD!
Model classes are iterable, so you can iterate over them step-by-step from one state to the next.
This is useful for gradient descent, live digital twins, composing one model within another within a Policy Function...
Here is an example of using a Model to update a Plotly figure live:
from radcad import Model
import time
import plotly.graph_objects as go
# Live update of figure using Model as a generator
fig = go.FigureWidget()
fig.add_scatter()
fig.show()
# Create a generator from the Model iterator
model_generator = iter(Model(initial_state=initial_state, state_update_blocks=state_update_blocks, params=params))
timesteps = 100
results = []
for t in range(timesteps):
# Step to next state
model = next(model_generator)
# Get state and update figure
state = model.state
a = state['a']
results.append(a)
fig.data[0].y = results[:t]You have access to the more advanced engine options too, using the __call__() method:
model(raise_exceptions=False, deepcopy=True, drop_substeps=False)
_model = next(model)Current limitations:
- Only works for single subsets (no parameter sweeps)
By default radCAD uses the multiprocessing library and sets the number of parallel processes used by the Engine to the number of system CPUs less one, but this can be customized as follows:
from radcad import Engine
...
experiment.engine = Engine(processes=1)
result = experiment.run()Alternatively, select the SINGLE_PROCESS Backend option which doesn't use any parallelisation library:
from radcad import Engine, Backend
...
experiment.engine = Engine(backend=Backend.SINGLE_PROCESS)
result = experiment.run()To improve performance, at the cost of mutability, the Engine module has the deepcopy option which is True by default:
experiment.engine = Engine(deepcopy=False)The deepcopy method used for creating deepcopies of state during the simulation can be customised by creating a custom SimulationExecution class. This is useful when trying to optimise simulation performance, where certain types of state are better suited to specific deepcopy methods.
class CustomSimulationExecution(SimulationExecution):
def __init_subclass__(cls) -> None:
return super().__init_subclass__()
@staticmethod
def deepcopy_method(obj) -> Any:
"""
Use copy.deepcopy instead of default Pickle dumps/loads for deepcopy operations.
"""
return copy.deepcopy(obj)
experiment.engine.simulation_execution = CustomSimulationExecutionThis same technique can be used for overriding any of the other simulation methods.
If you don't need the substeps in post-processing, you can both improve simulation performance and save post-processing time and dataset size by dropping the substeps:
experiment.engine = Engine(drop_substeps=True)radCAD allows you to choose whether to raise exceptions, ending the simulation, or to continue with the remaining runs and return the results along with the exceptions. Failed runs are returned as partial results - the part of the simulation result up until the timestep where the simulation failed.
...
experiment.engine = Engine(raise_exceptions=False)
experiment.run()
results = experiment.results # e.g. [[{...}, {...}], ..., [{...}, {...}]]
exceptions = experiment.exceptions # A dataframe of exceptions, tracebacks, and simulations metadataThis also means you can run a specific simulation directly, and access the results later:
predator_prey_simulation.run()
...
results = predator_prey_simulation.resultsTo use the Ray backend, install radCAD with the extension-backend-ray dependencies:
pip install -e ".[extension-backend-ray]"
# Or
poetry install -E extension-backend-ray
# Or
pdm install -G extension-backend-rayExport the following AWS credentials (or see Ray documentation for alternative providers):
BACKEND=AWS
AWS_ACCESS_KEY_ID=***
AWS_SECRET_ACCESS_KEY=***Start a new cluster (or use existing):
# Cluster config: single m5.large EC2 instance, us-west-2 region
ray up cluster/aws/minimal.yaml
# Test the connection
ray exec cluster/aws/minimal.yaml 'echo "hello world"'Change the execution backend to RAY_REMOTE:
from radcad.engine import Engine, Backend
import ray
# Connect to cluster head
ray.init(address='***:6379', _redis_password='***')
...
experiment.engine = Engine(backend=Backend.RAY_REMOTE)
result = experiment.run()Finally, spin down the cluster:
ray down cluster/ray-aws.yamlHooks allow you to easily extend the functionality of radCAD with a stable API, and without having to manipulate the robust core.
experiment.before_experiment = lambda experiment: print(f"Before experiment with {len(experiment.simulations)} simulations")
experiment.after_experiment = lambda experiment: print(f"After experiment with {len(experiment.simulations)} simulations")
experiment.before_simulation = lambda simulation: print(f"Before simulation {simulation.index} with params {simulation.model.params}")
experiment.after_simulation = lambda simulation: print(f"After simulation {simulation.index} with params {simulation.model.params}")
experiment.before_run = lambda context: print(f"Before run {context}")
experiment.after_run = lambda context: print(f"After run {context}")
experiment.before_subset = lambda context: print(f"Before subset {context}")
experiment.after_subset = lambda context: print(f"After subset {context}")See tests/test_hooks.py for expected functionality.
import pandas as pd
import datetime
def save_to_HDF5(experiment, store_file_name, store_key):
now = datetime.datetime.now()
store = pd.HDFStore(store_file_name)
store.put(store_key, pd.DataFrame(experiment.results))
store.get_storer(store_key).attrs.metadata = {
'date': now.isoformat()
}
store.close()
print(f"Saved experiment results to HDF5 store file {store_file_name} with key {store_key}")
experiment.after_experiment = lambda experiment: save_to_HDF5(experiment, 'experiment_results.hdf5', 'experiment_0')The biggest performance bottleneck with radCAD, and cadCAD for that matter, is avoiding mutation of state variables by creating a deep copy of the state passed to the state update function. This avoids the state update function mutating state variables outside of the framework by creating a copy of it first - a deep copy creates a copy of the object itself, and the key value pairs, which gets expensive.
To avoid the additional overhead, mutation of state history is allowed, and left up to the developer to avoid using standard Python best practises, but mutation of the current state is disabled.
See https://stackoverflow.com/questions/24756712/deepcopy-is-extremely-slow for some performance benchmarks of different methods. radCAD uses cPickle, which is faster than using deepcopy, but less flexible about what types it can handle (Pickle depends on serialization) - these could be interchanged in future.
Set up and enter the Python radCAD development environment with the PDM package manager:
pdm --help
pdm use "python3.10"
# NB: Use pdm-py38.lock for Python3.8 and pdm.lock for Python3.9+
pdm install --lockfile pdm.lockERROR:: Could not find a local HDF5 installation.
You may need to explicitly state where your local HDF5 headers and
library can be found by setting the ``HDF5_DIR`` environment
variable or by using the ``--hdf5`` command-line option.
# 1. Update `pyproject.toml` package version using semantic versioning
# 2. Update CHANGELOG.md
# 3. Submit PR and run tests
# 4. Merge into master on success
# 5. Build and publish package
pdm publish
# Enter in PyPI package repository credentials
# 6. Tag master commit with version e.g. `v0.5.0` and pushExport requirements.txt using PDM:
pdm export -o requirements.txt --without-hashesNote: the root requirements.txt is used for Streamlit app in examples, and is not for development.
This project uses Nox for executing tests and benchmarks across different Python versions. Nox expects that each of the Python versions defined in the tests is available.
nox is a command-line tool that automates testing in multiple Python environments, similar to tox. Unlike tox, Nox uses a standard Python file for configuration.
To install Nox with pipx:
pipx install noxYou can also use pip in your global python:
python3 -m pip install noxTo run all the tests:
nox --session testsTo run the default radCAD benchmark in benchmarks/ across all Python versions:
nox --session benchmarksSee noxfile.py for other session types.
# Install kernel
pdm run ipykernel install --user --name python3-radcad
# Start Jupyter
pdm run jupyter labpdm run pytest benchmarks/benchmark_radcad.py
pdm run pytest benchmarks/benchmark_single_process.pypdm run mprof run benchmarks/benchmark_memory_radcad.py
pdm run mprof plot- @danlessa: Thank you for contribution of Predator-Prey example
- @rogervs: Thank you for contribution of Harmonic Oscillator example
- @abzaremba: Thank you for contribution to documentation
- @smngvlkz: Thank you for contribution to cadCAD compatibility mode and CI pipeline