[Python] scipy.integrate.quad taking too long to excute when being used inside odeint

The coupled differential equations The equations have form factors denoted by Fxx, Fxy, Fyy,etc which are integrated with rms beam size as input as shown below in the image.coefficients of the coupled differential equations The terms Fxx,Fxy, Fyy, etc are as follows, Form factor solved using quad The goal of the program is to find the eigenvalues of the transfer matrix of the coupled differential equations mentioned above. The form factor function is defined as follows:

def Fxx(t,x,y,z):
x=x*1000
y=y*1000
z=z*1000
fxx=((x**2+t)**(-5/2))*((y**2+t)**(-1/2))*((gamma**2*z**2+t)**(-1/2))
return fxx*1e15


This form factor is called inside the function Maxexix below which calculates elements of one of the columns of transfer matrix.

def Maxexix(a,t,cc,xenv,yenv,zenv):
index=int(((t%lattlen)/lattlen)*steps)
rx=xenv[index]
ry=yenv[index]
rz=zenv[index]
fxz=quad(Fxz,0,np.inf,args=(rx,ry,rz),limit=100,epsabs=1e-5)
fxz=cc*fxz[0]
fxy=quad(Fxy,0,np.inf,args=(rx,ry,rz),limit=100,epsabs=1e-5)
fxy=cc*fxy[0]
fyz=quad(Fyz,0,np.inf,args=(rx,ry,rz),limit=100,epsabs=1e-5)
fyz=cc*fyz[0]
fxx=quad(Fxx,0,np.inf,args=(rx,ry,rz),limit=100,epsabs=1e-5)
fxx=cc*fxx[0]
fyy=quad(Fyy,0,np.inf,args=(rx,ry,rz),limit=100,epsabs=1e-5)
fyy=cc*fyy[0]
fzz=quad(Fzz,0,np.inf,args=(rx,ry,rz),limit=100,epsabs=1e-5)
fzz=cc*fzz[0]
spx=quad(Ix,0,np.inf,args=(rx,ry,rz),limit=100,epsabs=1e-5)
spx=cc*spx[0]
spy=quad(Iy,0,np.inf,args=(rx,ry,rz),limit=100,epsabs=1e-5)
spy=cc*spy[0]
spz=quad(Iz,0,np.inf,args=(rx,ry,rz),limit=100,epsabs=1e-5)
spz=cc*spz[0]

a12=(rx*ry*fxy)
a1=(kmag(t)+((3*e**2)/rx**4)-spx+(3*rx**2*fxx))
a13=(rx*rz*fxz)
a2=(-kmag(t)+((3*e**2)/ry**4)-spy+(3*ry**2*fyy))
a23=(ry*rz*fyz)
a3=(kzmag(t)+(3*(elong/gamma**2)**2/rz**4)-spz+(3*gamma**2*rz**2*fzz))

ax=a[1]
axp=-a1*a[0]-a12*a[2]-(gamma**2*a13)*a[4]
ex=a[3]
exp=-a12*a[0]-a2*a[2]-(gamma**2*a23)*a[4]
ix=a[5]
ixp=-a13*a[0]-a23*a[2]-a3*a[4]

return [ax,axp,ex,exp,ix,ixp]


This function is integrated using odeint as follows:

axsol=odeint(Maxexix,ax0,s,args=(k,rx,ry,rz))
aysol=odeint(Mayeyiy,ay0,s,args=(k,rx,ry,rz))
azsol=odeint(Mizxizyaz,az0,s,args=(k,rx,ry,rz))
bxsol=odeint(Mbxfxjx,bx0,s,args=(k,rx,ry,rz))
bysol=odeint(Mbyfyjy,by0,s,args=(k,rx,ry,rz))
bzsol=odeint(Mjzxjzybz,bz0,s,args=(k,rx,ry,rz))


Each of these solutions gives 6 elements of a single column of a 6X6 transfer matrix. This odeint is taking too long as it is stuck in calculating the formfactors. Number of steps for odeint is 500 i.e. is s.

The arguments in quad were of the order 1e-4, so to help quad converge I multiplied and divided rx,ry,rz by 1000 and took out the denominator out from each bracket of Fxx and substituting 1e6*dt= dT, I was left with 1e15 factor which was multiplied in the return command. This helped converge the integral but it is taking too long to solve a simple integral. 9*500 quad commands are taking more than 15mins=900sec.

Continue reading...