round.F

Go to the documentation of this file.

      MODULE round_mod
!
!git $Id$
!====================================================== H. D. Knoble ===
!  Copyright (c) 2002-2025 The ROMS Group           Hernan G. Arango   !
!    Licensed under a MIT/X style license                              !
!    See License_ROMS.md                                               !
!=======================================================================
!                                                                      !
!  Floating point rounding function with a Fuzzy or Tolerant Floor     !
!  function.                                                           !
!                                                                      !
!  On Input:                                                           !
!                                                                      !
!     X         Double precision argument to be operated on. It is     !
!                 assumed that X is represented with m mantissa bits.  !
!                                                                      !
!     CT        Comparison Tolerance such that 0 < CT <= 3-SQRT(5)/2.  !
!                 If the relative difference between X and a whole     !
!                 number is less than CT, then TFLOOR is returned as   !
!                 this whole number. By treating the floating-point    !
!                 numbers as a finite ordered set, note that the       !
!                 heuristic EPS=2.**(-(m-1)) and CT=3*eps causes       !
!                 arguments of TFLOOR/TCEIL to be treated as whole     !
!                 numbers if they are exactly whole numbers or are     !
!                 immediately adjacent to whole number representations.!
!                 Since EPS, the  "distance"  between  floating-point  !
!                 numbers on the unit interval, and m, the number of   !
!                 bits in X mantissa, exist on every  floating-point   !
!                 computer, TFLOOR/TCEIL are consistently definable    !
!                 on every floating-point computer.                    !
!                                                                      !
!  Usage:                                                              !
!                                                                      !
!    CT = 3 * EPSILON(X)     That is, CT is about 1 bit on either      !
!                            side of X mantissa bits.                  !
!     Y = round(X, CT)                                                 !
!                                                                      !
!  References:                                                         !
!                                                                      !
!    P. E. Hagerty, 1978: More on Fuzzy Floor and Ceiling, APL QUOTE   !
!      QUAD 8(4):20-24. (The TFLOOR=FL5 took five years of refereed    !
!      evolution publication).                                         !
!                                                                      !
!    L. M. Breed, 1978: Definitions for Fuzzy Floor and Ceiling, APL   !
!        QUOTE QUAD 8(3):16-23.                                        !
!                                                                      !
!  Adapted from H.D. Knoble code (Penn State University).              !
!                                                                      !
!=======================================================================
!
      USE mod_kinds
!
      implicit none
!
      PUBLIC  :: round           ! Tolerant round function
      PUBLIC  :: tceil           ! Tolerant ceiling function
      PUBLIC  :: tfloor          ! Tolerant floor function
      PRIVATE :: ufloor          ! Unfuzzy floor function
!
      CONTAINS
!
!***********************************************************************
      FUNCTION round (X, CT) RESULT (Y)
!***********************************************************************
!
!  Imported variable declarations.
!
      real (dp), intent(in) :: x, ct
!
!  Local variable declarations.
!
      real(dp) :: y
!
!------------------------------------------------------------------------
!  Compute tolerant round function.
!------------------------------------------------------------------------
!
      y=tfloor(x+0.5_dp,ct)
!
      RETURN
      END FUNCTION round
!
!***********************************************************************
      FUNCTION tceil (X,CT) RESULT (Y)
!***********************************************************************
!
!  Imported variable declarations.
!
      real (dp), intent(in) :: x, ct
!
!  Local variable declarations.
!
      real(dp) :: y
!
!------------------------------------------------------------------------
!  Compute tolerant ceiling function.
!------------------------------------------------------------------------
!
      y=-tfloor(-x,ct)
!
      RETURN
      END FUNCTION tceil
!
!***********************************************************************
      FUNCTION tfloor (X, CT) RESULT (Y)
!***********************************************************************
!
!  Imported variable declarations.
!
      real (dp), intent(in) :: x, ct
!
!  Local variable declarations.
!
      real (dp) :: eps5, q, rmax, y
!
!------------------------------------------------------------------------
!  Compute tolerant floor function.
!------------------------------------------------------------------------
!
!  Hagerty FL5 function
!
      q=1.0_dp
      IF (x.lt.0.0_dp) q=1.0_dp-ct
      rmax=q/(2.0_dp-ct)
      eps5=ct/q
      y=ufloor(x+max(ct,min(rmax,eps5*abs(1.0_dp+ufloor(x)))))
      IF ((x.le.0.0_dp).or.(y-x).lt.rmax) RETURN
      y=y-1.0_dp
!
      RETURN
      END FUNCTION tfloor
!
!***********************************************************************
      FUNCTION ufloor (X) RESULT (Y)
!***********************************************************************
!
!  Imported variable declarations.
!
      real (dp), intent(in) :: x
!
!  Local variable declarations.
!
      real(dp) :: y
!
!-----------------------------------------------------------------------
!  Compute the largest integer algebraically less than or equal to X;
!  that is, the unfuzzy Floor Function.
!-----------------------------------------------------------------------
!
      y=x-mod(x,1.0_dp)-mod(2.0_dp+sign(1.0_dp,x),3.0_dp)
!
      RETURN
      END FUNCTION ufloor
!
      END MODULE round_mod