1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103 | with Ada.Text_IO; use Ada.Text_IO;
with Ada.Numerics.Generic_Real_Arrays;
with U_Rand;
with Copulas;
procedure Test_Copulas is
subtype Real is Long_Float;
type Integer_Vector is array (Integer range <>) of Integer;
package RIO is new Float_IO (Real);
package GRA is new Ada.Numerics.Generic_Real_Arrays (Real);
package Real_U_Rand is new U_Rand (Real);
-- *** Choice of a random generator: A.N.F_R, or U_Rand (faster), or...:
package RRand renames
-- Ada.Numerics.Float_Random;
Real_U_Rand;
package RCopulas is new Copulas (
Real,
RRand.Uniformly_Distributed,
RRand.Generator,
RRand.Random,
GRA,
Integer_Vector
);
use GRA, RIO;
procedure Put (A : Real_Matrix) is
begin
for i in A'Range(1) loop
for j in A'Range(2) loop
Put (A(i,j), 4, Real'Digits, 0);
end loop;
New_Line;
end loop;
end Put;
procedure Test (A : Real_Matrix; res_name: String) is
dim_dep : constant Integer := A'Length(1);
dim_indep : constant := 3;
U01 : Real_Vector (1..dim_dep + dim_indep);
gen : RRand.Generator;
use RCopulas;
g_copula : Copula_access;
f : File_Type;
L, B, Z : Real_Matrix (A'Range(1), A'Range(2));
res : Real:= 0.0;
begin
RRand.Reset (gen, 1);
Construct_as_Gauss (g_copula, U01'Length, A);
Create (f, Out_File, res_name & ".csv");
for i in 1 .. dim_dep loop
Put (f, "dep;");
end loop;
for i in 1 .. dim_indep loop
Put (f, "indep;");
end loop;
New_Line(f);
for n in 1 .. 60_000 loop
U01 := Simulate (g_copula.all, gen);
for i in U01'Range loop
Put (f, U01(i));
Put (f, ';');
end loop;
New_Line(f);
end loop;
Put_Line("Results are stored in this file: " & Name(f));
L := Get_Cholesky_Matrix (Gauss_Copula (g_copula.all));
B := L * Transpose(L);
Put_Line("B := L*Lt.");
Put_Line("""A = B"" Matrix equality test (bitwise) " & Boolean'Image(A = B));
Z := A - B;
for i in Z'Range(1) loop
for j in Z'Range(2) loop
res:= res + abs Z(i,j);
end loop;
end loop;
Put_Line("""A - B = 0"" Matrix equality test (L1) " & Real'Image(res));
Put_Line("A =");
Put(A);
New_Line;
Put_Line("B =");
Put(B);
New_Line;
Put_Line("L =");
Put(L);
end Test;
A33: constant Real_Matrix(1..3, 1..3):=
((1.0 , 0.5 , 0.25),
(0.5 , 1.0 , 0.06),
(0.25, 0.06, 1.0));
begin
Test( A33, "Gaussian_3x3" );
-- Skip_Line;
end Test_Copulas;
|