/**
  * Tömörített ritkamátrix reprezentációban megadott mátrix rendes mátrix-á alakítása.
  * A ritka mátrix tárolása MatrixMarket(MM) formátumban legyen.
  * A program csak adott típusú mátrixra működik:
  * %%MatrixMarket matrix coordinate real symmetric
  * ahol
  * %%MatrixMarket:   MM-ben van megadva
  * matrix:           mátrixos alak.
  * coordinate:       koordináták vannak megadva.
  * real:             a mátrix elemei valós számok.
  * symmetric:        a mátrix szimmetrikus.
  *
  */
  
  /**
    * A program indítása:
    * ./mm2matrix param1 param2
    * ahol
    * param1: a tömörített ritka mátrix, ebben a fájban van: *.mtx
    * param2: a tömörített mátrixból az eredeti mátrixot ebbe a fájlba írja.
    */
typedef struct sor
{
       int oszlop;
       double ertek;
       struct sor *kov;
} SOR;
  
#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#include <string.h>
#include <time.h>

int hiba=0;

FILE *fin, *fout;
void beolvas2 ( char a [], char b [], int mit, int maxit, double e, double w);
double ELEM ( SOR ** f, int j );
void TOROL_LISTA ( SOR ** f );
void BESZUR_VEGE ( SOR ** f, int o, double e );
void rendez( SOR ** f, int n);
void Beszuro_rendezes(double *A, double *oszlop, int hossz);
void beolvas ( char a [], char b [], int mit, int maxit, double e, double w);
void Inkomplett_Cholesky( double *A, int n );
int sorokSzama ( char a[], int mit );

int main ( int argc, char *argv [] )
{
  int mit=0, reprezentacio=0, maxit=0;
  double e=0.0, w=0.0;
  
  printf("Valassz, hogy milyen formaban akarod kiiratni a(z) %s -be a mátrixot/mátrixokat.\n", argv [2]);
  puts("Valaszz [1|2]\n");
  while ( reprezentacio<1 || reprezentacio>2 ){
    puts("1. A teljes matrixot kiiratasa.");
    puts("2. Ritka matrix reprezentacioban.");
    scanf("%d", &reprezentacio);
    puts("");
  }
  
  printf("Milyen módszerhez kell majd a bemenet!\n");
  puts("Valaszz 1-6 ig\n");
  while ( mit<1 || mit>7 ){
    puts("1. LDLT.");                                                                             /* Bemenet: n A b */
    puts("2. Jacobi iteracio / Jacobi Seidel.");                                                  /* Bemenet: n A b x0 maxit e */ 
    puts("3. SOR iteracio.");                                                                     /* Bemenet: n A b x0 maxit e w */ 
    puts("4. Steepest Descent / Csokkeno iranyok modszere.");                                     /* Bemenet: n A b x0 maxit */ 
    puts("5. Konjugalt gradiens.");                                                               /* Bemenet: n A b x0 maxit e */ 
    puts("6. Prek. konj. grad. A prek. mátrix a mátrix főátlója.");                               /* Bemenet: n A H b x0 maxit e */
    puts("7. Prek. konj. grad. A prek. mátrix az inkomlett Cholesky módzserrel meghatározva.");   /* Bemenet: n A H b x0 maxit e */
    scanf("%d", &mit);
    puts("");
  }
  if( mit == 2 || mit == 3 || mit == 4 || mit == 5 || mit == 6 || mit == 7 ){
    puts("Adja meg, hogy mennyi legyen a maximális iterációszám.");
    scanf("%d", &maxit);
  }
  if ( mit == 2 || mit == 3 || mit == 5 || mit == 6 || mit == 7 ){
    puts("Adja meg, hogy mekkora legyen az e pontosság.");
    scanf("%lf", &e);
  }
  if ( mit == 3 ){
    puts("Adja meg, hogy mekkora legyen a w érétke.");
    scanf("%lf", &w);
  }
  if ( reprezentacio==1 )
    beolvas ( argv[1], argv [2], mit, maxit, e, w );
  else if ( reprezentacio==2 )
    beolvas2 ( argv[1], argv [2], mit, maxit, e, w );

return 0;
}
/*______________________________________________________________________________________*/
void beolvas2 ( char a [], char b [], int mit, int maxit, double e, double w){
    puts("\nA kiiratas formaja:");
    if (mit==1){
      /* LDLT. */      
      printf("n\nA\nb\n");
    }
    /* Jacobi iteracio / Jacobi Seidel. */
    else if(mit==2){
      printf("n\nA\nb\nx0\nmaxit\ne\n");     
    }
    /* SOR iteracio. */
    else if(mit==3){
      printf("n\nA\nb\nx0\nmaxit\ne\nw\n");
    }
    /* Steepest Descent / Csokkeno iranyok modszere. */
    else if(mit==4){
      printf("n\nA\nb\nx0\nmaxit\n");      
    }
    /* Konjugalt gradiens. */
    else if(mit==5){
      printf("n\nA\nb\nx0\nmaxit\ne\n");
    }
    /* Prekondicionalt konjugalt gradiens. */
    else if(mit==6 || mit==7){
      printf("n\nA\nH\nb\nx0\nmaxit\ne\n");
    }
    
    char  Data[100];
	  char *ReturnCode;
    int hanyadiksor = 1, megjSSzama=sorokSzama(a, 1), osszsor=sorokSzama(a, 2);
    int i, j, n, nnz;
    
    i=j=n=0;
        
    fin=fopen ( a , "r" );
    /* Átugrom a megjegyzéseket. */
    while ( (ReturnCode = fgets(Data, 100, fin)) != NULL && hanyadiksor < megjSSzama ){ hanyadiksor++; /*printf("%s",ReturnCode);*/ }
    
    /* Belovasom, hogy hányszor hányas */
    char * tok = strtok(ReturnCode, " ");   n = atoi(tok);
    tok = strtok(NULL, " ");
    tok = strtok(NULL, " ");  nnz = atoi(tok);
    hanyadiksor++;
    
    double (*Matrix)[n] = NULL;
    Matrix = calloc (n, sizeof (*Matrix));
    
    SOR *sorok[n];
    for(i=0;i<n;i++)
      sorok[i] = NULL;
    
    fout=fopen(b, "w");
    fprintf (fout, "%d\n%d", n, nnz);    fputs("\n", fout);
    /* Beolvasom a sor oszlop érték hármasokat. */
    /* akkor úgy ahogy van, kiirom. */
    /*if ( mit != 6 && mit != 7 ){ */
      double ertek;
      while ( (ReturnCode = fgets(Data, 100, fin)) != NULL && hanyadiksor < osszsor ){
        tok = strtok(ReturnCode, " ");    i = atoi(tok);
        tok = strtok(NULL, " ");          j = atoi(tok);
        tok = strtok(NULL, " ");
        ertek = atof(tok);
        Matrix[i-1][j-1] = ertek;     Matrix[j-1][i-1] = ertek;
        fprintf (fout, "%d %d %.15lf\n", i-1, j-1, ertek);
        hanyadiksor++;
      }
    /*}*/
    /* le kell mentenem ritka mátrix reprezentációban és ki kell irnom. */
    /*else{
      double ertek;
      while ( (ReturnCode = fgets(Data, 100, fin)) != NULL && hanyadiksor < osszsor ){
        tok = strtok(ReturnCode, " ");    i = atoi(tok);
        tok = strtok(NULL, " ");          j = atoi(tok);
        tok = strtok(NULL, " ");
        ertek = atof(tok);
        fprintf (fout, "%d %d %.15lf\n", i-1, j-1, ertek);
        printf ("%d %d %.15lf\n", i-1, j-1, ertek);
        Matrix[i-1][j-1] = ertek;     Matrix[j-1][i-1] = ertek;
        hanyadiksor++;
      }
      }*/

    
    fflush(fout);  
        
    /* Kiirom a param2-be a mátrixot + ami kell még. */
    /* LDLT. */
    if (mit==1){
      /* b */
      double temp;
      for ( i=0; i<n; i++ ){
        temp = 0.0;
        for ( j=0; j<n; j++ )
          temp += /*ELEM(&sorok[j], i);*/Matrix[j][i];
        fprintf (fout, "%lf ", temp);
      }
      fputs("\n", fout);
      
      
    }
    /* Jacobi iteracio / Jacobi Seidel. */
    else if(mit==2){
      /* b */
      double temp;
      for ( i=0; i<n; i++ ){
        temp = 0.0;
        for ( j=0; j<n; j++ )
          temp += /*ELEM(&sorok[j], i);*/Matrix[j][i];
        fprintf (fout, "%lf ", temp);
      }
      fputs("\n", fout);
      
      fflush(fout);
      /* x0 */
      for (i=0; i<n; i++ ){
        fprintf (fout, "%d", i);
        if ( (i+1) < n ) fprintf (fout, " ");
      }
      fputs("\n", fout);
      
      fflush(fout);
      /* maxit */
      fprintf (fout, "%d", maxit);
      fputs("\n", fout);
      
      /* e */
      fprintf (fout, "%lf", e);
      fputs("\n", fout);
      
    }
    /* SOR iteracio. */
    else if(mit==3){
      /* b */
      double temp;
      for ( i=0; i<n; i++ ){
        temp = 0.0;
        for ( j=0; j<n; j++ )
          temp += /*ELEM(&sorok[j], i);*/Matrix[j][i];
        fprintf (fout, "%lf ", temp);
      }
      fputs("\n", fout);
      
      fflush(fout);
      /* x0 */
      for (i=0; i<n; i++ ){
        fprintf (fout, "%d", i);
        if ( (i+1) < n ) fprintf (fout, " ");
      }
      fputs("\n", fout);
      
      fflush(fout);
      /* maxit */
      fprintf (fout, "%d", maxit);
      fputs("\n", fout);
      
      /* e */
      fprintf (fout, "%lf", e);
      fputs("\n", fout);
      
      /* w */
      fprintf (fout, "%lf", w);
      fputs("\n", fout);
      
    }
    /* Steepest Descent / Csokkeno iranyok modszere. */
    else if(mit==4){
      /* b */
      double temp;
      for ( i=0; i<n; i++ ){
        temp = 0.0;
        for ( j=0; j<n; j++ )
          temp += /*ELEM(&sorok[j], i);*/Matrix[j][i];
        fprintf (fout, "%lf ", temp);
      }
      fputs("\n", fout);
      
      fflush(fout);      
      /* x0 */
      for (i=0; i<n; i++ ){
        fprintf (fout, "%d", i);
        if ( (i+1) < n ) fprintf (fout, " ");
      }
      fputs("\n", fout);
      
      /* maxit */
      fprintf (fout, "%d", maxit);
      fputs("\n", fout);
      
    }
    /* Konjugalt gradiens. */
    else if(mit==5){
      /* b */
      double temp;
      for ( i=0; i<n; i++ ){
        temp = 0.0;
        for ( j=0; j<n; j++ )
          temp += /*ELEM(&sorok[j], i);*/Matrix[j][i];
        fprintf (fout, "%lf ", temp);
      }
      fputs("\n", fout);
      
      fflush(fout);
      /* x0 */
      for (i=0; i<n; i++ ){
        fprintf (fout, "%d", i);
        if ( (i+1) < n ) fprintf (fout, " ");
      }
      fputs("\n", fout);
      
      /* maxit */
      fprintf (fout, "%d", maxit);
      fputs("\n", fout);
      
      /* e */
      fprintf (fout, "%lf", e);
      fputs("\n", fout);
      
    }
    /* Prekondicionalt konjugalt gradiens. */
    else if(mit==6 || mit==7){      
      /* H */
      if ( mit==6 ){
        int nemnulladb=0;
        for (i=0; i<n;i++){
          for(j=0; j<n;j++){
            if ( i == j && fabs(Matrix[i][j]) > pow(10,-15))
              nemnulladb++;
          }
        }
        fprintf (fout, "%d\n", nemnulladb );
        for (i=0; i<n;i++){
          for(j=0; j<n;j++){
            if ( i == j )
              fprintf (fout, "%d %d %.15lf\n", i, j, Matrix[i][j] );
          }
        }
      }
      else if ( mit == 7 ){
        double (*Temp)[n] = NULL;
        Temp = calloc (n, sizeof (*Matrix));
        
        for (i=0; i<n;i++){
          for(j=0; j<n;j++){
            Temp[i][j] = Matrix[i][j];
          }
        }
        Inkomplett_Cholesky ( (double *)Temp, n );
        if ( hiba == 1 ) { puts("hiba"); return; }
        int nemnulladb=0;
        for (i=0; i<n;i++){
          for(j=0; j<n;j++){
            if ( fabs(Temp[i][j]) > pow(10,-15) )
              nemnulladb++;
          }
        }
        fprintf (fout, "%d\n", nemnulladb );
        for (i=0; i<n;i++){
          for(j=0; j<n;j++){
            if ( fabs(Temp[i][j]) > pow(10,-15) )
              fprintf (fout, "%d %d %.15lf\n", i, j, Temp[i][j]);
          }
        }
      }
      
      fflush(fout);
      /* b */
      double temp;
      for ( i=0; i<n; i++ ){
        temp = 0.0;
        for ( j=0; j<n; j++ )
          temp += /*ELEM(&sorok[j], i);*/Matrix[j][i];
        fprintf (fout, "%lf ", temp);
      }
      fputs("\n", fout);
      
      fflush(fout);
      /* x0 */
      for (i=0; i<n; i++ ){
        fprintf (fout, "%d", i);
        if ( (i+1) < n ) fprintf (fout, " ");
      }
      fputs("\n", fout);
      
      /* maxit */
      fprintf (fout, "%d", maxit);
      fputs("\n", fout);
      
      fflush(fout);
      /* e */
      fprintf (fout, "%lf", e);
      fputs("\n", fout);
      
    }      
    
    fclose(fin);
    fflush(fout);
    fclose(fout);
    for( i=0; i<n; i++ ){ TOROL_LISTA(&sorok[i]); }
return ;
}
/*____________________________________________________________________________*/
double ELEM ( SOR ** f, int j ){
  if( *f==NULL ) return 0.0;
  
  SOR *x;
  x = *f;
  if( x->oszlop > j ){  return 0.0;}
  while ( x != NULL ){
    if ( x->oszlop == j ){
      return x->ertek;}
    else if( x->oszlop >j ){
      return 0.0;}
    x=x->kov;
  }
  return 0.0;
}
/*____________________________________________________________________________*/
void TOROL_LISTA ( SOR ** f ){
  SOR *temp;
  
  while ( *f != NULL)
  {  
     temp = (*f) -> kov;
     free (*f);
     *f = temp;
  }
  return;
}
/*____________________________________________________________________________*/
void BESZUR_VEGE ( SOR ** f, int o, double e ){
  SOR *uj;
  uj = (SOR *) malloc (sizeof (SOR));
  uj->oszlop=o;
  uj->ertek=e;
    
  uj->kov=*f;
  *f=uj;
  return ;
}
/*______________________________________________________________________________________*/
void rendez( SOR ** f, int n){
  int i;
  
  for( i=0; i<n; i++){
    SOR *x=f[i];
    int sorhossz;
    for( sorhossz=0; x!=NULL; sorhossz++ ) x=x->kov;
    if ( sorhossz > 0 ){
      double *rendezett_oszlop = (double *) calloc( sorhossz, ( sizeof ( double ) ) );
      double *rendezett_ertek = (double *) calloc( sorhossz, ( sizeof ( double ) ) );
      for( x=f[i], sorhossz=0; x!=NULL; sorhossz++, x=x->kov ){
        rendezett_ertek[sorhossz]=x->ertek;
        rendezett_oszlop[sorhossz]=x->oszlop;
      }
      Beszuro_rendezes(rendezett_ertek, rendezett_oszlop, sorhossz);
      for( x=f[i], sorhossz=0; x!=NULL; sorhossz++, x=x->kov ){
        x->ertek=rendezett_ertek[sorhossz];
        x->oszlop=rendezett_oszlop[sorhossz];
      }
      free(rendezett_ertek);
      free(rendezett_oszlop);
    }    
  }
}
/*______________________________________________________________________________________*/
void Beszuro_rendezes(double *A, double *oszlop, int hossz){
  int i,j;
  double kulcs_e, kulcs_o;
  for( j=1; hossz; j++ ){
    kulcs_e=A[j];
    kulcs_o=oszlop[j];
    i=j-1;
    while ( i>0 && A[i]>kulcs_e ){
      A[i+1]=A[i];
      oszlop[i+1]=oszlop[i];
      i=i-1;
    }
    A[i+1]=kulcs_e;
    oszlop[i+1]=kulcs_o;
  }
}
/*______________________________________________________________________________________*/
void beolvas ( char a [], char b [], int mit, int maxit, double e, double w){
    char  Data[100];
	  char *ReturnCode;
    int hanyadiksor = 1, megjSSzama=sorokSzama(a, 1), osszsor=sorokSzama(a, 2);
    int i, j, n;
    
    i=j=n=0;
    
    fin=fopen ( a , "r" );
    /* Átugrom a megjegyzéseket. */
    while ( (ReturnCode = fgets(Data, 100, fin)) != NULL && hanyadiksor < megjSSzama ){ hanyadiksor++; /*printf("%s",ReturnCode);*/ }
    
    /* Belovasom, hogy hányszor hányas */
    char * tok = strtok(ReturnCode, " ");   n = atoi(tok);
    hanyadiksor++;

    double (*Matrix)[n] = NULL;
    Matrix = calloc (n, sizeof (*Matrix));
    
    /* Beolvasom a sor oszlop érték hármasokat. */
    while ( (ReturnCode = fgets(Data, 100, fin)) != NULL && hanyadiksor < osszsor){
      tok = strtok(ReturnCode, " ");    i = atoi(tok);
      tok = strtok(NULL, " ");          j = atoi(tok);
      tok = strtok(NULL, " ");          /*puts("tok vége");*/   Matrix[i-1][j-1] = atof(tok);     Matrix[j-1][i-1] = atof(tok);
            
      hanyadiksor++;
    }
    
    fout=fopen(b, "w");
    
    /* Kiirom a param2-be a mátrixot + ami kell még. */
    /* LDLT. */
    if (mit==1){
      puts("\nA kiiratas formaja:");
      printf("n\nA\nb\n");
      
      /* n */
      fprintf (fout, "%d", n);    fputs("\n", fout);
      
      /* A */
      for (i=0; i<n;i++){
        for(j=0; j<n;j++){
          fprintf (fout, "%.15lf",Matrix[i][j]);
          if ( (j+1) < n ) fprintf (fout, " ");
        }
        fputs("\n", fout);
      }
      
      fflush(fout);
      /* b */
      double temp;
      for ( i=0; i<n; i++ ){
        temp = 0.0;
        for ( j=0; j<n; j++ )
          temp += Matrix[j][i];
        fprintf (fout, "%lf ", temp);
      }
      fputs("\n", fout);
      
      
    }
    /* Jacobi iteracio / Jacobi Seidel. */
    else if(mit==2){
      puts("\nA kiiratas formaja:");
      printf("n\nA\nb\nx0\nmaxit\ne\n");
      
      /* n */
      fprintf (fout, "%d", n);    fputs("\n", fout);
      
      /* A */
      for (i=0; i<n;i++){
        for(j=0; j<n;j++){
          fprintf (fout, "%.15lf",Matrix[i][j]);
          if ( (j+1) < n ) fprintf (fout, " ");
        }
        fputs("\n", fout);
      }
      
      fflush(fout);
      /* b */
      double temp;
      for ( i=0; i<n; i++ ){
        temp = 0.0;
        for ( j=0; j<n; j++ )
          temp += Matrix[j][i];
        fprintf (fout, "%lf ", temp);
      }
      fputs("\n", fout);
      
      fflush(fout);
      /* x0 */
      for (i=0; i<n; i++ ){
        fprintf (fout, "%d", i);
        if ( (i+1) < n ) fprintf (fout, " ");
      }
      fputs("\n", fout);
      
      fflush(fout);
      /* maxit */
      fprintf (fout, "%d", maxit);
      fputs("\n", fout);
      
      /* e */
      fprintf (fout, "%lf", e);
      fputs("\n", fout);
      
    }
    /* SOR iteracio. */
    else if(mit==3){
      puts("\nA kiiratas formaja:");
      printf("n\nA\nb\nx0\nmaxit\ne\nw\n");
      
      /* n */
      fprintf (fout, "%d", n);    fputs("\n", fout);
      
      /* A */
      for (i=0; i<n;i++){
        for(j=0; j<n;j++){
          fprintf (fout, "%.15lf",Matrix[i][j]);
          if ( (j+1) < n ) fprintf (fout, " ");
        }
        fputs("\n", fout);
      }
      
      fflush(fout);
      /* b */
      double temp;
      for ( i=0; i<n; i++ ){
        temp = 0.0;
        for ( j=0; j<n; j++ )
          temp += Matrix[j][i];
        fprintf (fout, "%lf ", temp);
      }
      fputs("\n", fout);
      
      fflush(fout);
      /* x0 */
      for (i=0; i<n; i++ ){
        fprintf (fout, "%d", i);
        if ( (i+1) < n ) fprintf (fout, " ");
      }
      fputs("\n", fout);
      
      fflush(fout);
      /* maxit */
      fprintf (fout, "%d", maxit);
      fputs("\n", fout);
      
      /* e */
      fprintf (fout, "%lf", e);
      fputs("\n", fout);
      
      /* w */
      fprintf (fout, "%lf", w);
      fputs("\n", fout);
      
    }
    /* Steepest Descent / Csokkeno iranyok modszere. */
    else if(mit==4){
      puts("\nA kiiratas formaja:");
      printf("n\nA\nb\nx0\nmaxit\n");
      
      /* n */
      fprintf (fout, "%d", n);    fputs("\n", fout);
      
      /* A */
      for (i=0; i<n;i++){
        for(j=0; j<n;j++){
          fprintf (fout, "%.15lf",Matrix[i][j]);
          if ( (j+1) < n ) fprintf (fout, " ");
        }
        fputs("\n", fout);
      }
      
      fflush(fout);
      /* b */
      double temp;
      for ( i=0; i<n; i++ ){
        temp = 0.0;
        for ( j=0; j<n; j++ )
          temp += Matrix[j][i];
        fprintf (fout, "%lf ", temp);
      }
      fputs("\n", fout);
      
      fflush(fout);      
      /* x0 */
      for (i=0; i<n; i++ ){
        fprintf (fout, "%d", i);
        if ( (i+1) < n ) fprintf (fout, " ");
      }
      fputs("\n", fout);
      
      /* maxit */
      fprintf (fout, "%d", maxit);
      fputs("\n", fout);
      
    }
    /* Konjugalt gradiens. */
    else if(mit==5){
      puts("\nA kiiratas formaja:");
      printf("n\nA\nb\nx0\nmaxit\ne\n");
      
      /* n */
      fprintf (fout, "%d", n);    fputs("\n", fout);
      
      /* A */
      for (i=0; i<n;i++){
        for(j=0; j<n;j++){
          fprintf (fout, "%.15lf",Matrix[i][j]);
          if ( (j+1) < n ) fprintf (fout, " ");
        }
        fputs("\n", fout);
      }
      
      fflush(fout);
      /* b */
      double temp;
      for ( i=0; i<n; i++ ){
        temp = 0.0;
        for ( j=0; j<n; j++ )
          temp += Matrix[j][i];
        fprintf (fout, "%lf ", temp);
      }
      fputs("\n", fout);
      
      fflush(fout);
      /* x0 */
      for (i=0; i<n; i++ ){
        fprintf (fout, "%d", i);
        if ( (i+1) < n ) fprintf (fout, " ");
      }
      fputs("\n", fout);
      
      /* maxit */
      fprintf (fout, "%d", maxit);
      fputs("\n", fout);
      
      /* e */
      fprintf (fout, "%lf", e);
      fputs("\n", fout);
      
    }
    /* Prekondicionalt konjugalt gradiens. */
    else if(mit==6 || mit==7){
      puts("\nA kiiratas formaja:");
      printf("n\nA\nH\nb\nx0\nmaxit\ne\n");
      
      /* n */
      fprintf (fout, "%d", n);    fputs("\n", fout);
      
      /* A */
      for (i=0; i<n;i++){
        for(j=0; j<n;j++){
          fprintf (fout, "%.15lf",Matrix[i][j]);
          if ( (j+1) < n ) fprintf (fout, " ");
        }
        fputs("\n", fout);
      }
      
      puts("\nfflush:");
      fflush(fout);
      
      /* H */
      if ( mit==6 ){
        for (i=0; i<n;i++){
          for(j=0; j<n;j++){
            if ( i == j )
              fprintf (fout, "%.15lf",Matrix[i][j]);
            else fprintf (fout, "%.15lf", 0.0);
            if ( (j+1) < n ) fprintf (fout, " ");
          }
          fputs("\n", fout);
        }
      }
      else if ( mit == 7 ){
        double (*Temp)[n] = NULL;
        Temp = calloc (n, sizeof (*Matrix));
        
        for (i=0; i<n;i++){
          for(j=0; j<n;j++){
            Temp[i][j] = Matrix[i][j];
          }
        }
        Inkomplett_Cholesky ( (double *)Temp, n );
        if ( hiba == 1 ) { puts("hiba"); return; }
        for (i=0; i<n;i++){
          for(j=0; j<n;j++){
              fprintf (fout, "%.15lf",Temp[i][j]);
            if ( (j+1) < n ) fprintf (fout, " ");
          }
          fputs("\n", fout);
        }
      }
      
      fflush(fout);
      /* b */
      double temp;
      for ( i=0; i<n; i++ ){
        temp = 0.0;
        for ( j=0; j<n; j++ )
          temp += Matrix[j][i];
        fprintf (fout, "%lf ", temp);
      }
      fputs("\n", fout);
      
      fflush(fout);
      /* x0 */
      for (i=0; i<n; i++ ){
        fprintf (fout, "%d", i);
        if ( (i+1) < n ) fprintf (fout, " ");
      }
      fputs("\n", fout);
      
      /* maxit */
      fprintf (fout, "%d", maxit);
      fputs("\n", fout);
      
      fflush(fout);
      /* e */
      fprintf (fout, "%lf", e);
      fputs("\n", fout);
      
    }     
    
    fclose(fin);
    fflush(fout);
    fclose(fout);
    free(Matrix);
return ;
}
/*____________________________________________________________________________*/
void Inkomplett_Cholesky( double *A, int n ){
    int i, j, k;
    i = j = k = 0;

    /* Az algoritmus */
    for ( k=0; k<n; k++){
      A[k*n+k] = sqrt( A[k*n+k] );
      for( i=k+1; i<n; i++){
        if ( fabs( A[i*n+k] ) > pow ( 10,-15 ) ){
          A[i*n+k] = A[i*n+k] / A[k*n+k];
        }
      }
      for( j=k+1; j<n; j++ ){
        for( i=j; i<n; i++ ){
          if ( fabs( A[i*n+j] ) > pow ( 10,-15 ) ){
            A[i*n+j] = A[i*n+j] - ( A[i*n+k] * A[j*n+k] );
          }
        }
      }
    }
    
    for ( i=0; i<n; i++){
      for( j=i+1; j<n; j++ ){
        A[i*n+j] = 0.0;
      }
    }
    
return;
}
/*______________________________________________________________________________________*/
int sorokSzama ( char a[], int mit ){
    char  Data[200];
	  char *ReturnCode;
    int sorokSzama = 1;

    fin=fopen ( a , "r" );
    
    if ( fin == NULL )
    {
         fprintf( stderr,"A bemenet megnyitasa sikertelen\n" );
         exit (1);
    }
    
    /* megkegyzések */
    if ( mit == 1 ){
      while ( (ReturnCode = fgets(Data, 200, fin)) != NULL ){
          if ( ReturnCode[0] != '%' )
            return sorokSzama;
          sorokSzama++;
      }
    }
    /* minden sor */
    else if ( mit == 2 ) {
       while ( (ReturnCode = fgets(Data, 200, fin)) != NULL ){
          sorokSzama++;
      }
       return sorokSzama;
    }
    
    fclose(fin);
 
 return sorokSzama;
}
