mymkl.c

#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#include "mymkl.h"
/**********************************************/

void rcopy (int *n, double *x, int *incx, double *y, int *incy) {
  int i, dim = *n, ix = *incx, iy = *incy;
  double *px = x, *py = y;

  for (i=0; i<dim; i++) {
    *py = *px; px += ix; py += iy;
  }
}

void rscal (int *n, double *alpha, double *x, int *incx) {
  int i, dim = *n, ix = *incx;
  double *px = x, a = *alpha;

  for (i=0; i<dim; i++) {
    *px *= a; px += ix;
  }
}

void raxpy (int *n, double *alpha, double *x, int *incx, double *y, int *incy) {
  int i, dim = *n, ix = *incx, iy = *incy;
  double *px = x, *py = y, a = *alpha;

  for (i=0; i<dim; i++) {
    *py += *px * a; px += ix; py += iy;
  }
}

double rdot (int *n, double *x, int *incx, double *y, int *incy) {
  int i, dim = *n, ix = *incx, iy = *incy;
  double aux = 0.0, *px = x, *py = y;

  for (i=0; i<dim; i++) {
    aux += *py * *px; px += ix; py += iy;
  }

  return aux;
}

double rnrm2(int *n, double *x, int *incx) {
  int i, dim = *n, ix = *incx;
  double *px = x, sum = 0.0;

  for (i = 0; i < dim; i++) {
    sum += *px * *px; px += ix;
  }
  return sqrt(sum);
}

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