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 `MatType` 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: [](ch_matrices), `Mat`, `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;
33: PetscFunctionBegin;
37: PetscCheck(col >= 0, PetscObjectComm((PetscObject)A), PETSC_ERR_ARG_OUTOFRANGE, "Requested negative column: %" PetscInt_FMT, col);
38: PetscCall(MatGetSize(A, NULL, &N));
39: PetscCheck(col < N, PetscObjectComm((PetscObject)A), PETSC_ERR_ARG_OUTOFRANGE, "Requested column %" PetscInt_FMT " larger than number columns in matrix %" PetscInt_FMT, col, N);
40: PetscCall(MatGetOwnershipRange(A, &Rs, &Re));
41: PetscCall(VecGetOwnershipRange(yy, &rs, &re));
42: PetscCheck(Rs == rs && Re == re, PETSC_COMM_SELF, PETSC_ERR_ARG_INCOMP, "Matrix %" PetscInt_FMT " %" PetscInt_FMT " does not have same ownership range (size) as vector %" PetscInt_FMT " %" PetscInt_FMT, Rs, Re, rs, re);
44: if (A->ops->getcolumnvector) PetscUseTypeMethod(A, getcolumnvector, yy, col);
45: else {
46: PetscCall(VecSet(yy, 0.0));
47: PetscCall(VecGetArray(yy, &y));
48: /* TODO for general matrices */
49: for (i = Rs; i < Re; i++) {
50: PetscCall(MatGetRow(A, i, &nz, &idx, &v));
51: if (nz && idx[0] <= col) {
52: /*
53: Should use faster search here
54: */
55: for (j = 0; j < nz; j++) {
56: if (idx[j] >= col) {
57: if (idx[j] == col) y[i - rs] = v[j];
58: break;
59: }
60: }
61: }
62: PetscCall(MatRestoreRow(A, i, &nz, &idx, &v));
63: }
64: PetscCall(VecRestoreArray(yy, &y));
65: }
66: PetscFunctionReturn(PETSC_SUCCESS);
67: }
69: /*@
70: MatGetColumnNorms - Gets the norms of each column of a sparse or dense matrix.
72: Input Parameters:
73: + A - the matrix
74: - type - `NORM_2`, `NORM_1` or `NORM_INFINITY`
76: Output Parameter:
77: . norms - an array as large as the TOTAL number of columns in the matrix
79: Level: intermediate
81: Note:
82: Each process has ALL the column norms after the call. Because of the way this is computed each process gets all the values,
83: if each process wants only some of the values it should extract the ones it wants from the array.
85: .seealso: [](ch_matrices), `Mat`, `NormType`, `MatNorm()`
86: @*/
87: PetscErrorCode MatGetColumnNorms(Mat A, NormType type, PetscReal norms[])
88: {
89: /* NOTE: MatGetColumnNorms() could simply be a macro that calls MatGetColumnReductions().
90: * I've kept this as a function because it allows slightly more in the way of error checking,
91: * erroring out if MatGetColumnNorms() is not called with a valid NormType. */
93: PetscFunctionBegin;
94: if (type == NORM_2 || type == NORM_1 || type == NORM_FROBENIUS || type == NORM_INFINITY || type == NORM_1_AND_2) {
95: PetscCall(MatGetColumnReductions(A, type, norms));
96: } else SETERRQ(PetscObjectComm((PetscObject)A), PETSC_ERR_ARG_WRONG, "Unknown NormType");
97: PetscFunctionReturn(PETSC_SUCCESS);
98: }
100: /*@
101: MatGetColumnSumsRealPart - Gets the sums of the real part of each column of a sparse or dense matrix.
103: Input Parameter:
104: . A - the matrix
106: Output Parameter:
107: . sums - an array as large as the TOTAL number of columns in the matrix
109: Level: intermediate
111: Note:
112: Each process has ALL the column sums after the call. Because of the way this is computed each process gets all the values,
113: if each process wants only some of the values it should extract the ones it wants from the array.
115: .seealso: [](ch_matrices), `Mat`, `MatGetColumnSumsImaginaryPart()`, `VecSum()`, `MatGetColumnMeans()`, `MatGetColumnNorms()`, `MatGetColumnReductions()`
116: @*/
117: PetscErrorCode MatGetColumnSumsRealPart(Mat A, PetscReal sums[])
118: {
119: PetscFunctionBegin;
120: PetscCall(MatGetColumnReductions(A, REDUCTION_SUM_REALPART, sums));
121: PetscFunctionReturn(PETSC_SUCCESS);
122: }
124: /*@
125: MatGetColumnSumsImaginaryPart - Gets the sums of the imaginary part of each column of a sparse or dense matrix.
127: Input Parameter:
128: . A - the matrix
130: Output Parameter:
131: . sums - an array as large as the TOTAL number of columns in the matrix
133: Level: intermediate
135: Note:
136: Each process has ALL the column sums after the call. Because of the way this is computed each process gets all the values,
137: if each process wants only some of the values it should extract the ones it wants from the array.
139: .seealso: [](ch_matrices), `Mat`, `MatGetColumnSumsRealPart()`, `VecSum()`, `MatGetColumnMeans()`, `MatGetColumnNorms()`, `MatGetColumnReductions()`
140: @*/
141: PetscErrorCode MatGetColumnSumsImaginaryPart(Mat A, PetscReal sums[])
142: {
143: PetscFunctionBegin;
144: PetscCall(MatGetColumnReductions(A, REDUCTION_SUM_IMAGINARYPART, sums));
145: PetscFunctionReturn(PETSC_SUCCESS);
146: }
148: /*@
149: MatGetColumnSums - Gets the sums of each column of a sparse or dense matrix.
151: Input Parameter:
152: . A - the matrix
154: Output Parameter:
155: . sums - an array as large as the TOTAL number of columns in the matrix
157: Level: intermediate
159: Note:
160: Each process has ALL the column sums after the call. Because of the way this is computed each process gets all the values,
161: if each process wants only some of the values it should extract the ones it wants from the array.
163: .seealso: [](ch_matrices), `Mat`, `VecSum()`, `MatGetColumnMeans()`, `MatGetColumnNorms()`, `MatGetColumnReductions()`
164: @*/
165: PetscErrorCode MatGetColumnSums(Mat A, PetscScalar sums[])
166: {
167: #if defined(PETSC_USE_COMPLEX)
168: PetscInt i, n;
169: PetscReal *work;
170: #endif
172: PetscFunctionBegin;
174: #if !defined(PETSC_USE_COMPLEX)
175: PetscCall(MatGetColumnSumsRealPart(A, sums));
176: #else
177: PetscCall(MatGetSize(A, NULL, &n));
178: PetscCall(PetscArrayzero(sums, n));
179: PetscCall(PetscCalloc1(n, &work));
180: PetscCall(MatGetColumnSumsRealPart(A, work));
181: for (i = 0; i < n; i++) sums[i] = work[i];
182: PetscCall(MatGetColumnSumsImaginaryPart(A, work));
183: for (i = 0; i < n; i++) sums[i] += work[i] * PETSC_i;
184: PetscCall(PetscFree(work));
185: #endif
186: PetscFunctionReturn(PETSC_SUCCESS);
187: }
189: /*@
190: MatGetColumnMeansRealPart - Gets the arithmetic means of the real part of each column of a sparse or dense matrix.
192: Input Parameter:
193: . A - the matrix
195: Output Parameter:
196: . sums - an array as large as the TOTAL number of columns in the matrix
198: Level: intermediate
200: Note:
201: Each process has ALL the column means after the call. Because of the way this is computed each process gets all the values,
202: if each process wants only some of the values it should extract the ones it wants from the array.
204: .seealso: [](ch_matrices), `Mat`, `MatGetColumnMeansImaginaryPart()`, `VecSum()`, `MatGetColumnSums()`, `MatGetColumnNorms()`, `MatGetColumnReductions()`
205: @*/
206: PetscErrorCode MatGetColumnMeansRealPart(Mat A, PetscReal means[])
207: {
208: PetscFunctionBegin;
209: PetscCall(MatGetColumnReductions(A, REDUCTION_MEAN_REALPART, means));
210: PetscFunctionReturn(PETSC_SUCCESS);
211: }
213: /*@
214: MatGetColumnMeansImaginaryPart - Gets the arithmetic means of the imaginary part of each column of a sparse or dense matrix.
216: Input Parameter:
217: . A - the matrix
219: Output Parameter:
220: . sums - an array as large as the TOTAL number of columns in the matrix
222: Level: intermediate
224: Note:
225: Each process has ALL the column means after the call. Because of the way this is computed each process gets all the values,
226: if each process wants only some of the values it should extract the ones it wants from the array.
228: .seealso: [](ch_matrices), `Mat`, `MatGetColumnMeansRealPart()`, `VecSum()`, `MatGetColumnSums()`, `MatGetColumnNorms()`, `MatGetColumnReductions()`
229: @*/
230: PetscErrorCode MatGetColumnMeansImaginaryPart(Mat A, PetscReal means[])
231: {
232: PetscFunctionBegin;
233: PetscCall(MatGetColumnReductions(A, REDUCTION_MEAN_IMAGINARYPART, means));
234: PetscFunctionReturn(PETSC_SUCCESS);
235: }
237: /*@
238: MatGetColumnMeans - Gets the arithmetic means of each column of a sparse or dense matrix.
240: Input Parameter:
241: . A - the matrix
243: Output Parameter:
244: . means - an array as large as the TOTAL number of columns in the matrix
246: Level: intermediate
248: Note:
249: Each process has ALL the column means after the call. Because of the way this is computed each process gets all the values,
250: if each process wants only some of the values it should extract the ones it wants from the array.
252: .seealso: [](ch_matrices), `Mat`, `VecSum()`, `MatGetColumnSums()`, `MatGetColumnNorms()`, `MatGetColumnReductions()`
253: @*/
254: PetscErrorCode MatGetColumnMeans(Mat A, PetscScalar means[])
255: {
256: #if defined(PETSC_USE_COMPLEX)
257: PetscInt i, n;
258: PetscReal *work;
259: #endif
261: PetscFunctionBegin;
263: #if !defined(PETSC_USE_COMPLEX)
264: PetscCall(MatGetColumnMeansRealPart(A, means));
265: #else
266: PetscCall(MatGetSize(A, NULL, &n));
267: PetscCall(PetscArrayzero(means, n));
268: PetscCall(PetscCalloc1(n, &work));
269: PetscCall(MatGetColumnMeansRealPart(A, work));
270: for (i = 0; i < n; i++) means[i] = work[i];
271: PetscCall(MatGetColumnMeansImaginaryPart(A, work));
272: for (i = 0; i < n; i++) means[i] += work[i] * PETSC_i;
273: PetscCall(PetscFree(work));
274: #endif
275: PetscFunctionReturn(PETSC_SUCCESS);
276: }
278: /*@
279: MatGetColumnReductions - Gets the reductions of each column of a sparse or dense matrix.
281: Input Parameters:
282: + A - the matrix
283: - type - A constant defined in `NormType` or `ReductionType`: `NORM_2`, `NORM_1`, `NORM_INFINITY`, `REDUCTION_SUM_REALPART`,
284: `REDUCTION_SUM_IMAGINARYPART`, `REDUCTION_MEAN_REALPART`, `REDUCTION_MEAN_IMAGINARYPART`
286: Output Parameter:
287: . reductions - an array as large as the TOTAL number of columns in the matrix
289: Level: developer
291: Note:
292: Each process has ALL the column reductions after the call. Because of the way this is computed each process gets all the values,
293: if each process wants only some of the values it should extract the ones it wants from the array.
295: Developer Note:
296: This routine is primarily intended as a back-end.
297: `MatGetColumnNorms()`, `MatGetColumnSums()`, and `MatGetColumnMeans()` are implemented using this routine.
299: .seealso: [](ch_matrices), `Mat`, `ReductionType`, `NormType`, `MatGetColumnNorms()`, `MatGetColumnSums()`, `MatGetColumnMeans()`
300: @*/
301: PetscErrorCode MatGetColumnReductions(Mat A, PetscInt type, PetscReal reductions[])
302: {
303: PetscFunctionBegin;
305: PetscUseTypeMethod(A, getcolumnreductions, type, reductions);
306: PetscFunctionReturn(PETSC_SUCCESS);
307: }