crtbp_derivs Subroutine

    subroutine crtbp_derivs(mu,x,dx)

    implicit none

    real(wp),intent(in)               :: mu   !! CRTBP parameter (See [[compute_crtpb_parameter]])
    real(wp),dimension(6),intent(in)  :: x    !! normalized state \([\mathbf{r},\mathbf{v}]\)
    real(wp),dimension(6),intent(out) :: dx   !! normalized state derivative \([\dot{\mathbf{r}},\dot{\mathbf{v}}]\)

    !local variables:
    real(wp),dimension(3) :: r1,r2,rb1,rb2,r,v,g
    real(wp) :: r13,r23,omm,c1,c2

    !extract variables from x vector:
    r = x(1:3) ! position
    v = x(4:6) ! velocity

    !other parameters:
    omm = one - mu
    rb1 = [-mu,zero,zero] ! location of body 1
    rb2 = [omm,zero,zero] ! location of body 2
    r1  = r - rb1         ! body1 -> sc vector
    r2  = r - rb2         ! body2 -> sc vector
    r13 = norm2(r1)**3
    r23 = norm2(r2)**3
    c1  = omm/r13
    c2  = mu/r23

    !normalized gravity from both bodies:
    g(1) = -c1*(r(1) + mu) - c2*(r(1)-one+mu)
    g(2) = -c1*r(2)        - c2*r(2)
    g(3) = -c1*r(3)        - c2*r(3)

    ! derivative of x:
    dx(1:3) = v                       ! rdot
    dx(4)   =  two*v(2) + r(1) + g(1) ! vdot
    dx(5)   = -two*v(1) + r(2) + g(2) !
    dx(6)   =                    g(3) !

    end subroutine crtbp_derivs