Main Page
Namespaces
Classes
Files
File List
File Members
contrib
contact
clib
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
}
arith.h
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 Fri Feb 21 2020 03:20:54