src/gromacs/linearalgebra/gmx_blas/dger.c

   1 #include <math.h>
   2
   3 #include "gromacs/utility/real.h"
   4
   5 #include "../gmx_blas.h"
   6
   7 void
   8 F77_FUNC(dger,DGER)(int *m__,
   9                     int *n__,
  10                     double *alpha__,
  11                     double *x,
  12                     int *incx__,
  13                     double *y,
  14                     int *incy__,
  15                     double *a,
  16                     int *lda__)
  17 {
  18     int ix,kx,jy;
  19     int i,j;
  20     double temp;
  21
  22
  23     int m = *m__;
  24     int n = *n__;
  25     int incx = *incx__;
  26     int incy = *incy__;
  27     int lda = *lda__;
  28     double alpha = *alpha__;
  29
  30     if(m<=0 || n<=0 || fabs(alpha)<GMX_DOUBLE_MIN)
  31         return;
  32
  33     if(incy>0)
  34         jy = 0;
  35     else
  36         jy = incy * (1 - n);
  37
  38     if(incx==1) {
  39         for(j=0;j<n;j++,jy+=incy)
  40             if(fabs(y[jy])>GMX_DOUBLE_MIN) {
  41                 temp = alpha * y[jy];
  42                 for(i=0;i<m;i++)
  43                     a[j*(lda)+i] += temp*x[i];
  44             }
  45     } else {
  46         /* non-unit incx */
  47         if(incx>0)
  48             kx = 0;
  49         else
  50             kx = incx * (1 - m);
  51
  52         for(j=0;j<n;j++,jy+=incy) {
  53             if(fabs(y[jy])>GMX_DOUBLE_MIN) {
  54                 temp = alpha * y[jy];
  55                 ix = kx;
  56                 for(i=0;i<m;i++,ix+=incx)
  57                     a[j*(lda)+i] += temp*x[ix];
  58             }
  59         }
  60     }
  61         return;
  62 }