neurophysik / jitcdde

Just-in-time compilation for delay differential equations
Other
57 stars 14 forks source link

How to pass variable parameter to function #33

Closed Ziaeemehr closed 3 years ago

Ziaeemehr commented 3 years ago

It's Kuramoto model, I defined coupling g=Symbole("g") as parameter. I wand to pass the coupling to kuramotos_f(). if else statement at main function just check for presence of compiled file.

one more question is that is it possible to pass the adjacency matrix also as a parameter?

import os
import numpy as np
import pylab as plt
from jitcode import jitcode, y
from numpy.random import choice
from numpy.random import uniform
from symengine import sin, Symbol

def kuramotos_f():  # kuramotos_f(g, A)??
    for i in range(n):
        coupling_sum = sum(sin(y(j)-y(i))
                           for j in range(n)
                           if A[i, j])
        yield omega[i] + g / (n-1)*coupling_sum

def main(coupling, time_simulation):
        cfilename = "data/jitced.so"
        if not os.path.exists(cfilename):

            I = jitcode(kuramotos_f, n=n, control_pars=[g])
            I.set_parameters(coupling)
            I.generate_f_C()
            I.compile_C()
            I.save_compiled(overwrite=True, destination=cfilename)
        else:
            I = jitcode(n=n, module_location=cfilename)
            I.set_parameters(coupling)

        I.set_integrator("dopri5", atol=1e-6, rtol=1e-5)
        I.set_initial_value(initial_state, time=0.0)

        times = range(0, time_simulation)
        phases = np.empty((len(times), n))

        for i in range(len(times)):
            phases[i, :] = I.integrate(times[i]) % (2*np.pi)

        return times, phases

n = 100
q = 0.2
g = Symbol("g")
time_simulation = 2001

if __name__ == "__main__":

    np.random.seed(1)

    A = choice([1, 0], size=(n, n), p=[q, 1-q])
    initial_state = uniform(0, 2*np.pi, n)
    omega = uniform(-0.5, 0.5, n)
    omega.sort()

    coupling = 3.0
    times, phases = main(coupling, time_simulation)

Thank you for any guide.

Ziaeemehr commented 3 years ago

Oh, I should have put this in Jitcode.