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

Return a random positive-definite power-law-distributed (correlation) matrix.
More...

Detailed Description

Return a random positive-definite power-law-distributed (correlation) matrix.

See the documentation of pm_distCov for details.

Parameters
[in,out]rng: The input/output scalar that can be an object of,
  1. type rngf_type, implying the use of intrinsic Fortran uniform RNG.
  2. type xoshiro256ssw_type, implying the use of xoshiro256** uniform RNG.
[out]rand: The output 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 random (optionally power-law-distributed determinant) positive-definite matrix.
The output rand can of complex type if and only if the optional input argument method is missing.
[in,out]method: The input/output scalar constant that can be one of the following:
  1. The scalar input constant dvine implying the use of the Dvine algorithm for generating random covariance matrices whose determinants are power-law distributed with exponent eta.
    In this case, the argument method has intent(in).
  2. A scalar output variable of type onion_type such as onion implying the use of the Onion algorithm for generating random covariance matrices whose determinants are power-law distributed with exponent eta.
    In this case, the argument method has intent(out).
    If the Cholesky factorization within the Onion algorithm fails, methodinfo will be 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.
    Otherwise, the info component of the onion method is set to 0.
The resulting matrix distribution from dvine and onion are identically distributed but onion method tends to have slightly faster runtime.
The larger eta is, the more the output random matrix looks like the Identity matrix.
Setting eta = 0. corresponds to a uniform distribution of the output matrix over the space of positive-definite correlation matrices.
See the description of the output argument rand for more information on the effects of eta on the off-diagonal elements of the output positive-definite matrix.
(optional. If missing the Gram method is used for random matrix generation. It must be missing for output rand of type complex.)
[in]eta: The input non-negative scalar of type real of the same kind as the output argument rand.
The larger eta is, the more the output random matrix looks like the Identity matrix.
Setting eta = 0. corresponds to a uniform distribution of the output matrix over the space of positive-definite correlation matrices.
See the description of the output argument rand for more information on the effects of eta on the off-diagonal elements of the output positive-definite matrix.
(optional. It must be present if and only if the input argument method is also present.)
[in]scale: The input scalar or contiguous vector of size ndim of type real of the same kind as the output argument rand, representing the scale of the matrix (e.g., the standard deviation of a covariance matrix) along each dimension.
(optional. default = 1.)


Possible calling interfaces

! Default (Gram) method.
call setCovRand(rng, rand(1:ndim, 1:ndim))
call setCovRand(rng, rand(1:ndim, 1:ndim), scale)
call setCovRand(rng, rand(1:ndim, 1:ndim), scale(1:ndim))
! Other methods.
call setCovRand(rng, rand(1:ndim, 1:ndim), method, eta)
call setCovRand(rng, rand(1:ndim, 1:ndim), method, eta, scale)
call setCovRand(rng, rand(1:ndim, 1:ndim), method, eta, scale(1:ndim))
Return a random positive-definite power-law-distributed (correlation) matrix.
Definition: pm_distCov.F90:787
This module contains classes and procedures for generating random matrices distributed on the space o...
Definition: pm_distCov.F90:72
Warning
The condition 0 <= eta must hold for the corresponding input arguments.
The condition all([0 < scale]) must hold for the corresponding input arguments.
The condition size(rand, 1) == size(rand, 2) must hold for the corresponding input arguments.
The condition rank(scale) == 0 .or. all(size(scale) == shape(rand)) must hold for the corresponding input arguments.
These conditions are verified only if the library is built with the preprocessor macro CHECK_ENABLED=1.
Beware that when the input argument scale is missing, the diagonal elements of the output correlation matrix are not enforced to match 1.
As such, numerical matrix multiplication errors may lead to diagonal matrix values slightly deviating from 1. If you need such a guarantee on the diagonal elements of the output random correlation matrix, use getCovRand.


Example usage

1program example
2
3 use pm_kind, only: SK
4 use pm_kind, only: IK, LK
5 use pm_io, only: field_type
6 use pm_io, only: display_type
7 use pm_matrixChol, only: setChoLow
10 use pm_distCov, only: setCovRand, dvine, onion
12 use pm_arrayResize, only: setResized
13 use pm_matrixDet, only: getMatDet
14
15 implicit none
16
17 integer(IK) :: itry, ndim
18
19 type(display_type) :: disp
20 disp = display_type(file = SK_"main.out.F90", format = field_type(complex = SK_"math"))
21
22 call disp%skip()
23 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
24 call disp%show("!Gram method for real covariance.")
25 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
26 call disp%skip()
27
28 block
29
30 use pm_kind, only: TKG => RKS ! all kinds are supported.
31 real(TKG), allocatable :: scale(:)
32 real(TKG), allocatable :: rand(:,:)
33 do itry = 1, 5
34
35 call disp%skip()
36 call disp%show("ndim = getUnifRand(3, 9)")
37 ndim = getUnifRand(3, 9)
38 call disp%show("ndim")
39 call disp%show( ndim )
40 call disp%show("call setResized(rand, [ndim, ndim])")
41 call setResized(rand, [ndim, ndim])
42 call disp%show("scale = getUnifRand(1, 10, ndim)")
43 scale = getUnifRand(1, 10, ndim)
44 call disp%show("scale")
45 call disp%show( scale )
46 call disp%skip()
47
48 call disp%show("call setCovRand(rngf, rand)")
49 call setCovRand(rngf, rand)
50 call disp%show("rand")
51 call disp%show( rand )
52 call disp%show("isMatClass(rand, posdefmat)")
53 call disp%show( isMatClass(rand, posdefmat) )
54 call disp%show("isMatClass(rand, hermitian)")
55 call disp%show( isMatClass(rand, hermitian) )
56 call disp%show("getMatDet(rand)")
57 call disp%show( getMatDet(rand) )
58
59 call disp%show("call setCovRand(rngf, rand, scale(1))")
60 call setCovRand(rngf, rand, scale(1))
61 call disp%show("rand")
62 call disp%show( rand )
63 call disp%show("isMatClass(rand, posdefmat)")
64 call disp%show( isMatClass(rand, posdefmat) )
65 call disp%show("isMatClass(rand, hermitian)")
66 call disp%show( isMatClass(rand, hermitian) )
67 call disp%show("getMatDet(rand)")
68 call disp%show( getMatDet(rand) )
69 call disp%skip()
70
71 call disp%show("call setCovRand(rngf, rand, scale)")
72 call setCovRand(rngf, rand, scale)
73 call disp%show("rand")
74 call disp%show( rand )
75 call disp%show("isMatClass(rand, posdefmat)")
76 call disp%show( isMatClass(rand, posdefmat) )
77 call disp%show("isMatClass(rand, hermitian)")
78 call disp%show( isMatClass(rand, hermitian) )
79 call disp%show("getMatDet(rand)")
80 call disp%show( getMatDet(rand) )
81 call disp%skip()
82
83 end do
84
85 end block
86
87 call disp%skip()
88 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
89 call disp%show("!Gram method for complex covariance.")
90 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
91 call disp%skip()
92
93 block
94
95 use pm_kind, only: TKG => RKD ! all kinds are supported.
96 real(TKG), allocatable :: scale(:)
97 complex(TKG), allocatable :: rand(:,:)
98 do itry = 1, 5
99
100 call disp%skip()
101 call disp%show("ndim = getUnifRand(3, 9)")
102 ndim = getUnifRand(3, 9)
103 call disp%show("ndim")
104 call disp%show( ndim )
105 call disp%show("call setResized(rand, [ndim, ndim])")
106 call setResized(rand, [ndim, ndim])
107 call disp%show("scale = getUnifRand(1, 10, ndim)")
108 scale = getUnifRand(1, 10, ndim)
109 call disp%show("scale")
110 call disp%show( scale )
111 call disp%skip()
112
113 call disp%show("call setCovRand(rngf, rand)")
114 call setCovRand(rngf, rand)
115 call disp%show("rand")
116 call disp%show( rand )
117 call disp%show("isMatClass(rand, posdefmat)")
118 call disp%show( isMatClass(rand, posdefmat) )
119 call disp%show("isMatClass(rand, hermitian)")
120 call disp%show( isMatClass(rand, hermitian) )
121 call disp%show("getMatDet(rand)")
122 call disp%show( getMatDet(rand) )
123
124 call disp%show("call setCovRand(rngf, rand, scale(1))")
125 call setCovRand(rngf, rand, scale(1))
126 call disp%show("rand")
127 call disp%show( rand )
128 call disp%show("isMatClass(rand, posdefmat)")
129 call disp%show( isMatClass(rand, posdefmat) )
130 call disp%show("isMatClass(rand, hermitian)")
131 call disp%show( isMatClass(rand, hermitian) )
132 call disp%show("getMatDet(rand)")
133 call disp%show( getMatDet(rand) )
134 call disp%skip()
135
136 call disp%show("call setCovRand(rngf, rand, scale)")
137 call setCovRand(rngf, rand, scale)
138 call disp%show("rand")
139 call disp%show( rand )
140 call disp%show("isMatClass(rand, posdefmat)")
141 call disp%show( isMatClass(rand, posdefmat) )
142 call disp%show("isMatClass(rand, hermitian)")
143 call disp%show( isMatClass(rand, hermitian) )
144 call disp%show("getMatDet(rand)")
145 call disp%show( getMatDet(rand) )
146 call disp%skip()
147
148 end do
149
150 end block
151
152 call disp%skip()
153 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%")
154 call disp%show("!Dvine and Onion methods.")
155 call disp%show("!%%%%%%%%%%%%%%%%%%%%%%%%")
156 call disp%skip()
157
158 block
159
160 use pm_kind, only: TKG => RKS ! all kinds are supported.
161 real(TKG), allocatable :: rand(:,:)
162 real(TKG) :: eta, scale
163 do itry = 1, 10
164
165 call disp%skip()
166 call disp%show("eta = getUnifRand(1, 10)")
167 eta = getUnifRand(1, 10)
168 call disp%show("ndim = getUnifRand(2, 5)")
169 ndim = getUnifRand(2, 5)
170 call disp%show("call setResized(rand, [ndim, ndim])")
171 call setResized(rand, [ndim, ndim])
172 call disp%show("scale = getUnifRand(1, 10)")
173 scale = getUnifRand(1, 10)
174 call disp%skip()
175
176 call disp%show("call setCovRand(rngf, rand, dvine, eta)")
177 call setCovRand(rngf, rand, dvine, eta)
178 call disp%show("onion%info")
179 call disp%show( onion%info )
180 call disp%show("rand")
181 call disp%show( rand )
182 call disp%show("isMatClass(rand, posdefmat)")
183 call disp%show( isMatClass(rand, posdefmat) )
184 call disp%show("isMatClass(rand, hermitian)")
185 call disp%show( isMatClass(rand, hermitian) )
186 call disp%show("getMatDet(rand)")
187 call disp%show( getMatDet(rand) )
188
189 call disp%show("call setCovRand(rngf, rand, onion, eta)")
190 call setCovRand(rngf, rand, onion, eta)
191 call disp%show("onion%info")
192 call disp%show( onion%info )
193 call disp%show("rand")
194 call disp%show( rand )
195 call disp%show("isMatClass(rand, posdefmat)")
196 call disp%show( isMatClass(rand, posdefmat) )
197 call disp%show("isMatClass(rand, hermitian)")
198 call disp%show( isMatClass(rand, hermitian) )
199 call disp%show("getMatDet(rand)")
200 call disp%show( getMatDet(rand) )
201
202 call disp%show("call setCovRand(rngf, rand, dvine, eta, scale)")
203 call setCovRand(rngf, rand, dvine, eta, scale)
204 call disp%show("rand")
205 call disp%show( rand )
206 call disp%show("isMatClass(rand, posdefmat)")
207 call disp%show( isMatClass(rand, posdefmat) )
208 call disp%show("isMatClass(rand, hermitian)")
209 call disp%show( isMatClass(rand, hermitian) )
210 call disp%show("getMatDet(rand)")
211 call disp%show( getMatDet(rand) )
212 call disp%skip()
213
214 call disp%show("call setCovRand(rngf, rand, onion, eta, scale)")
215 call setCovRand(rngf, rand, onion, eta, scale)
216 call disp%show("onion%info")
217 call disp%show( onion%info )
218 call disp%show("rand")
219 call disp%show( rand )
220 call disp%show("isMatClass(rand, posdefmat)")
221 call disp%show( isMatClass(rand, posdefmat) )
222 call disp%show("isMatClass(rand, hermitian)")
223 call disp%show( isMatClass(rand, hermitian) )
224 call disp%show("getMatDet(rand)")
225 call disp%show( getMatDet(rand) )
226 call disp%skip()
227
228 end do
229
230 call disp%skip()
231 call disp%show("ndim = getUnifRand(2, 10)")
232 ndim = getUnifRand(2, 10)
233 call disp%show("call setResized(rand, [ndim, ndim])")
234 call setResized(rand, [ndim, ndim])
235 call disp%skip()
236
237 call disp%show("call setCovRand(rngf, rand, dvine, eta = 0._TKG, scale = [(real(itry, TKG), itry = 1, ndim)])")
238 call setCovRand(rngf, rand, dvine, eta = 0._TKG, scale = [(real(itry, TKG), itry = 1, ndim)])
239 call disp%show("rand")
240 call disp%show( rand )
241 call disp%show("isMatClass(rand, posdefmat)")
242 call disp%show( isMatClass(rand, posdefmat) )
243 call disp%skip()
244
245 call disp%show("call setCovRand(rngf, rand, onion, eta = 0._TKG, scale = [(real(itry, TKG), itry = 1, ndim)])")
246 call setCovRand(rngf, rand, onion, eta = 0._TKG, scale = [(real(itry, TKG), itry = 1, ndim)])
247 call disp%show("onion%info")
248 call disp%show( onion%info )
249 call disp%show("rand")
250 call disp%show( rand )
251 call disp%show("isMatClass(rand, posdefmat)")
252 call disp%show( isMatClass(rand, posdefmat) )
253 call disp%skip()
254
255 end block
256
257end program example
Allocate or resize (shrink or expand) an input allocatable scalar string or array of rank 1....
Generate and return a (collection) of random vector(s) of size ndim from the ndim-dimensional MultiVa...
Generate and return a scalar or a contiguous array of rank 1 of length s1 of randomly uniformly distr...
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
[LEGACY code] Return the lower-triangle of the Cholesky factorization of the symmetric positive-def...
Generate and return .true. if and only if the input matrix is of the specified input class.
Generate and return the determinant of the input general square matrix.
This module contains procedures and generic interfaces for resizing allocatable arrays of various typ...
type(onion_type) onion
The scalar module variable object of type onion_type implying the use of the Onion algorithm for gene...
Definition: pm_distCov.F90:331
type(dvine_type), parameter dvine
The scalar constant of type dvine_type implying the use of the Dvine algorithm for generating random ...
Definition: pm_distCov.F90:238
This module contains classes and procedures for computing various statistical quantities related to t...
This module contains classes and procedures for computing various statistical quantities related to t...
type(rngf_type) rngf
The scalar constant object of type rngf_type whose presence signified the use of the Fortran intrinsi...
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 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 RKD
The double precision real kind in Fortran mode. On most platforms, this is an 64-bit real kind.
Definition: pm_kind.F90:568
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 abstract and concrete derived types that are required for compile-time resolutio...
type(posdefmat_type), parameter posdefmat
This is a scalar parameter object of type hermitian_type that is exclusively used to signify the Herm...
type(hermitian_type), parameter hermitian
This is a scalar parameter object of type hermitian_type that is exclusively used to signify the Herm...
This module contains procedures and generic interfaces relevant to the computation of the determinant...
Generate and return an object of type display_type.
Definition: pm_io.F90:10282
The derived type that can be used for constructing containers of format or left and right delimiters ...
Definition: pm_io.F90:482

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!Gram method for real covariance.
4!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
5
6
7ndim = getUnifRand(3, 9)
8ndim
9+4
10call setResized(rand, [ndim, ndim])
11scale = getUnifRand(1, 10, ndim)
12scale
13+6.00000000, +3.00000000, +10.0000000, +3.00000000
14
15call setCovRand(rngf, rand)
16rand
17+1.00000000, -0.607202053, -0.637740970, -0.547712743
18-0.607202053, +0.999999762, +0.999207258, -0.156479537
19-0.637740970, +0.999207258, +0.999999881, -0.120625295
20-0.547712743, -0.156479537, -0.120625295, +1.00000012
22T
24T
25getMatDet(rand)
26+0.137529719E-5
27call setCovRand(rngf, rand, scale(1))
28rand
29+36.0000000, +13.1712132, -3.95407987, +5.32567883
30+13.1712132, +35.9999924, -21.2181053, +17.6303234
31-3.95407987, -21.2181053, +36.0000000, +13.3233232
32+5.32567883, +17.6303234, +13.3233232, +36.0000038
34T
36T
37getMatDet(rand)
38+67756.2734
39
40call setCovRand(rngf, rand, scale)
41rand
42+36.0000000, +6.08032608, -8.81618404, -10.6031647
43+6.08032608, +9.00000000, +24.1664600, -8.61673069
44-8.81618404, +24.1664600, +99.9999847, -20.0050793
45-10.6031647, -8.61673069, -20.0050793, +8.99999905
47T
49T
50getMatDet(rand)
51+29.2561016
52
53
54ndim = getUnifRand(3, 9)
55ndim
56+3
57call setResized(rand, [ndim, ndim])
58scale = getUnifRand(1, 10, ndim)
59scale
60+9.00000000, +10.0000000, +8.00000000
61
62call setCovRand(rngf, rand)
63rand
64+1.00000000, +0.783989906, -0.674711525
65+0.783989906, +0.999999881, -0.969199538
66-0.674711525, -0.969199538, +1.00000012
68T
70T
71getMatDet(rand)
72+0.161256138E-1
73call setCovRand(rngf, rand, scale(1))
74rand
75+81.0000000, -74.5279846, -57.5895309
76-74.5279846, +80.9999924, +44.2866096
77-57.5895309, +44.2866096, +80.9999924
79T
81T
82getMatDet(rand)
83+34185.4648
84
85call setCovRand(rngf, rand, scale)
86rand
87+81.0000000, +65.3768387, -36.0589676
88+65.3768387, +100.000000, +18.2492867
89-36.0589676, +18.2492867, +64.0000076
91T
93T
94getMatDet(rand)
95+1812.26831
96
97
98ndim = getUnifRand(3, 9)
99ndim
100+6
101call setResized(rand, [ndim, ndim])
102scale = getUnifRand(1, 10, ndim)
103scale
104+6.00000000, +3.00000000, +10.0000000, +7.00000000, +2.00000000, +5.00000000
105
106call setCovRand(rngf, rand)
107rand
108+1.00000000, +0.670397580, -0.647957921, -0.164709136, +0.583028853, -0.301006526
109+0.670397580, +0.999999881, -0.874956191, -0.481517136, +0.549441814, +0.725783557E-1
110-0.647957921, -0.874956191, +0.999999821, +0.741656959, -0.504804969, -0.238775432
111-0.164709136, -0.481517136, +0.741656959, +0.999999881, +0.103976876, -0.139477789
112+0.583028853, +0.549441814, -0.504804969, +0.103976876, +1.00000000, +0.374837965
113-0.301006526, +0.725783557E-1, -0.238775432, -0.139477789, +0.374837965, +0.999999881
115T
117T
118getMatDet(rand)
119+0.612580974E-3
120call setCovRand(rngf, rand, scale(1))
121rand
122+36.0000000, -15.9641104, +9.85854435, -21.9038925, +14.5169621, +5.42454338
123-15.9641104, +35.9999962, -28.5266743, +25.9438381, +12.5339375, -22.6445885
124+9.85854435, -28.5266743, +35.9999962, -7.68231773, -20.3564091, +11.8755932
125-21.9038925, +25.9438381, -7.68231773, +35.9999962, -1.92234778, -26.2596760
126+14.5169621, +12.5339375, -20.3564091, -1.92234778, +36.0000000, -21.0867615
127+5.42454338, -22.6445885, +11.8755932, -26.2596760, -21.0867615, +35.9999962
129T
131T
132getMatDet(rand)
133+114327.844
134
135call setCovRand(rngf, rand, scale)
136rand
137+36.0000000, -9.82194901, +12.7345095, -3.56191635, +4.12076521, -17.2265530
138-9.82194901, +9.00000000, +16.0429077, -14.0670309, +0.381992877, +9.73865604
139+12.7345095, +16.0429077, +99.9999924, -56.8835449, +0.894868374E-1, +25.2015533
140-3.56191635, -14.0670309, -56.8835449, +49.0000000, -3.31724262, -12.7640724
141+4.12076521, +0.381992877, +0.894868374E-1, -3.31724262, +4.00000048, -6.64766788
142-17.2265530, +9.73865604, +25.2015533, -12.7640724, -6.64766788, +25.0000019
144T
146T
147getMatDet(rand)
148+165501.781
149
150
151ndim = getUnifRand(3, 9)
152ndim
153+8
154call setResized(rand, [ndim, ndim])
155scale = getUnifRand(1, 10, ndim)
156scale
157+7.00000000, +10.0000000, +3.00000000, +10.0000000, +7.00000000, +2.00000000, +7.00000000, +1.00000000
158
159call setCovRand(rngf, rand)
160rand
161+0.999999881, +0.932833791, -0.694727182, +0.340731204, -0.330846459, +0.174946308, +0.695621789, -0.293543994
162+0.932833791, +1.00000000, -0.801351666, +0.500382304, -0.390431345, +0.745074525E-1, +0.504800558, -0.468121380
163-0.694727182, -0.801351666, +1.00000012, -0.622935295, +0.578649223, -0.384351820, -0.283424228, +0.330920368
164+0.340731204, +0.500382304, -0.622935295, +0.999999940, -0.822480142, +0.449057817, +0.338483930, +0.425164104E-1
165-0.330846459, -0.390431345, +0.578649223, -0.822480142, +1.00000000, -0.784671664, -0.281936318, -0.806260258E-1
166+0.174946308, +0.745074525E-1, -0.384351820, +0.449057817, -0.784671664, +1.00000012, +0.386945695, +0.259479702
167+0.695621789, +0.504800558, -0.283424228, +0.338483930, -0.281936318, +0.386945695, +1.00000000, +0.161374629
168-0.293543994, -0.468121380, +0.330920368, +0.425164104E-1, -0.806260258E-1, +0.259479702, +0.161374629, +1.00000000
170T
172T
173getMatDet(rand)
174+0.936835931E-6
175call setCovRand(rngf, rand, scale(1))
176rand
177+49.0000000, -46.6558113, -12.1458063, +20.8278484, +13.2235899, +10.7070847, -36.0320396, +18.4112339
178-46.6558113, +49.0000076, +1.22345459, -20.7341843, -3.48749733, -10.7935467, +38.2470398, -27.0036316
179-12.1458063, +1.22345459, +48.9999924, +3.15711951, -14.3523932, -23.9516144, +8.48469543, +16.4751625
180+20.8278484, -20.7341843, +3.15711951, +49.0000076, -23.2696800, +12.9882441, +9.29045486, +26.1748333
181+13.2235899, -3.48749733, -14.3523932, -23.2696800, +49.0000000, -19.9384289, -15.6480350, -25.5257282
182+10.7070847, -10.7935467, -23.9516144, +12.9882441, -19.9384289, +49.0000076, -10.3639507, +1.67178953
183-36.0320396, +38.2470398, +8.48469543, +9.29045486, -15.6480350, -10.3639507, +48.9999924, -8.03328037
184+18.4112339, -27.0036316, +16.4751625, +26.1748333, -25.5257282, +1.67178953, -8.03328037, +49.0000076
186T
188T
189getMatDet(rand)
190+761912.625
191
192call setCovRand(rngf, rand, scale)
193rand
194+48.9999924, -47.8065071, +0.787280917, -4.64001751, -14.9500217, -0.775342286, -22.4217110, +3.22581267
195-47.8065071, +99.9999924, -15.4403381, +50.0532722, +4.83650017, -0.363408089, +47.8149719, -5.91230011
196+0.787280917, -15.4403381, +9.00000191, -4.87562656, +3.24086332, +3.39702082, -5.30184412, +0.528273046
197-4.64001751, +50.0532722, -4.87562656, +100.000000, -50.8122292, +5.55411339, -4.63661957, -6.18737793
198-14.9500217, +4.83650017, +3.24086332, -50.8122292, +48.9999847, +0.519130707, +32.6832962, +2.50101399
199-0.775342286, -0.363408089, +3.39702082, +5.55411339, +0.519130707, +4.00000000, -0.222875535, +0.155383393
200-22.4217110, +47.8149719, -5.30184412, -4.63661957, +32.6832962, -0.222875535, +48.9999962, -1.02455914
201+3.22581267, -5.91230011, +0.528273046, -6.18737793, +2.50101399, +0.155383393, -1.02455914, +0.999999940
203T
205T
206getMatDet(rand)
207+15420.9053
208
209
210ndim = getUnifRand(3, 9)
211ndim
212+9
213call setResized(rand, [ndim, ndim])
214scale = getUnifRand(1, 10, ndim)
215scale
216+2.00000000, +7.00000000, +4.00000000, +4.00000000, +4.00000000, +8.00000000, +9.00000000, +10.0000000, +7.00000000
217
218call setCovRand(rngf, rand)
219rand
220+1.00000000, -0.863053083, -0.499363065, +0.590936303, +0.503084123, -0.476946533, +0.574725643E-1, +0.540428758E-1, +0.632998049
221-0.863053083, +1.00000000, +0.191976726E-1, -0.491667122, -0.239261240, +0.374010772, -0.213920832, -0.333160043, -0.466541767
222-0.499363065, +0.191976726E-1, +1.00000000, -0.561094284, -0.724674821, +0.199216068, +0.926321894E-1, +0.507102728, -0.333975345
223+0.590936303, -0.491667122, -0.561094284, +1.00000000, +0.744778633, -0.777084287E-2, +0.421764702, -0.182073280, +0.756997317E-1
224+0.503084123, -0.239261240, -0.724674821, +0.744778633, +0.999999881, -0.328844547, -0.226961255, -0.580158114, +0.116216071
225-0.476946533, +0.374010772, +0.199216068, -0.777084287E-2, -0.328844547, +1.00000012, +0.682225108, +0.399040133, -0.269291937
226+0.574725643E-1, -0.213920832, +0.926321894E-1, +0.421764702, -0.226961255, +0.682225108, +0.999999940, +0.621668637, -0.144103855
227+0.540428758E-1, -0.333160043, +0.507102728, -0.182073280, -0.580158114, +0.399040133, +0.621668637, +0.999999881, +0.329894833E-1
228+0.632998049, -0.466541767, -0.333975345, +0.756997317E-1, +0.116216071, -0.269291937, -0.144103855, +0.329894833E-1, +1.00000000
230T
232T
233getMatDet(rand)
234+0.546071218E-12
235call setCovRand(rngf, rand, scale(1))
236rand
237+4.00000000, +3.69166160, -0.676362813, +1.53341615, -0.277488917, -0.997197926, -1.04324114, -1.26370764, +0.340714186
238+3.69166160, +4.00000000, -0.351941437, +1.86836648, -0.197398067, -0.377988219, -0.474917322, -1.57892191, +1.06760609
239-0.676362813, -0.351941437, +4.00000000, +3.28583193, +1.17186558, -0.895448446, +2.05779529, -2.02394605, +2.08349657
240+1.53341615, +1.86836648, +3.28583193, +3.99999976, +1.39138305, -1.58286810, +1.69477510, -2.38440514, +2.48896408
241-0.277488917, -0.197398067, +1.17186558, +1.39138305, +4.00000000, -1.19381618, +2.67669034, +0.617098391, +0.423802376
242-0.997197926, -0.377988219, -0.895448446, -1.58286810, -1.19381618, +3.99999976, -0.201324925, -0.277451247, -1.61344695
243-1.04324114, -0.474917322, +2.05779529, +1.69477510, +2.67669034, -0.201324925, +4.00000000, +0.339222074, +0.660300970
244-1.26370764, -1.57892191, -2.02394605, -2.38440514, +0.617098391, -0.277451247, +0.339222074, +3.99999928, -1.82959795
245+0.340714186, +1.06760609, +2.08349657, +2.48896408, +0.423802376, -1.61344695, +0.660300970, -1.82959795, +4.00000000
247T
249T
250getMatDet(rand)
251+0.318054197E-3
252
253call setCovRand(rngf, rand, scale)
254rand
255+3.99999952, +6.94543982, +0.155016966E-1, +5.63472509, -1.67644095, +10.0074263, -0.995428264, -3.18645835, +4.56106091
256+6.94543982, +48.9999924, -22.4333954, +0.734390974, -15.4319572, -3.70864463, +24.1159592, +10.1568604, +8.75650787
257+0.155016966E-1, -22.4333954, +16.0000038, +7.82718945, +7.80707312, +17.5217705, -23.1779633, -12.6105022, +0.744934618
258+5.63472509, +0.734390974, +7.82718945, +16.0000000, +6.62431145, +24.4801483, -15.7340050, -1.98638141, +13.9765568
259-1.67644095, -15.4319572, +7.80707312, +6.62431145, +15.9999981, +9.31474018, -4.02922249, +15.0052624, +5.39749098
260+10.0074263, -3.70864463, +17.5217705, +24.4801483, +9.31474018, +63.9999924, -12.7945099, -20.9939957, +2.50078368
261-0.995428264, +24.1159592, -23.1779633, -15.7340050, -4.02922249, -12.7945099, +80.9999924, +45.3658104, -34.1223450
262-3.18645835, +10.1568604, -12.6105022, -1.98638141, +15.0052624, -20.9939957, +45.3658104, +100.000000, -3.06978846
263+4.56106091, +8.75650787, +0.744934618, +13.9765568, +5.39749098, +2.50078368, -34.1223450, -3.06978846, +49.0000076
265T
267T
268getMatDet(rand)
269+2297580.50
270
271
272!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
273!Gram method for complex covariance.
274!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
275
276
277ndim = getUnifRand(3, 9)
278ndim
279+4
280call setResized(rand, [ndim, ndim])
281scale = getUnifRand(1, 10, ndim)
282scale
283+9.0000000000000000, +8.0000000000000000, +9.0000000000000000, +4.0000000000000000
284
285call setCovRand(rngf, rand)
286rand
287+1.0000000000000000+0.0000000000000000i, -0.20180659324403502-0.11953182405609890i, +0.13963947477913755+0.19403983734931152i, -0.26265226608130338-0.41784568122692017i
288-0.20180659324403502+0.11953182405609890i, +1.0000000000000000+0.0000000000000000i, -0.47864999392567709+0.45720922367635281i, +0.54686010568054022-0.37809564990600630i
289+0.13963947477913755-0.19403983734931152i, -0.47864999392567709-0.45720922367635281i, +0.99999999999999978+0.0000000000000000i, -0.62682295321059955+0.22057380660122833i
290-0.26265226608130338+0.41784568122692017i, +0.54686010568054022+0.37809564990600630i, -0.62682295321059955-0.22057380660122833i, +1.0000000000000002+0.0000000000000000i
292T
294T
295getMatDet(rand)
296+0.97356066846538122E-1-0.53384497744928176E-17i
297call setCovRand(rngf, rand, scale(1))
298rand
299+81.000000000000000+0.0000000000000000i, -6.3753381749735549-18.031016080224877i, +16.396633025246381+33.695497127611922i, +15.407174067172694-3.1668577056836025i
300-6.3753381749735549+18.031016080224877i, +80.999999999999986+0.0000000000000000i, +29.886291914149687+21.596827325757506i, -21.050354740651571-0.23722782815243670i
301+16.396633025246381-33.695497127611922i, +29.886291914149687-21.596827325757506i, +81.000000000000000+0.0000000000000000i, -45.301861006938559-16.415315750486545i
302+15.407174067172694+3.1668577056836025i, -21.050354740651571+0.23722782815243670i, -45.301861006938559+16.415315750486545i, +81.000000000000014+0.0000000000000000i
304T
306T
307getMatDet(rand)
308+8341061.7250382919+0.26519867572506582E-10i
309
310call setCovRand(rngf, rand, scale)
311rand
312+81.000000000000000+0.0000000000000000i, -28.402866228462457+49.932829787513299i, +30.736198035990288+23.505227614558549i, -18.040085965453649-20.739528602106446i
313-28.402866228462457-49.932829787513299i, +64.000000000000000+0.0000000000000000i, +13.050352167414442-52.344520156064242i, -3.8206038657962247+25.618257226775576i
314+30.736198035990288-23.505227614558549i, +13.050352167414442+52.344520156064242i, +80.999999999999986+0.0000000000000000i, -20.135178734597865-7.3114796721506226i
315-18.040085965453649+20.739528602106446i, -3.8206038657962247-25.618257226775576i, -20.135178734597865+7.3114796721506226i, +15.999999999999995+0.0000000000000000i
317T
319T
320getMatDet(rand)
321+39983.573824705236-0.10913936421275139E-9i
322
323
324ndim = getUnifRand(3, 9)
325ndim
326+7
327call setResized(rand, [ndim, ndim])
328scale = getUnifRand(1, 10, ndim)
329scale
330+10.000000000000000, +4.0000000000000000, +7.0000000000000000, +2.0000000000000000, +9.0000000000000000, +2.0000000000000000, +7.0000000000000000
331
332call setCovRand(rngf, rand)
333rand
334+0.99999999999999978+0.0000000000000000i, +0.68046953576093561+0.43531934389214638i, -0.46718586718140898-0.33948626130008835i, +0.51566563754673411+0.39908145941168965i, +0.47545567707218012-0.45712340958796788i, -0.29488799607859278+0.10217746471386757i, -0.72904271920285230E-1+0.25801323994817812E-1i
335+0.68046953576093561-0.43531934389214638i, +1.0000000000000002+0.0000000000000000i, -0.63048107190345226-0.36061339193778108i, +0.71016390135277507+0.13083933479229992i, +0.60937439773577476E-1-0.55731896401903747i, -0.19605081251852999+0.45252801857048414i, +0.19118963900371011-0.60299568060957612E-1i
336-0.46718586718140898+0.33948626130008835i, -0.63048107190345226+0.36061339193778108i, +1.0000000000000000+0.0000000000000000i, -0.37557280907407631-0.13928434569303152i, +0.16179647126022528+0.89042520832421479E-1i, +0.60418256725510350E-1-0.48902122133717324i, -0.12070351585010886+0.30993857406430447i
337+0.51566563754673411-0.39908145941168965i, +0.71016390135277507-0.13083933479229992i, -0.37557280907407631+0.13928434569303152i, +1.0000000000000002+0.0000000000000000i, +0.37454448127400652-0.45175534522265104i, +0.22107459904707305+0.45742949491680118i, -0.21928783651564210-0.68212088406256394E-1i
338+0.47545567707218012+0.45712340958796788i, +0.60937439773577476E-1+0.55731896401903747i, +0.16179647126022528-0.89042520832421479E-1i, +0.37454448127400652+0.45175534522265104i, +0.99999999999999978+0.0000000000000000i, -0.29076463053532356E-1+0.16688379900271474i, -0.13248049925411851+0.45900189634859756E-1i
339-0.29488799607859278-0.10217746471386757i, -0.19605081251852999-0.45252801857048414i, +0.60418256725510350E-1+0.48902122133717324i, +0.22107459904707305-0.45742949491680118i, -0.29076463053532356E-1-0.16688379900271474i, +0.99999999999999989+0.0000000000000000i, -0.16356374398343226-0.29404297613369768i
340-0.72904271920285230E-1-0.25801323994817812E-1i, +0.19118963900371011+0.60299568060957612E-1i, -0.12070351585010886-0.30993857406430447i, -0.21928783651564210+0.68212088406256394E-1i, -0.13248049925411851-0.45900189634859756E-1i, -0.16356374398343226+0.29404297613369768i, +1.0000000000000000+0.0000000000000000i
342T
344T
345getMatDet(rand)
346+0.33270956084922134E-4+0.41674021004911577E-18i
347call setCovRand(rngf, rand, scale(1))
348rand
349+100.00000000000000+0.0000000000000000i, +40.899330441071257+68.891386627842095i, -19.222587293714227-6.4214966601012478i, +27.728231005108782+45.154863922083528i, -43.688214132173606-54.443913281842669i, -18.417611238372690+37.533816202182201i, +19.101241208510388-34.310158510470117i
350+40.899330441071257-68.891386627842095i, +100.00000000000000+0.0000000000000000i, -3.9471300906014681+52.351555919443186i, +69.191682142886620+16.711162233235395i, -79.792223294103053-7.8621464289633316i, +39.517262753925152+24.158427241157433i, -25.718409773641870-48.800893578847010i
351-19.222587293714227+6.4214966601012478i, -3.9471300906014681-52.351555919443186i, +100.00000000000000+0.0000000000000000i, +52.894551007006484-23.691319661312015i, -18.769164549487794+26.157448188232529i, +2.9077054398517044-4.8239413870796168i, -12.354043313472271-7.2774712440426885i
352+27.728231005108782-45.154863922083528i, +69.191682142886620-16.711162233235395i, +52.894551007006484+23.691319661312015i, +100.00000000000000+0.0000000000000000i, -78.614127788466277-12.197739238687998i, +43.104487948233967+19.905925220879389i, -28.944293477574778-59.680468569800802i
353-43.688214132173606+54.443913281842669i, -79.792223294103053+7.8621464289633316i, -18.769164549487794-26.157448188232529i, -78.614127788466277+12.197739238687998i, +99.999999999999986+0.0000000000000000i, -38.776373059058130-4.7067947156534684i, +52.509075339097237+30.292652641268237i
354-18.417611238372690-37.533816202182201i, +39.517262753925152-24.158427241157433i, +2.9077054398517044+4.8239413870796168i, +43.104487948233967-19.905925220879389i, -38.776373059058130+4.7067947156534684i, +100.00000000000004+0.0000000000000000i, -33.407666039024875-30.466146758069890i
355+19.101241208510388+34.310158510470117i, -25.718409773641870+48.800893578847010i, -12.354043313472271+7.2774712440426885i, -28.944293477574778+59.680468569800802i, +52.509075339097237-30.292652641268237i, -33.407666039024875+30.466146758069890i, +100.00000000000007+0.0000000000000000i
357T
359T
360getMatDet(rand)
361+74374284.110985488-0.27000904083251953E-4i
362
363call setCovRand(rngf, rand, scale)
364rand
365+100.00000000000000+0.0000000000000000i, +25.535011517098887+30.697741719986084i, -6.9385436298740846+16.529168070282942i, -5.6698805742247949-4.1871202075067506i, +31.509447354801733+15.839516662254947i, -0.76070462547688511+7.9944917313211095i, +23.712661525493985+32.789694189764496i
366+25.535011517098887-30.697741719986084i, +16.000000000000000+0.0000000000000000i, +4.2479403851125053+5.7655939813598112i, -2.7350415238710060+0.40480090696672294i, +13.692951103648968-5.0234327449212781i, +2.4697904995334916+2.0652582216117090i, +15.717489869168439+0.82143638790737139i
367-6.9385436298740846-16.529168070282942i, +4.2479403851125053-5.7655939813598112i, +48.999999999999986+0.0000000000000000i, +6.1559273571152175+2.4966381780194129i, +9.4538890772132564+10.522779715791833i, +8.5992010973895141+0.20133628328452136i, -1.3906630494092307-23.055342860749469i
368-5.6698805742247949+4.1871202075067506i, -2.7350415238710060-0.40480090696672294i, +6.1559273571152175-2.4966381780194129i, +4.0000000000000009+0.0000000000000000i, -7.4499219603040618+1.3135929842057950i, +1.1550838651833200+0.78177425318592819i, -3.7836673298247931-3.4112616061734644i
369+31.509447354801733-15.839516662254947i, +13.692951103648968+5.0234327449212781i, +9.4538890772132564-10.522779715791833i, -7.4499219603040618-1.3135929842057950i, +81.000000000000000+0.0000000000000000i, +2.0088865667177331-10.232614750211493i, +14.264844788177886+9.7418656012634894i
370-0.76070462547688511-7.9944917313211095i, +2.4697904995334916-2.0652582216117090i, +8.5992010973895141-0.20133628328452136i, +1.1550838651833200-0.78177425318592819i, +2.0088865667177331+10.232614750211493i, +3.9999999999999996+0.0000000000000000i, +2.3351404252123547-3.4235223228760461i
371+23.712661525493985-32.789694189764496i, +15.717489869168439-0.82143638790737139i, -1.3906630494092307+23.055342860749469i, -3.7836673298247931+3.4112616061734644i, +14.264844788177886-9.7418656012634894i, +2.3351404252123547+3.4235223228760461i, +48.999999999999993+0.0000000000000000i
373T
375T
376getMatDet(rand)
377+261.45747863464351-0.27088731258118059E-9i
378
379
380ndim = getUnifRand(3, 9)
381ndim
382+3
383call setResized(rand, [ndim, ndim])
384scale = getUnifRand(1, 10, ndim)
385scale
386+4.0000000000000000, +1.0000000000000000, +9.0000000000000000
387
388call setCovRand(rngf, rand)
389rand
390+1.0000000000000000+0.0000000000000000i, +0.66664517195900164E-1-0.55778315856007166i, -0.64780309429153726+0.39454616811992582i
391+0.66664517195900164E-1+0.55778315856007166i, +1.0000000000000000+0.0000000000000000i, -0.44605605646738794E-1-0.60678178512269298i
392-0.64780309429153726-0.39454616811992582i, -0.44605605646738794E-1+0.60678178512269298i, +0.99999999999999978+0.0000000000000000i
394T
396T
397getMatDet(rand)
398+0.16901184944484912+0.0000000000000000i
399call setCovRand(rngf, rand, scale(1))
400rand
401+15.999999999999996+0.0000000000000000i, +5.1707898972597981+10.842471616178694i, +5.9776473663490997-8.2955295288357309i
402+5.1707898972597981-10.842471616178694i, +16.000000000000007+0.0000000000000000i, +0.46715959982477173-4.4907430196746914i
403+5.9776473663490997+8.2955295288357309i, +0.46715959982477173+4.4907430196746914i, +16.000000000000000+0.0000000000000000i
405T
407T
408getMatDet(rand)
409+700.54113261805458-0.99212852236579138E-13i
410
411call setCovRand(rngf, rand, scale)
412rand
413+16.000000000000000+0.0000000000000000i, +1.7831872289489683+0.83434874776835777i, +5.8001978907778291+19.881042585706773i
414+1.7831872289489683-0.83434874776835777i, +1.0000000000000000+0.0000000000000000i, +5.4205514609353829+3.4794510027154617i
415+5.8001978907778291-19.881042585706773i, +5.4205514609353829-3.4794510027154617i, +80.999999999999986+0.0000000000000000i
417T
419T
420getMatDet(rand)
421+394.31559919753153-0.10768368839829166E-13i
422
423
424ndim = getUnifRand(3, 9)
425ndim
426+4
427call setResized(rand, [ndim, ndim])
428scale = getUnifRand(1, 10, ndim)
429scale
430+8.0000000000000000, +1.0000000000000000, +3.0000000000000000, +3.0000000000000000
431
432call setCovRand(rngf, rand)
433rand
434+1.0000000000000000+0.0000000000000000i, -0.94516336795558251E-1-0.57383921974168950i, -0.47001006542591955+0.45183765337598714i, -0.18074635604780748-0.45113853065796528i
435-0.94516336795558251E-1+0.57383921974168950i, +1.0000000000000000+0.0000000000000000i, -0.71953449818763082-0.32825633004439136E-1i, +0.68985842774589323-0.38623817305365948i
436-0.47001006542591955-0.45183765337598714i, -0.71953449818763082+0.32825633004439136E-1i, +1.0000000000000000+0.0000000000000000i, -0.44956558064339674+0.32223017862235753i
437-0.18074635604780748+0.45113853065796528i, +0.68985842774589323+0.38623817305365948i, -0.44956558064339674-0.32223017862235753i, +0.99999999999999989+0.0000000000000000i
439T
441T
442getMatDet(rand)
443+0.52134341571258881E-2+0.65052130349130266E-18i
444call setCovRand(rngf, rand, scale(1))
445rand
446+64.000000000000000+0.0000000000000000i, -40.781145924801905-36.250797309460403i, +14.055436440533210+21.806493532370155i, -3.6892320899289559+33.487352423784593i
447-40.781145924801905+36.250797309460403i, +64.000000000000000+0.0000000000000000i, -32.896248614590618-20.384104778880960i, -8.8030618267047629-21.458306164399602i
448+14.055436440533210-21.806493532370155i, -32.896248614590618+20.384104778880960i, +63.999999999999993+0.0000000000000000i, +27.049577130388464+37.100946122958256i
449-3.6892320899289559-33.487352423784593i, -8.8030618267047629+21.458306164399602i, +27.049577130388464-37.100946122958256i, +64.000000000000014+0.0000000000000000i
451T
453T
454getMatDet(rand)
455+397268.41333152971+0.58207660913467407E-10i
456
457call setCovRand(rngf, rand, scale)
458rand
459+64.000000000000000+0.0000000000000000i, +3.8652973055184190-0.23702498443956815i, +1.7126150808827960-10.112924958315379i, -8.3961811162232216+6.7317157616441268i
460+3.8652973055184190+0.23702498443956815i, +0.99999999999999978+0.0000000000000000i, -0.79730815480236661+1.0550933889463234i, -1.9626870294092149-1.0129391660130835i
461+1.7126150808827960+10.112924958315379i, -0.79730815480236661-1.0550933889463234i, +9.0000000000000036+0.0000000000000000i, -0.98035866841549635+2.7148712894722822i
462-8.3961811162232216-6.7317157616441268i, -1.9626870294092149+1.0129391660130835i, -0.98035866841549635-2.7148712894722822i, +8.9999999999999964+0.0000000000000000i
464T
466T
467getMatDet(rand)
468+92.467225883435475-0.74606987254810520E-13i
469
470
471ndim = getUnifRand(3, 9)
472ndim
473+6
474call setResized(rand, [ndim, ndim])
475scale = getUnifRand(1, 10, ndim)
476scale
477+4.0000000000000000, +10.000000000000000, +2.0000000000000000, +8.0000000000000000, +7.0000000000000000, +2.0000000000000000
478
479call setCovRand(rngf, rand)
480rand
481+1.0000000000000000+0.0000000000000000i, +0.65964687443127590+0.31088497916767077i, -0.39080504712913522-0.18591567379581071E-1i, -0.46565175065257902-0.22742869946037800i, +0.36926877589192963-0.20335310376346005i, -0.29661361762439042+0.33882554249795560i
482+0.65964687443127590-0.31088497916767077i, +1.0000000000000002+0.0000000000000000i, -0.59281283384174388+0.18485576505202841i, -0.60067445116543039+0.35385024648504138i, +0.35923611845257164E-2+0.43427879550218706E-1i, -0.11978012873194639+0.36841043826427167i
483-0.39080504712913522+0.18591567379581071E-1i, -0.59281283384174388-0.18485576505202841i, +1.0000000000000000+0.0000000000000000i, +0.60725240668799818-0.49268453266290790i, -0.36388745922922505+0.19975503855183763i, -0.67400130618396004E-1-0.25666346714670540i
484-0.46565175065257902+0.22742869946037800i, -0.60067445116543039-0.35385024648504138i, +0.60725240668799818+0.49268453266290790i, +0.99999999999999989+0.0000000000000000i, -0.91787930499120152E-1+0.30445813310613308E-1i, +0.19417273463534621-0.44601325685135640i
485+0.36926877589192963+0.20335310376346005i, +0.35923611845257164E-2-0.43427879550218706E-1i, -0.36388745922922505-0.19975503855183763i, -0.91787930499120152E-1-0.30445813310613308E-1i, +1.0000000000000002+0.0000000000000000i, -0.11320973281560223+0.61971398297198610E-1i
486-0.29661361762439042-0.33882554249795560i, -0.11978012873194639-0.36841043826427167i, -0.67400130618396004E-1+0.25666346714670540i, +0.19417273463534621+0.44601325685135640i, -0.11320973281560223-0.61971398297198610E-1i, +0.99999999999999978+0.0000000000000000i
488T
490T
491getMatDet(rand)
492+0.82925978877962121E-3-0.54210108624275222E-19i
493call setCovRand(rngf, rand, scale(1))
494rand
495+16.000000000000000+0.0000000000000000i, -10.348762221587977+2.5524219954850018i, +0.42616165173157861+10.640875734585187i, +6.5093832522711814-0.61276083344164267i, -8.8987998708555942+0.71992075522027943i, +2.8144874584537010+6.8349501771515238i
496-10.348762221587977-2.5524219954850018i, +16.000000000000000+0.0000000000000000i, +8.7879451118454721-11.360611601865509i, -4.3920212070438049-0.16725498225133278i, +8.4797382769817364+0.26632781552027374i, +4.2439234701161510-7.1963545235056152i
497+0.42616165173157861-10.640875734585187i, +8.7879451118454721+11.360611601865509i, +16.000000000000000+0.0000000000000000i, -0.97946319566062390-3.7167708358642026i, +2.5297087557502138+6.1907878689366758i, +10.702817921998566-1.4552255066713875i
498+6.5093832522711814+0.61276083344164267i, -4.3920212070438049+0.16725498225133278i, -0.97946319566062390+3.7167708358642026i, +16.000000000000000+0.0000000000000000i, -1.6016197543800987-1.3348680496088103i, -3.5042666229842458+1.8422781278030413i
499-8.8987998708555942-0.71992075522027943i, +8.4797382769817364-0.26632781552027374i, +2.5297087557502138-6.1907878689366758i, -1.6016197543800987+1.3348680496088103i, +16.000000000000000+0.0000000000000000i, -4.5759180264462582-7.1175310483602541i
500+2.8144874584537010-6.8349501771515238i, +4.2439234701161510+7.1963545235056152i, +10.702817921998566+1.4552255066713875i, -3.5042666229842458-1.8422781278030413i, -4.5759180264462582+7.1175310483602541i, +16.000000000000000+0.0000000000000000i
502T
504T
505getMatDet(rand)
506+41399.733899019506-0.14551915228366852E-10i
507
508call setCovRand(rngf, rand, scale)
509rand
510+16.000000000000000+0.0000000000000000i, -24.280675864037686+25.311252563548273i, +2.9595525316084905-1.1032735768371642i, -4.9347464410394624-10.574370302343128i, -4.8311962199095646-13.482214953191610i, +2.8938391588772050+2.0713332124388448i
511-24.280675864037686-25.311252563548273i, +100.00000000000000+0.0000000000000000i, -0.59675973922481362-4.7945053595619349i, +7.5301699216422673+6.6075302661073687i, -17.381433658340207+11.924944877477234i, -3.6295367195711661-9.7886743944100196i
512+2.9595525316084905+1.1032735768371642i, -0.59675973922481362+4.7945053595619349i, +3.9999999999999991+0.0000000000000000i, +7.8783965735078940-9.8523955104707781i, +3.3497455655402053-4.0387076883922450i, -0.77746042221751155+0.62745166570608957i
513-4.9347464410394624+10.574370302343128i, +7.5301699216422673-6.6075302661073687i, +7.8783965735078940+9.8523955104707781i, +63.999999999999986+0.0000000000000000i, +15.653262104748773+4.3075731829953146i, -4.2639364241000592-0.63727945021334742i
514-4.8311962199095646+13.482214953191610i, -17.381433658340207-11.924944877477234i, +3.3497455655402053+4.0387076883922450i, +15.653262104748773-4.3075731829953146i, +48.999999999999993+0.0000000000000000i, -0.27844176755848066-0.71095320933782369i
515+2.8938391588772050-2.0713332124388448i, -3.6295367195711661+9.7886743944100196i, -0.77746042221751155-0.62745166570608957i, -4.2639364241000592+0.63727945021334742i, -0.27844176755848066+0.71095320933782369i, +4.0000000000000000+0.0000000000000000i
517T
519T
520getMatDet(rand)
521+149325.42574507976-0.36379788070917130E-10i
522
523
524!%%%%%%%%%%%%%%%%%%%%%%%%
525!Dvine and Onion methods.
526!%%%%%%%%%%%%%%%%%%%%%%%%
527
528
529eta = getUnifRand(1, 10)
530ndim = getUnifRand(2, 5)
531call setResized(rand, [ndim, ndim])
532scale = getUnifRand(1, 10)
533
534call setCovRand(rngf, rand, dvine, eta)
535onion%info
536+0
537rand
538+1.00000000, +0.124258757, +0.144787431, +0.270395279E-1
539+0.124258757, +1.00000000, -0.132766470, +0.294473708
540+0.144787431, -0.132766470, +1.00000000, -0.575901754E-2
541+0.270395279E-1, +0.294473708, -0.575901754E-2, +1.00000000
543T
545T
546getMatDet(rand)
547+0.858290315
548call setCovRand(rngf, rand, onion, eta)
549onion%info
550+0
551rand
552+1.00000000, +0.948636532E-1, +0.570246316E-1, -0.456157833
553+0.948636532E-1, +1.00000000, +0.191861242, +0.492983200E-1
554+0.570246316E-1, +0.191861242, +1.00000000, +0.850529820E-1
555-0.456157833, +0.492983200E-1, +0.850529820E-1, +1.00000000
557T
559T
560getMatDet(rand)
561+0.737778962
562call setCovRand(rngf, rand, dvine, eta, scale)
563rand
564+4.00000000, -0.264009476, +0.205822945, -0.644260406
565-0.264009476, +4.00000000, +0.956239343, +0.737017691E-1
566+0.205822945, +0.956239343, +4.00000000, +1.05530357
567-0.644260406, +0.737017691E-1, +1.05530357, +4.00000000
569T
571T
572getMatDet(rand)
573+214.331314
574
575call setCovRand(rngf, rand, onion, eta, scale)
576onion%info
577+0
578rand
579+4.00000000, -0.709191561, +1.41572726, -0.699726880
580-0.709191561, +4.00000000, -0.804665089, +0.867180824E-1
581+1.41572726, -0.804665089, +4.00000000, +0.165075362E-1
582-0.699726880, +0.867180824E-1, +0.165075362E-1, +4.00000000
584T
586T
587getMatDet(rand)
588+204.443359
589
590
591eta = getUnifRand(1, 10)
592ndim = getUnifRand(2, 5)
593call setResized(rand, [ndim, ndim])
594scale = getUnifRand(1, 10)
595
596call setCovRand(rngf, rand, dvine, eta)
597onion%info
598+0
599rand
600+1.00000000, -0.285848260, +0.718290806E-1
601-0.285848260, +1.00000000, -0.295854419
602+0.718290806E-1, -0.295854419, +1.00000000
604T
606T
607getMatDet(rand)
608+0.837750614
609call setCovRand(rngf, rand, onion, eta)
610onion%info
611+0
612rand
613+1.00000000, +0.765621662E-2, +0.500656962
614+0.765621662E-2, +1.00000000, -0.993103758E-1
615+0.500656962, -0.993103758E-1, +1.00000000
617T
619T
620getMatDet(rand)
621+0.738660157
622call setCovRand(rngf, rand, dvine, eta, scale)
623rand
624+49.0000000, +6.16029263, -7.43757391
625+6.16029263, +49.0000000, -2.02561140
626-7.43757391, -2.02561140, +49.0000000
628T
630T
631getMatDet(rand)
632+113063.492
633
634call setCovRand(rngf, rand, onion, eta, scale)
635onion%info
636+0
637rand
638+49.0000000, -1.60309744, -6.53829575
639-1.60309744, +49.0000000, -58.2667770
640-6.53829575, -58.2667770, +49.0000000
642F
644T
645getMatDet(rand)
646-52148.9492
647
648
649eta = getUnifRand(1, 10)
650ndim = getUnifRand(2, 5)
651call setResized(rand, [ndim, ndim])
652scale = getUnifRand(1, 10)
653
654call setCovRand(rngf, rand, dvine, eta)
655onion%info
656+0
657rand
658+1.00000000, -0.806381702E-1
659-0.806381702E-1, +1.00000000
661T
663T
664getMatDet(rand)
665+0.993497491
666call setCovRand(rngf, rand, onion, eta)
667onion%info
668+0
669rand
670+1.00000000, +0.158905029
671+0.158905029, +1.00000000
673T
675T
676getMatDet(rand)
677+0.974749207
678call setCovRand(rngf, rand, dvine, eta, scale)
679rand
680+100.000000, -19.3109341
681-19.3109341, +100.000000
683T
685T
686getMatDet(rand)
687+9627.08789
688
689call setCovRand(rngf, rand, onion, eta, scale)
690onion%info
691+0
692rand
693+100.000000, +19.3522091
694+19.3522091, +100.000000
696T
698T
699getMatDet(rand)
700+9625.49219
701
702
703eta = getUnifRand(1, 10)
704ndim = getUnifRand(2, 5)
705call setResized(rand, [ndim, ndim])
706scale = getUnifRand(1, 10)
707
708call setCovRand(rngf, rand, dvine, eta)
709onion%info
710+0
711rand
712+1.00000000, -0.223070204, +0.108405948, -0.207901001E-3
713-0.223070204, +1.00000000, +0.203596994, -0.117064409
714+0.108405948, +0.203596994, +1.00000000, -0.104648732
715-0.207901001E-3, -0.117064409, -0.104648732, +1.00000000
717T
719T
720getMatDet(rand)
721+0.868815422
722call setCovRand(rngf, rand, onion, eta)
723onion%info
724+0
725rand
726+1.00000000, +0.136029482, +0.177782625, +0.251700938
727+0.136029482, +1.00000000, -0.134319574, -0.296457887
728+0.177782625, -0.134319574, +1.00000000, +0.931486487E-1
729+0.251700938, -0.296457887, +0.931486487E-1, +1.00000000
731T
733T
734getMatDet(rand)
735+0.763597786
736call setCovRand(rngf, rand, dvine, eta, scale)
737rand
738+36.0000000, +3.93651295, -12.7299852, -6.42633438
739+3.93651295, +36.0000000, +0.119677752, -5.01387882
740-12.7299852, +0.119677752, +36.0000000, +19.1610661
741-6.42633438, -5.01387882, +19.1610661, +36.0000000
743T
745T
746getMatDet(rand)
747+1008653.50
748
749call setCovRand(rngf, rand, onion, eta, scale)
750onion%info
751+0
752rand
753+36.0000000, -2.53202677, -5.23935986, +4.53819752
754-2.53202677, +36.0000000, -32.1713676, +18.0697021
755-5.23935986, -32.1713676, +36.0000000, -57.1106682
756+4.53819752, +18.0697021, -57.1106682, +36.0000000
758F
760T
761getMatDet(rand)
762-1846907.62
763
764
765eta = getUnifRand(1, 10)
766ndim = getUnifRand(2, 5)
767call setResized(rand, [ndim, ndim])
768scale = getUnifRand(1, 10)
769
770call setCovRand(rngf, rand, dvine, eta)
771onion%info
772+0
773rand
774+1.00000000, -0.330212712E-1, +0.210489154, +0.707058430
775-0.330212712E-1, +1.00000000, +0.259490460, -0.626730695E-1
776+0.210489154, +0.259490460, +1.00000000, +0.437535524
777+0.707058430, -0.626730695E-1, +0.437535524, +1.00000000
779T
781T
782getMatDet(rand)
783+0.351102054
784call setCovRand(rngf, rand, onion, eta)
785onion%info
786+0
787rand
788+1.00000000, +0.574183464E-1, +0.171530053, +0.283750981
789+0.574183464E-1, +1.00000000, +0.137752131, +0.192851096
790+0.171530053, +0.137752131, +1.00000000, -0.599809587E-1
791+0.283750981, +0.192851096, -0.599809587E-1, +1.00000000
793T
795T
796getMatDet(rand)
797+0.827517927
798call setCovRand(rngf, rand, dvine, eta, scale)
799rand
800+9.00000000, -2.35516238, +2.65871930, +2.31648350
801-2.35516238, +9.00000000, -0.586642146, -6.11838436
802+2.65871930, -0.586642146, +9.00000000, +3.13148808
803+2.31648350, -6.11838436, +3.13148808, +9.00000000
805T
807T
808getMatDet(rand)
809+2483.19019
810
811call setCovRand(rngf, rand, onion, eta, scale)
812onion%info
813+0
814rand
815+9.00000000, -1.80684328, +1.74620855, -1.47584176
816-1.80684328, +9.00000000, -4.73778629, +7.40137577
817+1.74620855, -4.73778629, +9.00000000, -7.45849466
818-1.47584176, +7.40137577, -7.45849466, +9.00000000
820T
822T
823getMatDet(rand)
824+469.842041
825
826
827eta = getUnifRand(1, 10)
828ndim = getUnifRand(2, 5)
829call setResized(rand, [ndim, ndim])
830scale = getUnifRand(1, 10)
831
832call setCovRand(rngf, rand, dvine, eta)
833onion%info
834+0
835rand
836+1.00000000, +0.321278811
837+0.321278811, +1.00000000
839T
841T
842getMatDet(rand)
843+0.896779895
844call setCovRand(rngf, rand, onion, eta)
845onion%info
846+0
847rand
848+1.00000000, -0.214562118
849-0.214562118, +1.00000000
851T
853T
854getMatDet(rand)
855+0.953963101
856call setCovRand(rngf, rand, dvine, eta, scale)
857rand
858+81.0000000, +15.8622971
859+15.8622971, +81.0000000
861T
863T
864getMatDet(rand)
865+6309.38770
866
867call setCovRand(rngf, rand, onion, eta, scale)
868onion%info
869+0
870rand
871+81.0000000, +17.2423534
872+17.2423534, +81.0000000
874T
876T
877getMatDet(rand)
878+6263.70117
879
880
881eta = getUnifRand(1, 10)
882ndim = getUnifRand(2, 5)
883call setResized(rand, [ndim, ndim])
884scale = getUnifRand(1, 10)
885
886call setCovRand(rngf, rand, dvine, eta)
887onion%info
888+0
889rand
890+1.00000000, +0.357404232, +0.787377357E-2, +0.447191715
891+0.357404232, +1.00000000, +0.458437949, +0.224996239
892+0.787377357E-2, +0.458437949, +1.00000000, -0.184157062E-1
893+0.447191715, +0.224996239, -0.184157062E-1, +1.00000000
895T
897T
898getMatDet(rand)
899+0.525736034
900call setCovRand(rngf, rand, onion, eta)
901onion%info
902+0
903rand
904+1.00000000, +0.361596346E-1, -0.679276943, -0.543750823
905+0.361596346E-1, +1.00000000, -0.278235435, -0.283562422
906-0.679276943, -0.278235435, +1.00000000, -0.390929103
907-0.543750823, -0.283562422, -0.390929103, +1.00000000
909F
911T
912getMatDet(rand)
913-0.383059710
914call setCovRand(rngf, rand, dvine, eta, scale)
915rand
916+49.0000000, +12.7469683, +5.59396553, +13.6033316
917+12.7469683, +49.0000000, +19.2898960, -22.1901150
918+5.59396553, +19.2898960, +49.0000000, -26.6197014
919+13.6033316, -22.1901150, -26.6197014, +49.0000000
921T
923T
924getMatDet(rand)
925+2112440.00
926
927call setCovRand(rngf, rand, onion, eta, scale)
928onion%info
929+3
930rand
931+1.00000000, +0.794275284, -0.115199409, +13.6033316
932+12.7469683, +1.00000000, -1.34799528, -22.1901150
933+5.59396553, +19.2898960, +1.00000000, -26.6197014
934+13.6033316, -22.1901150, -26.6197014, +49.0000000
936F
938F
939getMatDet(rand)
940+88311.1797
941
942
943eta = getUnifRand(1, 10)
944ndim = getUnifRand(2, 5)
945call setResized(rand, [ndim, ndim])
946scale = getUnifRand(1, 10)
947
948call setCovRand(rngf, rand, dvine, eta)
949onion%info
950+3
951rand
952+1.00000000, -0.107081652, -0.364922762
953-0.107081652, +1.00000000, -0.187360838
954-0.364922762, -0.187360838, +1.00000000
956T
958T
959getMatDet(rand)
960+0.805617988
961call setCovRand(rngf, rand, onion, eta)
962onion%info
963+0
964rand
965+1.00000000, +0.494294643, +0.243051872
966+0.494294643, +1.00000000, +0.537788868
967+0.243051872, +0.537788868, +1.00000000
969T
971T
972getMatDet(rand)
973+0.536600828
974call setCovRand(rngf, rand, dvine, eta, scale)
975rand
976+9.00000000, +1.98896849, -0.411219120
977+1.98896849, +9.00000000, -2.47881413
978-0.411219120, -2.47881413, +9.00000000
980T
982T
983getMatDet(rand)
984+640.628296
985
986call setCovRand(rngf, rand, onion, eta, scale)
987onion%info
988+0
989rand
990+9.00000000, -4.67568970, +0.231513917
991-4.67568970, +9.00000000, +0.992481351
992+0.231513917, +0.992481351, +9.00000000
994T
996T
997getMatDet(rand)
998+520.745056
999
1000
1001eta = getUnifRand(1, 10)
1002ndim = getUnifRand(2, 5)
1003call setResized(rand, [ndim, ndim])
1004scale = getUnifRand(1, 10)
1005
1006call setCovRand(rngf, rand, dvine, eta)
1007onion%info
1008+0
1009rand
1010+1.00000000, +0.273714900, -0.183587849, -0.410694301, +0.237391114
1011+0.273714900, +1.00000000, +0.531440616, -0.328188419, +0.157563150
1012-0.183587849, +0.531440616, +1.00000000, +0.283617914, -0.799388885E-1
1013-0.410694301, -0.328188419, +0.283617914, +1.00000000, -0.548463315E-2
1014+0.237391114, +0.157563150, -0.799388885E-1, -0.548463315E-2, +1.00000000
1015isMatClass(rand, posdefmat)
1016T
1017isMatClass(rand, hermitian)
1018T
1019getMatDet(rand)
1020+0.282707572
1021call setCovRand(rngf, rand, onion, eta)
1022onion%info
1023+0
1024rand
1025+1.00000000, +0.713528395E-1, +0.236614347E-1, +0.151054919, +0.133788111E-1
1026+0.713528395E-1, +1.00000000, +0.312997550, +0.387103409E-1, +0.391311079
1027+0.236614347E-1, +0.312997550, +1.00000000, +0.166263849, -0.761416703E-1
1028+0.151054919, +0.387103409E-1, +0.166263849, +1.00000000, -0.382541597
1029+0.133788111E-1, +0.391311079, -0.761416703E-1, -0.382541597, +1.00000000
1030isMatClass(rand, posdefmat)
1031T
1032isMatClass(rand, hermitian)
1033T
1034getMatDet(rand)
1035+0.567284465
1036call setCovRand(rngf, rand, dvine, eta, scale)
1037rand
1038+1.00000000, -0.290929794, +0.153805614, +0.435251594, +0.157630444E-2
1039-0.290929794, +1.00000000, -0.297143459, -0.531557202E-2, -0.125118464
1040+0.153805614, -0.297143459, +1.00000000, +0.262270272, -0.209394787E-1
1041+0.435251594, -0.531557202E-2, +0.262270272, +1.00000000, +0.261243820
1042+0.157630444E-2, -0.125118464, -0.209394787E-1, +0.261243820, +1.00000000
1043isMatClass(rand, posdefmat)
1044T
1045isMatClass(rand, hermitian)
1046T
1047getMatDet(rand)
1048+0.528713942
1049
1050call setCovRand(rngf, rand, onion, eta, scale)
1051onion%info
1052+0
1053rand
1054+1.00000000, +0.255149603E-1, -0.402773947, -0.152608141, +0.140496776
1055+0.255149603E-1, +1.00000000, +0.299819916, -0.201248109, -0.922721252E-2
1056-0.402773947, +0.299819916, +1.00000000, +0.250953913, +0.244573921
1057-0.152608141, -0.201248109, +0.250953913, +1.00000000, +0.232021123
1058+0.140496776, -0.922721252E-2, +0.244573921, +0.232021123, +1.00000000
1059isMatClass(rand, posdefmat)
1060T
1061isMatClass(rand, hermitian)
1062T
1063getMatDet(rand)
1064+0.523975492
1065
1066
1067eta = getUnifRand(1, 10)
1068ndim = getUnifRand(2, 5)
1069call setResized(rand, [ndim, ndim])
1070scale = getUnifRand(1, 10)
1071
1072call setCovRand(rngf, rand, dvine, eta)
1073onion%info
1074+0
1075rand
1076+1.00000000, -0.706497252, +0.156208038
1077-0.706497252, +1.00000000, +0.259474277
1078+0.156208038, +0.259474277, +1.00000000
1079isMatClass(rand, posdefmat)
1080T
1081isMatClass(rand, hermitian)
1082T
1083getMatDet(rand)
1084+0.351862341
1085call setCovRand(rngf, rand, onion, eta)
1086onion%info
1087+0
1088rand
1089+1.00000000, -0.419878364E-1, -0.431069285
1090-0.419878364E-1, +1.00000000, +0.454514861
1091-0.431069285, +0.454514861, +1.00000000
1092isMatClass(rand, posdefmat)
1093T
1094isMatClass(rand, hermitian)
1095T
1096getMatDet(rand)
1097+0.622285724
1098call setCovRand(rngf, rand, dvine, eta, scale)
1099rand
1100+64.0000000, -2.47647858, +50.4180679
1101-2.47647858, +64.0000000, +17.8578033
1102+50.4180679, +17.8578033, +64.0000000
1103isMatClass(rand, posdefmat)
1104T
1105isMatClass(rand, hermitian)
1106T
1107getMatDet(rand)
1108+74195.5703
1109
1110call setCovRand(rngf, rand, onion, eta, scale)
1111onion%info
1112+0
1113rand
1114+64.0000000, +31.3292770, +11.6796083
1115+31.3292770, +64.0000000, -6.68511963
1116+11.6796083, -6.68511963, +64.0000000
1117isMatClass(rand, posdefmat)
1118T
1119isMatClass(rand, hermitian)
1120T
1121getMatDet(rand)
1122+182843.469
1123
1124
1125ndim = getUnifRand(2, 10)
1126call setResized(rand, [ndim, ndim])
1127
1128call setCovRand(rngf, rand, dvine, eta = 0._TKG, scale = [(real(itry, TKG), itry = 1, ndim)])
1129rand
1130+1.00000000, +0.246458769, -0.243736982, -1.52473831, -2.87061977, -3.05276108, -1.00652838, -2.64686298, +3.13220429, -5.18655205
1131+0.246458769, +4.00000000, -2.37410259, -1.36218476, -0.974722803, +1.87992275, -4.82264233, -11.9923048, +1.16202521, -8.33574104
1132-0.243736982, -2.37410259, +9.00000000, +3.55996442, -1.53313172, +2.43992090, +11.0935087, +13.1541252, -5.98270512, +0.921773911
1133-1.52473831, -1.36218476, +3.55996442, +16.0000000, +2.11172652, +8.31253242, +3.03086877, +4.49136877, -5.22442341, +10.3387661
1134-2.87061977, -0.974722803, -1.53313172, +2.11172652, +25.0000000, -1.01309991, -0.312538445E-1, +4.16990089, +16.7755089, +18.7293224
1135-3.05276108, +1.87992275, +2.43992090, +8.31253242, -1.01309991, +36.0000000, +6.75722027, -0.115621805, -19.9659119, +12.3576193
1136-1.00652838, -4.82264233, +11.0935087, +3.03086877, -0.312538445E-1, +6.75722027, +49.0000000, +6.15420055, +5.00263977, +1.66495907
1137-2.64686298, -11.9923048, +13.1541252, +4.49136877, +4.16990089, -0.115621805, +6.15420055, +64.0000000, -19.6195259, +31.3703384
1138+3.13220429, +1.16202521, -5.98270512, -5.22442341, +16.7755089, -19.9659119, +5.00263977, -19.6195259, +81.0000000, -43.7680244
1139-5.18655205, -8.33574104, +0.921773911, +10.3387661, +18.7293224, +12.3576193, +1.66495907, +31.3703384, -43.7680244, +100.000000
1140isMatClass(rand, posdefmat)
1141T
1142
1143call setCovRand(rngf, rand, onion, eta = 0._TKG, scale = [(real(itry, TKG), itry = 1, ndim)])
1144onion%info
1145+5
1146rand
1147+1.00000000, -0.655725002E-1, -0.919323862E-1, +0.184340984, -0.286550045, -3.05276108, -1.00652838, -2.64686298, +3.13220429, -5.18655205
1148+0.246458769, +1.00000000, -0.193759486, -0.383268856E-1, -0.226552099, +1.87992275, -4.82264233, -11.9923048, +1.16202521, -8.33574104
1149-0.243736982, -2.37410259, +1.00000000, -0.649544179, +0.293768376, +2.43992090, +11.0935087, +13.1541252, -5.98270512, +0.921773911
1150-1.52473831, -1.36218476, +3.55996442, +1.00000000, +1.72349942, +8.31253242, +3.03086877, +4.49136877, -5.22442341, +10.3387661
1151-2.87061977, -0.974722803, -1.53313172, +2.11172652, +1.00000000, -1.01309991, -0.312538445E-1, +4.16990089, +16.7755089, +18.7293224
1152-3.05276108, +1.87992275, +2.43992090, +8.31253242, -1.01309991, +36.0000000, +6.75722027, -0.115621805, -19.9659119, +12.3576193
1153-1.00652838, -4.82264233, +11.0935087, +3.03086877, -0.312538445E-1, +6.75722027, +49.0000000, +6.15420055, +5.00263977, +1.66495907
1154-2.64686298, -11.9923048, +13.1541252, +4.49136877, +4.16990089, -0.115621805, +6.15420055, +64.0000000, -19.6195259, +31.3703384
1155+3.13220429, +1.16202521, -5.98270512, -5.22442341, +16.7755089, -19.9659119, +5.00263977, -19.6195259, +81.0000000, -43.7680244
1156-5.18655205, -8.33574104, +0.921773911, +10.3387661, +18.7293224, +12.3576193, +1.66495907, +31.3703384, -43.7680244, +100.000000
1157isMatClass(rand, posdefmat)
1158F
1159
1160
Test:
test_pm_distCov
Todo:
High Priority: The current implementation of this generic interface uses a naive method of computing the Cholesky factorization with a default matrix packing for the Onion method.
The RFP packing format must be also implemented for this generic interface.
Todo:
High Priority: The current implementation of the Gram method can be significantly improved, both computationally and functionally.


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, Monday March 6, 2017, 3:22 pm, Institute for Computational Engineering and Sciences (ICES), The University of Texas at Austin.

Definition at line 787 of file pm_distCov.F90.


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