Add procedure diff.

Author: zoziha

CHANGELOG.md

@@ -16,6 +16,8 @@ Features available from the latest git source
 - new module `stdlib_io_npy`
   [#581](https://github.com/fortran-lang/stdlib/pull/581)
   - new procedures `save_npy`, `load_npy`
+- update module `stdlib_math`
+  - new procedure `diff`
 
 Changes to existing modules
 

doc/specs/stdlib_math.md

@@ -528,3 +528,72 @@ program demo_math_all_close
     
 end program demo_math_all_close
 ```
+
+### `diff`
+
+#### Description
+
+Computes differences between adjacent elements of the input array.
+
+#### Syntax
+
+```fortran
+!> For rank-1 array
+Y = [[stdlib_math(module):diff(interface)]](X [, n])
+!> and for rank-2 array
+Y = [[stdlib_math(module):diff(interface)]](X [, n, dim])
+```
+
+#### Status
+
+Experimental.
+
+#### Class
+
+Pure function.
+
+#### Arguments
+
+Note: If the value of `dim` is not equal to `1` or `2`,
+`1` will be used by the internal process of `diff`. (Not recommended)
+
+`X`: Shall be a `real/integer` and `rank-1/rank-2` array.
+This argument is `intent(in)`.
+
+`n`: Shall be an `integer` scalar.
+This argument is `intent(in)` and `optional`, which is `1` by default.
+It represents to calculate the n-th order difference.
+
+`dim`: Shall be an `integer` scalar.
+This argument is `intent(in)` and `optional`, which is `1` by default.
+It represents to calculate the difference along which dimension.
+
+#### Result value
+
+Note: That `Y` has one fewer element than `X`.
+
+Returns a `real/integer` and `rank-1/rank-2` array.
+
+#### Example
+
+```fortran
+program demo_diff
+
+    use stdlib_math, only: diff
+    implicit none
+    
+    integer :: i(7) = [1, 1, 2, 3, 5, 8, 13]
+    real    :: x(6) = [0, 5, 15, 30, 50, 75]
+    integer :: A(3, 3) = reshape([1, 7, 17, 3, 11, 19, 5, 13, 23], [3, 3])
+    integer :: Y(3, 2)
+    
+    print *, diff(i)        !! [0, 1, 1, 2, 3, 5]
+    print *, diff(x, 2)     !! [5.0, 5.0, 5.0, 5.0]
+    
+    Y = diff(A, n=1, dim=2)
+    print *, Y(1, :)        !! [2, 2]
+    print *, Y(2, :)        !! [4, 2]
+    print *, Y(3, :)        !! [2, 4]
+
+end program demo_diff
+```

src/CMakeLists.txt

@@ -42,6 +42,7 @@ set(fppFiles
     stdlib_math_arange.fypp
     stdlib_math_is_close.fypp
     stdlib_math_all_close.fypp
+    stdlib_math_diff.fypp
     stdlib_string_type.fypp
     stdlib_string_type_constructor.fypp
     stdlib_strings_to_string.fypp

src/Makefile.manual

@@ -39,6 +39,7 @@ SRCFYPP = \
         stdlib_math_logspace.fypp \
         stdlib_math_is_close.fypp \
         stdlib_math_all_close.fypp \
+        stdlib_math_diff.fypp \
         stdlib_string_type.fypp \
         stdlib_string_type_constructor.fypp \
         stdlib_strings.fypp \
@@ -205,6 +206,7 @@ stdlib_math_is_close.o: \
 stdlib_math_all_close.o: \
     stdlib_math.o \
     stdlib_math_is_close.o
+stdlib_math_diff.o: stdlib_math.o
 stdlib_stringlist_type.o: stdlib_string_type.o \
     stdlib_math.o \
     stdlib_optval.o

src/stdlib_math.fypp

@@ -15,6 +15,7 @@ module stdlib_math
 #:endif
     public :: DEFAULT_LINSPACE_LENGTH, DEFAULT_LOGSPACE_BASE, DEFAULT_LOGSPACE_LENGTH
     public :: arange, is_close, all_close
+    public :: diff
 
     integer, parameter :: DEFAULT_LINSPACE_LENGTH = 100
     integer, parameter :: DEFAULT_LOGSPACE_LENGTH = 50
@@ -326,6 +327,26 @@ module stdlib_math
         #:endfor
         #:endfor
     end interface all_close
+    
+    !> Version: experimental
+    !>
+    !> Computes differences between adjacent elements of the input array.
+    !> ([Specification](../page/specs/stdlib_math.html#diff))
+    interface diff
+        #:set RI_KINDS_TYPES = REAL_KINDS_TYPES + INT_KINDS_TYPES
+        #:for k1, t1 in RI_KINDS_TYPES
+        pure module function diff_1_${k1}$(x, n) result(Y)
+            ${t1}$, intent(in) :: x(:)
+            integer, intent(in), optional :: n
+            ${t1}$, allocatable :: Y(:)
+        end function diff_1_${k1}$
+        pure module function diff_2_${k1}$(X, n, dim) result(Y)
+            ${t1}$, intent(in) :: X(:, :)
+            integer, intent(in), optional :: n, dim
+            ${t1}$, allocatable :: Y(:,:)
+        end function diff_2_${k1}$
+        #:endfor
+    end interface diff
 
 contains
 

src/stdlib_math_diff.fypp

@@ -0,0 +1,61 @@
+#! Originally written in 2016 by Keurfon Luu (keurfonluu@outlook.com)
+#! https://github.com/keurfonluu/Forlab/blob/master/src/lib/forlab.f90#L2673
+
+#:include "common.fypp"
+#:set RI_KINDS_TYPES = REAL_KINDS_TYPES + INT_KINDS_TYPES
+submodule (stdlib_math) stdlib_math_diff
+
+    implicit none
+
+contains
+
+    #! `diff` computes differences of arrays of the ${t1}$ type.
+    
+    #:for k1, t1 in RI_KINDS_TYPES
+    pure module function diff_1_${k1}$(X, n) result(Y)
+        ${t1}$, intent(in) :: X(:)
+        integer, intent(in), optional :: n
+        ${t1}$, allocatable :: Y(:)
+        integer :: n_, i
+
+        n_ = optval(n, 1)
+
+        Y = X
+        do i = 1, n_
+            Y = Y(2:) - Y(:size(Y) - 1)
+        end do
+
+    end function diff_1_${k1}$
+
+    pure module function diff_2_${k1}$(X, n, dim) result(Y)
+        ${t1}$, intent(in) :: X(:, :)
+        integer, intent(in), optional :: n, dim
+        ${t1}$, allocatable :: Y(:, :)
+        integer :: n_, dim_, i
+
+        n_ = optval(n, 1)
+        if (present(dim)) then
+            if (dim == 1 .or. dim == 2) then
+                dim_ = dim
+            else
+                dim_ = 1
+            end if
+        else
+            dim_ = 1
+        end if
+        
+        Y = X
+        if (dim_ == 1) then
+            do i = 1, n_
+                Y = Y(2:, :) - Y(:size(Y, 1) - 1, :)
+            end do
+        elseif (dim == 2) then
+            do i = 1, n_
+                Y = Y(:, 2:) - Y(:, :size(Y, 2) - 1)
+            end do
+        end if
+
+    end function diff_2_${k1}$
+    #:endfor
+
+end submodule stdlib_math_diff

src/tests/math/test_stdlib_math.fypp

@@ -4,7 +4,7 @@
 
 module test_stdlib_math
     use testdrive, only : new_unittest, unittest_type, error_type, check, skip_test
-    use stdlib_math, only: clip, is_close, all_close
+    use stdlib_math, only: clip, is_close, all_close, diff
     use stdlib_kinds, only: int8, int16, int32, int64, sp, dp, xdp, qp
     implicit none
 
@@ -40,6 +40,14 @@ contains
             , new_unittest("all_close-real-${k1}$", test_all_close_real_${k1}$) &
             , new_unittest("all_close-cmplx-${k1}$", test_all_close_cmplx_${k1}$) &
             #:endfor
+            
+            !> Tests for `diff`
+            #:for k1 in REAL_KINDS
+            , new_unittest("diff-real-${k1}$", test_diff_real_${k1}$) &
+            #:endfor
+            #:for k1 in INT_KINDS
+            , new_unittest("diff-int-${k1}$", test_diff_int_${k1}$) &
+            #:endfor
             ]
 
     end subroutine collect_stdlib_math
@@ -294,6 +302,48 @@ contains
 
     end subroutine test_all_close_cmplx_${k1}$
     #:endfor
+    
+    #:for k1, t1 in REAL_KINDS_TYPES
+    subroutine test_diff_real_${k1}$(error)
+        type(error_type), allocatable, intent(out) :: error
+        ${t1}$ :: x(6) = [${t1}$ :: 0, 5, 15, 30, 50, 75]
+        ${t1}$ :: A(1, 3) = reshape([${t1}$ :: 1, 3, 5], [1, 3])
+        
+        call check(error, all_close(diff(x), [${t1}$ :: 5, 10, 15, 20, 25]), &
+            "diff(x) in test_diff_real_${k1}$ failed")
+        if (allocated(error)) return
+        
+        call check(error, all_close(diff(x, n=2), [${t1}$ :: 5, 5, 5, 5]), &
+            "diff(x, n=2) in test_diff_real_${k1}$ failed")
+        if (allocated(error)) return
+        
+        call check(error, all_close(diff(A, n=1, dim=2), reshape([${t1}$ :: 2, 2], [1, 2])), &
+            "diff(x, n=1, dim=2) in test_diff_real_${k1}$ failed")
+        if (allocated(error)) return
+        
+    end subroutine test_diff_real_${k1}$
+    #:endfor
+    
+    #:for k1, t1 in INT_KINDS_TYPES
+    subroutine test_diff_int_${k1}$(error)
+        type(error_type), allocatable, intent(out) :: error
+        ${t1}$ :: x(6) = [${t1}$ :: 0, 5, 15, 30, 50, 75]
+        ${t1}$ :: A(1, 3) = reshape([${t1}$ :: 1, 3, 5], [1, 3])
+        
+        call check(error, all(diff(x) == [${t1}$ :: 5, 10, 15, 20, 25]), &
+            "diff(x) in test_diff_int_${k1}$ failed")
+        if (allocated(error)) return
+        
+        call check(error, all(diff(x, n=2) == [${t1}$ :: 5, 5, 5, 5]), &
+            "diff(x, n=2) in test_diff_int_${k1}$ failed")
+        if (allocated(error)) return
+        
+        call check(error, all(diff(A, n=1, dim=2) == reshape([${t1}$ :: 2, 2], [1, 2])), &
+            "diff(A, n=1, dim=2) in test_diff_int_${k1}$ failed")
+        if (allocated(error)) return
+        
+    end subroutine test_diff_int_${k1}$
+    #:endfor
 
 end module test_stdlib_math