Actual source code: getcolv.c
2: #include <petsc/private/matimpl.h>
4: /*@
5: MatGetColumnVector - Gets the values from a given column of a matrix.
7: Not Collective
9: Input Parameters:
10: + A - the matrix
11: . yy - the vector
12: - col - the column requested (in global numbering)
14: Level: advanced
16: Notes:
17: If a Mat type does not implement the operation, each processor for which this is called
18: gets the values for its rows using `MatGetRow()`.
20: The vector must have the same parallel row layout as the matrix.
22: Contributed by: Denis Vanderstraeten
24: .seealso: `MatGetRow()`, `MatGetDiagonal()`, `MatMult()`
25: @*/
26: PetscErrorCode MatGetColumnVector(Mat A, Vec yy, PetscInt col)
27: {
28: PetscScalar *y;
29: const PetscScalar *v;
30: PetscInt i, j, nz, N, Rs, Re, rs, re;
31: const PetscInt *idx;
37: MatGetSize(A, NULL, &N);
39: MatGetOwnershipRange(A, &Rs, &Re);
40: VecGetOwnershipRange(yy, &rs, &re);
43: if (A->ops->getcolumnvector) PetscUseTypeMethod(A, getcolumnvector, yy, col);
44: else {
45: VecSet(yy, 0.0);
46: VecGetArray(yy, &y);
47: /* TODO for general matrices */
48: for (i = Rs; i < Re; i++) {
49: MatGetRow(A, i, &nz, &idx, &v);
50: if (nz && idx[0] <= col) {
51: /*
52: Should use faster search here
53: */
54: for (j = 0; j < nz; j++) {
55: if (idx[j] >= col) {
56: if (idx[j] == col) y[i - rs] = v[j];
57: break;
58: }
59: }
60: }
61: MatRestoreRow(A, i, &nz, &idx, &v);
62: }
63: VecRestoreArray(yy, &y);
64: }
65: return 0;
66: }
68: /*@
69: MatGetColumnNorms - Gets the norms of each column of a sparse or dense matrix.
71: Input Parameters:
72: + A - the matrix
73: - type - `NORM_2`, `NORM_1` or `NORM_INFINITY`
75: Output Parameter:
76: . norms - an array as large as the TOTAL number of columns in the matrix
78: Level: intermediate
80: Note:
81: Each process has ALL the column norms after the call. Because of the way this is computed each process gets all the values,
82: if each process wants only some of the values it should extract the ones it wants from the array.
84: .seealso: `NormType`, `MatNorm()`
85: @*/
86: PetscErrorCode MatGetColumnNorms(Mat A, NormType type, PetscReal norms[])
87: {
88: /* NOTE: MatGetColumnNorms() could simply be a macro that calls MatGetColumnReductions().
89: * I've kept this as a function because it allows slightly more in the way of error checking,
90: * erroring out if MatGetColumnNorms() is not called with a valid NormType. */
92: if (type == NORM_2 || type == NORM_1 || type == NORM_FROBENIUS || type == NORM_INFINITY || type == NORM_1_AND_2) {
93: MatGetColumnReductions(A, type, norms);
94: } else SETERRQ(PetscObjectComm((PetscObject)A), PETSC_ERR_ARG_WRONG, "Unknown NormType");
95: return 0;
96: }
98: /*@
99: MatGetColumnSumsRealPart - Gets the sums of the real part of each column of a sparse or dense matrix.
101: Input Parameter:
102: . A - the matrix
104: Output Parameter:
105: . sums - an array as large as the TOTAL number of columns in the matrix
107: Level: intermediate
109: Note:
110: Each process has ALL the column sums after the call. Because of the way this is computed each process gets all the values,
111: if each process wants only some of the values it should extract the ones it wants from the array.
113: .seealso: `MatGetColumnSumsImaginaryPart()`, `VecSum()`, `MatGetColumnMeans()`, `MatGetColumnNorms()`, `MatGetColumnReductions()`
114: @*/
115: PetscErrorCode MatGetColumnSumsRealPart(Mat A, PetscReal sums[])
116: {
117: MatGetColumnReductions(A, REDUCTION_SUM_REALPART, sums);
118: return 0;
119: }
121: /*@
122: MatGetColumnSumsImaginaryPart - Gets the sums of the imaginary part of each column of a sparse or dense matrix.
124: Input Parameter:
125: . A - the matrix
127: Output Parameter:
128: . sums - an array as large as the TOTAL number of columns in the matrix
130: Level: intermediate
132: Note:
133: Each process has ALL the column sums after the call. Because of the way this is computed each process gets all the values,
134: if each process wants only some of the values it should extract the ones it wants from the array.
136: .seealso: `MatGetColumnSumsRealPart()`, `VecSum()`, `MatGetColumnMeans()`, `MatGetColumnNorms()`, `MatGetColumnReductions()`
137: @*/
138: PetscErrorCode MatGetColumnSumsImaginaryPart(Mat A, PetscReal sums[])
139: {
140: MatGetColumnReductions(A, REDUCTION_SUM_IMAGINARYPART, sums);
141: return 0;
142: }
144: /*@
145: MatGetColumnSums - Gets the sums of each column of a sparse or dense matrix.
147: Input Parameter:
148: . A - the matrix
150: Output Parameter:
151: . sums - an array as large as the TOTAL number of columns in the matrix
153: Level: intermediate
155: Note:
156: Each process has ALL the column sums after the call. Because of the way this is computed each process gets all the values,
157: if each process wants only some of the values it should extract the ones it wants from the array.
159: .seealso: `VecSum()`, `MatGetColumnMeans()`, `MatGetColumnNorms()`, `MatGetColumnReductions()`
160: @*/
161: PetscErrorCode MatGetColumnSums(Mat A, PetscScalar sums[])
162: {
163: #if defined(PETSC_USE_COMPLEX)
164: PetscInt i, n;
165: PetscReal *work;
166: #endif
169: #if !defined(PETSC_USE_COMPLEX)
170: MatGetColumnSumsRealPart(A, sums);
171: #else
172: MatGetSize(A, NULL, &n);
173: PetscArrayzero(sums, n);
174: PetscCalloc1(n, &work);
175: MatGetColumnSumsRealPart(A, work);
176: for (i = 0; i < n; i++) sums[i] = work[i];
177: MatGetColumnSumsImaginaryPart(A, work);
178: for (i = 0; i < n; i++) sums[i] += work[i] * PETSC_i;
179: PetscFree(work);
180: #endif
181: return 0;
182: }
184: /*@
185: MatGetColumnMeansRealPart - Gets the arithmetic means of the real part of each column of a sparse or dense matrix.
187: Input Parameter:
188: . A - the matrix
190: Output Parameter:
191: . sums - an array as large as the TOTAL number of columns in the matrix
193: Level: intermediate
195: Note:
196: Each process has ALL the column means after the call. Because of the way this is computed each process gets all the values,
197: if each process wants only some of the values it should extract the ones it wants from the array.
199: .seealso: `MatGetColumnMeansImaginaryPart()`, `VecSum()`, `MatGetColumnSums()`, `MatGetColumnNorms()`, `MatGetColumnReductions()`
200: @*/
201: PetscErrorCode MatGetColumnMeansRealPart(Mat A, PetscReal means[])
202: {
203: MatGetColumnReductions(A, REDUCTION_MEAN_REALPART, means);
204: return 0;
205: }
207: /*@
208: MatGetColumnMeansImaginaryPart - Gets the arithmetic means of the imaginary part of each column of a sparse or dense matrix.
210: Input Parameter:
211: . A - the matrix
213: Output Parameter:
214: . sums - an array as large as the TOTAL number of columns in the matrix
216: Level: intermediate
218: Note:
219: Each process has ALL the column means after the call. Because of the way this is computed each process gets all the values,
220: if each process wants only some of the values it should extract the ones it wants from the array.
222: .seealso: `MatGetColumnMeansRealPart()`, `VecSum()`, `MatGetColumnSums()`, `MatGetColumnNorms()`, `MatGetColumnReductions()`
223: @*/
224: PetscErrorCode MatGetColumnMeansImaginaryPart(Mat A, PetscReal means[])
225: {
226: MatGetColumnReductions(A, REDUCTION_MEAN_IMAGINARYPART, means);
227: return 0;
228: }
230: /*@
231: MatGetColumnMeans - Gets the arithmetic means of each column of a sparse or dense matrix.
233: Input Parameter:
234: . A - the matrix
236: Output Parameter:
237: . means - an array as large as the TOTAL number of columns in the matrix
239: Level: intermediate
241: Note:
242: Each process has ALL the column means after the call. Because of the way this is computed each process gets all the values,
243: if each process wants only some of the values it should extract the ones it wants from the array.
245: .seealso: `VecSum()`, `MatGetColumnSums()`, `MatGetColumnNorms()`, `MatGetColumnReductions()`
246: @*/
247: PetscErrorCode MatGetColumnMeans(Mat A, PetscScalar means[])
248: {
249: #if defined(PETSC_USE_COMPLEX)
250: PetscInt i, n;
251: PetscReal *work;
252: #endif
255: #if !defined(PETSC_USE_COMPLEX)
256: MatGetColumnMeansRealPart(A, means);
257: #else
258: MatGetSize(A, NULL, &n);
259: PetscArrayzero(means, n);
260: PetscCalloc1(n, &work);
261: MatGetColumnMeansRealPart(A, work);
262: for (i = 0; i < n; i++) means[i] = work[i];
263: MatGetColumnMeansImaginaryPart(A, work);
264: for (i = 0; i < n; i++) means[i] += work[i] * PETSC_i;
265: PetscFree(work);
266: #endif
267: return 0;
268: }
270: /*@
271: MatGetColumnReductions - Gets the reductions of each column of a sparse or dense matrix.
273: Input Parameters:
274: + A - the matrix
275: - type - A constant defined in `NormType` or `ReductionType`: `NORM_2`, `NORM_1`, `NORM_INFINITY`, `REDUCTION_SUM_REALPART`,
276: `REDUCTION_SUM_IMAGINARYPART`, `REDUCTION_MEAN_REALPART`, `REDUCTION_MEAN_IMAGINARYPART`
278: Output Parameter:
279: . reductions - an array as large as the TOTAL number of columns in the matrix
281: Level: developer
283: Note:
284: Each process has ALL the column reductions after the call. Because of the way this is computed each process gets all the values,
285: if each process wants only some of the values it should extract the ones it wants from the array.
287: Developer Note:
288: This routine is primarily intended as a back-end.
289: `MatGetColumnNorms()`, `MatGetColumnSums()`, and `MatGetColumnMeans()` are implemented using this routine.
291: .seealso: `ReductionType`, `NormType`, `MatGetColumnNorms()`, `MatGetColumnSums()`, `MatGetColumnMeans()`
292: @*/
293: PetscErrorCode MatGetColumnReductions(Mat A, PetscInt type, PetscReal reductions[])
294: {
296: PetscUseTypeMethod(A, getcolumnreductions, type, reductions);
297: return 0;
298: }