xref: /petsc/src/mat/utils/getcolv.c (revision 7854c7f8112cae921e4370fc5d1e189095717a21)
1 
2 #include "src/mat/matimpl.h"  /*I   "petscmat.h"  I*/
3 
4 #undef __FUNCT__
5 #define __FUNCT__ "MatGetColumnVector"
6 /*@
7    MatGetColumnVector - Gets the values from a given column of a matrix.
8 
9    Not Collective
10 
11    Input Parameters:
12 +  A - the matrix
13 .  yy - the vector
14 -  c - the column requested (in global numbering)
15 
16    Level: advanced
17 
18    Notes:
19    Each processor for which this is called gets the values for its rows.
20 
21    Since PETSc matrices are usually stored in compressed row format, this routine
22    will generally be slow.
23 
24    The vector must have the same parallel row layout as the matrix.
25 
26    Contributed by: Denis Vanderstraeten
27 
28 .keywords: matrix, column, get
29 
30 .seealso: MatGetRow(), MatGetDiagonal()
31 
32 @*/
33 PetscErrorCode MatGetColumnVector(Mat A,Vec yy,PetscInt col)
34 {
35   PetscScalar        *y,zero = 0.0;
36   const PetscScalar  *v;
37   PetscErrorCode     ierr;
38   PetscInt           i,j,nz,N,Rs,Re,rs,re;
39   const PetscInt     *idx;
40   MPI_Comm           comm;
41 
42   PetscFunctionBegin;
43   PetscValidHeaderSpecific(A,MAT_COOKIE,1);
44   PetscValidHeaderSpecific(yy,VEC_COOKIE,2);
45 
46   if (col < 0)  SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,"Requested negative column: %D",col);
47   ierr = MatGetSize(A,PETSC_NULL,&N);CHKERRQ(ierr);
48   if (col >= N)  SETERRQ2(PETSC_ERR_ARG_OUTOFRANGE,"Requested column %D larger than number columns in matrix %D",col,N);
49 
50   ierr = MatGetOwnershipRange(A,&Rs,&Re);CHKERRQ(ierr);
51 
52   ierr = PetscObjectGetComm((PetscObject)yy,&comm);CHKERRQ(ierr);
53   ierr = VecGetOwnershipRange(yy,&rs,&re);CHKERRQ(ierr);
54   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);
55 
56   ierr = VecSet(&zero,yy);CHKERRQ(ierr);
57   ierr = VecGetArray(yy,&y);CHKERRQ(ierr);
58 
59   for (i=Rs; i<Re; i++) {
60     ierr = MatGetRow(A,i,&nz,&idx,&v);CHKERRQ(ierr);
61     if (nz && idx[0] <= col) {
62       /*
63           Should use faster search here
64       */
65       for (j=0; j<nz; j++) {
66         if (idx[j] >= col) {
67           if (idx[j] == col) y[i-rs] = v[j];
68           break;
69         }
70       }
71     }
72     ierr = MatRestoreRow(A,i,&nz,&idx,&v);CHKERRQ(ierr);
73   }
74 
75   ierr = VecRestoreArray(yy,&y);CHKERRQ(ierr);
76   PetscFunctionReturn(0);
77 }
78