#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#include <string.h>
#include <ScalarVectors.h>

/*********************************************************************************/

void CreateInts (int **vint, int dim) {
	if ((*vint = (int *) malloc (sizeof(int)*dim)) == NULL)
		{ printf ("Memory Error (CreateInts(%d))\n", dim); exit (1); }
}

void RemoveInts (int **vint) { 
	if (*vint != NULL) { free (*vint); *vint = NULL; }  
}

void InitInts (int *vint, int dim, int frst, int incr) {
	int i, *p1 = vint, num = frst;

	for (i=0; i<dim; i++) 
		{ *(p1++) = num; num += incr; }
}
		
void CopyInts (int *src, int *dst, int dim) { 
	memmove (dst, src, sizeof(int) * dim);
}

void CopyShiftInts (int *src, int *dst, int dim, int shft) {
	int i, *p1 = src, *p2 = dst;

	if (shft == 0)
		CopyInts (src, dst, dim);
	else
		for (i=0; i<dim; i++)
			*(p2++) = *(p1++) + shft;
}
	
void TransformLengthtoHeader (int *vint, int dim) {
	int i, *pi = vint; 

	for (i=0; i<dim; i++) { *(pi+1) += *pi; pi++; }
}

void TransformHeadertoLength (int *vint, int dim) {
	int i, *pi = vint+dim; 

	for (i=dim; i>0; i--) { *(pi) -= *(pi-1); pi--; }
}

void ComputeHeaderfromLength (int *len, int *head, int dim) {
	int i, *pi1 = len, *pi2 = head; 

	for (i=0; i<dim; i++) { *(pi2+1) = (*pi2) +(*(pi1++)); pi2++; }
}

void ComputeLengthfromHeader (int *head, int *len, int dim) {
	int i, *pi1 = head, *pi2 = len; 

	for (i=0; i<dim; i++) { *(pi2++) = (*(pi1+1)) -(*pi1); pi1++; }
}

int AddInts (int *vint, int dim) {
	int i, *pi = vint, aux = 0;

	for (i=0; i<dim; i++) { 
		aux += *pi; pi++; 
	}

	return aux;
}

// The permutation defined by perm is applied on vec, whose size is dim. 
void PermuteInts (int *vec, int *perm, int dim) {
  int i, *pi = vec;

  for (i=0; i<dim; i++) { *pi = perm[*pi]; pi++; }
}

// Apply the inverse of perm, and store it on iperm, whose size is dim. 
void ComputeInvPermutation (int *perm, int *iperm, int dim) {
  int i, *pi1 = perm;

  for (i=0; i<dim; i++) { iperm[*(pi1++)] = i; }
}

// Scale by scal the elements of vint, whose size is dim. 
void ScaleInts (int *vint, int scal, int dim) {
  int i;
  int *pi = vint;

  for (i=0; i<dim; i++)
    *(pi++) *= scal;
}

void GetIntFromString (char *string, int *pnum, int numC, int shft) {
  int j = 0, num = 0, neg = 0;
  char *pchar = string;

  while ((j < numC) && ((*pchar < '0') || (*pchar > '9')) &&
         (*pchar != '+') && (*pchar != '-')) { j++; pchar++; }
  if (j < numC) {
    if ((*pchar == '+') || (*pchar == '-'))
      { neg = (*pchar == '-'); j++; pchar++; }
    while ((j < numC) && (*pchar >= '0') && (*pchar <= '9'))
      { num = num * 10 + (*pchar - 48); j++; pchar++; }
  }
  if (neg) num = -num;
  *pnum = num + shft;
}

void GetIntsFromString (char *string, int *vec, int numN, int numC, int shft) {
  int i, *pint = vec;
  char *pchar = string;

  for (i=0; i<numN; i++)
    { GetIntFromString (pchar, (pint++), numC, shft); pchar += numC; }
}

void GetFormatsFromString (char *string, int *vec, int numN, int numC) {
  int i, k = 0;
  int *pint = vec;
  char *pchar = string, *pch = NULL, c = ' ', c2;

  for (i=0; i<numN; i++) {
    pch = pchar;
    while (*pch == ' ') pch++;
    sscanf (pch, "(%i%c", pint, &c);
    if ((c == 'P') || (c == 'p')) {
      sscanf (pch, "(%i%c%i%c%i.%i)", &k, &c2, pint, &c, pint+1, pint+2);
      pint += 3;
    }
    else if ((c == 'E') || (c == 'e') || (c == 'D') || (c == 'd') ||
             (c == 'F') || (c == 'f') || (c == 'G') || (c == 'g')) {
      sscanf (pch, "(%i%c%i.%i)", pint, &c, pint+1, pint+2);
      pint += 3;
    }
    else
      { sscanf (pch, "(%i%c%i)", pint, &c, pint+1); pint += 2; }
    pchar += numC;
  }
}

/*********************************************************************************/

void CreateDoubles (double **vdbl, int dim) {
	if ((*vdbl = (double *) malloc (sizeof(double)*dim)) == NULL)
		{ printf ("Memory Error (CreateDoubles(%d))\n", dim); exit (1); }
}

void RemoveDoubles (double **vdbl) { 
	if (*vdbl != NULL) { free (*vdbl); *vdbl = NULL; }
}

void InitDoubles (double *vdbl, int dim, double frst, double incr) {
	int i; 
	double *pd = vdbl, num = frst;

	for (i=0; i<dim; i++) 
		{ *(pd++) = num; num += incr; }
}
		
void InitRandDoubles (double *vdbl, int dim, double frst, double last) {
	int i; 
	double *pd = vdbl, size = last - frst;

	for (i=0; i<dim; i++) 
		{ *(pd++) = frst + (size * (rand() / (RAND_MAX + 1.0))); }
}
		
void CopyDoubles (double *src, double *dst, int dim) { 
	memmove (dst, src, sizeof(double) * dim);
}

void ScaleDoubles (double *vdbl, double scal, int dim) {
	int i; 
	double *pd = vdbl;

	for (i=0; i<dim; i++) 
		*(pd++) *= scal;
}

double DotDoubles (double *vdbl1, double *vdbl2, int dim) {
	int i; 
	double *pd1 = vdbl1, *pd2 = vdbl2, res = 0.0;

	for (i=0; i<dim; i++) 
		res += (*(pd2++)) * (*(pd1++));

	return res;
}

void VvecDoubles (double alfa, double *src1, double *src2, double beta, double *dst, int dim) {
    int i;

    for (i = 0; i < dim; i++) {
        //dst[i] = (beta * dst[i]) + (alfa * src1[i] * src2[i]); 
        double tmp = alfa * src1[i] * src2[i];
        dst[i] = fma(beta, dst[i], tmp);
    }
}

void GetDoubleFromString (char *string, double *pdbl, int numC) {
  int j, k, exp, neg;
  double num, frac;
  char *pchar = string;

  j = 0; exp = 0; neg = 0; num = 0.0; frac = 1.0;
  while ((j < numC) && ((*pchar < '0') || (*pchar > '9')) &&
         (*pchar != '+') && (*pchar != '-') && (*pchar != '.')) { j++; pchar++; }
  if (j < numC) {
    if ((*pchar == '+') || (*pchar == '-'))
      { neg = (*pchar == '-'); j++; pchar++; }
    if (j < numC) {
      if (*pchar != '.')
        while ((j < numC) && (*pchar >= '0') && (*pchar <= '9'))
          { num = num * 10 + (*pchar - 48); j++; pchar++; }
      if (j < numC) {
        if (*pchar == '.') {
          j++; pchar++;
          while ((j < numC) && (*pchar >= '0') && (*pchar <= '9'))
            { frac /= 10; num += (*pchar-48) * frac; j++; pchar++; }
        }
        if (neg) num = -num;
        if (j < numC) {
          if ((*pchar == 'e') || (*pchar == 'E') || (*pchar == 'd') || (*pchar == 'D')) {
            neg = 0; j++; pchar++;
            if (j < numC) {
              if ((*pchar == '+') || (*pchar == '-'))
                { neg = (*pchar == '-'); j++; pchar++; }
              if (j < numC) {
                while ((j < numC) && (*pchar >= '0') && (*pchar <= '9'))
                  { exp = exp*10 + (*pchar-48); j++; pchar++; }
                if (neg) exp = -exp;
                for (k=0; k<exp; k++) num *= 10;
                for (k=0; k>exp; k--) num /= 10;
              }
            }
          }
        }
      } else
        if (neg) num = -num;
    }
  }
  *pdbl = num;
}

void GetDoublesFromString (char *string, double *vec, int numN, int numC) {
  int i;
  double *paux = vec;
  char *pchar = string;

  for (i=0; i<numN; i++)
    { GetDoubleFromString (pchar, (paux++), numC); pchar += numC; }
}

/*********************************************************************************/