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