Решение системы нелинейных уравнений. - C (СИ)

Формулировка задачи:

Нужно решить систему уравнений, методом ньютона. Решаю задачу на СИ, все по блок схеме. Пересел с паскаля буквально месяц. Возможно где то накосячил, помогите с ошибкой. При запуске программы пишет Floating Point error: Domain. Прикрепил блок схему и саму задачу. Неделю парюсь(

                    Листинг программы
                
 
                            
                        
#include <stdio.h>
#include <stdlib.h>
#include <math.h>
 
double fun1(double x1,double x2)
{
double f;
f=x1+3*(log(x1)/log(10))-x2*x2;
return f;
}
 
double fun2(double x1,double x2)
{
double f;
f=2*x1*2*x1-x1*x2-5*x1+1;
return f;
}
 
double fun1x1(double x1)
{
double f;
f=1+3/x1;
return f;
}
 
double fun1x2(double x2)
{
double f;
f=2*x2;
return f;
}
 
double fun2x1(double x1, double x2)
{
double f;
f=4*x1-x2-5;
return f;
}
 
double fun2x2(double x1)
{
double f;
f=-x1;
return f;
}
 
void main()
{
    double x1,x2,e,x10,x20,a[3][3],b[3],dx1,dx2;
        clrscr();
        //  printf ("vvedite x1 ");
        //scanf ("%f", &x1);
        //printf ("vvedite x2 ");
        //scanf ("%f", &x2);
        x1=3.4;
        x2=2.2;
        e=0.00001;
        do
        {
            x10=x1;
                x20=x2;
                b[1]=-fun1(x10,x20);
                b[2]=-fun2(x10,x20);
                a[1][1]=fun1x1(x10); a[1][2]=fun1x2(x20);
                a[2][1]=fun2x1(x10,x20); a[2][2]=fun2x2(x10);
                dx1=(b[1]*a[2][2]-a[1][2]*b[2])/(a[1][1]*a[2][2]-a[1][2]*a[2][1]);
                dx2=(a[1][1]*b[2]-b[1]*a[2][1])/(a[1][1]*a[2][2]-a[1][2]*a[2][1]);
                x1=x10+dx1;
                x2=x20+dx2;
        } while ((fabs(dx1)>=e) || (fabs(dx2)>=e));
        printf ("x1 = %f\n",x1);
        printf ("x2 = %f\n",x2);
        getch();
}

Решение задачи: «Решение системы нелинейных уравнений.»

Листинг программы

#include <cstdlib>
#include <iostream>
#include <math.h>
 
 
#if 0  // код в этом блоке компилятор не воспринимает, так как препроцессор его сьедает, можно здесь писать всякую чепуху не по делу и ничего не бояться
xy-tg(x-y)=0
0,5x^2+2y^2-1=0
 
x=0,975247
y=0,512078     
#endif
 
double fun1(double x,double y)
{
    double f;
    f=x*y-tan(x-y); //здесь твоя функция №1
    return f;
}
 
double fun2(double x,double y)
{
    double f;
    f=0.5*x*x+2*y*y-1; //здесь твоя функция №2
    return f;
}
 
double dfun1_dx(double x,double y)
{
    double f;
    f=(y-1)*x-(tan(x-y)*tan(x-y)+1); //здесь твоя производная от функции №1 по х
    return f;
}
 
double dfun1_dy(double x,double y)
{
    double f;
    f=(x-1)*y+(tan(x-y)*tan(x-y)+1);  //здесь твоя производная от функции №1 по y
    return f;
}
 
double dfun2_dx(double x, double y)
{
    double f;
    f=x;   //здесь твоя производная от функции №2 по х
    return f;
}
 
double dfun2_dy(double x,double y)
{
    double f;
    f=4*y;   //здесь твоя производная от функции №2 по y
    return f;
}
 
using namespace std;
 
int main()
{
    double x,y,e,x0,y0,a,b,c,d,lick,suck,dx,dy; // переменные, все те же что и у тебя только без массивов ( массивов в с++ нету как таковых, но можно сделать но это надо отдельно обьяснять)
    
    x=0.9; // (приближенное х)
    y=0.5; // (приближенное у)
    e=0.000001; //точность
 
    do
    {
        x0=x;
        y0=y;
        lick=-fun1(x0,y0);
        suck=-fun2(x0,y0);
        a=dfun1_dx(x0,y0);
        b=dfun1_dy(x0,y0);
        c=dfun2_dx(x0,y0);
        d=dfun2_dy(x0,y0);
        dx=(lick*d-b*suck)/(a*d-b*c);
        dy=(a*suck-lick*c)/(a*d-b*c);
        x=x0+dx;
        y=y0+dy;
    } while (((dx<e && dx>(-e)) || (dy<e && dy>(-e)))==false);  // в этом цикле все абсолютно в точности как у тебя, ничего менять не нужно (кстати в твоем коде ошибки, здесь исправлены)
 
    printf ("x = %f\n",x);
    printf ("y = %f\n",y);
    system("PAUSE");
    return 1;
}

Объяснение кода листинга программы

Код решает систему нелинейных уравнений методом итераций. Вот список действий, которые код выполняет по порядку:

Объявляет функции fun1 и fun2, которые представляют собой левую и правую части уравнений системы.
Объявляет функции dfun1_dx и dfun1_dy, которые представляют собой производные левой части уравнений системы по переменным x и y.
Объявляет функции dfun2_dx и dfun2_dy, которые представляют собой производные правой части уравнений системы по переменным x и y.
В функции main задаёт начальные приближения для переменных x и y (в данном случае 0.9 и 0.5 соответственно).
Задаёт точность e, которая определяет, когда следует остановить итерационный процесс.
Запускает итерационный процесс с помощью цикла do-while, в котором на каждой итерации обновляет значения переменных x и y с помощью формулы итерационного процесса.
Проверяет условие остановки: если изменение x и y на текущей итерации меньше e, то процесс продолжается.
После выхода из цикла выводит значения x и y на экран.
Вызывает функцию system(PAUSE), чтобы программа не закрылась сразу после вывода, а дождалась нажатия клавиши.
Возвращает 1, чтобы программа не закрылась сразу после вызова функции main.