From 7c3005670c4263b060ac0e7b21ae69eaf5c93dbe Mon Sep 17 00:00:00 2001 From: Juan Fiol Date: Tue, 4 Feb 2020 10:58:45 -0300 Subject: [PATCH 1/3] Added support and test for mean of complex arrays --- src/stdlib_experimental_stats.fypp | 73 +++++++++++++++--------- src/stdlib_experimental_stats.md | 4 +- src/stdlib_experimental_stats_mean.fypp | 75 ++++++++++++++++--------- src/tests/stats/test_mean.f90 | 33 +++++++++++ src/tests/stats/test_mean_f03.f90 | 15 +++++ 5 files changed, 143 insertions(+), 57 deletions(-) diff --git a/src/stdlib_experimental_stats.fypp b/src/stdlib_experimental_stats.fypp index 0edb2ba9e..48beeedf8 100644 --- a/src/stdlib_experimental_stats.fypp +++ b/src/stdlib_experimental_stats.fypp @@ -1,8 +1,6 @@ #:include "common.fypp" - #:set RANKS = range(1, MAXRANK + 1) - - +#:set RC_KINDS_TYPES = REAL_KINDS_TYPES + CMPLX_KINDS_TYPES module stdlib_experimental_stats use stdlib_experimental_kinds, only: sp, dp, qp, & int8, int16, int32, int64 @@ -12,92 +10,113 @@ module stdlib_experimental_stats public :: mean interface mean - - #:for k1, t1 in REAL_KINDS_TYPES +#:def name(Rank, Type, Kind) +$:"mean_{0}_all_{1}{2}_{1}{2}".format(Rank, Type[0], Kind) +#:enddef + #:for k1, t1 in RC_KINDS_TYPES #:for rank in RANKS - module function mean_${rank}$_all_${k1}$_${k1}$(x, mask) result(res) + module function ${name(rank, t1, k1)}$ (x, mask) result(res) ${t1}$, intent(in) :: x${ranksuffix(rank)}$ logical, intent(in), optional :: mask ${t1}$ :: res - end function mean_${rank}$_all_${k1}$_${k1}$ + end function ${name(rank, t1, k1)}$ #:endfor #:endfor +#:def name(Rank, Type, Kind) +$:"mean_{0}_all_{1}{2}_rdp".format(Rank, Type[0], Kind) +#:enddef #:for k1, t1 in INT_KINDS_TYPES #:for rank in RANKS - module function mean_${rank}$_all_${k1}$_dp(x, mask) result(res) + module function ${name(rank, t1, k1)}$ (x, mask) result(res) ${t1}$, intent(in) :: x${ranksuffix(rank)}$ logical, intent(in), optional :: mask real(dp) :: res - end function mean_${rank}$_all_${k1}$_dp + end function ${name(rank, t1, k1)}$ #:endfor #:endfor - #:for k1, t1 in REAL_KINDS_TYPES +#:def name(Rank, Type, Kind) +$:"mean_{0}_{1}{2}_{1}{2}".format(Rank, Type[0], Kind) +#:enddef + #:for k1, t1 in RC_KINDS_TYPES #:for rank in RANKS - module function mean_${rank}$_${k1}$_${k1}$(x, dim, mask) result(res) + module function ${name(rank, t1, k1)}$(x, dim, mask) result(res) ${t1}$, intent(in) :: x${ranksuffix(rank)}$ integer, intent(in) :: dim logical, intent(in), optional :: mask ${t1}$ :: res${reduced_shape('x', rank, 'dim')}$ - end function mean_${rank}$_${k1}$_${k1}$ + end function ${name(rank, t1, k1)}$ #:endfor #:endfor +#:def name(Rank, Type, Kind) +$:"mean_{0}_{1}{2}_rdp".format(Rank, Type[0], Kind) +#:enddef #:for k1, t1 in INT_KINDS_TYPES #:for rank in RANKS - module function mean_${rank}$_${k1}$_dp(x, dim, mask) result(res) + module function ${name(rank, t1, k1)}$(x, dim, mask) result(res) ${t1}$, intent(in) :: x${ranksuffix(rank)}$ integer, intent(in) :: dim logical, intent(in), optional :: mask real(dp) :: res${reduced_shape('x', rank, 'dim')}$ - end function mean_${rank}$_${k1}$_dp + end function ${name(rank, t1, k1)}$ #:endfor #:endfor - - #:for k1, t1 in REAL_KINDS_TYPES +#:def name(Rank, Type, Kind) +$:"mean_{0}_mask_all_{1}{2}_{1}{2}".format(Rank, Type[0], Kind) +#:enddef + #:for k1, t1 in RC_KINDS_TYPES #:for rank in RANKS - module function mean_${rank}$_mask_all_${k1}$_${k1}$(x, mask) result(res) + module function ${name(rank, t1, k1)}$(x, mask) result(res) ${t1}$, intent(in) :: x${ranksuffix(rank)}$ logical, intent(in) :: mask${ranksuffix(rank)}$ ${t1}$ :: res - end function mean_${rank}$_mask_all_${k1}$_${k1}$ + end function ${name(rank, t1, k1)}$ + #:endfor #:endfor - +#:def name(Rank, Type, Kind) +$:"mean_{0}_mask_all_{1}{2}_rdp".format(Rank, Type[0], Kind) +#:enddef #:for k1, t1 in INT_KINDS_TYPES #:for rank in RANKS - module function mean_${rank}$_mask_all_${k1}$_dp(x, mask) result(res) + module function ${name(rank, t1, k1)}$(x, mask) result(res) ${t1}$, intent(in) :: x${ranksuffix(rank)}$ logical, intent(in) :: mask${ranksuffix(rank)}$ real(dp) :: res - end function mean_${rank}$_mask_all_${k1}$_dp + end function ${name(rank, t1, k1)}$ #:endfor #:endfor - #:for k1, t1 in REAL_KINDS_TYPES +#:def name(Rank, Type, Kind) +$:"mean_{0}_mask_{1}{2}_{1}{2}".format(Rank, Type[0], Kind) +#:enddef + #:for k1, t1 in RC_KINDS_TYPES #:for rank in RANKS - module function mean_${rank}$_mask_${k1}$_${k1}$(x, dim, mask) result(res) + module function ${name(rank, t1, k1)}$(x, dim, mask) result(res) ${t1}$, intent(in) :: x${ranksuffix(rank)}$ integer, intent(in) :: dim logical, intent(in) :: mask${ranksuffix(rank)}$ ${t1}$ :: res${reduced_shape('x', rank, 'dim')}$ - end function mean_${rank}$_mask_${k1}$_${k1}$ + end function ${name(rank, t1, k1)}$ #:endfor #:endfor - +#:def name(Rank, Type, Kind) +$:"mean_{0}_mask_{1}{2}_rdp".format(Rank, Type[0], Kind) +#:enddef #:for k1, t1 in INT_KINDS_TYPES #:for rank in RANKS - module function mean_${rank}$_mask_${k1}$_dp(x, dim, mask) result(res) + module function ${name(rank, t1, k1)}$(x, dim, mask) result(res) ${t1}$, intent(in) :: x${ranksuffix(rank)}$ integer, intent(in) :: dim logical, intent(in) :: mask${ranksuffix(rank)}$ real(dp) :: res${reduced_shape('x', rank, 'dim')}$ - end function mean_${rank}$_mask_${k1}$_dp + end function ${name(rank, t1, k1)}$ #:endfor #:endfor diff --git a/src/stdlib_experimental_stats.md b/src/stdlib_experimental_stats.md index 1059ef65c..0c11c6822 100644 --- a/src/stdlib_experimental_stats.md +++ b/src/stdlib_experimental_stats.md @@ -18,7 +18,7 @@ Returns the mean of all the elements of `array`, or of the elements of `array` a ### Arguments -`array`: Shall be an array of type `integer`, or `real`. +`array`: Shall be an array of type `integer`, `real`, or `complex`. `dim`: Shall be a scalar of type `integer` with a value in the range from 1 to n, where n is the rank of `array`. @@ -26,7 +26,7 @@ Returns the mean of all the elements of `array`, or of the elements of `array` a ### Return value -If `array` is of type `real`, the result is of the same type as `array`. +If `array` is of type `real` or `complex`, the result is of the same type as `array`. If `array` is of type `integer`, the result is of type `double precision`. If `dim` is absent, a scalar with the mean of all elements in `array` is returned. Otherwise, an array of rank n-1, where n equals the rank of `array`, and a shape similar to that of `array` with dimension `dim` dropped is returned. diff --git a/src/stdlib_experimental_stats_mean.fypp b/src/stdlib_experimental_stats_mean.fypp index 26b5d47ba..1c79f9853 100644 --- a/src/stdlib_experimental_stats_mean.fypp +++ b/src/stdlib_experimental_stats_mean.fypp @@ -1,8 +1,6 @@ #:include "common.fypp" - #:set RANKS = range(1, MAXRANK + 1) - - +#:set RC_KINDS_TYPES = REAL_KINDS_TYPES + CMPLX_KINDS_TYPES submodule (stdlib_experimental_stats) stdlib_experimental_stats_mean use, intrinsic:: ieee_arithmetic, only: ieee_value, ieee_quiet_nan @@ -12,28 +10,33 @@ submodule (stdlib_experimental_stats) stdlib_experimental_stats_mean contains - #:for k1, t1 in REAL_KINDS_TYPES +#:def name(Rank, Type, Kind) +$:"mean_{0}_all_{1}{2}_{1}{2}".format(Rank, Type[0], Kind) +#:enddef + #:for k1, t1 in RC_KINDS_TYPES #:for rank in RANKS - module function mean_${rank}$_all_${k1}$_${k1}$(x, mask) result(res) + module function ${name(rank, t1, k1)}$ (x, mask) result(res) ${t1}$, intent(in) :: x${ranksuffix(rank)}$ logical, intent(in), optional :: mask ${t1}$ :: res if (.not.optval(mask, .true.)) then - res = ieee_value(res, ieee_quiet_nan) + res = ieee_value(real(res, kind=${k1}$), ieee_quiet_nan) return end if res = sum(x) / real(size(x, kind = int64), ${k1}$) - end function mean_${rank}$_all_${k1}$_${k1}$ + end function ${name(rank, t1, k1)}$ #:endfor #:endfor - +#:def name(Rank, Type, Kind) +$:"mean_{0}_all_{1}{2}_rdp".format(Rank, Type[0], Kind) +#:enddef #:for k1, t1 in INT_KINDS_TYPES #:for rank in RANKS - module function mean_${rank}$_all_${k1}$_dp(x, mask) result(res) + module function ${name(rank, t1, k1)}$(x, mask) result(res) ${t1}$, intent(in) :: x${ranksuffix(rank)}$ logical, intent(in), optional :: mask real(dp) :: res @@ -45,21 +48,24 @@ contains res = sum(real(x, dp)) / real(size(x, kind = int64), dp) - end function mean_${rank}$_all_${k1}$_dp + end function ${name(rank, t1, k1)}$ #:endfor #:endfor - #:for k1, t1 in REAL_KINDS_TYPES +#:def name(Rank, Type, Kind) +$:"mean_{0}_{1}{2}_{1}{2}".format(Rank, Type[0], Kind) +#:enddef + #:for k1, t1 in RC_KINDS_TYPES #:for rank in RANKS - module function mean_${rank}$_${k1}$_${k1}$(x, dim, mask) result(res) + module function ${name(rank, t1, k1)}$(x, dim, mask) result(res) ${t1}$, intent(in) :: x${ranksuffix(rank)}$ integer, intent(in) :: dim logical, intent(in), optional :: mask ${t1}$ :: res${reduced_shape('x', rank, 'dim')}$ if (.not.optval(mask, .true.)) then - res = ieee_value(res, ieee_quiet_nan) + res = ieee_value(real(res, kind=${k1}$), ieee_quiet_nan) return end if @@ -69,14 +75,17 @@ contains call error_stop("ERROR (mean): wrong dimension") end if - end function mean_${rank}$_${k1}$_${k1}$ + end function ${name(rank, t1, k1)}$ #:endfor #:endfor +#:def name(Rank, Type, Kind) +$:"mean_{0}_{1}{2}_rdp".format(Rank, Type[0], Kind) +#:enddef #:for k1, t1 in INT_KINDS_TYPES #:for rank in RANKS - module function mean_${rank}$_${k1}$_dp(x, dim, mask) result(res) + module function ${name(rank, t1, k1)}$(x, dim, mask) result(res) ${t1}$, intent(in) :: x${ranksuffix(rank)}$ integer, intent(in) :: dim logical, intent(in), optional :: mask @@ -93,42 +102,49 @@ contains call error_stop("ERROR (mean): wrong dimension") end if - end function mean_${rank}$_${k1}$_dp + end function ${name(rank, t1, k1)}$ #:endfor #:endfor - - #:for k1, t1 in REAL_KINDS_TYPES +#:def name(Rank, Type, Kind) +$:"mean_{0}_mask_all_{1}{2}_{1}{2}".format(Rank, Type[0], Kind) +#:enddef + #:for k1, t1 in RC_KINDS_TYPES #:for rank in RANKS - module function mean_${rank}$_mask_all_${k1}$_${k1}$(x, mask) result(res) + module function ${name(rank, t1, k1)}$(x, mask) result(res) ${t1}$, intent(in) :: x${ranksuffix(rank)}$ logical, intent(in) :: mask${ranksuffix(rank)}$ ${t1}$ :: res res = sum(x, mask) / real(count(mask, kind = int64), ${k1}$) - end function mean_${rank}$_mask_all_${k1}$_${k1}$ + end function ${name(rank, t1, k1)}$ #:endfor #:endfor +#:def name(Rank, Type, Kind) +$:"mean_{0}_mask_all_{1}{2}_rdp".format(Rank, Type[0], Kind) +#:enddef #:for k1, t1 in INT_KINDS_TYPES #:for rank in RANKS - module function mean_${rank}$_mask_all_${k1}$_dp(x, mask) result(res) + module function ${name(rank, t1, k1)}$(x, mask) result(res) ${t1}$, intent(in) :: x${ranksuffix(rank)}$ logical, intent(in) :: mask${ranksuffix(rank)}$ real(dp) :: res res = sum(real(x, dp), mask) / real(count(mask, kind = int64), dp) - end function mean_${rank}$_mask_all_${k1}$_dp + end function ${name(rank, t1, k1)}$ #:endfor #:endfor - - #:for k1, t1 in REAL_KINDS_TYPES +#:def name(Rank, Type, Kind) +$:"mean_{0}_mask_{1}{2}_{1}{2}".format(Rank, Type[0], Kind) +#:enddef + #:for k1, t1 in RC_KINDS_TYPES #:for rank in RANKS - module function mean_${rank}$_mask_${k1}$_${k1}$(x, dim, mask) result(res) + module function ${name(rank, t1, k1)}$(x, dim, mask) result(res) ${t1}$, intent(in) :: x${ranksuffix(rank)}$ integer, intent(in) :: dim logical, intent(in) :: mask${ranksuffix(rank)}$ @@ -140,14 +156,17 @@ contains call error_stop("ERROR (mean): wrong dimension") end if - end function mean_${rank}$_mask_${k1}$_${k1}$ + end function ${name(rank, t1, k1)}$ #:endfor #:endfor +#:def name(Rank, Type, Kind) +$:"mean_{0}_mask_{1}{2}_rdp".format(Rank, Type[0], Kind) +#:enddef #:for k1, t1 in INT_KINDS_TYPES #:for rank in RANKS - module function mean_${rank}$_mask_${k1}$_dp(x, dim, mask) result(res) + module function ${name(rank, t1, k1)}$(x, dim, mask) result(res) ${t1}$, intent(in) :: x${ranksuffix(rank)}$ integer, intent(in) :: dim logical, intent(in) :: mask${ranksuffix(rank)}$ @@ -159,7 +178,7 @@ contains call error_stop("ERROR (mean): wrong dimension") end if - end function mean_${rank}$_mask_${k1}$_dp + end function ${name(rank, t1, k1)}$ #:endfor #:endfor diff --git a/src/tests/stats/test_mean.f90 b/src/tests/stats/test_mean.f90 index 49e5a98fc..e2f165124 100644 --- a/src/tests/stats/test_mean.f90 +++ b/src/tests/stats/test_mean.f90 @@ -13,9 +13,13 @@ program test_mean real(sp), allocatable :: s(:, :) real(dp), allocatable :: d(:, :) +complex(dp), allocatable :: cs(:, :) +complex(dp), allocatable :: cd(:, :) + real(dp), allocatable :: d3(:, :, :) real(dp), allocatable :: d4(:, :, :, :) +complex(dp), allocatable :: cd3(:, :, :) !sp call loadtxt("array3.dat", s) @@ -36,6 +40,23 @@ program test_mean call assert( sum( abs( mean(d,1) - sum(d,1)/real(size(d,1), dp) )) < dptol) call assert( sum( abs( mean(d,2) - sum(d,2)/real(size(d,2), dp) )) < dptol) +!csp + +call loadtxt("array3.dat", d) +cs = cmplx(1._sp, 1._sp)*d + +call assert( abs(mean(cs) - sum(cs)/real(size(cs), sp)) < sptol) +call assert( sum( abs( mean(cs,1) - sum(cs,1)/real(size(cs,1), sp) )) < sptol) +call assert( sum( abs( mean(cs,2) - sum(cs,2)/real(size(cs,2), sp) )) < sptol) + +!cdp + +call loadtxt("array3.dat", d) +cd = cmplx(1._dp, 1._dp)*d + +call assert( abs(mean(cd) - sum(cd)/real(size(cd), dp)) < dptol) +call assert( sum( abs( mean(cd,1) - sum(cd,1)/real(size(cd,1), dp) )) < dptol) +call assert( sum( abs( mean(cd,2) - sum(cd,2)/real(size(cd,2), dp) )) < dptol) ! check mask = .false. @@ -75,6 +96,18 @@ program test_mean call assert( sum( abs( mean(d3,2) - sum(d3,2)/real(size(d3,2), dp) )) < dptol) call assert( sum( abs( mean(d3,3) - sum(d3,3)/real(size(d3,3), dp) )) < dptol) +!cdp rank 3 +allocate(cd3(size(d,1),size(d,2),3)) +cd3(:,:,1)=d; +cd3(:,:,2)=d*1.5; +cd3(:,:,3)=d*4; +cd3 = cmplx(1._sp, 1._sp)*cd3 + +call assert( abs(mean(cd3) - sum(cd3)/real(size(cd3), dp)) < dptol) +call assert( sum( abs( mean(cd3,1) - sum(cd3,1)/real(size(cd3,1), dp) )) < dptol) +call assert( sum( abs( mean(cd3,2) - sum(cd3,2)/real(size(cd3,2), dp) )) < dptol) +call assert( sum( abs( mean(cd3,3) - sum(cd3,3)/real(size(cd3,3), dp) )) < dptol) + !dp rank 4 allocate(d4(size(d,1),size(d,2),3,9)) diff --git a/src/tests/stats/test_mean_f03.f90 b/src/tests/stats/test_mean_f03.f90 index fbb09b7ba..ce6650d4a 100644 --- a/src/tests/stats/test_mean_f03.f90 +++ b/src/tests/stats/test_mean_f03.f90 @@ -7,6 +7,7 @@ program test_mean real(dp), allocatable :: d(:, :) real(dp), allocatable :: d8(:, :, :, :, :, :, :, :) +complex(dp), allocatable :: cd8(:, :, :, :, :, :, :, :) !dp @@ -33,6 +34,20 @@ program test_mean call assert( sum( abs( mean(d8,7) - sum(d8,7)/real(size(d8,7), dp) )) == 0.0_dp) call assert( sum( abs( mean(d8,8) - sum(d8,8)/real(size(d8,8), dp) )) == 0.0_dp) +!cdp rank 8 +allocate(cd8(size(d,1), size(d,2), 3, 4, 5, 6, 7, 8)) +cd8 = cmplx(1._dp, 1._dp, kind=dp)*d8 + +call assert( mean(cd8) - sum(cd8)/real(size(cd8), dp) == 0.0_dp) + +call assert( sum( abs( mean(cd8,1) - sum(cd8,1)/real(size(cd8,1), dp) )) == 0.0_dp) +call assert( sum( abs( mean(cd8,2) - sum(cd8,2)/real(size(cd8,2), dp) )) == 0.0_dp) +call assert( sum( abs( mean(cd8,3) - sum(cd8,3)/real(size(cd8,3), dp) )) == 0.0_dp) +call assert( sum( abs( mean(cd8,4) - sum(cd8,4)/real(size(cd8,4), dp) )) == 0.0_dp) +call assert( sum( abs( mean(cd8,5) - sum(cd8,5)/real(size(cd8,5), dp) )) == 0.0_dp) +call assert( sum( abs( mean(cd8,6) - sum(cd8,6)/real(size(cd8,6), dp) )) == 0.0_dp) +call assert( sum( abs( mean(cd8,7) - sum(cd8,7)/real(size(cd8,7), dp) )) == 0.0_dp) +call assert( sum( abs( mean(cd8,8) - sum(cd8,8)/real(size(cd8,8), dp) )) == 0.0_dp) contains end program From bb852d232d4dfa0c82f401db15570e93674482aa Mon Sep 17 00:00:00 2001 From: Juan Fiol Date: Tue, 4 Feb 2020 14:30:15 -0300 Subject: [PATCH 2/3] Added macro to create routine names --- src/common.fypp | 15 ++++++ src/stdlib_experimental_stats.fypp | 66 +++++++++---------------- src/stdlib_experimental_stats_mean.fypp | 64 +++++++++--------------- src/tests/stats/test_mean_f03.f90 | 44 +++++++++-------- 4 files changed, 86 insertions(+), 103 deletions(-) diff --git a/src/common.fypp b/src/common.fypp index 85d82704c..99fed2250 100644 --- a/src/common.fypp +++ b/src/common.fypp @@ -97,5 +97,20 @@ ${prefix + joinstr.join([line.strip() for line in txt.split("\n")]) + suffix}$ #:endif #:enddef +#! Generates a routine name from a generic name, rank, type and kind +#! +#! Args: +#! gname (str): Generic name +#! rank (integer): Rank if exist +#! type (str): Type of the input +#! kind (str): kind of inputs variable +#! suffix (str): other identifier (could be used for output type/kind) +#! +#! Returns: +#! A string with a new name +#! +#:def rname(gname, rank, type, kind, suffix='') + $:"{0}_{1}_{2}{3}_{2}{3}".format(gname, rank, type[0], kind) if suffix == '' else "{0}_{1}_{2}{3}_{4}".format(gname, rank, type[0], kind, suffix) +#:enddef #:endmute diff --git a/src/stdlib_experimental_stats.fypp b/src/stdlib_experimental_stats.fypp index 48beeedf8..f16fa0a92 100644 --- a/src/stdlib_experimental_stats.fypp +++ b/src/stdlib_experimental_stats.fypp @@ -10,113 +10,95 @@ module stdlib_experimental_stats public :: mean interface mean -#:def name(Rank, Type, Kind) -$:"mean_{0}_all_{1}{2}_{1}{2}".format(Rank, Type[0], Kind) -#:enddef #:for k1, t1 in RC_KINDS_TYPES #:for rank in RANKS - module function ${name(rank, t1, k1)}$ (x, mask) result(res) + #:set RName = rname("mean_all",rank, t1, k1) + module function ${RName}$ (x, mask) result(res) ${t1}$, intent(in) :: x${ranksuffix(rank)}$ logical, intent(in), optional :: mask ${t1}$ :: res - end function ${name(rank, t1, k1)}$ + end function ${RName}$ #:endfor #:endfor -#:def name(Rank, Type, Kind) -$:"mean_{0}_all_{1}{2}_rdp".format(Rank, Type[0], Kind) -#:enddef #:for k1, t1 in INT_KINDS_TYPES #:for rank in RANKS - module function ${name(rank, t1, k1)}$ (x, mask) result(res) + #:set RName = rname('mean_all', rank, t1, k1,'dp') + module function ${RName}$(x, mask) result(res) ${t1}$, intent(in) :: x${ranksuffix(rank)}$ logical, intent(in), optional :: mask real(dp) :: res - end function ${name(rank, t1, k1)}$ + end function ${RName}$ #:endfor #:endfor -#:def name(Rank, Type, Kind) -$:"mean_{0}_{1}{2}_{1}{2}".format(Rank, Type[0], Kind) -#:enddef #:for k1, t1 in RC_KINDS_TYPES #:for rank in RANKS - module function ${name(rank, t1, k1)}$(x, dim, mask) result(res) + #:set RName = rname("mean",rank, t1, k1) + module function ${RName}$(x, dim, mask) result(res) ${t1}$, intent(in) :: x${ranksuffix(rank)}$ integer, intent(in) :: dim logical, intent(in), optional :: mask ${t1}$ :: res${reduced_shape('x', rank, 'dim')}$ - end function ${name(rank, t1, k1)}$ + end function ${RName}$ #:endfor #:endfor -#:def name(Rank, Type, Kind) -$:"mean_{0}_{1}{2}_rdp".format(Rank, Type[0], Kind) -#:enddef #:for k1, t1 in INT_KINDS_TYPES #:for rank in RANKS - module function ${name(rank, t1, k1)}$(x, dim, mask) result(res) + #:set RName = rname("mean",rank, t1, k1,'dp') + module function ${RName}$(x, dim, mask) result(res) ${t1}$, intent(in) :: x${ranksuffix(rank)}$ integer, intent(in) :: dim logical, intent(in), optional :: mask real(dp) :: res${reduced_shape('x', rank, 'dim')}$ - end function ${name(rank, t1, k1)}$ + end function ${RName}$ #:endfor #:endfor -#:def name(Rank, Type, Kind) -$:"mean_{0}_mask_all_{1}{2}_{1}{2}".format(Rank, Type[0], Kind) -#:enddef #:for k1, t1 in RC_KINDS_TYPES #:for rank in RANKS - module function ${name(rank, t1, k1)}$(x, mask) result(res) + #:set RName = rname('mean_mask_all',rank, t1, k1) + module function ${RName}$(x, mask) result(res) ${t1}$, intent(in) :: x${ranksuffix(rank)}$ logical, intent(in) :: mask${ranksuffix(rank)}$ ${t1}$ :: res - end function ${name(rank, t1, k1)}$ - + end function ${RName}$ #:endfor #:endfor -#:def name(Rank, Type, Kind) -$:"mean_{0}_mask_all_{1}{2}_rdp".format(Rank, Type[0], Kind) -#:enddef #:for k1, t1 in INT_KINDS_TYPES #:for rank in RANKS - module function ${name(rank, t1, k1)}$(x, mask) result(res) + #:set RName = rname('mean_mask_all',rank, t1, k1, 'dp') + module function ${RName}$(x, mask) result(res) ${t1}$, intent(in) :: x${ranksuffix(rank)}$ logical, intent(in) :: mask${ranksuffix(rank)}$ real(dp) :: res - end function ${name(rank, t1, k1)}$ + end function ${RName}$ #:endfor #:endfor - -#:def name(Rank, Type, Kind) -$:"mean_{0}_mask_{1}{2}_{1}{2}".format(Rank, Type[0], Kind) -#:enddef #:for k1, t1 in RC_KINDS_TYPES #:for rank in RANKS - module function ${name(rank, t1, k1)}$(x, dim, mask) result(res) + #:set RName = rname('mean_mask',rank, t1, k1) + module function ${RName}$(x, dim, mask) result(res) ${t1}$, intent(in) :: x${ranksuffix(rank)}$ integer, intent(in) :: dim logical, intent(in) :: mask${ranksuffix(rank)}$ ${t1}$ :: res${reduced_shape('x', rank, 'dim')}$ - end function ${name(rank, t1, k1)}$ + end function ${RName}$ #:endfor #:endfor -#:def name(Rank, Type, Kind) -$:"mean_{0}_mask_{1}{2}_rdp".format(Rank, Type[0], Kind) -#:enddef #:for k1, t1 in INT_KINDS_TYPES #:for rank in RANKS - module function ${name(rank, t1, k1)}$(x, dim, mask) result(res) + #:set RName = rname('mean_mask',rank, t1, k1, 'dp') + module function ${RName}$(x, dim, mask) result(res) ${t1}$, intent(in) :: x${ranksuffix(rank)}$ integer, intent(in) :: dim logical, intent(in) :: mask${ranksuffix(rank)}$ real(dp) :: res${reduced_shape('x', rank, 'dim')}$ - end function ${name(rank, t1, k1)}$ + end function ${RName}$ #:endfor #:endfor diff --git a/src/stdlib_experimental_stats_mean.fypp b/src/stdlib_experimental_stats_mean.fypp index 1c79f9853..f12e4e336 100644 --- a/src/stdlib_experimental_stats_mean.fypp +++ b/src/stdlib_experimental_stats_mean.fypp @@ -10,12 +10,10 @@ submodule (stdlib_experimental_stats) stdlib_experimental_stats_mean contains -#:def name(Rank, Type, Kind) -$:"mean_{0}_all_{1}{2}_{1}{2}".format(Rank, Type[0], Kind) -#:enddef #:for k1, t1 in RC_KINDS_TYPES #:for rank in RANKS - module function ${name(rank, t1, k1)}$ (x, mask) result(res) + #:set RName = rname("mean_all",rank, t1, k1) + module function ${RName}$ (x, mask) result(res) ${t1}$, intent(in) :: x${ranksuffix(rank)}$ logical, intent(in), optional :: mask ${t1}$ :: res @@ -27,16 +25,14 @@ $:"mean_{0}_all_{1}{2}_{1}{2}".format(Rank, Type[0], Kind) res = sum(x) / real(size(x, kind = int64), ${k1}$) - end function ${name(rank, t1, k1)}$ + end function ${RName}$ #:endfor #:endfor -#:def name(Rank, Type, Kind) -$:"mean_{0}_all_{1}{2}_rdp".format(Rank, Type[0], Kind) -#:enddef #:for k1, t1 in INT_KINDS_TYPES #:for rank in RANKS - module function ${name(rank, t1, k1)}$(x, mask) result(res) + #:set RName = rname('mean_all', rank, t1, k1,'dp') + module function ${RName}$(x, mask) result(res) ${t1}$, intent(in) :: x${ranksuffix(rank)}$ logical, intent(in), optional :: mask real(dp) :: res @@ -48,17 +44,15 @@ $:"mean_{0}_all_{1}{2}_rdp".format(Rank, Type[0], Kind) res = sum(real(x, dp)) / real(size(x, kind = int64), dp) - end function ${name(rank, t1, k1)}$ + end function ${RName}$ #:endfor #:endfor -#:def name(Rank, Type, Kind) -$:"mean_{0}_{1}{2}_{1}{2}".format(Rank, Type[0], Kind) -#:enddef #:for k1, t1 in RC_KINDS_TYPES #:for rank in RANKS - module function ${name(rank, t1, k1)}$(x, dim, mask) result(res) + #:set RName = rname("mean",rank, t1, k1) + module function ${RName}$(x, dim, mask) result(res) ${t1}$, intent(in) :: x${ranksuffix(rank)}$ integer, intent(in) :: dim logical, intent(in), optional :: mask @@ -75,17 +69,15 @@ $:"mean_{0}_{1}{2}_{1}{2}".format(Rank, Type[0], Kind) call error_stop("ERROR (mean): wrong dimension") end if - end function ${name(rank, t1, k1)}$ + end function ${RName}$ #:endfor #:endfor -#:def name(Rank, Type, Kind) -$:"mean_{0}_{1}{2}_rdp".format(Rank, Type[0], Kind) -#:enddef #:for k1, t1 in INT_KINDS_TYPES #:for rank in RANKS - module function ${name(rank, t1, k1)}$(x, dim, mask) result(res) + #:set RName = rname("mean",rank, t1, k1,'dp') + module function ${RName}$(x, dim, mask) result(res) ${t1}$, intent(in) :: x${ranksuffix(rank)}$ integer, intent(in) :: dim logical, intent(in), optional :: mask @@ -102,49 +94,43 @@ $:"mean_{0}_{1}{2}_rdp".format(Rank, Type[0], Kind) call error_stop("ERROR (mean): wrong dimension") end if - end function ${name(rank, t1, k1)}$ + end function ${RName}$ #:endfor #:endfor -#:def name(Rank, Type, Kind) -$:"mean_{0}_mask_all_{1}{2}_{1}{2}".format(Rank, Type[0], Kind) -#:enddef #:for k1, t1 in RC_KINDS_TYPES #:for rank in RANKS - module function ${name(rank, t1, k1)}$(x, mask) result(res) + #:set RName = rname('mean_mask_all',rank, t1, k1) + module function ${RName}$(x, mask) result(res) ${t1}$, intent(in) :: x${ranksuffix(rank)}$ logical, intent(in) :: mask${ranksuffix(rank)}$ ${t1}$ :: res res = sum(x, mask) / real(count(mask, kind = int64), ${k1}$) - end function ${name(rank, t1, k1)}$ + end function ${RName}$ #:endfor #:endfor -#:def name(Rank, Type, Kind) -$:"mean_{0}_mask_all_{1}{2}_rdp".format(Rank, Type[0], Kind) -#:enddef #:for k1, t1 in INT_KINDS_TYPES #:for rank in RANKS - module function ${name(rank, t1, k1)}$(x, mask) result(res) + #:set RName = rname('mean_mask_all',rank, t1, k1, 'dp') + module function ${RName}$(x, mask) result(res) ${t1}$, intent(in) :: x${ranksuffix(rank)}$ logical, intent(in) :: mask${ranksuffix(rank)}$ real(dp) :: res res = sum(real(x, dp), mask) / real(count(mask, kind = int64), dp) - end function ${name(rank, t1, k1)}$ + end function ${RName}$ #:endfor #:endfor -#:def name(Rank, Type, Kind) -$:"mean_{0}_mask_{1}{2}_{1}{2}".format(Rank, Type[0], Kind) -#:enddef #:for k1, t1 in RC_KINDS_TYPES #:for rank in RANKS - module function ${name(rank, t1, k1)}$(x, dim, mask) result(res) + #:set RName = rname('mean_mask',rank, t1, k1) + module function ${RName}$(x, dim, mask) result(res) ${t1}$, intent(in) :: x${ranksuffix(rank)}$ integer, intent(in) :: dim logical, intent(in) :: mask${ranksuffix(rank)}$ @@ -156,17 +142,15 @@ $:"mean_{0}_mask_{1}{2}_{1}{2}".format(Rank, Type[0], Kind) call error_stop("ERROR (mean): wrong dimension") end if - end function ${name(rank, t1, k1)}$ + end function ${RName}$ #:endfor #:endfor -#:def name(Rank, Type, Kind) -$:"mean_{0}_mask_{1}{2}_rdp".format(Rank, Type[0], Kind) -#:enddef #:for k1, t1 in INT_KINDS_TYPES #:for rank in RANKS - module function ${name(rank, t1, k1)}$(x, dim, mask) result(res) + #:set RName = rname('mean_mask',rank, t1, k1, 'dp') + module function ${RName}$(x, dim, mask) result(res) ${t1}$, intent(in) :: x${ranksuffix(rank)}$ integer, intent(in) :: dim logical, intent(in) :: mask${ranksuffix(rank)}$ @@ -178,7 +162,7 @@ $:"mean_{0}_mask_{1}{2}_rdp".format(Rank, Type[0], Kind) call error_stop("ERROR (mean): wrong dimension") end if - end function ${name(rank, t1, k1)}$ + end function ${RName}$ #:endfor #:endfor diff --git a/src/tests/stats/test_mean_f03.f90 b/src/tests/stats/test_mean_f03.f90 index ce6650d4a..01434ffcb 100644 --- a/src/tests/stats/test_mean_f03.f90 +++ b/src/tests/stats/test_mean_f03.f90 @@ -5,6 +5,8 @@ program test_mean use stdlib_experimental_stats, only: mean implicit none +real(dp), parameter :: dptol = 2.2e-15_dp + real(dp), allocatable :: d(:, :) real(dp), allocatable :: d8(:, :, :, :, :, :, :, :) complex(dp), allocatable :: cd8(:, :, :, :, :, :, :, :) @@ -13,9 +15,9 @@ program test_mean !dp call loadtxt("array3.dat", d) -call assert( mean(d) - sum(d)/real(size(d), dp) == 0.0_dp) -call assert( sum( abs( mean(d,1) - sum(d,1)/real(size(d,1), dp) )) == 0.0_dp) -call assert( sum( abs( mean(d,2) - sum(d,2)/real(size(d,2), dp) )) == 0.0_dp) +call assert( mean(d) - sum(d)/real(size(d), dp) < dptol) +call assert( sum( abs( mean(d,1) - sum(d,1)/real(size(d,1), dp) )) < dptol) +call assert( sum( abs( mean(d,2) - sum(d,2)/real(size(d,2), dp) )) < dptol) !dp rank 8 allocate(d8(size(d,1), size(d,2), 3, 4, 5, 6, 7, 8)) @@ -23,31 +25,31 @@ program test_mean d8(:, :, 2, 4, 5 ,6 ,7 ,8)=d * 1.5_dp; d8(:, :, 3, 4, 5 ,6 ,7 ,8)=d * 4._dp; -call assert( mean(d8) - sum(d8)/real(size(d8), dp) == 0.0_dp) +call assert( mean(d8) - sum(d8)/real(size(d8), dp) < dptol) -call assert( sum( abs( mean(d8,1) - sum(d8,1)/real(size(d8,1), dp) )) == 0.0_dp) -call assert( sum( abs( mean(d8,2) - sum(d8,2)/real(size(d8,2), dp) )) == 0.0_dp) -call assert( sum( abs( mean(d8,3) - sum(d8,3)/real(size(d8,3), dp) )) == 0.0_dp) -call assert( sum( abs( mean(d8,4) - sum(d8,4)/real(size(d8,4), dp) )) == 0.0_dp) -call assert( sum( abs( mean(d8,5) - sum(d8,5)/real(size(d8,5), dp) )) == 0.0_dp) -call assert( sum( abs( mean(d8,6) - sum(d8,6)/real(size(d8,6), dp) )) == 0.0_dp) -call assert( sum( abs( mean(d8,7) - sum(d8,7)/real(size(d8,7), dp) )) == 0.0_dp) -call assert( sum( abs( mean(d8,8) - sum(d8,8)/real(size(d8,8), dp) )) == 0.0_dp) +call assert( sum( abs( mean(d8,1) - sum(d8,1)/real(size(d8,1), dp) )) < dptol) +call assert( sum( abs( mean(d8,2) - sum(d8,2)/real(size(d8,2), dp) )) < dptol) +call assert( sum( abs( mean(d8,3) - sum(d8,3)/real(size(d8,3), dp) )) < dptol) +call assert( sum( abs( mean(d8,4) - sum(d8,4)/real(size(d8,4), dp) )) < dptol) +call assert( sum( abs( mean(d8,5) - sum(d8,5)/real(size(d8,5), dp) )) < dptol) +call assert( sum( abs( mean(d8,6) - sum(d8,6)/real(size(d8,6), dp) )) < dptol) +call assert( sum( abs( mean(d8,7) - sum(d8,7)/real(size(d8,7), dp) )) < dptol) +call assert( sum( abs( mean(d8,8) - sum(d8,8)/real(size(d8,8), dp) )) < dptol) !cdp rank 8 allocate(cd8(size(d,1), size(d,2), 3, 4, 5, 6, 7, 8)) cd8 = cmplx(1._dp, 1._dp, kind=dp)*d8 -call assert( mean(cd8) - sum(cd8)/real(size(cd8), dp) == 0.0_dp) +call assert( abs(mean(cd8) - sum(cd8)/real(size(cd8), dp)) < dptol) -call assert( sum( abs( mean(cd8,1) - sum(cd8,1)/real(size(cd8,1), dp) )) == 0.0_dp) -call assert( sum( abs( mean(cd8,2) - sum(cd8,2)/real(size(cd8,2), dp) )) == 0.0_dp) -call assert( sum( abs( mean(cd8,3) - sum(cd8,3)/real(size(cd8,3), dp) )) == 0.0_dp) -call assert( sum( abs( mean(cd8,4) - sum(cd8,4)/real(size(cd8,4), dp) )) == 0.0_dp) -call assert( sum( abs( mean(cd8,5) - sum(cd8,5)/real(size(cd8,5), dp) )) == 0.0_dp) -call assert( sum( abs( mean(cd8,6) - sum(cd8,6)/real(size(cd8,6), dp) )) == 0.0_dp) -call assert( sum( abs( mean(cd8,7) - sum(cd8,7)/real(size(cd8,7), dp) )) == 0.0_dp) -call assert( sum( abs( mean(cd8,8) - sum(cd8,8)/real(size(cd8,8), dp) )) == 0.0_dp) +call assert( sum( abs( mean(cd8,1) - sum(cd8,1)/real(size(cd8,1), dp) )) < dptol) +call assert( sum( abs( mean(cd8,2) - sum(cd8,2)/real(size(cd8,2), dp) )) < dptol) +call assert( sum( abs( mean(cd8,3) - sum(cd8,3)/real(size(cd8,3), dp) )) < dptol) +call assert( sum( abs( mean(cd8,4) - sum(cd8,4)/real(size(cd8,4), dp) )) < dptol) +call assert( sum( abs( mean(cd8,5) - sum(cd8,5)/real(size(cd8,5), dp) )) < dptol) +call assert( sum( abs( mean(cd8,6) - sum(cd8,6)/real(size(cd8,6), dp) )) < dptol) +call assert( sum( abs( mean(cd8,7) - sum(cd8,7)/real(size(cd8,7), dp) )) < dptol) +call assert( sum( abs( mean(cd8,8) - sum(cd8,8)/real(size(cd8,8), dp) )) < dptol) contains end program From 828cda0133c2271eba784dd4d3588eec315e6b27 Mon Sep 17 00:00:00 2001 From: Juan Fiol Date: Tue, 4 Feb 2020 14:54:57 -0300 Subject: [PATCH 3/3] Changed parameters as suggested by @certik --- src/tests/stats/test_mean.f90 | 4 ++-- src/tests/stats/test_mean_f03.f90 | 2 +- 2 files changed, 3 insertions(+), 3 deletions(-) diff --git a/src/tests/stats/test_mean.f90 b/src/tests/stats/test_mean.f90 index e2f165124..69e816710 100644 --- a/src/tests/stats/test_mean.f90 +++ b/src/tests/stats/test_mean.f90 @@ -5,8 +5,8 @@ program test_mean use stdlib_experimental_stats, only: mean implicit none -real(sp), parameter :: sptol = 1.2e-06_sp -real(dp), parameter :: dptol = 2.2e-15_dp +real(sp), parameter :: sptol =10 * epsilon(1._sp) +real(dp), parameter :: dptol =10 * epsilon(1._dp) real(sp) :: s1(3) = [1.0_sp, 2.0_sp, 3.0_sp] diff --git a/src/tests/stats/test_mean_f03.f90 b/src/tests/stats/test_mean_f03.f90 index 01434ffcb..c1d9cd2e1 100644 --- a/src/tests/stats/test_mean_f03.f90 +++ b/src/tests/stats/test_mean_f03.f90 @@ -5,7 +5,7 @@ program test_mean use stdlib_experimental_stats, only: mean implicit none -real(dp), parameter :: dptol = 2.2e-15_dp +real(dp), parameter :: dptol =10 * epsilon(1._dp) real(dp), allocatable :: d(:, :) real(dp), allocatable :: d8(:, :, :, :, :, :, :, :)