The ParaMonte Documentation Website
Current view: top level - kernel/tests - Test_Cosmology_mod.f90 (source / functions) Hit Total Coverage
Test: ParaMonte 1.5.1 :: Serial Kernel - Code Coverage Report Lines: 58 58 100.0 %
Date: 2021-01-08 13:03:42 Functions: 7 7 100.0 %
Legend: Lines: hit not hit

          Line data    Source code
       1             : !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
       2             : !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
       3             : !!!!
       4             : !!!!   MIT License
       5             : !!!!
       6             : !!!!   ParaMonte: plain powerful parallel Monte Carlo library.
       7             : !!!!
       8             : !!!!   Copyright (C) 2012-present, The Computational Data Science Lab
       9             : !!!!
      10             : !!!!   This file is part of the ParaMonte library.
      11             : !!!!
      12             : !!!!   Permission is hereby granted, free of charge, to any person obtaining a
      13             : !!!!   copy of this software and associated documentation files (the "Software"),
      14             : !!!!   to deal in the Software without restriction, including without limitation
      15             : !!!!   the rights to use, copy, modify, merge, publish, distribute, sublicense,
      16             : !!!!   and/or sell copies of the Software, and to permit persons to whom the
      17             : !!!!   Software is furnished to do so, subject to the following conditions:
      18             : !!!!
      19             : !!!!   The above copyright notice and this permission notice shall be
      20             : !!!!   included in all copies or substantial portions of the Software.
      21             : !!!!
      22             : !!!!   THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
      23             : !!!!   EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF
      24             : !!!!   MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT.
      25             : !!!!   IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM,
      26             : !!!!   DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR
      27             : !!!!   OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE
      28             : !!!!   OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.
      29             : !!!!
      30             : !!!!   ACKNOWLEDGMENT
      31             : !!!!
      32             : !!!!   ParaMonte is an honor-ware and its currency is acknowledgment and citations.
      33             : !!!!   As per the ParaMonte library license agreement terms, if you use any parts of
      34             : !!!!   this library for any purposes, kindly acknowledge the use of ParaMonte in your
      35             : !!!!   work (education/research/industry/development/...) by citing the ParaMonte
      36             : !!!!   library as described on this page:
      37             : !!!!
      38             : !!!!       https://github.com/cdslaborg/paramonte/blob/main/ACKNOWLEDGMENT.md
      39             : !!!!
      40             : !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      41             : !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      42             : 
      43             : !>  \brief This module contains tests of the module [Cosmology_mod](@ref cosmology_mod).
      44             : !>  \author Amir Shahmoradi
      45             : 
      46             : module Test_Cosmology_mod
      47             : 
      48             :     use Cosmology_mod
      49             :     use Err_mod, only: Err_type
      50             :     use Test_mod, only: Test_type
      51             :     implicit none
      52             : 
      53             :     private
      54             :     public :: test_Cosmology
      55             : 
      56             :     type(Test_type) :: Test
      57             : 
      58             : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
      59             : 
      60             : contains
      61             : 
      62             : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
      63             : 
      64           1 :     subroutine test_Cosmology()
      65             :         implicit none
      66           1 :         Test = Test_type(moduleName=MODULE_NAME)
      67           1 :         call Test%run(test_getlogdvdz, "test_getlogdvdz")
      68           1 :         call Test%run(test_ldiswickram, "test_ldiswickram")
      69           1 :         call Test%run(test_getLogLumDisWicMpc, "test_getLogLumDisWicMpc")
      70           1 :         call Test%run(test_getLookBackTime_1, "test_getLookBackTime_1") ! The internal function passing as actual argument causes segfault with Gfortran (any version) on Windows subsystem for Linux.
      71           1 :         call Test%run(test_getLookBackTime_2, "test_getLookBackTime_2") ! The internal function passing as actual argument causes segfault with Gfortran (any version) on Windows subsystem for Linux.
      72           1 :         call Test%run(test_getUniverseAgeDerivative, "test_getUniverseAgeDerivative")
      73           1 :         call Test%finalize()
      74           1 :     end subroutine test_Cosmology
      75             : 
      76             : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
      77             : 
      78           1 :     function test_getlogdvdz() result(assertion)
      79             : 
      80           1 :         use Constants_mod, only: RK, IK
      81             : 
      82             :         implicit none
      83             : 
      84             :         logical :: assertion
      85             :         real(RK), parameter :: zplus1 = 1.e1_RK
      86             :         real(RK), parameter :: LOGzplus1 = log(zplus1)
      87             :         real(RK), parameter :: twiceLogLumDisMpc = 0._RK
      88             :         real(RK), parameter :: logdvdz_ref = 3.421655392580980_RK
      89             :         real(RK), parameter :: tolerance = 1.e-10_RK
      90           1 :         real(RK)            :: difference
      91           1 :         real(RK)            :: logdvdz
      92             : 
      93           1 :         logdvdz = getlogdvdz(zplus1,LOGzplus1,twiceLogLumDisMpc)
      94           1 :         difference = abs( (logdvdz - logdvdz_ref) / logdvdz_ref )
      95           1 :         assertion = difference < tolerance
      96             : 
      97             :         ! LCOV_EXCL_START
      98             :         if (Test%isDebugMode .and. .not. assertion) then
      99             :             write(Test%outputUnit,"(*(g0.15))")
     100             :             write(Test%outputUnit,"(*(g0.15))") "logdvdz_ref   = ", logdvdz_ref
     101             :             write(Test%outputUnit,"(*(g0.15))") "logdvdz       = ", logdvdz
     102             :             write(Test%outputUnit,"(*(g0.15))") "difference    = ", difference
     103             :             write(Test%outputUnit,"(*(g0.15))")
     104             :         end if
     105             :         ! LCOV_EXCL_STOP
     106             : 
     107           1 :     end function test_getlogdvdz
     108             : 
     109             : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
     110             : 
     111           1 :     function test_ldiswickram() result(assertion)
     112             : 
     113           1 :         use Constants_mod, only: RK, IK
     114             : 
     115             :         implicit none
     116             : 
     117             :         logical :: assertion
     118             :         real(RK), parameter :: zplus1 = 1.e1_RK
     119             :         real(RK), parameter :: tolerance = 1.e-8_RK
     120             :         real(RK), parameter :: lumDisWicMpc_ref = 90887.79805125466_RK
     121           1 :         real(RK)            :: lumDisWicMpc
     122           1 :         real(RK)            :: difference
     123             : 
     124           1 :         lumDisWicMpc = ldiswickram(zplus1)
     125           1 :         difference = abs( (lumDisWicMpc - lumDisWicMpc_ref) / lumDisWicMpc_ref )
     126           1 :         assertion = difference < tolerance
     127             : 
     128             :         ! LCOV_EXCL_START
     129             :         if (Test%isDebugMode .and. .not. assertion) then
     130             :             write(Test%outputUnit,"(*(g0.15))")
     131             :             write(Test%outputUnit,"(*(g0.15))") "lumDisWicMpc_ref  = ", lumDisWicMpc_ref
     132             :             write(Test%outputUnit,"(*(g0.15))") "lumDisWicMpc      = ", lumDisWicMpc
     133             :             write(Test%outputUnit,"(*(g0.15))") "difference        = ", difference
     134             :             write(Test%outputUnit,"(*(g0.15))")
     135             :         end if
     136             :         ! LCOV_EXCL_STOP
     137             : 
     138           1 :     end function test_ldiswickram
     139             : 
     140             : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
     141             : 
     142           1 :     function test_getLogLumDisWicMpc() result(assertion)
     143             : 
     144           1 :         use Constants_mod, only: RK, IK
     145             : 
     146             :         implicit none
     147             : 
     148             :         logical :: assertion
     149             :         real(RK), parameter :: zplus1 = 1.e1_RK
     150             :         real(RK), parameter :: tolerance = 1.e-8_RK
     151             :         real(RK), parameter :: logLumDisWicMpc_ref = 11.417381027586867_RK
     152           1 :         real(RK)            :: logLumDisWicMpc
     153           1 :         real(RK)            :: difference
     154             : 
     155           1 :         logLumDisWicMpc = getLogLumDisWicMpc(zplus1)
     156           1 :         difference = abs( (logLumDisWicMpc - logLumDisWicMpc_ref) / logLumDisWicMpc_ref )
     157           1 :         assertion = difference < tolerance
     158             : 
     159             :         ! LCOV_EXCL_START
     160             :         if (Test%isDebugMode .and. .not. assertion) then
     161             :             write(Test%outputUnit,"(*(g0.15))")
     162             :             write(Test%outputUnit,"(*(g0.15))") "logLumDisWicMpc_ref   = ", logLumDisWicMpc_ref
     163             :             write(Test%outputUnit,"(*(g0.15))") "logLumDisWicMpc       = ", logLumDisWicMpc
     164             :             write(Test%outputUnit,"(*(g0.15))") "difference            = ", difference
     165             :             write(Test%outputUnit,"(*(g0.15))")
     166             :         end if
     167             :         ! LCOV_EXCL_STOP
     168             : 
     169           1 :     end function test_getLogLumDisWicMpc
     170             : 
     171             : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
     172             : 
     173           1 :     function test_getLookBackTime_1() result(assertion)
     174             : 
     175           1 :         use Constants_mod, only: RK, IK
     176             : 
     177             :         implicit none
     178             : 
     179             :         logical                 :: assertion
     180             :         real(RK)    , parameter :: zplus1 = 1.e1_RK
     181             :         real(RK)    , parameter :: tolerance = 1.e-12_RK
     182             :         real(RK)    , parameter :: lookBackTime_ref = 12.7740854156781_RK
     183             :         real(RK)    , parameter :: maxRelativeError = 1.e-5_RK
     184             :         integer(IK) , parameter :: nRefinement = 6
     185           1 :         real(RK)                :: lookBackTime
     186           1 :         real(RK)                :: difference
     187             : 
     188           2 :         lookBackTime = getLookBackTime(zplus1 = zplus1, maxRelativeError = maxRelativeError, nRefinement = nRefinement)
     189           1 :         difference = abs( (lookBackTime - lookBackTime_ref) / lookBackTime_ref )
     190           1 :         assertion = difference < tolerance
     191             : 
     192             :         ! LCOV_EXCL_START
     193             :         if (Test%isDebugMode .and. .not. assertion) then
     194             :             write(Test%outputUnit,"(*(g0.15))")
     195             :             write(Test%outputUnit,"(*(g0.15))") "lookBackTime_ref  = ", lookBackTime_ref
     196             :             write(Test%outputUnit,"(*(g0.15))") "lookBackTime      = ", lookBackTime
     197             :             write(Test%outputUnit,"(*(g0.15))") "difference        = ", difference
     198             :             write(Test%outputUnit,"(*(g0.15))")
     199             :         end if
     200             :         ! LCOV_EXCL_STOP
     201             : 
     202           1 :     end function test_getLookBackTime_1
     203             : 
     204             : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
     205             : 
     206           1 :     function test_getLookBackTime_2() result(assertion)
     207             : 
     208           1 :         use Constants_mod, only: RK, IK
     209             : 
     210             :         implicit none
     211             : 
     212             :         logical                 :: assertion
     213             :         real(RK)    , parameter :: zplus1 = 1.e1_RK
     214             :         real(RK)    , parameter :: tolerance = 1.e-12_RK
     215             :         real(RK)    , parameter :: lookBackTime_ref = 12.7736258974511_RK
     216           1 :         real(RK)                :: lookBackTime
     217           1 :         real(RK)                :: difference
     218             : 
     219           2 :         lookBackTime = getLookBackTime(zplus1 = zplus1)
     220           1 :         difference = abs( (lookBackTime - lookBackTime_ref) / lookBackTime_ref )
     221           1 :         assertion = difference < tolerance
     222             : 
     223             :         ! LCOV_EXCL_START
     224             :         if (Test%isDebugMode .and. .not. assertion) then
     225             :             write(Test%outputUnit,"(*(g0.15))")
     226             :             write(Test%outputUnit,"(*(g0.15))") "lookBackTime_ref  = ", lookBackTime_ref
     227             :             write(Test%outputUnit,"(*(g0.15))") "lookBackTime      = ", lookBackTime
     228             :             write(Test%outputUnit,"(*(g0.15))") "difference        = ", difference
     229             :             write(Test%outputUnit,"(*(g0.15))")
     230             :         end if
     231             :         ! LCOV_EXCL_STOP
     232             : 
     233           1 :     end function test_getLookBackTime_2
     234             : 
     235             : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
     236             : 
     237           1 :     function test_getUniverseAgeDerivative() result(assertion)
     238             : 
     239           1 :         use Constants_mod, only: RK, IK
     240             : 
     241             :         implicit none
     242             : 
     243             :         logical                 :: assertion
     244             :         real(RK)    , parameter :: zplus1 = 1.e1_RK
     245             :         real(RK)    , parameter :: tolerance = 1.e-12_RK
     246             :         real(RK)    , parameter :: universeAgeDerivative_ref = 8.122223525986105e-05_RK
     247           1 :         real(RK)                :: universeAgeDerivative
     248           1 :         real(RK)                :: difference
     249             : 
     250           1 :         universeAgeDerivative = getUniverseAgeDerivative(zplus1 = zplus1)
     251           1 :         difference = abs( (universeAgeDerivative - universeAgeDerivative_ref) / universeAgeDerivative_ref )
     252           1 :         assertion = difference < tolerance
     253             : 
     254             :         ! LCOV_EXCL_START
     255             :         if (Test%isDebugMode .and. .not. assertion) then
     256             :             write(Test%outputUnit,"(*(g0.15))")
     257             :             write(Test%outputUnit,"(*(g0.15))") "universeAgeDerivative_ref = ", universeAgeDerivative_ref
     258             :             write(Test%outputUnit,"(*(g0.15))") "universeAgeDerivative     = ", universeAgeDerivative
     259             :             write(Test%outputUnit,"(*(g0.15))") "difference                = ", difference
     260             :             write(Test%outputUnit,"(*(g0.15))")
     261             :         end if
     262             :         ! LCOV_EXCL_STOP
     263             : 
     264           1 :     end function test_getUniverseAgeDerivative
     265             : 
     266             : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
     267             : 
     268             : end module Test_Cosmology_mod ! LCOV_EXCL_LINE

ParaMonte: Plain Powerful Parallel Monte Carlo Library 
The Computational Data Science Lab
© Copyright 2012 - 2021