anchor.diff:     Eref = 0.0D0
   
      DO i=1,niter
  -!WRITE(*,*) "N-R Iteration #", i
          StifK = 0.0D0
          RHS = 0.0D0
          MassM = 0.0D0
  
positive.code: 
     CALL BD_CopyInput(u, u_tmp, MESH_NEWCOPY, ErrStat2, ErrMsg2)
        CALL SetErrStat( ErrStat2, ErrMsg2, ErrStat, ErrMsg, RoutineName )
     CALL BD_InputGlobalLocal(p,u_tmp,0,ErrStat2,ErrMsg2)
        CALL SetErrStat( ErrStat2, ErrMsg2, ErrStat, ErrMsg, RoutineName )
     !Initialize displacements and rotations
     DO i=1,p%node_total
  
negative.code:    call BD_Input_extrapinterp( u, utimes, u_interp, t+p%dt, ErrStat2, ErrMsg2 )
        call SetErrStat(ErrStat2,ErrMsg2,ErrStat,ErrMsg,RoutineName)
           
     CALL BD_InputGlobalLocal(p,u_interp,ErrStat2,ErrMsg2)
        call SetErrStat(ErrStat2,ErrMsg2,ErrStat, ErrMsg, RoutineName)
     ! Incorporate boundary conditions
     CALL BD_BoundaryGA2(x,p,u_interp,t,OtherState,ErrStat2,ErrMsg2)
  