From 286ddf624cc4e010d9cd067ff82a7fd9bc603cb1 Mon Sep 17 00:00:00 2001 From: zoziha Date: Sat, 28 Aug 2021 18:36:42 +0800 Subject: [PATCH 1/6] add `arg/argd/argpi` for `stdlib_math` --- doc/specs/stdlib_math.md | 131 ++++++++++++++++++++++++++++++- src/stdlib_math.fypp | 64 ++++++++++++++- src/tests/math/CMakeLists.txt | 3 +- src/tests/math/Makefile.manual | 3 +- src/tests/math/test_math_arg.f90 | 41 ++++++++++ 5 files changed, 238 insertions(+), 4 deletions(-) create mode 100644 src/tests/math/test_math_arg.f90 diff --git a/doc/specs/stdlib_math.md b/doc/specs/stdlib_math.md index 8d7fb516b..18908d6fd 100644 --- a/doc/specs/stdlib_math.md +++ b/doc/specs/stdlib_math.md @@ -342,4 +342,133 @@ program demo_math_arange print *, arange(0.0,2.0,0.0) !! [0.0,1.0,2.0]. Not recommended: `step` argument is zero! end program demo_math_arange -``` \ No newline at end of file +``` + +## `arg` + +### Status + +Experimental + +### Class + +Elemental function. + +### Description + +`arg` computes the phase angle (radian version) of `complex` scalar in the interval [-π,π]. +The angles in `theta` are such that `z = abs(z)*exp((0.0, theta))`. + +### Syntax + +`result = [[stdlib_math(module):arg(interface)]](z)` + +### Arguments + +`z`: Shall be a `complex` scalar/array. +This is an `intent(in)` argument. + +### Return value + +Returns the `real` type phase angle (radian version) of the `complex` argument `z`. + +#### Notes + +Although the angle of the complex number `0` is undefined, `arg((0,0))` returns the value `0`. + +### Example + +```fortran +program demo_math_arg + use stdlib_math, only: arg + print *, arg((0.0, 0.0)) !! 0.0 + print *, arg((3.0, 4.0)) !! 0.927 + print *, arg(2.0*exp((0.0, 0.5))) !! 0.5 +end program demo_math_arg +``` + +## `argd` + +### Status + +Experimental + +### Class + +Elemental function. + +### Description + +`argd` computes the phase angle (degree version) of `complex` scalar in the interval [-180.0,180.0]. +The angles in `theta` are such that `z = abs(z)*exp((0.0, theta*π/180.0))`. + +### Syntax + +`result = [[stdlib_math(module):argd(interface)]](z)` + +### Arguments + +`z`: Shall be a `complex` scalar/array. +This is an `intent(in)` argument. + +### Return value + +Returns the `real` type phase angle (degree version) of the `complex` argument `z`. + +#### Notes + +Although the angle of the complex number `0` is undefined, `argd((0,0))` returns the value `0`. + +### Example + +```fortran +program demo_math_argd + use stdlib_math, only: argd + print *, argd((0.0, 0.0)) !! 0.0 + print *, argd((3.0, 4.0)) !! 53.1° + print *, argd(2.0*exp((0.0, 0.5))) !! 28.64° +end program demo_math_argd +``` + +## `argpi` + +### Status + +Experimental + +### Class + +Elemental function. + +### Description + +`argpi` computes the phase angle (circular version) of `complex` scalar in the interval [-1.0,1.0]. +The angles in `theta` are such that `z = abs(z)*exp((0.0, theta*π))`. + +### Syntax + +`result = [[stdlib_math(module):argpi(interface)]](z)` + +### Arguments + +`z`: Shall be a `complex` scalar/array. +This is an `intent(in)` argument. + +### Return value + +Returns the `real` type phase angle (circular version) of the `complex` argument `z`. + +#### Notes + +Although the angle of the complex number `0` is undefined, `argpi((0,0))` returns the value `0`. + +### Example + +```fortran +program demo_math_argpi + use stdlib_math, only: argpi + print *, argpi((0.0, 0.0)) !! 0.0 + print *, argpi((3.0, 4.0)) !! 0.295 + print *, argpi(2.0*exp((0.0, 0.5))) !! 0.159 +end program demo_math_argpi +``` diff --git a/src/stdlib_math.fypp b/src/stdlib_math.fypp index 2ca0f543d..793d3081a 100644 --- a/src/stdlib_math.fypp +++ b/src/stdlib_math.fypp @@ -11,7 +11,7 @@ module stdlib_math public :: clip, linspace, logspace public :: EULERS_NUMBER_SP, EULERS_NUMBER_DP, EULERS_NUMBER_QP public :: DEFAULT_LINSPACE_LENGTH, DEFAULT_LOGSPACE_BASE, DEFAULT_LOGSPACE_LENGTH - public :: arange + public :: arange, arg, argd, argpi integer, parameter :: DEFAULT_LINSPACE_LENGTH = 100 integer, parameter :: DEFAULT_LOGSPACE_LENGTH = 50 @@ -22,6 +22,11 @@ module stdlib_math real(dp), parameter :: EULERS_NUMBER_DP = exp(1.0_dp) real(qp), parameter :: EULERS_NUMBER_QP = exp(1.0_qp) + !> Useful constants `PI` for `argd/argpi` + #:for k1 in REAL_KINDS + real(kind=${k1}$), parameter :: PI_${k1}$ = 4.0_${k1}$*atan(1.0_${k1}$) + #:endfor + interface clip #:for k1, t1 in IR_KINDS_TYPES module procedure clip_${k1}$ @@ -279,6 +284,36 @@ module stdlib_math #:endfor end interface arange + !> Version: experimental + !> + !> `arg` computes the phase angle in the interval [-π,π]. + !> ([Specification](../page/specs/stdlib_math.html#arg)) + interface arg + #:for k1, t1 in CMPLX_KINDS_TYPES + procedure :: arg_${k1}$ + #:endfor + end interface arg + + !> Version: experimental + !> + !> `argd` computes the phase angle of degree version in the interval [-180.0,180.0]. + !> ([Specification](../page/specs/stdlib_math.html#argd)) + interface argd + #:for k1, t1 in CMPLX_KINDS_TYPES + procedure :: argd_${k1}$ + #:endfor + end interface argd + + !> Version: experimental + !> + !> `argpi` computes the phase angle of IEEE circular version in the interval [-1.0,1.0]. + !> ([Specification](../page/specs/stdlib_math.html#argpi)) + interface argpi + #:for k1, t1 in CMPLX_KINDS_TYPES + procedure :: argpi_${k1}$ + #:endfor + end interface argpi + contains #:for k1, t1 in IR_KINDS_TYPES @@ -292,4 +327,31 @@ contains end function clip_${k1}$ #:endfor + + #:for k1, t1 in CMPLX_KINDS_TYPES + elemental function arg_${k1}$(z) result(result) + ${t1}$, intent(in) :: z + real(${k1}$) :: result + + result = aimag(log(z)) + + end function arg_${k1}$ + + elemental function argd_${k1}$(z) result(result) + ${t1}$, intent(in) :: z + real(${k1}$) :: result + + result = aimag(log(z))*180.0_${k1}$/PI_${k1}$ + + end function argd_${k1}$ + + elemental function argpi_${k1}$(z) result(result) + ${t1}$, intent(in) :: z + real(${k1}$) :: result + + result = aimag(log(z))/PI_${k1}$ + + end function argpi_${k1}$ + #:endfor + end module stdlib_math diff --git a/src/tests/math/CMakeLists.txt b/src/tests/math/CMakeLists.txt index 0bcaf1e1a..98cb9a762 100644 --- a/src/tests/math/CMakeLists.txt +++ b/src/tests/math/CMakeLists.txt @@ -1,4 +1,5 @@ ADDTEST(stdlib_math) ADDTEST(linspace) ADDTEST(logspace) -ADDTEST(math_arange) \ No newline at end of file +ADDTEST(math_arange) +ADDTEST(math_arg) \ No newline at end of file diff --git a/src/tests/math/Makefile.manual b/src/tests/math/Makefile.manual index f11cbf7a4..0796545d6 100644 --- a/src/tests/math/Makefile.manual +++ b/src/tests/math/Makefile.manual @@ -1,5 +1,6 @@ PROGS_SRC = test_stdlib_math.f90 test_linspace.f90 test_logspace.f90 \ - test_math_arange.f90 + test_math_arange.f90 \ + test_math_arg.f90 include ../Makefile.manual.test.mk diff --git a/src/tests/math/test_math_arg.f90 b/src/tests/math/test_math_arg.f90 new file mode 100644 index 000000000..c4836c9a7 --- /dev/null +++ b/src/tests/math/test_math_arg.f90 @@ -0,0 +1,41 @@ +program tester + + real :: tol = 1.0e-3 + call test_math_arg_complex + call test_math_argd_complex + call test_math_argpi_complex + print *, "All tests in `test_math_arg` passed." + +contains + + subroutine test_math_arg_complex + use stdlib_math, only: arg + use stdlib_error, only: check + + call check(abs(arg(2.0*exp((0.0e0, 0.5))) - 0.5) < tol, msg="arg(2.0*exp((0.0, 0.5))) failed.") + call check(abs(arg((1.7552, 0.9589)) - 0.5) < tol, msg="arg((1.7552, 0.9589)) failed.") + call check(abs(arg((0.0, 0.0)) - 0.0) < tol, msg="arg((0.0, 0.0)) failed.") + + end subroutine test_math_arg_complex + + subroutine test_math_argd_complex + use stdlib_math, only: argd + use stdlib_error, only: check + + call check(abs(argd(2.0*exp((0.0, 0.5))) - 28.648) < tol, msg="argd(2.0*exp((0.0, 0.5))) failed.") + call check(abs(argd((1.7552, 0.9589)) - 28.648) < tol, msg="argd((1.7552, 0.9589)) failed.") + call check(abs(argd((0.0, 0.0)) - 0.0) < tol, msg="argd((0.0, 0.0)) failed.") + + end subroutine test_math_argd_complex + + subroutine test_math_argpi_complex + use stdlib_math, only: argpi + use stdlib_error, only: check + + call check(abs(argpi(2.0*exp((0.0, 0.5))) - 0.159) < tol, msg="argpi(2.0*exp((0.0, 0.5))) failed.") + call check(abs(argpi((1.7552, 0.9589)) - 0.159) < tol, msg="argpi((1.7552, 0.9589)) failed.") + call check(abs(argpi((0.0, 0.0)) - 0.0) < tol, msg="argpi((0.0, 0.0)) failed.") + + end subroutine test_math_argpi_complex + +end program tester From feca677eba8d88a2c5a14b020878330c6672380c Mon Sep 17 00:00:00 2001 From: zoziha Date: Sun, 19 Sep 2021 10:16:08 +0800 Subject: [PATCH 2/6] Tighter tests for `arg/argd/argpi`. --- doc/specs/stdlib_math.md | 6 +++--- src/tests/math/test_math_arg.f90 | 36 ++++++++++++++++++++------------ 2 files changed, 26 insertions(+), 16 deletions(-) diff --git a/doc/specs/stdlib_math.md b/doc/specs/stdlib_math.md index 18908d6fd..d3157ba3b 100644 --- a/doc/specs/stdlib_math.md +++ b/doc/specs/stdlib_math.md @@ -356,7 +356,7 @@ Elemental function. ### Description -`arg` computes the phase angle (radian version) of `complex` scalar in the interval [-π,π]. +`arg` computes the phase angle (radian version) of `complex` scalar in the interval (-π,π]. The angles in `theta` are such that `z = abs(z)*exp((0.0, theta))`. ### Syntax @@ -399,7 +399,7 @@ Elemental function. ### Description -`argd` computes the phase angle (degree version) of `complex` scalar in the interval [-180.0,180.0]. +`argd` computes the phase angle (degree version) of `complex` scalar in the interval (-180.0,180.0]. The angles in `theta` are such that `z = abs(z)*exp((0.0, theta*π/180.0))`. ### Syntax @@ -442,7 +442,7 @@ Elemental function. ### Description -`argpi` computes the phase angle (circular version) of `complex` scalar in the interval [-1.0,1.0]. +`argpi` computes the phase angle (circular version) of `complex` scalar in the interval (-1.0,1.0]. The angles in `theta` are such that `z = abs(z)*exp((0.0, theta*π))`. ### Syntax diff --git a/src/tests/math/test_math_arg.f90 b/src/tests/math/test_math_arg.f90 index c4836c9a7..8649bae66 100644 --- a/src/tests/math/test_math_arg.f90 +++ b/src/tests/math/test_math_arg.f90 @@ -1,6 +1,6 @@ program tester - real :: tol = 1.0e-3 + real :: tol = epsilon(1.0) call test_math_arg_complex call test_math_argd_complex call test_math_argpi_complex @@ -9,32 +9,42 @@ program tester contains subroutine test_math_arg_complex - use stdlib_math, only: arg + use stdlib_math, only: arg use stdlib_error, only: check - call check(abs(arg(2.0*exp((0.0e0, 0.5))) - 0.5) < tol, msg="arg(2.0*exp((0.0, 0.5))) failed.") - call check(abs(arg((1.7552, 0.9589)) - 0.5) < tol, msg="arg((1.7552, 0.9589)) failed.") - call check(abs(arg((0.0, 0.0)) - 0.0) < tol, msg="arg((0.0, 0.0)) failed.") + print *, 2.0*exp((0.0e0, 0.5)) !! (1.75516510, 0.958851099) + + call check(abs(arg(2.0*exp((0.0e0, 0.5))) - 0.5) < tol, msg="arg(2.0*exp((0.0e0, 0.5))) - 0.5) failed.") + call check(abs(arg((1.75516510, 0.958851099)) - 0.5) < tol, msg="arg((1.75516510, 0.958851099)) - 0.5) failed.") + call check(abs(arg((0.0, 0.0)) - 0.0) < tol, msg="arg((0.0, 0.0)) failed.") end subroutine test_math_arg_complex subroutine test_math_argd_complex - use stdlib_math, only: argd + use stdlib_math, only: argd use stdlib_error, only: check - call check(abs(argd(2.0*exp((0.0, 0.5))) - 28.648) < tol, msg="argd(2.0*exp((0.0, 0.5))) failed.") - call check(abs(argd((1.7552, 0.9589)) - 28.648) < tol, msg="argd((1.7552, 0.9589)) failed.") - call check(abs(argd((0.0, 0.0)) - 0.0) < tol, msg="argd((0.0, 0.0)) failed.") + call check(abs(argd(2.0*exp((0.0, 0.5))) - 28.6478882) < tol, msg="argd(2.0*exp((0.0, 0.5))) - 28.6478882) failed.") + call check(abs(argd((1.75516510, 0.958851099)) - 28.6478882) < tol, & + msg="argd((1.75516510, 0.958851099)) - 28.6478882) failed.") + call check(abs(argd((0.0, 0.0)) - 0.0) < tol, msg="argd((0.0, 0.0)) failed.") end subroutine test_math_argd_complex subroutine test_math_argpi_complex - use stdlib_math, only: argpi + use stdlib_math, only: argpi use stdlib_error, only: check + real, parameter :: PI = 4*atan(1.0) + + !> Power exponent calculation will introduce calculation errors: 2.0*exp((0.0, PI)) /= (-2.0, 0.0). + print *, 2.0*exp((0.0, PI)), 2.0*exp(cmplx(0.0, -PI)) !! (-2.00000000,-1.748455531E-07), (-2.00000000,1.748455531E-07) + print *, argpi(2.0*exp((0.0, PI))), argpi(2.0*exp(cmplx(0.0, -PI))) !! -0.999999940, 0.999999940 + + call check(abs(argpi(2.0*exp((0.0, PI))) - (- 1.0)) < tol, msg="argpi(2.0*exp((0.0, PI))) - (- 1.0)) failed.") + call check(abs(argpi(2.0*exp(cmplx(0.0, -PI))) - 1.0) < tol, msg="argpi(2.0*exp(cmplx(0.0, -PI))) - 1.0) failed.") - call check(abs(argpi(2.0*exp((0.0, 0.5))) - 0.159) < tol, msg="argpi(2.0*exp((0.0, 0.5))) failed.") - call check(abs(argpi((1.7552, 0.9589)) - 0.159) < tol, msg="argpi((1.7552, 0.9589)) failed.") - call check(abs(argpi((0.0, 0.0)) - 0.0) < tol, msg="argpi((0.0, 0.0)) failed.") + call check(abs(argpi((-2.0, 0.0)) - 1.0) < tol, msg="argpi((-2.0, 0.0)) - 1.0) failed.") !! (-1.0, 1.0] + call check(abs(argpi(( 0.0, 0.0)) - 0.0) < tol, msg="argpi((0.0, 0.0)) failed.") end subroutine test_math_argpi_complex From 6e00619677ea1c1fff60542707e3a160a44f033f Mon Sep 17 00:00:00 2001 From: zoziha Date: Fri, 10 Dec 2021 16:16:17 +0800 Subject: [PATCH 3/6] Update `test-drive` style tests for arg* routines. --- src/stdlib_math.fypp | 12 ++--- src/tests/math/test_stdlib_math.fypp | 71 +++++++++++++++++++++++++++- 2 files changed, 76 insertions(+), 7 deletions(-) diff --git a/src/stdlib_math.fypp b/src/stdlib_math.fypp index 60ee81701..d7d01fd35 100644 --- a/src/stdlib_math.fypp +++ b/src/stdlib_math.fypp @@ -301,30 +301,30 @@ module stdlib_math !> Version: experimental !> - !> `arg` computes the phase angle in the interval [-π,π]. + !> `arg` computes the phase angle in the interval (-π,π]. !> ([Specification](../page/specs/stdlib_math.html#arg)) interface arg - #:for k1, t1 in CMPLX_KINDS_TYPES + #:for k1 in CMPLX_KINDS procedure :: arg_${k1}$ #:endfor end interface arg !> Version: experimental !> - !> `argd` computes the phase angle of degree version in the interval [-180.0,180.0]. + !> `argd` computes the phase angle of degree version in the interval (-180.0,180.0]. !> ([Specification](../page/specs/stdlib_math.html#argd)) interface argd - #:for k1, t1 in CMPLX_KINDS_TYPES + #:for k1 in CMPLX_KINDS procedure :: argd_${k1}$ #:endfor end interface argd !> Version: experimental !> - !> `argpi` computes the phase angle of IEEE circular version in the interval [-1.0,1.0]. + !> `argpi` computes the phase angle of IEEE circular version in the interval (-1.0,1.0]. !> ([Specification](../page/specs/stdlib_math.html#argpi)) interface argpi - #:for k1, t1 in CMPLX_KINDS_TYPES + #:for k1 in CMPLX_KINDS procedure :: argpi_${k1}$ #:endfor end interface argpi diff --git a/src/tests/math/test_stdlib_math.fypp b/src/tests/math/test_stdlib_math.fypp index 651fea610..82ced0cab 100644 --- a/src/tests/math/test_stdlib_math.fypp +++ b/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 + use stdlib_math, only: clip, arg, argd, argpi use stdlib_kinds, only: int8, int16, int32, int64, sp, dp, xdp, qp implicit none @@ -32,6 +32,13 @@ contains new_unittest("clip-real-double-bounds", test_clip_rdp_bounds), & new_unittest("clip-real-quad", test_clip_rqp), & new_unittest("clip-real-quad-bounds", test_clip_rqp_bounds) & + + !> Tests for arg/argd/argpi + #:for k1 in CMPLX_KINDS + , new_unittest("arg-cmplx-${k1}$", test_arg_${k1}$) & + , new_unittest("argd-cmplx-${k1}$", test_argd_${k1}$) & + , new_unittest("argpi-cmplx-${k1}$", test_argpi_${k1}$) & + #:endfor ] end subroutine collect_stdlib_math @@ -203,6 +210,68 @@ contains #:endif end subroutine test_clip_rqp_bounds + + #:for k1 in CMPLX_KINDS + subroutine test_arg_${k1}$(error) + type(error_type), allocatable, intent(out) :: error + real(${k1}$), parameter :: tol = epsilon(1.0_${k1}$) + + #! For scalar + call check(error, abs(arg(2*exp((0.0_${k1}$, 0.5_${k1}$))) - 0.5_${k1}$) < tol, & + "test_nonzero_scalar") + if (allocated(error)) return + call check(error, abs(arg((0.0_${k1}$, 0.0_${k1}$)) - 0.0_${k1}$) < tol, & + "test_zero_scalar") + + #! and for array + call check(error, all(abs(arg([2*exp((0.0_${k1}$, 0.5_${k1}$))]) - [0.5_${k1}$]) < [tol]), & + "test_nonzero_array") + if (allocated(error)) return + call check(error, all(abs(arg([(0.0_${k1}$, 0.0_${k1}$)]) - [0.0_${k1}$]) < [tol]), & + "test_zero_array") + + end subroutine test_arg_${k1}$ + + subroutine test_argd_${k1}$(error) + type(error_type), allocatable, intent(out) :: error + real(${k1}$), parameter :: tol = epsilon(1.0_${k1}$) + + #! For scalar + call check(error, abs(argd((-1.0_${k1}$, 0.0_${k1}$)) - 180.0_${k1}$) < tol, & + "test_nonzero_scalar") + if (allocated(error)) return + call check(error, abs(argd((0.0_${k1}$, 0.0_${k1}$)) - 0.0_${k1}$) < tol, & + "test_zero_scalar") + + #! and for array + call check(error, all(abs(argd([(-1.0_${k1}$, 0.0_${k1}$)]) - [180.0_${k1}$]) < [tol]), & + "test_nonzero_array") + if (allocated(error)) return + call check(error, all(abs(argd([(0.0_${k1}$, 0.0_${k1}$)]) - [0.0_${k1}$]) < [tol]), & + "test_zero_array") + + end subroutine test_argd_${k1}$ + + subroutine test_argpi_${k1}$(error) + type(error_type), allocatable, intent(out) :: error + real(${k1}$), parameter :: tol = epsilon(1.0_${k1}$) + + #! For scalar + call check(error, abs(argpi((-1.0_${k1}$, 0.0_${k1}$)) - 1.0_${k1}$) < tol, & + "test_nonzero_scalar") + if (allocated(error)) return + call check(error, abs(argpi((0.0_${k1}$, 0.0_${k1}$)) - 0.0_${k1}$) < tol, & + "test_zero_scalar") + + #! and for array + call check(error, all(abs(argpi([(-1.0_${k1}$, 0.0_${k1}$)]) - [1.0_${k1}$]) < [tol]), & + "test_nonzero_array") + if (allocated(error)) return + call check(error, all(abs(argpi([(0.0_${k1}$, 0.0_${k1}$)]) - [0.0_${k1}$]) < [tol]), & + "test_zero_array") + + end subroutine test_argpi_${k1}$ + #:endfor end module test_stdlib_math From ec04eca8660588a5411bccf9d9cf3a628c4f867e Mon Sep 17 00:00:00 2001 From: zoziha Date: Fri, 10 Dec 2021 16:24:34 +0800 Subject: [PATCH 4/6] Remove old tests for arg* routines. --- src/tests/math/CMakeLists.txt | 1 - src/tests/math/Makefile.manual | 1 - src/tests/math/test_math_arg.f90 | 51 -------------------------------- 3 files changed, 53 deletions(-) delete mode 100644 src/tests/math/test_math_arg.f90 diff --git a/src/tests/math/CMakeLists.txt b/src/tests/math/CMakeLists.txt index ef03af7f3..9d11bf765 100644 --- a/src/tests/math/CMakeLists.txt +++ b/src/tests/math/CMakeLists.txt @@ -8,4 +8,3 @@ ADDTEST(stdlib_math) ADDTEST(linspace) ADDTEST(logspace) ADDTEST(math_arange) -ADDTEST(math_arg) diff --git a/src/tests/math/Makefile.manual b/src/tests/math/Makefile.manual index 49a71edb6..5da73c5da 100644 --- a/src/tests/math/Makefile.manual +++ b/src/tests/math/Makefile.manual @@ -4,7 +4,6 @@ SRCGEN = $(SRCFYPP:.fypp=.f90) PROGS_SRC = test_linspace.f90 test_logspace.f90 \ test_math_arange.f90 \ - test_math_arg.f90 \ $(SRCGEN) $(SRCGEN): %.f90: %.fypp ../../common.fypp diff --git a/src/tests/math/test_math_arg.f90 b/src/tests/math/test_math_arg.f90 deleted file mode 100644 index 8649bae66..000000000 --- a/src/tests/math/test_math_arg.f90 +++ /dev/null @@ -1,51 +0,0 @@ -program tester - - real :: tol = epsilon(1.0) - call test_math_arg_complex - call test_math_argd_complex - call test_math_argpi_complex - print *, "All tests in `test_math_arg` passed." - -contains - - subroutine test_math_arg_complex - use stdlib_math, only: arg - use stdlib_error, only: check - - print *, 2.0*exp((0.0e0, 0.5)) !! (1.75516510, 0.958851099) - - call check(abs(arg(2.0*exp((0.0e0, 0.5))) - 0.5) < tol, msg="arg(2.0*exp((0.0e0, 0.5))) - 0.5) failed.") - call check(abs(arg((1.75516510, 0.958851099)) - 0.5) < tol, msg="arg((1.75516510, 0.958851099)) - 0.5) failed.") - call check(abs(arg((0.0, 0.0)) - 0.0) < tol, msg="arg((0.0, 0.0)) failed.") - - end subroutine test_math_arg_complex - - subroutine test_math_argd_complex - use stdlib_math, only: argd - use stdlib_error, only: check - - call check(abs(argd(2.0*exp((0.0, 0.5))) - 28.6478882) < tol, msg="argd(2.0*exp((0.0, 0.5))) - 28.6478882) failed.") - call check(abs(argd((1.75516510, 0.958851099)) - 28.6478882) < tol, & - msg="argd((1.75516510, 0.958851099)) - 28.6478882) failed.") - call check(abs(argd((0.0, 0.0)) - 0.0) < tol, msg="argd((0.0, 0.0)) failed.") - - end subroutine test_math_argd_complex - - subroutine test_math_argpi_complex - use stdlib_math, only: argpi - use stdlib_error, only: check - real, parameter :: PI = 4*atan(1.0) - - !> Power exponent calculation will introduce calculation errors: 2.0*exp((0.0, PI)) /= (-2.0, 0.0). - print *, 2.0*exp((0.0, PI)), 2.0*exp(cmplx(0.0, -PI)) !! (-2.00000000,-1.748455531E-07), (-2.00000000,1.748455531E-07) - print *, argpi(2.0*exp((0.0, PI))), argpi(2.0*exp(cmplx(0.0, -PI))) !! -0.999999940, 0.999999940 - - call check(abs(argpi(2.0*exp((0.0, PI))) - (- 1.0)) < tol, msg="argpi(2.0*exp((0.0, PI))) - (- 1.0)) failed.") - call check(abs(argpi(2.0*exp(cmplx(0.0, -PI))) - 1.0) < tol, msg="argpi(2.0*exp(cmplx(0.0, -PI))) - 1.0) failed.") - - call check(abs(argpi((-2.0, 0.0)) - 1.0) < tol, msg="argpi((-2.0, 0.0)) - 1.0) failed.") !! (-1.0, 1.0] - call check(abs(argpi(( 0.0, 0.0)) - 0.0) < tol, msg="argpi((0.0, 0.0)) failed.") - - end subroutine test_math_argpi_complex - -end program tester From c37c0b2b75c5d0f485c83d9b753629c0b592d02c Mon Sep 17 00:00:00 2001 From: zoziha Date: Tue, 14 Dec 2021 20:28:57 +0800 Subject: [PATCH 5/6] arg**: `aimag(log(..)) -> atan2(..)`, update tests. --- doc/specs/stdlib_math.md | 26 +++++++--------- src/stdlib_math.fypp | 12 ++++---- src/tests/math/test_stdlib_math.fypp | 45 ++++++++++++++-------------- 3 files changed, 40 insertions(+), 43 deletions(-) diff --git a/doc/specs/stdlib_math.md b/doc/specs/stdlib_math.md index 46b3f1bf5..d5ba100b2 100644 --- a/doc/specs/stdlib_math.md +++ b/doc/specs/stdlib_math.md @@ -389,7 +389,7 @@ program demo_math_arange end program demo_math_arange ``` -## `arg` +## `arg` - Computes the phase angle in radian of a complex scalar ### Status @@ -402,7 +402,7 @@ Elemental function. ### Description `arg` computes the phase angle (radian version) of `complex` scalar in the interval (-π,π]. -The angles in `theta` are such that `z = abs(z)*exp((0.0, theta))`. +The angles in `θ` are such that `z = abs(z)*exp((0.0, θ))`. ### Syntax @@ -417,9 +417,7 @@ This is an `intent(in)` argument. Returns the `real` type phase angle (radian version) of the `complex` argument `z`. -#### Notes - -Although the angle of the complex number `0` is undefined, `arg((0,0))` returns the value `0`. +Notes: Although the angle of the complex number `0` is undefined, `arg((0,0))` returns the value `0`. ### Example @@ -432,7 +430,7 @@ program demo_math_arg end program demo_math_arg ``` -## `argd` +## `argd` - Computes the phase angle in degree of a complex scalar ### Status @@ -445,7 +443,7 @@ Elemental function. ### Description `argd` computes the phase angle (degree version) of `complex` scalar in the interval (-180.0,180.0]. -The angles in `theta` are such that `z = abs(z)*exp((0.0, theta*π/180.0))`. +The angles in `θ` are such that `z = abs(z)*exp((0.0, θ*π/180.0))`. ### Syntax @@ -460,9 +458,7 @@ This is an `intent(in)` argument. Returns the `real` type phase angle (degree version) of the `complex` argument `z`. -#### Notes - -Although the angle of the complex number `0` is undefined, `argd((0,0))` returns the value `0`. +Notes: Although the angle of the complex number `0` is undefined, `argd((0,0))` returns the value `0`. ### Example @@ -475,7 +471,7 @@ program demo_math_argd end program demo_math_argd ``` -## `argpi` +## `argpi` - Computes the phase angle in circular of a complex scalar ### Status @@ -487,8 +483,8 @@ Elemental function. ### Description -`argpi` computes the phase angle (circular version) of `complex` scalar in the interval (-1.0,1.0]. -The angles in `theta` are such that `z = abs(z)*exp((0.0, theta*π))`. +`argpi` computes the phase angle (IEEE circular version) of `complex` scalar in the interval (-1.0,1.0]. +The angles in `θ` are such that `z = abs(z)*exp((0.0, θ*π))`. ### Syntax @@ -503,9 +499,7 @@ This is an `intent(in)` argument. Returns the `real` type phase angle (circular version) of the `complex` argument `z`. -#### Notes - -Although the angle of the complex number `0` is undefined, `argpi((0,0))` returns the value `0`. +Notes: Although the angle of the complex number `0` is undefined, `argpi((0,0))` returns the value `0`. ### Example diff --git a/src/stdlib_math.fypp b/src/stdlib_math.fypp index d7d01fd35..3b11771d1 100644 --- a/src/stdlib_math.fypp +++ b/src/stdlib_math.fypp @@ -29,7 +29,7 @@ module stdlib_math !> Useful constants `PI` for `argd/argpi` #:for k1 in REAL_KINDS - real(kind=${k1}$), parameter :: PI_${k1}$ = 4.0_${k1}$*atan(1.0_${k1}$) + real(kind=${k1}$), parameter :: PI_${k1}$ = acos(-1.0_${k1}$) #:endfor interface clip @@ -321,7 +321,7 @@ module stdlib_math !> Version: experimental !> - !> `argpi` computes the phase angle of IEEE circular version in the interval (-1.0,1.0]. + !> `argpi` computes the phase angle of circular version in the interval (-1.0,1.0]. !> ([Specification](../page/specs/stdlib_math.html#argpi)) interface argpi #:for k1 in CMPLX_KINDS @@ -348,7 +348,7 @@ contains ${t1}$, intent(in) :: z real(${k1}$) :: result - result = aimag(log(z)) + result = merge(0.0_${k1}$, atan2(z%im, z%re), z == (0.0_${k1}$, 0.0_${k1}$)) end function arg_${k1}$ @@ -356,7 +356,8 @@ contains ${t1}$, intent(in) :: z real(${k1}$) :: result - result = aimag(log(z))*180.0_${k1}$/PI_${k1}$ + result = merge(0.0_${k1}$, atan2(z%im, z%re), z == (0.0_${k1}$, 0.0_${k1}$)) & + *180.0_${k1}$/PI_${k1}$ end function argd_${k1}$ @@ -364,7 +365,8 @@ contains ${t1}$, intent(in) :: z real(${k1}$) :: result - result = aimag(log(z))/PI_${k1}$ + result = merge(0.0_${k1}$, atan2(z%im, z%re), z == (0.0_${k1}$, 0.0_${k1}$)) & + /PI_${k1}$ end function argpi_${k1}$ #:endfor diff --git a/src/tests/math/test_stdlib_math.fypp b/src/tests/math/test_stdlib_math.fypp index 82ced0cab..04b0ec5c3 100644 --- a/src/tests/math/test_stdlib_math.fypp +++ b/src/tests/math/test_stdlib_math.fypp @@ -4,11 +4,15 @@ module test_stdlib_math use testdrive, only : new_unittest, unittest_type, error_type, check, skip_test - use stdlib_math, only: clip, arg, argd, argpi + use stdlib_math, only: clip, arg, argd, argpi, arange use stdlib_kinds, only: int8, int16, int32, int64, sp, dp, xdp, qp implicit none public :: collect_stdlib_math + + #:for k1 in REAL_KINDS + real(kind=${k1}$), parameter :: PI_${k1}$ = acos(-1.0_${k1}$) + #:endfor contains @@ -214,7 +218,8 @@ contains #:for k1 in CMPLX_KINDS subroutine test_arg_${k1}$(error) type(error_type), allocatable, intent(out) :: error - real(${k1}$), parameter :: tol = epsilon(1.0_${k1}$) + real(${k1}$), parameter :: tol = sqrt(epsilon(1.0_${k1}$)) + real(${k1}$), allocatable :: theta(:) #! For scalar call check(error, abs(arg(2*exp((0.0_${k1}$, 0.5_${k1}$))) - 0.5_${k1}$) < tol, & @@ -223,18 +228,17 @@ contains call check(error, abs(arg((0.0_${k1}$, 0.0_${k1}$)) - 0.0_${k1}$) < tol, & "test_zero_scalar") - #! and for array - call check(error, all(abs(arg([2*exp((0.0_${k1}$, 0.5_${k1}$))]) - [0.5_${k1}$]) < [tol]), & - "test_nonzero_array") - if (allocated(error)) return - call check(error, all(abs(arg([(0.0_${k1}$, 0.0_${k1}$)]) - [0.0_${k1}$]) < [tol]), & - "test_zero_array") + #! and for array (180.0° see scalar version) + theta = arange(-179.0_${k1}$, 179.0_${k1}$, 3.58_${k1}$) + call check(error, all(abs(arg(exp(cmplx(0.0_${k1}$, theta/180*PI_${k1}$, ${k1}$))) - theta/180*PI_${k1}$) < tol), & + "test_array") end subroutine test_arg_${k1}$ subroutine test_argd_${k1}$(error) type(error_type), allocatable, intent(out) :: error - real(${k1}$), parameter :: tol = epsilon(1.0_${k1}$) + real(${k1}$), parameter :: tol = sqrt(epsilon(1.0_${k1}$)) + real(${k1}$), allocatable :: theta(:) #! For scalar call check(error, abs(argd((-1.0_${k1}$, 0.0_${k1}$)) - 180.0_${k1}$) < tol, & @@ -243,18 +247,17 @@ contains call check(error, abs(argd((0.0_${k1}$, 0.0_${k1}$)) - 0.0_${k1}$) < tol, & "test_zero_scalar") - #! and for array - call check(error, all(abs(argd([(-1.0_${k1}$, 0.0_${k1}$)]) - [180.0_${k1}$]) < [tol]), & - "test_nonzero_array") - if (allocated(error)) return - call check(error, all(abs(argd([(0.0_${k1}$, 0.0_${k1}$)]) - [0.0_${k1}$]) < [tol]), & - "test_zero_array") + #! and for array (180.0° see scalar version) + theta = arange(-179.0_${k1}$, 179.0_${k1}$, 3.58_${k1}$) + call check(error, all(abs(argd(exp(cmplx(0.0_${k1}$, theta/180*PI_${k1}$, ${k1}$))) - theta) < tol), & + "test_array") end subroutine test_argd_${k1}$ subroutine test_argpi_${k1}$(error) type(error_type), allocatable, intent(out) :: error - real(${k1}$), parameter :: tol = epsilon(1.0_${k1}$) + real(${k1}$), parameter :: tol = sqrt(epsilon(1.0_${k1}$)) + real(${k1}$), allocatable :: theta(:) #! For scalar call check(error, abs(argpi((-1.0_${k1}$, 0.0_${k1}$)) - 1.0_${k1}$) < tol, & @@ -263,12 +266,10 @@ contains call check(error, abs(argpi((0.0_${k1}$, 0.0_${k1}$)) - 0.0_${k1}$) < tol, & "test_zero_scalar") - #! and for array - call check(error, all(abs(argpi([(-1.0_${k1}$, 0.0_${k1}$)]) - [1.0_${k1}$]) < [tol]), & - "test_nonzero_array") - if (allocated(error)) return - call check(error, all(abs(argpi([(0.0_${k1}$, 0.0_${k1}$)]) - [0.0_${k1}$]) < [tol]), & - "test_zero_array") + #! and for array (180.0° see scalar version) + theta = arange(-179.0_${k1}$, 179.0_${k1}$, 3.58_${k1}$) + call check(error, all(abs(argpi(exp(cmplx(0.0_${k1}$, theta/180*PI_${k1}$, ${k1}$))) - theta/180) < tol), & + "test_array") end subroutine test_argpi_${k1}$ #:endfor From 7e1e4fa5fdb46f4247dc4318406a62202fbaf018 Mon Sep 17 00:00:00 2001 From: zoziha Date: Sun, 19 Dec 2021 22:04:34 +0800 Subject: [PATCH 6/6] Update CHANGELOG.md --- CHANGELOG.md | 5 +++++ doc/specs/stdlib_math.md | 2 -- 2 files changed, 5 insertions(+), 2 deletions(-) diff --git a/CHANGELOG.md b/CHANGELOG.md index da425d130..491a0e74e 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -13,6 +13,11 @@ 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 procedures `is_close` and `all_close` + [#488](https://github.com/fortran-lang/stdlib/pull/488) + - new procedures `arg`, `argd` and `argpi` + [#498](https://github.com/fortran-lang/stdlib/pull/498) Changes to existing modules diff --git a/doc/specs/stdlib_math.md b/doc/specs/stdlib_math.md index 4d2e6d402..665d6ee5f 100644 --- a/doc/specs/stdlib_math.md +++ b/doc/specs/stdlib_math.md @@ -572,7 +572,6 @@ Returns a `logical` scalar/array. 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 @@ -637,7 +636,6 @@ Returns a `logical` scalar. 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)