Creating User-defined Functions in SD Models
system dynamics, systemdynamics, sd dsl, bptk, bptk-py, python, business simulation
Creating User-defined Functions in SD Models
One of the benefits of creating System Dynamics models in Python is that we can use the full power of Python to create our own functions, which we can then use in our models.
This how to illustrates how to do this.
First of all, lets set up our model:
Now let’s define a function we would like to use in our model. A user defined function can have as many arguments as you like, but it must accept at least a model and time parameter (you don’t need to use the parameters if you don’t want to).
How you define your function is up to you - you can use any of the methods available in Python, such as class methods, using def, or lambda functions.
The example below uses a lambda function which simply multiplies the current time t with 5.
As you can see, much like with stocks and converters, we associate our function with the model by calling the model’s function method.
Next we set up a converter whose equation calls the function, and test it at t = 5:
25.0
Let’s plot the function over time:
We can also create a stock that has the converter as an inflow:
You can do all the usual arithmetic in an equation. Here the same converter is divided by time before it flows in - a second stock rather than a second equation on the first one, because an element whose equation is set in two cells has no defined value outside the cell that just set it:
The function we created above was just dependent on time and not on other model variables. Let’s create a function that takes more arguments, e.g. one that multiplies a model variable with time.
You can add as many arguments as you like, but they must come after the model and t arguments.
The function needs two inputs, so we define a converter for each of them and then a third converter that applies another_model_function to both:
Of course functions defined in this way can also be used within scenarios. The quickest way to set up a scenario manager for a given model is register_model, which creates a scenario manager named after the model with the prefix “sm” - the name is normalized to start with a capital letter, so TestModel becomes smTestmodel - together with a “base” scenario that runs the model as-is.
Three more scenarios are added below, each with a different multiplier, and list_scenarios shows what the manager now holds:
*** smTestmodel ***
base
multiplier5
multiplier10
multiplier15
And plotting the four scenarios against each other shows what the multiplier does:
Functions and arrays
An argument can be an arrayed element, and a function says once - when it is registered - what that should mean.
By default a function is applied per index, which is the rule every operator in the DSL follows: the result is an array of the same shape, and a named array keeps its labels.
A function declared with elementwise=False receives the whole array instead and answers once. An unnamed array arrives as a list, a named one as a dict keyed by its labels, and a matrix arrives nested. The result is a single value, so anything built on top of it is scalar too.
250.0
Which of the two you want is a property of the function, not of the call, so it is stated where the function is registered. The default is the element-wise one, because that is what every other operator does with an array.
The aggregations the DSL already carries - arr_sum, arr_mean, arr_max and the rest - need no function at all, and they run on the Rust engine. Reach for elementwise=False when the reduction is one the DSL does not have.
Where a user-defined function may stand
A function computes the value of an element, so it goes where an equation goes: a converter, a flow, a biflow, and the equation of a stock. It can also be a stock’s initial value, since an initial value is an expression like any other:
starting_point = model.function("starting_point", lambda model, t: 42.0)
stock = model.stock("stock")
stock.initial_value = starting_point()An initial value is evaluated before the run rather than during it, so a function standing there is called once, with t at the start time.
Which functions run on the Rust engine
Since 3.2.0 a model with user-defined functions runs on the Rust engine. The engine evaluates the model and calls back into Python at those nodes, once per node and timestep, so the function you wrote is the function that runs.
The run says so, because it is not a pure engine run:
[WARN] This model calls 1 user-defined Python function(s) ('largest_region'). Those
nodes are evaluated in Python while the rest of the model runs on the Rust engine.Nothing is wrong there - the numbers are the same either way - but each call crosses between two runtimes, and what that costs is measured.
Two shapes stay on the Python engine, and asking for backend="rust" with them raises RustBackendError rather than quietly computing the answer elsewhere:
elementwise=False. Such a function is handed a whole array and answers once, where a callback answers with one number per node.- A hybrid model - agents and a System Dynamics side - whose function reads what the agents produced. The two halves advance one step at a time, and the engine computes every step of a run at once, so the function would read a step that has not happened yet. A model that merely has agents is fine.
On the server
/execute refuses a model that carries user-defined functions, with HTTP 400 naming them. The model arrived over HTTP, and running the Python that goes with it is a different decision from running the model. Serve such a model from a process that holds the functions - a bptk instance with the scenario manager registered - rather than by posting the model itself.