Dear Chris and Dawbarton, your post are very helpful, thank you. I use Python 3+ on daily basis, but would like to transfer my code to Julia. Very recently I tried to implement a short code for one of my student in Julia dealing with the fiiting of experimental data based on simple differential equation, namely dc/dt = -k1*c^k2. Please have look on the following code.

How can I start with parametrized ODE´s and data fitting in Julia?! Any suggestions will be very appreciated

I had great troubles installing (still not working) DifferentialEquations on Windows 10 64-bit (JuliaPro 6.0.1 with mingw -32)

Thank you,

Pavel Dytrych

import numpy as np

import pandas as pd

from scipy import integrate, interpolate

from scipy import optimize

import matplotlib.pylab as plt

import matplotlib

matplotlib.rcParams[‘text.usetex’] = True

matplotlib.rcParams[‘text.latex.unicode’] = True

“”" Old - data given by hand

x_data = np.array([

0,

60,

120,

240,

360,

600,

1080,

2040,

3960,

7800

])

y_data = np.array([

970,

739,

699,

654,

609,

556,

482,

393,

254,

52

])

"""

file = 'data_kinetics.xlsx’

df = pd.read_excel(file, sheetname=‘List1’)

print(df)

print(“The list of row indicies”)

print(df.index)

print(“The column headings”)

print(df.columns)

x_data = df[‘t’]

y_data = df[‘y’]

x, y = x_data, y_data

def f_pseudo(y, t, k):

dy0 = -k[0]*y[0]**k[1]

return dy0

def f_pseudo2(y, t, k):

dy0 = -k[0]*(y[0]-k[1])**2

return dy0

def my_ls_func(x,teta):

""“definition of function for LS fit

x gives evaluation points,

teta is an array of parameters to be varied for fit”""

# create an alias to f which passes the optional params

f2 = lambda y,t: f_pseudo(y, t, teta)

# calculate ode solution, retuen values for each entry of "x"

r = integrate.odeint(f2,y0,x)

#in this case, we only need one of the dependent variable values

return r[:,0]

def f_resid§:

""" function to pass to optimize.leastsq

The routine will square and sum the values returned by

this function"""

return y_data-my_ls_func(x_data,p)

#solve the system - the solution is in variable c

guess = [0,1] #initial guess for params

y0 = [970] #inital conditions for ODEs

(c,kvg) = optimize.leastsq(f_resid,guess) #get params

popt, pcov = optimize.leastsq(f_resid,guess)

# Parameters for pseudo 1. order

print (c[0], c[1])

pf11=c[0]

pf12=c[1]

def my_ls_func2(x,teta):

""“definition of function for LS fit

x gives evaluation points,

teta is an array of parameters to be varied for fit”""

# create an alias to f which passes the optional params

f3 = lambda y,t: f_pseudo2(y, t, teta)

# calculate ode solution, retuen values for each entry of "x"

r = integrate.odeint(f3,y02,x)

#in this case, we only need one of the dependent variable values

return r[:,0]

def f_resid2§:

""" function to pass to optimize.leastsq

The routine will square and sum the values returned by

this function"""

return y_data-my_ls_func2(x_data,p)

#solve the system - the solution is in variable c

guess2 = [0,1] #initial guess for params

y02 = [970] #inital conditions for ODEs

(c2,kvg) = optimize.leastsq(f_resid2,guess2) #get params

popt2, pcov2 = optimize.leastsq(f_resid2,guess2)

# Parameters for pseudo 2. order

print (c2[0], c2[1])

ps21=c2[0]

ps22=c2[1]

xeval=np.linspace(min(x_data), max(x_data),500)

gls = interpolate.UnivariateSpline(xeval, my_ls_func(xeval,c), k=5, s=0)

gls2 = interpolate.UnivariateSpline(xeval, my_ls_func2(xeval,c2), k=5, s=0)

f_experimental = open(“Data_experimental.csv”, “w”)

for i in range(len(x_data)):

f_experimental.write("{}, {} \n".format(x_data[i], y_data[i]))

f_experimental.close()

f_data = open(“Data_fitting.csv”, “w”)

for i in range(len(xeval)):

f_data.write("{}, {}, {} \n".format(xeval[i], my_ls_func(xeval,c)[i], my_ls_func2(xeval,c2)[i]))

f_data.close()

f_L = open(“Fitting_parameters_pseudo_first_order.csv”, “w”)

f_L.write("{}, {} \n".format(pf11, pf12))

f_L.close()

f_F = open(“Fitting_parameters_pseudo_second_order.csv”, “w”)

f_F.write("{}, {} \n".format(ps21, ps22))

f_F.close()

plt.figure(1, figsize=(8, 6))

plt.xlabel(r’t’, fontsize=18)

plt.ylabel(r’q_t’, fontsize=18)

plt.title(r"Experimental data adsorption",fontsize=18)

plt.xticks(fontsize=18)

plt.yticks(fontsize=18)

plt.plot(x_data,y_data,’*’, label=‘Experimental data’)

plt.plot(xeval, gls(xeval), ‘b-’, label=‘Pseudo-first order’)

plt.plot(xeval, gls2(xeval), ‘r-’, label=‘Pseudo-second order’)

legend = plt.legend(loc=‘best’, shadow=True, fontsize=‘x-large’)

legend.get_frame().set_facecolor(’#FFFFFF’)

plt.savefig(“Rozklad_silanu.png”, dpi=300)

plt.show()