anchor.diff:     ALLOCATE(Nrrr(rot_elem),STAT = allo_stat)
      IF(allo_stat/=0) GOTO 9999
      Nrrr = 0.0D0
  -   
  +
      ALLOCATE(Nvvv(dof_elem),STAT = allo_stat)
      IF(allo_stat/=0) GOTO 9999
      Nvvv = 0.0D0
  
positive.code:    INTEGER(IntKi),    INTENT(IN   ):: ngp
  
     REAL(ReKi),        INTENT(  OUT):: StifK(:,:)
     REAL(ReKi),        INTENT(  OUT):: RHS(:) 
     REAL(ReKi),        INTENT(  OUT):: MassM(:,:)
     REAL(ReKi),        INTENT(  OUT):: DampG(:,:)
  
  
negative.code:    END SUBROUTINE BD_CrvMatrixH
  
  !**********************************************************************************************************************************
  !   This subroutine calculates the rotation parameters for each element using the Wiener-Milenkovic method.
  !   This method is detailed in the paper: Bauchau, O.A., 2008, "Interpolation of finite rotations in flexible
  !   multi-body dynamics simulations", IMechE, Equation (9).
  !**********************************************************************************************************************************
     SUBROUTINE BD_CrvCompose(rr,pp,qq,flag)
  
  