📄 nlequations.java
字号:
dblCoef[i]=dblCoef[i]/q;
q=q*p;
}
x=0.0001;
x1=x;
y=0.2;
y1=y;
dx=1.0;
g=1.0e+37;
while (true)
{
u=dblCoef[0]; v=0.0;
for (i=1; i<=k; i++)
{
p=u*x1;
q=v*y1;
pq=(u+v)*(x1+y1);
u=p-q+dblCoef[i];
v=pq-p-q;
}
g1=u*u+v*v;
if (g1>=g)
{
if (is!=0)
{
it=1;
Real xtmp = new Real(x);
Real ytmp = new Real(y);
Real x1tmp = new Real(x1);
Real y1tmp = new Real(y1);
Real dxtmp = new Real(dx);
Real dytmp = new Real(dy);
Real ddtmp = new Real(dd);
Real dctmp = new Real(dc);
Real ctmp = new Real(c);
Real ktmp = new Real(k);
Real istmp = new Real(is);
Real ittmp = new Real(it);
g65(xtmp, ytmp, x1tmp, y1tmp, dxtmp, dytmp, ddtmp, dctmp, ctmp, ktmp, istmp, ittmp);
x=xtmp.doubleValue();
y=ytmp.doubleValue();
x1=x1tmp.doubleValue();
y1=y1tmp.doubleValue();
dx=dxtmp.doubleValue();
dy=dytmp.doubleValue();
dd=ddtmp.doubleValue();
dc=dctmp.doubleValue();
c=ctmp.doubleValue();
k=ktmp.intValue();
is=istmp.intValue();
it=ittmp.intValue();
if (it==0)
continue;
}
else
{
Real ttmp = new Real(t);
Real xtmp = new Real(x);
Real ytmp = new Real(y);
Real x1tmp = new Real(x1);
Real y1tmp = new Real(y1);
Real dxtmp = new Real(dx);
Real dytmp = new Real(dy);
Real ptmp = new Real(p);
Real qtmp = new Real(q);
Real ktmp = new Real(k);
Real ittmp = new Real(it);
g60(ttmp,xtmp,ytmp,x1tmp,y1tmp,dxtmp,dytmp,ptmp,qtmp,ktmp,ittmp);
t=ttmp.doubleValue();
x=xtmp.doubleValue();
y=ytmp.doubleValue();
x1=x1tmp.doubleValue();
y1=y1tmp.doubleValue();
dx=dxtmp.doubleValue();
dy=dytmp.doubleValue();
p=ptmp.doubleValue();
q=qtmp.doubleValue();
k=ktmp.intValue();
it=ittmp.intValue();
if (t>=1.0e-03)
continue;
if (g>1.0e-18)
{
it=0;
Real xtmp1 = new Real(x);
Real ytmp1 = new Real(y);
Real x1tmp1 = new Real(x1);
Real y1tmp1 = new Real(y1);
Real dxtmp1 = new Real(dx);
Real dytmp1 = new Real(dy);
Real ddtmp1 = new Real(dd);
Real dctmp1 = new Real(dc);
Real ctmp1 = new Real(c);
Real ktmp1 = new Real(k);
Real istmp1 = new Real(is);
Real ittmp1 = new Real(it);
g65(xtmp1, ytmp1, x1tmp1, y1tmp1, dxtmp1, dytmp1, ddtmp1, dctmp1, ctmp1, ktmp1, istmp1, ittmp1);
x=xtmp1.doubleValue();
y=ytmp1.doubleValue();
x1=x1tmp1.doubleValue();
y1=y1tmp1.doubleValue();
dx=dxtmp1.doubleValue();
dy=dytmp1.doubleValue();
dd=ddtmp1.doubleValue();
dc=dctmp1.doubleValue();
c=ctmp1.doubleValue();
k=ktmp1.intValue();
is=istmp1.intValue();
it=ittmp1.intValue();
if (it==0)
continue;
}
}
Real xtmp = new Real(x);
Real ytmp = new Real(y);
Real ptmp = new Real(p);
Real qtmp = new Real(q);
Real wtmp = new Real(w);
Real ktmp = new Real(k);
g90(xr,xi,dblCoef,xtmp,ytmp,ptmp,qtmp,wtmp,ktmp);
x=xtmp.doubleValue();
y=ytmp.doubleValue();
p=ptmp.doubleValue();
q=qtmp.doubleValue();
w=wtmp.doubleValue();
k=ktmp.intValue();
break;
}
else
{
g=g1;
x=x1;
y=y1;
is=0;
if (g<=1.0e-22)
{
Real xtmp = new Real(x);
Real ytmp = new Real(y);
Real ptmp = new Real(p);
Real qtmp = new Real(q);
Real wtmp = new Real(w);
Real ktmp = new Real(k);
g90(xr,xi,dblCoef,xtmp,ytmp,ptmp,qtmp,wtmp,ktmp);
x=xtmp.doubleValue();
y=ytmp.doubleValue();
p=ptmp.doubleValue();
q=qtmp.doubleValue();
w=wtmp.doubleValue();
k=ktmp.intValue();
}
else
{
u1=k*dblCoef[0];
v1=0.0;
for (i=2; i<=k; i++)
{
p=u1*x;
q=v1*y;
pq=(u1+v1)*(x+y);
u1=p-q+(k-i+1)*dblCoef[i-1];
v1=pq-p-q;
}
p=u1*u1+v1*v1;
if (p<=1.0e-20)
{
it=0;
Real xtmp = new Real(x);
Real ytmp = new Real(y);
Real x1tmp = new Real(x1);
Real y1tmp = new Real(y1);
Real dxtmp = new Real(dx);
Real dytmp = new Real(dy);
Real ddtmp = new Real(dd);
Real dctmp = new Real(dc);
Real ctmp = new Real(c);
Real ktmp = new Real(k);
Real istmp = new Real(is);
Real ittmp = new Real(it);
g65(xtmp, ytmp, x1tmp, y1tmp, dxtmp, dytmp, ddtmp, dctmp, ctmp, ktmp, istmp, ittmp);
x=xtmp.doubleValue();
y=ytmp.doubleValue();
x1=x1tmp.doubleValue();
y1=y1tmp.doubleValue();
dx=dxtmp.doubleValue();
dy=dytmp.doubleValue();
dd=ddtmp.doubleValue();
dc=dctmp.doubleValue();
c=ctmp.doubleValue();
k=ktmp.intValue();
is=istmp.intValue();
it=ittmp.intValue();
if (it==0)
continue;
Real xtmp1 = new Real(x);
Real ytmp1 = new Real(y);
Real ptmp1 = new Real(p);
Real qtmp1 = new Real(q);
Real wtmp1 = new Real(w);
Real ktmp1 = new Real(k);
g90(xr,xi,dblCoef,xtmp1,ytmp1,ptmp1,qtmp1,wtmp1,ktmp1);
x=xtmp1.doubleValue();
y=ytmp1.doubleValue();
p=ptmp1.doubleValue();
q=qtmp1.doubleValue();
w=wtmp1.doubleValue();
k=ktmp1.intValue();
}
else
{
dx=(u*u1+v*v1)/p;
dy=(u1*v-v1*u)/p;
t=1.0+4.0/k;
Real ttmp = new Real(t);
Real xtmp = new Real(x);
Real ytmp = new Real(y);
Real x1tmp = new Real(x1);
Real y1tmp = new Real(y1);
Real dxtmp = new Real(dx);
Real dytmp = new Real(dy);
Real ptmp = new Real(p);
Real qtmp = new Real(q);
Real ktmp = new Real(k);
Real ittmp = new Real(it);
g60(ttmp,xtmp,ytmp,x1tmp,y1tmp,dxtmp,dytmp,ptmp,qtmp,ktmp,ittmp);
t=ttmp.doubleValue();
x=xtmp.doubleValue();
y=ytmp.doubleValue();
x1=x1tmp.doubleValue();
y1=y1tmp.doubleValue();
dx=dxtmp.doubleValue();
dy=dytmp.doubleValue();
p=ptmp.doubleValue();
q=qtmp.doubleValue();
k=ktmp.intValue();
it=ittmp.intValue();
if (t>=1.0e-03)
continue;
if (g>1.0e-18)
{
it=0;
Real xtmp1 = new Real(x);
Real ytmp1 = new Real(y);
Real x1tmp1 = new Real(x1);
Real y1tmp1 = new Real(y1);
Real dxtmp1 = new Real(dx);
Real dytmp1 = new Real(dy);
Real ddtmp1 = new Real(dd);
Real dctmp1 = new Real(dc);
Real ctmp1 = new Real(c);
Real ktmp1 = new Real(k);
Real istmp1 = new Real(is);
Real ittmp1 = new Real(it);
g65(xtmp1, ytmp1, x1tmp1, y1tmp1, dxtmp1, dytmp1, ddtmp1, dctmp1, ctmp1, ktmp1, istmp1, ittmp1);
x=xtmp1.doubleValue();
y=ytmp1.doubleValue();
x1=x1tmp1.doubleValue();
y1=y1tmp1.doubleValue();
dx=dxtmp1.doubleValue();
dy=dytmp1.doubleValue();
dd=ddtmp1.doubleValue();
dc=dctmp1.doubleValue();
c=ctmp1.doubleValue();
k=ktmp1.intValue();
is=istmp1.intValue();
it=ittmp1.intValue();
if (it==0)
continue;
}
Real xtmp1 = new Real(x);
Real ytmp1 = new Real(y);
Real ptmp1 = new Real(p);
Real qtmp1 = new Real(q);
Real wtmp1 = new Real(w);
Real ktmp1 = new Real(k);
g90(xr,xi,dblCoef,xtmp1,ytmp1,ptmp1,qtmp1,wtmp1,ktmp1);
x=xtmp1.doubleValue();
y=ytmp1.doubleValue();
p=ptmp1.doubleValue();
q=qtmp1.doubleValue();
w=wtmp1.doubleValue();
k=ktmp1.intValue();
}
}
break;
}
}
if (k==1)
jt=0;
else
jt=1;
}
return true;
}
/**
* 内部函数
*/
private void g60(Real t, Real x, Real y, Real x1, Real y1,
Real dx, Real dy, Real p, Real q, Real k, Real it)
{
it.setValue(1);
while (it.intValue()==1)
{
t.setValue(t.doubleValue()/1.67);
it.setValue(0);
x1.setValue(x.doubleValue()-t.doubleValue()*dx.doubleValue());
y1.setValue(y.doubleValue()-t.doubleValue()*dy.doubleValue());
if (k.intValue()>=50)
{
p.setValue(Math.sqrt(x1.doubleValue()*x1.doubleValue()+y1.doubleValue()*y1.doubleValue()));
q.setValue(Math.exp(85.0/k.doubleValue()));
if (p.doubleValue()>=q.doubleValue())
it.setValue(1);
}
}
}
/**
* 内部函数
*/
private void g90(double[] xr, double[] xi, double[] dblCoef,
Real x, Real y, Real p, Real q, Real w, Real k)
{
int i;
if (Math.abs(y.doubleValue())<=1.0e-06)
{
p.setValue(-x.doubleValue());
y.setValue(0.0);
q.setValue(0.0);
}
else
{
p.setValue(-2.0*x.doubleValue());
q.setValue(x.doubleValue()*x.doubleValue() + y.doubleValue()*y.doubleValue());
xr[k.intValue()-1]=x.doubleValue()*w.doubleValue();
xi[k.intValue()-1]=-y.doubleValue()*w.doubleValue();
k.setValue(k.intValue()-1);
}
for (i=1; i<=k.intValue(); i++)
{
dblCoef[i]=dblCoef[i]-dblCoef[i-1]*p.doubleValue();
dblCoef[i+1]=dblCoef[i+1]-dblCoef[i-1]*q.doubleValue();
}
xr[k.intValue()-1]=x.doubleValue()*w.doubleValue();
xi[k.intValue()-1]=y.doubleValue()*w.doubleValue();
k.setValue(k.intValue()-1);
if (k.intValue()==1)
{
xr[0]=-dblCoef[1]*w.doubleValue()/dblCoef[0];
xi[0]=0.0;
}
}
/**
* 内部函数
*/
private void g65(Real x,Real y,Real x1,Real y1,Real dx,Real dy,Real dd,Real dc,Real c,
Real k,Real is,Real it)
{
if (it.intValue()==0)
{
is.setValue(1);
dd.setValue(Math.sqrt(dx.doubleValue()*dx.doubleValue()+dy.doubleValue()*dy.doubleValue()));
if (dd.doubleValue()>1.0)
dd.setValue(1.0);
dc.setValue(6.28/(4.5*k.doubleValue()));
c.setValue(0.0);
}
while(true)
{
c.setValue(c.doubleValue()+dc.doubleValue());
dx.setValue(dd.doubleValue()*Math.cos(c.doubleValue()));
dy.setValue(dd.doubleValue()*Math.sin(c.doubleValue()));
x1.setValue(x.doubleValue()+dx.doubleValue());
y1.setValue(y.doubleValue()+dy.doubleValue());
if (c.doubleValue()<=6.29)
{
it.setValue(0);
return;
}
dd.setValue(dd.doubleValue()/1.67);
if (dd.doubleValue()<=1.0e-07)
{
it.setValue(1);
return;
}
c.setValue(0.0);
}
}
/**
* 求复系数代数方程全部根的牛顿下山法
*
* @param n - 多项式方程的次数
* @param ar - 一维数组，长度为n+1，按降幂次序依次存放n次多项式方程的
* n+1个系数的实部
* @param ai - 一维数组，长度为n+1，按降幂次序依次存放n次多项式方程的
* n+1个系数的虚部
* @param xr - 一维数组，长度为n，返回n个根的实部
* @param xi - 一维数组，长度为n，返回n个根的虚部
* @return boolean 型，求解是否成功
*/
public boolean getRootNewtonDownHill(int n, double[] ar, double[] ai, double[] xr, double[] xi)
{
int m=0,i=0,jt=0,k=0,is=0,it=0;
double t=0,x=0,y=0,x1=0,y1=0,dx=0,dy=0,p=0,q=0,w=0,dd=0,dc=0,c=0;
double g=0,u=0,v=0,pq=0,g1=0,u1=0,v1=0;
// 初始判断
m=n;
p=Math.sqrt(ar[m]*ar[m]+ai[m]*ai[m]);
while ((m>0)&&(p+1.0==1.0))
{
m=m-1;
p=Math.sqrt(ar[m]*ar[m]+ai[m]*ai[m]);
}
// 求解失败
if (m<=0)
return false;
for (i=0; i<=m; i++)
{
ar[i]=ar[i]/p;
ai[i]=ai[i]/p;
}
for (i=0; i<=m/2; i++)
{
w=ar[i];
ar[i]=ar[m-i];
ar[m-i]=w;
w=ai[i];
ai[i]=ai[m-i];
ai[m-i]=w;
}
// 迭代求解
k=m;
is=0;
w=1.0;
jt=1;
while (jt==1)
{
⌨️ 快捷键说明
复制代码
Ctrl + C
搜索代码
Ctrl + F
全屏模式
F11
切换主题
Ctrl + Shift + D
显示快捷键
?
增大字号
Ctrl + =
减小字号
Ctrl + -