kepler Function

private pure function kepler(ma, ec)

solve kepler's equation ma = ea + ec*sin(ea)

Arguments

Type IntentOptional Attributes Name
real(kind=wp), intent(in) :: ma

eccentricity
real(kind=wp), intent(in) :: ec

mean anomaly in rad

Return Value real(kind=wp)


Called by

proc~~kepler~~CalledByGraph proc~kepler standish_module::kepler proc~calcelements standish_module::calcelements proc~calcelements->proc~kepler proc~helio standish_module::helio proc~helio->proc~calcelements proc~standish_module_test standish_module::standish_module_test proc~standish_module_test->proc~helio proc~standish_rv_func standish_module::standish_ephemeris%standish_rv_func proc~standish_module_test->proc~standish_rv_func proc~standish_rv_func->proc~helio

Source Code

    pure real(wp) function kepler (ma, ec)

    implicit none

    real(wp), intent (in) :: ma !! eccentricity
    real(wp), intent (in) :: ec !! mean anomaly in rad

    real(wp) :: r, ea
    integer :: i

    ! acceptable accuracy for this calculation
    real(wp),parameter :: tol = 1.0e-08_wp  !! max error in eccentric anomaly `ea` in rad
    integer,parameter :: maxit = 12         !! max iterations (1-4 typical for `ec<0.3`)

    ea = ma + ec * sin (ma)! starting value
    do i = 1, maxit ! newton(-raphson) iterations
        r = (ma-ea+ec*sin(ea)) / (one-ec*cos(ea))
        ea = ea + r
        if (abs(r) <= tol) exit
    end do
    kepler = modulo (ea, twopi) ! eccentric anomaly adjusted 0-2pi

    end function kepler