Math improvements

This commit is contained in:
Martin Diehl 2023-08-11 08:01:55 +00:00 committed by Franz Roters
parent 651d0d989b
commit a1eb7829ee
3 changed files with 63 additions and 41 deletions

View File

@ -27,8 +27,10 @@ module math
#if __INTEL_COMPILER >= 1900 #if __INTEL_COMPILER >= 1900
! do not make use of associated entities available to other modules ! do not make use of associated entities available to other modules
private :: & private :: &
misc, &
IO, & IO, &
config config, &
parallelization
#endif #endif
real(pREAL), parameter :: & real(pREAL), parameter :: &
@ -38,11 +40,11 @@ module math
INRAD = TAU/360.0_pREAL !< conversion from degree to radian INRAD = TAU/360.0_pREAL !< conversion from degree to radian
real(pREAL), dimension(3,3), parameter :: & real(pREAL), dimension(3,3), parameter :: &
math_I3 = reshape([& math_I3 = real(reshape([&
1.0_pREAL,0.0_pREAL,0.0_pREAL, & 1, 0, 0, &
0.0_pREAL,1.0_pREAL,0.0_pREAL, & 0, 1, 0, &
0.0_pREAL,0.0_pREAL,1.0_pREAL & 0, 0, 1 &
],shape(math_I3)) !< 3x3 Identity ],shape(math_I3)),pREAL) !< 3x3 Identity
real(pREAL), dimension(*), parameter, private :: & real(pREAL), dimension(*), parameter, private :: &
NRMMANDEL = [1.0_pREAL, 1.0_pREAL,1.0_pREAL, sqrt(2.0_pREAL), sqrt(2.0_pREAL), sqrt(2.0_pREAL)] !< forward weighting for Mandel notation NRMMANDEL = [1.0_pREAL, 1.0_pREAL,1.0_pREAL, sqrt(2.0_pREAL), sqrt(2.0_pREAL), sqrt(2.0_pREAL)] !< forward weighting for Mandel notation
@ -83,9 +85,6 @@ module math
3,3 & 3,3 &
],shape(MAPPLAIN)) !< arrangement in Plain notation ],shape(MAPPLAIN)) !< arrangement in Plain notation
!---------------------------------------------------------------------------------------------------
private :: &
selfTest
contains contains
@ -109,20 +108,21 @@ subroutine math_init()
allocate(seed(randSize)) allocate(seed(randSize))
if (num_generic%contains('random_seed')) then if (num_generic%contains('random_seed')) then
seed = num_generic%get_as1dInt('random_seed',requiredSize=randSize) seed = num_generic%get_as1dInt('random_seed',requiredSize=randSize) &
+ worldrank*42_MPI_INTEGER_KIND
else else
call random_seed() call random_seed()
call random_seed(get = seed) call random_seed(get = seed)
end if end if
call random_seed(put = seed + worldrank*42_MPI_INTEGER_KIND) call random_seed(put = seed)
call random_number(randTest) call random_number(randTest)
print'(/,a,i2)', ' size of random seed: ', randSize print'(/,a,i2)', ' size of random seed: ', randSize
print*, 'value of random seed: ', seed print*, 'value of random seed: ', seed
print'( a,4(/,26x,f17.14))', ' start of random sequence: ', randTest print'( a,4(/,26x,f17.14))', ' start of random sequence: ', randTest
call selfTest() call math_selfTest()
end subroutine math_init end subroutine math_init
@ -1275,7 +1275,7 @@ end function math_clip
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
!> @brief Check correctness of some math functions. !> @brief Check correctness of some math functions.
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
subroutine selfTest() subroutine math_selfTest()
integer, dimension(2,4) :: & integer, dimension(2,4) :: &
sort_in_ = reshape([+1,+5, +5,+6, -1,-1, +3,-2],[2,4]) sort_in_ = reshape([+1,+5, +5,+6, -1,-1, +3,-2],[2,4])
@ -1447,6 +1447,6 @@ subroutine selfTest()
error stop 'math_normal(sigma)' error stop 'math_normal(sigma)'
end block normal_distribution end block normal_distribution
end subroutine selfTest end subroutine math_selfTest
end module math end module math

View File

@ -5,51 +5,56 @@ program DAMASK_test
use IO use IO
use test_prec use test_prec
use test_crystal
use test_IO
use test_rotations
use test_misc use test_misc
use test_math
use test_tables use test_tables
use test_crystal
use test_rotations
use test_IO
use test_HDF5_utilities use test_HDF5_utilities
external :: quit external :: quit
character(len=*), parameter :: & character(len=*), parameter :: &
tab = '19', & ok = achar(27)//'[32mok'//achar(27)//'[0m', &
ok = achar(27)//'[32mok'//achar(27)//'[0m' fmt = '(3x,a,T19,a,1x)'
call parallelization_init() call parallelization_init()
call HDF5_utilities_init() call HDF5_utilities_init()
write(IO_STDOUT,fmt='(/,1x,a,/)') achar(27)//'[1m'//'testing'//achar(27)//'[0m' write(IO_STDOUT,fmt='(/,1x,a,/)') achar(27)//'[1m'//'testing'//achar(27)//'[0m'
write(IO_STDOUT,fmt='(3x,a,T'//tab//',a)', advance='no') 'prec','...' write(IO_STDOUT,fmt=fmt, advance='no') 'prec','...'
call test_prec_run() call test_prec_run()
write(IO_STDOUT,fmt='(1x,a)') ok write(IO_STDOUT,fmt='(a)') ok
write(IO_STDOUT,fmt='(3x,a,T'//tab//',a)', advance='no') 'tables','...' write(IO_STDOUT,fmt=fmt, advance='no') 'misc','...'
call test_tables_run()
write(IO_STDOUT,fmt='(1x,a)') ok
write(IO_STDOUT,fmt='(3x,a,T'//tab//',a)', advance='no') 'crystal','...'
call test_crystal_run()
write(IO_STDOUT,fmt='(1x,a)') ok
write(IO_STDOUT,fmt='(3x,a,T'//tab//',a)', advance='no') 'rotations','...'
call test_rotations_run()
write(IO_STDOUT,fmt='(1x,a)') ok
write(IO_STDOUT,fmt='(3x,a,T'//tab//',a)', advance='no') 'IO','...'
call test_IO_run()
write(IO_STDOUT,fmt='(1x,a)') ok
write(IO_STDOUT,fmt='(3x,a,T'//tab//',a)', advance='no') 'misc','...'
call test_misc_run() call test_misc_run()
write(IO_STDOUT,fmt='(1x,a)') ok write(IO_STDOUT,fmt='(a)') ok
write(IO_STDOUT,fmt='(3x,a,T'//tab//',a)', advance='no') 'HDF5_utilities','...' write(IO_STDOUT,fmt=fmt, advance='no') 'math','...'
call test_math_run()
write(IO_STDOUT,fmt='(a)') ok
write(IO_STDOUT,fmt=fmt, advance='no') 'tables','...'
call test_tables_run()
write(IO_STDOUT,fmt='(a)') ok
write(IO_STDOUT,fmt=fmt, advance='no') 'crystal','...'
call test_crystal_run()
write(IO_STDOUT,fmt='(a)') ok
write(IO_STDOUT,fmt=fmt, advance='no') 'rotations','...'
call test_rotations_run()
write(IO_STDOUT,fmt='(a)') ok
write(IO_STDOUT,fmt=fmt, advance='no') 'IO','...'
call test_IO_run()
write(IO_STDOUT,fmt='(a)') ok
write(IO_STDOUT,fmt=fmt, advance='no') 'HDF5_utilities','...'
call test_HDF5_utilities_run() call test_HDF5_utilities_run()
write(IO_STDOUT,fmt='(1x,a)') ok write(IO_STDOUT,fmt='(a)') ok
call quit(0) call quit(0)

17
src/test/test_math.f90 Normal file
View File

@ -0,0 +1,17 @@
module test_math
use math
implicit none(type,external)
private
public :: test_math_run
contains
subroutine test_math_run()
call math_selfTest()
end subroutine test_math_run
end module test_math