stat_dev: init

Author: jvdp1

src/CMakeLists.txt

@@ -5,6 +5,7 @@ set(SRC
     stdlib_experimental_kinds.f90
     stdlib_experimental_optval.f90
     stdlib_experimental_system.F90
+    stdlib_experimental_stat.f90
 )
 
 add_library(fortran_stdlib ${SRC})

src/stdlib_experimental_stat.f90

@@ -0,0 +1,50 @@
+module stdlib_experimental_stat
+use stdlib_experimental_kinds, only: sp, dp, qp
+use stdlib_experimental_error, only: error_stop
+use stdlib_experimental_optval, only: optval
+implicit none
+private
+! Public API
+public :: mean
+
+
+interface mean
+    module procedure mean_1_dp_dp
+    module procedure mean_2_dp_dp
+end interface
+
+contains
+
+pure function mean_1_dp_dp(mat) result(res)
+    real(dp), intent(in) :: mat(:)
+    real(dp) ::res
+
+    res = sum(mat) / real(size(mat), dp)
+
+end function mean_1_dp_dp
+
+function mean_2_dp_dp(mat, dim) result(res)
+    real(dp), intent(in) :: mat(:,:)
+    integer, intent(in), optional :: dim
+    real(dp), allocatable ::res(:)
+
+    integer :: i
+    integer :: dim_
+
+    dim_ = optval(dim, 1)
+
+    allocate(res(size(mat, dim_)))
+
+    if (dim_ == 1) then
+        do i=1, size(mat, dim_)
+            res(i) = mean_1_dp_dp(mat(i,:))
+        end do
+    else if (dim_ == 2) then
+        do i=1, size(mat, dim_)
+            res(i) = mean_1_dp_dp(mat(:,i))
+        end do
+    end if
+
+end function mean_2_dp_dp
+
+end module

src/tests/CMakeLists.txt

@@ -9,6 +9,7 @@ endmacro(ADDTEST)
 add_subdirectory(ascii)
 add_subdirectory(io)
 add_subdirectory(optval)
+add_subdirectory(stat)
 add_subdirectory(system)
 
 ADDTEST(always_skip)

src/tests/stat/CMakeLists.txt

@@ -0,0 +1 @@
+ADDTEST(mean)

src/tests/stat/array1.dat

@@ -0,0 +1,4 @@
+1 2
+3 4
+5 6
+7 8

src/tests/stat/array2.dat

@@ -0,0 +1,4 @@
+1 2 9
+3 4 10
+5 6 11
+7 8 12

src/tests/stat/array3.dat

@@ -0,0 +1,16 @@
+1.000000000000000021e-08 9.199998759392489944e+01
+1.024113254885563425e-08 9.199998731474968849e+01
+1.048233721895820948e-08 9.199998703587728244e+01
+1.072361403187881949e-08 9.199998675729767683e+01
+1.096496300919481796e-08 9.199998647900135040e+01
+1.120638417249036630e-08 9.199998620097916557e+01
+1.144787754335570897e-08 9.199998592322251056e+01
+1.168944314338753750e-08 9.199998564572304360e+01
+1.193108099418952317e-08 9.199998536847290609e+01
+1.217279111737088596e-08 9.199998509146449521e+01
+1.241457353454836993e-08 9.199998481469057765e+01
+1.265642826734443823e-08 9.199998453814424693e+01
+1.289835533738818635e-08 9.199998426181879552e+01
+1.314035476631514857e-08 9.199998398570787117e+01
+1.338242657576766519e-08 9.199998370980536322e+01
+1.362457078739434161e-08 9.199998343410533153e+01

src/tests/stat/array4.dat

@@ -0,0 +1,3 @@
+  1.56367173122998851E-010  4.51568171776229776E-007  4.96568621780730290E-006  5.01068666781180638E-005  5.01518671281225327E-004  5.01763629287519872E-003  5.58487648776459511E-002  0.32618374746711520        1.7639051761733842        9.4101331514118236
+  8.23481961129666271E-010  4.58239319656296504E-007  5.03239769660796763E-006  5.07739814661247314E-005  5.08189819161291786E-004  5.09287863145356859E-003  5.62489258981838380E-002  0.32831192218075922        1.7752234390209392        9.4703270222745211
+  2.02201163784892633E-009  4.70224616423489051E-007  5.15225066427989480E-006  5.19725111428439625E-005  5.20175115928484585E-004  5.22805802989171828E-003  5.69678499382489378E-002  0.33213537295325257        1.7955576815764616        9.5784705410250410

src/tests/stat/test_mean.f90

@@ -0,0 +1,59 @@
+program test_mean
+use stdlib_experimental_error, only: assert
+use stdlib_experimental_kinds, only: sp, dp
+use stdlib_experimental_io, only: loadtxt
+use stdlib_experimental_stat, only: mean
+use stdlib_experimental_error, only: error_stop
+implicit none
+
+real(sp), allocatable :: s(:, :)
+real(dp), allocatable :: d(:, :)
+real(dp), allocatable :: res(:)
+
+!call loadtxt("array1.dat", s)
+!call print_array(s)
+
+call loadtxt("array1.dat", d)
+
+res = mean(d)
+call print_array(d)
+print *,'Mean = ', res
+call assert(sum( res - [1.5_dp, 3.5_dp, 5.5_dp, 7.5_dp] ) == 0.0_dp)
+
+res = mean(d, dim = 2)
+call print_array(d)
+print *,'Mean = ', res
+call assert(sum( res - [4.0_dp, 5.0_dp] ) == 0.0_dp)
+
+!call loadtxt("array2.dat", d)
+!call print_array(d)
+!
+!call loadtxt("array3.dat", d)
+!call print_array(d)
+!
+!call loadtxt("array4.dat", d)
+!call print_array(d)
+
+contains
+
+subroutine print_array(a)
+class(*),intent(in) :: a(:, :)
+integer :: i
+print *, "Array, shape=(", size(a, 1), ",", size(a, 2), ")"
+
+ select type(a)
+  type is(real(sp))
+   do i = 1, size(a, 1)
+    print *, a(i, :)
+   end do
+  type is(real(dp))
+   do i = 1, size(a, 1)
+    print *, a(i, :)
+   end do
+  class default
+   call error_stop('The proposed type is not supported')
+ end select
+
+end subroutine
+
+end program