src/gromacs/linearalgebra/gmx_blas/dnrm2.c

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