anchor.diff:      OVz = OVRR*OVRR-OVL*OVL
   
       ! Finish calculating the distance from the intersection line to the outer edge of the wake
  -    if (OVz > 0.0_dp + tol) then
  +    !if (OVz > 0.0_dp + tol) then
  +    if (OVz > 0.0_dp) then
           OVz = sqrt(OVz)
       else
           OVz = 0.0_dp
       end if
  +    
  +    !print *, "OVRR, OVL, OVRR, OVr, OVdYd, OVz ", OVRR, OVL, OVRR, OVr, OVdYd, OVz
  +    
  +    
   
       if (OVdYd < (OVr+OVRR)) then ! if the rotor overlaps the wake
           !print *, "OVL: ", OVL
  -        !if (OVL < OVRR .and. (OVdYd-OVL) < OVr) then
  -        if (OVdYd > 0.0_dp + tol) then
  -            wake_overlap = OVRR*OVRR*dacos(OVL/OVRR) + OVr*OVr*dacos((OVdYd-OVL)/OVr) - OVdYd*OVz
  +        if (OVL < OVRR .and. (OVdYd-OVL) < OVr) then
  +!         if (OVdYd > 0.0_dp + tol) then
  +!         if ((OVdYd > 0.0_dp) .and. (OVdYd > (OVRR - OVr))) then
  +            ! print *, "acos(OVL/OVRR), acos((OVdYd-OVL)/OVr), OVRR, OVL, OVr, OVdYd, OVL/OVRR, (OVdYd-OVL)/OVr ", &
  +!     & acos(OVL/OVRR), acos((OVdYd-OVL)/OVr), OVRR, OVL, OVr, OVdYd, OVL/OVRR, (OVdYd-OVL)/OVr
  +            wake_overlap = OVRR*OVRR*acos(OVL/OVRR) + OVr*OVr*acos((OVdYd-OVL)/OVr) - OVdYd*OVz
           else if (OVRR > OVr) then
               wake_overlap = pi*OVr*OVr
               !print *, "wake ovl: ", wake_overlap
  
positive.code:         
          ! Check if this is the max and use it if it is
          if (TI_tmp > TI_dst) then
              TI_dst = TI_tmp
          end if
      !print *, "sigma: ", sigma
      ! TODO add other TI calculation methods
          
      ! wake combination method error 
      else
          !print *, "Invalid added TI calculation method. Must be one of [1,2,3]."
          stop 1
      end if            
      
  
negative.code:     
      ! calculate the distance from the wake center to the line perpendicular to the 
      ! line between the two circle intersection points
      !if (OVdYd >= 0.0_dp + tol) then ! check case to avoid division by zero
      if (OVdYd >= 0.0_dp) then ! check case to avoid division by zero
          OVL = (-OVr*OVr+OVRR*OVRR+OVdYd*OVdYd)/(2.0_dp*OVdYd)
      else
          OVL = 0.0_dp
  