/*
 $Id: utils.cc,v 1.3 1996/10/11 11:30:11 roitzsch Exp $
 (C)opyright 1996 by Konrad-Zuse-Center, Berlin
 All rights reserved.
 Part of the Kaskade distribution
*/

#include <stdarg.h>

#include "utils.h"

static double z_seconds();

MultiCompNode MCNode;

//-------------------------------------------------------------------------

void strToLower(char* s)  { while(*s != '\0') { *s = tolower(*s); ++s; } } 
void strToUpper(char* s)  { while(*s != '\0') { *s = toupper(*s); ++s; } }
 
//-------------------------------------------------------------------------

Timer:: Timer()
{ 
    t0 = z_seconds();
    time(&abst0);
}

float Timer:: cpu(int print)
{ 
    double seconds = z_seconds()-t0;

    if (print) 
    {
        cout << " " << Form("%1.2f", seconds) << " cpu sec.";

        int  minutes = int(seconds)/60;
        float restSeconds = seconds - minutes*60;
        cout << " ("  << minutes << ":" << Form("%1.2f",restSeconds) << ")\n";
    }

    return seconds;
}

long Timer:: total(int print)
{ 
    time_t seconds; 
    time(&seconds);
    seconds = seconds - abst0;
    
    if (print) cout << " " << seconds << " sec. (" << seconds/60 << ":" 
    		    << seconds%60 << ")\n";
    return seconds;
}
//-------------------------------------------------------------------------
//-------------------------------------------------------------------------


Parser:: Parser(const char* file, const char* endToken0, const char commentFlag0)

	: commentFlag(commentFlag0), fp(0)

{
	fp = fopen(file, "r");
	if (fp == 0) 
	{ 
	    cout << "\n*** Parser: " << file << " not found\n";
	    cout.flush(); abort();
	}
    endToken = new char[strlen(endToken0)+1];
    strcpy(endToken, endToken0);

    fileName = new char[strlen(file)+1];
    strcpy(fileName, file);
}
//-------------------------------------------------------------------------

Parser:: ~Parser() { delete fileName;  delete endToken;  if (fp) fclose(fp); }

//-------------------------------------------------------------------------

void Parser:: readError(const char* function, char* word)
{ 
    cout << "\n*** Parser." << function << ": read error on word '" << word
         << "' in file " << fileName << "\n"; 
    cout.flush(); abort(); 
}
//-------------------------------------------------------------------------

int Parser:: skipRestOfLine()
{
    int c;
    while ((c=getc(fp)) != '\n') if (c==EOF) {status = eof; return status;}
    status = ok; 
    return status;
}
//-------------------------------------------------------------------------

int Parser:: getString(char* word)
{
    while (fscanf(fp,"%s",word) != EOF)
    {
	strToLower(word);
	if (word[0] == commentFlag)
	{
	    if (skipRestOfLine() == eof) { status = eof; return status; }
	}
	else
	{
	    if (!strcmp(word,endToken)) status = end;
	    else  status = ok;
	    return status;
	}
    }
    status = eof;
    return status;
}
//-------------------------------------------------------------------------


int Parser:: findWord(const char* target, int completeMatch)
{
    char word[WORDLENGTH], s[WORDLENGTH];
    strcpy(s, target);
    strToLower(s);

    if (completeMatch) 
    {
	while (getString(word) != eof)
	    if (!strcmp(word,s)) { status = ok; return status; }
    }
    else
    {
	int length = strlen(target);

	while (getString(word) != eof)
	    if (!strncmp(word,s,length)) { status = ok; return status; }
    }
    return status;
}
//-------------------------------------------------------------------------


int Parser:: findChar(const char target)
{
    char ch;
    int c;

    ch = getChar(); 
    c = ch;

    while (c != EOF) 
    {
	  if (ch == target) { status = ok; return status; }
          ch = getChar(); 
	  c = ch;
    }
    return status;
}
//-------------------------------------------------------------------------

int Parser:: getValue(char* c)
{
    char word[WORDLENGTH];

    if (getString(word) != ok)  return status;

    if (strlen(word) != 1)        readError("getValue(char)", word);
    if (sscanf(word,"%c",c) != 1) readError("getValue(char)", word);
    status = ok;
    return status;
}

int Parser:: getValue(int* i)
{
    char word[WORDLENGTH];

    if (getString(word) != ok)  return status;
    if (sscanf(word,"%d",i) != 1) readError("getValue(int)", word);
    status = ok;
    return status;
}

int Parser:: getValue(float* x)
{
    char word[WORDLENGTH];

    if (getString(word) != ok)  return status;
    if (sscanf(word,"%f",x) != 1) readError("getValue(float)", word);
    status = ok;
    return status;
}

int Parser:: getValue(double* x)
{
    char word[WORDLENGTH];

    if (getString(word) != ok)  return status;
    if (sscanf(word,"%lf",x) != 1) readError("getValue(double)", word);
    status = ok;
    return status;
}
//-------------------------------------------------------------------------

char Parser:: getChar()
{
    int c;
    if ((c=getc(fp)) == EOF) status = eof;
    else status = ok; 
    return char(c);
}

int Parser:: getLine(char* line)
{
    int count = 0;
    int c;

    while ((c=getc(fp)) != '\n') 
    {
	if (c==EOF) { line[count++] = '\n'; status = eof; return status; }
	line[count++] = char(c);
    }
    line[count++] = '\n';
    status = ok; 
    return status;
}
//-------------------------------------------------------------------------

void Parser:: rewind() { ::rewind(fp); }
long Parser:: ftell () { return ::ftell(fp); }
int  Parser:: fseek (long offset, int mode) { return ::fseek(fp,offset,mode); }

//-------------------------------------------------------------------------


BufferedParser:: BufferedParser(const char* file, const char* endToken0, 
				const char commentFlag0)
	: Parser(file, endToken0, commentFlag0), buffer(0), pos(0)
{
    if (::fseek(fp,0,endOfFile) != 0) 
    {
	cout << "\n*** BufferedParser: fseek failed on " << fileName << "\n"; 
	fclose(fp);
	cout.flush(); abort();
    }

    bufferLength = ::ftell(fp);
    ::rewind(fp);
   
    ++bufferLength;
    buffer = new char[bufferLength];

    if (buffer == 0)
    {
	cout << "\n*** BufferedParser: not enough memory to store " 
	     << bufferLength << " bytes for " << fileName << "\n";
	fclose(fp);
	cout.flush(); abort();
    }

    int size = 1;
    fread(buffer, size, bufferLength-1, fp);

    buffer[bufferLength-1] = '\0';
    pos = buffer;

    fclose(fp); 
    fp = 0;
}
//-------------------------------------------------------------------------

BufferedParser:: ~BufferedParser()  { delete buffer; }

//-------------------------------------------------------------------------

//-------------------------------------------------------------------------

int BufferedParser:: getValue(char* c)
{
    char word[WORDLENGTH];

    posBefore = pos;

    if (getString(word) != ok)  return status;

    if (strlen(word) != 1)     
    {
	status=end; pos=posBefore; return status;
    }
    if (sscanf(word,"%c",c) != 1) 
    {
	status=end; pos=posBefore; return status;
    }
    status = ok;
    return status;
}

int BufferedParser:: getSpecValue(char* c)
{
    char word[WORDLENGTH];

    posBefore = pos;

    if (getString(word) != ok)  {
	pos=posBefore; return status;
    }

    if (strlen(word) != 1)     
    {
	status=end; pos=posBefore; return status;
    }
    if (sscanf(word,"%c",c) != 1) 
    {
	status=end; pos=posBefore; return status;
    }
    if (   (*c != 'c') && (*c != 's') && (*c != 'p') 
        && (*c != 'C') && (*c != 'S') && (*c != 'P') )
    {
	status=end; pos=posBefore; return status;
    }
    status = ok;
    return status;
}



int BufferedParser:: getValue(int* i)
{
    char word[WORDLENGTH];


    if (getString(word) != ok)  return status;
    if (sscanf(word,"%d",i) != 1) readError("getValue(int)", word);
    status = ok;
    return status;
}

int BufferedParser:: getValue(float* x)
{
    char word[WORDLENGTH];

    if (getString(word) != ok)  return status;
    if (sscanf(word,"%f",x) != 1) readError("getValue(float)", word);
    status = ok;
    return status;
}

int BufferedParser:: getValue(double* x)
{
    char word[WORDLENGTH];

    if (getString(word) != ok)  return status;
    if (sscanf(word,"%lf",x) != 1) readError("getValue(double)", word);
    status = ok;
    return status;
}
//-------------------------------------------------------------------------

int BufferedParser:: getString(char* word)
{
  //    while (sscanf(pos,"%s",word) == 1)		// beware of EOF
  //    {
  //	while (isspace(*pos)) ++pos;		// skip spaces
  //	pos += strlen(word);
  //
  // 	strToLower(word);
  //
  //	if (word[0] == commentFlag)
  //	{
  //	    if (skipRestOfLine() == eof) { status = eof; return status; }
  //	}
  //	else 
  //	{
  //	    if (!strcmp(word,endToken)) status = end;
  //	    else  status = ok;
  //	    return status;
  //	}
  //    }
  //    status = eof;
  //    return status;

  // This code snippet is a performance hack to replace the
  // code above and gains a factor 60 speedup in some, if not
  // many, situations.
  while (1) {
    char* w = word;
    while (isspace(*pos)) ++pos;                  // skip leading whitespace
    while (*pos && !isspace(*pos)) *w++ = *pos++; // copy string
    *w = 0;                                       // terminate copied string

    // copied string empty? Must have been eof... (rather eos here ;-)
    if (0 == w-word) { status = eof; return status; }
    
    
    if (word[0] == commentFlag) { // a comment! skip it.
      if (skipRestOfLine()==eof) { status = eof; return status; }
    } else { // A token! Wow!
      strToLower(word);
      if (!strcmp(word,endToken)) status = end;
      else status = ok;
      return status;
    }
  }
  
  // never get here...
  assert(0);
  return 0;

}
//-------------------------------------------------------------------------

int BufferedParser:: skipRestOfLine()
{
    while (*pos != '\n') 
    {
	if (*pos=='\0') { status = eof; return status; }
	else ++pos;
    }
    ++pos;
    status = ok; 
    return status;
}

char BufferedParser:: getChar()
{
    char c;
    if ((c = *pos) == '\0') status = eof;
    else {
	++pos;
	status = ok; 
    }
    return c;
}

int BufferedParser:: getLine(char* line)
{
    int count = 0;
    while (*pos != '\n')
    {
	if (*pos == '\0') { line[++count] = '\n'; status = eof; return status; }
	line[++count] = *pos;
	++pos;
    }
    line[++count] = '\n';
    ++pos;
    status = ok; 
    return status;
}
//-------------------------------------------------------------------------

void BufferedParser:: rewind()  { pos = buffer; }
long BufferedParser:: ftell()   { return long(pos-buffer); }

int  BufferedParser:: fseek(long offset, int mode) 
{  
    if      (mode == beginning)  pos = buffer + offset;
    else if (mode == current)    pos = pos + offset;
    else if (mode == endOfFile)  pos = buffer + bufferLength + offset;
    else
    {
	cout << "\n*** BufferedParser: fseek called on unknown mode "
	<< mode <<"\n"; 
	cout.flush(); abort();
    }

    if (pos < buffer || pos > buffer+bufferLength)
    {
	cout << "\n*** BufferedParser: fseek requested position outside file "
	<< "\n";
	cout.flush(); abort();
    }
    return 1;
}
//-------------------------------------------------------------------------
//-------------------------------------------------------------------------


// 	invert matrix A (which is destroyed in LUDecompose)

void invert(Matrix<Real>& A, Matrix<Real>& AInv)
{
    int i,j;
    const int dim = A.cHigh();

    if (dim == 1) 
    {
	AInv(1,1) = 1.0/A(1,1); 
    }
    else if (dim == 2)
    {
	AInv(1,1) =  A(2,2);   AInv(1,2) = -A(1,2); 
	AInv(2,1) = -A(2,1);   AInv(2,2) =  A(1,1);
	
	Real det = A(1,1)*A(2,2) - A(1,2)*A(2,1);

	for (i=1; i<=dim; ++i)  
	for (j=1; j<=dim; ++j) AInv(i,j) /= det;
    }
    else
    {
	Vector<int>  index(dim);
	Vector<Real> col(dim);

	Decompose(A, index);
	
	for (j=1; j<=dim; ++j)
	{
	    for (i=1; i<=dim; ++i) col[i] = 0.0;
	    col[j] = 1.0;
	    
	    FBSubst(A, index, col);
	    for (i=1; i<=dim; ++i) AInv(i,j) = col[i];
	}
    }
};
//-------------------------------------------------------------------------

// --LU decomposition with partial pivoting
//			(~ Crout's algorithm taken from Numerical Recipes)


void Decompose(Matrix<Real>& A, Vector<int>& index)
{
    int   i, imax=0, j, k;
    Real  big, dum, sum, temp;
    Real  d = 1.0;
    const Real tiny = machMin(Real(0.0));

    const int dim = index.high();
    Vector<Real> vv(dim);
    
    for (i=1; i<=dim; i++) 
    {
	big=0.0;
	for (j=1; j<=dim; j++)
	{
	    if (Abs(temp=A(i,j)) > Abs(big))  big = temp;
	}
	if (big == 0.0) {
	    cout << "\n*** Decompose: Singular Matrix\n";  cout.flush(); abort();
	}
	vv[i] = 1.0/big;
    }
    
    for (j=1; j<=dim; j++) 
    {
	for (i=1; i<j; i++) 
	{
	    sum = A(i,j);
	    for (k=1; k<i; k++) sum -= A(i,k)*A(k,j);
	    A(i,j) = sum;
	}
	
	big=0.0;

	for (i=j; i<=dim; i++) 
	{
	    sum = A(i,j);
	    for (k=1; k<j; k++) sum -= A(i,k)*A(k,j);
	    A(i,j) = sum;

	    if (Abs(dum=vv[i]*sum) >= Abs(big)) 
	    {
		big=dum;
		imax=i;
	    }
	} 
	
	if (j != imax) 
	{
	    for (k=1; k<=dim; k++) 
	    {
		dum = A(imax,k);
		A(imax,k) = A(j,k);
		A(j,k) = dum;
	    }
	    d = -d;
	    vv[imax]=vv[j];
	}
	
	index[j] = imax;
	if (A(j,j) == 0.0) A(j,j) = tiny;
	if (j != dim) 
	{
	    dum=1.0/A(j,j);
	    for (i=j+1; i<=dim; i++) A(i,j) *= dum;
	}
    }
};
//-------------------------------------------------------------------------

// 	LU-FB-Substitution

// 	solve  lhs = A**-1 * rhs;  A = L*U
//	FB Substitution may be done in place, i.e. FBSubst(x,x) 


void FBSubst(Matrix<Real>& A, Vector<int>& index, Vector<Real>& x)
{
    int i, ii = 0, ip, j;
    Real sum;

    const int dim = index.high();

    for (i=1; i<=dim; i++) 
    {
	ip    = index[i];
	sum   = x[ip];
	x[ip] = x[i];
	if (ii) for (j=ii; j<=i-1; j++) sum -= A(i,j) * x[j];
	else if (Abs(sum)) ii = i;
	x[i] = sum;
    }
    for (i=dim; i>=1; i--) 
    {
	sum = x[i];
	for (j=i+1; j<=dim; j++) sum -= A(i,j) * x[j];
	x[i] = sum / A(i,i);
    }
};
//-------------------------------------------------------------------------
//-------------------------------------------------------------------------


// 	invert matrix A (which is destroyed in LUDecompose)

void invert(Matrix<Complex>& A, Matrix<Complex>& AInv)
{
    int i,j;
    const int dim = A.cHigh();

    if (dim == 1) 
    {
	AInv(1,1) = 1.0/A(1,1); 
    }
    else if (dim == 2)
    {
	AInv(1,1) =  A(2,2);   AInv(1,2) = -A(1,2); 
	AInv(2,1) = -A(2,1);   AInv(2,2) =  A(1,1);
	
	Complex det = A(1,1)*A(2,2) - A(1,2)*A(2,1);

	for (i=1; i<=dim; ++i)  
	for (j=1; j<=dim; ++j) AInv(i,j) /= det;
    }
    else
    {
	Vector<int>  index(dim);
	Vector<Complex> col(dim);

	Decompose(A, index);
	
	for (j=1; j<=dim; ++j)
	{
	    for (i=1; i<=dim; ++i) col[i] = 0.0;
	    col[j] = 1.0;
	    
	    FBSubst(A, index, col);
	    for (i=1; i<=dim; ++i) AInv(i,j) = col[i];
	}
    }
};
//-------------------------------------------------------------------------

// --LU decomposition with partial pivoting 
//			(~ Crout's algorithm taken from Numerical Recipes)


void Decompose(Matrix<Complex>& A, Vector<int>& index)
{
    int   i, imax=0, j, k;
    Complex big, dum, sum, temp;
    Complex d = 1.0;
    const Real tiny = machMin(Real(0.0));

    const int dim = index.high();
    Vector<Complex> vv(dim);
    
    for (i=1; i<=dim; i++) 
    {
	big=0.0;
	for (j=1; j<=dim; j++)
	{
	    if (Abs(temp=A(i,j)) > Abs(big))  big = temp;
	}
	if (big == 0.0) {
	    cout << "\n*** Decompose: Singular Matrix\n";  cout.flush(); abort();
	}
	vv[i] = 1.0/big;
    }
    
    for (j=1; j<=dim; j++) 
    {
	for (i=1; i<j; i++) 
	{
	    sum = A(i,j);
	    for (k=1; k<i; k++) sum -= A(i,k)*A(k,j);
	    A(i,j) = sum;
	}
	
	big=0.0;

	for (i=j; i<=dim; i++) 
	{
	    sum = A(i,j);
	    for (k=1; k<j; k++) sum -= A(i,k)*A(k,j);
	    A(i,j) = sum;

	    if (Abs(dum=vv[i]*sum) >= Abs(big)) 
	    {
		big=dum;
		imax=i;
	    }
	} 
	
	if (j != imax) 
	{
	    for (k=1; k<=dim; k++) 
	    {
		dum = A(imax,k);
		A(imax,k) = A(j,k);
		A(j,k) = dum;
	    }
	    d = -d;
	    vv[imax]=vv[j];
	}
	
	index[j] = imax;
	if (A(j,j) == 0.0) A(j,j) = tiny;
	if (j != dim) 
	{
	    dum=1.0/A(j,j);
	    for (i=j+1; i<=dim; i++) A(i,j) *= dum;
	}
    }
};
//-------------------------------------------------------------------------

// 	LU-FB-Substitution

// 	solve  lhs = A**-1 * rhs;  A = L*U
//	FB Substitution may be done in place, i.e. FBSubst(x,x) 


void FBSubst(Matrix<Complex>& A, Vector<int>& index, Vector<Complex>& x)
{
    int i, ii = 0, ip, j;
    Complex sum;

    const int dim = index.high();

    for (i=1; i<=dim; i++) 
    {
	ip    = index[i];
	sum   = x[ip];
	x[ip] = x[i];
	if (ii) for (j=ii; j<=i-1; j++) sum -= A(i,j) * x[j];
	else if (Abs(sum)) ii = i;
	x[i] = sum;
    }
    for (i=dim; i>=1; i--) 
    {
	sum = x[i];
	for (j=i+1; j<=dim; j++) sum -= A(i,j) * x[j];
	x[i] = sum / A(i,i);
    }
};
//-------------------------------------------------------------------------
//-------------------------------------------------------------------------

static const int stringSize = 2000;
static char String[stringSize];

char *Form(const char *format ...) 
{
    va_list va;
    va_start(va, format);
    if (strlen(format) > 100)
      cout << 
	"\n--- Warning --- Form: standard output string may be too small\n";

    for (int i=0; i<stringSize; ++i) String[i] = '\0'; 
    vsprintf(String, format, va);
    va_end(va);
    return String;
}
//-------------------------------------------------------------------------
//-------------------------------------------------------------------------

#if GETHRVTIME 
/*
 * 	It uses the system call gethrtime(3C), which is accurate to 
 *	nanoseconds. 
*/
#include <sys/time.h>
 
static double z_seconds() {
    return ( gethrvtime() / 1e9 );
}

#elif MACOS

#include <time.h>

static double z_seconds() {
    double t = clock(), cps = CLOCKS_PER_SEC;
    return t/cps;
}

#else

#include <sys/times.h>
#include <unistd.h>

#ifndef CLK_TCK
#define CLK_TCK 60
#endif

static double z_seconds()
{
    struct tms rusage;
    times(&rusage);
    return (double)rusage.tms_utime / CLK_TCK;
}

#endif