Tips and Tricks¶

Model Initialisation¶

Ensure the neworder.Model base class is properly initialised:

class MyModel(neworder.Model):
  def __init__(self, args...) -> None:
    timeline = ... # initialise an appropriate timeline
    super().__init__(timeline) # (1)!
    ... # now initialise the subclass

1. this line is essential

Custom Seeding Strategies¶

neworder provides three basic seeding functions which initialise the model's random stream so that they are either non-reproducible (neworder.MonteCarlo.nondeterministic_stream), or reproducible and either identical (neworder.MonteCarlo.deterministic_identical_stream) or independent across parallel runs (neworder.MonteCarlo.deterministic_independent_stream). Typically, a user would select identical streams (and perturbed inputs) for sensitivity analysis, and independent streams (with identical inputs) for convergence analysis.

If necessary, you can supply your own seeding strategy, for instance if you required half the processes to have identical streams:

import neworder
def hybrid_seeder() -> int:
    return (neworder.mpi.RANK % 2) + 12345

which returns the same seed for all odd-ranked processes and a different seed for the even-ranked ones. You can define your seeder inline when you instantiate the Model, e.g.

class MyModel(neworder.Model):
    def __init__(self, timeline: neworder.Timeline) -> None:
        super().__init__(timeline, lambda: (neworder.mpi.RANK % 2) + 12345)
        ...

If there was a requirement for multiple processes to all have the same nondeterministic stream, you could implement a seeding strategy like so:

def nondeterministic_identical_stream() -> int:
    # only process 0 gets a seed
    seed = neworder.MonteCarlo.nondeterministic_stream() if neworder.mpi.RANK == 0 else None
    # then broadcasts it to the other processes
    seed = neworder.mpi.COMM.bcast(seed, root=0)
    return seed

Identical Streams¶

The "option" example relies on parallel processes with identical random streams to reduce noise when computing differences for sensitivity analysis. It implements a check step that compares the internal states of the random stream in each process and fails if any are different (see the example code).

Ultimate Reproducibility¶

The MonteCarlo engine is a sequential stream: every draw advances its internal state, so the value you get depends on how many draws have been taken before it. This couples reproducibility to draw order - if agents are added, removed, or processed in a different sequence, the stream diverges.

SplitMix64 eliminates this coupling. Each variate is computed by hashing a set of integer keys (e.g. person ID, process ID, timestep) together with the seed. There is no state to advance, so:

• the draw for person i is the same whether you compute the full population or just person i in isolation,
• draws can be computed in any order, on any thread, without coordination.

Basic usage¶

Construct within your model class, passing a seeder (the same callables used by MonteCarlo):

self.rng = neworder.SplitMix64(neworder.MonteCarlo.deterministic_identical_stream)

Call uarray with any mix of scalar integers (used as context, adding no output dimension) and 1-D integer arrays (each adding one output dimension). Multi-dimensional arrays are not supported and will raise a TypeError.

# 1-D: one variate per person
draws = self.rng.uarray(person_ids, process_id, self.timeline.index)

# 2-D: one variate per (person, draw_index) pair
draws = self.rng.uarray(person_ids, process_id, self.timeline.index, draw_indices)

String keys (e.g. the name of a stochastic process) can be converted to stable integers with SplitMix64.hash64:

draws = self.rng.uarray(person_ids, neworder.SplitMix64.hash64("mortality"), self.timeline.index)

Repeated calls with the same arguments¶

Because SplitMix64 has no state, two calls with identical arguments return identical values. To get independent draws across repeated calls, either:

• use a non-deterministic seeder (this is called each time uarray is called), or
• construct with use_counter=True, which mixes an auto-incrementing counter into each call. In this case calling reset() rewinds the counter and will then replay the same sequence.

External Sources of Randomness¶

Other libraries, such as numpy, contain a much broader selection of random number functionality than neworder does, and it makes no sense to reimplement such functionality. If you are using a specific seeding strategy within neworder, and are also using an external random generator, it is important to ensure they are also following the same strategy, otherwise reproducibility may be compromised.

In your model constructor, you can seed the numpy generator like so

ext_seed = self.mc.raw()
self.nprand = np.random.Generator(np.random.MT19937(ext_seed))
# ...get some values
x = self.nprand.normal(size=5)

If you've chosen a deterministic seeding strategy, then ext_seed will be reproducible, and if you've chosen an independent strategy, then ext_seed will be different for each process, thus propagating your chosen seeding strategy to the external generator.

Using neworder's random generator with numpy¶

It is now possible to use the the neworder model's Monte-Carlo engine as a numpy generator. In this way all of numpy's functionality is available with neworder's MonteCarlo RNG. To achieve this use the adapter function as_np. Similarly to the example above, in your model constructor create the numpy generator, then:

self.nprand = no.as_np(self.mc)
# ...get some values
x = self.nprand.normal(size=5)

NB as there is only one RNG state, you can safely get independent variates when calling both the RNG directly and via numpy.

Ending the model run¶

Models will continue to run until the end of their timeline is reached, unless explicitly told otherwise (see next section).

Open-ended timelines and Conditional Halting¶

In some models, rather than (or as well as) evolving the population over a fixed timeline, it may make more sense to iterate timesteps until some condition is met. The "Schelling" example illustrates this - it runs until all agents are in a satisfied state. Currently, the inbuilt LinearTimeline and CalendarTimeline classes support both fixed and open-ended timelines. In other cases it may be useful to temporarily exit the model for later resumption.

The model's halt method can be used to stop the model run. In these situations, the step method should have some logic to (conditionally) call the halt method.

Overriding the halt method should not be necessary and is not recommended. The finalise method, if needed, must be called explicitly for models with open-ended timelines that have been halted.

Deadlocks¶

Blocking communications between processes will deadlock if, for instance, the receiving process has ended due to an error. This will cause the entire run to hang (and may impact your HPC bill). The option example, as described above, has a check for random stream synchronisation that looks like this:

    def check(self) -> bool:
        # check the rng streams are still in sync by sampling from each one,
        # comparing, and broadcasting the result. If one process fails the
        # check and exits without notifying the others, deadlocks can result.
        # send the state representation to process 0 (others will get None)
        states = neworder.mpi.COMM.gather(self.mc.state(), 0)
        # process 0 checks the values
        ok = all(s == states[0] for s in states) if states else True
        # broadcast process 0's ok to all processes
        ok = neworder.mpi.COMM.bcast(ok, root=0)
        return ok

The key here is that there is only one result, shared between all processes. In this case only one process is performing the check and broadcasting the result to the others.

Time Comparison¶

neworder uses 64-bit floating-point numbers to represent time, and the values -inf, +inf and nan respectively to represent the concepts of the distant past, the far future and never. This allows users to define, or compare against, values that are:

• unequal to any time value, including itself (neworder.time.NEVER),
• before any other (non-never) time value (neworder.time.DISTANT_PAST) , or
• after any other (non-never) time value (neworder.time.FAR_FUTURE)

To compare time values with "never", use the supplied function isnever():

import neworder
n = neworder.time.NEVER
neworder.log(n == n) # False!
neworder.log(neworder.time.isnever(n)) # True

Data Types¶

If an argument to a neworder method or function is not the correct type, it will fail immediately (as opposed to python, which will fail only if an invalid operation for the given type is attempted (a.k.a. "duck typing")). This applies to contained types (numpy's dtype) too. In the example below, the function is expecting an integer, and will complain if you pass it a floating-point argument:

>>> import neworder
>>> neworder.df.unique_index(3.0)
Traceback (most recent call last):
  File "<stdin>", line 1, in <module>
TypeError: unique_index(): incompatible function arguments. The following argument types are supported:
    1. (n: int) -> numpy.ndarray[int64]

Invoked with: 3.0

Project Structure¶

Although obvious to many users, in order to promote reusability, it is recommended to separate out functionality into logical units, for example:

• model definition - the actual model implementation
• model data - loading and preprocessing of input data
• model execution - defining the parameters of the model and running it
• result postprocessing and visualisation

This makes life much easier when you want to:

• use the same model with different parameters and/or input data,
• run the model on different plaforms without modification (think desktop vs HPC cluster vs web service).
• have visualisations tailored to the platform you are working on.
• run multiple models from one script.

The examples use canned (i.e. already preprocessed) data but otherwise largely adhere to this pattern.