anchor.diff:        IF ( theta > 3.1_ReKi ) THEN  ! theta/(2*sin(theta)) blows up quickly as theta approaches pi, 
            ! so I'm putting a pretty large tolerance on pi here, and using a different equation to find the solution near pi
                        
  -         !d11 = 1 - (1-cos(theta))/theta^2 * (logMap3^2 + logMap2^2)
  -         !d22 = 1 - (1-cos(theta))/theta^2 * (logMap3^2 + logMap1^2)
  -         !d33 = 1 - (1-cos(theta))/theta^2 * (logMap2^2 + logMap1^2)
  -         
  -         logMap(1) = theta * sqrt(abs( 0.5_ReKi * ( 1.0_ReKi + DCM(1,1) - DCM(2,2) - DCM(3,3) ) / (1.0_ReKi-cosTheta) ))
  -         logMap(2) = theta * sqrt(abs( 0.5_ReKi * ( 1.0_ReKi - DCM(1,1) + DCM(2,2) - DCM(3,3) ) / (1.0_ReKi-cosTheta) ))
  -         logMap(3) = theta * sqrt(abs( 0.5_ReKi * ( 1.0_ReKi - DCM(1,1) - DCM(2,2) + DCM(3,3) ) / (1.0_ReKi-cosTheta) ))
  -               
  -         ! we choose logMap1 positive then we get the signs for logMap2 and logMap3:
  -         if ( .not. EqualRealNos( logMap(1), 0.0_ReKi ) ) then
  -            !d12+d21=2*(1-cos(theta))/theta**2 * logMap(1)*logMap(2); 2*(1-cos(theta))/theta**2 * logMap(1)>0 so logMap(2) is sign(logMap(2),d12+d21)
  -            !d13+d31=2*(1-cos(theta))/theta**2 * logMap(1)*logMap(3); 2*(1-cos(theta))/theta**2 * logMap(1)>0 so logMap(3) is sign(logMap(3),d13+d31)
  -            
  -            logMap(2) = sign( logMap(2), DCM(1,2)+DCM(2,1) )
  -            logMap(3) = sign( logMap(3), DCM(1,3)+DCM(3,1) )            
  +         logMap(1) = 1.0_ReKi + DCM(1,1) - DCM(2,2) - DCM(3,3);
  +         logMap(2) = 1.0_ReKi - DCM(1,1) + DCM(2,2) - DCM(3,3);
  +         logMap(3) = 1.0_ReKi - DCM(1,1) - DCM(2,2) + DCM(3,3);
  +             
  +         indx_max = maxloc( abs(logMap), 1 )
  +             
  +         divisor = sqrt(abs( logMap(indx_max) *  2.0_ReKi*(1.0_ReKi - cosTheta)  )) / theta  ! 2*(1-cosTheta)/theta^2 * abs(lambda(indx))
  +         if (indx_max == 1) then
  +           !logMap(1) = 1.0 + DCM(1,1) - DCM(2,2) - DCM(3,3)               ! 2*(1-cosTheta)/theta^2 * lambda(1) * lambda(1)
  +            logMap(2) = DCM(1,2) + DCM(2,1)                                ! 2*(1-cosTheta)/theta^2 * lambda(1) * lambda(2)
  +            logMap(3) = DCM(1,3) + DCM(3,1)                                ! 2*(1-cosTheta)/theta^2 * lambda(1) * lambda(3)
  +         elseif (indx_max == 2) then
  +            logMap(1) = DCM(1,2) + DCM(2,1)                                ! 2*(1-cosTheta)/theta^2 * lambda(2) * lambda(1)
  +           !logMap(2) = 1.0 - DCM(1,1) + DCM(2,2) - DCM(3,3)               ! 2*(1-cosTheta)/theta^2 * lambda(2) * lambda(2)
  +            logMap(3) = DCM(2,3) + DCM(3,2)                                ! 2*(1-cosTheta)/theta^2 * lambda(2) * lambda(3)
            else
  -            ! because logMap1 is zero, we can choose logMap2 positive:
  -            
  -            !d23+d32=2*(1-cos(theta))/theta**2 * logMap(2)*logMap(3); 2*(1-cos(theta))/theta**2 * logMap(2)>0 so logMap(3) is sign(logMap(3),d23+d32)
  -            logMap(3) = sign( logMap(3), DCM(2,3)+DCM(3,2) )            
  -            
  +            logMap(1) = DCM(1,3) + DCM(3,1)                                ! 2*(1-cosTheta)/theta^2 * lambda(3) * lambda(1)
  +            logMap(2) = DCM(2,3) + DCM(3,2)                                ! 2*(1-cosTheta)/theta^2 * lambda(3) * lambda(2)
  +           !logMap(3) = 1.0 - DCM(1,1) - DCM(2,2) + DCM(3,3)               ! 2*(1-cosTheta)/theta^2 * lambda(3) * lambda(3)
            end if
  -                    
  +         logMap = logMap / divisor                                         ! lambda * sign(lambda(indx_max))
  +         
            ! at this point we may have the wrong sign for logMap (though if theta==pi, it doesn't matter because we can change it in the DCM_setLogMapforInterp() routines)
            ! we'll do a little checking to see if we should change the sign:
            
  
positive.code:          thetaOut = theta
        END IF      
        
        
        !IF ( EqualRealNos( pi_D, theta )  ) THEN
        IF ( theta > 3.1_DbKi ) THEN  ! theta/(2*sin(theta)) blows up quickly as theta approaches pi, 
           ! so I'm putting a pretty large tolerance on pi here, and using a different equation to find the solution near pi
                       
           !d11 = 1 - (1-cos(theta))/theta^2 * (logMap3^2 + logMap2^2)
           !d22 = 1 - (1-cos(theta))/theta^2 * (logMap3^2 + logMap1^2)
           !d33 = 1 - (1-cos(theta))/theta^2 * (logMap2^2 + logMap1^2)
                    
           logMap(1) = theta * sqrt(abs( 0.5_DbKi * ( 1.0_DbKi + DCM(1,1) - DCM(2,2) - DCM(3,3) ) / (1.0_DbKi-cosTheta) ))
           logMap(2) = theta * sqrt(abs( 0.5_DbKi * ( 1.0_DbKi - DCM(1,1) + DCM(2,2) - DCM(3,3) ) / (1.0_DbKi-cosTheta) ))
           logMap(3) = theta * sqrt(abs( 0.5_DbKi * ( 1.0_DbKi - DCM(1,1) - DCM(2,2) + DCM(3,3) ) / (1.0_DbKi-cosTheta) ))
                 
           ! we choose logMap1 positive then we get the signs for logMap2 and logMap3:
           if ( .not. EqualRealNos( logMap(1), 0.0_DbKi ) ) then
              !d12+d21=2*(1-cos(theta))/theta**2 * logMap(1)*logMap(2); 2*(1-cos(theta))/theta**2 * logMap(1)>0 so logMap(2) is sign(logMap(2),d12+d21)
              !d13+d31=2*(1-cos(theta))/theta**2 * logMap(1)*logMap(3); 2*(1-cos(theta))/theta**2 * logMap(1)>0 so logMap(3) is sign(logMap(3),d13+d31)
              
              logMap(2) = sign( logMap(2), DCM(1,2)+DCM(2,1) )
              logMap(3) = sign( logMap(3), DCM(1,3)+DCM(3,1) )            
           else
              ! because logMap1 is zero, we can choose logMap2 positive:
              
              !d23+d32=2*(1-cos(theta))/theta**2 * logMap(2)*logMap(3); 2*(1-cos(theta))/theta**2 * logMap(2)>0 so logMap(3) is sign(logMap(3),d23+d32)
              logMap(3) = sign( logMap(3), DCM(2,3)+DCM(3,2) )            
              
           end if
                      
           ! at this point we may have the wrong sign for logMap (though if theta==pi, it doesn't matter because we can change it in the DCM_setLogMapforInterp() routines)
           ! we'll do a little checking to see if we should change the sign:
           
           IF ( EqualRealNos( pi_D, theta )  ) RETURN
           
           v(1) = -DCM(3,2) + DCM(2,3) !-skewSym(3,2)
           v(2) =  DCM(3,1) - DCM(1,3) ! skewSym(3,1)
           v(3) = -DCM(2,1) + DCM(1,2) !-skewSym(2,1)
   
           indx_max = 1
           do i=2,3
              if ( abs(v(i)) > abs(v(indx_max)) ) indx_max = i
           end do
           
           if ( .not. EqualRealNos( sign(1.0_DbKi,v(indx_max)), sign(1.0_DbKi,logMap(indx_max)) )) logMap = -logMap
           
  
negative.code: !! \f}
  !! where \f$v\f$ is the skew-symmetric matrix associated with the unit-length eigenvector of \f$\Lambda\f$ associated with the eigenvalue 1.
  !! However, this equation has numerical issues near \f$\theta = \pi\f$, so for \f$\theta > 3.1\f$  we instead implement
  !! a separate equation to find lambda * sign(lambda(indx_max))
  !! and use \f$\Lambda - \Lambda^T\f$ to choose the appropriate signs. 
  !!   
  !! This routine is the inverse of DCM_exp (nwtc_num::dcm_exp). \n
  