rotmg(3F)
SROTMG, DROTMG - Constructs a modified Givens plane rotation
As shipped in IRIX 6.5.7. Last changed in IRIX 6.5.5.
NAME SROTMG, DROTMG - Constructs a modified Givens plane rotation SYNOPSIS Real CALL SROTMG (d , d , b , b , rparam) 1 2 1 2 Double precision CALL DROTMG (d , d , b , b , rparam) 1 2 1 2 IMPLEMENTATION IRIX systems DESCRIPTION These routines compute the elements of a modified Givens plane rotation matrix. These routines have the following arguments: d First diagonal element. (input and output) 1 SROTMG: Real. DROTMG: Double precision. On input, this value is the first diagonal element of the scaling matrix D. On the first call to the routine, this value is typically 1.0. Subsequent calls typically use the value from the previous call. On output, this value is the first diagonal element of the updated scaling matrix D'. d First diagonal element. (input and output) 2 SROTMG: Real. DROTMG: Double precision. On input, this is the first diagonal element of the scaling matrix D. On the first call to the routine, this value is typically 1.0. Subsequent calls typically use the value from the previous call. On output, this value is the first diagonal element of the updated scaling matrix D'. b x-coordinate. (input and output) 1 SROTMG: Real. DROTMG: Double precision. On input, this value is the x-coordinate of the vector used to define the angle of rotation, before scaling (multiplying by the matrix D). On output, this value is the x-coordinate of the rotated vector, before scaling (multiplying by the matrix D'). b y-coordinate. (input) 2 SROTMG: Real. DROTMG: Double precision. On input, this value is the y-coordinate of the vector used to define the angle of rotation, before scaling (multiplying by the matrix D). It is unchanged on output. rparam Real array of dimension 5. (output) SROTMG: Real array. DROTMG: Double precision array. This array contains rotation matrix information. The routine sets up the computed elements in rparam from inputs d , d , b , and b . 1 2 1 2 Standard Givens Rotation A standard Givens rotation (see SROTG(3F)) is based on an orthogonal matrix G that rotates points on a Cartesian xy-coordinate plane. To calculate the rotation matrix, you must provide the angle of rotation desired, or, equivalently, a vector (point) that lies along the angle of rotation. For a given planar point (x , y ), G is formed so that: r r _ _ _ _ _ _ _ _ | x' | | c s | | x | G | x | | | | | | r | | r | | 0 | = |-s c | | y | = | y | | | | | | r | | r | - - - - - - - - 2 2 where x' = sqrt (x + y ) . r r With this rotation matrix G, you can then convert any number of existing planar points to the new (rotated) xy-coordinate system. For n points, the rotations would be as follows: _ _ _ _ _ _ | x | | c s | | x | | i | | | | i | | 0 |<- |-s c | | y | | i | | | | i | - - - - - - for i = 1, 2, ..., n Modified Givens Rotation The algorithm for these routines is based on the following observation. The rotation matrix G can be factored into a scaling matrix (diagonal matrix) and modified rotation matrix H, for which either the diagonal or the off-diagonal elements are units (that is, +1 or -1). Thus, to perform m modified (scaled) rotations on n planar points, requires only 2nm, rather than 4nm multiplications for the standard rotation. Because you may want to perform several successive rotations, this routine assumes that you have leftover scaling factors from your previous modified Givens rotation; that is, the routine requires you to input not only a planar rotation vector (b , b ) but also the squares of the diagonal elements of the scali1g m2trix, d and d . The actual rotation vector is specified as follows: 1 2 _ _ _ _ _ _ 1/2_ _ | xr | | sqrt(d1) 0 | | b1 | D | b1 | | yr | = | 0 sqrt(d2) | | b2 | = | b2 | - - - - - - - - where d1 and d2 are the input scaling factors. Given these inputs, the routine generates a new modified rotation matrix H with units for either the diagonal or off-diagonal elements, and new elements d1' and d2' for the new scaling factors, and rotated but unscaled vector (b1', 0), with the following results: _ _ 1/2_ _ 1/2 _ _ 1/2_ _ _ _ G | xr | G D | b1 | D' H | b1 | = D' | h11 h12 | | b1 | | yr | = | b2 | = | b2 | | h21 h22 | | b2 | - - - - - - - - - - 1/2_ _ _ _ D' | b1' | | x' | = | 0 | = | 0 | - - - - where: 1/2 D' ( = diag {sqrt(d1'), sqrt(d2')} ) uses the updated scaling factors d ' and d ', which are d and d on output. 1 2 1 2 H is stored in the output array argument rparam; b1' is stored as b1 on output. 1/2 1/2 D' H equals G D , not G, as implied earlier. You must account for the old scaling factors when calculating the new scaling factors. After calculating the matrix H by using SROTMG/DROTMG, you can then use it in SROTM/DROTM to convert points to the new coordinate system. NOTES Meaning of the Output Values The output values are returned through arguments d , d , b , and rparam. 1 2 1 The scaling factors d and d are updated with each call to the routine. Although SR1TM/DRO2M does not need the updated factors, they are needed in two other important contexts: * As input for subsequent calls to this routine. * As scaling factors for rotated but unscaled points (x , y ), which are output from SROTM./DROTM. i i In this second usage, the actual (scaled) points would be given by (sqrt(d1)*x(i), sqrt(d2)*y(i)). Doing this operation frequently on all your points is counterproductive. The main advantage of the modified rotation algorithm is to reduce the number of operations. If you fold in the scaling factors after each rotation, you are performing the same number of operations as in the standard Givens rotation. These two uses for the scaling factors are mutually exclusive; that is, if you fold the scaling factors back into all your points, you no longer need those factors for these routines. After folding in the scaling factors, and before the next call to SROTMG/DROTMG, reset d and d to 1.0. 1 2 On output, b represents the new x-coordinate after rotating (but before scali1g) the rotation vector. Although the y-coordinate of this vector is 0.0 (see the the previous discussion), the corresponding value b is unchanged on output. 2 The output array argument rparam specifies the format of matrix H, and it holds the nonunit values of H. This is the only output of these routines that SROTM/DROTM requires. Each element of rparam has a specific meaning, as follows: rparam(1) a flag parameter that specifies how the matrix is stored = 0.0 Off-diagonal elements of H are units. = 1.0 Diagonal elements of H are units. = -1.0 Rescaling case (see the following subsection). = -2.0 H is the identity matrix; no rotation needed. rparam(2) = h if needed 1,1 rparam(3) = h if needed 2,1 rparam(4) = h if needed 1,2 rparam(5) = h if needed 2,2 Calculating the Output Values The following presents the algorithm for calculating the output values d , d , b , and rparam, based on the input values d , d , b , and b2. T1is 2lgo1ithm is presented without explanation or 1roo2. 1or a more complete discussion of the modified Givens algorithm, see the papers listed in the SEE ALSO section. Case 1: b = 0 (trivial case) 2 In this case, no rotation is needed. The flag value, rparam(1), is set to -2.0. When passed to routine SROTM/DROTM, this flag indicates that it should not do any rotations. On output from SROTMG/DROTMG, d1, d2, and b1 are unchanged. Case 2: | sqrt(d1)*b1 | > | sqrt(d2)*b2 | (|xr| > |yr|) In this case, the diagonal elements of H (h and h ) are set to 1.0 . Thus, the rparam values set on outpu1,1re as 2o2lows: rparam(1) <-- 0.0 rparam(3) <-- h21 = -b2/b1 rparam(4) <-- h1,2 = (d2*b2)/(d1*b1) The output values of the scaling factors are as follows: d1 <-- d1' = d1/u d2 <-- d2' = d2/u where u = det(H) = 1 + (d2*b2*b2)/(d1*b1*b1) . The output value of b is as follows: 1 b1 <-- b1' = b1*u If rescaling is needed, SROTMG/DROTMG will further modify these output values before the end of the routine. See case 4 later in this subsection. Case 3: | sqrt(d1)*b1 | <= | sqrt(d2)*b2 | (|xr| <= |yr|) In this case, the off-diagonal elements of H are units (to be specific, h2,1 = -1 and h1,2 = 1). Thus, the rparam values set on output are as follows: rparam(1) <-- 1.0 rparam(2) <-- h1,1 = (d1*b1)/(d2*b2) rparam(5) <-- h2,2 = b1/b2 The output values of the scaling factors are as follows: d1 <-- d1' = d2/u d2 <-- d2' = d1/u where u = det(H) = 1 + (d1*b1*b1)/(d2*b2*b2) . The output value of b is as follows: 1 b1 <-- b1' = b2*u If rescaling is needed, SROTMG/DROTMG will further modify these output values before the end of the routine. See case 4. Case 4: Rescaling If the scaling factors become either very large or very small, the scaling and rotation operations may lose a lot of accuracy; therefore, each scaling factor from case 2 or case 3 is kept within the range: gamma**(-2) <= | d(i)' | <= gamma**2 , for i = 1, 2 where d(1)' = d1', d(2)' = d2', and gamma = 2.0**12 = 4096.0. At the end of case 2 or case 3 assignments, if either of the scaling factors falls outside this range, this routine must rescale that factor (and the corresponding elements of H) to bring its size back within the specified range. 48 - log2(gamma) = 48 - 12 = 36 bits (~ 10 decimal digits) Rescaling is performed as follows: If either d1 or d2 is 0, no rescaling is done; otherwise, let q(i) = int(log_base_gamma(sqrt(|d(i)'|))) = int(log2(|d(i)'|)/24) , for i = 1 , 2. Then the following is true: q(i) < 0 , if |d(i)'| < gamma**2 q(i) = 0 , if gamma**(-2) <= |d(i)'| <= gamma**2 q(i) > 0 , if |d(i)'| > gamma**2 Furthermore, |q(i)| represents the number of times d(i)' must be multiplied (or divided) by gamma**2 to return it to the proper range of values. In this case, the rparam values set on output are as follows: rparam(1) <-- -1.0 rparam(2) <-- h1,1' = h1,1*gamma**q(1) rparam(3) <-- h2,1' = h2,1*gamma**q(2) rparam(4) <-- h1,2' = h1,2*gamma**q(1) rparam(5) <-- h2,2' = h2,2*gamma**q(2) The output values of the scaling factors are as follows: d1 <-- d1'' = d1'*gamma**(-2*q(1)) d2 <-- d2'' = d2'*gamma**(-2*q(2)) The output value of b1 is as follows: b1 <-- b1'' = b1'*gamma**q(1) SEE ALSO ROT(3F), ROTG(3F), ROTM(3F) Gentleman, W. M., "Least Squares Computations by Givens Transformations Without Square Roots," Journal of the Institute for Mathematical Applications 12 (1973), pp. 329 - 336. Lawson, C., Hanson, R., Kincaid, D., and Krogh, F., "Basic Linear Algebra Subprograms for Fortran Usage," ACM Transactions on Mathematical Software, 5 (1979), pp. 308 - 325. This man page is available only online.