anchor.diff:  
      ! local variables
      INTEGER(IntKi)                          :: k                   ! generic loop counters
  -   INTEGER(IntKi)                          :: i                   ! generic loop counters
  +   INTEGER(IntKi)                          :: i, iStart           ! generic loop counters
      INTEGER(IntKi)                          :: iBody               ! WAMIT body loop counter
      INTEGER(IntKi)                          :: j                   ! generic loop counters
      INTEGER(IntKi)                          :: indx                ! generic loop counters
  
positive.code:    !!!IF ( p_FAST%CompServo == Module_SrvD ) THEN
     !!!END IF
        
        
     ! HydroDyn: copy op to actual states and inputs
     IF ( p_FAST%CompHydro == Module_HD ) THEN
        ThisModule = Module_HD
        if (allocated(y_FAST%Lin%Modules(ThisModule)%Instance(1)%op_x_eig_mag)) then
           ! TODO: WAMIT parameter and continuous states are now an arrays of length NBody, so we need to modify this following code to handle that
           ! We will try to loop over NBody and add each to the state array
           do iBody = 1, HD%p%NBody
              nStates = HD%p%WAMIT(iBody)%SS_Exctn%numStates
              if (nStates > 0) then
                 call GetStateAry(p_FAST, iMode, t, HD%x( STATE_CURR)%WAMIT(iBody)%SS_Exctn%x, y_FAST%Lin%Modules(ThisModule)%Instance(1)%op_x_eig_mag(         :nStates), y_FAST%Lin%Modules(ThisModule)%Instance(1)%op_x_eig_phase(         :nStates))
              end if
           end do
           do iBody = 1, HD%p%NBody
              nStates = HD%p%WAMIT(iBody)%SS_Rdtn%numStates
              if (nStates < size(y_FAST%Lin%Modules(ThisModule)%Instance(1)%op_x_eig_mag)) then
                 call GetStateAry(p_FAST, iMode, t, HD%x( STATE_CURR)%WAMIT(iBody)%SS_Rdtn%x,  y_FAST%Lin%Modules(ThisModule)%Instance(1)%op_x_eig_mag(1+nStates:       ), y_FAST%Lin%Modules(ThisModule)%Instance(1)%op_x_eig_phase(1+nStates:       ))
              end if
           end do
        end if
     END IF
              
     !LIN-TODO: we need this for SD GJH          
     !!!! SubDyn: copy final predictions to actual states
     !!!IF ( p_FAST%CompSub == Module_SD ) THEN
        ! SubDyn:
        ThisModule = Module_SD
        if (allocated(y_FAST%Lin%Modules(ThisModule)%Instance(1)%op_x_eig_mag)) then
           nStates = size(y_FAST%Lin%Modules(ThisModule)%Instance(1)%op_x_eig_mag)/2
           
           call GetStateAry(p_FAST, iMode, t, SD%x( STATE_CURR)%qm,    y_FAST%Lin%Modules(ThisModule)%Instance(1)%op_x_eig_mag(         :nStates), y_FAST%Lin%Modules(ThisModule)%Instance(1)%op_x_eig_phase(         :nStates))
           call GetStateAry(p_FAST, iMode, t, SD%x( STATE_CURR)%qmdot, y_FAST%Lin%Modules(ThisModule)%Instance(1)%op_x_eig_mag(1+nStates:       ), y_FAST%Lin%Modules(ThisModule)%Instance(1)%op_x_eig_phase(1+nStates:       ))
        end if
  
     !!!ELSE IF ( p_FAST%CompSub == Module_ExtPtfm ) THEN
     !!!END IF   
  
  END SUBROUTINE PerturbOP
  !----------------------------------------------------------------------------------------------------------------------------------
  
negative.code:       call SD_JacobianPInput( t_global, SD%Input(1), SD%p, SD%x(STATE_CURR), SD%xd(STATE_CURR), &
                SD%z(STATE_CURR), SD%OtherSt(STATE_CURR),  SD%y, SD%m, ErrStat2, ErrMsg2, &
                dYdu=y_FAST%Lin%Modules(Module_SD)%Instance(1)%D, dXdu=y_FAST%Lin%Modules(Module_SD)%Instance(1)%B )
        if(Failed()) return;
  
        call SD_JacobianPContState( t_global, SD%Input(1), SD%p, SD%x(STATE_CURR), SD%xd(STATE_CURR), &
            SD%z(STATE_CURR), SD%OtherSt(STATE_CURR), SD%y, SD%m, ErrStat2, ErrMsg2,&
            dYdx=y_FAST%Lin%Modules(Module_SD)%Instance(1)%C, dXdx=y_FAST%Lin%Modules(Module_SD)%Instance(1)%A )
        if(Failed()) return;
  
        ! get the operating point
        call SD_GetOP(t_global, SD%Input(1), SD%p, SD%x(STATE_CURR), SD%xd(STATE_CURR), SD%z(STATE_CURR),&
            SD%OtherSt(STATE_CURR), SD%y, SD%m, ErrStat2, ErrMsg2, u_op=y_FAST%Lin%Modules(Module_SD)%Instance(1)%op_u,&
            y_op=y_FAST%Lin%Modules(Module_SD)%Instance(1)%op_y, &
            x_op=y_FAST%Lin%Modules(Module_SD)%Instance(1)%op_x, dx_op=y_FAST%Lin%Modules(Module_SD)%Instance(1)%op_dx)
        if(Failed()) return;
  
        ! write the module matrices:
        if (p_FAST%LinOutMod) then
           OutFileName = trim(LinRootName)//'.'//TRIM(y_FAST%Module_Abrev(Module_SD))      
           call WrLinFile_txt_Head(t_global, p_FAST, y_FAST, y_FAST%Lin%Modules(Module_SD)%Instance(1), OutFileName, Un, ErrStat2, ErrMsg2)
           if(Failed()) return;
              
           if (p_FAST%LinOutJac) then
              ! Jacobians
  