GRASS GIS 8 Programmer's Manual 8.3.2(2024)-exported
Loading...
Searching...
No Matches
symmetric_band_matrix.c
Go to the documentation of this file.
1#include <stdlib.h>
2#include <stdio.h>
3#include <math.h>
4#include <grass/gis.h>
5#include <grass/gmath.h>
6#include <grass/glocale.h>
7
8/**
9 * \brief Convert a symmetrix matrix into a symmetric band matrix
10 *
11 * \verbatim
12 Symmetric matrix with bandwidth of 3
13
14 5 2 1 0
15 2 5 2 1
16 1 2 5 2
17 0 1 2 5
18
19 will be converted into the symmetric band matrix
20
21 5 2 1
22 5 2 1
23 5 2 0
24 5 0 0
25
26 \endverbatim
27
28 * \param A (double**) the symmetric matrix
29 * \param rows (int)
30 * \param bandwidth (int)
31 * \return B (double**) new created symmetric band matrix
32 * */
33
34double **G_math_matrix_to_sband_matrix(double **A, int rows, int bandwidth)
35{
36 int i, j;
37 double **B = G_alloc_matrix(rows, bandwidth);
38 double tmp;
39
40 for (i = 0; i < rows; i++) {
41 for (j = 0; j < bandwidth; j++) {
42 if (i + j < rows) {
43 tmp = A[i][i + j];
44 B[i][j] = tmp;
45 }
46 else {
47 B[i][j] = 0.0;
48 }
49 }
50 }
51
52 return B;
53}
54
55/**
56 * \brief Convert a symmetric band matrix into a symmetric matrix
57 *
58 * \verbatim
59 Such a symmetric band matrix with banwidth 3
60
61 5 2 1
62 5 2 1
63 5 2 0
64 5 0 0
65
66 Will be converted into this symmetric matrix
67
68 5 2 1 0
69 2 5 2 1
70 1 2 5 2
71 0 1 2 5
72
73 \endverbatim
74 * \param A (double**) the symmetric band matrix
75 * \param rows (int)
76 * \param bandwidth (int)
77 * \return B (double**) new created symmetric matrix
78 * */
79
80double **G_math_sband_matrix_to_matrix(double **A, int rows, int bandwidth)
81{
82 int i, j;
83 double **B = G_alloc_matrix(rows, rows);
84 double tmp;
85
86 for (i = 0; i < rows; i++) {
87 for (j = 0; j < bandwidth; j++) {
88 tmp = A[i][j];
89 if (i + j < rows) {
90 B[i][i + j] = tmp;
91 }
92 }
93 }
94
95 /*Symmetry */
96 for (i = 0; i < rows; i++) {
97 for (j = i; j < rows; j++) {
98 B[j][i] = B[i][j];
99 }
100 }
101
102 return B;
103}
104
105/*!
106 * \brief Compute the matrix - vector product
107 * of symmetric band matrix A and vector x.
108 *
109 * This function is multi-threaded with OpenMP and can be called within a
110 * parallel OpenMP region.
111 *
112 * y = A * x
113 *
114 *
115 * \param A (double **)
116 * \param x (double) *)
117 * \param y (double * )
118 * \param rows (int)
119 * \param bandwidth (int)
120 * \return (void)
121 *
122 * */
123void G_math_Ax_sband(double **A, double *x, double *y, int rows, int bandwidth)
124{
125 int i, j;
126 double tmp;
127
128#pragma omp for schedule(static) private(i, j, tmp)
129 for (i = 0; i < rows; i++) {
130 tmp = 0;
131 for (j = 0; j < bandwidth; j++) {
132 if ((i + j) < rows)
133 tmp += A[i][j] * x[i + j];
134 }
135 y[i] = tmp;
136 }
137
138#pragma omp single
139 {
140 for (i = 0; i < rows; i++) {
141 tmp = 0;
142 for (j = 1; j < bandwidth; j++) {
143 if (i + j < rows)
144 y[i + j] += A[i][j] * x[i];
145 }
146 }
147 }
148 return;
149}
double ** G_alloc_matrix(int rows, int cols)
Matrix memory allocation.
Definition dalloc.c:57
void G_math_Ax_sband(double **A, double *x, double *y, int rows, int bandwidth)
Compute the matrix - vector product of symmetric band matrix A and vector x.
double ** G_math_matrix_to_sband_matrix(double **A, int rows, int bandwidth)
Convert a symmetrix matrix into a symmetric band matrix.
double ** G_math_sband_matrix_to_matrix(double **A, int rows, int bandwidth)
Convert a symmetric band matrix into a symmetric matrix.
#define x