/usr/share/doc/libsuperlu-dev/tests/MATGEN/clartg.c is in libsuperlu3-dev 3.0+20070106-3.
This file is owned by root:root, with mode 0o644.
The actual contents of the file can be viewed below.
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 | #include "f2c.h"
/* Subroutine */ int clartg_(complex *f, complex *g, real *cs, complex *sn,
complex *r)
{
/* -- LAPACK auxiliary routine (version 2.0) --
Univ. of Tennessee, Univ. of California Berkeley, NAG Ltd.,
Courant Institute, Argonne National Lab, and Rice University
September 30, 1994
Purpose
=======
CLARTG generates a plane rotation so that
[ CS SN ] [ F ] [ R ]
[ __ ] . [ ] = [ ] where CS**2 + |SN|**2 = 1.
[ -SN CS ] [ G ] [ 0 ]
This is a faster version of the BLAS1 routine CROTG, except for
the following differences:
F and G are unchanged on return.
If G=0, then CS=1 and SN=0.
If F=0 and (G .ne. 0), then CS=0 and SN=1 without doing any
floating point operations.
Arguments
=========
F (input) COMPLEX
The first component of vector to be rotated.
G (input) COMPLEX
The second component of vector to be rotated.
CS (output) REAL
The cosine of the rotation.
SN (output) COMPLEX
The sine of the rotation.
R (output) COMPLEX
The nonzero component of the rotated vector.
=====================================================================
[ 25 or 38 ops for main paths ] */
/* System generated locals */
real r__1, r__2;
doublereal d__1;
complex q__1, q__2, q__3;
/* Builtin functions */
void r_cnjg(complex *, complex *);
double c_abs(complex *), r_imag(complex *), sqrt(doublereal);
/* Local variables */
static real d, f1, f2, g1, g2, fa, ga, di;
static complex fs, gs, ss;
if (g->r == 0.f && g->i == 0.f) {
*cs = 1.f;
sn->r = 0.f, sn->i = 0.f;
r->r = f->r, r->i = f->i;
} else if (f->r == 0.f && f->i == 0.f) {
*cs = 0.f;
r_cnjg(&q__2, g);
d__1 = c_abs(g);
q__1.r = q__2.r / d__1, q__1.i = q__2.i / d__1;
sn->r = q__1.r, sn->i = q__1.i;
d__1 = c_abs(g);
r->r = d__1, r->i = 0.f;
/* SN = ONE
R = G */
} else {
f1 = (r__1 = f->r, dabs(r__1)) + (r__2 = r_imag(f), dabs(r__2));
g1 = (r__1 = g->r, dabs(r__1)) + (r__2 = r_imag(g), dabs(r__2));
if (f1 >= g1) {
q__1.r = g->r / f1, q__1.i = g->i / f1;
gs.r = q__1.r, gs.i = q__1.i;
/* Computing 2nd power */
r__1 = gs.r;
/* Computing 2nd power */
r__2 = r_imag(&gs);
g2 = r__1 * r__1 + r__2 * r__2;
q__1.r = f->r / f1, q__1.i = f->i / f1;
fs.r = q__1.r, fs.i = q__1.i;
/* Computing 2nd power */
r__1 = fs.r;
/* Computing 2nd power */
r__2 = r_imag(&fs);
f2 = r__1 * r__1 + r__2 * r__2;
d = sqrt(g2 / f2 + 1.f);
*cs = 1.f / d;
r_cnjg(&q__3, &gs);
q__2.r = q__3.r * fs.r - q__3.i * fs.i, q__2.i = q__3.r * fs.i +
q__3.i * fs.r;
d__1 = *cs / f2;
q__1.r = d__1 * q__2.r, q__1.i = d__1 * q__2.i;
sn->r = q__1.r, sn->i = q__1.i;
q__1.r = d * f->r, q__1.i = d * f->i;
r->r = q__1.r, r->i = q__1.i;
} else {
q__1.r = f->r / g1, q__1.i = f->i / g1;
fs.r = q__1.r, fs.i = q__1.i;
/* Computing 2nd power */
r__1 = fs.r;
/* Computing 2nd power */
r__2 = r_imag(&fs);
f2 = r__1 * r__1 + r__2 * r__2;
fa = sqrt(f2);
q__1.r = g->r / g1, q__1.i = g->i / g1;
gs.r = q__1.r, gs.i = q__1.i;
/* Computing 2nd power */
r__1 = gs.r;
/* Computing 2nd power */
r__2 = r_imag(&gs);
g2 = r__1 * r__1 + r__2 * r__2;
ga = sqrt(g2);
d = sqrt(f2 / g2 + 1.f);
di = 1.f / d;
*cs = fa / ga * di;
r_cnjg(&q__3, &gs);
q__2.r = q__3.r * fs.r - q__3.i * fs.i, q__2.i = q__3.r * fs.i +
q__3.i * fs.r;
d__1 = fa * ga;
q__1.r = q__2.r / d__1, q__1.i = q__2.i / d__1;
ss.r = q__1.r, ss.i = q__1.i;
q__1.r = di * ss.r, q__1.i = di * ss.i;
sn->r = q__1.r, sn->i = q__1.i;
q__2.r = g->r * ss.r - g->i * ss.i, q__2.i = g->r * ss.i + g->i *
ss.r;
q__1.r = d * q__2.r, q__1.i = d * q__2.i;
r->r = q__1.r, r->i = q__1.i;
}
}
return 0;
/* End of CLARTG */
} /* clartg_ */
|