/**
   * Az Ax=b egyenletrendszer megoldasahoz a Jacobi iteráció.
   * x* az Ax*=b megoldása.
   */
   
   /*
   * Bemenet:
   * n
   * A
   * b
   * x0
   * maxit
   * e
   */
   
   /* Ahol
   * n       : nxn -es mátrix.
   * A       : mátrix
   * b       : Ax=b egyenletrendszer jobb oldala.
   * x0      : kezdővektor
   * maxit   : maximális iteráció száma.
   * e       : küszöbszám. xk és xk+1 közötti különbség ha ettől kisebb, akkor vége.
   */

void Jacobi_iteracio(char in[], char out[]);
void Eredmeny_kiiras( char out[], int n, double * vektor, int it );
double Ellenorzes_legnagyobb_kulonbseg( double * x, int n);
double Euklideszi_norma_vektorra( double * vektor, int o );
void Ellenorzes(double *A, int n, double * b, double * x);

#include <stdio.h>
#include <math.h>
#include <stdlib.h>
#include <string.h>
#include <ctype.h>
/*#include <sys\timeb.h>*/
/*#include "leak_detector_c.h"*/

int main (int argc, char* argv[]){
    if (argc<3){
      puts("Hiba! Kevés argumentumot adott meg.");
      puts("A program helyes futtatása:\nprogramnév bemeneti_fájl_neve_kiterjesztéssel kimeneti_fájl_neve_kiterjesztéssel\n");
      return 1;
     }
     else if (argc>3){
       puts("Hiba! Túl sok argumentumot adott meg.");
       puts("A program helyes futtatása:\nprogramnév bemeneti_fájl_neve_kiterjesztéssel kimeneti_fájl_neve_kiterjesztéssel\n");
       return 1;
     }

    /*struct timeb t_start, t_current;
    int t_diff_in_sec, t_diff_in_ms;
    ftime(&t_start);*/
    
    Jacobi_iteracio(argv[1], argv[2]);
    
    /*ftime(&t_current);
      t_diff_in_sec = (int)  (t_current.time - t_start.time);
      t_diff_in_ms = (int ) (t_current.millitm - t_start.millitm);

    FILE *fout;
    fout = fopen(argv[2], "a");    
    printf("\nEltelt idő: %d ms",t_diff_in_sec*1000+t_diff_in_ms);
    fprintf(fout, "\nEltelt idő: %d ms",t_diff_in_sec*1000+t_diff_in_ms);
    fclose(fout);*/
    
    /*atexit(report_mem_leak);*/
    return 0;
}
/*____________________________________________________________________________*/
void Jacobi_iteracio(char in[], char out[]){
    int n, i, j, k, l, maxit;
    n = i = j = k = l = maxit = 0;
    
    FILE *fin, *fout;
    
    fin = fopen(in, "r");
    fout = fopen(out, "w");
    
    fscanf(fin, "%d\n",&n);
    double (*A)[n] = NULL;
    A = calloc (n, sizeof (*A));
    
    double *b = (double *) calloc( n, ( sizeof ( double ) ) );
    double *x0 = (double *) calloc( n, ( sizeof ( double ) ) );
    double *xk = (double *) calloc( n, ( sizeof ( double ) ) );
    double e;
    
    /*----------------------------------------*/
    /* A matrix beolvasasa */
    for ( i=0; i<n; i++ ){
        for ( j=0; j<n; j++ ){
            fscanf(fin, "%lf", &A[i][j]);
        }
    }
    /* b vektor beolvasása */
    for ( i=0; i<n; i++ ){
    	fscanf(fin, "%lf", &b[i] );
    }
	/* x0 vektor beolvasása */
	for ( i=0; i<n; i++ ){
    	fscanf(fin, "%lf", &x0[i] );
    	xk[i] = x0[i];
	}
	fscanf(fin, "%d",&maxit);
	fscanf(fin, "%lf",&e);
	fclose(fin);
			
	/*----------------------------------------*/
	for ( k=0; k<maxit; k++ ){
        
        for ( i=0; i<n; i++ ){
            double sum=0;
                for ( j=0; j<n; j++ ){
                    if ( j!=i ) sum += A[i][j]*xk[j];
                }
            sum *= -1;
            if ( fabs( A[i][i] ) <= pow ( 10,-15 ) ) {
              puts("hiba"); fputs("hiba", fout);
              fflush(fout);
              fclose(fout);
              free(A); free(b); free(x0); free(xk);
              return;
            }
            xk[i] = (1/A[i][i])*( sum + b[i] );
        }
        double y[n];
        for( i=0; i<n; i++ ){
          y[i]=x0[i]-xk[i];
        }
        double temp2 = Euklideszi_norma_vektorra( y, n );
        temp2 = sqrt(temp2);
        if ( temp2 <=e  ){
          fflush(fout);
          Eredmeny_kiiras( out, n, x0, k );
          /*Ellenorzes( (double *)A, n, b, x0 );*/
          fflush(fout);
          fclose(fout);
          free(A); free(b); free(x0); free(xk);
          return;
        }
        
        for ( j=0; j<n; j++ ) x0[j] = xk[j];
        
        printf("%d.\n", k);
  
        if( k+1 ==maxit ){
    puts("A következő iterációs lépésben véget ér a program!\n");
    printf("A legnagyobb külünbség: kezdővektor-kapottmegoldás: \n");
    puts("Szeretné folytatni az iterációt, annak érdkében, hogy jobb eredményt kapjon?\n");
    int jo=0;
    char valasz[10];
    
    while( jo == 0 ){
      puts("Parancsok: [i|n], [igen|nem], [y|n], [yes|no]");
      scanf("%s", valasz);
      if( strncmp(valasz, "i", strlen("i"))==0 || strncmp(valasz, "igen", strlen("igen"))==0 || 
          strncmp(valasz, "y", strlen("y"))==0 || strncmp(valasz, "yes", strlen("yes"))==0){
        int ujmaxit=0;
        printf("uj maxit: "); scanf("%d", &ujmaxit);
        while( ujmaxit<=maxit ){
          puts("HIBA! A megadott új maximális iterációszám kisebb, mint a jelenlegi!\n");
          printf("Adjon %d -tól/-től nagyobb számot!\n", maxit);
          printf("uj maxit: "); scanf("%d", &ujmaxit);
        }
        maxit=ujmaxit;
        jo=1;
      }
      else if( strncmp(valasz, "n", strlen("n"))==0 || strncmp(valasz, "nem", strlen("nem"))==0 || 
                strncmp(valasz, "no", strlen("no"))==0)
                jo=1;
      else      
        puts("HIBA! Rossz parancsot adott meg!");
                
    }
  }
    }
    
    fflush(fout);
    Eredmeny_kiiras( out, n, x0, k );
    /*Ellenorzes( (double *)A, n, b, x0 );*/
    fflush(fout);
    fclose(fout);
    free(A); free(b); free(x0); free(xk);
     /*Ellenorzes( (double *)A, n, b, x0 );*/
     
return ;
}
/*____________________________________________________________________________*/
void Eredmeny_kiiras( char out[], int n, double * vektor, int it ){
  
  FILE *fout;

  fout = fopen(out, "a");
  
  int i;
  for ( i=0; i<n; i++ ){
    if (i<n-1)
      fprintf(fout, "%.5lf ", vektor[i] );
     else
      fprintf(fout, "%.5lf", vektor[i] );
  }
  fputs("", fout);
  printf("A legnagyobb külünbség: ||kezdővektor-kapottmegoldás||= %.16lf\n", Ellenorzes_legnagyobb_kulonbseg(vektor, n) );
  fprintf(fout, "\n\nA legnagyobb külünbség: ||kezdővektor-kapottmegoldás||= %.16lf\n", Ellenorzes_legnagyobb_kulonbseg(vektor, n) );
  printf("Lépések száma: %d", it);
  fprintf(fout, "Lépések száma: %d", it);
  fflush(fout);
  fclose(fout);
return;
}
/*____________________________________________________________________________*/
double Ellenorzes_legnagyobb_kulonbseg( double * x, int n){
  double *temp=(double *) calloc( n, ( sizeof ( double ) ) );
  int i=0;
  for( ; i<n; i++ )
    temp[i] = 1-x[i];
    
  return sqrt( Euklideszi_norma_vektorra(temp, n) );
  
}
/*____________________________________________________________________________*/
double Euklideszi_norma_vektorra( double * vektor, int o ){
       double temp=0.0;
       int szamlal;
       for ( szamlal=0; szamlal<o; szamlal++ ){
           temp = temp + (vektor[szamlal]*vektor[szamlal]);
       }
       return temp;
}
/*____________________________________________________________________________*/
void Ellenorzes(double *A, int n, double * b, double * x){
  int i, j;
  double temp[n];
  puts("\n\nEllenorzes: Ax=b");    
  puts("b:\t\t\tAx:");
    for ( i=0; i<n; i++ ){
      temp[i]=0;
      for ( j=0; j<n; j++ ){
        temp[i] += A[i*n+j]*x[j];
        /*printf("\n%lf += %lf*%lf ", temp[i], A[i*n+j], x[j]);*/
      }
      printf( "b[%d]=%lf\t\t%lf\t\t\n",i, b[i], temp[i] );
    }
}
