anchor.diff:  
        end do newton
   
  +     !----------------------------------------------------------------!
        ! print the summary of the newton iteration
  +     !----------------------------------------------------------------!
        write(*,'(f15.2,e15.6,e15.6,e15.6,e15.6,e15.6,xi4,xi8)') &
             &time, update_norm, res_norm, body1%KE, body1%PE, &
             & body1%KE + body1%PE, newton_cnt, fcnt
   
  -     
  +!!$     !----------------------------------------------------------------!     
  +!!$     ! extrapolate to next time step
  +!!$     !-----------------------------------------------------------------!
  +!!$     dq     = get_approx_q(q,q_dot)
  +!!$     q_dot  = get_approx_q_dot(q,1)
  +
     end do time_march
   
     call print_body(body1)
  
positive.code: 
    ! The following defaults are already set in global_varaibles.f90.
    ! The user is free to change here too.
    dT         = 0.01_dp
    start_time = 0.0_dp
    end_time   = 1.0_dp
  
  
negative.code:         !dq = matmul(jac, -res)
          dq = direct_solve(jac, -res, TOT_NDOF)
  
          !-------------------------------------------------------------!
          ! Calcualte residual tolratances                    
          !-------------------------------------------------------------!
  
          update_norm  =  maxval(abs(dq)); res_norm  =  maxval(abs(res));
  
          !-------------------------------------------------------------!
          ! Update the state variables and then update the body
          !-------------------------------------------------------------!
  
          q            = get_updated_q(q, dq)
          q_dot        = get_updated_q_dot(q_dot, dq)
          ! q_double_dot = get_updated_q_double_dot(q_double_dot, dq) 
  