https://www.cdslab.org/paramonte/fortran/2
Current view: top level - test - test_pm_bench.F90 (source / functions) Hit Total Coverage
Test: ParaMonte 2.0.0 :: Serial Fortran - Code Coverage Report Lines: 138 138 100.0 %
Date: 2024-04-08 03:18:57 Functions: 12 12 100.0 %
Legend: Lines: hit not hit

          Line data    Source code
       1             : !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
       2             : !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
       3             : !!!!                                                                                                                            !!!!
       4             : !!!!    ParaMonte: Parallel Monte Carlo and Machine Learning Library.                                                           !!!!
       5             : !!!!                                                                                                                            !!!!
       6             : !!!!    Copyright (C) 2012-present, The Computational Data Science Lab                                                          !!!!
       7             : !!!!                                                                                                                            !!!!
       8             : !!!!    This file is part of the ParaMonte library.                                                                             !!!!
       9             : !!!!                                                                                                                            !!!!
      10             : !!!!    LICENSE                                                                                                                 !!!!
      11             : !!!!                                                                                                                            !!!!
      12             : !!!!       https://github.com/cdslaborg/paramonte/blob/main/LICENSE.md                                                          !!!!
      13             : !!!!                                                                                                                            !!!!
      14             : !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      15             : !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      16             : 
      17             : !>  \brief
      18             : !>  This module contains tests of the module [pm_bench](@ref pm_bench).
      19             : !>
      20             : !>  \todo
      21             : !>  \phigh 
      22             : !>  The following tests should be extended to include allocatable arrays of arbitrary bounds.
      23             : !>
      24             : !>  \fintest
      25             : !>
      26             : !>  \author
      27             : !>  \FatemehBagheri, Wednesday 12:20 AM, October 13, 2021, Dallas, TX
      28             : 
      29             : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
      30             : 
      31             : module test_pm_bench
      32             : 
      33             :     use pm_bench
      34             :     use pm_test, only: test_type, LK
      35             :     use pm_timer, only: timerCPU_type
      36             :     use pm_timer, only: timerDAT_type
      37             :     use pm_timer, only: timerSYS_type
      38             :     use pm_kind, only: IK, RK, RKS, RKD, RKQ
      39             : 
      40             :     implicit none
      41             : 
      42             :     private
      43             :     public :: setTest
      44             :     type(test_type) :: test
      45             : 
      46             :     real(RKS)       :: unifrnd_RKS
      47             :     real(RKS)       :: unifsum_RKS = 0._RKS ! dummy calculation to prevent aggressive compiler optimizations.
      48             :     real(RKD)       :: unifrnd_RKD
      49             :     real(RKD)       :: unifsum_RKD = 0._RKD ! dummy calculation to prevent aggressive compiler optimizations.
      50             :     real(RKQ)       :: unifrnd_RKQ
      51             :     real(RKQ)       :: unifsum_RKQ = 0._RKQ ! dummy calculation to prevent aggressive compiler optimizations.
      52             : 
      53             : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
      54             : 
      55             : contains
      56             : 
      57             : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
      58             : 
      59           1 :     subroutine setTest()
      60             : 
      61           1 :         test = test_type(MODULE_NAME)
      62           1 :         call test%run(test_constructBench, SK_"test_constructBench")
      63           1 :         call test%run(test_constructBenchMulti, SK_"test_constructBenchMulti")
      64           1 :         call test%summarize()
      65             : 
      66           1 :     end subroutine setTest
      67             : 
      68             : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
      69             : 
      70           1 :     function test_constructBench() result(assertion)
      71             :         use pm_timer, only: timerCPU_type
      72             :         use pm_timer, only: timerDAT_type
      73             :         use pm_timer, only: timerSYS_type
      74             :         logical(LK) :: assertion
      75             :         type(timerCPU_type) :: TimerCPU
      76             :         type(timerDAT_type) :: TimerDAT
      77             :         type(timerSYS_type) :: TimerSYS
      78           1 :         type(bench_type) :: bench
      79             :         integer(IK) :: miniter_def
      80             :         integer(IK) :: miniter
      81             :         real(RKD) :: minsec_def
      82             :         real(RKD) :: minsec
      83             :         
      84           1 :         miniter = 100_IK
      85           1 :         minsec = 0.02_RKD
      86             :         miniter_def = 1
      87             :         minsec_def = 0.05_RKD
      88           1 :         assertion = .true._LK
      89             : 
      90           1 :         TimerCPU = timerCPU_type()
      91           1 :         TimerDAT = timerDAT_type()
      92           1 :         TimerSYS = timerSYS_type()
      93             : 
      94           1 :         bench = bench_type("testing", exec)
      95           1 :         assertion = assertion .and. logical(bench%name == SK_"testing", LK)
      96           1 :         call test%assert(assertion, line = int(__LINE__, IK))
      97             : 
      98           1 :         TimerSYS%start = TimerSYS%time()
      99           1 :         call bench%setTiming()
     100           1 :         TimerSYS%delta = TimerSYS%time(since = TimerSYS%start)
     101             :         !assertion = assertion .and. logical(TimerSYS%delta >= minsec_def, LK)
     102           1 :         call test%assert(assertion, line = int(__LINE__, IK))
     103           1 :         assertion = assertion .and. logical(size(bench%timing%values, kind = IK) >= miniter_def, LK)
     104           1 :         call test%assert(assertion, line = int(__LINE__, IK))
     105             : 
     106           1 :         bench = bench_type("testing", exec, overhead)
     107           1 :         assertion = assertion .and. logical(bench%name == SK_"testing", LK)
     108           1 :         call test%assert(assertion, line = int(__LINE__, IK))
     109             : 
     110           1 :         call bench%setTiming()
     111           1 :         assertion = assertion .and. logical(size(bench%timing%values, kind = IK) >= miniter_def, LK)
     112           1 :         call test%assert(assertion, line = int(__LINE__, IK))
     113             : 
     114           1 :         bench = bench_type("testing", exec, overhead, minsec = minsec)
     115           1 :         assertion = assertion .and. logical(bench%name == SK_"testing", LK)
     116           1 :         call test%assert(assertion, line = int(__LINE__, IK))
     117           1 :         assertion = assertion .and. logical(bench%minsec == minsec, LK)
     118           1 :         call test%assert(assertion, line = int(__LINE__, IK))
     119             : 
     120           1 :         call bench%setTiming()
     121           1 :         assertion = assertion .and. logical(size(bench%timing%values, kind = IK) >= miniter_def, LK)
     122           1 :         call test%assert(assertion, line = int(__LINE__, IK))
     123             : 
     124           1 :         bench = bench_type("testing", exec, overhead, minsec = minsec, miniter = miniter)
     125           1 :         assertion = assertion .and. logical(bench%name == SK_"testing", LK)
     126           1 :         call test%assert(assertion, line = int(__LINE__, IK))
     127           1 :         assertion = assertion .and. logical(bench%minsec == minsec, LK)
     128           1 :         call test%assert(assertion, line = int(__LINE__, IK))
     129           1 :         assertion = assertion .and. logical(bench%miniter == miniter, LK)
     130           1 :         call test%assert(assertion, line = int(__LINE__, IK))
     131             : 
     132           1 :         bench = bench_type("testing", exec, overhead, minsec = minsec, miniter = miniter, timer = timerCPU_type())
     133           1 :         assertion = assertion .and. logical(bench%name == SK_"testing", LK)
     134           1 :         call test%assert(assertion, line = int(__LINE__, IK))
     135           1 :         assertion = assertion .and. logical(bench%minsec == minsec, LK)
     136           1 :         call test%assert(assertion, line = int(__LINE__, IK))
     137           1 :         assertion = assertion .and. logical(bench%miniter == miniter, LK)
     138           1 :         call test%assert(assertion, line = int(__LINE__, IK))
     139             : 
     140           1 :         TimerCPU%start = TimerCPU%time()
     141           1 :         call bench%setTiming()
     142           1 :         TimerCPU%delta = TimerCPU%time(since = TimerCPU%start)
     143             :         !assertion = assertion .and. logical(TimerCPU%delta >= minsec, LK)
     144           1 :         call test%assert(assertion, line = int(__LINE__, IK))
     145             : 
     146           1 :         bench = bench_type("testing", exec, overhead, minsec = minsec, miniter = miniter, timer = timerDAT_type())
     147           1 :         assertion = assertion .and. logical(bench%name == SK_"testing", LK)
     148           1 :         call test%assert(assertion, line = int(__LINE__, IK))
     149           1 :         assertion = assertion .and. logical(bench%minsec == minsec, LK)
     150           1 :         call test%assert(assertion, line = int(__LINE__, IK))
     151           1 :         assertion = assertion .and. logical(bench%miniter == miniter, LK)
     152           1 :         call test%assert(assertion, line = int(__LINE__, IK))
     153             : 
     154           1 :         TimerDAT%start = TimerDAT%time()
     155           1 :         call bench%setTiming(miniter = miniter)
     156           1 :         TimerDAT%delta = TimerDAT%time(since = TimerDAT%start)
     157             :         !assertion = assertion .and. logical(TimerDAT%delta >= minsec, LK)
     158             :         !print *, TimerDAT%delta, minsec
     159             :         !print *, TimerDAT%resol, bench%timer%resol
     160           1 :         call test%assert(assertion, line = int(__LINE__, IK))
     161           1 :         assertion = assertion .and. logical(size(bench%timing%values, kind = IK) >= miniter, LK)
     162           1 :         call test%assert(assertion, line = int(__LINE__, IK))
     163             : 
     164           1 :         bench = bench_type("testing", exec, overhead, minsec = minsec, miniter = miniter, timer = timerSYS_type())
     165           1 :         assertion = assertion .and. logical(bench%name == SK_"testing", LK)
     166           1 :         call test%assert(assertion, line = int(__LINE__, IK))
     167           1 :         assertion = assertion .and. logical(bench%minsec == minsec, LK)
     168           1 :         call test%assert(assertion, line = int(__LINE__, IK))
     169           1 :         assertion = assertion .and. logical(bench%miniter == miniter, LK)
     170           1 :         call test%assert(assertion, line = int(__LINE__, IK))
     171             : 
     172           1 :         TimerSYS%start = TimerSYS%time()
     173           1 :         call bench%setTiming(miniter = miniter)
     174           1 :         TimerSYS%delta = TimerSYS%time(since = TimerSYS%start)
     175             :         !assertion = assertion .and. logical(TimerSYS%delta >= minsec_def, LK)
     176           1 :         call test%assert(assertion, line = int(__LINE__, IK))
     177           1 :         assertion = assertion .and. logical(size(bench%timing%values, kind = IK) >= miniter, LK)
     178           1 :         call test%assert(assertion, line = int(__LINE__, IK))
     179             : 
     180           1 :         bench%timing = bench%getTiming(miniter = miniter)
     181           1 :         assertion = assertion .and. logical(size(bench%timing%values, kind = IK) >= miniter, LK)
     182           1 :         call test%assert(assertion, line = int(__LINE__, IK))
     183             : 
     184           1 :         bench%timing = bench%getTiming(miniter = miniter)
     185           1 :         assertion = assertion .and. logical(size(bench%timing%values, kind = IK) >= miniter, LK)
     186           1 :         call test%assert(assertion, line = int(__LINE__, IK))
     187             : 
     188           1 :         bench%timing = bench%getTiming(minsec = minsec, miniter = miniter)
     189           1 :         assertion = assertion .and. logical(size(bench%timing%values, kind = IK) >= miniter, LK)
     190           1 :         call test%assert(assertion, line = int(__LINE__, IK))
     191             : 
     192           1 :         bench%timing = bench%getTiming(miniter = miniter)
     193           1 :         assertion = assertion .and. logical(size(bench%timing%values, kind = IK) >= miniter, LK)
     194           1 :         call test%assert(assertion, line = int(__LINE__, IK))
     195             : 
     196           1 :         call bench%setTiming(miniter = miniter)
     197           1 :         assertion = assertion .and. logical(size(bench%timing%values, kind = IK) >= miniter, LK)
     198           1 :         call test%assert(assertion, line = int(__LINE__, IK))
     199             : 
     200           1 :         call bench%setTiming(minsec = minsec)
     201           1 :         call test%assert(assertion, line = int(__LINE__, IK))
     202             : 
     203           1 :         call bench%setTiming(minsec = minsec, miniter = miniter)
     204           1 :         assertion = assertion .and. logical(size(bench%timing%values, kind = IK) >= miniter, LK)
     205           1 :         call test%assert(assertion, line = int(__LINE__, IK))
     206             : 
     207           1 :         miniter = 10001_IK
     208           1 :         call bench%setTiming(miniter = miniter)
     209           1 :         assertion = assertion .and. logical(size(bench%timing%values, kind = IK) >= miniter, LK)
     210           1 :         call test%assert(assertion, line = int(__LINE__, IK))
     211             : 
     212           1 :     end function
     213             : 
     214             : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
     215             : 
     216           1 :     function test_constructBenchMulti() result(assertion)
     217             :         use pm_io, only: display_type
     218             :         logical(LK)             :: assertion
     219           1 :         type(display_type)      :: disp
     220           1 :         type(benchMulti_type)   :: benchMulti
     221             :      
     222           1 :         assertion = .true._LK
     223             : 
     224           1 :         disp = display_type(file = "test_constructBenchMulti.tmp")
     225             : 
     226             :         benchMulti= benchMulti_type([ bench_type(SK_'rng_single', procwrap_RKS, overhead_RKS) & ! LCOV_EXCL_LINE
     227             :                                     , bench_type(SK_'rng_double', procwrap_RKD, overhead_RKD) & ! LCOV_EXCL_LINE
     228             :                                     , bench_type(SK_'rng_quadro', procwrap_RKQ, overhead_RKQ) & ! LCOV_EXCL_LINE
     229           7 :                                     ])
     230           1 :         call test%assert(assertion, line = int(__LINE__, IK))
     231           1 :         call benchMulti%showsum(unit = disp%unit) ! Display the results in a nice tabular format.
     232           1 :         call test%assert(assertion, line = int(__LINE__, IK))
     233           1 :         call benchMulti%showsum(unit = disp%unit, tabular = .true._LK) ! Display the results in a nice tabular format.
     234           1 :         call test%assert(assertion, line = int(__LINE__, IK))
     235           1 :         call benchMulti%showsum(unit = disp%unit, tabular = .false._LK) ! Display the results in a nice tabular format.
     236           1 :         call test%assert(assertion, line = int(__LINE__, IK))
     237             :         benchMulti= benchMulti_type([ bench_type(SK_'rng_single', procwrap_RKS, overhead_RKS) & ! LCOV_EXCL_LINE
     238             :                                     , bench_type(SK_'rng_double', procwrap_RKD, overhead_RKD) & ! LCOV_EXCL_LINE
     239             :                                     , bench_type(SK_'rng_quadro', procwrap_RKQ, overhead_RKQ) & ! LCOV_EXCL_LINE
     240             :                                     ] & ! LCOV_EXCL_LINE
     241             :                                     , sorted = .false. & ! LCOV_EXCL_LINE
     242           7 :                                     )
     243             :         benchMulti= benchMulti_type([ bench_type(SK_'rng_single', procwrap_RKS, overhead_RKS) & ! LCOV_EXCL_LINE
     244             :                                     , bench_type(SK_'rng_double', procwrap_RKD, overhead_RKD) & ! LCOV_EXCL_LINE
     245             :                                     , bench_type(SK_'rng_quadro', procwrap_RKQ, overhead_RKQ) & ! LCOV_EXCL_LINE
     246             :                                     ] & ! LCOV_EXCL_LINE
     247             :                                     , sorted = .true. & ! LCOV_EXCL_LINE
     248           7 :                                     )
     249             :         benchMulti= benchMulti_type([ bench_type(SK_'rng_single', procwrap_RKS, overhead_RKS) & ! LCOV_EXCL_LINE
     250             :                                     , bench_type(SK_'rng_double', procwrap_RKD, overhead_RKD) & ! LCOV_EXCL_LINE
     251             :                                     , bench_type(SK_'rng_quadro', procwrap_RKQ, overhead_RKQ) & ! LCOV_EXCL_LINE
     252             :                                     ] & ! LCOV_EXCL_LINE
     253             :                                     , repeat = 1_IK & ! LCOV_EXCL_LINE
     254           7 :                                     )
     255             : 
     256           1 :     end function
     257             : 
     258             : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
     259             : 
     260      590162 :     subroutine overhead_RKS(); call random_number(unifrnd_RKS); unifsum_RKS = unifsum_RKS + unifrnd_RKS; end
     261      577455 :     subroutine overhead_RKD(); call random_number(unifrnd_RKD); unifsum_RKD = unifsum_RKD + unifrnd_RKD; end
     262      180633 :     subroutine overhead_RKQ(); call random_number(unifrnd_RKQ); unifsum_RKQ = unifsum_RKQ + unifrnd_RKQ; end
     263      590169 :     subroutine procwrap_RKS(); call random_number(unifrnd_RKS); unifsum_RKS = unifsum_RKS + 2._RKS**unifrnd_RKS; end
     264      577462 :     subroutine procwrap_RKD(); call random_number(unifrnd_RKD); unifsum_RKD = unifsum_RKD + 2._RKD**unifrnd_RKD; end
     265      180640 :     subroutine procwrap_RKQ(); call random_number(unifrnd_RKQ); unifsum_RKQ = unifsum_RKQ + 2._RKQ**unifrnd_RKQ; end
     266             : 
     267             : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
     268             : 
     269       15892 :     subroutine exec()
     270             :         real :: summ
     271       15892 :         summ = wasteTime()
     272       15892 :     end subroutine
     273             : 
     274       14915 :     subroutine overhead()
     275       14915 :     end subroutine
     276             : 
     277             : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
     278             : 
     279       15892 :     function wasteTime() result(summ)
     280             :         real :: summ, Matrix(100,100)
     281       15892 :         call random_number(Matrix)
     282   160525092 :         summ = sum(Matrix)
     283       15892 :     end function
     284             : 
     285             : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
     286             : 
     287             : end module test_pm_bench

ParaMonte: Parallel Monte Carlo and Machine Learning Library 
The Computational Data Science Lab
© Copyright 2012 - 2024