|
AbstractLinAlgPack: C++ Interfaces For Vectors, Matrices And Related Linear Algebra Objects Version of the Day
|
00001 /* translated by f2c from dchud.f and hand modified. */ 00002 00003 #include "Moocho_Config.h" 00004 00005 00006 #if !defined(HAVE_MOOCHO_FORTRAN) 00007 00008 00009 #include "Teuchos_BLAS_wrappers.hpp" 00010 #include <math.h> 00011 00012 00013 void dchud_c(double r__[], int *ldr, int *p, 00014 double x[], double z__[], int *ldz, int *nz, 00015 double y[], double *rho, double c__[], double s[] 00016 ) 00017 /* 00018 double *r__; 00019 int *ldr, *p; 00020 double *x, *z__; 00021 int *ldz, *nz; 00022 double *y, *rho, *c__, *s; 00023 */ 00024 { 00025 /* System generated locals */ 00026 int r_dim1, r_offset, z_dim1, z_offset, i__1, i__2; 00027 double d__1, d__2; 00028 00029 /* Local variables */ 00030 double zeta; 00031 int i__, j; 00032 double t, scale, azeta; 00033 double xj; 00034 int jm1; 00035 00036 /* ***FIRST EXECUTABLE STATEMENT DCHUD */ 00037 /* Parameter adjustments */ 00038 r_dim1 = *ldr; 00039 r_offset = 1 + r_dim1 * 1; 00040 r__ -= r_offset; 00041 --x; 00042 z_dim1 = *ldz; 00043 z_offset = 1 + z_dim1 * 1; 00044 z__ -= z_offset; 00045 --y; 00046 --rho; 00047 --c__; 00048 --s; 00049 00050 /* Function Body */ 00051 i__1 = *p; 00052 for (j = 1; j <= i__1; ++j) { 00053 xj = x[j]; 00054 00055 /* APPLY THE PREVIOUS ROTATIONS. */ 00056 00057 jm1 = j - 1; 00058 if (jm1 < 1) { 00059 goto L20; 00060 } 00061 i__2 = jm1; 00062 for (i__ = 1; i__ <= i__2; ++i__) { 00063 t = c__[i__] * r__[i__ + j * r_dim1] + s[i__] * xj; 00064 xj = c__[i__] * xj - s[i__] * r__[i__ + j * r_dim1]; 00065 r__[i__ + j * r_dim1] = t; 00066 /* L10: */ 00067 } 00068 L20: 00069 00070 /* COMPUTE THE NEXT ROTATION. */ 00071 00072 DROTG_F77(&r__[j + j * r_dim1], &xj, &c__[j], &s[j]); 00073 /* L30: */ 00074 } 00075 00076 /* IF REQUIRED, UPDATE Z AND RHO. */ 00077 00078 if (*nz < 1) { 00079 goto L70; 00080 } 00081 i__1 = *nz; 00082 for (j = 1; j <= i__1; ++j) { 00083 zeta = y[j]; 00084 i__2 = *p; 00085 for (i__ = 1; i__ <= i__2; ++i__) { 00086 t = c__[i__] * z__[i__ + j * z_dim1] + s[i__] * zeta; 00087 zeta = c__[i__] * zeta - s[i__] * z__[i__ + j * z_dim1]; 00088 z__[i__ + j * z_dim1] = t; 00089 /* L40: */ 00090 } 00091 azeta = fabs(zeta); 00092 if (azeta == 0. || rho[j] < 0.) { 00093 goto L50; 00094 } 00095 scale = azeta + rho[j]; 00096 /* Computing 2nd power */ 00097 d__1 = azeta / scale; 00098 /* Computing 2nd power */ 00099 d__2 = rho[j] / scale; 00100 rho[j] = scale * sqrt(d__1 * d__1 + d__2 * d__2); 00101 L50: 00102 /* L60: */ 00103 ; 00104 } 00105 L70: 00106 00107 return; 00108 00109 } /* dchud_ */ 00110 00111 00112 #endif // !defined(HAVE_MOOCHO_FORTRAN)
1.7.4