smooth_cutoff Subroutine

pure subroutine smooth_cutoff(r, cutoff, width, sw, dswdr)

   !> Interatomic distance
   real(wp), intent(in) :: r

   !> Real space cutoff
   real(wp), intent(in) :: cutoff

   !> Width of smooth cutoff
   real(wp), intent(in) :: width

   !> Switching function value
   real(wp), intent(out) :: sw

   !> Derivative of the switching function with respect to distance
   real(wp), intent(out) :: dswdr

   real(wp) :: inner, effective_width, x

   if (width <= 0.0_wp .or. cutoff <= 0.0_wp) then
      sw = 1.0_wp
      dswdr = 0.0_wp
   else
      ! Keep the switching interval within the physical range 0 <= r <= cutoff.
      effective_width = min(width, cutoff)
      inner = cutoff - effective_width
      if (r <= inner) then
         sw = 1.0_wp
         dswdr = 0.0_wp
      else if (r >= cutoff) then
         sw = 0.0_wp
         dswdr = 0.0_wp
      else
         x = (cutoff - r) / effective_width
         ! Quintic Hermite switch with zero first derivatives at both boundaries.
         sw = x**3 * (smoothstep3 + x*(smoothstep4 + smoothstep5*x))
         dswdr = -smoothstep_deriv * x**2 * (1.0_wp - x)**2 / effective_width
      end if
   end if

end subroutine smooth_cutoff