10925cdddSBarry Smith #ifdef PETSC_RCS_HEADER 2*82bf6240SBarry Smith static char vcid[] = "$Id: getcolv.c,v 1.1 1998/01/28 21:03:00 bsmith Exp bsmith $"; 30925cdddSBarry Smith #endif 40925cdddSBarry Smith 50925cdddSBarry Smith #include "src/mat/matimpl.h" /*I "mat.h" I*/ 60925cdddSBarry Smith 70925cdddSBarry Smith #undef __FUNC__ 8*82bf6240SBarry Smith #define __FUNC__ "MatGetColumnVector" 90925cdddSBarry Smith /*@ 10*82bf6240SBarry Smith MatGetColumnVector - Gets the values from a given column of a matrix. 110925cdddSBarry Smith 120925cdddSBarry Smith Input Parameters: 13*82bf6240SBarry Smith . X - the matrix 14*82bf6240SBarry Smith . v - the vector 15*82bf6240SBarry Smith . c - the column requested 160925cdddSBarry Smith 17*82bf6240SBarry Smith Contributed by: Denis Vanderstraeten 180925cdddSBarry Smith 19*82bf6240SBarry Smith .keywords: matrix, column, get 200925cdddSBarry Smith 210925cdddSBarry Smith @*/ 22*82bf6240SBarry Smith int MatGetColumnVector(Mat A, Vec yy, int col) 230925cdddSBarry Smith { 24*82bf6240SBarry Smith Scalar *y,*v,zero = 0.0; 25*82bf6240SBarry Smith int ierr,i,j,nz,*idx,M,N,Mv,Rs,Re,rs,re; 260925cdddSBarry Smith 270925cdddSBarry Smith PetscFunctionBegin; 28*82bf6240SBarry Smith PetscValidHeaderSpecific(A,MAT_COOKIE); 29*82bf6240SBarry Smith PetscValidHeaderSpecific(yy,VEC_COOKIE); 300925cdddSBarry Smith 31*82bf6240SBarry Smith if (col < 0) SETERRQ(1,1,"Requested negative column"); 32*82bf6240SBarry Smith ierr = MatGetSize(A,&M,&N); CHKERRQ(ierr); 33*82bf6240SBarry Smith if (col >= N) SETERRQ(1,1,"Requested column larger than number columns in matrix"); 340925cdddSBarry Smith 35*82bf6240SBarry Smith ierr = VecGetSize(yy,&Mv); CHKERRQ(ierr); 36*82bf6240SBarry Smith if (M != Mv) SETERRQ(1,1,"Matrix does not have same number of columns as vector"); 37*82bf6240SBarry Smith 38*82bf6240SBarry Smith ierr = MatGetOwnershipRange(A,&Rs,&Re);CHKERRQ(ierr); 39*82bf6240SBarry Smith ierr = VecGetOwnershipRange(yy,&rs,&re);CHKERRQ(ierr); 40*82bf6240SBarry Smith if (Rs != rs || Re != re) SETERRQ(1,1,"Matrix does not have same ownership range as vector"); 41*82bf6240SBarry Smith 42*82bf6240SBarry Smith ierr = VecSet(&zero,yy);CHKERRQ(ierr); 43*82bf6240SBarry Smith ierr = VecGetArray(yy,&y);CHKERRQ(ierr); 44*82bf6240SBarry Smith 45*82bf6240SBarry Smith for ( i=Rs; i<Re; i++ ) { 46*82bf6240SBarry Smith ierr = MatGetRow(A,i,&nz,&idx,&v);CHKERRQ(ierr); 47*82bf6240SBarry Smith if (nz && idx[0] <= col) { 48*82bf6240SBarry Smith /* 49*82bf6240SBarry Smith Should use faster search here 50*82bf6240SBarry Smith */ 51*82bf6240SBarry Smith for ( j=0; j<nz; j++ ) { 52*82bf6240SBarry Smith if (idx[j] >= col) { 53*82bf6240SBarry Smith if (idx[j] == col) y[i-rs] = v[j]; 54*82bf6240SBarry Smith break; 550925cdddSBarry Smith } 560925cdddSBarry Smith } 570925cdddSBarry Smith } 58*82bf6240SBarry Smith ierr = MatRestoreRow(A,i,&nz,&idx,&v);CHKERRQ(ierr); 590925cdddSBarry Smith } 600925cdddSBarry Smith 61*82bf6240SBarry Smith ierr = VecRestoreArray(yy,&y);CHKERRQ(ierr); 620925cdddSBarry Smith PetscFunctionReturn(0); 630925cdddSBarry Smith } 64