The computation uses the formulas |x| = sqrt( Re(x)2 + Im(x)2 ) sgn(x) = x / |x| if x /= 0 = 1 if x = 0 c = |a| / sqrt(|a|2 + |b|2) s = sgn(a) * conjg(b) / sqrt(|a|2 + |b|2) When a and b are real and r /= 0, the formulas simplify to r = sgn(a)sqrt(|a|2 + |b|*2) c = a / r s = b / r the same as in DROTG when |a| > |b|. When |b| >= |a|, the sign of c and s will be different from those computed by DROTG if the signs of a and b are not the same.
Type | Intent | Optional | Attributes | Name | ||
---|---|---|---|---|---|---|
complex(kind=dp), | intent(inout) | :: | a | |||
complex(kind=dp), | intent(in) | :: | b | |||
real(kind=dp), | intent(out) | :: | c | |||
complex(kind=dp), | intent(out) | :: | s |
pure subroutine stdlib_zrotg( a, b, c, s ) !! The computation uses the formulas !! |x| = sqrt( Re(x)**2 + Im(x)**2 ) !! sgn(x) = x / |x| if x /= 0 !! = 1 if x = 0 !! c = |a| / sqrt(|a|**2 + |b|**2) !! s = sgn(a) * conjg(b) / sqrt(|a|**2 + |b|**2) !! When a and b are real and r /= 0, the formulas simplify to !! r = sgn(a)*sqrt(|a|**2 + |b|**2) !! c = a / r !! s = b / r !! the same as in DROTG when |a| > |b|. When |b| >= |a|, the !! sign of c and s will be different from those computed by DROTG !! if the signs of a and b are not the same. ! -- reference blas level1 routine -- ! -- reference blas is a software package provided by univ. of tennessee, -- ! -- univ. of california berkeley, univ. of colorado denver and nag ltd..-- ! Scaling Constants ! Scalar Arguments real(dp), intent(out) :: c complex(dp), intent(inout) :: a complex(dp), intent(in) :: b complex(dp), intent(out) :: s ! Local Scalars real(dp) :: d, f1, f2, g1, g2, h2, p, u, uu, v, vv, w complex(dp) :: f, fs, g, gs, r, t ! Intrinsic Functions intrinsic :: abs,aimag,conjg,max,min,real,sqrt ! Statement Functions real(dp) :: abssq ! Statement Function Definitions abssq( t ) = real( t,KIND=dp)**2 + aimag( t )**2 ! Executable Statements f = a g = b if( g == czero ) then c = one s = czero r = f else if( f == czero ) then c = zero g1 = max( abs(real(g,KIND=dp)), abs(aimag(g)) ) if( g1 > rtmin .and. g1 < rtmax ) then ! use unscaled algorithm g2 = abssq( g ) d = sqrt( g2 ) s = conjg( g ) / d r = d else ! use scaled algorithm u = min( safmax, max( safmin, g1 ) ) uu = one / u gs = g*uu g2 = abssq( gs ) d = sqrt( g2 ) s = conjg( gs ) / d r = d*u end if else f1 = max( abs(real(f,KIND=dp)), abs(aimag(f)) ) g1 = max( abs(real(g,KIND=dp)), abs(aimag(g)) ) if( f1 > rtmin .and. f1 < rtmax .and. g1 > rtmin .and. g1 < rtmax ) then ! use unscaled algorithm f2 = abssq( f ) g2 = abssq( g ) h2 = f2 + g2 if( f2 > rtmin .and. h2 < rtmax ) then d = sqrt( f2*h2 ) else d = sqrt( f2 )*sqrt( h2 ) end if p = 1 / d c = f2*p s = conjg( g )*( f*p ) r = f*( h2*p ) else ! use scaled algorithm u = min( safmax, max( safmin, f1, g1 ) ) uu = one / u gs = g*uu g2 = abssq( gs ) if( f1*uu < rtmin ) then ! f is not well-scaled when scaled by g1. ! use a different scaling for f. v = min( safmax, max( safmin, f1 ) ) vv = one / v w = v * uu fs = f*vv f2 = abssq( fs ) h2 = f2*w**2 + g2 else ! otherwise use the same scaling for f and g. w = one fs = f*uu f2 = abssq( fs ) h2 = f2 + g2 end if if( f2 > rtmin .and. h2 < rtmax ) then d = sqrt( f2*h2 ) else d = sqrt( f2 )*sqrt( h2 ) end if p = 1 / d c = ( f2*p )*w s = conjg( gs )*( fs*p ) r = ( fs*( h2*p ) )*u end if end if a = r return end subroutine stdlib_zrotg