Mercurial > hg > octave-nkf
annotate libcruft/lapack/slasd5.f @ 9665:1dba57e9d08d
use blas_trans_type for xgemm
author | Jaroslav Hajek <highegg@gmail.com> |
---|---|
date | Sat, 26 Sep 2009 10:41:07 +0200 |
parents | 82be108cc558 |
children |
rev | line source |
---|---|
7789
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
1 SUBROUTINE SLASD5( I, D, Z, DELTA, RHO, DSIGMA, WORK ) |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
2 * |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
3 * -- LAPACK auxiliary routine (version 3.1) -- |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
4 * Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd.. |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
5 * November 2006 |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
6 * |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
7 * .. Scalar Arguments .. |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
8 INTEGER I |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
9 REAL DSIGMA, RHO |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
10 * .. |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
11 * .. Array Arguments .. |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
12 REAL D( 2 ), DELTA( 2 ), WORK( 2 ), Z( 2 ) |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
13 * .. |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
14 * |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
15 * Purpose |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
16 * ======= |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
17 * |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
18 * This subroutine computes the square root of the I-th eigenvalue |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
19 * of a positive symmetric rank-one modification of a 2-by-2 diagonal |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
20 * matrix |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
21 * |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
22 * diag( D ) * diag( D ) + RHO * Z * transpose(Z) . |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
23 * |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
24 * The diagonal entries in the array D are assumed to satisfy |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
25 * |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
26 * 0 <= D(i) < D(j) for i < j . |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
27 * |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
28 * We also assume RHO > 0 and that the Euclidean norm of the vector |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
29 * Z is one. |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
30 * |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
31 * Arguments |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
32 * ========= |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
33 * |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
34 * I (input) INTEGER |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
35 * The index of the eigenvalue to be computed. I = 1 or I = 2. |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
36 * |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
37 * D (input) REAL array, dimension (2) |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
38 * The original eigenvalues. We assume 0 <= D(1) < D(2). |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
39 * |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
40 * Z (input) REAL array, dimension (2) |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
41 * The components of the updating vector. |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
42 * |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
43 * DELTA (output) REAL array, dimension (2) |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
44 * Contains (D(j) - sigma_I) in its j-th component. |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
45 * The vector DELTA contains the information necessary |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
46 * to construct the eigenvectors. |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
47 * |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
48 * RHO (input) REAL |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
49 * The scalar in the symmetric updating formula. |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
50 * |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
51 * DSIGMA (output) REAL |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
52 * The computed sigma_I, the I-th updated eigenvalue. |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
53 * |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
54 * WORK (workspace) REAL array, dimension (2) |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
55 * WORK contains (D(j) + sigma_I) in its j-th component. |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
56 * |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
57 * Further Details |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
58 * =============== |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
59 * |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
60 * Based on contributions by |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
61 * Ren-Cang Li, Computer Science Division, University of California |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
62 * at Berkeley, USA |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
63 * |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
64 * ===================================================================== |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
65 * |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
66 * .. Parameters .. |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
67 REAL ZERO, ONE, TWO, THREE, FOUR |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
68 PARAMETER ( ZERO = 0.0E+0, ONE = 1.0E+0, TWO = 2.0E+0, |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
69 $ THREE = 3.0E+0, FOUR = 4.0E+0 ) |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
70 * .. |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
71 * .. Local Scalars .. |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
72 REAL B, C, DEL, DELSQ, TAU, W |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
73 * .. |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
74 * .. Intrinsic Functions .. |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
75 INTRINSIC ABS, SQRT |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
76 * .. |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
77 * .. Executable Statements .. |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
78 * |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
79 DEL = D( 2 ) - D( 1 ) |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
80 DELSQ = DEL*( D( 2 )+D( 1 ) ) |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
81 IF( I.EQ.1 ) THEN |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
82 W = ONE + FOUR*RHO*( Z( 2 )*Z( 2 ) / ( D( 1 )+THREE*D( 2 ) )- |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
83 $ Z( 1 )*Z( 1 ) / ( THREE*D( 1 )+D( 2 ) ) ) / DEL |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
84 IF( W.GT.ZERO ) THEN |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
85 B = DELSQ + RHO*( Z( 1 )*Z( 1 )+Z( 2 )*Z( 2 ) ) |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
86 C = RHO*Z( 1 )*Z( 1 )*DELSQ |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
87 * |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
88 * B > ZERO, always |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
89 * |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
90 * The following TAU is DSIGMA * DSIGMA - D( 1 ) * D( 1 ) |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
91 * |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
92 TAU = TWO*C / ( B+SQRT( ABS( B*B-FOUR*C ) ) ) |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
93 * |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
94 * The following TAU is DSIGMA - D( 1 ) |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
95 * |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
96 TAU = TAU / ( D( 1 )+SQRT( D( 1 )*D( 1 )+TAU ) ) |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
97 DSIGMA = D( 1 ) + TAU |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
98 DELTA( 1 ) = -TAU |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
99 DELTA( 2 ) = DEL - TAU |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
100 WORK( 1 ) = TWO*D( 1 ) + TAU |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
101 WORK( 2 ) = ( D( 1 )+TAU ) + D( 2 ) |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
102 * DELTA( 1 ) = -Z( 1 ) / TAU |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
103 * DELTA( 2 ) = Z( 2 ) / ( DEL-TAU ) |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
104 ELSE |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
105 B = -DELSQ + RHO*( Z( 1 )*Z( 1 )+Z( 2 )*Z( 2 ) ) |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
106 C = RHO*Z( 2 )*Z( 2 )*DELSQ |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
107 * |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
108 * The following TAU is DSIGMA * DSIGMA - D( 2 ) * D( 2 ) |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
109 * |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
110 IF( B.GT.ZERO ) THEN |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
111 TAU = -TWO*C / ( B+SQRT( B*B+FOUR*C ) ) |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
112 ELSE |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
113 TAU = ( B-SQRT( B*B+FOUR*C ) ) / TWO |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
114 END IF |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
115 * |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
116 * The following TAU is DSIGMA - D( 2 ) |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
117 * |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
118 TAU = TAU / ( D( 2 )+SQRT( ABS( D( 2 )*D( 2 )+TAU ) ) ) |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
119 DSIGMA = D( 2 ) + TAU |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
120 DELTA( 1 ) = -( DEL+TAU ) |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
121 DELTA( 2 ) = -TAU |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
122 WORK( 1 ) = D( 1 ) + TAU + D( 2 ) |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
123 WORK( 2 ) = TWO*D( 2 ) + TAU |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
124 * DELTA( 1 ) = -Z( 1 ) / ( DEL+TAU ) |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
125 * DELTA( 2 ) = -Z( 2 ) / TAU |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
126 END IF |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
127 * TEMP = SQRT( DELTA( 1 )*DELTA( 1 )+DELTA( 2 )*DELTA( 2 ) ) |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
128 * DELTA( 1 ) = DELTA( 1 ) / TEMP |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
129 * DELTA( 2 ) = DELTA( 2 ) / TEMP |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
130 ELSE |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
131 * |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
132 * Now I=2 |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
133 * |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
134 B = -DELSQ + RHO*( Z( 1 )*Z( 1 )+Z( 2 )*Z( 2 ) ) |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
135 C = RHO*Z( 2 )*Z( 2 )*DELSQ |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
136 * |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
137 * The following TAU is DSIGMA * DSIGMA - D( 2 ) * D( 2 ) |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
138 * |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
139 IF( B.GT.ZERO ) THEN |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
140 TAU = ( B+SQRT( B*B+FOUR*C ) ) / TWO |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
141 ELSE |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
142 TAU = TWO*C / ( -B+SQRT( B*B+FOUR*C ) ) |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
143 END IF |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
144 * |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
145 * The following TAU is DSIGMA - D( 2 ) |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
146 * |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
147 TAU = TAU / ( D( 2 )+SQRT( D( 2 )*D( 2 )+TAU ) ) |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
148 DSIGMA = D( 2 ) + TAU |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
149 DELTA( 1 ) = -( DEL+TAU ) |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
150 DELTA( 2 ) = -TAU |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
151 WORK( 1 ) = D( 1 ) + TAU + D( 2 ) |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
152 WORK( 2 ) = TWO*D( 2 ) + TAU |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
153 * DELTA( 1 ) = -Z( 1 ) / ( DEL+TAU ) |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
154 * DELTA( 2 ) = -Z( 2 ) / TAU |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
155 * TEMP = SQRT( DELTA( 1 )*DELTA( 1 )+DELTA( 2 )*DELTA( 2 ) ) |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
156 * DELTA( 1 ) = DELTA( 1 ) / TEMP |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
157 * DELTA( 2 ) = DELTA( 2 ) / TEMP |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
158 END IF |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
159 RETURN |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
160 * |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
161 * End of SLASD5 |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
162 * |
82be108cc558
First attempt at single precision tyeps
David Bateman <dbateman@free.fr>
parents:
diff
changeset
|
163 END |