bisection_special_3 Subroutine

private subroutine bisection_special_3(cx, cy, Xo, Yo, Zo, tol, n, m, Hm)

Arguments

Type IntentOptional Attributes Name
real(kind=wp), intent(in) :: cx
real(kind=wp), intent(in) :: cy
real(kind=wp), intent(in) :: Xo
real(kind=wp), intent(in) :: Yo
real(kind=wp), intent(in) :: Zo
real(kind=wp), intent(in) :: tol
integer, intent(out) :: n
real(kind=wp), intent(out) :: m
real(kind=wp), intent(out) :: Hm

Called by

proc~~bisection_special_3~~CalledByGraph proc~bisection_special_3 bisection_special_3 proc~cartesian_to_geodetic_triaxial cartesian_to_geodetic_triaxial proc~cartesian_to_geodetic_triaxial->proc~bisection_special_3

Source Code

subroutine bisection_special_3(cx, cy, Xo, Yo, Zo, tol, n, m, Hm)

    real(wp),intent(in) :: cx, cy, Xo, Yo, Zo, tol
    integer,intent(out) :: n
    real(wp),intent(out) :: m, Hm

    real(wp) :: d1,Hd1,d2,d,MM

    d1 = -one+Zo
    Hd1 = ((cx*Xo*cx*Xo)/((cx+d1)*(cx+d1)))+((cy*Yo*cy*Yo)/((cy+d1)*(cy+d1)))
    d2 = -one+sqrt(cx*Xo*cx*Xo+cy*Yo*cy*Yo+Zo*Zo)
    d = (d2 - d1)/two
    n = 0
    m = -two

    do while (d > tol)
        n = n + 1
        MM = m
        m = d1 + d
        Hm = ((cx*Xo*cx*Xo)/((cx+m)*(cx+m)))+((cy*Yo*cy*Yo)/&
             ((cy+m)*(cy+m)))+((Zo*Zo)/((one+m)**2))-one
        if (MM == m + tol .or. Hm == zero) return
        if (sign(one,Hm) == sign(one,Hd1)) then
            d1 = m
            Hd1 = Hm
        else
            d2 = m
        end if
        d = (d2 - d1)/two
    end do

    n = n + 1
    m = d1 + d
    Hm = ((cx*Xo*cx*Xo)/((cx+m)*(cx+m)))+((cy*Yo*cy*Yo)/&
         ((cy+m)*(cy+m)))+((Zo*Zo)/((one+m)**2))-one

end subroutine bisection_special_3