Class ODE
- All Implemented Interfaces:
Discrete
,CLOProvider
- Direct Known Subclasses:
PDE
,RungeKutta
,SDE
Note: Currently differential equation models are restricted to discrete strategy sets.
- Author:
- Christoph Hauert
-
Nested Class Summary
Nested ClassesModifier and TypeClassDescriptionstatic interface
Methods that everyModule
must implement, which advertises numerical solutions based on differential equations.static interface
Additional methods that must be implemented byModule
s that advertise numerical solutions based on ordinary differential equations.static enum
Types of initial configurations. -
Field Summary
FieldsModifier and TypeFieldDescription(package private) double
Desired accuracy to determine whether numerical integration has converged or a population has become monomorphic.final CLOption
Command line option to activate adjusted replicator dynamics (if possible) in ODE models.final CLOption
Command line option to set the desired accuracy of ODE calculations as criteria for convergence and whether population is monomorphic.final CLOption
Command line option to set the time incrementdt
in DE models.final CLOption
Command line option to set the initial configuration in ODE/SDE models.final CLOption
Command line option to set the number of generations between reports forEvoLudo.modelNext()
.(package private) int[]
Array containing the indices of the dependent trait in each species.(package private) double[]
Pointer to temporarily store a reference to the current stateyt
.(package private) double
Discretization of time increment for continuous time models.(package private) double
The size of step taken by integrator.(package private) double
The attempted size of next step to take.(package private) double[]
(package private) boolean
The flag istrue
if the numerical integration moves forward in time.(package private) double[]
The current fitness of all traits in the ODE.(package private) int[]
Array containing indices that delimit individual species in the ODE state variables.protected ODE.InitType[]
Type of initial configuration for each species.(package private) double[]
Array containing the inverse of the fitness range:1.0/(maxFitness - minFitness)
for each species.(package private) boolean
true
if the adjusted replicator dynamics should be applied.protected boolean
Set totrue
to requests that the numerical integration stops once the population reaches a monomorphic state.(package private) Mutation.Discrete[]
Array containing the mutation rates/probabilities for each species.(package private) String[]
The names of the traits (dynamical variables) in the ODE.(package private) int
Total number of traits (dynamical variables) in the ODE.(package private) double[]
Array containing the update rates for each species.(package private) double[]
ForModule
s with static fitness, e.g.(package private) double[]
The initial frequencies/densities of the ODE.(package private) double[]
The previous frequencies/densities of the ODE.(package private) double[]
The current frequencies/densities of the ODE.Fields inherited from class Model
cloSamples, cloTimeRelax, cloTimeStep, cloTimeStop, connect, converged, engine, fixData, isMultispecies, isRelaxing, logger, mode, nSamples, nSpecies, nStatisticsFailed, nStatisticsSamples, rng, species, statisticsSampleNew, time, timeRelax, timeStep, timeStop, type
-
Constructor Summary
Constructors -
Method Summary
Modifier and TypeMethodDescriptionboolean
check()
Check consistency of parameters and adjust if necessary (and possible).protected boolean
checkConvergence
(double dist2) Helper method to check whether the squared distancedist2
qualifies to signal convergence.void
collectCLO
(CLOParser parser) All providers of command line options must implement this method to collect their options.protected double
deStep
(double step) Attempts a numerical integration step of sizestep
.(package private) void
encodeFitness
(StringBuilder plist) Encodes the fitness of the model in the form of aplist
string.void
encodeState
(StringBuilder plist) Encode the state of the model in aplist
inspiredXML
string.(package private) void
encodeStrategies
(StringBuilder plist) Encodes state of the model in the form of aplist
string.double
Gets the numerical accuracy for determining convergence.protected void
getDerivatives
(double t, double[] state, double[] fitness, double[] change) Calculate the rates of change for all species instate
attime
given the fitnessfitness
and returned inchange
.double
getDt()
Gets the discretization of time increments in continuous time models.<T> void
getFitnessData
(int id, T[] colors, ColorMap.Gradient1D<T> colorMap) Translates fitness data into colors using ColorMapcolorMap
.void
getFitnessHistogramData
(int id, double[][] bins) Generates a histogram of fitness data and returns the result in the provided arraybins
.void
getInitialTraits
(double[] init) Collect and return initial trait values for all species.void
getInitialTraits
(int id, double[] init) Return initial trait values for species with IDid
.double
getMaxFitness
(int id) Calculates and returns the absolute fitness maximum.double
getMaxScore
(int id) Returns the maximum score that individuals of species with IDid
can achieve in this model.boolean
getMeanFitness
(double[] mean) Gets the mean fitness values for traits in all species.boolean
getMeanFitness
(int id, double[] mean) Gets the mean fitness values for species with IDid
.boolean
getMeanTraits
(double[] mean) Collect and return mean trait values for all species.boolean
getMeanTraits
(int id, double[] mean) Return mean trait values for species with IDid
.double
getMinFitness
(int id) Calculates and returns the absolute fitness minimum.double
getMinScore
(int id) Returns the minimum score that individuals of species with IDid
can achieve in this model.double
getMonoScore
(int id, int idx) Calculate and return the payoff/score of individuals in monomorphic populations with trait/strategyidx
but also deals with payoff accounting (averaged versus accumulated).int
getNMean()
Return the number of mean values for this model including all species (for traits or fitness).int
getNMean
(int id) Return the number of mean trait values for species with IDid
.Returns status message from model.<T> void
getTraitData
(int id, T[] colors, ColorMap<T> colorMap) Unused interface method.void
init()
Milestone: Initialize this modelprivate void
init
(boolean doRandom) Helper method to set the intital state of the model.boolean
Return whether this DE model tracks frequencies or densities.boolean
Check if population is monomorphic.private boolean
isMonomorphic
(int start, int stop, int dep, int vac) Helper method to check if the species with trait indicesstart
throughstop
is monomorphic.boolean
Checks if time is reversed.void
load()
Milestone: Load this model and allocate resources (if applicable).boolean
next()
Advance model by one step.protected void
normalizeState
(double[] state) Convenience method to normalize state in frequency based models.boolean
Parse initializer stringarg
.boolean
Checks if time reversal is permitted.void
reset()
Milestone: Reset this model(package private) boolean
restoreFitness
(Plist plist) Restores the fitness encoded inplist
.boolean
restoreState
(Plist plist) Restore the state encoded in theplist
inspiredmap
ofkey, value
-pairs.(package private) boolean
restoreStrategies
(Plist plist) Restores the state encoded inplist
.void
setAccuracy
(double acc) Sets the desired accuracy for determining convergence.void
setAdjustedDynamics
(boolean adjusted) Sets whether adjusted dynamics is requested.void
setDt
(double deltat) Sets the discretization of time increments in continuous time models.boolean
setInitialTraits
(double[] init) Set initial traits in one or all species.boolean
setInitialTraits
(int id, double[] init) Set initial trait values for species with IDid
.void
setTimeReversed
(boolean reversed) Request time reversal ifreversed==true
.void
unload()
Milestone: Unload this model and free resources (if applicable).void
update()
Update this model.protected double
updateBest
(Module mod, double[] state, double[] fitness, int nGroup, int index, double[] change) Implementation of player updatePlayerUpdate.Type.BEST
.protected double
updateBestResponse
(Module mod, double[] state, double[] fitness, int nGroup, int index, double[] change) Implementation of the player updatePlayerUpdate.Type.BEST_RESPONSE
.protected double
updateEcology
(Module mod, double[] state, double[] fitness, int nGroup, int index, double[] change) Implementation of the player updates for modules with populations of variable size (density based or with vacant 'space', i.e.protected double
updateImitate
(Module mod, double[] state, double[] fitness, int nGroup, int index, double[] change) Implementation of the player updatePlayerUpdate.Type.IMITATE
.protected double
updateImitateBetter
(Module mod, double[] state, double[] fitness, int nGroup, int index, double[] change) Implementation of the player updatePlayerUpdate.Type.IMITATE_BETTER
.protected double
updateProportional
(Module mod, double[] state, double[] fitness, int nGroup, int index, double[] change) Implementation of the player updatePlayerUpdate.Type.PROPORTIONAL
.private double
updateReplicate
(Module mod, double[] state, double[] fitness, int nGroup, int index, double[] change, double noise) Helper method to calculate the rate of change for the standard replicator dynamics with different amounts of noise arising from the microscopic update rule.protected double
updateThermal
(Module mod, double[] state, double[] fitness, int nGroup, int index, double[] change) Implementation of the player updatePlayerUpdate.Type.THERMAL
.Methods inherited from class Model
debugStep, getCounter, getFitnessNameAt, getFixationData, getLogger, getMeanFitnessAt, getMeanName, getMeanNames, getMeanTraitAt, getMode, getNextHalt, getNSamples, getNSpecies, getNStatisticsFailed, getNStatisticsSamples, getScoreNameAt, getSpecies, getTime, getTimeRelax, getTimeStep, getTimeStop, getTraitNameAt, getType, hasConverged, initStatisticsFailed, initStatisticsSample, isConnected, isContinuous, isDE, isIBS, isODE, isPDE, isRelaxing, isSDE, isType, permitsDebugStep, permitsMode, permitsSampleStatistics, permitsUpdateStatistics, readStatisticsSample, relax, requestMode, resetStatisticsSample, setMode, setNSamples, setTimeRelax, setTimeStep, setTimeStop, useScheduling
Methods inherited from class Object
clone, equals, finalize, getClass, hashCode, notify, notifyAll, toString, wait, wait, wait
Methods inherited from interface CLOProvider
adjustCLO
-
Field Details
-
dt
double dtDiscretization of time increment for continuous time models. This is the attempted size of next step to take.Important: always positive regardless of direction of integration.
- See Also:
-
dtTry
double dtTryThe attempted size of next step to take. InODE
this is always equal todt
but for more sophisticated integrators,RungeKutta
for example, the attempted step size may get adjusted.Important: always positive regardless of direction of integration.
- See Also:
-
dtTaken
double dtTakenThe size of step taken by integrator. This can be less than the step size attempted,dtTaken ≤ dtTry
.Important: always positive regardless of direction of integration.
- See Also:
-
forward
boolean forwardThe flag istrue
if the numerical integration moves forward in time.- See Also:
-
nDim
int nDimTotal number of traits (dynamical variables) in the ODE. -
y0
double[] y0The initial frequencies/densities of the ODE. In multi-species modules the initial states of each species are concatenated into this single array. The initial frequencies/densities of the first species are stored inyt[0]
throughyt[n1]
wheren1
denotes the number of traits in the first species. The initial frequencies/densities of the second species are stored inyt[n1+1]
throughyt[n1+n2]
wheren2
denotes the number of traits in the second species, etc. -
yt
double[] ytThe current frequencies/densities of the ODE. In multi-species modules the current states of each species are concatenated into this single array. The current frequencies/densities of the first species are stored inyt[0]
throughyt[n1]
wheren1
denotes the number of traits in the first species. The current frequencies/densities of the second species are stored inyt[n1+1]
throughyt[n1+n2]
wheren2
denotes the number of traits in the second species, etc. -
ft
double[] ftThe current fitness of all traits in the ODE. In multi-species modules the current fitness of each species are concatenated into this single array. The current fitness of the first species are stored inft[0]
throughft[n1]
wheren1
denotes the number of traits in the first species. The current fitness of the second species are stored inft[n1+1]
throughft[n1+n2]
wheren2
denotes the number of traits in the second species, etc. -
dyt
double[] dytThe frequency/density increments of the ODE between the current stateyt
and the previous stateyout
. In multi-species modules the increments of each species are concatenated into this single array. The increments in the first species are stored indyt[0]
throughdyt[n1]
wheren1
denotes the number of traits in the first species. The increments in the second species are stored indyt[n1+1]
throughdyt[n1+n2]
wheren2
denotes the number of traits in the second species, etc. -
yout
double[] youtThe previous frequencies/densities of the ODE. In multi-species modules the previous states of each species are concatenated into this single array. The previous frequencies/densities of the first species are stored inyout[0]
throughyout[n1]
wheren1
denotes the number of traits in the first species. The current frequencies/densities of the second species are stored inyout[n1+1]
throughyout[n1+n2]
wheren2
denotes the number of traits in the second species, etc. -
dstate
double[] dstatePointer to temporarily store a reference to the current stateyt
. This affects only multi-species modules because fitness calculations for the derivatives may need access to the state of the other populations. To allow thisgetMeanTraits(int, double[])
usesdstate
while the derivative calculations are in progress (dstate != null
). Otherwiseyt
is used. -
staticfit
double[] staticfitForModule
s with static fitness, e.g. the original Moran processMoran
, this array stores the fitness of all types for efficiency reasons. In multi-species modules the static fitness of each species are concatenated into this single array. The fitness of the first species are stored instaticfit[0]
throughstaticfit[n1]
wheren1
denotes the number of traits in the first species. The static fitness of the second species are stored instaticfit[n1+1]
throughstaticfit[n1+n2]
wheren2
denotes the number of traits in the second species, etc. If only some species have static fitness only those are stored here. All other entries in the array are unused but the size remains unchanged and is always the same as e.g.yt
ornull
if no species has static fitness. -
initType
Type of initial configuration for each species.- See Also:
-
idxSpecies
int[] idxSpeciesArray containing indices that delimit individual species in the ODE state variables. In a module withN
species the array hasN+1
entries. For example, for a single species module withn1
traitsidxSpecies = int[] { 0, n1 }
holds, whereas for a two species module withn1
andn2
traits, respectively,idxSpecies = int[] { 0, n1, n1+n2 }
holds, and for three speciesidxSpecies = int[] { 0, n1, n1+n2, n1+n2+n3 }
etc. In other words, the first entry is always zero and the last always contains the total number of dynamical variables in this ODE,idxSpecies[N]==nDim
. This means that, e.g. the frequencies/densities referring to speciesi
are stored inyt[idxSpecies[i]]
throughyt[idxSpecies[i+1]]
.- See Also:
-
rates
double[] ratesArray containing the update rates for each species. Determines the relative update rate between species, Only relevant in multi-species modules.- See Also:
-
mutation
Mutation.Discrete[] mutationArray containing the mutation rates/probabilities for each species.- See Also:
-
dependents
int[] dependentsArray containing the indices of the dependent trait in each species. Dependent traits make only sense in replicator type model with frequencies, i.e. normalized states. The dependent trait is the one that gets derived from the changes in all others in order to maintain normalization.- See Also:
-
invFitRange
double[] invFitRangeArray containing the inverse of the fitness range:1.0/(maxFitness - minFitness)
for each species. This is used to normalize imitation rules of players.- See Also:
-
names
String[] namesThe names of the traits (dynamical variables) in the ODE. In multi-species modules all names are concatenated into this single array. The names of the first species are stored innames[0]
throughnames[n1]
wheren1
denotes the number of traits in the first species. The names of the second species are stored innames[n1+1]
throughnames[n1+n2]
wheren2
denotes the number of traits in the second species, etc. -
monoStop
protected boolean monoStopSet totrue
to requests that the numerical integration stops once the population reaches a monomorphic state.- See Also:
-
isAdjustedDynamics
boolean isAdjustedDynamicstrue
if the adjusted replicator dynamics should be applied. Compared to the standard replicator equation nothing changes in terms of trajectories, location of fixed points or their stability. However, the speed of evolutionary changes depend inversely proportional on the overall fitness of the population(s) such that changes happen rapidly in a population with low fitness but slowly in populations with high fitness. The adjusted replicator dynamics requires that all payoffs are positive.- See Also:
-
accuracy
double accuracyDesired accuracy to determine whether numerical integration has converged or a population has become monomorphic.- See Also:
-
cloDEdt
Command line option to set the time incrementdt
in DE models. If the requesteddt
turns out to be too big, e.g. due to diffusion inPDE
, it is reduced appropriately. ForRungeKutta
this represents the initial step.- See Also:
-
cloInit
Command line option to set the initial configuration in ODE/SDE models.Note: PDE models use their own set of initial configurations.
- See Also:
-
cloAdjustedDynamics
Command line option to activate adjusted replicator dynamics (if possible) in ODE models.- See Also:
-
cloDEAccuracy
Command line option to set the desired accuracy of ODE calculations as criteria for convergence and whether population is monomorphic.- See Also:
-
cloTimeReversed
Command line option to set the number of generations between reports forEvoLudo.modelNext()
.
-
-
Constructor Details
-
ODE
Constructs a new model for the numerical integration of the system of ordinary differential equations representing the dynamics specified by theModule
module
using theEvoLudo
pacemakerengine
to control the numerical evaluations. The integrator implements Euler's method (fixed step size).- Parameters:
engine
- the pacemaker for running the model
-
-
Method Details
-
load
public void load()Description copied from class:Model
Milestone: Load this model and allocate resources (if applicable). -
unload
public void unload()Description copied from class:Model
Milestone: Unload this model and free resources (if applicable). -
check
public boolean check()Description copied from class:Model
Check consistency of parameters and adjust if necessary (and possible). All issues and modifications should be reported throughlogger
. Some parameters can be adjusted while the model remains active or even while running, whereas others require a reset. An example of the former category is in general simple adjustments of payoffs, while an example of the latter category is a change of the population structure. -
reset
public void reset()Description copied from class:Model
Milestone: Reset this model -
init
public void init()Description copied from class:Model
Milestone: Initialize this model -
setDt
public void setDt(double deltat) Sets the discretization of time increments in continuous time models.Note: Some models may need to adjust, i.e. reduce,
deltat
(seePDE.checkDt()
) or choose a variable step size in which casedeltat
is ignored (seeRungeKutta.getAutoDt()
).- Parameters:
deltat
- the time increments in continuous time models.
-
getDt
public double getDt()Gets the discretization of time increments in continuous time models.Note: This may be different from
dt
set throughsetDt(double)
if the model required adjustments.- Returns:
- the time increment in continuous time models.
-
isDensity
public boolean isDensity()Return whether this DE model tracks frequencies or densities. Returnsfalse
(i.e. frequency based model) by default.- Returns:
true
if state refers to densities.
-
setInitialTraits
public boolean setInitialTraits(double[] init) Description copied from class:Model
Set initial traits in one or all species.NOTE: this is a convenience method for multi-species modules to set inital states efficiently for interactions with GUI.
- Overrides:
setInitialTraits
in classModel
- Parameters:
init
- the array with the initial trait values- Returns:
true
if initial traits successfully set
-
getInitialTraits
public void getInitialTraits(double[] init) Description copied from class:Model
Collect and return initial trait values for all species.NOTE: this is a convenience method for multi-species modules to retrieve states efficiently for further processing or visualization.
- Specified by:
getInitialTraits
in classModel
- Parameters:
init
- the array for storing the initial trait values
-
setInitialTraits
public boolean setInitialTraits(int id, double[] init) Description copied from class:Model
Set initial trait values for species with IDid
.- Overrides:
setInitialTraits
in classModel
- Parameters:
id
- the species identifierinit
- the array with the initial trait values- Returns:
true
if initial traits successfully set
-
getInitialTraits
public void getInitialTraits(int id, double[] init) Description copied from class:Model
Return initial trait values for species with IDid
.- Specified by:
getInitialTraits
in classModel
- Parameters:
id
- the species identifierinit
- the array for storing the initial trait values
-
getMonoScore
public double getMonoScore(int id, int idx) Description copied from interface:Discrete
Calculate and return the payoff/score of individuals in monomorphic populations with trait/strategyidx
but also deals with payoff accounting (averaged versus accumulated).- Specified by:
getMonoScore
in interfaceDiscrete
- Parameters:
id
- the id of the population for multi-species modelsidx
- trait/strategy- Returns:
- payoff/score in monomorphic population with trait/strategy
idx
. ReturnsNaN
if scores ill defined - See Also:
-
getNMean
public int getNMean()Description copied from class:Model
Return the number of mean values for this model including all species (for traits or fitness). By default this returns the number of traits in the module. Models that report a different number of mean traits must override this method -
getNMean
public int getNMean(int id) Description copied from class:Model
Return the number of mean trait values for species with IDid
. -
getMeanTraits
public boolean getMeanTraits(int id, double[] mean) Description copied from class:Model
Return mean trait values for species with IDid
.- Specified by:
getMeanTraits
in classModel
- Parameters:
id
- the species identifiermean
- the array for storing the mean trait values- Returns:
true
if this and the previous data point should be connected, i.e. no reset had been requested in the mean time.
-
getMeanTraits
public boolean getMeanTraits(double[] mean) Description copied from class:Model
Collect and return mean trait values for all species.NOTE: this is a convenience method for multi-species modules to retrieve states efficiently for further processing or visualization.
- Specified by:
getMeanTraits
in classModel
- Parameters:
mean
- the array for storing the mean trait values- Returns:
true
if this and previous data point should be connected, i.e. no reset had been requested in the mean time.
-
getTraitData
Unused interface method.- Specified by:
getTraitData
in classModel
- Type Parameters:
T
- color data type.Color
forPopGraph2D
andPopGraph2D
as well asMeshLambertMaterial
forPopGraph3D
- Parameters:
id
- the species identifiercolors
- the array for storing the colors for individualscolorMap
- the map for translating individual traits into colors
-
getMinScore
public double getMinScore(int id) Description copied from class:Model
Returns the minimum score that individuals of species with IDid
can achieve in this model. Takes into account potential adjustments due to population structure and payoff accounting.- Specified by:
getMinScore
in classModel
- Parameters:
id
- the id of the population for multi-species models- Returns:
- the minimum score
- See Also:
-
getMaxScore
public double getMaxScore(int id) Description copied from class:Model
Returns the maximum score that individuals of species with IDid
can achieve in this model. Takes into account potential adjustments due to population structure and payoff accounting.- Specified by:
getMaxScore
in classModel
- Parameters:
id
- the id of the population for multi-species models- Returns:
- the maximum score
- See Also:
-
getMinFitness
public double getMinFitness(int id) Description copied from class:Model
Calculates and returns the absolute fitness minimum. This is important to- determine probabilities or rates for adopting the strategy of another player,
- optimize fitness based picking of individuals, and
- scaling graphical output.
- Specified by:
getMinFitness
in classModel
- Parameters:
id
- the id of the population for multi-species models- Returns:
- the minimum fitness
- See Also:
-
getMaxFitness
public double getMaxFitness(int id) Description copied from class:Model
Calculates and returns the absolute fitness maximum. This is important to- determine probabilities or rates for adopting the strategy of another player,
- optimize fitness based picking of individuals, and
- scaling graphical output.
- Specified by:
getMaxFitness
in classModel
- Parameters:
id
- the id of the population for multi-species models- Returns:
- the maximum fitness
- See Also:
-
getMeanFitness
public boolean getMeanFitness(int id, double[] mean) Description copied from class:Model
Gets the mean fitness values for species with IDid
. The result is stored and returned inmean
. Used by GUI to visualize local dynamics atidx
.- Specified by:
getMeanFitness
in classModel
- Parameters:
id
- the species identifiermean
- the array for storing the mean fitness values- Returns:
true
if this and the previous data point should be connected, i.e. no reset had been requested in the mean time.
-
getMeanFitness
public boolean getMeanFitness(double[] mean) Description copied from class:Model
Gets the mean fitness values for traits in all species. The result is stored and returned inmean
. Used by GUI to visualize the current the state of the model.Note: this is a convenience method for multi-species modules to retrieve states efficiently for further processing or visualization.
- Specified by:
getMeanFitness
in classModel
- Parameters:
mean
- the array for storing the mean fitness values- Returns:
true
if this and the previous data point should be connected, i.e. no reset had been requested in the mean time.
-
getFitnessData
Description copied from class:Model
Translates fitness data into colors using ColorMapcolorMap
. Used by GUI to visualize current state of model.- Specified by:
getFitnessData
in classModel
- Type Parameters:
T
- color data type.Color
forPopGraph2D
andPopGraph2D
as well asMeshLambertMaterial
forPopGraph3D
- Parameters:
id
- the species identifiercolors
- the array for storing color valuescolorMap
- the map to use for translating traits to colors
-
getFitnessHistogramData
public void getFitnessHistogramData(int id, double[][] bins) Description copied from class:Model
Generates a histogram of fitness data and returns the result in the provided arraybins
. For Discrete modules a fitness histogram is returned for each trait separately. For Continuous modules there is, in general, only a single fitness dimension.- Specified by:
getFitnessHistogramData
in classModel
- Parameters:
id
- the species identifierbins
- the array for storing histogram. For Discrete modules this is always one dimensional
-
next
public boolean next()Advance model by one step. The details of what happens during one step depends on the modelsType
as well as itsMode
.Implementation Notes:
Before:yt
current stateft
current fitnessdyt
current derivatives fitnessdtTry
size of time step attempting to make
yt
next stateft
next fitnessdyt
next derivativesyout
contains previous state (i.e.yt
when entering method)stepTaken
the time betweenyt
andyout
- Specified by:
next
in classModel
- Returns:
true
ifnext()
can be called again. Typicallyfalse
is returned if the model requires attention, such as the following conditions:- the model has converged
- the model turned monomorphic (stops only if requested)
- a statistics sample is available
- a preset time has been reached
- See Also:
-
checkConvergence
protected boolean checkConvergence(double dist2) Helper method to check whether the squared distancedist2
qualifies to signal convergence.Implementation Notes:
- ODE's can converge even with non-zero mutation rates
- step size may keep decreasing without ever reaching the accuracy
threshold (or reach step). This scenario requires an emergency brake
(otherwise the browser becomes completely unresponsive - nasty!). Emergency
brake triggered if
stepTaken
is so small that more than1/(accuracy * stepTaken)
updates are required to completedt
. stepTaken < 0
may hold but irrelevant because only the squared value is used.
- Parameters:
dist2
- the squared distance between the current and previous states(y[t+dt]-y[t])2
.- Returns:
true
if convergence criteria passed
-
isMonomorphic
public boolean isMonomorphic()Check if population is monomorphic. Note, in multi-species modules all species need to be monomorphic.- Returns:
true
if criteria for monomorphic state passed
-
isMonomorphic
private boolean isMonomorphic(int start, int stop, int dep, int vac) Helper method to check if the species with trait indicesstart
throughstop
is monomorphic. The species may have dependent traitdep
and allow for vacant sites (with indexvac
).- Parameters:
start
- the start index of the speciesstop
- the stop index of the speciesdep
- the index of the dependent traitvac
- the index of the vacant trait- Returns:
true
if the species is monomorphic
-
deStep
protected double deStep(double step) Attempts a numerical integration step of sizestep
. The baseline are steps of fixed size following Euler's method.- Parameters:
step
- the time step to attempt- Returns:
- squared distance between this state and previous one,
(yt-yout)2
.
-
normalizeState
protected void normalizeState(double[] state) Convenience method to normalize state in frequency based models. For density models nothing needs to be done.- Parameters:
state
- the array that needs to be normalized if appropriate
-
getDerivatives
protected void getDerivatives(double t, double[] state, double[] fitness, double[] change) Calculate the rates of change for all species instate
attime
given the fitnessfitness
and returned inchange
. For replicator models the dynamics depends on the selected type of player updating, seePlayerUpdate.Type
, while for modules with variable population sizes (density based or with vaccant 'space' (reproductive opportunities)) the fitness denotes the rate of reproduction moderated by the available 'space'.IMPORTANT: parallel processing for PDE's in JRE requires this method to be thread safe.
Note: adding
synchronized
to this method would seem to make sense but this results in a deadlock for multiple threads in PDE's and JRE. However, apparently not needed as javascript and java yield the exact same results (after 15k generations!)- Parameters:
t
- the time at which to calculate the derivativestate
- the array of frequencies/densities denoting the state populationfitness
- the array of fitness values of types in populationchange
- the array to return the rate of change of each type in the population- See Also:
-
updateEcology
protected double updateEcology(Module mod, double[] state, double[] fitness, int nGroup, int index, double[] change) Implementation of the player updates for modules with populations of variable size (density based or with vacant 'space', i.e. reproductive opportunities). In models with varying reproductive opportunities the rate of reproduction is given by the product of the individual's fitness moderated by the available 'space'. This results in the dynamical equation given by where denotes vacant space, the fitness of type and is the death rate. This works nicely for replicator type equations but may not fit more general dynamical equations.- Parameters:
mod
- the module representing the current speciesstate
- array of frequencies/densities denoting the state populationfitness
- array of fitness values of types in populationnGroup
- the interaction group sizeindex
- the index of the modulemod
in multi-species moduleschange
- array to return the rate of change of each type in the population- Returns:
- the total change (should be zero in theory)
-
updateThermal
protected double updateThermal(Module mod, double[] state, double[] fitness, int nGroup, int index, double[] change) Implementation of the player updatePlayerUpdate.Type.THERMAL
. This calculates the rates of change for each type in speciesmod
for the popular choice for 'pairwise comparisons' where the focal player and one of its neighbours are randomly chosen. The focal player adopts the strategy of player with a probability where denotes the selection strength and the payoffs of individuals and , respectively. The resulting dynamics for the frequencies of the different strategic types is then given by- Parameters:
mod
- the module representing the current speciesstate
- array of frequencies/densities denoting the state populationfitness
- array of fitness values of types in populationnGroup
- the interaction group sizeindex
- the index of the modulemod
in multi-species moduleschange
- array to return the rate of change of each type in the population- Returns:
- the total change (should be zero in theory)
- See Also:
-
updateBest
protected double updateBest(Module mod, double[] state, double[] fitness, int nGroup, int index, double[] change) Implementation of player updatePlayerUpdate.Type.BEST
. This calculates the rate of change individuals adopt the strategy of better performing individuals with certainty (and never those of worse performing individuals). This calculates the rates of change for each type in speciesmod
for the popular choice for 'pairwise comparisons' where the focal player and one of its neighbours are randomly chosen. The focal player adopts the strategy of player if but never adopts those of a player with , where denote the fitness of players and , respectively. In case of a tie the individual sticks to its strategy. More specifically, the probability to adopt the strategy of is given by where denotes the Heaviside step function with for , for . In principle but assuming that players need at least a marginal incentive to switch strategies we set . In most cases this choice is inconsequential. The resulting dynamics for the frequencies of the different strategic types is then given byNote:In the limit of vanishing noise the updates
updateThermal(Module, double[], double[], int, int, double[])
andupdateImitate(Module, double[], double[], int, int, double[])
recover thePlayerUpdate.Type.BEST
updating type as well.- Parameters:
mod
- the module representing the current speciesstate
- array of frequencies/densities denoting the state populationfitness
- array of fitness values of types in populationnGroup
- the interaction group sizeindex
- the index of the modulemod
in multi-species moduleschange
- array to return the rate of change of each type in the population- Returns:
- the total change (should be zero in theory)
-
updateImitate
protected double updateImitate(Module mod, double[] state, double[] fitness, int nGroup, int index, double[] change) Implementation of the player updatePlayerUpdate.Type.IMITATE
. This calculates the rates of change for each type in speciesmod
for the popular choice for 'pairwise comparisons' where the focal player and one of its neighbours are randomly chosen. The focal player adopts the strategy of with a probability proportional to the payoff difference : where denote the fitness of players and , respectively. The resulting dynamics for the frequencies of the different strategic types is then given by the standard replicator equation where denotes the the average population payoff.Note: in multi-species models all rates of change are scaled by the maximum fitness difference for all species to preserve their relative time scales.
- Parameters:
mod
- the module representing the current speciesstate
- array of frequencies/densities denoting the state populationfitness
- array of fitness values of types in populationnGroup
- the interaction group sizeindex
- the index of the modulemod
in multi-species moduleschange
- array to return the rate of change of each type in the population- Returns:
- the total change (should be zero in theory)
-
updateImitateBetter
protected double updateImitateBetter(Module mod, double[] state, double[] fitness, int nGroup, int index, double[] change) Implementation of the player updatePlayerUpdate.Type.IMITATE_BETTER
. This calculates the rates of change for each type in speciesmod
for the popular choice for 'pairwise comparisons' where the focal player and one of its neighbours are randomly chosen. The focal player adopts the strategy of a better performing player with a probability proportional to the payoff difference : Incidentally and somewhat surprisingly this update also reduces to the standard replicator dynamics albeit with a constant rescaling of time such that the rate of change are half of the standard replicator equation obtained fromPlayerUpdate.Type#IMITATE
.Note: in multi-species models all rates of change are scaled by the maximum fitness difference for all species to preserve their relative time scales.
- Parameters:
mod
- the module representing the current speciesstate
- array of frequencies/densities denoting the state populationfitness
- array of fitness values of types in populationnGroup
- the interaction group sizeindex
- the index of the modulemod
in multi-species moduleschange
- array to return the rate of change of each type in the population- Returns:
- the total change (should be zero in theory)
-
updateReplicate
private double updateReplicate(Module mod, double[] state, double[] fitness, int nGroup, int index, double[] change, double noise) Helper method to calculate the rate of change for the standard replicator dynamics with different amounts of noise arising from the microscopic update rule. In the microscopic implementation of the replicator dynamics the noise arises from focal individuals that adopt the strategy of a neighbour with a probability proportional to the payoff difference , where refer to the respective payoffs of and . However, the noise is cut in half if the focal imitates only better performing neighbours, i.e. proportional to .- Parameters:
mod
- the module representing the current speciesstate
- array of frequencies/densities denoting the state populationfitness
- array of fitness values of types in populationnGroup
- the interaction group sizeindex
- the index of the modulemod
in multi-species moduleschange
- array to return the rate of change of each type in the populationnoise
- the noise arising from probabilistical updates- Returns:
- the total change (should be zero in theory)
-
updateProportional
protected double updateProportional(Module mod, double[] state, double[] fitness, int nGroup, int index, double[] change) Implementation of the player updatePlayerUpdate.Type.PROPORTIONAL
. This calculates the rates of change for each type in speciesmod
for a 'pairwise comparison' where the focal player and one of its neighbours are randomly chosen. The focal player adopts the strategy of a player with probability: where refer to the respective payoffs of and . In the continuum limit this yields- Parameters:
mod
- the module representing the current speciesstate
- array of frequencies/densities denoting the state populationfitness
- array of fitness values of types in populationnGroup
- the interaction group sizeindex
- the index of the modulemod
in multi-species moduleschange
- array to return the rate of change of each type in the population- Returns:
- the total change (should be zero in theory)
-
updateBestResponse
protected double updateBestResponse(Module mod, double[] state, double[] fitness, int nGroup, int index, double[] change) Implementation of the player updatePlayerUpdate.Type.BEST_RESPONSE
. This calculates the rates of change for each type in speciesmod
for the best-response dynamic where the focal player switches to the strategy that has the highest payoff given the current state of the population. Note, in contrast to other player updates, such asPlayerUpdate.Type.THERMAL
orPlayerUpdate.Type.IMITATE
the best-response is an innovative update rule, which means that strategic types that are currently not present in the population can get introduced. Consequently, homogeneous states may not be absorbing.Note: inactive strategies are excluded and do not qualify as a best response.
IMPORTANT: needs to be thread safe (must supply memory for calculations and results)
- Parameters:
mod
- the module representing the current speciesstate
- array of frequencies/densities denoting the state populationfitness
- array of fitness values of types in populationnGroup
- the interaction group sizeindex
- the index of the modulemod
in multi-species moduleschange
- array to return the rate of change of each type in the population- Returns:
- the total change (should be zero in theory)
-
update
public void update()Description copied from class:Model
Update this model. For example called after initialization and when parameters changed. -
getStatus
Description copied from class:Model
Returns status message from model. Typically this is a string summarizing the current state of the simulation. For example, models with discrete strategy sets (such as 2x2 games, seeTBT
) return the average frequencies of each strategy type in the population(s), seeIBSDPopulation
. Similarly, models with continuous strategies (such as continuous snowdrift games, seeCSD
) return the mean, minimum and maximum trait value(s) in the population(s), seeIBSMCPopulation
. The status message is displayed along the bottom of the GUI.Note: if the model runs into difficulties, problems should be reported through the logging mechanism. Messages with severity
Level.WARNING
or higher are displayed in the status of the GUI and override status messages returned here. -
permitsTimeReversal
public boolean permitsTimeReversal()Checks if time reversal is permitted. By default returnsfalse
. Only few models are capable of time reversal.ODE and SDE models return
true
by default.- Overrides:
permitsTimeReversal
in classModel
- Returns:
true
if time reversal permissible.- See Also:
-
isTimeReversed
public boolean isTimeReversed()Description copied from class:Model
Checks if time is reversed. By default returnsfalse
. Only few models are capable of time reversal.- Overrides:
isTimeReversed
in classModel
- Returns:
true
if time is reversed- See Also:
-
setTimeReversed
public void setTimeReversed(boolean reversed) Description copied from class:Model
Request time reversal ifreversed==true
. The model may not be able to honour the request. However, some models allow to travel back in time. In general, this is only possible for ODE and SDE models and even for those it may not be feasible due to details of the dynamics, such as dissipative terms in the differential equations. For example, ininvalid reference
org.evoludo.simulator.modules.CG
- Overrides:
setTimeReversed
in classModel
- Parameters:
reversed
- the request whether time should be reversed.- See Also:
-
setAdjustedDynamics
public void setAdjustedDynamics(boolean adjusted) Sets whether adjusted dynamics is requested. It may not be possible to honour the request, e.g. if payoffs can be negative.- Parameters:
adjusted
- iftrue
requests adjusted dynamics
-
setAccuracy
public void setAccuracy(double acc) Sets the desired accuracy for determining convergence. Ify(t+dt)-y(t)<acc dt
holds, wherey(t+dt)
denotes the new state andy(t)
the previous state, the numerical integration is reported as having converged and stops.- Parameters:
acc
- the numerical accuracy
-
getAccuracy
public double getAccuracy()Gets the numerical accuracy for determining convergence.- Returns:
- the numerical accuracy
-
init
private void init(boolean doRandom) Helper method to set the intital state of the model. The method initializes the state vectoryt
and the initial state vectory0
. IfdoRandom
istrue
random initial frequencies are set for species withInitType.RANDOM
.- Parameters:
doRandom
- iftrue
use random initial frequencies as requested
-
parse
Parse initializer stringarg
. Determine type of initialization and process its arguments as appropriate.Note: Not possible to perform parsing in
CLODelegate
ofcloInit
because PDE model provide their ownPDE.InitType
s.- Parameters:
arg
- the arguments to parse- Returns:
true
if parsing successful- See Also:
-
collectCLO
Description copied from interface:CLOProvider
All providers of command line options must implement this method to collect their options.Each command line option is (uniquely) identified by it's name (see
CLOption.getName()
), which corresponds to the long version of the option. If an attempt is made to add an option with a name that already exists, theparser
issues a warning and ignores the option. Thus, in general, implementing subclasses should first register their options and callsuper.collectCLO(CLOParser)
at the end such that subclasses are able to override command line options specified in a parental class.Override this method in subclasses to add further command line options. Subclasses must make sure that they include a call to super.
- Specified by:
collectCLO
in interfaceCLOProvider
- Overrides:
collectCLO
in classModel
- Parameters:
parser
- the reference to parser that manages command line options- See Also:
-
encodeState
Description copied from class:Model
Encode the state of the model in aplist
inspiredXML
string. This allows to save the state and restore later with the exact same results as when continuing to run the model. This even allows to switch from JRE to GWT or back and obtain identical results!- Specified by:
encodeState
in classModel
- Parameters:
plist
- theStringBuilder
to write the encoded state to- See Also:
-
restoreState
Description copied from class:Model
Restore the state encoded in theplist
inspiredmap
ofkey, value
-pairs.- Specified by:
restoreState
in classModel
- Parameters:
plist
- the map ofkey, value
-pairs- Returns:
true
if successful- See Also:
-
encodeStrategies
Encodes state of the model in the form of aplist
string.- Parameters:
plist
- the string builder for the encoded state
-
restoreStrategies
Restores the state encoded inplist
.- Parameters:
plist
- the encoded state- Returns:
true
if successful
-
encodeFitness
Encodes the fitness of the model in the form of aplist
string.- Parameters:
plist
- the string builder for the encoded state
-
restoreFitness
Restores the fitness encoded inplist
.- Parameters:
plist
- the encoded state- Returns:
true
if successful
-