using System;
using System.Collections.Generic;
using System.Text;

namespace ZI_22_KMDI_6
{
    class RyngeRytta4
    {
        const int n = 2; // Порядок системи
        const double h = 0.1; // Крок
        double x, a = 0, b = 1; // Проміжок
        global double[] y = new double[n]; //если небольшой масив
        global double[] f = new double[n]; //тогда можна без global

        void dy() // Праві частини рівнянь
        {
            f[0] = y[1];
            f[1] = (8 * y[1] - 5 * y[0] - 4 * Math.Exp(x) * Math.Sin(x / 2)) / 4;
        }
        public void Calculate()
        {
            double[] K1 = new double[n];
            double[] K2 = new double[n];
            double[] K3 = new double[n];
            double[] K4 = new double[n];
            double[] z = new double[n];

            // Початкові значення
            y[0] = 1; y[1] = 1;
            x = a;

            do
            {
                dy();
                for (int i = 0; i < n; i++)
                {
                    K1[i] = f[i];
                    z[i] = y[i];
                    y[i] = z[i] + 0.5 * h * K1[i];
                }
                x = x + 0.5 * h;
                dy();
                for (int i = 0; i < n; i++)
                {
                    K2[i] = f[i];
                    y[i] = z[i] + 0.5 * h * K2[i];
                }
                dy();
                for (int i = 0; i < n; i++)
                {
                    K3[i] = f[i];
                    y[i] = z[i] + 0.5 * h * K3[i];
                }
                x = x + 0.5 * h;
                dy();
                for (int i = 0; i < n; i++)
                {
                    K4[i] = f[i];
                    y[i] = z[i] + h * (K1[i] + K4[i] + 2 * (K2[i] + K3[i])) / 6;
                }
                Console.WriteLine("x = " + x + "; y = " + y[0] + "; y' = " + y[1]);
            }
            while (x < b);
        }
    }
    class Program
    {
        static void Main(string[] args)
        {
            RyngeRytta4 Lab6 = new RyngeRytta4();
            Lab6.Calculate();
            Console.Read();
        }
    }
}

#include <stdio.h>
#include <math.h>
#include <conio.h>
#define n 2
void main(void)
{
clrscr();
int i, j, L, k, m;
const double E=1e-5;
double x[n], Xprev[n], F[n], J[n][n], a[n][n+1], deltaX[n];
double f(int i,double x[n]);
double Jacobi(int i,int j,double x[n]);
double max(double x[n],double Xprev[n]);

x[0]=0.00000001;     //pochatkove
x[1]=0.00000001;     //nablyzhennya
m=0;

do
{
  for(i=0;i<=n-1;i++)
    Xprev[i]=x[i];
  //----------FORMUYEMO F TA J------------------------
  for(i=0;i<=n-1;i++)
  {
   F[i]=f(i,x);
   for(j=0;j<=n-1;j++)
    J[i][j]=Jacobi(i,j,x);
  }
  //--------------------------------------------------
  for(i=0;i<=n-1;i++)
   for(j=0;j<=n-1;j++)
    a[i][j]=J[i][j];
  a[0][n]=F[0];
  a[1][n]=F[1];

   //------------- GAUSS PO STOVPTSYAH----------------------
  for(L=0;L<=n-2;L++)
  {
    for(k=L+1;k<=n;k++)
    {
      a[L][k]=-a[L][k]/a[L][L];
      for(i=L+1;i<=n-1;i++)
        a[i][k]=a[i][k]+a[i][L]*a[L][k];
    }
  }

    // ------ZVOROTNIY HID-------
  deltaX[n-1]=-a[n-1][n]/a[n-1][n-1];
  for(i=n-2;i>=0;i--)
  {
    deltaX[i]=a[i][n];
    for(k=i+1;k<=n-1;k++)
      deltaX[i]+=a[i][k]*deltaX[k];
  }

  for(i=0;i<=n-1;i++)
    x[i]+=deltaX[i];    //--UTOCHNENE ZNACHENNYA x[i]--
  m++;
}
while(max(x,Xprev)>E);

printf("\t ---ROZVYAZOK SYSTEMY---\n\t|   (ZA %d ITERATSII):   |\n",m);
for(i=0;i<=n-1;i++)
  printf("\t|    x[%d]= %lf    |\n",i,x[i]);
printf("\t -----------------------");
printf("\nPEREVIRKA_1!!!  %lf\n",-x[0]+x[0]*x[0]-x[1]*x[1]-0.1);
printf("PEREVIRKA_2!!!  %lf\n",-x[1]+2*x[0]*x[1]-0.1);

getch();
}

double f(int i,double x[n])
{
 if(i==0)
  return -x[0]+x[0]*x[0]-x[1]*x[1]-0.1;
 else
  return -x[1]+2*x[0]*x[1]-0.1;
}
double Jacobi(int i,int j,double x[n])
{
 const double h=1e-9;
 double prev=f(i,x);
 x[j]+=h;
 return (f(i,x)-prev)/h;
}
double max(double x[n],double Xprev[n])
{
int i;
double z,q;
z=fabs((x[0]-Xprev[0])/Xprev[0]);
for(i=1;i<=n-1;i++)
{
 q=fabs((x[i]-Xprev[i])/Xprev[i]);
 if(q>z)
  z=q;
}
return z;
}

#include <stdio.h>
#include <math.h>
#include <conio.h>
#define n 2
#define pi 3.14159265
void main(void)
{
int i;
double x,a,b,h,y[n],yn[n],yz[n];
double dy(int i,double x,double y[n]);
double modyficEYLER(int i,double h,double x,double y[n],double yz[n]);
a=pi/4;
b=pi/3;
h=0.026;
x=a;
y[0]=0.0;
y[1]=2.0;
printf("x= %lf\ty[1]= %lf\ty[2]= %lf\n",x,y[0],y[1]);
do
{
 printf("x= %lf\t",x+h);
 for(i=0;i<=n-1;i++)
 {
  yz[i]=y[i];
  yn[i]=modyficEYLER(i,h,x,y,yz);
  printf("y[%d]= %lf\t",i+1,yn[i]);
 }
 printf("\n");
 for(i=0;i<=n-1;i++)
  y[i]=yn[i];
 x+=h;
}
while(x<=(b-h));
printf("\n");

getch();
}

double dy(int i,double x,double y[n])
{
 if(i==0)
  return y[1];
 else
  return -4*y[0]+8*(1/tan(2*x));
}

double modyficEYLER(int i,double h,double x,double y[n],double yz[n])
{
 yz[i]=yz[i]+h*dy(i,x,yz);
 return y[i]+0.5*h*( dy(i,x,y)+dy(i,x+h,yz) );
}

program lab4_11v;
uses Crt;
const a = 0;
      b = 3;
      E = 0.001;

 function arcsin(x:real):real;
 begin
  if x=1.0 then
     arcsin := Pi/2.0
  else
     arcsin := arctan(x/sqrt(1-x*x));
 end;

 function f(x:double):double;
 begin
  f:= arcsin(sqrt(x/(1+x)));
 end;

 var
  S, h,n,M, x,Xn,I : double;

begin
 ClrScr;

 M := 0.05413; {max|f''(x)|}
 h := sqrt( (24*E) / ((b-a)*M) ); {krok}
 n := (b-a)/h; {4aCTuH}

 S := 0;
 x := a;
 repeat
 begin
     Xn := x+h;
      S := S + f((x+Xn)/2);
      x := Xn;
 end;
 until (Xn <= (b-h));

 I:= h * S;
 writeln('I = ',I:5:5);

 readln;
end.