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+9
10call setResized(rand, [ndim, ndim])
11scale = getUnifRand(1, 10, ndim)
12scale
13+8.00000000, +3.00000000, +4.00000000, +10.0000000, +6.00000000, +4.00000000, +10.0000000, +10.0000000, +3.00000000
14
15call setCovRand(rngf, rand)
16rand
17+1.00000000, +0.957517326, +0.568198740, +0.542673826, +0.489191800, -0.356931627, -0.962954666E-2, +0.654745340, -0.282974422
18+0.957517326, +1.00000012, +0.755919039, +0.337433100, +0.334102571, -0.458285511, -0.631121844E-1, +0.538206458, -0.418572485
19+0.568198740, +0.755919039, +1.00000012, -0.280721724, -0.250233740, -0.308606625, -0.270803452, +0.208552212, -0.540053070
20+0.542673826, +0.337433100, -0.280721724, +0.999999881, +0.564207017, -0.326139718, +0.437263787, +0.371268362, -0.449737161E-1
21+0.489191800, +0.334102571, -0.250233740, +0.564207017, +1.00000000, -0.233650595, -0.141660780, +0.511791646, +0.288811713
22-0.356931627, -0.458285511, -0.308606625, -0.326139718, -0.233650595, +1.00000000, -0.117863551, +0.504920855E-1, +0.496963412
23-0.962954666E-2, -0.631121844E-1, -0.270803452, +0.437263787, -0.141660780, -0.117863551, +1.00000000, +0.739401132E-1, +0.155677482
24+0.654745340, +0.538206458, +0.208552212, +0.371268362, +0.511791646, +0.504920855E-1, +0.739401132E-1, +1.00000000, +0.362090945
25-0.282974422, -0.418572485, -0.540053070, -0.449737161E-1, +0.288811713, +0.496963412, +0.155677482, +0.362090945, +1.00000024
27T
29T
30getMatDet(rand)
31+0.128518263E-7
32call setCovRand(rngf, rand, scale(1))
33rand
34+64.0000000, +50.4022102, +32.4713364, -21.6564598, +40.4578247, +27.3008213, -5.24287510, -15.2494507, +15.4966755
35+50.4022102, +64.0000076, +4.51029778, +5.01695251, +6.92782593, +29.0213833, -22.3060226, -17.9950905, +15.9546604
36+32.4713364, +4.51029778, +64.0000000, -4.66325188, +42.9756355, +19.6075058, +35.9222946, -10.8078575, +4.54169369
37-21.6564598, +5.01695251, -4.66325188, +63.9999962, -27.7601547, +27.1944046, -3.22864532, -22.4624348, +12.6389408
38+40.4578247, +6.92782593, +42.9756355, -27.7601547, +63.9999924, +32.8835220, -0.750291824, -15.6556320, +21.1342888
39+27.3008213, +29.0213833, +19.6075058, +27.1944046, +32.8835220, +63.9999886, -24.8746185, -32.8423042, +37.0633850
40-5.24287510, -22.3060226, +35.9222946, -3.22864532, -0.750291824, -24.8746185, +63.9999886, +16.4446545, -28.3457317
41-15.2494507, -17.9950905, -10.8078575, -22.4624348, -15.6556320, -32.8423042, +16.4446545, +64.0000000, -0.239372253E-1
42+15.4966755, +15.9546604, +4.54169369, +12.6389408, +21.1342888, +37.0633850, -28.3457317, -0.239372253E-1, +64.0000000
44T
46T
47getMatDet(rand)
48+0.855332506E+10
49
50call setCovRand(rngf, rand, scale)
51rand
52+64.0000000, +5.92940044, +14.2702913, +67.8983994, +6.26989841, -9.10739803, -30.5183334, +42.8082962, -4.16588497
53+5.92940044, +8.99999809, -1.28502917, +3.38465500, -1.56195664, -7.32579088, -5.83927631, -10.7693233, -0.629105568
54+14.2702913, -1.28502917, +15.9999981, +2.35771537, -12.0460091, -3.79404116, -11.1901169, +17.9625244, +3.52691841
55+67.8983994, +3.38465500, +2.35771537, +99.9999924, +10.3877478, +6.49343443, -22.6729660, +50.6163979, -12.6810360
56+6.26989841, -1.56195664, -12.0460091, +10.3877478, +36.0000000, -5.13571167, -13.6847868, -7.86860085, -2.44178915
57-9.10739803, -7.32579088, -3.79404116, +6.49343443, -5.13571167, +16.0000019, +9.04471588, +8.30318069, -4.32467318
58-30.5183334, -5.83927631, -11.1901169, -22.6729660, -13.6847868, +9.04471588, +100.000000, +28.7007694, -3.43619967
59+42.8082962, -10.7693233, +17.9625244, +50.6163979, -7.86860085, +8.30318069, +28.7007694, +99.9999847, -10.3450241
60-4.16588497, -0.629105568, +3.52691841, -12.6810360, -2.44178915, -4.32467318, -3.43619967, -10.3450241, +9.00000000
62T
64T
65getMatDet(rand)
66+1882902.62
67
68
69ndim = getUnifRand(3, 9)
70ndim
71+7
72call setResized(rand, [ndim, ndim])
73scale = getUnifRand(1, 10, ndim)
74scale
75+5.00000000, +8.00000000, +5.00000000, +2.00000000, +5.00000000, +1.00000000, +3.00000000
76
77call setCovRand(rngf, rand)
78rand
79+1.00000000, -0.968285501, -0.649314642, -0.912330821E-1, -0.155782253, -0.288999975, +0.498902321
80-0.968285501, +1.00000000, +0.440165997, +0.235607222, +0.324412793, +0.168063670, -0.375256956
81-0.649314642, +0.440165997, +0.999999881, -0.451365829, -0.476461202, +0.544231296, -0.678153098
82-0.912330821E-1, +0.235607222, -0.451365829, +0.999999881, +0.682478368, -0.156975895, +0.272120059
83-0.155782253, +0.324412793, -0.476461202, +0.682478368, +1.00000012, -0.472894877, +0.637292922
84-0.288999975, +0.168063670, +0.544231296, -0.156975895, -0.472894877, +1.00000000, -0.464530885
85+0.498902321, -0.375256956, -0.678153098, +0.272120059, +0.637292922, -0.464530885, +1.00000000
87T
89T
90getMatDet(rand)
91+0.588229021E-8
92call setCovRand(rngf, rand, scale(1))
93rand
94+25.0000000, -19.9442482, +24.4878025, -16.2620773, -17.8176136, -12.7186184, +12.6411219
95-19.9442482, +25.0000000, -19.2004433, +23.3937302, +21.6040516, +7.32900810, -1.02460456
96+24.4878025, -19.2004433, +24.9999981, -16.1335182, -18.4336433, -12.1916275, +11.3499441
97-16.2620773, +23.3937302, -16.1335182, +25.0000000, +23.6711597, +7.10887766, +2.22400022
98-17.8176136, +21.6040516, -18.4336433, +23.6711597, +25.0000000, +6.44416952, +0.361280888
99-12.7186184, +7.32900810, -12.1916275, +7.10887766, +6.44416952, +25.0000038, -19.8210487
100+12.6411219, -1.02460456, +11.3499441, +2.22400022, +0.361280888, -19.8210487, +24.9999943
102T
104T
105getMatDet(rand)
106+1528.03247
107
108call setCovRand(rngf, rand, scale)
109rand
110+25.0000000, +4.43194056, +7.39055586, -2.01091743, -14.0758848, +1.15525091, +2.65164709
111+4.43194056, +63.9999847, -35.8767853, +6.31929636, +19.9159508, +5.34445381, +14.0770884
112+7.39055586, -35.8767853, +25.0000000, -5.95658064, -20.0584373, -3.00715637, -6.84581470
113-2.01091743, +6.31929636, -5.95658064, +3.99999976, +8.52338696, +1.40324390, -0.838441014
114-14.0758848, +19.9159508, -20.0584373, +8.52338696, +25.0000000, +2.24430227, +1.21620166
115+1.15525091, +5.34445381, -3.00715637, +1.40324390, +2.24430227, +1.00000012, +0.244440839
116+2.65164709, +14.0770884, -6.84581470, -0.838441014, +1.21620166, +0.244440839, +9.00000000
118T
120T
121getMatDet(rand)
122+16.7363205
123
124
125ndim = getUnifRand(3, 9)
126ndim
127+4
128call setResized(rand, [ndim, ndim])
129scale = getUnifRand(1, 10, ndim)
130scale
131+1.00000000, +10.0000000, +6.00000000, +10.0000000
132
133call setCovRand(rngf, rand)
134rand
135+0.999999881, -0.983515084, +0.607585371, -0.476067483
136-0.983515084, +1.00000000, -0.696168244, +0.466376722
137+0.607585371, -0.696168244, +1.00000000, -0.625689209
138-0.476067483, +0.466376722, -0.625689209, +1.00000000
140T
142T
143getMatDet(rand)
144+0.460834661E-2
145call setCovRand(rngf, rand, scale(1))
146rand
147+1.00000000, +0.859935403, +0.795688748, -0.468120366
148+0.859935403, +0.999999881, +0.630496800, -0.167468503
149+0.795688748, +0.630496800, +0.999999940, -0.199580088
150-0.468120366, -0.167468503, -0.199580088, +1.00000000
152T
154T
155getMatDet(rand)
156+0.399442911E-1
157
158call setCovRand(rngf, rand, scale)
159rand
160+1.00000000, -5.46995926, +1.94590950, -0.435469180
161-5.46995926, +99.9999924, -33.1928864, -12.4557476
162+1.94590950, -33.1928864, +36.0000000, -27.3139687
163-0.435469180, -12.4557476, -27.3139687, +100.000008
165T
167T
168getMatDet(rand)
169+100682.133
170
171
172ndim = getUnifRand(3, 9)
173ndim
174+7
175call setResized(rand, [ndim, ndim])
176scale = getUnifRand(1, 10, ndim)
177scale
178+3.00000000, +1.00000000, +4.00000000, +2.00000000, +6.00000000, +5.00000000, +5.00000000
179
180call setCovRand(rngf, rand)
181rand
182+1.00000000, +0.225868165, +0.473733872, +0.146691009, -0.453077018, +0.238284275, +0.318854988
183+0.225868165, +1.00000000, +0.954240263, -0.536291957, +0.403288901, -0.504550576, -0.154705748
184+0.473733872, +0.954240263, +1.00000012, -0.355038047, +0.256449938, -0.469558716, -0.499768704E-1
185+0.146691009, -0.536291957, -0.355038047, +1.00000012, -0.136459470E-1, -0.248969764, +0.285035342
186-0.453077018, +0.403288901, +0.256449938, -0.136459470E-1, +0.999999821, -0.695822537, -0.326630354
187+0.238284275, -0.504550576, -0.469558716, -0.248969764, -0.695822537, +0.999999940, +0.113475487
188+0.318854988, -0.154705748, -0.499768704E-1, +0.285035342, -0.326630354, +0.113475487, +1.00000000
190T
192T
193getMatDet(rand)
194+0.947085994E-6
195call setCovRand(rngf, rand, scale(1))
196rand
197+9.00000000, -6.79219151, -5.17774534, +3.17225838, -2.92784262, -4.81202221, -3.20212984
198-6.79219151, +8.99999905, +4.87215614, +0.535456896, +3.96795011, +0.948569655, +2.97115159
199-5.17774534, +4.87215614, +9.00000000, -5.03748512, +6.97629642, +3.51789188, +3.23459983
200+3.17225838, +0.535456896, -5.03748512, +8.99999905, -5.91726732, -4.40062666, -4.47510433
201-2.92784262, +3.96795011, +6.97629642, -5.91726732, +9.00000000, +2.28185582, +6.06398773
202-4.81202221, +0.948569655, +3.51789188, -4.40062666, +2.28185582, +8.99999905, +3.17725801
203-3.20212984, +2.97115159, +3.23459983, -4.47510433, +6.06398773, +3.17725801, +9.00000095
205T
207T
208getMatDet(rand)
209+393.917938
210
211call setCovRand(rngf, rand, scale)
212rand
213+9.00000000, -2.43390346, +4.59102201, +0.553656518, +4.19657230, +7.46858788, -0.108692795
214-2.43390346, +1.00000000, -3.28771925, -0.103933565, -0.692990482, -2.00443339, -0.934959888
215+4.59102201, -3.28771925, +16.0000000, +1.94677031, +4.36853123, +2.08205390, +2.78367543
216+0.553656518, -0.103933565, +1.94677031, +4.00000048, +4.28786755, -4.73784924, -2.57987452
217+4.19657230, -0.692990482, +4.36853123, +4.28786755, +36.0000000, -2.45145988, -6.10959625
218+7.46858788, -2.00443339, +2.08205390, -4.73784924, -2.45145988, +25.0000000, +0.482200831
219-0.108692795, -0.934959888, +2.78367543, -2.57987452, -6.10959625, +0.482200831, +25.0000019
221T
223T
224getMatDet(rand)
225+450.903381
226
227
228ndim = getUnifRand(3, 9)
229ndim
230+4
231call setResized(rand, [ndim, ndim])
232scale = getUnifRand(1, 10, ndim)
233scale
234+4.00000000, +4.00000000, +5.00000000, +6.00000000
235
236call setCovRand(rngf, rand)
237rand
238+0.999999881, +0.136644214, -0.149096757, +0.489646077
239+0.136644214, +1.00000000, +0.599548578, -0.440750450
240-0.149096757, +0.599548578, +1.00000000, +0.448246002E-1
241+0.489646077, -0.440750450, +0.448246002E-1, +0.999999881
243T
245T
246getMatDet(rand)
247+0.975263268E-1
248call setCovRand(rngf, rand, scale(1))
249rand
250+16.0000000, +5.57715702, +4.32213497, -15.2615728
251+5.57715702, +16.0000000, -10.2397299, -2.76209116
252+4.32213497, -10.2397299, +16.0000019, -5.60533381
253-15.2615728, -2.76209116, -5.60533381, +16.0000019
255T
257T
258getMatDet(rand)
259+1006.20703
260
261call setCovRand(rngf, rand, scale)
262rand
263+16.0000000, -12.6560812, +12.2601757, -10.4007177
264-12.6560812, +16.0000038, -19.1436672, +14.0240307
265+12.2601757, -19.1436672, +25.0000000, -20.5133114
266-10.4007177, +14.0240307, -20.5133114, +36.0000076
268T
270T
271getMatDet(rand)
272+493.964905
273
274
275!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
276!Gram method for complex covariance.
277!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
278
279
280ndim = getUnifRand(3, 9)
281ndim
282+8
283call setResized(rand, [ndim, ndim])
284scale = getUnifRand(1, 10, ndim)
285scale
286+7.0000000000000000, +1.0000000000000000, +5.0000000000000000, +8.0000000000000000, +10.000000000000000, +5.0000000000000000, +2.0000000000000000, +5.0000000000000000
287
288call setCovRand(rngf, rand)
289rand
290+1.0000000000000000+0.0000000000000000i, -0.49832483866328869-0.69447989667501520i, -0.47972002420122223-0.25278382851354031E-1i, -0.53992853404372335E-1-0.51257778051012259i, -0.92321424482674219E-2+0.28110491412977828i, -0.24292094570834227+0.15185500597686385i, +0.33521229221315535E-1+0.27040602525048235i, +0.20526230740400536+0.11428534027238189i
291-0.49832483866328869+0.69447989667501520i, +1.0000000000000002+0.0000000000000000i, +0.45037142274296599-0.39477567341640390i, +0.27162568505049800+0.34678278472359092i, +0.73480786140718041E-1-0.35637139544124496i, +0.19978295535539997-0.11622916662629768E-1i, -0.30314452271131270-0.72443359474813895E-1i, -0.19465409279869628-0.54120136934354579E-1i
292-0.47972002420122223+0.25278382851354031E-1i, +0.45037142274296599+0.39477567341640390i, +1.0000000000000002+0.0000000000000000i, -0.54394560430452499+0.63726998826841597i, +0.45293024053587305-0.58937766906215372i, +0.13057382952964777-0.20506064785636452i, -0.23978566018174419-0.74026106587421592E-1i, -0.26263264488359572-0.44388202677625299i
293-0.53992853404372335E-1+0.51257778051012259i, +0.27162568505049800-0.34678278472359092i, -0.54394560430452499-0.63726998826841597i, +1.0000000000000002+0.0000000000000000i, -0.74949790603658284+0.74179737503371521E-1i, -0.96334293293891468E-1-0.48203868591324292E-1i, -0.90882879842327041E-1+0.69199766857724679E-2i, -0.14666654332803933+0.50214428042729975i
294-0.92321424482674219E-2-0.28110491412977828i, +0.73480786140718041E-1+0.35637139544124496i, +0.45293024053587305+0.58937766906215372i, -0.74949790603658284-0.74179737503371521E-1i, +0.99999999999999956+0.0000000000000000i, +0.13251080006424554+0.17457757700008039i, -0.11398602463933540+0.26245745466804554E-1i, +0.26172855995687933-0.36951618262055330i
295-0.24292094570834227-0.15185500597686385i, +0.19978295535539997+0.11622916662629768E-1i, +0.13057382952964777+0.20506064785636452i, -0.96334293293891468E-1+0.48203868591324292E-1i, +0.13251080006424554-0.17457757700008039i, +1.0000000000000000+0.0000000000000000i, -0.16269451299502213-0.24434748541414592i, -0.19823304207800527-0.33681520940505671i
296+0.33521229221315535E-1-0.27040602525048235i, -0.30314452271131270+0.72443359474813895E-1i, -0.23978566018174419+0.74026106587421592E-1i, -0.90882879842327041E-1-0.69199766857724679E-2i, -0.11398602463933540-0.26245745466804554E-1i, -0.16269451299502213+0.24434748541414592i, +1.0000000000000002+0.0000000000000000i, +0.43609467885263070E-1-0.90703647025568795E-1i
297+0.20526230740400536-0.11428534027238189i, -0.19465409279869628+0.54120136934354579E-1i, -0.26263264488359572+0.44388202677625299i, -0.14666654332803933-0.50214428042729975i, +0.26172855995687933+0.36951618262055330i, -0.19823304207800527+0.33681520940505671i, +0.43609467885263070E-1+0.90703647025568795E-1i, +1.0000000000000000+0.0000000000000000i
299T
301T
302getMatDet(rand)
303+0.82399401886198851E-5-0.72844833463869829E-19i
304call setCovRand(rngf, rand, scale(1))
305rand
306+49.000000000000000+0.0000000000000000i, -4.5174804212418893-16.560085840717001i, +8.3600036428003737+5.5151436149788973i, -25.200402790443800+19.134777030244816i, -7.0022733877475867+24.643040667601944i, -6.1701016171056358+10.930932151132515i, +0.62071493896050745+19.843955916149049i, -13.769720711096404+18.317169097622067i
307-4.5174804212418893+16.560085840717001i, +49.000000000000000+0.0000000000000000i, -20.826785484140576-27.922514394669847i, -14.387897615190655-6.2383779058941622i, +14.828008380132662-3.5668237566063525i, -0.57257417561588697+2.4855103773212326i, -1.6405232157644636-19.860471709635693i, -8.6800421247209520+10.726845697064210i
308+8.3600036428003737-5.5151436149788973i, -20.826785484140576+27.922514394669847i, +49.000000000000000+0.0000000000000000i, -8.0948367601173672+1.1167928356200152i, -3.9512782479219837+6.9814797479110826i, +8.5405082421495191+13.967914642988502i, +17.981220441529441+18.475265127863040i, -6.9526446853609940+4.4284695969111763i
309-25.200402790443800-19.134777030244816i, -14.387897615190655+6.2383779058941622i, -8.0948367601173672-1.1167928356200152i, +48.999999999999993+0.0000000000000000i, +6.4403388654424258-5.8807620307221979i, -8.0523843369438275-12.284023381745669i, -6.1903642863591983-20.205428237467494i, +19.223215471338357-20.738553148365156i
310-7.0022733877475867-24.643040667601944i, +14.828008380132662+3.5668237566063525i, -3.9512782479219837-6.9814797479110826i, +6.4403388654424258+5.8807620307221979i, +49.000000000000014+0.0000000000000000i, -10.006258478171125+17.277269646371959i, +9.0995302792721020-19.655343882269822i, -2.6535601603737549+11.103063434703525i
311-6.1701016171056358-10.930932151132515i, -0.57257417561588697-2.4855103773212326i, +8.5405082421495191-13.967914642988502i, -8.0523843369438275+12.284023381745669i, -10.006258478171125-17.277269646371959i, +48.999999999999986+0.0000000000000000i, +17.529022089657385+4.4699847317828887i, +28.431876618007319+10.326940013205563i
312+0.62071493896050745-19.843955916149049i, -1.6405232157644636+19.860471709635693i, +17.981220441529441-18.475265127863040i, -6.1903642863591983+20.205428237467494i, +9.0995302792721020+19.655343882269822i, +17.529022089657385-4.4699847317828887i, +49.000000000000000+0.0000000000000000i, +14.817728086381226+4.9252500160084738i
313-13.769720711096404-18.317169097622067i, -8.6800421247209520-10.726845697064210i, -6.9526446853609940-4.4284695969111763i, +19.223215471338357+20.738553148365156i, -2.6535601603737549-11.103063434703525i, +28.431876618007319-10.326940013205563i, +14.817728086381226-4.9252500160084738i, +49.000000000000014+0.0000000000000000i
315T
317T
318getMatDet(rand)
319+1671059209.2868035-0.11974756705805618E-4i
320
321call setCovRand(rngf, rand, scale)
322rand
323+49.000000000000000+0.0000000000000000i, -6.1612860414941073+2.0012408162165567i, +2.2240426525562991+7.5546154655638489i, +3.8586660822987309+30.114708684502581i, -28.989864511941768-42.209195888284015i, -1.1378722493831885-4.0838402059374586i, +2.8740306224584269-6.1151019226214238i, +3.7497833200578716-13.681979753815369i
324-6.1612860414941073-2.0012408162165567i, +1.0000000000000002+0.0000000000000000i, -0.74338902742851332-1.7144305494668579i, +1.6484892711513575-5.0079752464818243i, +2.4492175802020539+6.1461612372403476i, -0.83266894114766343+1.3453860400162176i, -0.58692732470136544+0.53347264646785209i, -0.66377140200104745+1.8121418261847968i
325+2.2240426525562991-7.5546154655638489i, -0.74338902742851332+1.7144305494668579i, +24.999999999999996+0.0000000000000000i, +9.0881451119497054+14.583996495037482i, -10.072020975530780+19.857294917630224i, -3.3430579786446146-3.7023698708437225i, -2.5485469814848569-1.3580112870658756i, +2.4021913759205669-4.6943359843703290i
326+3.8586660822987309-30.114708684502581i, +1.6484892711513575+5.0079752464818243i, +9.0881451119497054-14.583996495037482i, +64.000000000000000+0.0000000000000000i, -16.077786058362733+26.500416474628455i, -0.93342013892648179+6.6910242635485693i, -5.8900767835004935-2.2303202911730722i, -2.4978749791202359-6.3382448954961106i
327-28.989864511941768+42.209195888284015i, +2.4492175802020539-6.1461612372403476i, -10.072020975530780-19.857294917630224i, -16.077786058362733-26.500416474628455i, +100.00000000000000+0.0000000000000000i, +8.9119508309388333+12.899656157225131i, +2.8269335531165192+3.0969899885391161i, -1.8908095855220450-1.8916326785560891i
328-1.1378722493831885+4.0838402059374586i, -0.83266894114766343-1.3453860400162176i, -3.3430579786446146+3.7023698708437225i, -0.93342013892648179-6.6910242635485693i, +8.9119508309388333-12.899656157225131i, +24.999999999999996+0.0000000000000000i, -1.0638876925864840+2.2761744335602927i, -7.8049825794282448-8.3068224205299224i
329+2.8740306224584269+6.1151019226214238i, -0.58692732470136544-0.53347264646785209i, -2.5485469814848569+1.3580112870658756i, -5.8900767835004935+2.2303202911730722i, +2.8269335531165192-3.0969899885391161i, -1.0638876925864840-2.2761744335602927i, +4.0000000000000000+0.0000000000000000i, -0.12198309460948249-0.63703794298224503i
330+3.7497833200578716+13.681979753815369i, -0.66377140200104745-1.8121418261847968i, +2.4021913759205669+4.6943359843703290i, -2.4978749791202359+6.3382448954961106i, -1.8908095855220450+1.8916326785560891i, -7.8049825794282448+8.3068224205299224i, -0.12198309460948249+0.63703794298224503i, +25.000000000000004+0.0000000000000000i
332T
334T
335getMatDet(rand)
336+1086057.1147289227+0.17113052308559418E-7i
337
338
339ndim = getUnifRand(3, 9)
340ndim
341+9
342call setResized(rand, [ndim, ndim])
343scale = getUnifRand(1, 10, ndim)
344scale
345+6.0000000000000000, +5.0000000000000000, +9.0000000000000000, +6.0000000000000000, +10.000000000000000, +4.0000000000000000, +9.0000000000000000, +2.0000000000000000, +1.0000000000000000
346
347call setCovRand(rngf, rand)
348rand
349+1.0000000000000000+0.0000000000000000i, +0.32900590699352944-0.60215539007590091i, -0.34878222193818620+0.57382377724911748i, +0.16293905335421063-0.65573657212641386E-2i, -0.14541883389026303-0.48321938761512639i, -0.47500268105014903-0.18021057147627817i, -0.37803009397803139-0.40150551543415125i, +0.18789855007849834-0.26353962103767120i, +0.24189566659532263+0.14379234458328108i
350+0.32900590699352944+0.60215539007590091i, +1.0000000000000000+0.0000000000000000i, -0.51170406346552189-0.21284472942399368i, -0.32265236111794438E-1+0.26336784042009720i, +0.37617455313999915-0.62275982071497066E-1i, -0.25626420579401843-0.31428706749708324i, +0.14305014065680322-0.53174028392126993i, +0.49360796123049644+0.28407613520212727i, -0.18960586081750747+0.19854237796870111i
351-0.34878222193818620-0.57382377724911748i, -0.51170406346552189+0.21284472942399368i, +1.0000000000000000+0.0000000000000000i, -0.13902369578998863-0.47763544720401130i, -0.56687492425499486+0.31223469307791435i, -0.24538379829727000+0.14697527918102854i, +0.16674605509478441+0.39617620500194328i, -0.47170141579206670+0.41802140358943321E-1i, +0.92931037383273060E-1-0.16908208691741061i
352+0.16293905335421063+0.65573657212641386E-2i, -0.32265236111794438E-1-0.26336784042009720i, -0.13902369578998863+0.47763544720401130i, +1.0000000000000000+0.0000000000000000i, -0.58837852051103624E-1-0.59695899004357433i, +0.33380470567800824-0.15919370972745603i, -0.65821374869493784E-1-0.28741285216199147E-1i, -0.13073184136097038-0.24863831264263264i, -0.25491833717523438-0.15626919859592545i
353-0.14541883389026303+0.48321938761512639i, +0.37617455313999915+0.62275982071497066E-1i, -0.56687492425499486-0.31223469307791435i, -0.58837852051103624E-1+0.59695899004357433i, +1.0000000000000000+0.0000000000000000i, +0.40323704722156989E-1-0.14685274532269429i, +0.33779508856120566+0.47965506313200929E-1i, +0.40255567534359560-0.14571434722310361i, +0.18707174670290055+0.12702426507958525i
354-0.47500268105014903+0.18021057147627817i, -0.25626420579401843+0.31428706749708324i, -0.24538379829727000-0.14697527918102854i, +0.33380470567800824+0.15919370972745603i, +0.40323704722156989E-1+0.14685274532269429i, +1.0000000000000000+0.0000000000000000i, -0.40377668535180232E-1+0.30879236302780994i, +0.79268772662593204E-1+0.17606108947028576i, -0.42154880556331853-0.77186227855060047E-1i
355-0.37803009397803139+0.40150551543415125i, +0.14305014065680322+0.53174028392126993i, +0.16674605509478441-0.39617620500194328i, -0.65821374869493784E-1+0.28741285216199147E-1i, +0.33779508856120566-0.47965506313200929E-1i, -0.40377668535180232E-1-0.30879236302780994i, +1.0000000000000000+0.0000000000000000i, -0.36513437792295343E-1-0.15876156076685616i, +0.18942734694968388-0.25970121603848861i
356+0.18789855007849834+0.26353962103767120i, +0.49360796123049644-0.28407613520212727i, -0.47170141579206670-0.41802140358943321E-1i, -0.13073184136097038+0.24863831264263264i, +0.40255567534359560+0.14571434722310361i, +0.79268772662593204E-1-0.17606108947028576i, -0.36513437792295343E-1+0.15876156076685616i, +0.99999999999999989+0.0000000000000000i, +0.24380083098597743+0.49039827326152130i
357+0.24189566659532263-0.14379234458328108i, -0.18960586081750747-0.19854237796870111i, +0.92931037383273060E-1+0.16908208691741061i, -0.25491833717523438+0.15626919859592545i, +0.18707174670290055-0.12702426507958525i, -0.42154880556331853+0.77186227855060047E-1i, +0.18942734694968388+0.25970121603848861i, +0.24380083098597743-0.49039827326152130i, +0.99999999999999989+0.0000000000000000i
359T
361T
362getMatDet(rand)
363+0.38269802321705537E-5-0.17999450129153882E-19i
364call setCovRand(rngf, rand, scale(1))
365rand
366+36.000000000000000+0.0000000000000000i, +19.363829600610593+17.376844187893692i, +17.206759786725925+7.4597724196007045i, -17.821842480593190-17.494392823175971i, -0.77625016634315835-11.567041412422100i, +10.401533332524508+18.312324578761611i, -5.3978274985189421-9.5615459830395366i, -15.361305345619648+2.0940482816834938i, +0.49798982146565740+6.4597558567529214i
367+19.363829600610593-17.376844187893692i, +36.000000000000000+0.0000000000000000i, +26.523070245720216+8.5029359723746776i, -26.194663018313044-7.1927095942124764i, -15.423227910174313-15.206547553577025i, +15.295731624149660+5.7558519944592534i, -4.9222984200229014+5.3181319156497269i, -13.501024493142154+8.9278351933838351i, -1.8720409568027889+1.7302137445461248i
368+17.206759786725925-7.4597724196007045i, +26.523070245720216-8.5029359723746776i, +36.000000000000007+0.0000000000000000i, -25.505262998954109+0.48042367214056847i, -13.374683512559727-5.1771658445461810i, +15.779197784070799+2.1009045328073728i, +9.3878434491285567+1.7676819804291128i, -8.4210639671366874+5.4643958829855164i, -0.12814696131537628E-1+3.3831916342645290i
369-17.821842480593190+17.494392823175971i, -26.194663018313044+7.1927095942124764i, -25.505262998954109-0.48042367214056847i, +35.999999999999993+0.0000000000000000i, +19.