anchor.diff:     INTEGER(IntKi)                  :: j
      INTEGER(IntKi)                  :: temp_id
      INTEGER(IntKi)                  :: allo_stat
  -   INTEGER(IntKi)          :: ErrStat2                     ! Temporary Error status
  -   CHARACTER(ErrMsgLen)    :: ErrMsg2                      ! Temporary Error message
  -   CHARACTER(*), PARAMETER :: RoutineName = 'BD_GenerateStaticElement'
  +   INTEGER(IntKi)                  :: ErrStat2                     ! Temporary Error status
  +   CHARACTER(ErrMsgLen)            :: ErrMsg2                      ! Temporary Error message
  +   CHARACTER(*), PARAMETER         :: RoutineName = 'BD_GenerateStaticElement'
   
      dof_elem = dof_node * node_elem
      rot_elem = (dof_node/2) * node_elem
  
positive.code:    REAL(ReKi)                  :: Oe(6,6)
     REAL(ReKi)                  :: Pe(6,6)
     REAL(ReKi)                  :: Qe(6,6)
     REAL(ReKi)                  :: Nvvv(dof_node*node_elem)
     REAL(ReKi)                  :: temp_Naaa(dof_node*node_elem)
     REAL(ReKi)                  :: temp_aaa(6)
     INTEGER(IntKi)              :: igp
     INTEGER(IntKi)              :: i
     INTEGER(IntKi)              :: j
  
negative.code:    CHARACTER(ErrMsgLen)    :: ErrMsg2                      ! Temporary Error message
     CHARACTER(*), PARAMETER :: RoutineName = 'BD_StaticElementMatrix'
  
     ErrStat  = ErrID_None
     ErrMsg   = ""
     elk(:,:) = 0.0D0
     elf(:)   = 0.0D0
  
     CALL AllocAry(gp,ngp,'Gauss piont array',ErrStat2,ErrMsg2)
        CALL SetErrStat( ErrStat2, ErrMsg2, ErrStat, ErrMsg, RoutineName )
     CALL AllocAry(gw,ngp,'Gauss piont weight function array',ErrStat2,ErrMsg2)
        CALL SetErrStat( ErrStat2, ErrMsg2, ErrStat, ErrMsg, RoutineName )
     CALL AllocAry(hhx,node_elem,'Shape function array',ErrStat2,ErrMsg2)
        CALL SetErrStat( ErrStat2, ErrMsg2, ErrStat, ErrMsg, RoutineName )
     CALL AllocAry(hpx,node_elem,'Derivative of shape function array',ErrStat2,ErrMsg2)
        CALL SetErrStat( ErrStat2, ErrMsg2, ErrStat, ErrMsg, RoutineName )
     CALL AllocAry(GLL_temp,node_elem,'Gauss-Lobatto-Legendre (GLL) point array',ErrStat2,ErrMsg2)
        CALL SetErrStat( ErrStat2, ErrMsg2, ErrStat, ErrMsg, RoutineName )
     CALL AllocAry(w_temp,node_elem,'GLL weight function array',ErrStat2,ErrMsg2)
        CALL SetErrStat( ErrStat2, ErrMsg2, ErrStat, ErrMsg, RoutineName )
     CALL AllocAry(temp_Nvvv,dof_node*node_elem,'Temporary elemental velocity array',ErrStat2,ErrMsg2)
        CALL SetErrStat( ErrStat2, ErrMsg2, ErrStat, ErrMsg, RoutineName )
     CALL AllocAry(temp_Naaa,dof_node*node_elem,'Temporary elemental acceleration array',ErrStat2,ErrMsg2)
        CALL SetErrStat( ErrStat2, ErrMsg2, ErrStat, ErrMsg, RoutineName )
     if (ErrStat >= AbortErrLev) then
        call Cleanup()
        return
     end if
  
     temp_Nvvv(:)  = 0.0D0
     temp_Naaa(:)  = 0.0D0
  
     CALL BD_GenerateGLL(node_elem-1,GLL_temp,w_temp,ErrStat2,ErrMsg2)
        CALL SetErrStat( ErrStat2, ErrMsg2, ErrStat, ErrMsg, RoutineName )
     CALL BD_GaussPointWeight(ngp,gp,gw,ErrStat2,ErrMsg2)
        CALL SetErrStat( ErrStat2, ErrMsg2, ErrStat, ErrMsg, RoutineName )
     DO igp=1,ngp
         gpr = gp(igp)
         CALL BD_ComputeJacobian(gpr,Nuu0,node_elem,dof_node,gp,GLL_temp,ngp,igp,hhx,hpx,Jacobian,ErrStat2,ErrMsg2)
            CALL SetErrStat( ErrStat2, ErrMsg2, ErrStat, ErrMsg, RoutineName )
         CALL BD_GaussPointDataAt0(hhx,hpx,Nuu0,Nrr0,node_elem,dof_node,uu0,E10,ErrStat2,ErrMsg2)
            CALL SetErrStat( ErrStat2, ErrMsg2, ErrStat, ErrMsg, RoutineName )
         Stif(1:6,1:6) = EStif0_GL(1:6,1:6,igp)
         CALL BD_GaussPointData(hhx,hpx,Nuuu,Nrrr,uu0,E10,node_elem,dof_node,&
                                uuu,uup,E1,RR0,kapa,Stif,cet,ErrStat2,ErrMsg2)
            CALL SetErrStat( ErrStat2, ErrMsg2, ErrStat, ErrMsg, RoutineName )
         CALL BD_ElasticForce(E1,RR0,kapa,Stif,cet,Fc,Fd,Oe,Pe,Qe,ErrStat2,ErrMsg2)
            CALL SetErrStat( ErrStat2, ErrMsg2, ErrStat, ErrMsg, RoutineName )
         mmm          = EMass0_GL(1,1,igp)
         mEta(2)      =-EMass0_GL(1,6,igp)
         mEta(3)      = EMass0_GL(1,5,igp)
         rho(1:3,1:3) = EMass0_GL(4:6,4:6,igp)
         CALL BD_GaussPointDataMass(hhx,hpx,temp_Nvvv,temp_Naaa,RR0,node_elem,dof_node,&
                                    vvv,aaa,vvp,mmm,mEta,rho,ErrStat2,ErrMsg2)
            CALL SetErrStat( ErrStat2, ErrMsg2, ErrStat, ErrMsg, RoutineName )
         CALL BD_GravityForce(mmm,mEta,gravity,Fg,ErrStat2,ErrMsg2)
            CALL SetErrStat( ErrStat2, ErrMsg2, ErrStat, ErrMsg, RoutineName )
         if (ErrStat >= AbortErrLev) then
             call Cleanup()
             return
         end if
         Fd(:) = Fd(:) - Fg(:) - Distr_GL(:,igp)
  
         DO i=1,node_elem
  