Class ODE

All Implemented Interfaces:
Discrete, CLOProvider
Direct Known Subclasses:
PDE, RungeKutta, SDE

public class ODE extends Model implements Discrete
Common base class for all differential equations models. Provides the basic Euler method for the numerical integration of systems of ordinary differential equations.

Note: Currently differential equation models are restricted to discrete strategy sets.

Author:
Christoph Hauert
  • Field Details Link icon

    • dt Link icon

      double dt
      Discretization 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 Link icon

      double dtTry
      The attempted size of next step to take. In ODE this is always equal to dt 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 Link icon

      double dtTaken
      The 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 Link icon

      boolean forward
      The flag is true if the numerical integration moves forward in time.
      See Also:
    • nDim Link icon

      int nDim
      Total number of traits (dynamical variables) in the ODE.
    • y0 Link icon

      double[] y0
      The 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 in yt[0] through yt[n1] where n1 denotes the number of traits in the first species. The initial frequencies/densities of the second species are stored in yt[n1+1] through yt[n1+n2] where n2 denotes the number of traits in the second species, etc.
    • yt Link icon

      double[] yt
      The 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 in yt[0] through yt[n1] where n1 denotes the number of traits in the first species. The current frequencies/densities of the second species are stored in yt[n1+1] through yt[n1+n2] where n2 denotes the number of traits in the second species, etc.
    • ft Link icon

      double[] ft
      The 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 in ft[0] through ft[n1] where n1 denotes the number of traits in the first species. The current fitness of the second species are stored in ft[n1+1] through ft[n1+n2] where n2 denotes the number of traits in the second species, etc.
    • dyt Link icon

      double[] dyt
      The frequency/density increments of the ODE between the current state yt and the previous state yout. In multi-species modules the increments of each species are concatenated into this single array. The increments in the first species are stored in dyt[0] through dyt[n1] where n1 denotes the number of traits in the first species. The increments in the second species are stored in dyt[n1+1] through dyt[n1+n2] where n2 denotes the number of traits in the second species, etc.
    • yout Link icon

      double[] yout
      The 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 in yout[0] through yout[n1] where n1 denotes the number of traits in the first species. The current frequencies/densities of the second species are stored in yout[n1+1] through yout[n1+n2] where n2 denotes the number of traits in the second species, etc.
    • dstate Link icon

      double[] dstate
      Pointer to temporarily store a reference to the current state yt. This affects only multi-species modules because fitness calculations for the derivatives may need access to the state of the other populations. To allow this getMeanTraits(int, double[]) uses dstate while the derivative calculations are in progress (dstate != null). Otherwise yt is used.
    • staticfit Link icon

      double[] staticfit
      For Modules with static fitness, e.g. the original Moran process Moran, 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 in staticfit[0] through staticfit[n1] where n1 denotes the number of traits in the first species. The static fitness of the second species are stored in staticfit[n1+1] through staticfit[n1+n2] where n2 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 or null if no species has static fitness.
    • initType Link icon

      protected ODE.InitType[] initType
      Type of initial configuration for each species.
      See Also:
    • idxSpecies Link icon

      int[] idxSpecies
      Array containing indices that delimit individual species in the ODE state variables. In a module with N species the array has N+1 entries. For example, for a single species module with n1 traits idxSpecies = int[] { 0, n1 } holds, whereas for a two species module with n1 and n2 traits, respectively, idxSpecies = int[] { 0, n1, n1+n2 } holds, and for three species idxSpecies = 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 species i are stored in yt[idxSpecies[i]] through yt[idxSpecies[i+1]].
      See Also:
    • rates Link icon

      double[] rates
      Array containing the update rates for each species. Determines the relative update rate between species, Only relevant in multi-species modules.
      See Also:
    • mutation Link icon

      Mutation.Discrete[] mutation
      Array containing the mutation rates/probabilities for each species.
      See Also:
    • dependents Link icon

      int[] dependents
      Array 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 Link icon

      double[] invFitRange
      Array 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 Link icon

      String[] names
      The 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 in names[0] through names[n1] where n1 denotes the number of traits in the first species. The names of the second species are stored in names[n1+1] through names[n1+n2] where n2 denotes the number of traits in the second species, etc.
    • monoStop Link icon

      protected boolean monoStop
      Set to true to requests that the numerical integration stops once the population reaches a monomorphic state.
      See Also:
    • isAdjustedDynamics Link icon

      boolean isAdjustedDynamics
      true 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 Link icon

      double accuracy
      Desired accuracy to determine whether numerical integration has converged or a population has become monomorphic.
      See Also:
    • cloDEdt Link icon

      public final CLOption cloDEdt
      Command line option to set the time increment dt in DE models. If the requested dt turns out to be too big, e.g. due to diffusion in PDE, it is reduced appropriately. For RungeKutta this represents the initial step.
      See Also:
    • cloInit Link icon

      public final CLOption 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 Link icon

      public final CLOption cloAdjustedDynamics
      Command line option to activate adjusted replicator dynamics (if possible) in ODE models.
      See Also:
    • cloDEAccuracy Link icon

      public final CLOption cloDEAccuracy
      Command line option to set the desired accuracy of ODE calculations as criteria for convergence and whether population is monomorphic.
      See Also:
    • cloTimeReversed Link icon

      public final CLOption cloTimeReversed
      Command line option to set the number of generations between reports for EvoLudo.modelNext().
  • Constructor Details Link icon

    • ODE Link icon

      public ODE(EvoLudo engine)
      Constructs a new model for the numerical integration of the system of ordinary differential equations representing the dynamics specified by the Module module using the EvoLudo pacemaker engine to control the numerical evaluations. The integrator implements Euler's method (fixed step size).
      Parameters:
      engine - the pacemaker for running the model
  • Method Details Link icon

    • load Link icon

      public void load()
      Description copied from class: Model
      Milestone: Load this model and allocate resources (if applicable).
      Overrides:
      load in class Model
      See Also:
    • unload Link icon

      public void unload()
      Description copied from class: Model
      Milestone: Unload this model and free resources (if applicable).
      Overrides:
      unload in class Model
      See Also:
    • check Link icon

      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 through logger. 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.
      Overrides:
      check in class Model
      Returns:
      true if reset required
      See Also:
    • reset Link icon

      public void reset()
      Description copied from class: Model
      Milestone: Reset this model
      Overrides:
      reset in class Model
      See Also:
    • init Link icon

      public void init()
      Description copied from class: Model
      Milestone: Initialize this model
      Overrides:
      init in class Model
      See Also:
    • setDt Link icon

      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 (see PDE.checkDt()) or choose a variable step size in which case deltat is ignored (see RungeKutta.getAutoDt()).

      Parameters:
      deltat - the time increments in continuous time models.
    • getDt Link icon

      public double getDt()
      Gets the discretization of time increments in continuous time models.

      Note: This may be different from dt set through setDt(double) if the model required adjustments.

      Returns:
      the time increment in continuous time models.
    • isDensity Link icon

      public boolean isDensity()
      Return whether this DE model tracks frequencies or densities. Returns false (i.e. frequency based model) by default.
      Returns:
      true if state refers to densities.
    • setInitialTraits Link icon

      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 class Model
      Parameters:
      init - the array with the initial trait values
      Returns:
      true if initial traits successfully set
    • getInitialTraits Link icon

      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 class Model
      Parameters:
      init - the array for storing the initial trait values
    • setInitialTraits Link icon

      public boolean setInitialTraits(int id, double[] init)
      Description copied from class: Model
      Set initial trait values for species with ID id.
      Overrides:
      setInitialTraits in class Model
      Parameters:
      id - the species identifier
      init - the array with the initial trait values
      Returns:
      true if initial traits successfully set
    • getInitialTraits Link icon

      public void getInitialTraits(int id, double[] init)
      Description copied from class: Model
      Return initial trait values for species with ID id.
      Specified by:
      getInitialTraits in class Model
      Parameters:
      id - the species identifier
      init - the array for storing the initial trait values
    • getMonoScore Link icon

      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/strategy idx but also deals with payoff accounting (averaged versus accumulated).
      Specified by:
      getMonoScore in interface Discrete
      Parameters:
      id - the id of the population for multi-species models
      idx - trait/strategy
      Returns:
      payoff/score in monomorphic population with trait/strategy idx. Returns NaN if scores ill defined
      See Also:
    • getNMean Link icon

      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
      Specified by:
      getNMean in class Model
      Returns:
      the number of mean values for all species
    • getNMean Link icon

      public int getNMean(int id)
      Description copied from class: Model
      Return the number of mean trait values for species with ID id.
      Specified by:
      getNMean in class Model
      Parameters:
      id - the species identifier
      Returns:
      the number of mean values for species id
    • getMeanTraits Link icon

      public boolean getMeanTraits(int id, double[] mean)
      Description copied from class: Model
      Return mean trait values for species with ID id.
      Specified by:
      getMeanTraits in class Model
      Parameters:
      id - the species identifier
      mean - 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 Link icon

      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 class Model
      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 Link icon

      public <T> void getTraitData(int id, T[] colors, ColorMap<T> colorMap)
      Unused interface method.
      Specified by:
      getTraitData in class Model
      Type Parameters:
      T - color data type. Color for PopGraph2D and PopGraph2D as well as MeshLambertMaterial for PopGraph3D
      Parameters:
      id - the species identifier
      colors - the array for storing the colors for individuals
      colorMap - the map for translating individual traits into colors
    • getMinScore Link icon

      public double getMinScore(int id)
      Description copied from class: Model
      Returns the minimum score that individuals of species with ID id can achieve in this model. Takes into account potential adjustments due to population structure and payoff accounting.
      Specified by:
      getMinScore in class Model
      Parameters:
      id - the id of the population for multi-species models
      Returns:
      the minimum score
      See Also:
    • getMaxScore Link icon

      public double getMaxScore(int id)
      Description copied from class: Model
      Returns the maximum score that individuals of species with ID id can achieve in this model. Takes into account potential adjustments due to population structure and payoff accounting.
      Specified by:
      getMaxScore in class Model
      Parameters:
      id - the id of the population for multi-species models
      Returns:
      the maximum score
      See Also:
    • getMinFitness Link icon

      public double getMinFitness(int id)
      Description copied from class: Model
      Calculates and returns the absolute fitness minimum. This is important to
      1. determine probabilities or rates for adopting the strategy of another player,
      2. optimize fitness based picking of individuals, and
      3. scaling graphical output.
      Specified by:
      getMinFitness in class Model
      Parameters:
      id - the id of the population for multi-species models
      Returns:
      the minimum fitness
      See Also:
    • getMaxFitness Link icon

      public double getMaxFitness(int id)
      Description copied from class: Model
      Calculates and returns the absolute fitness maximum. This is important to
      1. determine probabilities or rates for adopting the strategy of another player,
      2. optimize fitness based picking of individuals, and
      3. scaling graphical output.
      Specified by:
      getMaxFitness in class Model
      Parameters:
      id - the id of the population for multi-species models
      Returns:
      the maximum fitness
      See Also:
    • getMeanFitness Link icon

      public boolean getMeanFitness(int id, double[] mean)
      Description copied from class: Model
      Gets the mean fitness values for species with ID id. The result is stored and returned in mean. Used by GUI to visualize local dynamics at idx.
      Specified by:
      getMeanFitness in class Model
      Parameters:
      id - the species identifier
      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.
    • getMeanFitness Link icon

      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 in mean. 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 class Model
      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 Link icon

      public <T> void getFitnessData(int id, T[] colors, ColorMap.Gradient1D<T> colorMap)
      Description copied from class: Model
      Translates fitness data into colors using ColorMap colorMap. Used by GUI to visualize current state of model.
      Specified by:
      getFitnessData in class Model
      Type Parameters:
      T - color data type. Color for PopGraph2D and PopGraph2D as well as MeshLambertMaterial for PopGraph3D
      Parameters:
      id - the species identifier
      colors - the array for storing color values
      colorMap - the map to use for translating traits to colors
    • getFitnessHistogramData Link icon

      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 array bins. 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 class Model
      Parameters:
      id - the species identifier
      bins - the array for storing histogram. For Discrete modules this is always one dimensional
    • next Link icon

      public boolean next()
      Advance model by one step. The details of what happens during one step depends on the models Type as well as its Mode.

      Implementation Notes: Link icon

      Before:
      1. yt current state
      2. ft current fitness
      3. dyt current derivatives fitness
      4. dtTry size of time step attempting to make
      After:
      1. yt next state
      2. ft next fitness
      3. dyt next derivatives
      4. yout contains previous state (i.e. yt when entering method)
      5. stepTaken the time between yt and yout
      Specified by:
      next in class Model
      Returns:
      true if next() can be called again. Typically false 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 Link icon

      protected boolean checkConvergence(double dist2)
      Helper method to check whether the squared distance dist2 qualifies to signal convergence.

      Implementation Notes: Link icon

      1. ODE's can converge even with non-zero mutation rates
      2. 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 than 1/(accuracy * stepTaken) updates are required to complete dt.
      3. 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 Link icon

      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 Link icon

      private boolean isMonomorphic(int start, int stop, int dep, int vac)
      Helper method to check if the species with trait indices start through stop is monomorphic. The species may have dependent trait dep and allow for vacant sites (with index vac).
      Parameters:
      start - the start index of the species
      stop - the stop index of the species
      dep - the index of the dependent trait
      vac - the index of the vacant trait
      Returns:
      true if the species is monomorphic
    • deStep Link icon

      protected double deStep(double step)
      Attempts a numerical integration step of size step. 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 Link icon

      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 Link icon

      protected void getDerivatives(double t, double[] state, double[] fitness, double[] change)
      Calculate the rates of change for all species in state at time given the fitness fitness and returned in change. For replicator models the dynamics depends on the selected type of player updating, see PlayerUpdate.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 derivative
      state - the array of frequencies/densities denoting the state population
      fitness - the array of fitness values of types in population
      change - the array to return the rate of change of each type in the population
      See Also:
    • updateEcology Link icon

      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 x˙=xi(zfid) where z denotes vacant space, fi the fitness of type i and d 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 species
      state - array of frequencies/densities denoting the state population
      fitness - array of fitness values of types in population
      nGroup - the interaction group size
      index - the index of the module mod in multi-species modules
      change - array to return the rate of change of each type in the population
      Returns:
      the total change (should be zero in theory)
    • updateThermal Link icon

      protected double updateThermal(Module mod, double[] state, double[] fitness, int nGroup, int index, double[] change)
      Implementation of the player update PlayerUpdate.Type.THERMAL. This calculates the rates of change for each type in species mod for the popular choice for 'pairwise comparisons' where the focal player i and one of its neighbours j are randomly chosen. The focal player i adopts the strategy of player j with a probability pij=11+exp[w(fifj)], where w0 denotes the selection strength and fi,fj the payoffs of individuals i and j, respectively. The resulting dynamics for the frequencies of the different strategic types is then given by x˙i=j=1nxixj(11+exp[w(fifj)]11+exp[w(fjfi)])=xij=1ntanh[w(fifj)/2].
      Parameters:
      mod - the module representing the current species
      state - array of frequencies/densities denoting the state population
      fitness - array of fitness values of types in population
      nGroup - the interaction group size
      index - the index of the module mod in multi-species modules
      change - 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 Link icon

      protected double updateBest(Module mod, double[] state, double[] fitness, int nGroup, int index, double[] change)
      Implementation of player update PlayerUpdate.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 species mod for the popular choice for 'pairwise comparisons' where the focal player i and one of its neighbours j are randomly chosen. The focal player i adopts the strategy of player j if fj>fi but never adopts those of a player j with fj<fi, where fi,fj denote the fitness of players i and j, respectively. In case of a tie fj=fi the individual sticks to its strategy. More specifically, the probability to adopt the strategy of j is given by pij=θ(fjfi), where θ(x) denotes the Heaviside step function with θ(x)=0 for x<0, θ(x)=1 for x>0. In principle θ(0)=1/2 but assuming that players need at least a marginal incentive to switch strategies we set θ(0)=0. In most cases this choice is inconsequential. The resulting dynamics for the frequencies of the different strategic types is then given by x˙i=j=1nxixjθ[fjfi)].

      Note:In the limit of vanishing noise the updates updateThermal(Module, double[], double[], int, int, double[]) and updateImitate(Module, double[], double[], int, int, double[]) recover the PlayerUpdate.Type.BEST updating type as well.

      Parameters:
      mod - the module representing the current species
      state - array of frequencies/densities denoting the state population
      fitness - array of fitness values of types in population
      nGroup - the interaction group size
      index - the index of the module mod in multi-species modules
      change - array to return the rate of change of each type in the population
      Returns:
      the total change (should be zero in theory)
    • updateImitate Link icon

      protected double updateImitate(Module mod, double[] state, double[] fitness, int nGroup, int index, double[] change)
      Implementation of the player update PlayerUpdate.Type.IMITATE. This calculates the rates of change for each type in species mod for the popular choice for 'pairwise comparisons' where the focal player i and one of its neighbours j are randomly chosen. The focal player i adopts the strategy of j with a probability proportional to the payoff difference fjfi: pij=1/2(1+(fjfi)/(fj+fi)), where fi,fj denote the fitness of players i and j, respectively. The resulting dynamics for the frequencies of the different strategic types is then given by the standard replicator equation x˙i=xi(fif¯) where f¯ 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 species
      state - array of frequencies/densities denoting the state population
      fitness - array of fitness values of types in population
      nGroup - the interaction group size
      index - the index of the module mod in multi-species modules
      change - array to return the rate of change of each type in the population
      Returns:
      the total change (should be zero in theory)
    • updateImitateBetter Link icon

      protected double updateImitateBetter(Module mod, double[] state, double[] fitness, int nGroup, int index, double[] change)
      Implementation of the player update PlayerUpdate.Type.IMITATE_BETTER. This calculates the rates of change for each type in species mod for the popular choice for 'pairwise comparisons' where the focal player i and one of its neighbours j are randomly chosen. The focal player i adopts the strategy of a better performing player j with a probability proportional to the payoff difference fjfi: x˙i=j=1nxixj(fifj)+=j=1nxixj(fifj)j=1nxixj(fifj)=j=1nxixj(fifj)/2=xi(fif¯)/2. 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 from PlayerUpdate.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 species
      state - array of frequencies/densities denoting the state population
      fitness - array of fitness values of types in population
      nGroup - the interaction group size
      index - the index of the module mod in multi-species modules
      change - array to return the rate of change of each type in the population
      Returns:
      the total change (should be zero in theory)
    • updateReplicate Link icon

      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 i that adopt the strategy of a neighbour j with a probability proportional to the payoff difference fjfi, where fi,fj refer to the respective payoffs of i and j. However, the noise is cut in half if the focal imitates only better performing neighbours, i.e. proportional to (fjfi)+.
      Parameters:
      mod - the module representing the current species
      state - array of frequencies/densities denoting the state population
      fitness - array of fitness values of types in population
      nGroup - the interaction group size
      index - the index of the module mod in multi-species modules
      change - array to return the rate of change of each type in the population
      noise - the noise arising from probabilistical updates
      Returns:
      the total change (should be zero in theory)
    • updateProportional Link icon

      protected double updateProportional(Module mod, double[] state, double[] fitness, int nGroup, int index, double[] change)
      Implementation of the player update PlayerUpdate.Type.PROPORTIONAL. This calculates the rates of change for each type in species mod for a 'pairwise comparison' where the focal player i and one of its neighbours j are randomly chosen. The focal player i adopts the strategy of a player j with probability: pij=fj/(fi+fj), where fi,fj refer to the respective payoffs of i and j. In the continuum limit this yields x˙i=xij=1nxjfifjfi+fj.
      Parameters:
      mod - the module representing the current species
      state - array of frequencies/densities denoting the state population
      fitness - array of fitness values of types in population
      nGroup - the interaction group size
      index - the index of the module mod in multi-species modules
      change - array to return the rate of change of each type in the population
      Returns:
      the total change (should be zero in theory)
    • updateBestResponse Link icon

      protected double updateBestResponse(Module mod, double[] state, double[] fitness, int nGroup, int index, double[] change)
      Implementation of the player update PlayerUpdate.Type.BEST_RESPONSE. This calculates the rates of change for each type in species mod for the best-response dynamic where the focal player i switches to the strategy j that has the highest payoff given the current state of the population. Note, in contrast to other player updates, such as PlayerUpdate.Type.THERMAL or PlayerUpdate.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 species
      state - array of frequencies/densities denoting the state population
      fitness - array of fitness values of types in population
      nGroup - the interaction group size
      index - the index of the module mod in multi-species modules
      change - array to return the rate of change of each type in the population
      Returns:
      the total change (should be zero in theory)
    • update Link icon

      public void update()
      Description copied from class: Model
      Update this model. For example called after initialization and when parameters changed.
      Specified by:
      update in class Model
      See Also:
    • getStatus Link icon

      public String 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, see TBT) return the average frequencies of each strategy type in the population(s), see IBSDPopulation. Similarly, models with continuous strategies (such as continuous snowdrift games, see CSD) return the mean, minimum and maximum trait value(s) in the population(s), see IBSMCPopulation. 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.

      Specified by:
      getStatus in class Model
      Returns:
      status of active model
    • permitsTimeReversal Link icon

      public boolean permitsTimeReversal()
      Checks if time reversal is permitted. By default returns false. Only few models are capable of time reversal.

      ODE and SDE models return true by default.

      Overrides:
      permitsTimeReversal in class Model
      Returns:
      true if time reversal permissible.
      See Also:
    • isTimeReversed Link icon

      public boolean isTimeReversed()
      Description copied from class: Model
      Checks if time is reversed. By default returns false. Only few models are capable of time reversal.
      Overrides:
      isTimeReversed in class Model
      Returns:
      true if time is reversed
      See Also:
    • setTimeReversed Link icon

      public void setTimeReversed(boolean reversed)
      Description copied from class: Model
      Request time reversal if reversed==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, in
      invalid reference
      org.evoludo.simulator.modules.CG
      the ecological dynamics of the patch quality prevents time reversal, i.e. results are numerically unstable due to exponential amplification of deviations in dissipative systems.
      Overrides:
      setTimeReversed in class Model
      Parameters:
      reversed - the request whether time should be reversed.
      See Also:
    • setAdjustedDynamics Link icon

      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 - if true requests adjusted dynamics
    • setAccuracy Link icon

      public void setAccuracy(double acc)
      Sets the desired accuracy for determining convergence. If y(t+dt)-y(t)<acc dt holds, where y(t+dt) denotes the new state and y(t) the previous state, the numerical integration is reported as having converged and stops.
      Parameters:
      acc - the numerical accuracy
    • getAccuracy Link icon

      public double getAccuracy()
      Gets the numerical accuracy for determining convergence.
      Returns:
      the numerical accuracy
    • init Link icon

      private void init(boolean doRandom)
      Helper method to set the intital state of the model. The method initializes the state vector yt and the initial state vector y0. If doRandom is true random initial frequencies are set for species with InitType.RANDOM.
      Parameters:
      doRandom - if true use random initial frequencies as requested
    • parse Link icon

      public boolean parse(String arg)
      Parse initializer string arg. Determine type of initialization and process its arguments as appropriate.

      Note: Not possible to perform parsing in CLODelegate of cloInit because PDE model provide their own PDE.InitTypes.

      Parameters:
      arg - the arguments to parse
      Returns:
      true if parsing successful
      See Also:
    • collectCLO Link icon

      public void collectCLO(CLOParser parser)
      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, the parser issues a warning and ignores the option. Thus, in general, implementing subclasses should first register their options and call super.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 interface CLOProvider
      Overrides:
      collectCLO in class Model
      Parameters:
      parser - the reference to parser that manages command line options
      See Also:
    • encodeState Link icon

      public void encodeState(StringBuilder plist)
      Description copied from class: Model
      Encode the state of the model in a plist inspired XML 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 class Model
      Parameters:
      plist - the StringBuilder to write the encoded state to
      See Also:
    • restoreState Link icon

      public boolean restoreState(Plist plist)
      Description copied from class: Model
      Restore the state encoded in the plist inspired map of key, value-pairs.
      Specified by:
      restoreState in class Model
      Parameters:
      plist - the map of key, value-pairs
      Returns:
      true if successful
      See Also:
    • encodeStrategies Link icon

      void encodeStrategies(StringBuilder plist)
      Encodes state of the model in the form of a plist string.
      Parameters:
      plist - the string builder for the encoded state
    • restoreStrategies Link icon

      boolean restoreStrategies(Plist plist)
      Restores the state encoded in plist.
      Parameters:
      plist - the encoded state
      Returns:
      true if successful
    • encodeFitness Link icon

      void encodeFitness(StringBuilder plist)
      Encodes the fitness of the model in the form of a plist string.
      Parameters:
      plist - the string builder for the encoded state
    • restoreFitness Link icon

      boolean restoreFitness(Plist plist)
      Restores the fitness encoded in plist.
      Parameters:
      plist - the encoded state
      Returns:
      true if successful