#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#define MAXN 100000


int jacobi(int n, double *x, double *b, double *xpom, double eps);
void uloz_reseni(const char *fname, int n, double *x);

main()
{
    const maxiter = 15000;
    int n, i, k, num;
    double rez, h;
    double X[MAXN], XPOM[MAXN], B[MAXN];

   n = 50;  
   h = 1. / (n + 1);

   for(i = 1; i < n; i++)
       B[i] = h * h * sin( i * h * M_PI);

   for(i = 0; i <= n; i++)
       X[i] = 0;

   num = jacobi(n, X, B, XPOM, 1e-12);

   printf("Jac. met pro %d x %d rez 1e-12 behem %d iteraci\n", n, n, num);

   uloz_reseni("reseni.dat", n, X);
}


int jacobi(int n, double *x, double *b, double *xpom, double eps)
{
   int i, iter;
   double rez, rezi; 
   rez = 1 + eps; iter = 0;

   while (rez > eps)
   {
       rez = 0;
       for(i = 1; i < n; i++)
       {
           xpom[i] =  (b[i] + x[i - 1] + x[i + 1]) / 2;
           rezi = fabs(xpom[i] - x[i]);
           rez  = (rezi > rez) ? rezi : rez;
       }
       printf("iter %d, reziduum %g \n", iter, rez);
       for(i = 1; i < n; i++)
           x[i] = xpom[i];
       iter++;
   }
   return iter;
}

void uloz_reseni(const char *fname, int n, double *x)
{
    FILE *fid;
    int i;

    printf("Reseni ukladam do souboru: %s \n", fname);
    fid = fopen(fname, "w");
    if (fid == NULL)
        exit(1);
    for(i = 1; i < n; i++)
        fprintf(fid, "%20.10g\n", x[i]);
    fclose(fid);
    
}