間違い？jacobiを用いたnewton法
以下のプログラムを作ったんですが実行したときに答えがnanになります。どなたか間違いがわかる方よろしくおねがします。
const double e =1e-10;
enum {N=100};
const double x0 =1;
const double y0 =1;
double f(const double x,double y) {return x*x-y*y;}
double g(const double x,double y) {return 6*x*y-5;}
double dfdx(const double x,double y){return 2*x;}
double dfdy(const double x,double y){return -2*y;}
double dgdx(const double x,double y){return 6*y;}
double dgdy(const double x,double y){return 6*x;}
double det(double x,double y);
main()
{
double x[N+1],y[N+1];
int i;
double det2;
x[0]=x0;
y[0]=y0;
for(i=0;i<=20;i++)
{
det2=det(x[i],y[i]);
x[i+1]=x[i]-(dgdy(x[i],y[i])*f(x[i],y[i]))-(dfdy(x[i],y[i])*g(x[i],y[i]))/det2;
y[i+1]=y[i]-(dgdx(x[i],y[i])*f(x[i],y[i]))-(dfdx(x[i],y[i])*g(x[i],y[i]))/det2;
printf("x = %lf\n",x[i+1]);
printf("y = %lf\n",y[i+1]);
printf("\n");
}
exit(0);
}
double det(double a,double b)
{
double ans;
ans=1/(dfdx(a,b)*dgdy(a,b)-dfdy(a,b)*dgdx(a,b));
return ans;
}
補足
日本ではないのですか？よくわかりましたね・・・どうやって？