Actual source code: getcolv.c
1: #define PETSCMAT_DLL
3: #include include/private/matimpl.h
7: /*@
8: MatGetColumnVector - Gets the values from a given column of a matrix.
10: Not Collective
12: Input Parameters:
13: + A - the matrix
14: . yy - the vector
15: - c - the column requested (in global numbering)
17: Level: advanced
19: Notes:
20: Each processor for which this is called gets the values for its rows.
22: Since PETSc matrices are usually stored in compressed row format, this routine
23: will generally be slow.
25: The vector must have the same parallel row layout as the matrix.
27: Contributed by: Denis Vanderstraeten
29: .keywords: matrix, column, get
31: .seealso: MatGetRow(), MatGetDiagonal()
33: @*/
34: PetscErrorCode MatGetColumnVector(Mat A,Vec yy,PetscInt col)
35: {
36: PetscScalar *y,zero = 0.0;
37: const PetscScalar *v;
38: PetscErrorCode ierr;
39: PetscInt i,j,nz,N,Rs,Re,rs,re;
40: const PetscInt *idx;
41: MPI_Comm comm;
42:
47: if (col < 0) SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,"Requested negative column: %D",col);
48: MatGetSize(A,PETSC_NULL,&N);
49: if (col >= N) SETERRQ2(PETSC_ERR_ARG_OUTOFRANGE,"Requested column %D larger than number columns in matrix %D",col,N);
51: MatGetOwnershipRange(A,&Rs,&Re);
53: PetscObjectGetComm((PetscObject)yy,&comm);
54: VecGetOwnershipRange(yy,&rs,&re);
55: if (Rs != rs || Re != re) SETERRQ4(PETSC_ERR_ARG_INCOMP,"Matrix %D %D does not have same ownership range (size) as vector %D %D",Rs,Re,rs,re);
57: VecSet(yy,zero);
58: VecGetArray(yy,&y);
60: for (i=Rs; i<Re; i++) {
61: MatGetRow(A,i,&nz,&idx,&v);
62: if (nz && idx[0] <= col) {
63: /*
64: Should use faster search here
65: */
66: for (j=0; j<nz; j++) {
67: if (idx[j] >= col) {
68: if (idx[j] == col) y[i-rs] = v[j];
69: break;
70: }
71: }
72: }
73: MatRestoreRow(A,i,&nz,&idx,&v);
74: }
76: VecRestoreArray(yy,&y);
77: return(0);
78: }