In [ ]:

    

import numpy as np
import matplotlib.pyplot as plt
from pyodesys.symbolic import SymbolicSys
%matplotlib inline


    

In [ ]:

    

def f(t, y, p):
    k = t**p['e']
    return {
        'a': -k*y['a'],
        'b': +k*y['a']
    }

def analytic(tout, init_y, p):
    y0ref = init_y[0]*np.exp(-tout**(p[0]+1)/(p[0]+1))
    return np.array([y0ref, init_y[0] - y0ref + init_y[1]]).T


    

In [ ]:

    

def integrate_and_plot(odesys, axes=None):
    result = odesys.integrate(3, {'a': 2, 'b': 1}, {'e': 2}, integrator='cvode')
    if axes is None:
        _fig, axes = plt.subplots(1, 3, figsize=(14, 4))
    result.plot(ax=axes[0])
    result.plot(ax=axes[1], y=result.yout - analytic(result.xout, result.yout[0, :], result.params))
    result.plot_invariant_violations(ax=axes[2])
    print({k: v for k, v in result.info.items() if not k.startswith('internal')})


    

In [ ]:

    

odes = SymbolicSys.from_callback(f, names='ab', param_names='e', dep_by_name=True, par_by_name=True,
                                     linear_invariants=[{'a': 1, 'b': 1}])


    

In [ ]:

    

integrate_and_plot(odes)


    

In [ ]:

    

odes_auto = odes.as_autonomous()


    

In [ ]:

    

integrate_and_plot(odes_auto)