ParaMonte Fortran 2.0.0
Parallel Monte Carlo and Machine Learning Library
See the latest version documentation.
pm_matrixDet::setMatDetSqrt Interface Reference

Return the determinant of the input positive-definite square matrix. More...

Detailed Description

Return the determinant of the input positive-definite square matrix.

The positive-definiteness guarantee by the user enables a fast method of computing the square root of the determinant of the input positive-definite square matrix via its Cholesky Factorization.
For high rank matrices, there is an overflow possibility for the output of this generic interface.
The overflow risk can be avoided by calling setMatDetSqrtLog instead.

Note
This generic interface is meant to facilitate simultaneous computation of the Cholesky factorization and determinant of a input positive definite matrix.
If the Cholesky factorization is already computed, simply compute the determinant as the square of the multiplicative trace of the Cholesky factor.
Parameters
[in]mat: The input contiguous matrix of shape (1:ndim, 1:ndim) of
  1. type complex of kind any supported by the processor (e.g., CK, CK32, CK64, or CK128),
  2. type real of kind any supported by the processor (e.g., RK, RK32, RK64, or RK128),
containing a subset of the positive-definite matrix whose log(sqrt(determinant)) is to be computed.
[in]subset: The input scalar that can be either,
  1. the constant uppDia implying that the upper-diagonal triangular block of mat should be used for building the lower-diagonal Cholesky factor in the output chol.
  2. the constant lowDia implying that the lower-diagonal triangular block of mat should be used for building the upper-diagonal Cholesky factor in the output chol.
This argument is merely a convenience to differentiate the different procedure functionalities within this generic interface.
[out]detSqrt: The output scalar of type real of the same kind as the input mat representing the square root of its determinant.
[out]info: The output scalar integer of default kind IK that is 0 if and only if the Cholesky factorization succeeds.
Otherwise, it is set to the order of the leading minor of the specified input subset of mat that is not positive definite, indicating the occurrence of an error and that the factorization could not be completed.
A non-zero value implies failure of the Cholesky factorization.
A non-zero info implied the input mat is not positive definite.
[in,out]chol: The input/output matrix of the same type, kind, and shape as the input mat containing the Cholesky factorization in its triangular subset dictated by the input arguments subset and operation.
See the documentation of the input argument mat above for more information.
[in]operation: The input scalar that can be either,
  1. the constant nothing implying that the same subset of the output chol as the input subset must be populated with the Cholesky factorization.
  2. the constant transHerm implying that the complex transpose of the computed Cholesky factorization must be written to the output chol.


Possible calling interfaces

call setMatDetSqrtLog(mat, subset, detSqrt, info, chol, operation)
Return the natural logarithm of the square-root of the determinant of the input positive-definite squ...
This module contains procedures and generic interfaces relevant to the computation of the determinant...
Warning
The condition size(mat, 1) == size(mat, 2) must hold for the corresponding input arguments.
This condition is verified only if the library is built with the preprocessor macro CHECK_ENABLED=1.
The input mat is assumed to be positive-definite, otherwise the Cholesky Factorization within the procedure will fail, in which case, the procedure will return the control with info /= 0.
The pure procedure(s) documented herein become impure when the ParaMonte library is compiled with preprocessor macro CHECK_ENABLED=1.
By default, these procedures are pure in release build and impure in debug and testing builds.
See also
getMatDet
setMatDet
getMatDetSqrt
setMatDetSqrt
getMatDetSqrtLog
setMatDetSqrtLog
getMatMulTrace


Example usage

1program example
2
3 use pm_kind, only: SK, IK, LK
4 use pm_matrixCopy, only: rdpack
5 use pm_matrixCopy, only: transHerm
10 use pm_matrixChol, only: uppDia, lowDia
11 use pm_matrixChol, only: getMatChol
12 use pm_distCov, only: getCovRand
13 use pm_arrayResize, only: setResized
14 use pm_distUnif, only: getUnifRand
15 use pm_matrixInv, only: getMatInv
16 use pm_io, only: display_type
17 use pm_io, only: getFormat
18
19 implicit none
20
21 integer(IK) :: info, ndim, itry, ntry = 10
22 character(:, SK), allocatable :: format
23 type(display_type) :: disp
24
25 disp = display_type(file = "main.out.F90")
26
27 call disp%skip()
28 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
29 call disp%show("! Compute the sqrt of the determinant of the positive definite matrix.")
30 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
31 call disp%skip()
32
33 block
34 use pm_kind, only: TKG => RKS
35 real(TKG), allocatable :: mat(:,:)
36 real(TKG) :: detSqrtLog
37 format = getFormat(mold = [0._TKG], ed = SK_"e")
38 do itry = 1, ntry
39 call disp%skip()
40 call disp%show("ndim = getUnifRand(1, 5)")
41 ndim = getUnifRand(1, 5)
42 call disp%show("ndim")
43 call disp%show( ndim )
44 call disp%show("call setResized(mat, [ndim, ndim + 1])")
45 call setResized(mat, [ndim, ndim + 1])
46 call disp%show("mat(:,1:ndim) = getCovRand(mold = 1._TKG, ndim = ndim)")
47 mat(:,1:ndim) = getCovRand(mold = 1._TKG, ndim = ndim)
48 call disp%show("mat(:,1:ndim)")
49 call disp%show( mat(:,1:ndim) , format = format )
50 call disp%show("mat(:,1:ndim) = getMatCopy(rdpack, mat(:,1:ndim), rdpack, lowDia, init = 0._TKG) ! reset the upper.")
51 mat(:,1:ndim) = getMatCopy(rdpack, mat(:,1:ndim), rdpack, lowDia, init = 0._TKG)
52 call disp%show("mat(:,1:ndim)")
53 call disp%show( mat(:,1:ndim) , format = format )
54 call disp%show("call setMatDetSqrtLog(mat(:,1:ndim), lowDia, detSqrtLog, info, mat(:,2:ndim+1), transHerm)")
55 call setMatDetSqrtLog(mat(:,1:ndim), lowDia, detSqrtLog, info, mat(:,2:ndim+1), transHerm)
56 call disp%show("if (info /= 0) error stop 'Cholesky factorization failed.'")
57 if (info /= 0) error stop 'Cholesky factorization failed.'
58 call disp%show("detSqrtLog")
59 call disp%show( detSqrtLog )
60 call disp%show("mat")
61 call disp%show( mat , format = format )
62 call disp%show("getMatChol(mat(:,1:ndim), lowDia)")
63 call disp%show( getMatChol(mat(:,1:ndim), lowDia) , format = format )
64 call disp%show("getMatMulTraceLog(getMatChol(mat(:,1:ndim), lowDia)) ! for comparison.")
65 call disp%show( getMatMulTraceLog(getMatChol(mat(:,1:ndim), lowDia)) )
66 call disp%show("detSqrtLog - getMatMulTraceLog(getMatChol(mat(:,1:ndim), lowDia)) ! must be one.")
67 call disp%show( detSqrtLog - getMatMulTraceLog(getMatChol(mat(:,1:ndim), lowDia)) )
68 call disp%skip()
69 call disp%show("mat(:,2:ndim+1) = getCovRand(mold = 1._TKG, ndim = ndim)")
70 mat(:,2:ndim+1) = getCovRand(mold = 1._TKG, ndim = ndim)
71 call disp%show("mat(:,2:ndim+1) = getMatCopy(rdpack, mat(:,2:ndim+1), rdpack, uppDia, init = 0._TKG) ! reset the lower.")
72 mat(:,2:ndim+1) = getMatCopy(rdpack, mat(:,2:ndim+1), rdpack, uppDia, init = 0._TKG)
73 call disp%show("mat(:,2:ndim+1)")
74 call disp%show( mat(:,2:ndim+1) , format = format )
75 call disp%show("call setMatDetSqrtLog(mat(:,2:ndim+1), uppDia, detSqrtLog, info, mat(:,1:ndim), transHerm)")
76 call setMatDetSqrtLog(mat(:,2:ndim+1), uppDia, detSqrtLog, info, mat(:,1:ndim), transHerm)
77 call disp%show("if (info /= 0) error stop 'Cholesky factorization failed.'")
78 if (info /= 0) error stop 'Cholesky factorization failed.'
79 call disp%show("detSqrtLog")
80 call disp%show( detSqrtLog )
81 call disp%show("mat")
82 call disp%show( mat , format = format )
83 call disp%show("getMatChol(mat(:,2:ndim+1), uppDia)")
84 call disp%show( getMatChol(mat(:,2:ndim+1), uppDia) , format = format )
85 call disp%show("getMatMulTraceLog(getMatChol(mat(:,2:ndim+1), uppDia)) ! for comparison.")
86 call disp%show( getMatMulTraceLog(getMatChol(mat(:,2:ndim+1), uppDia)) )
87 call disp%show("detSqrtLog - getMatMulTraceLog(getMatChol(mat(:,2:ndim+1), uppDia)) ! must be one.")
88 call disp%show( detSqrtLog - getMatMulTraceLog(getMatChol(mat(:,2:ndim+1), uppDia)) )
89 call disp%skip()
90 end do
91 end block
92
93 block
94 use pm_kind, only: TKG => CKS
95 integer(IK), parameter :: ndim = 3
96 complex(TKG), allocatable :: tmp(:,:)
97 complex(TKG), parameter :: mat(*,*) = reshape( [ (9.0, 0.0), (3.0, 3.0), (3.0, -3.0) &
98 , (3.0, -3.0),(18.0, 0.0), (8.0, -6.0) &
99 , (3.0, 3.0), (8.0, 6.0),(43.0, 0.0) &
100 ], shape = [ndim, ndim], order = [2, 1])
101 real(TKG) :: detSqrtLog
102 format = getFormat(mold = [(0._TKG, 0._TKG)], ed = SK_"e")
103 call disp%skip()
104 call disp%show("mat")
105 call disp%show( mat , format = format )
106 call disp%show("call setResized(tmp, [ndim, ndim + 1])")
107 call setResized(tmp, [ndim, ndim + 1])
108 call disp%show("tmp(:,1:ndim) = getMatCopy(rdpack, mat, rdpack, lowDia, init = (0._TKG, 0._TKG)) ! reset the upper.")
109 tmp(:,1:ndim) = getMatCopy(rdpack, mat, rdpack, lowDia, init = (0._TKG, 0._TKG))
110 call disp%show("tmp(:,1:ndim)")
111 call disp%show( tmp(:,1:ndim) , format = format )
112 call disp%show("call setMatDetSqrtLog(tmp(:,1:ndim), lowDia, detSqrtLog, info, tmp(:,2:ndim+1), transHerm)")
113 call setMatDetSqrtLog(tmp(:,1:ndim), lowDia, detSqrtLog, info, tmp(:,2:ndim+1), transHerm)
114 call disp%show("if (info /= 0) error stop 'Cholesky factorization failed.'")
115 if (info /= 0) error stop 'Cholesky factorization failed.'
116 call disp%show("detSqrtLog")
117 call disp%show( detSqrtLog )
118 call disp%show("tmp")
119 call disp%show( tmp , format = format )
120 call disp%show("getMatMulTraceLog(getMatChol(tmp(:,1:ndim), lowDia)) ! for comparison.")
121 call disp%show( getMatMulTraceLog(getMatChol(tmp(:,1:ndim), lowDia)) )
122 call disp%show("detSqrtLog - getMatMulTraceLog(getMatChol(tmp(:,1:ndim), lowDia)) ! must be one")
123 call disp%show( detSqrtLog - getMatMulTraceLog(getMatChol(tmp(:,1:ndim), lowDia)) )
124 call disp%skip()
125 call disp%show("tmp(:,2:ndim+1) = getMatCopy(rdpack, mat(:,1:ndim), rdpack, uppDia, init = (0._TKG, 0._TKG)) ! reset the upper.")
126 tmp(:,2:ndim+1) = getMatCopy(rdpack, mat(:,1:ndim), rdpack, uppDia, init = (0._TKG, 0._TKG))
127 call disp%show("tmp(:,2:ndim+1)")
128 call disp%show( tmp(:,2:ndim+1) , format = format )
129 call disp%show("call setMatDetSqrtLog(tmp(:,2:ndim+1), uppDia, detSqrtLog, info, tmp(:,1:ndim), transHerm)")
130 call setMatDetSqrtLog(tmp(:,2:ndim+1), uppDia, detSqrtLog, info, tmp(:,1:ndim), transHerm)
131 call disp%show("if (info /= 0) error stop 'Cholesky factorization failed.'")
132 if (info /= 0) error stop 'Cholesky factorization failed.'
133 call disp%show("detSqrtLog")
134 call disp%show( detSqrtLog )
135 call disp%show("tmp")
136 call disp%show( tmp , format = format )
137 call disp%show("getMatMulTraceLog(getMatChol(tmp(:,2:ndim+1), uppDia)) ! for comparison.")
138 call disp%show( getMatMulTraceLog(getMatChol(tmp(:,2:ndim+1), uppDia)) )
139 call disp%show("detSqrtLog - getMatMulTraceLog(getMatChol(tmp(:,2:ndim+1), uppDia)) ! must be one.")
140 call disp%show( detSqrtLog - getMatMulTraceLog(getMatChol(tmp(:,2:ndim+1), uppDia)) )
141 call disp%skip()
142 end block
143
144 block
145 use pm_kind, only: TKG => CKS
146 integer(IK), parameter :: ndim = 3
147 complex(TKG), allocatable :: tmp(:,:)
148 complex(TKG), parameter :: mat(*,*) = reshape( [ (25.0, 0.0), (-5.0, -5.0), (10.0, 5.0) &
149 , (-5.0, 5.0), (51.0, 0.0), (4.0, -6.0) &
150 , (10.0, -5.0), (4.0, 6.0), (71.0, 0.0) &
151 ], shape = [ndim, ndim], order = [2, 1])
152 real(TKG) :: detSqrtLog
153 format = getFormat(mold = [(0._TKG, 0._TKG)], ed = SK_"e")
154 call disp%skip()
155 call disp%show("mat")
156 call disp%show( mat , format = format )
157 call disp%show("call setResized(tmp, [ndim, ndim + 1])")
158 call setResized(tmp, [ndim, ndim + 1])
159 call disp%show("tmp(:,1:ndim) = getMatCopy(rdpack, mat, rdpack, lowDia, init = (0._TKG, 0._TKG)) ! reset the upper.")
160 tmp(:,1:ndim) = getMatCopy(rdpack, mat, rdpack, lowDia, init = (0._TKG, 0._TKG))
161 call disp%show("tmp(:,1:ndim)")
162 call disp%show( tmp(:,1:ndim) , format = format )
163 call disp%show("call setMatDetSqrtLog(tmp(:,1:ndim), lowDia, detSqrtLog, info, tmp(:,2:ndim+1), transHerm)")
164 call setMatDetSqrtLog(tmp(:,1:ndim), lowDia, detSqrtLog, info, tmp(:,2:ndim+1), transHerm)
165 call disp%show("if (info /= 0) error stop 'Cholesky factorization failed.'")
166 if (info /= 0) error stop 'Cholesky factorization failed.'
167 call disp%show("detSqrtLog")
168 call disp%show( detSqrtLog )
169 call disp%show("tmp")
170 call disp%show( tmp , format = format )
171 call disp%show("getMatMulTraceLog(getMatChol(tmp(:,1:ndim), lowDia)) ! for comparison.")
172 call disp%show( getMatMulTraceLog(getMatChol(tmp(:,1:ndim), lowDia)) )
173 call disp%show("detSqrtLog - getMatMulTraceLog(getMatChol(tmp(:,1:ndim), lowDia)) ! must be one")
174 call disp%show( detSqrtLog - getMatMulTraceLog(getMatChol(tmp(:,1:ndim), lowDia)) )
175 call disp%skip()
176 call disp%show("tmp(:,2:ndim+1) = getMatCopy(rdpack, mat(:,1:ndim), rdpack, uppDia, init = (0._TKG, 0._TKG)) ! reset the upper.")
177 tmp(:,2:ndim+1) = getMatCopy(rdpack, mat(:,1:ndim), rdpack, uppDia, init = (0._TKG, 0._TKG))
178 call disp%show("tmp(:,2:ndim+1)")
179 call disp%show( tmp(:,2:ndim+1) , format = format )
180 call disp%show("call setMatDetSqrtLog(tmp(:,2:ndim+1), uppDia, detSqrtLog, info, tmp(:,1:ndim), transHerm)")
181 call setMatDetSqrtLog(tmp(:,2:ndim+1), uppDia, detSqrtLog, info, tmp(:,1:ndim), transHerm)
182 call disp%show("if (info /= 0) error stop 'Cholesky factorization failed.'")
183 if (info /= 0) error stop 'Cholesky factorization failed.'
184 call disp%show("detSqrtLog")
185 call disp%show( detSqrtLog )
186 call disp%show("tmp")
187 call disp%show( tmp , format = format )
188 call disp%show("getMatMulTraceLog(getMatChol(tmp(:,2:ndim+1), uppDia)) ! for comparison.")
189 call disp%show( getMatMulTraceLog(getMatChol(tmp(:,2:ndim+1), uppDia)) )
190 call disp%show("detSqrtLog - getMatMulTraceLog(getMatChol(tmp(:,2:ndim+1), uppDia)) ! must be one.")
191 call disp%show( detSqrtLog - getMatMulTraceLog(getMatChol(tmp(:,2:ndim+1), uppDia)) )
192 call disp%skip()
193 end block
194
195end program example
Allocate or resize (shrink or expand) an input allocatable scalar string or array of rank 1....
Generate and return a random positive-definite (correlation or covariance) matrix using the Gram meth...
Definition: pm_distCov.F90:394
Generate and return a scalar or a contiguous array of rank 1 of length s1 of randomly uniformly distr...
Generate and return a generic or type/kind-specific IO format with the requested specifications that ...
Definition: pm_io.F90:18485
This is a generic method of the derived type display_type with pass attribute.
Definition: pm_io.F90:11726
This is a generic method of the derived type display_type with pass attribute.
Definition: pm_io.F90:11508
Generate and return the upper or the lower Cholesky factorization of the input symmetric positive-def...
Generate and return a copy of a desired subset of the input source matrix of arbitrary shape (:) or (...
Generate and return the natural logarithm of the square-root of the determinant of the input positive...
Generate and return the full inverse of an input matrix of general or triangular form directly or thr...
Generate and return the natural logarithm of the multiplicative trace of an input square matrix of ty...
This module contains procedures and generic interfaces for resizing allocatable arrays of various typ...
This module contains classes and procedures for generating random matrices distributed on the space o...
Definition: pm_distCov.F90:72
This module contains classes and procedures for computing various statistical quantities related to t...
This module contains classes and procedures for input/output (IO) or generic display operations on st...
Definition: pm_io.F90:252
type(display_type) disp
This is a scalar module variable an object of type display_type for general display.
Definition: pm_io.F90:11393
This module defines the relevant Fortran kind type-parameters frequently used in the ParaMonte librar...
Definition: pm_kind.F90:268
integer, parameter LK
The default logical kind in the ParaMonte library: kind(.true.) in Fortran, kind(....
Definition: pm_kind.F90:541
integer, parameter CKS
The single-precision complex kind in Fortran mode. On most platforms, this is a 32-bit real kind.
Definition: pm_kind.F90:570
integer, parameter IK
The default integer kind in the ParaMonte library: int32 in Fortran, c_int32_t in C-Fortran Interoper...
Definition: pm_kind.F90:540
integer, parameter SK
The default character kind in the ParaMonte library: kind("a") in Fortran, c_char in C-Fortran Intero...
Definition: pm_kind.F90:539
integer, parameter RKS
The single-precision real kind in Fortran mode. On most platforms, this is an 32-bit real kind.
Definition: pm_kind.F90:567
This module contains procedures and generic interfaces for computing the Cholesky factorization of po...
This module contains procedures and generic interfaces relevant to copying (diagonal or upper/lower t...
This module contains abstract and concrete derived types and procedures related to the inversion of s...
This module contains procedures and generic interfaces for computing the additive or multiplicative t...
Generate and return an object of type display_type.
Definition: pm_io.F90:10282

Example Unix compile command via Intel ifort compiler
1#!/usr/bin/env sh
2rm main.exe
3ifort -fpp -standard-semantics -O3 -Wl,-rpath,../../../lib -I../../../inc main.F90 ../../../lib/libparamonte* -o main.exe
4./main.exe

Example Windows Batch compile command via Intel ifort compiler
1del main.exe
2set PATH=..\..\..\lib;%PATH%
3ifort /fpp /standard-semantics /O3 /I:..\..\..\include main.F90 ..\..\..\lib\libparamonte*.lib /exe:main.exe
4main.exe

Example Unix / MinGW compile command via GNU gfortran compiler
1#!/usr/bin/env sh
2rm main.exe
3gfortran -cpp -ffree-line-length-none -O3 -Wl,-rpath,../../../lib -I../../../inc main.F90 ../../../lib/libparamonte* -o main.exe
4./main.exe

Example output
1
2!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
3! Compute the sqrt of the determinant of the positive definite matrix.
4!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
5
6
7ndim = getUnifRand(1, 5)
8ndim
9+5
10call setResized(mat, [ndim, ndim + 1])
11mat(:,1:ndim) = getCovRand(mold = 1._TKG, ndim = ndim)
12mat(:,1:ndim)
13 0.100000E+01, -0.981024E+00, -0.807195E+00, -0.702849E+00, -0.282062E+00
14-0.981024E+00, 0.100000E+01, 0.765316E+00, 0.657443E+00, 0.398481E+00
15-0.807195E+00, 0.765316E+00, 0.100000E+01, 0.294366E+00, 0.402081E+00
16-0.702849E+00, 0.657443E+00, 0.294366E+00, 0.100000E+01, -0.280745E+00
17-0.282062E+00, 0.398481E+00, 0.402081E+00, -0.280745E+00, 0.100000E+01
18mat(:,1:ndim) = getMatCopy(rdpack, mat(:,1:ndim), rdpack, lowDia, init = 0._TKG) ! reset the upper.
19mat(:,1:ndim)
20 0.100000E+01, 0.000000E+00, 0.000000E+00, 0.000000E+00, 0.000000E+00
21-0.981024E+00, 0.100000E+01, 0.000000E+00, 0.000000E+00, 0.000000E+00
22-0.807195E+00, 0.765316E+00, 0.100000E+01, 0.000000E+00, 0.000000E+00
23-0.702849E+00, 0.657443E+00, 0.294366E+00, 0.100000E+01, 0.000000E+00
24-0.282062E+00, 0.398481E+00, 0.402081E+00, -0.280745E+00, 0.100000E+01
25call setMatDetSqrtLog(mat(:,1:ndim), lowDia, detSqrtLog, info, mat(:,2:ndim+1), transHerm)
26if (info /= 0) error stop 'Cholesky factorization failed.'
27detSqrtLog
28-3.71000409
29mat
30 0.100000E+01, 0.100000E+01, -0.981024E+00, -0.807195E+00, -0.702849E+00, -0.282062E+00
31-0.981024E+00, 0.100000E+01, 0.193888E+00, -0.136995E+00, -0.165394E+00, 0.628048E+00
32-0.807195E+00, 0.765316E+00, 0.100000E+01, 0.574167E+00, -0.514883E+00, 0.453598E+00
33-0.702849E+00, 0.657443E+00, 0.294366E+00, 0.100000E+01, 0.462108E+00, -0.306351E+00
34-0.282062E+00, 0.398481E+00, 0.402081E+00, -0.280745E+00, 0.100000E+01, 0.475809E+00
35getMatChol(mat(:,1:ndim), lowDia)
36 0.100000E+01, 0.000000E+00, 0.000000E+00, 0.000000E+00, 0.000000E+00
37-0.981024E+00, 0.193888E+00, 0.000000E+00, 0.000000E+00, 0.000000E+00
38-0.807195E+00, -0.136995E+00, 0.574167E+00, 0.000000E+00, 0.000000E+00
39-0.702849E+00, -0.165394E+00, -0.514883E+00, 0.462108E+00, 0.000000E+00
40-0.282062E+00, 0.628048E+00, 0.453598E+00, -0.306351E+00, 0.475809E+00
41getMatMulTraceLog(getMatChol(mat(:,1:ndim), lowDia)) ! for comparison.
42-3.71000409
43detSqrtLog - getMatMulTraceLog(getMatChol(mat(:,1:ndim), lowDia)) ! must be one.
44+0.00000000
45
46mat(:,2:ndim+1) = getCovRand(mold = 1._TKG, ndim = ndim)
47mat(:,2:ndim+1) = getMatCopy(rdpack, mat(:,2:ndim+1), rdpack, uppDia, init = 0._TKG) ! reset the lower.
48mat(:,2:ndim+1)
49 0.100000E+01, -0.173145E+00, -0.155464E+00, -0.278116E+00, 0.445956E+00
50 0.000000E+00, 0.100000E+01, -0.819725E+00, 0.750250E+00, -0.548509E+00
51 0.000000E+00, 0.000000E+00, 0.100000E+01, -0.848054E+00, 0.678673E-01
52 0.000000E+00, 0.000000E+00, 0.000000E+00, 0.100000E+01, -0.845793E-01
53 0.000000E+00, 0.000000E+00, 0.000000E+00, 0.000000E+00, 0.100000E+01
54call setMatDetSqrtLog(mat(:,2:ndim+1), uppDia, detSqrtLog, info, mat(:,1:ndim), transHerm)
55if (info /= 0) error stop 'Cholesky factorization failed.'
56detSqrtLog
57-2.72252059
58mat
59 0.100000E+01, 0.100000E+01, -0.173145E+00, -0.155464E+00, -0.278116E+00, 0.445956E+00
60-0.173145E+00, 0.984896E+00, 0.100000E+01, -0.819725E+00, 0.750250E+00, -0.548509E+00
61-0.155464E+00, -0.859626E+00, 0.486697E+00, 0.100000E+01, -0.848054E+00, 0.678673E-01
62-0.278116E+00, 0.712863E+00, -0.572216E+00, 0.295038E+00, 0.100000E+01, -0.845793E-01
63 0.445956E+00, -0.478522E+00, -0.563292E+00, 0.197410E+00, 0.464620E+00, 0.100000E+01
64getMatChol(mat(:,2:ndim+1), uppDia)
65 0.100000E+01, -0.173145E+00, -0.155464E+00, -0.278116E+00, 0.445956E+00
66 0.000000E+00, 0.984896E+00, -0.859626E+00, 0.712863E+00, -0.478522E+00
67 0.000000E+00, 0.000000E+00, 0.486697E+00, -0.572216E+00, -0.563292E+00
68 0.000000E+00, 0.000000E+00, 0.000000E+00, 0.295038E+00, 0.197410E+00
69 0.000000E+00, 0.000000E+00, 0.000000E+00, 0.000000E+00, 0.464620E+00
70getMatMulTraceLog(getMatChol(mat(:,2:ndim+1), uppDia)) ! for comparison.
71-2.72252059
72detSqrtLog - getMatMulTraceLog(getMatChol(mat(:,2:ndim+1), uppDia)) ! must be one.
73+0.00000000
74
75
76ndim = getUnifRand(1, 5)
77ndim
78+1
79call setResized(mat, [ndim, ndim + 1])
80mat(:,1:ndim) = getCovRand(mold = 1._TKG, ndim = ndim)
81mat(:,1:ndim)
82 0.100000E+01
83mat(:,1:ndim) = getMatCopy(rdpack, mat(:,1:ndim), rdpack, lowDia, init = 0._TKG) ! reset the upper.
84mat(:,1:ndim)
85 0.100000E+01
86call setMatDetSqrtLog(mat(:,1:ndim), lowDia, detSqrtLog, info, mat(:,2:ndim+1), transHerm)
87if (info /= 0) error stop 'Cholesky factorization failed.'
88detSqrtLog
89+0.00000000
90mat
91 0.100000E+01, 0.100000E+01
92getMatChol(mat(:,1:ndim), lowDia)
93 0.100000E+01
94getMatMulTraceLog(getMatChol(mat(:,1:ndim), lowDia)) ! for comparison.
95+0.00000000
96detSqrtLog - getMatMulTraceLog(getMatChol(mat(:,1:ndim), lowDia)) ! must be one.
97+0.00000000
98
99mat(:,2:ndim+1) = getCovRand(mold = 1._TKG, ndim = ndim)
100mat(:,2:ndim+1) = getMatCopy(rdpack, mat(:,2:ndim+1), rdpack, uppDia, init = 0._TKG) ! reset the lower.
101mat(:,2:ndim+1)
102 0.100000E+01
103call setMatDetSqrtLog(mat(:,2:ndim+1), uppDia, detSqrtLog, info, mat(:,1:ndim), transHerm)
104if (info /= 0) error stop 'Cholesky factorization failed.'
105detSqrtLog
106+0.00000000
107mat
108 0.100000E+01, 0.100000E+01
109getMatChol(mat(:,2:ndim+1), uppDia)
110 0.100000E+01
111getMatMulTraceLog(getMatChol(mat(:,2:ndim+1), uppDia)) ! for comparison.
112+0.00000000
113detSqrtLog - getMatMulTraceLog(getMatChol(mat(:,2:ndim+1), uppDia)) ! must be one.
114+0.00000000
115
116
117ndim = getUnifRand(1, 5)
118ndim
119+2
120call setResized(mat, [ndim, ndim + 1])
121mat(:,1:ndim) = getCovRand(mold = 1._TKG, ndim = ndim)
122mat(:,1:ndim)
123 0.100000E+01, 0.803510E+00
124 0.803510E+00, 0.100000E+01
125mat(:,1:ndim) = getMatCopy(rdpack, mat(:,1:ndim), rdpack, lowDia, init = 0._TKG) ! reset the upper.
126mat(:,1:ndim)
127 0.100000E+01, 0.000000E+00
128 0.803510E+00, 0.100000E+01
129call setMatDetSqrtLog(mat(:,1:ndim), lowDia, detSqrtLog, info, mat(:,2:ndim+1), transHerm)
130if (info /= 0) error stop 'Cholesky factorization failed.'
131detSqrtLog
132-0.518703997
133mat
134 0.100000E+01, 0.100000E+01, 0.803510E+00
135 0.803510E+00, 0.100000E+01, 0.595292E+00
136getMatChol(mat(:,1:ndim), lowDia)
137 0.100000E+01, 0.000000E+00
138 0.803510E+00, 0.595292E+00
139getMatMulTraceLog(getMatChol(mat(:,1:ndim), lowDia)) ! for comparison.
140-0.518703997
141detSqrtLog - getMatMulTraceLog(getMatChol(mat(:,1:ndim), lowDia)) ! must be one.
142+0.00000000
143
144mat(:,2:ndim+1) = getCovRand(mold = 1._TKG, ndim = ndim)
145mat(:,2:ndim+1) = getMatCopy(rdpack, mat(:,2:ndim+1), rdpack, uppDia, init = 0._TKG) ! reset the lower.
146mat(:,2:ndim+1)
147 0.100000E+01, 0.563573E+00
148 0.000000E+00, 0.100000E+01
149call setMatDetSqrtLog(mat(:,2:ndim+1), uppDia, detSqrtLog, info, mat(:,1:ndim), transHerm)
150if (info /= 0) error stop 'Cholesky factorization failed.'
151detSqrtLog
152-0.191079870
153mat
154 0.100000E+01, 0.100000E+01, 0.563573E+00
155 0.563573E+00, 0.826067E+00, 0.100000E+01
156getMatChol(mat(:,2:ndim+1), uppDia)
157 0.100000E+01, 0.563573E+00
158 0.000000E+00, 0.826067E+00
159getMatMulTraceLog(getMatChol(mat(:,2:ndim+1), uppDia)) ! for comparison.
160-0.191079870
161detSqrtLog - getMatMulTraceLog(getMatChol(mat(:,2:ndim+1), uppDia)) ! must be one.
162+0.00000000
163
164
165ndim = getUnifRand(1, 5)
166ndim
167+4
168call setResized(mat, [ndim, ndim + 1])
169mat(:,1:ndim) = getCovRand(mold = 1._TKG, ndim = ndim)
170mat(:,1:ndim)
171 0.100000E+01, -0.989069E+00, 0.632195E+00, 0.501644E+00
172-0.989069E+00, 0.100000E+01, -0.676797E+00, -0.586130E+00
173 0.632195E+00, -0.676797E+00, 0.100000E+01, 0.281726E+00
174 0.501644E+00, -0.586130E+00, 0.281726E+00, 0.100000E+01
175mat(:,1:ndim) = getMatCopy(rdpack, mat(:,1:ndim), rdpack, lowDia, init = 0._TKG) ! reset the upper.
176mat(:,1:ndim)
177 0.100000E+01, 0.000000E+00, 0.000000E+00, 0.000000E+00
178-0.989069E+00, 0.100000E+01, 0.000000E+00, 0.000000E+00
179 0.632195E+00, -0.676797E+00, 0.100000E+01, 0.000000E+00
180 0.501644E+00, -0.586130E+00, 0.281726E+00, 0.100000E+01
181call setMatDetSqrtLog(mat(:,1:ndim), lowDia, detSqrtLog, info, mat(:,2:ndim+1), transHerm)
182if (info /= 0) error stop 'Cholesky factorization failed.'
183detSqrtLog
184-2.98239970
185mat
186 0.100000E+01, 0.100000E+01, -0.989069E+00, 0.632195E+00, 0.501644E+00
187-0.989069E+00, 0.100000E+01, 0.147455E+00, -0.349342E+00, -0.610152E+00
188 0.632195E+00, -0.676797E+00, 0.100000E+01, 0.691585E+00, -0.359409E+00
189 0.501644E+00, -0.586130E+00, 0.281726E+00, 0.100000E+01, 0.496884E+00
190getMatChol(mat(:,1:ndim), lowDia)
191 0.100000E+01, 0.000000E+00, 0.000000E+00, 0.000000E+00
192-0.989069E+00, 0.147455E+00, 0.000000E+00, 0.000000E+00
193 0.632195E+00, -0.349342E+00, 0.691585E+00, 0.000000E+00
194 0.501644E+00, -0.610152E+00, -0.359409E+00, 0.496884E+00
195getMatMulTraceLog(getMatChol(mat(:,1:ndim), lowDia)) ! for comparison.
196-2.98239970
197detSqrtLog - getMatMulTraceLog(getMatChol(mat(:,1:ndim), lowDia)) ! must be one.
198+0.00000000
199
200mat(:,2:ndim+1) = getCovRand(mold = 1._TKG, ndim = ndim)
201mat(:,2:ndim+1) = getMatCopy(rdpack, mat(:,2:ndim+1), rdpack, uppDia, init = 0._TKG) ! reset the lower.
202mat(:,2:ndim+1)
203 0.100000E+01, -0.860287E+00, 0.724564E+00, -0.345541E+00
204 0.000000E+00, 0.100000E+01, -0.295646E+00, -0.353903E-02
205 0.000000E+00, 0.000000E+00, 0.100000E+01, -0.538205E+00
206 0.000000E+00, 0.000000E+00, 0.000000E+00, 0.100000E+01
207call setMatDetSqrtLog(mat(:,2:ndim+1), uppDia, detSqrtLog, info, mat(:,1:ndim), transHerm)
208if (info /= 0) error stop 'Cholesky factorization failed.'
209detSqrtLog
210-2.52654362
211mat
212 0.100000E+01, 0.100000E+01, -0.860287E+00, 0.724564E+00, -0.345541E+00
213-0.860287E+00, 0.509811E+00, 0.100000E+01, -0.295646E+00, -0.353903E-02
214 0.724564E+00, 0.642762E+00, 0.248723E+00, 0.100000E+01, -0.538205E+00
215-0.345541E+00, -0.590030E+00, 0.367522E+00, 0.630392E+00, 0.100000E+01
216getMatChol(mat(:,2:ndim+1), uppDia)
217 0.100000E+01, -0.860287E+00, 0.724564E+00, -0.345541E+00
218 0.000000E+00, 0.509811E+00, 0.642762E+00, -0.590030E+00
219 0.000000E+00, 0.000000E+00, 0.248723E+00, 0.367522E+00
220 0.000000E+00, 0.000000E+00, 0.000000E+00, 0.630392E+00
221getMatMulTraceLog(getMatChol(mat(:,2:ndim+1), uppDia)) ! for comparison.
222-2.52654362
223detSqrtLog - getMatMulTraceLog(getMatChol(mat(:,2:ndim+1), uppDia)) ! must be one.
224+0.00000000
225
226
227ndim = getUnifRand(1, 5)
228ndim
229+2
230call setResized(mat, [ndim, ndim + 1])
231mat(:,1:ndim) = getCovRand(mold = 1._TKG, ndim = ndim)
232mat(:,1:ndim)
233 0.100000E+01, -0.597432E+00
234-0.597432E+00, 0.100000E+01
235mat(:,1:ndim) = getMatCopy(rdpack, mat(:,1:ndim), rdpack, lowDia, init = 0._TKG) ! reset the upper.
236mat(:,1:ndim)
237 0.100000E+01, 0.000000E+00
238-0.597432E+00, 0.100000E+01
239call setMatDetSqrtLog(mat(:,1:ndim), lowDia, detSqrtLog, info, mat(:,2:ndim+1), transHerm)
240if (info /= 0) error stop 'Cholesky factorization failed.'
241detSqrtLog
242-0.220746502
243mat
244 0.100000E+01, 0.100000E+01, -0.597432E+00
245-0.597432E+00, 0.100000E+01, 0.801920E+00
246getMatChol(mat(:,1:ndim), lowDia)
247 0.100000E+01, 0.000000E+00
248-0.597432E+00, 0.801920E+00
249getMatMulTraceLog(getMatChol(mat(:,1:ndim), lowDia)) ! for comparison.
250-0.220746502
251detSqrtLog - getMatMulTraceLog(getMatChol(mat(:,1:ndim), lowDia)) ! must be one.
252+0.00000000
253
254mat(:,2:ndim+1) = getCovRand(mold = 1._TKG, ndim = ndim)
255mat(:,2:ndim+1) = getMatCopy(rdpack, mat(:,2:ndim+1), rdpack, uppDia, init = 0._TKG) ! reset the lower.
256mat(:,2:ndim+1)
257 0.100000E+01, -0.275754E+00
258 0.000000E+00, 0.100000E+01
259call setMatDetSqrtLog(mat(:,2:ndim+1), uppDia, detSqrtLog, info, mat(:,1:ndim), transHerm)
260if (info /= 0) error stop 'Cholesky factorization failed.'
261detSqrtLog
262-0.395432599E-1
263mat
264 0.100000E+01, 0.100000E+01, -0.275754E+00
265-0.275754E+00, 0.961228E+00, 0.100000E+01
266getMatChol(mat(:,2:ndim+1), uppDia)
267 0.100000E+01, -0.275754E+00
268 0.000000E+00, 0.961228E+00
269getMatMulTraceLog(getMatChol(mat(:,2:ndim+1), uppDia)) ! for comparison.
270-0.395432599E-1
271detSqrtLog - getMatMulTraceLog(getMatChol(mat(:,2:ndim+1), uppDia)) ! must be one.
272+0.00000000
273
274
275ndim = getUnifRand(1, 5)
276ndim
277+4
278call setResized(mat, [ndim, ndim + 1])
279mat(:,1:ndim) = getCovRand(mold = 1._TKG, ndim = ndim)
280mat(:,1:ndim)
281 0.100000E+01, -0.989510E+00, 0.386073E+00, 0.819919E+00
282-0.989510E+00, 0.100000E+01, -0.490753E+00, -0.847928E+00
283 0.386073E+00, -0.490753E+00, 0.100000E+01, 0.666666E+00
284 0.819919E+00, -0.847928E+00, 0.666666E+00, 0.100000E+01
285mat(:,1:ndim) = getMatCopy(rdpack, mat(:,1:ndim), rdpack, lowDia, init = 0._TKG) ! reset the upper.
286mat(:,1:ndim)
287 0.100000E+01, 0.000000E+00, 0.000000E+00, 0.000000E+00
288-0.989510E+00, 0.100000E+01, 0.000000E+00, 0.000000E+00
289 0.386073E+00, -0.490753E+00, 0.100000E+01, 0.000000E+00
290 0.819919E+00, -0.847928E+00, 0.666666E+00, 0.100000E+01
291call setMatDetSqrtLog(mat(:,1:ndim), lowDia, detSqrtLog, info, mat(:,2:ndim+1), transHerm)
292if (info /= 0) error stop 'Cholesky factorization failed.'
293detSqrtLog
294-3.43692827
295mat
296 0.100000E+01, 0.100000E+01, -0.989510E+00, 0.386073E+00, 0.819919E+00
297-0.989510E+00, 0.100000E+01, 0.144468E+00, -0.752624E+00, -0.253416E+00
298 0.386073E+00, -0.490753E+00, 0.100000E+01, 0.533389E+00, 0.298825E+00
299 0.819919E+00, -0.847928E+00, 0.666666E+00, 0.100000E+01, 0.417394E+00
300getMatChol(mat(:,1:ndim), lowDia)
301 0.100000E+01, 0.000000E+00, 0.000000E+00, 0.000000E+00
302-0.989510E+00, 0.144468E+00, 0.000000E+00, 0.000000E+00
303 0.386073E+00, -0.752624E+00, 0.533389E+00, 0.000000E+00
304 0.819919E+00, -0.253416E+00, 0.298825E+00, 0.417394E+00
305getMatMulTraceLog(getMatChol(mat(:,1:ndim), lowDia)) ! for comparison.
306-3.43692827
307detSqrtLog - getMatMulTraceLog(getMatChol(mat(:,1:ndim), lowDia)) ! must be one.
308+0.00000000
309
310mat(:,2:ndim+1) = getCovRand(mold = 1._TKG, ndim = ndim)
311mat(:,2:ndim+1) = getMatCopy(rdpack, mat(:,2:ndim+1), rdpack, uppDia, init = 0._TKG) ! reset the lower.
312mat(:,2:ndim+1)
313 0.100000E+01, -0.947685E+00, 0.503152E+00, 0.739516E+00
314 0.000000E+00, 0.100000E+01, -0.715148E+00, -0.533841E+00
315 0.000000E+00, 0.000000E+00, 0.100000E+01, 0.469112E-01
316 0.000000E+00, 0.000000E+00, 0.000000E+00, 0.100000E+01
317call setMatDetSqrtLog(mat(:,2:ndim+1), uppDia, detSqrtLog, info, mat(:,1:ndim), transHerm)
318if (info /= 0) error stop 'Cholesky factorization failed.'
319detSqrtLog
320-2.89995837
321mat
322 0.100000E+01, 0.100000E+01, -0.947685E+00, 0.503152E+00, 0.739516E+00
323-0.947685E+00, 0.319206E+00, 0.100000E+01, -0.715148E+00, -0.533841E+00
324 0.503152E+00, -0.746598E+00, 0.435235E+00, 0.100000E+01, 0.469112E-01
325 0.739516E+00, 0.523136E+00, 0.150250E+00, 0.396067E+00, 0.100000E+01
326getMatChol(mat(:,2:ndim+1), uppDia)
327 0.100000E+01, -0.947685E+00, 0.503152E+00, 0.739516E+00
328 0.000000E+00, 0.319206E+00, -0.746598E+00, 0.523136E+00
329 0.000000E+00, 0.000000E+00, 0.435235E+00, 0.150250E+00
330 0.000000E+00, 0.000000E+00, 0.000000E+00, 0.396067E+00
331getMatMulTraceLog(getMatChol(mat(:,2:ndim+1), uppDia)) ! for comparison.
332-2.89995837
333detSqrtLog - getMatMulTraceLog(getMatChol(mat(:,2:ndim+1), uppDia)) ! must be one.
334+0.00000000
335
336
337ndim = getUnifRand(1, 5)
338ndim
339+5
340call setResized(mat, [ndim, ndim + 1])
341mat(:,1:ndim) = getCovRand(mold = 1._TKG, ndim = ndim)
342mat(:,1:ndim)
343 0.100000E+01, 0.736853E+00, 0.918539E+00, -0.263893E-01, 0.770495E+00
344 0.736853E+00, 0.100000E+01, 0.891344E+00, -0.358091E+00, 0.302522E+00
345 0.918539E+00, 0.891344E+00, 0.100000E+01, -0.323577E+00, 0.532899E+00
346-0.263893E-01, -0.358091E+00, -0.323577E+00, 0.100000E+01, 0.415947E+00
347 0.770495E+00, 0.302522E+00, 0.532899E+00, 0.415947E+00, 0.100000E+01
348mat(:,1:ndim) = getMatCopy(rdpack, mat(:,1:ndim), rdpack, lowDia, init = 0._TKG) ! reset the upper.
349mat(:,1:ndim)
350 0.100000E+01, 0.000000E+00, 0.000000E+00, 0.000000E+00, 0.000000E+00
351 0.736853E+00, 0.100000E+01, 0.000000E+00, 0.000000E+00, 0.000000E+00
352 0.918539E+00, 0.891344E+00, 0.100000E+01, 0.000000E+00, 0.000000E+00
353-0.263893E-01, -0.358091E+00, -0.323577E+00, 0.100000E+01, 0.000000E+00
354 0.770495E+00, 0.302522E+00, 0.532899E+00, 0.415947E+00, 0.100000E+01
355call setMatDetSqrtLog(mat(:,1:ndim), lowDia, detSqrtLog, info, mat(:,2:ndim+1), transHerm)
356if (info /= 0) error stop 'Cholesky factorization failed.'
357detSqrtLog
358-3.17460322
359mat
360 0.100000E+01, 0.100000E+01, 0.736853E+00, 0.918539E+00, -0.263893E-01, 0.770495E+00
361 0.736853E+00, 0.100000E+01, 0.676053E+00, 0.317307E+00, -0.500916E+00, -0.392306E+00
362 0.918539E+00, 0.891344E+00, 0.100000E+01, 0.235803E+00, -0.595383E+00, -0.213521E+00
363-0.263893E-01, -0.358091E+00, -0.323577E+00, 0.100000E+01, 0.627619E+00, 0.179473E+00
364 0.770495E+00, 0.302522E+00, 0.532899E+00, 0.415947E+00, 0.100000E+01, 0.417889E+00
365getMatChol(mat(:,1:ndim), lowDia)
366 0.100000E+01, 0.000000E+00, 0.000000E+00, 0.000000E+00, 0.000000E+00
367 0.736853E+00, 0.676053E+00, 0.000000E+00, 0.000000E+00, 0.000000E+00
368 0.918539E+00, 0.317307E+00, 0.235803E+00, 0.000000E+00, 0.000000E+00
369-0.263893E-01, -0.500916E+00, -0.595383E+00, 0.627619E+00, 0.000000E+00
370 0.770495E+00, -0.392306E+00, -0.213521E+00, 0.179473E+00, 0.417889E+00
371getMatMulTraceLog(getMatChol(mat(:,1:ndim), lowDia)) ! for comparison.
372-3.17460322
373detSqrtLog - getMatMulTraceLog(getMatChol(mat(:,1:ndim), lowDia)) ! must be one.
374+0.00000000
375
376mat(:,2:ndim+1) = getCovRand(mold = 1._TKG, ndim = ndim)
377mat(:,2:ndim+1) = getMatCopy(rdpack, mat(:,2:ndim+1), rdpack, uppDia, init = 0._TKG) ! reset the lower.
378mat(:,2:ndim+1)
379 0.100000E+01, 0.930178E+00, -0.225584E+00, -0.490474E+00, -0.372734E+00
380 0.000000E+00, 0.100000E+01, 0.494317E-01, -0.768959E+00, -0.475243E+00
381 0.000000E+00, 0.000000E+00, 0.100000E+01, -0.600441E+00, -0.496376E+00
382 0.000000E+00, 0.000000E+00, 0.000000E+00, 0.100000E+01, 0.519754E+00
383 0.000000E+00, 0.000000E+00, 0.000000E+00, 0.000000E+00, 0.100000E+01
384call setMatDetSqrtLog(mat(:,2:ndim+1), uppDia, detSqrtLog, info, mat(:,1:ndim), transHerm)
385if (info /= 0) error stop 'Cholesky factorization failed.'
386detSqrtLog
387-4.58287621
388mat
389 0.100000E+01, 0.100000E+01, 0.930178E+00, -0.225584E+00, -0.490474E+00, -0.372734E+00
390 0.930178E+00, 0.367109E+00, 0.100000E+01, 0.494317E-01, -0.768959E+00, -0.475243E+00
391-0.225584E+00, 0.706234E+00, 0.671078E+00, 0.100000E+01, -0.600441E+00, -0.496376E+00
392-0.490474E+00, -0.851874E+00, -0.163113E+00, 0.845012E-01, 0.100000E+01, 0.519754E+00
393-0.372734E+00, -0.350126E+00, -0.496496E+00, -0.500703E+00, 0.491192E+00, 0.100000E+01
394getMatChol(mat(:,2:ndim+1), uppDia)
395 0.100000E+01, 0.930178E+00, -0.225584E+00, -0.490474E+00, -0.372734E+00
396 0.000000E+00, 0.367109E+00, 0.706234E+00, -0.851874E+00, -0.350126E+00
397 0.000000E+00, 0.000000E+00, 0.671078E+00, -0.163113E+00, -0.496496E+00
398 0.000000E+00, 0.000000E+00, 0.000000E+00, 0.845012E-01, -0.500703E+00
399 0.000000E+00, 0.000000E+00, 0.000000E+00, 0.000000E+00, 0.491192E+00
400getMatMulTraceLog(getMatChol(mat(:,2:ndim+1), uppDia)) ! for comparison.
401-4.58287621
402detSqrtLog - getMatMulTraceLog(getMatChol(mat(:,2:ndim+1), uppDia)) ! must be one.
403+0.00000000
404
405
406ndim = getUnifRand(1, 5)
407ndim
408+2
409call setResized(mat, [ndim, ndim + 1])
410mat(:,1:ndim) = getCovRand(mold = 1._TKG, ndim = ndim)
411mat(:,1:ndim)
412 0.100000E+01, 0.926814E+00
413 0.926814E+00, 0.100000E+01
414mat(:,1:ndim) = getMatCopy(rdpack, mat(:,1:ndim), rdpack, lowDia, init = 0._TKG) ! reset the upper.
415mat(:,1:ndim)
416 0.100000E+01, 0.000000E+00
417 0.926814E+00, 0.100000E+01
418call setMatDetSqrtLog(mat(:,1:ndim), lowDia, detSqrtLog, info, mat(:,2:ndim+1), transHerm)
419if (info /= 0) error stop 'Cholesky factorization failed.'
420detSqrtLog
421-0.979442835
422mat
423 0.100000E+01, 0.100000E+01, 0.926814E+00
424 0.926814E+00, 0.100000E+01, 0.375520E+00
425getMatChol(mat(:,1:ndim), lowDia)
426 0.100000E+01, 0.000000E+00
427 0.926814E+00, 0.375520E+00
428getMatMulTraceLog(getMatChol(mat(:,1:ndim), lowDia)) ! for comparison.
429-0.979442835
430detSqrtLog - getMatMulTraceLog(getMatChol(mat(:,1:ndim), lowDia)) ! must be one.
431+0.00000000
432
433mat(:,2:ndim+1) = getCovRand(mold = 1._TKG, ndim = ndim)
434mat(:,2:ndim+1) = getMatCopy(rdpack, mat(:,2:ndim+1), rdpack, uppDia, init = 0._TKG) ! reset the lower.
435mat(:,2:ndim+1)
436 0.100000E+01, -0.360922E+00
437 0.000000E+00, 0.100000E+01
438call setMatDetSqrtLog(mat(:,2:ndim+1), uppDia, detSqrtLog, info, mat(:,1:ndim), transHerm)
439if (info /= 0) error stop 'Cholesky factorization failed.'
440detSqrtLog
441-0.697829649E-1
442mat
443 0.100000E+01, 0.100000E+01, -0.360922E+00
444-0.360922E+00, 0.932596E+00, 0.100000E+01
445getMatChol(mat(:,2:ndim+1), uppDia)
446 0.100000E+01, -0.360922E+00
447 0.000000E+00, 0.932596E+00
448getMatMulTraceLog(getMatChol(mat(:,2:ndim+1), uppDia)) ! for comparison.
449-0.697829649E-1
450detSqrtLog - getMatMulTraceLog(getMatChol(mat(:,2:ndim+1), uppDia)) ! must be one.
451+0.00000000
452
453
454ndim = getUnifRand(1, 5)
455ndim
456+2
457call setResized(mat, [ndim, ndim + 1])
458mat(:,1:ndim) = getCovRand(mold = 1._TKG, ndim = ndim)
459mat(:,1:ndim)
460 0.100000E+01, 0.939151E+00
461 0.939151E+00, 0.100000E+01
462mat(:,1:ndim) = getMatCopy(rdpack, mat(:,1:ndim), rdpack, lowDia, init = 0._TKG) ! reset the upper.
463mat(:,1:ndim)
464 0.100000E+01, 0.000000E+00
465 0.939151E+00, 0.100000E+01
466call setMatDetSqrtLog(mat(:,1:ndim), lowDia, detSqrtLog, info, mat(:,2:ndim+1), transHerm)
467if (info /= 0) error stop 'Cholesky factorization failed.'
468detSqrtLog
469-1.06855559
470mat
471 0.100000E+01, 0.100000E+01, 0.939151E+00
472 0.939151E+00, 0.100000E+01, 0.343504E+00
473getMatChol(mat(:,1:ndim), lowDia)
474 0.100000E+01, 0.000000E+00
475 0.939151E+00, 0.343504E+00
476getMatMulTraceLog(getMatChol(mat(:,1:ndim), lowDia)) ! for comparison.
477-1.06855559
478detSqrtLog - getMatMulTraceLog(getMatChol(mat(:,1:ndim), lowDia)) ! must be one.
479+0.00000000
480
481mat(:,2:ndim+1) = getCovRand(mold = 1._TKG, ndim = ndim)
482mat(:,2:ndim+1) = getMatCopy(rdpack, mat(:,2:ndim+1), rdpack, uppDia, init = 0._TKG) ! reset the lower.
483mat(:,2:ndim+1)
484 0.100000E+01, -0.150104E+00
485 0.000000E+00, 0.100000E+01
486call setMatDetSqrtLog(mat(:,2:ndim+1), uppDia, detSqrtLog, info, mat(:,1:ndim), transHerm)
487if (info /= 0) error stop 'Cholesky factorization failed.'
488detSqrtLog
489-0.113944411E-1
490mat
491 0.100000E+01, 0.100000E+01, -0.150104E+00
492-0.150104E+00, 0.988670E+00, 0.100000E+01
493getMatChol(mat(:,2:ndim+1), uppDia)
494 0.100000E+01, -0.150104E+00
495 0.000000E+00, 0.988670E+00
496getMatMulTraceLog(getMatChol(mat(:,2:ndim+1), uppDia)) ! for comparison.
497-0.113944411E-1
498detSqrtLog - getMatMulTraceLog(getMatChol(mat(:,2:ndim+1), uppDia)) ! must be one.
499+0.00000000
500
501
502ndim = getUnifRand(1, 5)
503ndim
504+3
505call setResized(mat, [ndim, ndim + 1])
506mat(:,1:ndim) = getCovRand(mold = 1._TKG, ndim = ndim)
507mat(:,1:ndim)
508 0.100000E+01, -0.691536E+00, 0.430271E+00
509-0.691536E+00, 0.100000E+01, -0.160519E+00
510 0.430271E+00, -0.160519E+00, 0.100000E+01
511mat(:,1:ndim) = getMatCopy(rdpack, mat(:,1:ndim), rdpack, lowDia, init = 0._TKG) ! reset the upper.
512mat(:,1:ndim)
513 0.100000E+01, 0.000000E+00, 0.000000E+00
514-0.691536E+00, 0.100000E+01, 0.000000E+00
515 0.430271E+00, -0.160519E+00, 0.100000E+01
516call setMatDetSqrtLog(mat(:,1:ndim), lowDia, detSqrtLog, info, mat(:,2:ndim+1), transHerm)
517if (info /= 0) error stop 'Cholesky factorization failed.'
518detSqrtLog
519-0.450205237
520mat
521 0.100000E+01, 0.100000E+01, -0.691536E+00, 0.430271E+00
522-0.691536E+00, 0.100000E+01, 0.722342E+00, 0.189701E+00
523 0.430271E+00, -0.160519E+00, 0.100000E+01, 0.882542E+00
524getMatChol(mat(:,1:ndim), lowDia)
525 0.100000E+01, 0.000000E+00, 0.000000E+00
526-0.691536E+00, 0.722342E+00, 0.000000E+00
527 0.430271E+00, 0.189701E+00, 0.882542E+00
528getMatMulTraceLog(getMatChol(mat(:,1:ndim), lowDia)) ! for comparison.
529-0.450205237
530detSqrtLog - getMatMulTraceLog(getMatChol(mat(:,1:ndim), lowDia)) ! must be one.
531+0.00000000
532
533mat(:,2:ndim+1) = getCovRand(mold = 1._TKG, ndim = ndim)
534mat(:,2:ndim+1) = getMatCopy(rdpack, mat(:,2:ndim+1), rdpack, uppDia, init = 0._TKG) ! reset the lower.
535mat(:,2:ndim+1)
536 0.100000E+01, 0.739000E+00, 0.697827E+00
537 0.000000E+00, 0.100000E+01, 0.831346E+00
538 0.000000E+00, 0.000000E+00, 0.100000E+01
539call setMatDetSqrtLog(mat(:,2:ndim+1), uppDia, detSqrtLog, info, mat(:,1:ndim), transHerm)
540if (info /= 0) error stop 'Cholesky factorization failed.'
541detSqrtLog
542-1.00787342
543mat
544 0.100000E+01, 0.100000E+01, 0.739000E+00, 0.697827E+00
545 0.739000E+00, 0.673705E+00, 0.100000E+01, 0.831346E+00
546 0.697827E+00, 0.468531E+00, 0.541772E+00, 0.100000E+01
547getMatChol(mat(:,2:ndim+1), uppDia)
548 0.100000E+01, 0.739000E+00, 0.697827E+00
549 0.000000E+00, 0.673705E+00, 0.468531E+00
550 0.000000E+00, 0.000000E+00, 0.541772E+00
551getMatMulTraceLog(getMatChol(mat(:,2:ndim+1), uppDia)) ! for comparison.
552-1.00787342
553detSqrtLog - getMatMulTraceLog(getMatChol(mat(:,2:ndim+1), uppDia)) ! must be one.
554+0.00000000
555
556
557mat
558( 0.900000E+01, 0.000000E+00), ( 0.300000E+01, 0.300000E+01), ( 0.300000E+01, -0.300000E+01)
559( 0.300000E+01, -0.300000E+01), ( 0.180000E+02, 0.000000E+00), ( 0.800000E+01, -0.600000E+01)
560( 0.300000E+01, 0.300000E+01), ( 0.800000E+01, 0.600000E+01), ( 0.430000E+02, 0.000000E+00)
561call setResized(tmp, [ndim, ndim + 1])
562tmp(:,1:ndim) = getMatCopy(rdpack, mat, rdpack, lowDia, init = (0._TKG, 0._TKG)) ! reset the upper.
563tmp(:,1:ndim)
564( 0.900000E+01, 0.000000E+00), ( 0.000000E+00, 0.000000E+00), ( 0.000000E+00, 0.000000E+00)
565( 0.300000E+01, -0.300000E+01), ( 0.180000E+02, 0.000000E+00), ( 0.000000E+00, 0.000000E+00)
566( 0.300000E+01, 0.300000E+01), ( 0.800000E+01, 0.600000E+01), ( 0.430000E+02, 0.000000E+00)
567call setMatDetSqrtLog(tmp(:,1:ndim), lowDia, detSqrtLog, info, tmp(:,2:ndim+1), transHerm)
568if (info /= 0) error stop 'Cholesky factorization failed.'
569detSqrtLog
570+4.27666616
571tmp
572( 0.900000E+01, 0.000000E+00), ( 0.300000E+01, 0.000000E+00), ( 0.100000E+01, 0.100000E+01), ( 0.100000E+01, -0.100000E+01)
573( 0.300000E+01, -0.300000E+01), ( 0.180000E+02, 0.000000E+00), ( 0.400000E+01, 0.000000E+00), ( 0.200000E+01, -0.100000E+01)
574( 0.300000E+01, 0.300000E+01), ( 0.800000E+01, 0.600000E+01), ( 0.430000E+02, 0.000000E+00), ( 0.600000E+01, 0.000000E+00)
575getMatMulTraceLog(getMatChol(tmp(:,1:ndim), lowDia)) ! for comparison.
576(+4.27666616, +0.00000000)
577detSqrtLog - getMatMulTraceLog(getMatChol(tmp(:,1:ndim), lowDia)) ! must be one
578(+0.00000000, +0.00000000)
579
580tmp(:,2:ndim+1) = getMatCopy(rdpack, mat(:,1:ndim), rdpack, uppDia, init = (0._TKG, 0._TKG)) ! reset the upper.
581tmp(:,2:ndim+1)
582( 0.900000E+01, 0.000000E+00), ( 0.300000E+01, 0.300000E+01), ( 0.300000E+01, -0.300000E+01)
583( 0.000000E+00, 0.000000E+00), ( 0.180000E+02, 0.000000E+00), ( 0.800000E+01, -0.600000E+01)
584( 0.000000E+00, 0.000000E+00), ( 0.000000E+00, 0.000000E+00), ( 0.430000E+02, 0.000000E+00)
585call setMatDetSqrtLog(tmp(:,2:ndim+1), uppDia, detSqrtLog, info, tmp(:,1:ndim), transHerm)
586if (info /= 0) error stop 'Cholesky factorization failed.'
587detSqrtLog
588+4.27666616
589tmp
590( 0.300000E+01, 0.000000E+00), ( 0.900000E+01, 0.000000E+00), ( 0.300000E+01, 0.300000E+01), ( 0.300000E+01, -0.300000E+01)
591( 0.100000E+01, -0.100000E+01), ( 0.400000E+01, 0.000000E+00), ( 0.180000E+02, 0.000000E+00), ( 0.800000E+01, -0.600000E+01)
592( 0.100000E+01, 0.100000E+01), ( 0.200000E+01, 0.100000E+01), ( 0.600000E+01, 0.000000E+00), ( 0.430000E+02, 0.000000E+00)
593getMatMulTraceLog(getMatChol(tmp(:,2:ndim+1), uppDia)) ! for comparison.
594(+4.27666616, +0.00000000)
595detSqrtLog - getMatMulTraceLog(getMatChol(tmp(:,2:ndim+1), uppDia)) ! must be one.
596(+0.00000000, +0.00000000)
597
598
599mat
600( 0.250000E+02, 0.000000E+00), (-0.500000E+01, -0.500000E+01), ( 0.100000E+02, 0.500000E+01)
601(-0.500000E+01, 0.500000E+01), ( 0.510000E+02, 0.000000E+00), ( 0.400000E+01, -0.600000E+01)
602( 0.100000E+02, -0.500000E+01), ( 0.400000E+01, 0.600000E+01), ( 0.710000E+02, 0.000000E+00)
603call setResized(tmp, [ndim, ndim + 1])
604tmp(:,1:ndim) = getMatCopy(rdpack, mat, rdpack, lowDia, init = (0._TKG, 0._TKG)) ! reset the upper.
605tmp(:,1:ndim)
606( 0.250000E+02, 0.000000E+00), ( 0.000000E+00, 0.000000E+00), ( 0.000000E+00, 0.000000E+00)
607(-0.500000E+01, 0.500000E+01), ( 0.510000E+02, 0.000000E+00), ( 0.000000E+00, 0.000000E+00)
608( 0.100000E+02, -0.500000E+01), ( 0.400000E+01, 0.600000E+01), ( 0.710000E+02, 0.000000E+00)
609call setMatDetSqrtLog(tmp(:,1:ndim), lowDia, detSqrtLog, info, tmp(:,2:ndim+1), transHerm)
610if (info /= 0) error stop 'Cholesky factorization failed.'
611detSqrtLog
612+5.63478947
613tmp
614( 0.250000E+02, 0.000000E+00), ( 0.500000E+01, 0.000000E+00), (-0.100000E+01, -0.100000E+01), ( 0.200000E+01, 0.100000E+01)
615(-0.500000E+01, 0.500000E+01), ( 0.510000E+02, 0.000000E+00), ( 0.700000E+01, 0.000000E+00), ( 0.100000E+01, -0.100000E+01)
616( 0.100000E+02, -0.500000E+01), ( 0.400000E+01, 0.600000E+01), ( 0.710000E+02, 0.000000E+00), ( 0.800000E+01, 0.000000E+00)
617getMatMulTraceLog(getMatChol(tmp(:,1:ndim), lowDia)) ! for comparison.
618(+5.63478947, +0.00000000)
619detSqrtLog - getMatMulTraceLog(getMatChol(tmp(:,1:ndim), lowDia)) ! must be one
620(+0.00000000, +0.00000000)
621
622tmp(:,2:ndim+1) = getMatCopy(rdpack, mat(:,1:ndim), rdpack, uppDia, init = (0._TKG, 0._TKG)) ! reset the upper.
623tmp(:,2:ndim+1)
624( 0.250000E+02, 0.000000E+00), (-0.500000E+01, -0.500000E+01), ( 0.100000E+02, 0.500000E+01)
625( 0.000000E+00, 0.000000E+00), ( 0.510000E+02, 0.000000E+00), ( 0.400000E+01, -0.600000E+01)
626( 0.000000E+00, 0.000000E+00), ( 0.000000E+00, 0.000000E+00), ( 0.710000E+02, 0.000000E+00)
627call setMatDetSqrtLog(tmp(:,2:ndim+1), uppDia, detSqrtLog, info, tmp(:,1:ndim), transHerm)
628if (info /= 0) error stop 'Cholesky factorization failed.'
629detSqrtLog
630+5.63478947
631tmp
632( 0.500000E+01, 0.000000E+00), ( 0.250000E+02, 0.000000E+00), (-0.500000E+01, -0.500000E+01), ( 0.100000E+02, 0.500000E+01)
633(-0.100000E+01, 0.100000E+01), ( 0.700000E+01, 0.000000E+00), ( 0.510000E+02, 0.000000E+00), ( 0.400000E+01, -0.600000E+01)
634( 0.200000E+01, -0.100000E+01), ( 0.100000E+01, 0.100000E+01), ( 0.800000E+01, 0.000000E+00), ( 0.710000E+02, 0.000000E+00)
635getMatMulTraceLog(getMatChol(tmp(:,2:ndim+1), uppDia)) ! for comparison.
636(+5.63478947, +0.00000000)
637detSqrtLog - getMatMulTraceLog(getMatChol(tmp(:,2:ndim+1), uppDia)) ! must be one.
638(+0.00000000, +0.00000000)
639
640
Test:
test_pm_matrixDet


Final Remarks


If you believe this algorithm or its documentation can be improved, we appreciate your contribution and help to edit this page's documentation and source file on GitHub.
For details on the naming abbreviations, see this page.
For details on the naming conventions, see this page.
This software is distributed under the MIT license with additional terms outlined below.

  1. If you use any parts or concepts from this library to any extent, please acknowledge the usage by citing the relevant publications of the ParaMonte library.
  2. If you regenerate any parts/ideas from this library in a programming environment other than those currently supported by this ParaMonte library (i.e., other than C, C++, Fortran, MATLAB, Python, R), please also ask the end users to cite this original ParaMonte library.

This software is available to the public under a highly permissive license.
Help us justify its continued development and maintenance by acknowledging its benefit to society, distributing it, and contributing to it.

Author:
Amir Shahmoradi, Apr 21, 2017, 1:43 PM, Institute for Computational Engineering and Sciences (ICES), The University of Texas Austin

Definition at line 822 of file pm_matrixDet.F90.


The documentation for this interface was generated from the following file: