-
-
Notifications
You must be signed in to change notification settings - Fork 302
/
bench_solver_krylov.c
122 lines (91 loc) · 4.08 KB
/
bench_solver_krylov.c
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
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
/*****************************************************************************
*
* MODULE: Grass PDE Numerical Library
* AUTHOR(S): Soeren Gebbert, Berlin (GER) Dec 2006
* soerengebbert <at> gmx <dot> de
*
* PURPOSE: benchmarking the krylov subspace solvers
*
* COPYRIGHT: (C) 2000 by the GRASS Development Team
*
* This program is free software under the GNU General Public
* License (>=v2). Read the file COPYING that comes with GRASS
* for details.
*
*****************************************************************************/
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <grass/glocale.h>
#include <grass/gmath.h>
#include "test_gmath_lib.h"
#include <sys/time.h>
/* prototypes */
static int bench_solvers(int rows);
/* ************************************************************************* */
/* Performe the solver unit tests ****************************************** */
/* ************************************************************************* */
int bench_solvers_krylov(int rows) {
G_message(_("\n++ Running krylov solver benchmark ++"));
bench_solvers(rows);
return 1;
}
/* *************************************************************** */
/* Test all implemented solvers for sparse and normal matrix *** */
/* *************************************************************** */
int bench_solvers(int rows) {
G_math_les *les;
G_math_les *sples;
struct timeval tstart;
struct timeval tend;
G_message("\t * benchmarking pcg solver with symmetric matrix and preconditioner 1\n");
les = create_normal_symmetric_les(rows);
sples = create_sparse_symmetric_les(rows);
gettimeofday(&tstart, NULL);
G_math_solver_pcg(les->A, les->x, les->b, les->rows, 250, 0.1e-9, 1);
gettimeofday(&tend, NULL);
printf("Computation time pcg normal matrix: %g\n", compute_time_difference(tstart, tend));
gettimeofday(&tstart, NULL);
G_math_solver_sparse_pcg(sples->Asp, sples->x, sples->b, les->rows, 250,
0.1e-9, 1);
gettimeofday(&tend, NULL);
printf("Computation time pcg sparse matrix: %g\n", compute_time_difference(tstart, tend));
G_math_free_les(les);
G_math_free_les(sples);
G_message("\t * benchmark cg solver with symmetric matrix\n");
les = create_normal_symmetric_les(rows);
sples = create_sparse_symmetric_les(rows);
gettimeofday(&tstart, NULL);
G_math_solver_cg(les->A, les->x, les->b, les->rows, 250, 0.1e-9);
gettimeofday(&tend, NULL);
printf("Computation time cg normal matrix: %g\n", compute_time_difference(tstart, tend));
gettimeofday(&tstart, NULL);
G_math_solver_sparse_cg(sples->Asp, sples->x, sples->b, les->rows, 250,
0.1e-9);
gettimeofday(&tend, NULL);
printf("Computation time cg sparse matrix: %g\n", compute_time_difference(tstart, tend));
G_math_free_les(les);
G_math_free_les(sples);
G_message("\t * benchmark cg solver with symmetric band matrix\n");
les = create_symmetric_band_les(rows);
gettimeofday(&tstart, NULL);
G_math_solver_cg_sband(les->A, les->x, les->b, les->rows, les->rows, 250, 0.1e-9);
gettimeofday(&tend, NULL);
printf("Computation time cg symmetric band matrix: %g\n", compute_time_difference(tstart, tend));
G_math_free_les(les);
G_message("\t * benchmark bicgstab solver with unsymmetric matrix\n");
les = create_normal_unsymmetric_les(rows);
sples = create_sparse_unsymmetric_les(rows);
gettimeofday(&tstart, NULL);
G_math_solver_bicgstab(les->A, les->x, les->b, les->rows, 250, 0.1e-9);
gettimeofday(&tend, NULL);
printf("Computation time bicgstab normal matrix: %g\n", compute_time_difference(tstart, tend));
gettimeofday(&tstart, NULL);
G_math_solver_sparse_bicgstab(sples->Asp, sples->x, sples->b, les->rows,
250, 0.1e-9);
gettimeofday(&tend, NULL);
printf("Computation time bicgstab sparse matrix: %g\n", compute_time_difference(tstart, tend));
G_math_free_les(les);
G_math_free_les(sples);
return 1;
}