Compute the coordinates of the libration points (L1,L2,L3,L4,L5).
This is just an alternate version of compute_libration_points.
| Type | Intent | Optional | Attributes | Name | ||
|---|---|---|---|---|---|---|
| real(kind=wp), | intent(in) | :: | mu |
CRTBP parameter |
||
| real(kind=wp), | intent(out), | optional | :: | r1 |
L1 x coordinate |
|
| real(kind=wp), | intent(out), | optional | :: | r2 |
L2 x coordinate |
|
| real(kind=wp), | intent(out), | optional | :: | r3 |
L3 x coordinate |
|
| real(kind=wp), | intent(out), | optional, | dimension(2) | :: | r4 |
L4 [x,y] coordinates |
| real(kind=wp), | intent(out), | optional, | dimension(2) | :: | r5 |
L5 [x,y] coordinates |
subroutine compute_libration_points_v2(mu,r1,r2,r3,r4,r5) use math_module, only: cube_root implicit none real(wp),intent(in) :: mu !! CRTBP parameter real(wp),intent(out),optional :: r1 !! L1 x coordinate real(wp),intent(out),optional :: r2 !! L2 x coordinate real(wp),intent(out),optional :: r3 !! L3 x coordinate real(wp),dimension(2),intent(out),optional :: r4 !! L4 [x,y] coordinates real(wp),dimension(2),intent(out),optional :: r5 !! L5 [x,y] coordinates integer,parameter :: maxiter = 100 !! maximum number of iterations real(wp),parameter :: tol = 1.0e-12_wp !! convergence tolerance real(wp) :: gamma, gamma0 integer :: i !! counter if (present(r1)) then gamma0 = cube_root(mu * (one-mu) / three) gamma = gamma0 + one do i = 1,maxiter if (abs(gamma-gamma0)<=tol) exit gamma0 = gamma gamma = cube_root((mu*(gamma0-one)**2)/(three-two*mu-gamma0*(three-mu-gamma0))) end do r1 = one - mu - gamma end if if (present(r2)) then gamma0 = cube_root(mu * (one-mu) / three) gamma = gamma0 + one do i = 1,maxiter if (abs(gamma-gamma0)<=tol) exit gamma0 = gamma gamma = cube_root((mu*(gamma0+one)**2)/(three-two*mu+gamma0*(three-mu+gamma0))) end do r2 = one - mu + gamma end if if (present(r3)) then gamma0 = cube_root(mu * (one-mu) / three) gamma = gamma0 + one do i = 1,maxiter if (abs(gamma-gamma0)<=tol) exit gamma0 = gamma gamma = cube_root((one-mu)*(gamma0+one)**2/(one+two*mu+gamma0*(two+mu+gamma0))) end do r3 = - mu - gamma end if call compute_libration_points(mu,r4=r4,r5=r5) end subroutine compute_libration_points_v2