Understand observers with Lotka-Volterra equations

In this example, we present what can we do with observers in ponio, with Lotka-Voletrra equations. The equations are the following

\[\begin{split}\begin{cases} \frac{\mathrm{d}x}{\mathrm{d}t} = \alpha x - \beta xy\\ \frac{\mathrm{d}y}{\mathrm{d}t} = \delta xy - \gamma y\\ \end{cases}\end{split}\]

with \(\alpha = \frac{2}{3}\), \(\beta = \frac{4}{3}\) and \(\delta = \gamma = 1\). Initial condition \((x, y)(t=0) = (x_0, x_0)\), and \(x_0 = 1\).

We write the problem like

    // parameters
    double alpha = 2. / 3.;
    double beta  = 4. / 3.;
    double gamma = 1.;
    double delta = 1.;

    auto lotka_volterra_pb = [=]( double, auto const& u, state_t& du )
    {
        du[0] = alpha * u[0] - beta * u[0] * u[1];
        du[1] = delta * u[0] * u[1] - gamma * u[1];
    };

And launch the simulation between \(0\) and \(15\) with the time step \(\Delta t = 0.1\)

    ponio::time_span<double> const t_span = { 0., 15. }; // begin and end time
    double const dt                       = 0.1;         // time step

    state_t const u0 = { 1., 1. }; // initial condition

    ponio::solve( lotka_volterra_pb, ponio::runge_kutta::rk_33(), u0, t_span, dt, obs );

The next sections specify how to define the observer obs, to write into a file, in standard output, a stream, or a user-defined observer. The observer is called at each iteration in time and take the current time \(t^n\), the current state \(u^n\) and the current time step \(\Delta t^n\).

The file observer

In this example we export data into a file with ponio::observer::file_observer class. The ponio library provides this output for simple types and containers (thanks concepts). The output is formate as following

0 1 1  0.1000000000000000056
0.1000000000000000056 0.9356451028806584969 0.9967476543209876638  0.1000000000000000056
0.2000000000000000111 0.8761754465532282099 0.9873714302907439233  0.1000000000000000056
0.3000000000000000444 0.8218181271885944827 0.9725318841195829123  0.1000000000000000056
0.4000000000000000222 0.7726097837233855126 0.9529683892879063922  0.1000000000000000056

First column is the current time, the last one is the current time step, and between we get the multiple values of current state. Other observers provide by ponio have the same format.

The ponio::observer::file_observer class can be build with a std::string or a std::filesystem::path. The constructor makes necessary directory if needed.

    std::string filename = "lotka_volterra_fobs.txt";
    auto obs             = ponio::observer::file_observer( filename );

Note

If name is fixed at compile time you can use literal _fobs by

using namespace ponio::observer;

auto obs = "output.txt"_fobs;

Note

Data are not flushed into the output file, you have to wait the flush of the buffer.

The full example can be found in lotka_volterra_fobs.cpp.

Solution in phase space \((x, y)\)

The cout observer

In this example we export data into the standard output with ponio::observer::cout_observer class.

    auto obs = ponio::observer::cout_observer();

Note

Data are not flushed into the standard output, you have to wait the flush of the buffer.

The full example can be found in lotka_volterra_cobs.cpp.

Solution in phase space \((x, y)\)

The stream observer

In this example we export data into the standard output with ponio::observer::stream_observer class. This observer is the more generally one and can be plug with any other type of stream, for example a std::stringstream. The ponio::observer::file_observer is a specialization with a std::ofstream and ponio::observer::cout_observer with std::cout.

    std::stringstream buffer;
    auto obs = ponio::observer::stream_observer( buffer );

Next we get all informations into the observer or our buffer.

    std::cout << buffer.str() << "\n";

The full example can be found in lotka_volterra_sobs.cpp.

Solution in phase space \((x, y)\)

The user-defined observer

An observer is a invocable object (function, functor, lambda) with the following signature

void operator() ( value_t current_time, state_t & current_state, value_t current_time_step );

In our case, value_t is double and state_t is std::valarray<double>. We propose to export also the invariant of Lotka-Voletrra equations given by

\[V = \delta x - \ln(x) + \beta y - \alpha \ln(y)\]

We implement the observer into a class

struct LV_observer
{
    std::ofstream output;
    double a;
    double b;
    double d;
    double g;

    LV_observer( std::string&& filename, double alpha, double beta, double delta, double gamma )
        : output( filename )
        , a( alpha )
        , b( beta )
        , d( delta )
        , g( gamma )
    {
    }

    ~LV_observer()
    {
        output.close();
    }

    double
    V( std::valarray<double> const& u ) const
    {
        double x = u[0];
        double y = u[1];

        return d * x - g * std::log( x ) + b * y - a * std::log( y );
    }

    void
    operator()( double t, std::valarray<double>& u, double /* dt */ )
    {
        output << t << " " << u[0] << " " << u[1];
        output << " " << V( u ) << "\n";
    }
};

Now we just have to build an instance of this class

    auto obs = LV_observer( "lotka_volterra_uobs.txt", alpha, beta, delta, gamma );

The output of this observer looks like this

0 1 1 2.33333
0.1 0.935645 0.996748 2.33333
0.2 0.876175 0.987371 2.33333
0.3 0.821818 0.972532 2.33333
0.4 0.77261 0.952968 2.33333

The full example can be found in lotka_volterra_uobs.cpp.

Solution in phase space \((x, y)\)

We can compute from previous simulation the invariant \(V\), or display it with out user-defined observer. We don’t change the default precision in C++ which is set to 6 by default (see std::ios_base::precision), in ponio observers the precision is set to the maximum of precision of given type.

Relative error of invariant \(V\)