xref: /petsc/src/mat/utils/getcolv.c (revision 95d5f7c29374efcfd2ca44c2fe93981fbc2b4454)
1 /*$Id: getcolv.c,v 1.10 2000/04/09 04:36:57 bsmith Exp bsmith $*/
2 
3 #include "src/mat/matimpl.h"  /*I   "mat.h"  I*/
4 
5 #undef __FUNC__
6 #define __FUNC__ /*<a name=""></a>*/"MatGetColumnVector"
7 /*@
8    MatGetColumnVector - Gets the values from a given column of a matrix.
9 
10    Collective on Mat and Vec
11 
12    Input Parameters:
13 +  X - the matrix
14 .  v - the vector
15 -  c - the column requested
16 
17    Level: advanced
18 
19    Contributed by: Denis Vanderstraeten
20 
21 .keywords: matrix, column, get
22 
23 .seealso: MatGetRow(), MatGetDiagonal()
24 
25 @*/
26 int MatGetColumnVector(Mat A,Vec yy,int col)
27 {
28   Scalar *y,*v,zero = 0.0;
29   int    ierr,i,j,nz,*idx,M,N,Mv,Rs,Re,rs,re;
30 
31   PetscFunctionBegin;
32   PetscValidHeaderSpecific(A,MAT_COOKIE);
33   PetscValidHeaderSpecific(yy,VEC_COOKIE);
34 
35   if (col < 0)  SETERRQ1(1,1,"Requested negative column: %d",col);
36   ierr = MatGetSize(A,&M,&N);CHKERRQ(ierr);
37   if (col >= N)  SETERRQ2(1,1,"Requested column %d larger than number columns in matrix %d",col,N);
38 
39   ierr = VecGetSize(yy,&Mv);CHKERRQ(ierr);
40   if (M != Mv) SETERRQ2(1,1,"Matrix does not have same number of columns %d as vector %d",M,Mv);
41 
42   ierr = MatGetOwnershipRange(A,&Rs,&Re);CHKERRQ(ierr);
43   ierr = VecGetOwnershipRange(yy,&rs,&re);CHKERRQ(ierr);
44   if (Rs != rs || Re != re) SETERRQ4(1,1,"Matrix %d %d does not have same ownership range as vector %d %d",Rs,Re,rs,re);
45 
46   ierr = VecSet(&zero,yy);CHKERRQ(ierr);
47   ierr = VecGetArray(yy,&y);CHKERRQ(ierr);
48 
49   for (i=Rs; i<Re; i++) {
50     ierr = MatGetRow(A,i,&nz,&idx,&v);CHKERRQ(ierr);
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     ierr = MatRestoreRow(A,i,&nz,&idx,&v);CHKERRQ(ierr);
63   }
64 
65   ierr = VecRestoreArray(yy,&y);CHKERRQ(ierr);
66   PetscFunctionReturn(0);
67 }
68