The ParaMonte Documentation Website
Current view: top level - kernel/tests - Test_Matrix_mod.f90 (source / functions) Hit Total Coverage
Test: ParaMonte 1.5.1 :: Serial Kernel - Code Coverage Report Lines: 244 244 100.0 %
Date: 2021-01-08 13:03:42 Functions: 25 25 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 [Matrix_mod](@ref matrix_mod).
      44             : !>  \author Amir Shahmoradi
      45             : 
      46             : module Test_Matrix_mod
      47             : 
      48             :     use Matrix_mod ! LCOV_EXCL_LINE
      49             :     use Err_mod, only: Err_type
      50             :     use Test_mod, only: Test_type
      51             :     implicit none
      52             : 
      53             :     private
      54             :     public :: test_Matrix
      55             : 
      56             :     type(Test_type) :: Test
      57             : 
      58             : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
      59             : 
      60             : contains
      61             : 
      62             : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
      63             : 
      64           1 :     subroutine test_Matrix()
      65             : 
      66             :         implicit none
      67             : 
      68           1 :         Test = Test_type(moduleName=MODULE_NAME)
      69           1 :         call Test%run(test_eye_1, "test_eye_1")
      70           1 :         call Test%run(test_eye_2, "test_eye_2")
      71           1 :         call Test%run(test_isPosDef_1, "test_isPosDef_1")
      72           1 :         call Test%run(test_isPosDef_2, "test_isPosDef_2")
      73           1 :         call Test%run(test_getInvMat_1, "test_getInvMat_1")
      74           1 :         call Test%run(test_getOuterProd_1, "test_getOuterProd_1")
      75           1 :         call Test%run(test_getInvMatDet_1, "test_getInvMatDet_1")
      76           1 :         call Test%run(test_sortPosDefMat_1, "test_sortPosDefMat_1")
      77           1 :         call Test%run(test_getRegresCoef_1, "test_getRegresCoef_1")
      78           1 :         call Test%run(test_multiplyMatrix_1, "test_multiplyMatrix_1")
      79           1 :         call Test%run(test_getDeterminant_1, "test_getDeterminant_1")
      80           1 :         call Test%run(test_getInvPosDefMat_1, "test_getInvPosDefMat_1")
      81           1 :         call Test%run(test_getInvPosDefMat_2, "test_getInvPosDefMat_2")
      82           1 :         call Test%run(test_getCholeskyFactor_1, "test_getCholeskyFactor_1")
      83           1 :         call Test%run(test_getCholeskyFactor_2, "test_getCholeskyFactor_2")
      84           1 :         call Test%run(test_getInvMatFromCholFac, "test_getInvMatFromCholFac")
      85           1 :         call Test%run(test_getSqrtDetPosDefMat_1, "test_getSqrtDetPosDefMat_1")
      86           1 :         call Test%run(test_getSqrtDetPosDefMat_2, "test_getSqrtDetPosDefMat_2")
      87           1 :         call Test%run(test_getInvPosDefMatSqrtDet_1, "test_getInvPosDefMatSqrtDet_1")
      88           1 :         call Test%run(test_getInvPosDefMatSqrtDet_2, "test_getInvPosDefMatSqrtDet_2")
      89           1 :         call Test%run(test_getInvPosDefMatSqrtDet_3, "test_getInvPosDefMatSqrtDet_3")
      90           1 :         call Test%run(test_getLogSqrtDetPosDefMat_1, "test_getLogSqrtDetPosDefMat_1")
      91           1 :         call Test%run(test_getLogSqrtDetPosDefMat_2, "test_getLogSqrtDetPosDefMat_2")
      92           1 :         call Test%run(test_symmetrizeUpperSquareMatrix_1, "test_symmetrizeUpperSquareMatrix_1")
      93           1 :         call Test%finalize()
      94             : 
      95           1 :     end subroutine test_Matrix
      96             : 
      97             : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
      98             : 
      99             :     !> \brief
     100             :     !> Test whether [getEye()](@ref matrix_mod::geteye) return a correct output eye matrix.
     101           1 :     function test_eye_1() result(assertion)
     102             : 
     103           1 :         use Constants_mod, only: IK, RK
     104             :         implicit none
     105             :         logical                 :: assertion
     106             :         integer(IK) , parameter :: n = 3_IK, m = 4_IK
     107             :         real(RK)    , parameter :: Eye_ref(n, m) = reshape( [ 1._RK, 0._RK, 0._RK &
     108             :                                                             , 0._RK, 1._RK, 0._RK &
     109             :                                                             , 0._RK, 0._RK, 1._RK &
     110             :                                                             , 0._RK, 0._RK, 0._RK ], shape = shape(eye_ref) )
     111             :         real(RK), allocatable   :: Eye(:,:)
     112             :         real(RK), allocatable   :: Difference(:,:)
     113           1 :         Eye = getEye(n,m)
     114             : 
     115          17 :         Difference = abs(Eye - Eye_ref)
     116          17 :         assertion = all(Difference==0._RK)
     117             : 
     118           1 :         if (Test%isDebugMode .and. .not. assertion) then
     119             :         ! LCOV_EXCL_START
     120             :             write(Test%outputUnit,"(*(g0,:,', '))")
     121             :             write(Test%outputUnit,"(*(g0,:,', '))") "Eye_ref    = ", Eye_ref
     122             :             write(Test%outputUnit,"(*(g0,:,', '))") "Eye        = ", Eye
     123             :             write(Test%outputUnit,"(*(g0,:,', '))") "Difference = ", Difference
     124             :             write(Test%outputUnit,"(*(g0,:,', '))")
     125             :         end if
     126             :         ! LCOV_EXCL_STOP
     127             : 
     128           1 :     end function test_eye_1
     129             : 
     130             : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
     131             : 
     132             :     !> \brief
     133             :     !> Test whether [getEye()](@ref matrix_mod::geteye) return a correct output eye matrix with the input optional argument.
     134           1 :     function test_eye_2() result(assertion)
     135             : 
     136           1 :         use Constants_mod, only: IK, RK
     137             :         implicit none
     138             :         logical                 :: assertion
     139             :         integer(IK) , parameter :: n = 3_IK, m = 3_IK
     140             :         real(RK)    , parameter :: Eye_ref(n, m) = reshape( [ 2._RK, 0._RK, 0._RK &
     141             :                                                             , 0._RK, 2._RK, 0._RK &
     142             :                                                             , 0._RK, 0._RK, 2._RK ], shape = shape(eye_ref) )
     143             :         real(RK), allocatable   :: Eye(:,:)
     144             :         real(RK), allocatable   :: Difference(:,:)
     145           1 :         Eye = getEye(n,m, diag = 2._RK)
     146             : 
     147          13 :         Difference = abs(Eye - Eye_ref)
     148          13 :         assertion = all(Difference==0._RK)
     149             : 
     150           1 :         if (Test%isDebugMode .and. .not. assertion) then
     151             :         ! LCOV_EXCL_START
     152             :             write(Test%outputUnit,"(*(g0,:,', '))")
     153             :             write(Test%outputUnit,"(*(g0,:,', '))") "Eye_ref    = ", Eye_ref
     154             :             write(Test%outputUnit,"(*(g0,:,', '))") "Eye        = ", Eye
     155             :             write(Test%outputUnit,"(*(g0,:,', '))") "Difference = ", Difference
     156             :             write(Test%outputUnit,"(*(g0,:,', '))")
     157             :         end if
     158             :         ! LCOV_EXCL_STOP
     159             : 
     160           1 :     end function test_eye_2
     161             : 
     162             : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
     163             : 
     164           1 :     function test_getCholeskyFactor_1() result(assertion)
     165             : 
     166           1 :         use Constants_mod, only: IK, RK
     167             :         implicit none
     168             :         logical                 :: assertion
     169             :         integer(IK) , parameter :: nd = 3_IK
     170             :         real(RK)    , parameter :: tolerance = 1.e-12_RK
     171             :         real(RK)    , parameter :: PosDefMat(nd,nd) = reshape(  [ 1._RK, 0._RK, 1._RK &
     172             :                                                                 , 0._RK, 2._RK, 0._RK &
     173             :                                                                 , 1._RK, 0._RK, 3._RK ], shape = shape(PosDefMat) )
     174             :         real(RK)    , parameter :: CholeskyLower_ref(nd,nd) = reshape(  [ 1.000000000000000_RK, 0.000000000000000_RK, 1.000000000000000_RK &
     175             :                                                                         , 0.000000000000000_RK, 2.000000000000000_RK, 0.000000000000000_RK &
     176             :                                                                         , 1.000000000000000_RK, 0.000000000000000_RK, 3.000000000000000_RK ] &
     177             :                                                                         , shape = shape(CholeskyLower_ref) )
     178             :         real(RK)    , parameter :: CholeskyDiagonal_ref(nd) = [ 1.000000000000000_RK, 1.414213562373095_RK, 1.414213562373095_RK ]
     179             :         real(RK)                :: CholeskyLower(nd,nd), CholeskyDiagonal(nd)
     180             :         real(RK), allocatable   :: CholeskyLower_diff(:,:), CholeskyDiagonal_diff(:)
     181           1 :         CholeskyLower = PosDefMat
     182           1 :         call getCholeskyFactor(nd = nd, PosDefMat = CholeskyLower, Diagonal = CholeskyDiagonal)
     183             : 
     184             :         ! Gfortran 7.1 fails to automatically reallocate this array. This is not implemented in Gfortran 7.0.0
     185           1 :         if (allocated(CholeskyLower_diff)) deallocate(CholeskyLower_diff); allocate(CholeskyLower_diff, mold = PosDefMat)
     186          14 :         CholeskyLower_diff = abs(PosDefMat - CholeskyLower_ref)
     187             : 
     188             :         ! Gfortran 7.1 fails to automatically reallocate this array. This is not implemented in Gfortran 7.0.0
     189           1 :         if (allocated(CholeskyDiagonal_diff)) deallocate(CholeskyDiagonal_diff); allocate(CholeskyDiagonal_diff, mold = CholeskyDiagonal)
     190           5 :         CholeskyDiagonal_diff = abs(CholeskyDiagonal - CholeskyDiagonal_ref)
     191          16 :         assertion = all(CholeskyLower_diff < tolerance) .and. all(CholeskyDiagonal_diff < tolerance)
     192             : 
     193           1 :         if (Test%isDebugMode .and. .not. assertion) then
     194             :         ! LCOV_EXCL_START
     195             :             write(Test%outputUnit,"(*(g0,:,', '))")
     196             :             write(Test%outputUnit,"(*(g0,:,', '))") "CholeskyLower_ref  = ", CholeskyLower_ref
     197             :             write(Test%outputUnit,"(*(g0,:,', '))") "CholeskyLower      = ", CholeskyLower
     198             :             write(Test%outputUnit,"(*(g0,:,', '))") "CholeskyLower_diff = ", CholeskyLower_diff
     199             :             write(Test%outputUnit,"(*(g0,:,', '))")
     200             :             write(Test%outputUnit,"(*(g0,:,', '))") "CholeskyDiagonal_ref   = ", CholeskyDiagonal_ref
     201             :             write(Test%outputUnit,"(*(g0,:,', '))") "CholeskyDiagonal       = ", CholeskyDiagonal
     202             :             write(Test%outputUnit,"(*(g0,:,', '))") "CholeskyDiagonal_diff  = ", CholeskyDiagonal_diff
     203             :             write(Test%outputUnit,"(*(g0,:,', '))")
     204             :         end if
     205             :         ! LCOV_EXCL_STOP
     206             : 
     207           1 :     end function test_getCholeskyFactor_1
     208             : 
     209             : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
     210             : 
     211           1 :     function test_getCholeskyFactor_2() result(assertion)
     212           1 :         use Constants_mod, only: IK, RK
     213             :         implicit none
     214             :         logical                 :: assertion
     215             :         integer(IK) , parameter :: nd = 3_IK
     216             :         real(RK)    , parameter :: PosDefMat(nd,nd) = reshape(  [ 1._RK, 0._RK, -1._RK &
     217             :                                                                 , 0._RK, 2._RK, -0._RK &
     218             :                                                                 , 1._RK, 0._RK, -3._RK ], shape = shape(PosDefMat) )
     219             :         real(RK)                :: CholeskyLower(nd,nd), CholeskyDiagonal(nd)
     220           1 :         CholeskyLower = PosDefMat
     221           1 :         call getCholeskyFactor(nd = nd, PosDefMat = CholeskyLower, Diagonal = CholeskyDiagonal)
     222           1 :         assertion = CholeskyDiagonal(1) < 0._RK
     223           1 :         if (Test%isDebugMode .and. .not. assertion) then
     224             :         ! LCOV_EXCL_START
     225             :             write(Test%outputUnit,"(*(g0,:,', '))")
     226             :             write(Test%outputUnit,"(*(g0,:,', '))") "CholeskyLower      = ", CholeskyLower
     227             :             write(Test%outputUnit,"(*(g0,:,', '))") "CholeskyDiagonal   = ", CholeskyDiagonal
     228             :             write(Test%outputUnit,"(*(g0,:,', '))")
     229             :         end if
     230             :         ! LCOV_EXCL_STOP
     231           1 :     end function test_getCholeskyFactor_2
     232             : 
     233             : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
     234             : 
     235           1 :     function test_getInvPosDefMatSqrtDet_1() result(assertion)
     236             : 
     237           1 :         use Constants_mod, only: IK, RK
     238             :         implicit none
     239             : 
     240             :         logical                 :: assertion
     241             :         integer(IK) , parameter :: nd = 3_IK
     242             :         real(RK)    , parameter :: tolerance = 1.e-12_RK
     243             :         real(RK)    , parameter :: PosDefMat(nd,nd) = reshape(  [ 1._RK, 0._RK, 1._RK &
     244             :                                                                 , 0._RK, 2._RK, 0._RK &
     245             :                                                                 , 1._RK, 0._RK, 3._RK ], shape = shape(PosDefMat) )
     246             :         real(RK)    , parameter :: MatInvMat_ref(nd,nd) = reshape(  [ 1.500000000000000_RK, 0.000000000000000_RK, -0.50000000000000_RK &
     247             :                                                                     , 0.000000000000000_RK, 0.500000000000000_RK, 0.000000000000000_RK &
     248             :                                                                     , -0.50000000000000_RK, 0.000000000000000_RK, 0.500000000000000_RK ] &
     249             :                                                                     , shape = shape(MatInvMat_ref) )
     250             :         real(RK)    , parameter :: CholeskyDiagonal_ref(nd) = [ 1.000000000000000_RK, 1.414213562373095_RK, 1.414213562373095_RK ]
     251             :         real(RK)    , parameter :: sqrtDetInvPosDefMat_ref = 0.5_RK
     252           1 :         real(RK)                :: MatInvMat(nd,nd), sqrtDetInvPosDefMat
     253           1 :         real(RK), allocatable   :: MatInvMat_diff(:,:), sqrtDetInvPosDefMat_diff
     254             : 
     255           1 :         MatInvMat = PosDefMat
     256             : 
     257           1 :         call getInvPosDefMatSqrtDet(nd = nd, MatInvMat = MatInvMat, sqrtDetInvPosDefMat = sqrtDetInvPosDefMat)
     258             : 
     259             :         ! Gfortran 7.1 fails to automatically reallocate this array. This is not implemented in Gfortran 7.0.0
     260           1 :         if (allocated(MatInvMat_diff)) deallocate(MatInvMat_diff); allocate(MatInvMat_diff, mold = MatInvMat)
     261             : 
     262          14 :         MatInvMat_diff = abs(MatInvMat - MatInvMat_ref)
     263           1 :         sqrtDetInvPosDefMat_diff = abs(sqrtDetInvPosDefMat - sqrtDetInvPosDefMat_ref)
     264             : 
     265          13 :         assertion = all(MatInvMat_diff < tolerance) .and. sqrtDetInvPosDefMat_diff < tolerance
     266             : 
     267           1 :         if (Test%isDebugMode .and. .not. assertion) then
     268             :         ! LCOV_EXCL_START
     269             :             write(Test%outputUnit,"(*(g0,:,', '))")
     270             :             write(Test%outputUnit,"(*(g0,:,', '))") "MatInvMat_ref  = ", MatInvMat_ref
     271             :             write(Test%outputUnit,"(*(g0,:,', '))") "MatInvMat      = ", MatInvMat
     272             :             write(Test%outputUnit,"(*(g0,:,', '))") "MatInvMat_diff = ", MatInvMat_diff
     273             :             write(Test%outputUnit,"(*(g0,:,', '))")
     274             :             write(Test%outputUnit,"(*(g0,:,', '))") "sqrtDetInvPosDefMat_ref   = ", sqrtDetInvPosDefMat_ref
     275             :             write(Test%outputUnit,"(*(g0,:,', '))") "sqrtDetInvPosDefMat_diff  = ", sqrtDetInvPosDefMat
     276             :             write(Test%outputUnit,"(*(g0,:,', '))") "sqrtDetInvPosDefMat       = ", sqrtDetInvPosDefMat_diff
     277             :             write(Test%outputUnit,"(*(g0,:,', '))")
     278             :         end if
     279             :         ! LCOV_EXCL_STOP
     280             : 
     281           1 :     end function test_getInvPosDefMatSqrtDet_1
     282             : 
     283             : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
     284             : 
     285           1 :     function test_getInvPosDefMatSqrtDet_2() result(assertion)
     286           1 :         use Constants_mod, only: IK, RK
     287             :         implicit none
     288             :         logical                 :: assertion
     289             :         integer(IK) , parameter :: nd = 3_IK
     290             :         real(RK)    , parameter :: PosDefMat(nd,nd) = reshape(  [ 1._RK, 0._RK, -1._RK &
     291             :                                                                 , 0._RK, 2._RK, -0._RK &
     292             :                                                                 , 1._RK, 0._RK, -3._RK ], shape = shape(PosDefMat) )
     293           1 :         real(RK)                :: MatInvMat(nd,nd), sqrtDetInvPosDefMat
     294             : 
     295           1 :         MatInvMat = PosDefMat
     296             : 
     297           1 :         call getInvPosDefMatSqrtDet(nd = nd, MatInvMat = MatInvMat, sqrtDetInvPosDefMat = sqrtDetInvPosDefMat)
     298             : 
     299           1 :         assertion = sqrtDetInvPosDefMat < 0._RK
     300             : 
     301           1 :         if (Test%isDebugMode .and. .not. assertion) then
     302             :         ! LCOV_EXCL_START
     303             :             write(Test%outputUnit,"(*(g0,:,', '))")
     304             :             write(Test%outputUnit,"(*(g0,:,', '))") "MatInvMat              = ", MatInvMat
     305             :             write(Test%outputUnit,"(*(g0,:,', '))") "sqrtDetInvPosDefMat    = ", sqrtDetInvPosDefMat
     306             :             write(Test%outputUnit,"(*(g0,:,', '))")
     307             :         end if
     308             :         ! LCOV_EXCL_STOP
     309             : 
     310           1 :     end function test_getInvPosDefMatSqrtDet_2
     311             : 
     312             : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
     313             : 
     314             :     !> \brief
     315             :     !> Test with an 1-dimensional input matrix.
     316           1 :     function test_getInvPosDefMatSqrtDet_3() result(assertion)
     317           1 :         use Constants_mod, only: IK, RK
     318             :         implicit none
     319             :         logical                 :: assertion
     320             :         integer(IK) , parameter :: nd = 1_IK
     321             :         real(RK)    , parameter :: tolerance = 1.e-12_RK
     322             :         real(RK)    , parameter :: MatInvMat_ref(nd,nd) = reshape( [ 0.5_RK ], shape = shape(MatInvMat_ref) )
     323             :         real(RK)    , parameter :: PosDefMat(nd,nd) = reshape( [ 2._RK ], shape = shape(PosDefMat) )
     324             :         real(RK)    , parameter :: sqrtDetInvPosDefMat_ref = 0.5_RK
     325           3 :         real(RK)                :: MatInvMat(nd,nd), sqrtDetInvPosDefMat
     326           1 :         real(RK), allocatable   :: MatInvMat_diff(:,:), sqrtDetInvPosDefMat_diff
     327             : 
     328           1 :         MatInvMat = PosDefMat
     329             : 
     330           1 :         call getInvPosDefMatSqrtDet(nd = nd, MatInvMat = MatInvMat, sqrtDetInvPosDefMat = sqrtDetInvPosDefMat)
     331           1 :         if (allocated(MatInvMat_diff)) deallocate(MatInvMat_diff); allocate(MatInvMat_diff, mold = MatInvMat)
     332             : 
     333           4 :         MatInvMat_diff = abs(MatInvMat - MatInvMat_ref)
     334           1 :         sqrtDetInvPosDefMat_diff = abs(sqrtDetInvPosDefMat - sqrtDetInvPosDefMat_ref)
     335             : 
     336           3 :         assertion = all(MatInvMat_diff < tolerance) .and. sqrtDetInvPosDefMat_diff < tolerance
     337             : 
     338           1 :         if (Test%isDebugMode .and. .not. assertion) then
     339             :         ! LCOV_EXCL_START
     340             :             write(Test%outputUnit,"(*(g0,:,', '))")
     341             :             write(Test%outputUnit,"(*(g0,:,', '))") "MatInvMat_ref  = ", MatInvMat_ref
     342             :             write(Test%outputUnit,"(*(g0,:,', '))") "MatInvMat      = ", MatInvMat
     343             :             write(Test%outputUnit,"(*(g0,:,', '))") "MatInvMat_diff = ", MatInvMat_diff
     344             :             write(Test%outputUnit,"(*(g0,:,', '))")
     345             :             write(Test%outputUnit,"(*(g0,:,', '))") "sqrtDetInvPosDefMat_ref   = ", sqrtDetInvPosDefMat_ref
     346             :             write(Test%outputUnit,"(*(g0,:,', '))") "sqrtDetInvPosDefMat_diff  = ", sqrtDetInvPosDefMat
     347             :             write(Test%outputUnit,"(*(g0,:,', '))") "sqrtDetInvPosDefMat       = ", sqrtDetInvPosDefMat_diff
     348             :             write(Test%outputUnit,"(*(g0,:,', '))")
     349             :         end if
     350             :         ! LCOV_EXCL_STOP
     351             : 
     352             : 
     353           1 :     end function test_getInvPosDefMatSqrtDet_3
     354             : 
     355             : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
     356             : 
     357           1 :     function test_getInvMatFromCholFac() result(assertion)
     358             : 
     359           1 :         use Constants_mod, only: IK, RK
     360             :         implicit none
     361             : 
     362             :         logical                 :: assertion
     363             :         integer(IK) , parameter :: nd = 3_IK
     364             :         real(RK)    , parameter :: tolerance = 1.e-12_RK
     365             :         real(RK)    , parameter :: CholeskyLower(nd,nd) = reshape(  [ 1.000000000000000_RK, 0.000000000000000_RK, 1.000000000000000_RK &
     366             :                                                                     , 0.000000000000000_RK, 2.000000000000000_RK, 0.000000000000000_RK &
     367             :                                                                     , 1.000000000000000_RK, 0.000000000000000_RK, 3.000000000000000_RK ] &
     368             :                                                                     , shape = shape(CholeskyLower) )
     369             :         real(RK)    , parameter :: CholeskyDiagonal(nd) = [ 1.000000000000000_RK, 1.414213562373095_RK, 1.414213562373095_RK ]
     370             :         real(RK)    , parameter :: CholeskyDiagonal_ref(nd) = [ 1.000000000000000_RK, 1.414213562373095_RK, 1.414213562373095_RK ]
     371             :         real(RK)    , parameter :: InvMatFromCholFac_ref(nd,nd) = reshape(  [ 1.500000000000000_RK, 0.000000000000000_RK, -0.50000000000000_RK &
     372             :                                                                             , 0.000000000000000_RK, 0.500000000000000_RK, 0.000000000000000_RK &
     373             :                                                                             , -0.50000000000000_RK, 0.000000000000000_RK, 0.500000000000000_RK ] &
     374             :                                                                             , shape = shape(InvMatFromCholFac_ref) )
     375             :         real(RK)    , parameter :: sqrtDetInvPosDefMat_ref = 0.5_RK
     376             :         real(RK)                :: InvMatFromCholFac(nd,nd)
     377             :         real(RK), allocatable   :: InvMatFromCholFac_diff(:,:)
     378             : 
     379           1 :         InvMatFromCholFac = getInvMatFromCholFac(nd = nd, CholeskyLower = CholeskyLower, Diagonal = CholeskyDiagonal)
     380             : 
     381             :         ! Gfortran 7.1 fails to automatically reallocate this array. This is not implemented in Gfortran 7.0.0
     382           1 :         if (allocated(InvMatFromCholFac_diff)) deallocate(InvMatFromCholFac_diff); allocate(InvMatFromCholFac_diff, mold = InvMatFromCholFac)
     383             : 
     384          14 :         InvMatFromCholFac_diff = abs(InvMatFromCholFac - InvMatFromCholFac_ref)
     385             : 
     386          13 :         assertion = all(InvMatFromCholFac_diff < tolerance)
     387             : 
     388           1 :         if (Test%isDebugMode .and. .not. assertion) then
     389             :         ! LCOV_EXCL_START
     390             :             write(Test%outputUnit,"(*(g0,:,', '))")
     391             :             write(Test%outputUnit,"(*(g0,:,', '))") "InvMatFromCholFac_ref  = ", InvMatFromCholFac_ref
     392             :             write(Test%outputUnit,"(*(g0,:,', '))") "InvMatFromCholFac      = ", InvMatFromCholFac
     393             :             write(Test%outputUnit,"(*(g0,:,', '))") "InvMatFromCholFac_diff = ", InvMatFromCholFac_diff
     394             :             write(Test%outputUnit,"(*(g0,:,', '))")
     395             :         end if
     396             :         ! LCOV_EXCL_STOP
     397             : 
     398           1 :     end function test_getInvMatFromCholFac
     399             : 
     400             : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
     401             : 
     402           1 :     function test_getInvPosDefMat_1() result(assertion)
     403             : 
     404           1 :         use Constants_mod, only: IK, RK
     405             :         implicit none
     406             : 
     407             :         logical                 :: assertion
     408             :         integer(IK) , parameter :: nd = 3_IK
     409             :         real(RK)    , parameter :: tolerance = 1.e-12_RK
     410             :         real(RK)    , parameter :: PosDefMat(nd,nd) = reshape(  [ 1._RK, 0._RK, 1._RK &
     411             :                                                                 , 0._RK, 2._RK, 0._RK &
     412             :                                                                 , 1._RK, 0._RK, 3._RK ], shape = shape(PosDefMat) )
     413             :         real(RK)    , parameter :: MatInvMat_ref(nd,nd) = reshape(  [ 1.500000000000000_RK, 0.000000000000000_RK, -0.50000000000000_RK &
     414             :                                                                     , 0.000000000000000_RK, 0.500000000000000_RK, 0.000000000000000_RK &
     415             :                                                                     , -0.50000000000000_RK, 0.000000000000000_RK, 0.500000000000000_RK ] &
     416             :                                                                     , shape = shape(MatInvMat_ref) )
     417             :         real(RK)                :: MatInvMat(nd,nd)
     418             :         real(RK), allocatable   :: MatInvMat_diff(:,:)
     419             : 
     420           1 :         MatInvMat = getInvPosDefMat(nd = nd, PosDefMat = PosDefMat)
     421             : 
     422             :         ! Gfortran 7.1 fails to automatically reallocate this array. This is not implemented in Gfortran 7.0.0
     423           1 :         if (allocated(MatInvMat_diff)) deallocate(MatInvMat_diff); allocate(MatInvMat_diff, mold = MatInvMat)
     424          14 :         MatInvMat_diff = abs(MatInvMat - MatInvMat_ref)
     425             : 
     426          13 :         assertion = all(MatInvMat_diff < tolerance)
     427             : 
     428           1 :         if (Test%isDebugMode .and. .not. assertion) then
     429             :         ! LCOV_EXCL_START
     430             :             write(Test%outputUnit,"(*(g0,:,', '))")
     431             :             write(Test%outputUnit,"(*(g0,:,', '))") "MatInvMat_ref  = ", MatInvMat_ref
     432             :             write(Test%outputUnit,"(*(g0,:,', '))") "MatInvMat      = ", MatInvMat
     433             :             write(Test%outputUnit,"(*(g0,:,', '))") "MatInvMat_diff = ", MatInvMat_diff
     434             :             write(Test%outputUnit,"(*(g0,:,', '))")
     435             :         end if
     436             :         ! LCOV_EXCL_STOP
     437             : 
     438           1 :     end function test_getInvPosDefMat_1
     439             : 
     440             : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
     441             : 
     442           1 :     function test_getInvPosDefMat_2() result(assertion)
     443             : 
     444           1 :         use Constants_mod, only: IK, RK
     445             :         implicit none
     446             : 
     447             :         logical                 :: assertion
     448             :         integer(IK) , parameter :: nd = 3_IK
     449             :         real(RK)    , parameter :: PosDefMat(nd,nd) = reshape(  [ 1._RK, 0._RK, -1._RK &
     450             :                                                                 , 0._RK, 2._RK, -0._RK &
     451             :                                                                 , 1._RK, 0._RK, -3._RK ], shape = shape(PosDefMat) )
     452             :         real(RK)                :: MatInvMat(nd,nd)
     453             : 
     454           1 :         MatInvMat = getInvPosDefMat(nd = nd, PosDefMat = PosDefMat)
     455             : 
     456           1 :         assertion = MatInvMat(1,1) < 0._RK
     457             : 
     458           1 :         if (Test%isDebugMode .and. .not. assertion) then
     459             :         ! LCOV_EXCL_START
     460             :             write(Test%outputUnit,"(*(g0,:,', '))")
     461             :             write(Test%outputUnit,"(*(g0,:,', '))") "MatInvMat      = ", MatInvMat
     462             :             write(Test%outputUnit,"(*(g0,:,', '))")
     463             :         end if
     464             :         ! LCOV_EXCL_STOP
     465             : 
     466           1 :     end function test_getInvPosDefMat_2
     467             : 
     468             : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
     469             : 
     470           1 :     function test_getInvMatDet_1() result(assertion)
     471             : 
     472           1 :         use Constants_mod, only: IK, RK
     473             :         implicit none
     474             : 
     475             :         logical                 :: assertion
     476             :         integer(IK) , parameter :: nd = 3_IK
     477             :         real(RK)    , parameter :: tolerance = 1.e-12_RK
     478             :         real(RK)    , parameter :: PosDefMat(nd,nd) = reshape(  [ 1._RK, 0._RK, 1._RK &
     479             :                                                                 , 0._RK, 2._RK, 0._RK &
     480             :                                                                 , 1._RK, 0._RK, 3._RK ], shape = shape(PosDefMat) )
     481             :         real(RK)    , parameter :: MatInvMat_ref(nd,nd) = reshape(  [ 1.500000000000000_RK, 0.000000000000000_RK, -0.50000000000000_RK &
     482             :                                                                     , 0.000000000000000_RK, 0.500000000000000_RK, 0.000000000000000_RK &
     483             :                                                                     , -0.50000000000000_RK, 0.000000000000000_RK, 0.500000000000000_RK ] &
     484             :                                                                     , shape = shape(MatInvMat_ref) )
     485             :         real(RK)    , parameter :: MatrixLU_ref(nd,nd)  = reshape(  [ 1.000000000000000_RK, 0.000000000000000_RK, 1.000000000000000_RK &
     486             :                                                                     , 0.000000000000000_RK, 2.000000000000000_RK, 0.000000000000000_RK &
     487             :                                                                     , 1.000000000000000_RK, 0.000000000000000_RK, 2.000000000000000_RK ] &
     488             :                                                                     , shape = shape(MatInvMat_ref) )
     489             :         real(RK)    , parameter :: detInvMat_ref = 0.25_RK
     490           1 :         real(RK)                :: MatInvMat(nd,nd), detInvMat, detInvMat_diff
     491             :         real(RK), allocatable   :: MatrixLU(:,:), MatInvMat_diff(:,:), MatrixLU_diff(:,:)
     492             : 
     493          13 :         MatrixLU = PosDefMat
     494             : 
     495           1 :         call getInvMatDet(nd = nd, MatrixLU = MatrixLU, InverseMatrix = MatInvMat, detInvMat = detInvMat)
     496             : 
     497             :         ! Gfortran 7.1 fails to automatically reallocate this array. This is not implemented in Gfortran 7.0.0
     498           1 :         if (allocated(MatrixLU_diff)) deallocate(MatrixLU_diff); allocate(MatrixLU_diff, mold = MatrixLU)
     499          14 :         MatrixLU_diff = abs(MatrixLU - MatrixLU_ref)
     500             : 
     501             :         ! Gfortran 7.1 fails to automatically reallocate this array. This is not implemented in Gfortran 7.0.0
     502           1 :         if (allocated(MatInvMat_diff)) deallocate(MatInvMat_diff); allocate(MatInvMat_diff, mold = MatInvMat)
     503          14 :         MatInvMat_diff = abs(MatInvMat - MatInvMat_ref)
     504             : 
     505           1 :         detInvMat_diff = abs(detInvMat - detInvMat_ref)
     506             : 
     507          13 :         assertion = all(MatInvMat_diff < tolerance) .and. detInvMat_diff < tolerance
     508             : 
     509           1 :         if (Test%isDebugMode .and. .not. assertion) then
     510             :         ! LCOV_EXCL_START
     511             :             write(Test%outputUnit,"(*(g0,:,', '))")
     512             :             write(Test%outputUnit,"(*(g0,:,', '))") "MatrixLU_ref  = ", MatrixLU_ref
     513             :             write(Test%outputUnit,"(*(g0,:,', '))") "MatrixLU      = ", MatrixLU
     514             :             write(Test%outputUnit,"(*(g0,:,', '))") "MatrixLU_diff = ", MatrixLU_diff
     515             :             write(Test%outputUnit,"(*(g0,:,', '))")
     516             :             write(Test%outputUnit,"(*(g0,:,', '))") "MatInvMat_ref  = ", MatInvMat_ref
     517             :             write(Test%outputUnit,"(*(g0,:,', '))") "MatInvMat      = ", MatInvMat
     518             :             write(Test%outputUnit,"(*(g0,:,', '))") "MatInvMat_diff = ", MatInvMat_diff
     519             :             write(Test%outputUnit,"(*(g0,:,', '))")
     520             :             write(Test%outputUnit,"(*(g0,:,', '))") "detInvMat_ref  = ", detInvMat_ref
     521             :             write(Test%outputUnit,"(*(g0,:,', '))") "detInvMat      = ", detInvMat
     522             :             write(Test%outputUnit,"(*(g0,:,', '))") "detInvMat_diff = ", detInvMat_diff
     523             :             write(Test%outputUnit,"(*(g0,:,', '))")
     524             :         end if
     525             :         ! LCOV_EXCL_STOP
     526             : 
     527           1 :     end function test_getInvMatDet_1
     528             : 
     529             : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
     530             : 
     531           1 :     function test_getInvMat_1() result(assertion)
     532             : 
     533           1 :         use Constants_mod, only: IK, RK
     534             :         implicit none
     535             : 
     536             :         logical                 :: assertion
     537             :         integer(IK) , parameter :: nd = 3_IK
     538             :         real(RK)    , parameter :: tolerance = 1.e-12_RK
     539             :         real(RK)    , parameter :: PosDefMat(nd,nd) = reshape(  [ 1._RK, 0._RK, 1._RK &
     540             :                                                                 , 0._RK, 2._RK, 0._RK &
     541             :                                                                 , 1._RK, 0._RK, 3._RK ], shape = shape(PosDefMat) )
     542             :         real(RK)    , parameter :: InverseMatrix_ref(nd,nd) = reshape(  [ 1.500000000000000_RK, 0.000000000000000_RK, -0.50000000000000_RK &
     543             :                                                                         , 0.000000000000000_RK, 0.500000000000000_RK, 0.000000000000000_RK &
     544             :                                                                         , -0.50000000000000_RK, 0.000000000000000_RK, 0.500000000000000_RK ] &
     545             :                                                                         , shape = shape(InverseMatrix_ref) )
     546             :         real(RK)                :: InverseMatrix(nd,nd)
     547             :         real(RK), allocatable   :: InverseMatrix_diff(:,:)
     548             : 
     549           1 :         InverseMatrix = getInvMat(nd = nd, Matrix = PosDefMat)
     550             : 
     551             :         ! Gfortran 7.1 fails to automatically reallocate this array. This is not implemented in Gfortran 7.0.0
     552           1 :         if (allocated(InverseMatrix_diff)) deallocate(InverseMatrix_diff); allocate(InverseMatrix_diff, mold = InverseMatrix)
     553             : 
     554          14 :         InverseMatrix_diff = abs(InverseMatrix - InverseMatrix_ref)
     555             : 
     556          13 :         assertion = all(InverseMatrix_diff < tolerance)
     557             : 
     558           1 :         if (Test%isDebugMode .and. .not. assertion) then
     559             :         ! LCOV_EXCL_START
     560             :             write(Test%outputUnit,"(*(g0,:,', '))")
     561             :             write(Test%outputUnit,"(*(g0,:,', '))") "InverseMatrix_ref  = ", InverseMatrix_ref
     562             :             write(Test%outputUnit,"(*(g0,:,', '))") "InverseMatrix      = ", InverseMatrix
     563             :             write(Test%outputUnit,"(*(g0,:,', '))") "InverseMatrix_diff = ", InverseMatrix_diff
     564             :             write(Test%outputUnit,"(*(g0,:,', '))")
     565             :         end if
     566             :         ! LCOV_EXCL_STOP
     567             : 
     568           1 :     end function test_getInvMat_1
     569             : 
     570             : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
     571             : 
     572           1 :     function test_multiplyMatrix_1() result(assertion)
     573             : 
     574           1 :         use Constants_mod, only: IK, RK
     575             :         implicit none
     576             : 
     577             :         logical                 :: assertion
     578             :         integer(IK) , parameter :: nd = 3_IK
     579             :         real(RK)    , parameter :: tolerance = 1.e-12_RK
     580             :         real(RK)    , parameter :: PosDefMat(nd,nd) = reshape(  [ 1._RK, 0._RK, 1._RK &
     581             :                                                                 , 0._RK, 2._RK, 0._RK &
     582             :                                                                 , 1._RK, 0._RK, 3._RK ], shape = shape(PosDefMat) )
     583             :         real(RK)                :: MatrixProduct(nd,nd)
     584             :         real(RK), allocatable   :: MatrixProduct_diff(:,:), MatrixProduct_ref(:,:)
     585             : 
     586          13 :         MatrixProduct_ref = matmul(PosDefMat,PosDefMat)
     587             : 
     588           1 :         call multiplyMatrix(A = PosDefMat, rowsA = nd, colsA = nd, B = PosDefMat, rowsB = nd, colsB = nd, C = MatrixProduct)
     589             : 
     590             :         ! Gfortran 7.1 fails to automatically reallocate this array. This is not implemented in Gfortran 7.0.0
     591           1 :         if (allocated(MatrixProduct_diff)) deallocate(MatrixProduct_diff); allocate(MatrixProduct_diff, mold = MatrixProduct)
     592             : 
     593          14 :         MatrixProduct_diff = abs(MatrixProduct - MatrixProduct_ref)
     594             : 
     595          13 :         assertion = all(MatrixProduct_diff < tolerance)
     596             : 
     597           1 :         if (Test%isDebugMode .and. .not. assertion) then
     598             :         ! LCOV_EXCL_START
     599             :             write(Test%outputUnit,"(*(g0,:,', '))")
     600             :             write(Test%outputUnit,"(*(g0,:,', '))") "MatrixProduct_ref  = ", MatrixProduct_ref
     601             :             write(Test%outputUnit,"(*(g0,:,', '))") "MatrixProduct      = ", MatrixProduct
     602             :             write(Test%outputUnit,"(*(g0,:,', '))") "MatrixProduct_diff = ", MatrixProduct_diff
     603             :             write(Test%outputUnit,"(*(g0,:,', '))")
     604             :         end if
     605             :         ! LCOV_EXCL_STOP
     606             : 
     607           1 :     end function test_multiplyMatrix_1
     608             : 
     609             : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
     610             : 
     611           1 :     function test_getDeterminant_1() result(assertion)
     612             : 
     613           1 :         use Constants_mod, only: IK, RK
     614             :         implicit none
     615             : 
     616             :         logical                 :: assertion
     617             :         integer(IK) , parameter :: nd = 3_IK
     618             :         real(RK)    , parameter :: tolerance = 1.e-12_RK
     619             :         real(RK)    , parameter :: PosDefMat(nd,nd) = reshape(  [ -1._RK, -0._RK, -1._RK &
     620             :                                                                 , -0._RK, -2._RK, -0._RK &
     621             :                                                                 , -1._RK, -0._RK, -3._RK ], shape = shape(PosDefMat) )
     622             :         real(RK)    , parameter :: determinant_ref = -4._RK
     623           1 :         real(RK)                :: determinant, determinant_diff
     624             : 
     625           2 :         determinant = getDeterminant(nd = nd, Matrix = PosDefMat)
     626             : 
     627           1 :         determinant_diff = abs(determinant - determinant_ref)
     628             : 
     629           1 :         assertion = determinant_diff < tolerance
     630             : 
     631           1 :         if (Test%isDebugMode .and. .not. assertion) then
     632             :         ! LCOV_EXCL_START
     633             :             write(Test%outputUnit,"(*(g0,:,', '))")
     634             :             write(Test%outputUnit,"(*(g0,:,', '))") "determinant_ref  = ", determinant_ref
     635             :             write(Test%outputUnit,"(*(g0,:,', '))") "determinant      = ", determinant
     636             :             write(Test%outputUnit,"(*(g0,:,', '))") "determinant_diff = ", determinant_diff
     637             :             write(Test%outputUnit,"(*(g0,:,', '))")
     638             :         end if
     639             :         ! LCOV_EXCL_STOP
     640             : 
     641           1 :     end function test_getDeterminant_1
     642             : 
     643             : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
     644             : 
     645           1 :     function test_getSqrtDetPosDefMat_1() result(assertion)
     646             : 
     647           1 :         use Constants_mod, only: IK, RK
     648             :         implicit none
     649             : 
     650             :         logical                 :: assertion
     651             :         integer(IK) , parameter :: nd = 3_IK
     652             :         real(RK)    , parameter :: tolerance = 1.e-12_RK
     653             :         real(RK)    , parameter :: PosDefMat(nd,nd) = reshape(  [ 1._RK, 0._RK, 1._RK &
     654             :                                                                 , 0._RK, 2._RK, 0._RK &
     655             :                                                                 , 1._RK, 0._RK, 3._RK ], shape = shape(PosDefMat) )
     656             :         real(RK)    , parameter :: sqrtDetPosDefMat_ref = 2._RK
     657           1 :         real(RK)                :: sqrtDetPosDefMat, sqrtDetPosDefMat_diff
     658             : 
     659           1 :         sqrtDetPosDefMat = getSqrtDetPosDefMat(nd = nd, PosDefMat = PosDefMat)
     660             : 
     661           1 :         sqrtDetPosDefMat_diff = abs(sqrtDetPosDefMat - sqrtDetPosDefMat_ref)
     662             : 
     663           1 :         assertion = sqrtDetPosDefMat_diff < tolerance
     664             : 
     665           1 :         if (Test%isDebugMode .and. .not. assertion) then
     666             :         ! LCOV_EXCL_START
     667             :             write(Test%outputUnit,"(*(g0,:,', '))")
     668             :             write(Test%outputUnit,"(*(g0,:,', '))") "sqrtDetPosDefMat_ref  = ", sqrtDetPosDefMat_ref
     669             :             write(Test%outputUnit,"(*(g0,:,', '))") "sqrtDetPosDefMat      = ", sqrtDetPosDefMat
     670             :             write(Test%outputUnit,"(*(g0,:,', '))") "sqrtDetPosDefMat_diff = ", sqrtDetPosDefMat_diff
     671             :             write(Test%outputUnit,"(*(g0,:,', '))")
     672             :         end if
     673             :         ! LCOV_EXCL_STOP
     674             : 
     675           1 :     end function test_getSqrtDetPosDefMat_1
     676             : 
     677             : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
     678             : 
     679           1 :     function test_getSqrtDetPosDefMat_2() result(assertion)
     680             : 
     681           1 :         use Constants_mod, only: IK, RK
     682             :         implicit none
     683             : 
     684             :         logical                 :: assertion
     685             :         integer(IK) , parameter :: nd = 3_IK
     686             :         real(RK)    , parameter :: tolerance = 1.e-12_RK
     687             :         real(RK)    , parameter :: PosDefMat(nd,nd) = reshape(  [ -1._RK, -0._RK, -1._RK &
     688             :                                                                 , -0._RK, -2._RK, -0._RK &
     689             :                                                                 , -1._RK, -0._RK, -3._RK ], shape = shape(PosDefMat) )
     690           1 :         real(RK)                :: sqrtDetPosDefMat
     691             : 
     692           1 :         sqrtDetPosDefMat = getSqrtDetPosDefMat(nd = nd, PosDefMat = PosDefMat)
     693             : 
     694           1 :         assertion = sqrtDetPosDefMat < 0._RK
     695             : 
     696           1 :         if (Test%isDebugMode .and. .not. assertion) then
     697             :         ! LCOV_EXCL_START
     698             :             write(Test%outputUnit,"(*(g0,:,', '))")
     699             :             write(Test%outputUnit,"(*(g0,:,', '))") "sqrtDetPosDefMat      = ", sqrtDetPosDefMat
     700             :             write(Test%outputUnit,"(*(g0,:,', '))")
     701             :         end if
     702             :         ! LCOV_EXCL_STOP
     703             : 
     704           1 :     end function test_getSqrtDetPosDefMat_2
     705             : 
     706             : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
     707             : 
     708           1 :     function test_getLogSqrtDetPosDefMat_1() result(assertion)
     709             : 
     710           1 :         use Constants_mod, only: IK, RK
     711             :         implicit none
     712             : 
     713             :         logical                 :: assertion
     714             :         integer(IK) , parameter :: nd = 3_IK
     715             :         real(RK)    , parameter :: tolerance = 1.e-12_RK
     716             :         real(RK)    , parameter :: PosDefMat(nd,nd) = reshape(  [ 1._RK, 0._RK, 1._RK &
     717             :                                                                 , 0._RK, 2._RK, 0._RK &
     718             :                                                                 , 1._RK, 0._RK, 3._RK ], shape = shape(PosDefMat) )
     719             :         real(RK)    , parameter :: logSqrtDetPosDefMat_ref = log(2._RK)
     720           1 :         real(RK)                :: logSqrtDetPosDefMat, logSqrtDetPosDefMat_diff
     721             :         real(RK), allocatable   :: Matrix(:,:)
     722             :         logical                 :: failed
     723             : 
     724          13 :         Matrix = PosDefMat
     725             : 
     726           1 :         call getLogSqrtDetPosDefMat(nd = nd, PosDefMat = Matrix, logSqrtDetPosDefMat = logSqrtDetPosDefMat, failed = failed)
     727             : 
     728           1 :         assertion = .not. failed
     729           1 :         if (.not. assertion) return
     730             : 
     731           1 :         logSqrtDetPosDefMat_diff = abs(logSqrtDetPosDefMat - logSqrtDetPosDefMat_ref)
     732             : 
     733           1 :         assertion = logSqrtDetPosDefMat_diff < tolerance
     734             : 
     735           1 :         if (Test%isDebugMode .and. .not. assertion) then
     736             :         ! LCOV_EXCL_START
     737             :             write(Test%outputUnit,"(*(g0,:,', '))")
     738             :             write(Test%outputUnit,"(*(g0,:,', '))") "logSqrtDetPosDefMat_ref  = ", logSqrtDetPosDefMat_ref
     739             :             write(Test%outputUnit,"(*(g0,:,', '))") "logSqrtDetPosDefMat      = ", logSqrtDetPosDefMat
     740             :             write(Test%outputUnit,"(*(g0,:,', '))") "logSqrtDetPosDefMat_diff = ", logSqrtDetPosDefMat_diff
     741             :             write(Test%outputUnit,"(*(g0,:,', '))")
     742             :         end if
     743             :         ! LCOV_EXCL_STOP
     744             : 
     745           1 :     end function test_getLogSqrtDetPosDefMat_1
     746             : 
     747             : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
     748             : 
     749           1 :     function test_getLogSqrtDetPosDefMat_2() result(assertion)
     750             : 
     751           1 :         use Constants_mod, only: IK, RK
     752             :         implicit none
     753             : 
     754             :         logical                 :: assertion
     755             :         integer(IK) , parameter :: nd = 3_IK
     756             :         real(RK)    , parameter :: tolerance = 1.e-12_RK
     757             :         real(RK)    , parameter :: PosDefMat(nd,nd) = reshape(  [ 1._RK, 0._RK, 1._RK &
     758             :                                                                 , 0._RK, 2._RK, 0._RK &
     759             :                                                                 , 1._RK, 0._RK, 3._RK ], shape = shape(PosDefMat) )
     760             :         real(RK)    , parameter :: logSqrtDetPosDefMat_ref = log(2._RK)
     761           1 :         real(RK)                :: logSqrtDetPosDefMat
     762             :         real(RK), allocatable   :: Matrix(:,:)
     763             :         logical                 :: failed
     764             : 
     765          13 :         Matrix = -PosDefMat
     766             : 
     767           1 :         call getLogSqrtDetPosDefMat(nd = nd, PosDefMat = Matrix, logSqrtDetPosDefMat = logSqrtDetPosDefMat, failed = failed)
     768             : 
     769           1 :         assertion = failed
     770           1 :         if (.not. assertion) return
     771             : 
     772           1 :     end function test_getLogSqrtDetPosDefMat_2
     773             : 
     774             : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
     775             : 
     776           1 :     function test_isPosDef_1() result(assertion)
     777           1 :         use Constants_mod, only: IK, RK
     778             :         implicit none
     779             :         logical                 :: assertion
     780             :         integer(IK) , parameter :: nd = 3_IK
     781             :         real(RK)    , parameter :: PosDefMat(nd,nd) = reshape(  [ 1._RK, 0._RK, 1._RK &
     782             :                                                                 , 0._RK, 2._RK, 0._RK &
     783             :                                                                 , 1._RK, 0._RK, 3._RK ], shape = shape(PosDefMat) )
     784             : 
     785           1 :         assertion = isPosDef(nd = nd, Matrix = PosDefMat)
     786           2 :     end function test_isPosDef_1
     787             : 
     788             : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
     789             : 
     790           1 :     function test_isPosDef_2() result(assertion)
     791           1 :         use Constants_mod, only: IK, RK
     792             :         implicit none
     793             :         logical                 :: assertion
     794             :         integer(IK) , parameter :: nd = 3_IK
     795             :         real(RK)    , parameter :: PosDefMat(nd,nd) = reshape(  [ -1._RK, -0._RK, -1._RK &
     796             :                                                                 , -0._RK, -2._RK, -0._RK &
     797             :                                                                 , -1._RK, -0._RK, -3._RK ], shape = shape(PosDefMat) )
     798             : 
     799           1 :         assertion = .not. isPosDef(nd = nd, Matrix = PosDefMat)
     800           2 :     end function test_isPosDef_2
     801             : 
     802             : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
     803             : 
     804           1 :     function test_symmetrizeUpperSquareMatrix_1() result(assertion)
     805           1 :         use Constants_mod, only: IK, RK
     806             :         implicit none
     807             :         logical                 :: assertion
     808             :         integer(IK) , parameter :: nd = 3_IK
     809             :         real(RK)    , parameter :: UpperSquareMatrix_ref(nd,nd) = reshape(  [ -1._RK, -0._RK, -1._RK &
     810             :                                                                             , -0._RK, -2._RK, -0._RK &
     811             :                                                                             , -1._RK, -0._RK, -3._RK ], shape = shape(UpperSquareMatrix_ref) )
     812             :         real(RK)                :: UpperSquareMatrix(nd,nd)
     813             :         integer(IK)             :: i, j
     814             : 
     815           1 :         UpperSquareMatrix = 0._RK
     816           4 :         do j = 1, nd
     817          10 :             do i = 1, j
     818           9 :                 UpperSquareMatrix(i,j) = UpperSquareMatrix_ref(i,j)
     819             :             end do
     820             :         end do
     821             : 
     822           1 :         call symmetrizeUpperSquareMatrix(nd,UpperSquareMatrix)
     823          13 :         assertion = all( abs(UpperSquareMatrix-UpperSquareMatrix_ref) < 1.e-14_RK )
     824             : 
     825           1 :         if (Test%isDebugMode .and. .not. assertion) then
     826             :         ! LCOV_EXCL_START
     827             :             write(Test%outputUnit,"(*(g0,:,', '))")
     828             :             write(Test%outputUnit,"(*(g0,:,', '))") "UpperSquareMatrix_ref  = ", UpperSquareMatrix_ref
     829             :             write(Test%outputUnit,"(*(g0,:,', '))") "UpperSquareMatrix      = ", UpperSquareMatrix
     830             :             write(Test%outputUnit,"(*(g0,:,', '))")
     831             :         end if
     832             :         ! LCOV_EXCL_STOP
     833             : 
     834           1 :     end function test_symmetrizeUpperSquareMatrix_1
     835             : 
     836             : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
     837             : 
     838           1 :     function test_getOuterProd_1() result(assertion)
     839             : 
     840           1 :         use Constants_mod, only: IK, RK
     841             :         implicit none
     842             : 
     843             :         logical                 :: assertion
     844             :         real(RK)    , parameter :: tolerance = 1.e-12_RK
     845             :         real(RK)    , parameter :: Vector2(*) = [4._RK, 5._RK]
     846             :         real(RK)    , parameter :: Vector1(*) = [1._RK, 2._RK, 3._RK]
     847             :         real(RK)    , parameter :: OuterProduct_ref(3,2) = reshape( [ 4._RK,  8._RK, 12._RK &
     848             :                                                                     , 5._RK, 10._RK, 15._RK ], shape = shape(OuterProduct_ref) )
     849             :         real(RK), allocatable   :: OuterProduct(:,:), OuterProduct_diff(:,:)
     850             : 
     851           1 :         OuterProduct = getOuterProd(Vector1 = Vector1, Vector2 = Vector2)
     852             : 
     853             :         ! Gfortran 7.1 fails to automatically reallocate this array. This is not implemented in Gfortran 7.0.0
     854           1 :         if (allocated(OuterProduct_diff)) deallocate(OuterProduct_diff); allocate(OuterProduct_diff, mold = OuterProduct)
     855             : 
     856          10 :         OuterProduct_diff = abs(OuterProduct - OuterProduct_ref)
     857             : 
     858           9 :         assertion = all(OuterProduct_diff < tolerance)
     859             : 
     860           1 :         if (Test%isDebugMode .and. .not. assertion) then
     861             :         ! LCOV_EXCL_START
     862             :             write(Test%outputUnit,"(*(g0,:,', '))")
     863             :             write(Test%outputUnit,"(*(g0,:,', '))") "OuterProduct_ref  = ", OuterProduct_ref
     864             :             write(Test%outputUnit,"(*(g0,:,', '))") "OuterProduct      = ", OuterProduct
     865             :             write(Test%outputUnit,"(*(g0,:,', '))") "OuterProduct_diff = ", OuterProduct_diff
     866             :             write(Test%outputUnit,"(*(g0,:,', '))")
     867             :         end if
     868             :         ! LCOV_EXCL_STOP
     869             : 
     870           1 :     end function test_getOuterProd_1
     871             : 
     872             : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
     873             : 
     874           1 :     function test_sortPosDefMat_1() result(assertion)
     875             : 
     876           1 :         use Constants_mod, only: RK, IK
     877             : 
     878             :         implicit none
     879             :         logical                 :: assertion
     880             :         integer(IK), parameter  :: ColIndx(*) = [4]
     881             :         integer(IK), parameter  :: ColIndxMap(*) = [5]
     882             :         integer(IK)             :: rank
     883             :         real(RK), allocatable   :: PosDefMat(:,:), OutPosDefMat(:,:), OutPosDefMat_ref(:,:)
     884             :         integer(IK)             :: i, j
     885             : 
     886           1 :         assertion = .true.
     887             : 
     888           1 :         rank = 5
     889           1 :         allocate(PosDefMat(rank,rank))
     890           6 :         do j = 1,rank
     891          21 :             do i = 1,j
     892          20 :                 PosDefMat(i,j) = i*10 + j
     893             :             end do
     894             :         end do
     895             : 
     896             :         ! switch the variable 4 with 5, such that the output remains a positive-definite matrix.
     897             : 
     898           1 :         OutPosDefMat = sortPosDefMat(rank, PosDefMat, 1_IK, ColIndx, ColIndxMap)
     899             : 
     900             :         OutPosDefMat_ref = reshape(  [ 11._RK, 12._RK, 13._RK, 15._RK, 14._RK &
     901             :                                     ,  0._RK, 22._RK, 23._RK, 25._RK, 24._RK &
     902             :                                     ,  0._RK,  0._RK, 33._RK, 35._RK, 34._RK &
     903             :                                     ,  0._RK,  0._RK,  0._RK, 55._RK, 45._RK &
     904             :                                     ,  0._RK,  0._RK,  0._RK,  0._RK, 44._RK ] &
     905           3 :                                     , shape = [rank,rank] )
     906          62 :         OutPosDefMat_ref = transpose(OutPosDefMat_ref)
     907           6 :         do j = 1, rank
     908          21 :             do i = 1, j
     909          20 :                 assertion = assertion .and. OutPosDefMat_ref(i,j) == OutPosDefMat(i,j)
     910             :             end do
     911             :         end do
     912             : 
     913           1 :         if (Test%isDebugMode .and. .not. assertion) then
     914             :         ! LCOV_EXCL_START
     915             : 
     916             :             write(Test%outputUnit,"(*(g0))")
     917             :             write(Test%outputUnit,"(*(g0))") "OutPosDefMat_ref:"
     918             :             do i = 1,rank
     919             :                 write(Test%outputUnit,"(*(F7.1))") (OutPosDefMat_ref(i,j),j=1,rank)
     920             :             end do
     921             :             write(Test%outputUnit,"(*(g0))")
     922             : 
     923             :             write(Test%outputUnit,"(*(g0))")
     924             :             write(Test%outputUnit,"(*(g0))") "Output Positive-Definite Matrix with variables 4 and 5 swapped:"
     925             :             do i = 1,rank
     926             :                 write(Test%outputUnit,"(*(F7.1))") (OutPosDefMat(i,j),j=1,rank)
     927             :             end do
     928             :             write(Test%outputUnit,"(*(g0))")
     929             : 
     930             :             write(Test%outputUnit,"(*(g0))")
     931             :             write(Test%outputUnit,"(*(g0))") "Original Matrix:"
     932             :             do i = 1,rank
     933             :                 write(Test%outputUnit,"(*(F7.1))") (PosDefMat(i,j),j=1,rank)
     934             :             end do
     935             : 
     936             :         end if
     937             :         ! LCOV_EXCL_STOP
     938             : 
     939           1 :     end function test_sortPosDefMat_1
     940             : 
     941             : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
     942             : 
     943           1 :     function test_getRegresCoef_1() result(assertion)
     944             : 
     945           1 :         use Constants_mod, only: RK, IK
     946             : 
     947             :         implicit none
     948             :         logical                 ::  assertion
     949           1 :         real(RK)                ::  normalizedDifference
     950             :         integer , parameter     ::  rankPDM = 4, rankS11 = 3, rankS22 = 1
     951             :         real(RK), parameter     ::  PosDefMat(rankPDM,rankPDM) = reshape( &
     952             :                                     [ 4.414182620515998_RK, 1.173760167060120_RK, 0.757607629189287_RK, 5.075277296976230_RK &
     953             :                                     , 1.173760167060120_RK, 0.866956750091570_RK, 0.310654936099342_RK, 1.621274787164182_RK &
     954             :                                     , 0.757607629189287_RK, 0.310654936099342_RK, 0.955157221699132_RK, 1.254186231887444_RK &
     955             :                                     , 5.075277296976230_RK, 1.621274787164182_RK, 1.254186231887444_RK, 6.407791808961157_RK ] &
     956             :                                     , shape=shape(PosDefMat) )
     957             :         real(RK), parameter     ::  RegresCoefMatRef(rankS11,rankS22) = reshape( &
     958             :                                     [ 0.792047783119072 &
     959             :                                     , 0.253016145889269 &
     960             :                                     , 0.195728305363088 ] &
     961             :                                     , shape=shape(RegresCoefMatRef) )
     962             :         real(RK), parameter     ::  SchurComplementRef(rankS11,rankS11) = reshape( &
     963             :                                     [ 0.3943204887314210_RK, -0.110366933940115_RK, -0.235767795395625_RK &
     964             :                                     , -0.110366933940115_RK,  0.456748052015843_RK, -0.006674430520204_RK &
     965             :                                     , -0.235767795395625_RK, -0.006674430520204_RK,  0.709677475922085_RK ] &
     966             :                                     , shape=shape(SchurComplementRef) )
     967             :         real(RK)                ::  SchurComplement(rankS11,rankS11)
     968             :         real(RK)                ::  RegresCoefMat(rankS11,rankS22)
     969             :         integer                 ::  i,j
     970             : 
     971           1 :         assertion = .true.
     972             : 
     973             :         call getRegresCoef  ( rankPDM           = rankPDM           &
     974             :                             , rankS11           = rankS11           &
     975             :                             , rankS22           = rankS22           &
     976             :                             , PosDefMat         = PosDefMat         &
     977             :                             , RegresCoefMat     = RegresCoefMat     &
     978             :                             , SchurComplement   = SchurComplement   &
     979           1 :                             )
     980             : 
     981           4 :         do i = 1,rankS11
     982             :             !write(Test%outputUnit,"(*(F22.15))") (RegresCoefMat(i,j),j=1,rankS22), (RegresCoefMatRef(i,j),j=1,rankS22)
     983           7 :             do j = 1,rankS22
     984           3 :                 normalizedDifference = abs(RegresCoefMatRef(i,j) - RegresCoefMat(i,j)) / (RegresCoefMatRef(i,j) + RegresCoefMat(i,j))
     985           6 :                 assertion = assertion .and. normalizedDifference < 1.e-5_RK
     986             :             end do
     987             :         end do
     988             : 
     989           1 :         if (Test%isDebugMode .and. .not. assertion) then
     990             :         ! LCOV_EXCL_START
     991             :             write(Test%outputUnit,"(*(g0))")
     992             :             write(Test%outputUnit,"(*(g0))") "Original Covariance Matrix:"
     993             :             do i = 1,rankPDM
     994             :                 write(Test%outputUnit,"(*(F22.15))") (PosDefMat(i,j),j=1,rankPDM)
     995             :             end do
     996             :             write(Test%outputUnit,"(*(g0))")
     997             :             write(Test%outputUnit,"(*(g0))") "RegresCoefMat, RegresCoefMatRef, difference, normalizedDifference:"
     998             :             write(Test%outputUnit,"(*(g0))")
     999             :             write(Test%outputUnit,"(*(g0))") "SchurComplement, SchurComplementRef, difference, normalizedDifference:"
    1000             :             write(Test%outputUnit,"(*(g0))")
    1001             :             do i = 1,rankS11
    1002             :                 do j = 1,rankS22
    1003             :                     write(Test%outputUnit,"(*(F22.15))" ) RegresCoefMat(i,j) &
    1004             :                                                         , RegresCoefMatRef(i,j) &
    1005             :                                                         , abs(RegresCoefMat(i,j)-RegresCoefMatRef(i,j)) &
    1006             :                                                         , normalizedDifference
    1007             :                 end do
    1008             :             end do
    1009             :             write(Test%outputUnit,"(*(g0))")
    1010             :         end if
    1011             :         ! LCOV_EXCL_STOP
    1012             : 
    1013           4 :         do i = 1,rankS11
    1014          13 :             do j = 1,rankS11
    1015           9 :                 normalizedDifference = abs(SchurComplementRef(i,j) - SchurComplement(i,j)) / (SchurComplementRef(i,j) + SchurComplement(i,j))
    1016          12 :                 assertion = assertion .and. normalizedDifference < 1.e-6_RK
    1017             :             end do
    1018             :         end do
    1019             : 
    1020           1 :         if (Test%isDebugMode .and. .not. assertion) then
    1021             :         ! LCOV_EXCL_START
    1022             :             do i = 1,rankS11
    1023             :                 do j = 1,rankS11
    1024             :                     write(Test%outputUnit,"(*(F22.15))" ) SchurComplement(i,j) &
    1025             :                                                         , SchurComplementRef(i,j) &
    1026             :                                                         , abs(SchurComplement(i,j)-SchurComplementRef(i,j)) &
    1027             :                                                         , normalizedDifference
    1028             :                 end do
    1029             :             end do
    1030             :         end if
    1031             :         ! LCOV_EXCL_STOP
    1032             : 
    1033           1 :     end function test_getRegresCoef_1
    1034             : 
    1035             : !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    1036             : 
    1037             : end module Test_Matrix_mod ! LCOV_EXCL_LINE

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