Line data Source code
1 : /*
2 : * File: TestLDL.cpp
3 : * Author: Pantelis Sopasakis
4 : *
5 : * Created on Aug 4, 2015, 8:18:33 PM
6 : */
7 :
8 : #include "TestLDL.h"
9 : #include <cmath>
10 :
11 :
12 1 : CPPUNIT_TEST_SUITE_REGISTRATION(TestLDL);
13 :
14 4 : TestLDL::TestLDL() {
15 4 : }
16 :
17 8 : TestLDL::~TestLDL() {
18 8 : }
19 :
20 4 : void TestLDL::setUp() {
21 4 : }
22 :
23 4 : void TestLDL::tearDown() {
24 4 : }
25 :
26 1 : void TestLDL::testSolveDense() {
27 1 : const size_t n = 40;
28 1 : const double tol = 1e-7;
29 1 : const size_t repetitions = 20;
30 :
31 21 : for (size_t k = 0; k < repetitions; k++) {
32 : FactoredSolver *ldlSolver;
33 20 : Matrix A = MatrixFactory::MakeRandomMatrix(n, n, 0.0, 1.0, Matrix::MATRIX_DENSE);
34 40 : Matrix At(A);
35 20 : At.transpose();
36 20 : A += At; // this will ensure that A is symmetric
37 20 : A *= 0.5;
38 :
39 820 : for (size_t i = 0; i < n; i++) {
40 32800 : for (size_t j = 0; j < n; j++) {
41 32000 : _ASSERT_NUM_EQ(A.get(i, j), A.get(j, i), tol);
42 : }
43 : }
44 :
45 20 : ldlSolver = new LDLFactorization(A);
46 20 : _ASSERT_EQ(ForBESUtils::STATUS_OK, ldlSolver->factorize());
47 40 : Matrix b = MatrixFactory::MakeRandomMatrix(n, 1, 0.0, 1.0, Matrix::MATRIX_DENSE);
48 40 : Matrix x;
49 20 : _ASSERT_EQ(ForBESUtils::STATUS_OK, ldlSolver->solve(b, x));
50 :
51 40 : Matrix err = A * x - b;
52 :
53 820 : for (size_t i = 0; i < n; i++) {
54 800 : _ASSERT(std::abs(err.get(i, 0)) < tol);
55 : }
56 :
57 20 : delete ldlSolver;
58 20 : }
59 1 : }
60 :
61 1 : void TestLDL::testSolveSymmetric() {
62 1 : const size_t n = 10;
63 1 : const double tol = 1e-7;
64 1 : const size_t repetitions = 15;
65 1 : Matrix b = MatrixFactory::MakeRandomMatrix(n, 1, 0.0, 1.0, Matrix::MATRIX_DENSE);
66 2 : Matrix x;
67 16 : for (size_t k = 0; k < repetitions; k++) {
68 : FactoredSolver *ldlSolver;
69 15 : Matrix S = MatrixFactory::MakeRandomMatrix(n, n, 0.0, 1.0, Matrix::MATRIX_SYMMETRIC);
70 30 : Matrix S_copy(S);
71 15 : ldlSolver = new LDLFactorization(S);
72 15 : _ASSERT_EQ(ForBESUtils::STATUS_OK, ldlSolver->factorize());
73 : /*
74 : * MODIFY S to make sure it doesn't affect the solver...
75 : */
76 15 : S = MatrixFactory::MakeRandomMatrix(n, n, 0.0, 100.0, Matrix::MATRIX_SYMMETRIC);
77 15 : _ASSERT_EQ(ForBESUtils::STATUS_OK, ldlSolver->solve(b, x));
78 30 : Matrix err = S_copy * x - b; /* Make sure that the original S is used */
79 165 : for (size_t i = 0; i < n; i++) {
80 150 : _ASSERT(std::abs(err.get(i, 0)) < tol);
81 : }
82 15 : delete ldlSolver;
83 16 : }
84 1 : }
85 :
86 1 : void TestLDL::testSolveSparse() {
87 1 : const double tol = 1e-7;
88 1 : const size_t N = 10;
89 1 : const size_t NNZ = 19;
90 :
91 1 : Matrix A = MatrixFactory::MakeSparse(N, N, NNZ, Matrix::SPARSE_SYMMETRIC_L);
92 1 : A.set(0, 0, 1.7);
93 1 : A.set(0, 8, 0.13);
94 1 : A.set(1, 1, 1.0);
95 1 : A.set(1, 9, 0.01);
96 1 : A.set(2, 2, 1.5);
97 1 : A.set(3, 3, 1.1);
98 1 : A.set(4, 4, 2.6);
99 1 : A.set(5, 5, 1.2);
100 1 : A.set(6, 6, 1.3);
101 1 : A.set(6, 6, 1.3);
102 1 : A.set(7, 7, 1.6);
103 1 : A.set(8, 8, 1.4);
104 1 : A.set(9, 9, 3.1);
105 1 : A.set(1, 4, 0.02);
106 1 : A.set(4, 6, 0.16);
107 1 : A.set(4, 7, 0.09);
108 1 : A.set(4, 8, 0.52);
109 1 : A.set(4, 9, 0.53);
110 1 : A.set(6, 9, 0.56);
111 1 : A.set(7, 8, 0.11);
112 :
113 2 : Matrix A_copy(A);
114 :
115 1 : LDLFactorization * solver = new LDLFactorization(A);
116 : int status;
117 1 : _ASSERT_OK(status = solver->factorize());
118 :
119 :
120 1 : _ASSERT_EQ(ForBESUtils::STATUS_OK, status);
121 :
122 1 : double b[N] = {.287, .22, .45, .44, 2.486, .72, 1.55, 1.424, 1.621, 3.759};
123 2 : Matrix rhs(N, 1, b);
124 :
125 :
126 : /*
127 : * MODIFY A to make sure it doesn't affect the solver...
128 : */
129 1 : A = Matrix(N - 1, N - 1);
130 2 : Matrix sol;
131 1 : status = solver->solve(rhs, sol);
132 1 : _ASSERT_EQ(ForBESUtils::STATUS_OK, status);
133 :
134 :
135 11 : for (size_t i = 0; i < N; i++)
136 10 : _ASSERT_NUM_EQ((i + 1.0) / 10.0, sol.get(i, 0), tol);
137 :
138 2 : Matrix err = A_copy * sol - rhs;
139 11 : for (size_t i = 0; i < N; i++)
140 10 : _ASSERT(std::abs(err.get(i, 0)) < tol);
141 :
142 2 : delete solver;
143 1 : }
144 :
145 1 : void TestLDL::testSolveSparse2() {
146 1 : const double tol = 1e-7;
147 1 : const size_t N = 40;
148 1 : const size_t NNZ = 2 * N;
149 1 : Matrix A = MatrixFactory::MakeRandomSparse(N, N, NNZ, 0.0, 1.0);
150 41 : for (size_t i = 0; i < N; i++) {
151 40 : A.set(i, i, A.get(i, i) + 1.5 * N);
152 : }
153 2 : Matrix At(A);
154 1 : At.transpose();
155 1 : A += At;
156 :
157 2 : Matrix A_copy(A);
158 :
159 1 : FactoredSolver * solver = new LDLFactorization(A);
160 : int status;
161 1 : status = solver->factorize();
162 1 : _ASSERT_EQ(ForBESUtils::STATUS_OK, status);
163 :
164 2 : Matrix b = MatrixFactory::MakeRandomMatrix(N, 1, 0.0, 1.0, Matrix::MATRIX_DENSE);
165 2 : Matrix c = MatrixFactory::MakeRandomMatrix(N, 1, 0.0, 1.0, Matrix::MATRIX_DENSE);
166 :
167 2 : Matrix sol_b(N, 1);
168 2 : Matrix sol_c(N, 1);
169 :
170 1 : A = Matrix(1, 1);
171 1 : solver->solve(b, sol_b);
172 1 : solver->solve(c, sol_c);
173 :
174 2 : Matrix err_b = A_copy * sol_b - b;
175 2 : Matrix err_c = A_copy * sol_c - c;
176 41 : for (size_t i = 0; i < N; i++) {
177 40 : _ASSERT(std::abs(err_b.get(i, 0)) < tol);
178 40 : _ASSERT(std::abs(err_c.get(i, 0)) < tol);
179 : }
180 :
181 2 : delete solver;
182 4 : }
183 :
|
