anchor.diff:              ! NOTE: This will find the closest WaveTime index (wvIndx) which is has waveTime(wvIndx) > = Time.  If WaveDT = DT then waveTime(wvIndx) will equal Time
               ! For WaveMod = 6 or WaveMod = 5 WaveDT must equal DT for the returned value of elementWaterState to be meaningful, for other WaveMod, 
               ! elementWaterState is the same for all time for a given node, J.
  -        elementWaterState = REAL( InterpWrappedStpInt( REAL(Time, SiKi), p%WaveTime(:), p%elementWaterState(:,J), m%LastIndWave, p%NStepWave + 1 ), ReKi )
  -       
  +        elementWaterState = REAL(p%elementWaterState(IntWrapIndx,J), ReKi)
            
            ! Determine the dynamic pressure at the marker
  -         m%D_FDynP(J) = InterpWrappedStpReal ( REAL(Time, SiKi), p%WaveTime(:), p%WaveDynP(:,nodeIndx), &
  -                                    m%LastIndWave, p%NStepWave + 1 )
  +         m%D_FDynP(J) = InterpolateWithSlope(InterpolationSlope, m%LastIndWave, p%WaveDynP(:,nodeIndx))
            
               
            DO I=1,3
                  ! Determine the fluid acceleration and velocity at the marker
  -            m%D_FA(I,J) = InterpWrappedStpReal ( REAL(Time, SiKi), p%WaveTime(:), p%WaveAcc(:,nodeIndx,I), &
  -                                    m%LastIndWave, p%NStepWave + 1       )
  -            m%D_FV(I,J) = InterpWrappedStpReal ( REAL(Time, SiKi), p%WaveTime(:), p%WaveVel(:,nodeIndx,I), &
  -                                    m%LastIndWave, p%NStepWave + 1       )
  +            m%D_FA(I,J) = InterpolateWithSlope(InterpolationSlope, m%LastIndWave, p%WaveAcc(:,nodeIndx,I))
  +            m%D_FV(I,J) = InterpolateWithSlope(InterpolationSlope, m%LastIndWave, p%WaveVel(:,nodeIndx,I))
               
               vrel(I) =  m%D_FV(I,J) - u%DistribMesh%TranslationVel(I,J)
               
  -            m%D_F_I(I,J) = elementWaterState * InterpWrappedStpReal ( REAL(Time, SiKi), p%WaveTime(:), p%D_F_I(:,I,J), &
  -                                    m%LastIndWave, p%NStepWave + 1       )
  +            m%D_F_I(I,J) = elementWaterState * InterpolateWithSlope(InterpolationSlope, m%LastIndWave, p%D_F_I(:,I,J) )
            END DO
            
               ! (k x vrel x k)
  
positive.code: 
        
  END SUBROUTINE Morison_UpdateStates
  
  !----------------------------------------------------------------------------------------------------------------------------------
  !> Routine for computing outputs, used in both loose and tight coupling.
  SUBROUTINE Morison_CalcOutput( Time, u, p, x, xd, z, OtherState, y, m, ErrStat, ErrMsg )   
  
negative.code:          DO I=1,6
                          
              ! We are now combining the dynamic pressure term into the inertia term
              m%L_F_I(I,J) = InterpolateWithSlope(InterpolationSlope, m%LastIndWave, p%L_F_I(:,I,J))
              
              IF (I < 4 ) THEN
        
  