1

I am trying to do a global fit with shared parameters using objective functions that contain exponential integrals. If I were using, curve_fit, I would simply include the line:

from scipy.special import expi

I haven't been able to get the equivalent symbolic function in Symfit, 'symfit.Ei' to work.

import numpy as np
import scipy.optimize
from scipy.special import expi
import numpy as np
import symfit as sf

y1 = np.array(data1)

y3 = np.array(ydata3)

x1 = np.array(xdata1)

x3 = np.array(xdata3)

t1 = 50     # AOM 1st OFF Time
t2 = 100    # AOM 2nd ON Time
t3 = 150    # AOM 2nd OFF Time 
t4 = 450


def mod1(data, a_decay, b_decay, constant, on1_a_amplitude, on1_b_amplitude, ARE_A, ARE_B): 
    first_term_a =  sf.Ei(-data/a_decay) 
    second_term_a = sf.Ei(-data*(1+ARE_A*a_decay)/a_decay)) 
    first_term_b =  sf.Ei(-data/b_decay)) 
    second_term_b = (sf.Ei(-data*(1+ARE_B*b_decay)/b_decay))
 return constant+on1_a_amplitude*(first_term_a-second_term_a+sf.log(1+ARE_A*a_decay)) + on1_b_amplitude*(first_term_b-second_term_b+sf.log(1+ARE_B*b_decay))

def mod2(data, a_decay, b_decay, constant, on2_a_amplitude, on2_b_amplitude, ARE_A, ARE_B): # not all parameters are used here
    first_term_a =  (sf.Ei((t2-t1-data)/a_decay)) - (sf.Ei((t2-data)/a_decay))
    second_term_a = (sf.Ei((t2-data)*(1+ARE_A*a_decay)/a_decay)) - (sf.Ei((t2-t1-data)*(1+ARE_A*a_decay)/a_decay))
    third_term_a = (sf.Ei((t2-data)/a_decay)) -(sf.Ei((t2-data)*(1+ARE_A*a_decay)/a_decay)) + sf.log(ARE_A*a_decay)

    first_term_b =  (sf.Ei((t2-t1-data)/b_decay)) - (sf.Ei((t2-data)/b_decay))
    second_term_b = (sf.Ei((t2-data)*(1+ARE_B*b_decay)/b_decay)) - (sf.Ei((t2-t1-data)*(1+ARE_B*b_decay)/a_decay))
    third_term_b = (sf.Ei((t2-data)/b_decay)) - (sf.Ei((t2-data)*(1+ARE_B*b_decay)/b_decay)) + sf.log(ARE_B*b_decay)

return constant+on2_a_amplitude*(sf.exp((t1-t2)/a_decay)*(first_term_a+second_term_a)+third_term_a)+on2_b_amplitude*(sf.exp((t1-t2)/b_decay)*(first_term_b+second_term_b)+third_term_b)  

x_1, x_3, y_1, y_3 = sf.variables('x_1, x_3, y_1, y_3')

a_dec = sf.Parameter(value=64.32, min=0.0, max=100)
b_dec = sf.Parameter(value=53.88, min=0.0, max=100)
con = sf.Parameter(value=2.911E-4, min=0.0, max=1)
off1_a_amp = sf.Parameter(value=0.015, min=0.0, max=1)
off1_b_amp = sf.Parameter(value=0.015, min=0.0, max=1)
off2_a_amp = sf.Parameter(value=0.015, min=0.0, max=1)
off2_b_amp = sf.Parameter(value=0.015, min=0.0, max=1)
AR_A = sf.Parameter(value=0.032, min=0.0, max=1)
AR_B = sf.Parameter(value=0.1, min=0.0, max=1)


model = sf.Model({
y_1: mod1(x_1, a_dec, b_dec, con, off1_a_amp, off1_b_amp, AR_A, AR_B),
y_3: mod2(x_3, a_dec, b_dec, con, off2_a_amp, off2_b_amp, AR_A, AR_B ),
})


fit = sf.Fit(model, x_1=x1, x_3=x3, y_1=y1, y_3=y3)
fit_result = fit.execute()
print(fit_result)

For some reason 'sf.Ei' throws a name error. i.e. 'Name Error: name 'Ei' is not defined.

  • shouldn't it be expi()? Believe the compiler - the name is not defined. – duffymo May 28 '23 at 23:06
  • Thanks for the quick reply @duffymo. I tried your suggestion, using 'expi()' and got the following error: 'TypeError: ufunc 'expi' not supported for the input types, and the inputs could not be safely coerced to any supported types according to the casting rule ''safe''' – luiz_iniciante May 29 '23 at 00:13
  • My guess is that you should find the documentation and find out what types that function expects. You can't make progress on this problem until you have the function you want and the correct input. – duffymo May 29 '23 at 16:25
  • I suspect the problem is related to related to the function not being defined in with [Lambdify](https://stackoverflow.com/questions/45273827/some-function-is-not-defined-with-sympy-lambdify). This post here says that calling Model in symfit [means the symbolic expression is converted into a lambda function](https://stackoverflow.com/questions/54134273/attributeerror-when-evaluating-symfit-model-with-substraction-in-exponent). I can't see all the way to a solution, though. – luiz_iniciante May 29 '23 at 19:18

0 Answers0