//The test problem is explained in detail below at the printout

#include <iostream.h>    //these four files standard for GENPREC
#include <stdlib.h>     //abs (absolute value)
#include <math.h>       //log,sqrt,etc
#include <time.h>       //use when timing (benchmarking)

#pragma off (check_stack);  //This can save a few ticks

const int precision = 128;  
const int zmaxsize=51;      
const int cmaxsize=0;
int qqq;

#include "genprec.h"       //use this if you want assembly language
//#include "genpnoas.sem"  //use this is you don't want assembly language
#include "cgenprec.h"
#include "apps.h"

int n;     //size of system       

realmp f2,lam2; 
realmp r[zmaxsizepad];                
realmp m[zmaxsizepad][zmaxsizepad];
realmp solution[zmaxsizepad];

main()
{ int i,j,k,l;
  double fs,lams,is;           //user defined ordinary C variables

  zmemory(1,3,2,0,0,0);        //for genlind.h
  zmemory(2,1,1,0,0,0);
  clock_t start_time, end_time;   //declaration statement for times
  start_time=clock();             //set start time
  n=50;                           //size of system in this code
  for (i=1;i<n+1;i++)             //set m=i*j + 1 delta(i,j)
      {r[i]=i;
       for (j=1;j<n+1;j++)
         {m[i][j]=i*j;
          if (i==j) m[i][j]=m[i][j]+1;
         }
      }                           //end of set m

/*
   for (i=1;i<n+1;i++)           //print r and m
     {
      cout<<i<<" ";
      printmp(r[i]);
      cin>>qqq;
      for (j=1;j<n+1;j++)
        {
         cout<<i<<" "<<j<<" ";
         printmp(m[i][j]);
         cin>>qqq;}
     }
*/

  solution<solvlin(n,m,r);      //solves the problem.
                                //< is overloaded in overpoly.new

  end_time=clock();
  cout<<"clock ticks"<<'\t'<<end_time - start_time<<'\n';
  cout<<"ticks/second"<<'\t'<<CLOCKS_PER_SEC<<'\n';
  cout<<"execution time"<<'\t'<<(end_time - start_time)/CLOCKS_PER_SEC<<'\t'<<"seconds"<<'\n';
  zmessage="timing data";
  zbreak();

  fs=n;                              //set up constants required in
  lams=fs*fs*fs/3 + fs*fs/2 + fs/6;  //exact solution (see manual--
                                     //notes on the test linear system
                                     //at end of section on examples)
  lams=lams/(1+lams);                //lambda=sum(j^2)/(1+sum(j^2))
  f2=n;      //The above was done in C++. Now do the same in genprec
  lam2=f2/6;
  lam2=lam2 + f2*f2/2;
  lam2=lam2 + f2*f2*f2/3;
  lam2=lam2/(1+lam2);
  cout<<"  "<<'\n';
  cout<<"The linear system (i*j + delta(i,j))*x(j)=i,1<=i,j<=n can be"<<'\n';
  cout<<"solved exactly for the x(j).  This problem is used as a test of"<<'\n';
  cout<<"genlin.h. The solution found by genlin.h is compared to the exact"<<'\n'; 
  cout<<"solution below for n=50. The MP version and the floating point"<<'\n'; 
  cout<<"version of the exact solution are both given."<<'\n';
  cout<<"  "<<'\n';
  cout<<"After each x(i),i=1-50,is printed enter 1 from keyboard to continue"<<'\n';
  cout<<" "<<'\n';
  zmessage=" ";
  zbreak();

  for (i=1;i<n+1;i++)
     {printmp(solution[i]);    
      printmp(i*(1-lam2));
      is=i;
      cout << is*(1-lams) << '\n';
      zmessage="solution followed by exact solution followed by float";
      zbreak();
     }
   return 0;
}