Estoy intentando pasar éste código de Matlab a Python, es la solución numérica de la integral elíptica completa generalizada de tercer tipo, la cual es la siguiente:
Éste es el código en Matlab que tengo:
function res = integral_eliptica(kc,p,c,s)
if kc==0
res=NaN;
quit;
end
errtol=0.000001;
k=abs(kc);
pp=p;
cc=c;
ss=s;
em=1;
if p>0
pp=sqrt(p);
ss=s/pp;
else
f=kc^2;
q=1-f;
g=1-pp;
f=f-pp;
q=q*(ss-c*pp);
pp=sqrt(f/g);
cc=(c-ss)/g;
ss=-q/(g^2*pp)+cc*pp;
end
f=cc;
cc=cc+ss/pp;
g=k/pp;
ss=2*(ss+f*g);
pp=g+pp;
g=em;
em=k+em;
kk=k;
while abs(g-k)>g*errtol
k=2*sqrt(kk);
kk=k*em;
f=cc;
cc=cc+ss/pp;
g=kk/pp;
ss=2*(ss+f*g);
pp=g+pp;
g=em;
em=k+em;
end
res=(pi/2)*(ss+cc*em)/(em*(em+pp));
end
Y esto es lo que he intentado hacer en Python:
import numpy as np
from numpy import npymath
def cel(kc,p,c,s):
if kc==0:
return cel = NPY_NAN
errtol=0.000001
k=ABS(kc)
pp=p
cc=c
ss=s
em=1.0
if (p>0)
pp=(p)**(1/2)
ss=s/pp
else
f=kc*kc
q=1.0-f
g=1.0-pp
f=f-pp
q=q*(ss-c*pp)
pp=(f/g)**(1/2)
cc=(c-ss)/g
ss=-q/(g*g*pp)+cc*pp
f=cc
cc=cc+ss/pp
g=k/pp
ss=2*(ss+f*g)
pp=g+pp
g=em
em=k+em
kk=k
while (abs(g-k)>g*errtol)
k=2*(kk)**(1/2)
kk=k*em
f=cc
cc=cc+ss/pp
g=kk/pp
ss=2*(ss+f*g)
pp=g+pp
g=em
em=k+em
return cel=(pi/2)*(ss+cc*em)/(em*(em+pp))
Sin embargo, no logro hacer que funcione, no estoy muy familiarizado con la sintaxis de Python para lo que requiero. Otra opción que he probado para resolver la integral es usando su equivalente con R_F y R_J de scipy.special:
Pero me aparece el siguiente error:
AttributeError: module 'scipy.special' has no attribute 'elliprf'