src/gromacs/linearalgebra/gmx_blas/snrm2.c

   1 #include <math.h>
   2
   3
   4 #include "gromacs/utility/real.h"
   5 #include "../gmx_blas.h"
   6
   7 float
   8 F77_FUNC(snrm2,SNRM2)(int  *     n__,
   9                       float *    x,
  10                       int    *    incx__)
  11 {
  12     int ix,max_ix;
  13     float ssq,scale,absxi,t;
  14
  15     int n = *n__;
  16     int incx = *incx__;
  17
  18     if(n<1 || incx<1)
  19         return 0;
  20     else if (n==1) {
  21         t = x[0];
  22         if(t>=0)
  23             return t;
  24         else
  25             return -t;
  26     }
  27
  28     scale = 0.0;
  29     ssq   = 1.0;
  30
  31     max_ix = 1+(n-1)*(incx);
  32     for(ix=1;ix<=max_ix;ix+=incx) {
  33         t = x[ix-1];
  34         if(fabs(t)>GMX_FLOAT_MIN) {
  35             absxi = (t>=0) ? t : (-t);
  36             if(scale<absxi) {
  37                 t = scale/absxi;
  38                 t = t*t;
  39                 ssq = ssq*t + 1.0;
  40                 scale = absxi;
  41             } else {
  42                 t = absxi/scale;
  43                 ssq += t*t;
  44             }
  45         }
  46     }
  47     return scale*sqrt(ssq);
  48
  49 }
  50
  51
  52