snrm2.c
Go to the documentation of this file.
1 /* Euclidean norm of the n-vector stored in sx()
2  with storage -- Simplified Version
3  Ver.1.0, May,26,1988. */
4 
5 #include "arith.h"
6 #define CUTLO 4.441e-16
7 #define CUTHI 1.304e19
8 #define ZERO 0.0e0
9 #define ONE 1.0e0
10 #define SXI sx[isx+i][jsx]
11 
12 REAL snrm2(n,sx,isx,jsx,incx)
13 int n,isx,jsx,incx;
14 MATRIX sx;
15 {
16  int i,j,next,nn;
17  REAL hitest,sum,xmax,snrm;
18 
19  if(n <=0)
20  return(0.0);
21 
22  next=30;
23  sum= ZERO;
24  nn= n*incx;
25 
26  i=0;
27 
28  label20:
29  switch(next){
30  case 30:
31  goto label30;
32  break;
33 
34  case 50:
35  goto label50;
36  break;
37 
38  case 70:
39  goto label70;
40  break;
41 
42  case 110:
43  goto label110;
44  break;
45  }
46 
47  label30:
48  if(fabs(SXI) > CUTLO)
49  goto label85;
50  next=50;
51  xmax=ZERO;
52 
53  /* phase 1. sum is zero */
54 
55  label50:
56  if(SXI == ZERO)
57  goto label200;
58  if(fabs(SXI) > CUTLO)
59  goto label85;
60 
61  /* prepare for pahse 2. */
62 
63  next=70;
64  goto label105;
65 
66  /* prepare for phase 4. */
67 
68  label100:
69  i=j;
70  next=110;
71  sum=(sum/SXI)/SXI;
72 
73  label105:
74  xmax=fabs(SXI);
75  goto label115;
76 
77  /* phase 2. sum is small.
78  scale to avoid destructive underflow. */
79 
80  label70:
81  if(fabs(SXI) > CUTLO)
82  goto label75;
83 
84  /* common code for phases 2 and 4.
85  in phase 4 sum is large. scale to avoid overflow. */
86 
87  label110:
88  if(fabs(SXI) <= xmax)
89  goto label115;
90  sum= ONE + sum*(xmax/SXI)*(xmax/SXI);
91  xmax=fabs(SXI);
92  goto label200;
93 
94  label115:
95  sum= sum+(SXI/xmax)*(SXI/xmax);
96  goto label200;
97 
98  /* prepare for phase 3. */
99 
100  label75:
101  sum= (sum*xmax)*xmax;
102 
103  /* for real or D.P. set hitest= CUTHI/n */
104 
105  label85:
106  hitest= CUTHI/(double)n;
107 
108  /* phase 3. sum is mid-range. no scaling. */
109 
110  for(j=i; j<nn; j += incx){
111  if(fabs(SXI) >= hitest)
112  goto label100;
113  sum= sum+ sx[isx+j][jsx]*sx[isx+j][jsx];
114  }
115 
116  snrm= sqrt(sum);
117  return(snrm);
118 
119  label200:
120  i= i+incx;
121  if(i < nn)
122  goto label20;
123 
124  /* end of main loop
125  computer square root and adjust for scaling. */
126 
127  snrm= xmax*sqrt(sum);
128  return(snrm);
129 }
n
GLfloat n[6][3]
Definition: cube.c:15
REAL
double REAL
Definition: arith.h:25
ZERO
#define ZERO
Definition: snrm2.c:8
SXI
#define SXI
Definition: snrm2.c:10
sqrt
double sqrt()
ONE
#define ONE
Definition: snrm2.c:9
CUTHI
#define CUTHI
Definition: snrm2.c:7
MATRIX
VECTOR MATRIX[MAX]
Definition: arith.h:27
snrm2
REAL snrm2(int n, MATRIX sx, int isx, int jsx, int incx)
Definition: snrm2.c:12
CUTLO
#define CUTLO
Definition: snrm2.c:6

euslisp
Author(s): Toshihiro Matsui
autogenerated on Mon Feb 28 2022 22:18:28