[stdlib_math] add is_close routines.

doc/specs/stdlib_math.md

@@ -320,25 +320,25 @@ program demo_logspace_rstart_cbase
 
 end program demo_logspace_rstart_cbase
 ```
-## `arange`
+### `arange`
 
-### Status
+#### Status
 
 Experimental
 
-### Class
+#### Class
 
 Pure function.
 
-### Description
+#### Description
 
 Creates a one-dimensional `array` of the `integer/real` type with fixed-spaced values of given spacing, within a given interval.
 
-### Syntax
+#### Syntax
 
 `result = [[stdlib_math(module):arange(interface)]](start [, end, step])`
 
-### Arguments
+#### Arguments
 
 All arguments should be the same type and kind.
 
@@ -354,18 +354,18 @@ The default `end` value is the inputted `start` value.
 This is an `intent(in)` and `optional` argument. 
 The default `step` value is `1`.
 
-#### Warning
+##### Warning
 If `step = 0`, the `step` argument will be corrected to `1/1.0` by the internal process of the `arange` function. 
 If `step < 0`, the `step` argument will be corrected to `abs(step)` by the internal process of the `arange` function. 
 
-### Return value
+#### Return value
 
 Returns a one-dimensional `array` of fixed-spaced values.
 
 For `integer` type arguments, the length of the result vector is `(end - start)/step + 1`. 
 For `real` type arguments, the length of the result vector is `floor((end - start)/step) + 1`.
 
-### Example
+#### Example
 
 ```fortran
 program demo_math_arange
@@ -388,3 +388,143 @@ program demo_math_arange
 
 end program demo_math_arange
 ```
+
+### `is_close`
+
+#### Description
+
+Returns a boolean scalar/array where two scalars/arrays are element-wise equal within a tolerance.
+
+```fortran
+!> For `real` type
+is_close(a, b, rel_tol, abs_tol) = abs(a - b) <= max(rel_tol*(abs(a), abs(b)), abs_tol)
+
+!> and for `complex` type
+is_close(a, b, rel_tol, abs_tol) = is_close(a%re, b%re, rel_tol, abs_tol) .and. &
+ is_close(a%im, b%im, rel_tol, abs_tol)
+```
+
+#### Syntax
+
+`bool = [[stdlib_math(module):is_close(interface)]] (a, b [, rel_tol, abs_tol, equal_nan])`
+
+#### Status
+
+Experimental.
+
+#### Class
+
+Elemental function.
+
+#### Arguments
+
+Note: All `real/complex` arguments must have same `kind`. 
+If the value of `rel_tol/abs_tol` is negative (not recommended), 
+it will be corrected to `abs(rel_tol/abs_tol)` by the internal process of `is_close`.
+
+`a`: Shall be a `real/complex` scalar/array.
+This argument is `intent(in)`.
+
+`b`: Shall be a `real/complex` scalar/array.
+This argument is `intent(in)`.
+
+`rel_tol`: Shall be a `real` scalar/array.
+This argument is `intent(in)` and `optional`, which is `sqrt(epsilon(..))` by default.
+
+`abs_tol`: Shall be a `real` scalar/array.
+This argument is `intent(in)` and `optional`, which is `0.0` by default.
+
+`equal_nan`: Shall be a `logical` scalar/array.
+This argument is `intent(in)` and `optional`, which is `.false.` by default.
+Whether to compare `NaN` values as equal. If `.true.`, 
+`NaN` values in `a` will be considered equal to `NaN` values in `b`.
+
+#### Result value
+
+Returns a `logical` scalar/array.
+
+#### Example
+
+```fortran
+program demo_math_is_close
+
+ use stdlib_math, only: is_close
+ use stdlib_error, only: check
+ real :: x(2) = [1, 2], y, NAN
+
+ y = -3
+ NAN = sqrt(y)
+
+ print *, is_close(x,[real :: 1, 2.1]) !! [T, F]
+ print *, is_close(2.0, 2.1, abs_tol=0.1) !! T
+ print *, NAN, is_close(2.0, NAN), is_close(2.0, NAN, equal_nan=.true.) !! NAN, F, F
+ print *, is_close(NAN, NAN), is_close(NAN, NAN, equal_nan=.true.) !! F, T
+
+end program demo_math_is_close
+```
+
+### `all_close`
+
+#### Description
+
+Returns a boolean scalar where two arrays are element-wise equal within a tolerance.
+
+#### Syntax
+
+`bool = [[stdlib_math(module):all_close(interface)]] (a, b [, rel_tol, abs_tol, equal_nan])`
+
+#### Status
+
+Experimental.
+
+#### Class
+
+Pure function.
+
+#### Arguments
+
+Note: All `real/complex` arguments must have same `kind`. 
+If the value of `rel_tol/abs_tol` is negative (not recommended), 
+it will be corrected to `abs(rel_tol/abs_tol)` by the internal process of `all_close`.
+
+`a`: Shall be a `real/complex` array.
+This argument is `intent(in)`.
+
+`b`: Shall be a `real/complex` array.
+This argument is `intent(in)`.
+
+`rel_tol`: Shall be a `real` scalar.
+This argument is `intent(in)` and `optional`, which is `sqrt(epsilon(..))` by default.
+
+`abs_tol`: Shall be a `real` scalar.
+This argument is `intent(in)` and `optional`, which is `0.0` by default.
+
+`equal_nan`: Shall be a `logical` scalar.
+This argument is `intent(in)` and `optional`, which is `.false.` by default.
+Whether to compare `NaN` values as equal. If `.true.`, 
+`NaN` values in `a` will be considered equal to `NaN` values in `b`.
+
+#### Result value
+
+Returns a `logical` scalar.
+
+#### Example
+
+```fortran
+program demo_math_all_close
+
+ use stdlib_math, only: all_close
+ use stdlib_error, only: check
+ real :: x(2) = [1, 2], y, NAN
+ complex :: z(4, 4)
+
+ y = -3
+ NAN = sqrt(y)
+ z = (1.0, 1.0)
+
+ print *, all_close(z+cmplx(1.0e-11, 1.0e-11), z) !! T
+ print *, NAN, all_close([NAN], [NAN]), all_close([NAN], [NAN], equal_nan=.true.) 
+ !! NAN, F, T
+
+end program demo_math_all_close
+```

src/CMakeLists.txt

@@ -39,6 +39,8 @@ set(fppFiles
  stdlib_math_linspace.fypp
  stdlib_math_logspace.fypp
  stdlib_math_arange.fypp
+ stdlib_math_is_close.fypp
+ stdlib_math_all_close.fypp
  stdlib_string_type.fypp
  stdlib_string_type_constructor.fypp
  stdlib_strings_to_string.fypp

src/Makefile.manual

@@ -34,8 +34,10 @@ SRCFYPP = \
  stdlib_stats_distribution_uniform.fypp \
  stdlib_stats_var.fypp \
  stdlib_math.fypp \
-stdlib_math_linspace.fypp \
-stdlib_math_logspace.fypp \
+ stdlib_math_linspace.fypp \
+ stdlib_math_logspace.fypp \
+ stdlib_math_is_close.fypp \
+ stdlib_math_all_close.fypp \
  stdlib_string_type.fypp \
  stdlib_string_type_constructor.fypp \
  stdlib_strings.fypp \
@@ -80,7 +82,8 @@ $(SRCGEN): %.f90: %.fypp common.fypp
 # Fortran module dependencies
 f18estop.o: stdlib_error.o
 stdlib_ascii.o: stdlib_kinds.o
-stdlib_bitsets.o: stdlib_kinds.o
+stdlib_bitsets.o: stdlib_kinds.o \
+ stdlib_optval.o
 stdlib_bitsets_64.o: stdlib_bitsets.o
 stdlib_bitsets_large.o: stdlib_bitsets.o
 stdlib_error.o: stdlib_optval.o
@@ -103,7 +106,8 @@ stdlib_io_npy_save.o: \
  stdlib_error.o \
  stdlib_strings.o
 stdlib_linalg.o: \
- stdlib_kinds.o
+ stdlib_kinds.o \
+ stdlib_optval.o
 stdlib_linalg_diag.o: \
  stdlib_linalg.o \
  stdlib_kinds.o
@@ -189,6 +193,11 @@ stdlib_math_logspace.o: \
 stdlib_math_arange.o: \
  stdlib_math.o
 stdlib_linalg_outer_product.o: stdlib_linalg.o
+stdlib_math_is_close.o: \
+ stdlib_math.o
+stdlib_math_all_close.o: \
+ stdlib_math.o \
+ stdlib_math_is_close.o
 stdlib_stringlist_type.o: stdlib_string_type.o \
  stdlib_math.o \
-stdlib_optval.o
+ stdlib_optval.o

src/stdlib_math.fypp

@@ -14,7 +14,7 @@ module stdlib_math
  public :: EULERS_NUMBER_QP
 #:endif
  public :: DEFAULT_LINSPACE_LENGTH, DEFAULT_LOGSPACE_BASE, DEFAULT_LOGSPACE_LENGTH
- public :: arange
+ public :: arange, is_close, all_close
 
  integer, parameter :: DEFAULT_LINSPACE_LENGTH = 100
  integer, parameter :: DEFAULT_LOGSPACE_LENGTH = 50
@@ -294,6 +294,39 @@ module stdlib_math
  #:endfor
  end interface arange
 
+ !> Version: experimental
+ !>
+ !> Returns a boolean scalar/array where two scalar/arrays are element-wise equal within a tolerance.
+ !> ([Specification](../page/specs/stdlib_math.html#is_close))
+ interface is_close
+ #:set RC_KINDS_TYPES = REAL_KINDS_TYPES + CMPLX_KINDS_TYPES
+ #:for k1, t1 in RC_KINDS_TYPES
+ elemental module logical function is_close_${t1[0]}$${k1}$(a, b, rel_tol, abs_tol, equal_nan) result(close)
+ ${t1}$, intent(in) :: a, b
+ real(${k1}$), intent(in), optional :: rel_tol, abs_tol
+ logical, intent(in), optional :: equal_nan
+ end function is_close_${t1[0]}$${k1}$
+ #:endfor
+ end interface is_close
+
+ !> Version: experimental
+ !>
+ !> Returns a boolean scalar where two arrays are element-wise equal within a tolerance.
+ !> ([Specification](../page/specs/stdlib_math.html#all_close))
+ interface all_close
+ #:set RC_KINDS_TYPES = REAL_KINDS_TYPES + CMPLX_KINDS_TYPES
+ #:set RANKS = range(1, MAXRANK + 1)
+ #:for k1, t1 in RC_KINDS_TYPES
+ #:for r1 in RANKS
+ logical pure module function all_close_${r1}$_${t1[0]}$${k1}$(a, b, rel_tol, abs_tol, equal_nan) result(close)
+ ${t1}$, intent(in) :: a${ranksuffix(r1)}$, b${ranksuffix(r1)}$
+ real(${k1}$), intent(in), optional :: rel_tol, abs_tol
+ logical, intent(in), optional :: equal_nan
+ end function all_close_${r1}$_${t1[0]}$${k1}$
+ #:endfor
+ #:endfor
+ end interface all_close
+
 contains
 
  #:for k1, t1 in IR_KINDS_TYPES

src/stdlib_math_all_close.fypp

@@ -0,0 +1,25 @@
+#:include "common.fypp"
+#:set RANKS = range(1, MAXRANK + 1)
+#:set RC_KINDS_TYPES = REAL_KINDS_TYPES + CMPLX_KINDS_TYPES
+
+submodule (stdlib_math) stdlib_math_all_close
+
+ implicit none
+
+contains
+
+ #:for k1, t1 in RC_KINDS_TYPES
+ #:for r1 in RANKS
+ logical pure module function all_close_${r1}$_${t1[0]}$${k1}$(a, b, rel_tol, abs_tol, equal_nan) result(close)
+
+ ${t1}$, intent(in) :: a${ranksuffix(r1)}$, b${ranksuffix(r1)}$
+ real(${k1}$), intent(in), optional :: rel_tol, abs_tol
+ logical, intent(in), optional :: equal_nan
+
+ close = all(is_close(a, b, rel_tol, abs_tol, equal_nan))
+
+ end function all_close_${r1}$_${t1[0]}$${k1}$
+ #:endfor
+ #:endfor
+
+end submodule stdlib_math_all_close

src/stdlib_math_arange.fypp

@@ -44,7 +44,7 @@ contains
  step_ = optval(step, 1_${k1}$)
  step_ = sign(merge(step_, 1_${k1}$, step_ /= 0_${k1}$), end_ - start_)
 
- allocate(result((end_ - start_)/step_ + 1))
+ allocate(result((end_ - start_)/step_ + 1_${k1}$))
 
  result = [(i, i=start_, end_, step_)]
 

src/stdlib_math_is_close.fypp

@@ -0,0 +1,47 @@
+#:include "common.fypp"
+
+submodule(stdlib_math) stdlib_math_is_close
+
+ use, intrinsic :: ieee_arithmetic, only: ieee_is_nan
+ implicit none
+
+ #:for k1 in REAL_KINDS
+ real(${k1}$), parameter :: sqrt_eps_${k1}$ = sqrt(epsilon(1.0_${k1}$))
+ #:endfor
+
+contains
+
+ #! Determines whether the values of `a` and `b` are close.
+
+ #:for k1, t1 in REAL_KINDS_TYPES
+ elemental module logical function is_close_${t1[0]}$${k1}$(a, b, rel_tol, abs_tol, equal_nan) result(close)
+ ${t1}$, intent(in) :: a, b
+ real(${k1}$), intent(in), optional :: rel_tol, abs_tol
+ logical, intent(in), optional :: equal_nan
+ logical :: equal_nan_
+
+ equal_nan_ = optval(equal_nan, .false.)
+
+ if (ieee_is_nan(a) .or. ieee_is_nan(b)) then
+ close = merge(.true., .false., equal_nan_ .and. ieee_is_nan(a) .and. ieee_is_nan(b))
+ else
+ close = abs(a - b) <= max( abs(optval(rel_tol, sqrt_eps_${k1}$)*max(abs(a), abs(b))), &
+ abs(optval(abs_tol, 0.0_${k1}$)) )
+ end if 
+
+ end function is_close_${t1[0]}$${k1}$
+ #:endfor
+
+ #:for k1, t1 in CMPLX_KINDS_TYPES
+ elemental module logical function is_close_${t1[0]}$${k1}$(a, b, rel_tol, abs_tol, equal_nan) result(close)
+ ${t1}$, intent(in) :: a, b
+ real(${k1}$), intent(in), optional :: rel_tol, abs_tol
+ logical, intent(in), optional :: equal_nan
+
+ close = is_close_r${k1}$(a%re, b%re, rel_tol, abs_tol, equal_nan) .and. &
+ is_close_r${k1}$(a%im, b%im, rel_tol, abs_tol, equal_nan)
+
+ end function is_close_${t1[0]}$${k1}$
+ #:endfor
+
+end submodule stdlib_math_is_close