Integrales elípticas sn, cn, dn
Posición x de la partícula
x=A·cn(
A
l
0
k
m
t
)
public class Eliptica {
final double CA=0.0003; // The accuracy is the square of CA.
public double sn;
public double cn;
public double dn;
public Eliptica() {
sn=dn=cn=0.0;
}
void sncndn(double uu, double emmc){
//Returns the Jacobian elliptic functions sn(u; kc), //cn(u; kc), and dn(u; kc). Here uu = u, while
//emmc = k2 c .
double a,b,c=0.0,d=0.0,emc,u;
double[] em=new double[14];
double[] en=new double[14];
int i,ii,l=0;
boolean bo;
emc=emmc;
u=uu;
if (emc!=0.0) {
bo=(emc < 0.0);
if (bo) {
d=1.0-emc;
emc /= -1.0/d;
u *= (d=Math.sqrt(d));
}
a=1.0;
dn=1.0;
for (i=1;i<=13;i++) {
l=i;
em[i]=a;
en[i]=(emc=Math.sqrt(emc));
c=0.5*(a+emc);
if (Math.abs(a-emc) <= CA*a) break;
emc *= a;
a=c;
}
u *= c;
sn=Math.sin(u);
cn=Math.cos(u);
if (sn!=0.0) {
a=cn/sn;
c *= a;
for (ii=l;ii>=1;ii--) {
b=em[ii];
a *= c;
c *= dn;
dn=(en[ii]+a)/(b+a);
a=c/b;
}
a=1.0/Math.sqrt(c*c+1.0);
sn=(sn >= 0.0 ? a : -a);
cn=c*sn;
}
if (bo) {
a=dn;
dn=cn;
cn=a;
sn /= d;
}
} else {
cn=1.0/cosh(u);
dn=cn;
sn=tanh(u);
}
}
double cosh(double x){
return((Math.exp(x)+Math.exp(-x))/2);
}
double tanh(double x){
return((Math.exp(x)-Math.exp(-x))/(Math.exp(x)+Math.exp(-x)));
}
}
//objeto
Eliptica obj=new Eliptica();
//posición del oscilador
obj.sncndn((x0*Math.sqrt(cte)*t/lonMuelle),0.5);
x=x0*obj.cn;
t+=dt;
|
Referencias
Press W. H., Teukolsky S. A., Vetterling W. T., Flannery B.
P. Numerical Recipes in C, Second edition, Special functions. Elliptic
integrals and Jacobian Elliptic functions, Chapter 6º. pp. 269-270. Cambridge
University Press. Código en C adaptado por el
autor al lenguaje Java