Boiling Water Reactor Dynamics

Replicated from “Hybrid Simulation of Boiling Water Reactor Dynamics Using A University Research Reactor” by James A. Turso, Robert M. Edwards, Jose March-Leuba, Nuclear Technology vol. 110, Apr. 1995.

This is a simple 5th-order representation of the operation of a boiling water reactor around its normal operating point, which is subject to interesting limit cycle dynamics.

The original article documents the model well, with the exception of the bifurcation parameter K and a nonlinear term, for which I’ve identified plausible values by experiment.


Oscillation from a purely positive loop

Replicated by Mohammad Mojtahedzadeh from Alan Graham’s thesis, or created anew with the same inspiration. He created these models in the course of his thesis work on structural analysis through pathway participation matrices.

Alan Graham, 1977. Principles on the Relationship Between Structure and Behavior of Dynamic Systems. MIT Thesis. Page 76+

These models are pure positive feedback loops that don’t exhibit exponential growth (under the right initial conditions). See my blog post for a discussion of the details.

These are generic models, and therefore don’t have units. All should run with Vensim PLE, except the generic gain matrix version which uses arrays and therefore requires an advanced version or the Model Reader.

The original 4th order model, replicated from Alan’s thesis: PurePosOscill4.vpm – note that this includes a .cin file with an alternate stable initialization.

My slightly modified version, permitting initialization with different gains at each level: PurePosOscill4alt.vpm

Loops of different orders: 3.vpm 6.vpm 8.vpm 12.vpm (I haven’t spent much time with these. It appears that the high-order versions transition to growth rather quickly – my guess is that this is an artifact of numerical precision, i.e. any tiny imprecision in the initialization introduces a bit of the growth eigenvector, which quickly swamps the oscillatory signal. It would be interesting to try these in double precision Vensim to see if I’m right.)

Stable initializations: 2stab.vpm 12stab.vpm

A generic version, representing a system as a generic gain matrix, so you can use it to explore any linear unforced variant: Generic.vpm

Bifurcating Salmon

A nifty paper on nonlinear dynamics of salmon populations caught my eye on today. The math is straightforward and elegant, so I replicated the model in Vensim.

A three-species model explaining cyclic dominance of pacific salmon

Authors: Christian Guill, Barbara Drossel, Wolfram Just, Eddy Carmack

Abstract: The four-year oscillations of the number of spawning sockeye salmon (Oncorhynchus nerka) that return to their native stream within the Fraser River basin in Canada are a striking example of population oscillations. The period of the oscillation corresponds to the dominant generation time of these fish. Various – not fully convincing – explanations for these oscillations have been proposed, including stochastic influences, depensatory fishing, or genetic effects. Here, we show that the oscillations can be explained as a stable dynamical attractor of the population dynamics, resulting from a strong resonance near a Neimark Sacker bifurcation. This explains not only the long-term persistence of these oscillations, but also reproduces correctly the empirical sequence of salmon abundance within one period of the oscillations. Furthermore, it explains the observation that these oscillations occur only in sockeye stocks originating from large oligotrophic lakes, and that they are usually not observed in salmon species that have a longer generation time.

The paper does a nice job of connecting behavior to structure, and of relating the emergence of oscillations to eigenvalues in the linearized system.

Units balance, though I had to add a couple implicit scale factors to do so.

The general results are qualitatitively replicable. I haven’t tried to precisely reproduce the authors’ bifurcation diagram and other experiments, in part because I couldn’t find a precise specification of numerical methods used (time step, integration method), so I wouldn’t expect to succeed.

Unlike most SD models, this is a hybrid discrete-continuous system. Salmon, predator and zooplankton populations evolve continuously during a growing season, but with discrete transitions between seasons.

The model uses SAMPLE IF TRUE, so you need an advanced version of Vensim to run it, or the free Model Reader. (It should be possible to replace the SAMPLE IF TRUE if an enterprising person wanted a PLE version). It would also be a good candidate for an application of SHIFT IF TRUE if someone wanted to experiment with the cohort age structure.


For a more policy-oriented take on salmon, check out Andy Ford’s work on smolt migration.

Lorenz Attractor

This is an implementation of Lorenz’ groundbreaking model that exhibits continuous-time chaos.

A google search turns up lots of good information on this model. For more advanced material, try google scholar.

I didn’t replicate this from Lorenz’ original 1963 article, Deterministic Nonperiodic Flow, but you can find a copy here.




Logistic Chaos

This is an implementation of the logistic model – a very simple example of discrete time chaotic behavior. It’s sometimes used to illustrate chaotic dynamics of insect populations.

There’s a nice description here, and the other top links on google tend to be good.

Note that this version corrects an equation error in previous versions.

Logistic (Vensim .vpm)

Logistic (Vensim .vmf)

Ultradian Oscillations of Insulin and Glucose

Citation: Jeppe Sturis, Kenneth S. Polonsky, Erik Mokilde, and Eve van Cauter. Computer Model for Mechanisms Underlying Ultradian Oscillations of Insulin and Glucose. Am. J. Physiol. 260 (Endocrinol. Metab. 23): E801-E809, 1991.

Source: Replicated by Hank Taylor

Units: No Yes!

Format: Vensim

Ultradian Oscillations of Insulin and Glucose (Vensim .vpm)

Update, 10/2017:

Refreshed, with units defined (mathematically the same as before): ultradia2.vpm ultradia2.mdl

Further refined, for initialization in equilibrium (insulin by analytic expression; glucose by parameter). Glucose infusion turned on by default. Graphs added.

ultradia-enhanced-3.mdl ultradia-enhanced-3.vpm