10925cdddSBarry Smith #ifdef PETSC_RCS_HEADER 2*15091d37SBarry Smith static char vcid[] = "$Id: getcolv.c,v 1.5 1999/01/04 21:51:52 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__ 882bf6240SBarry Smith #define __FUNC__ "MatGetColumnVector" 90925cdddSBarry Smith /*@ 1082bf6240SBarry Smith MatGetColumnVector - Gets the values from a given column of a matrix. 110925cdddSBarry Smith 12fee21e36SBarry Smith Collective on Mat and Vec 13fee21e36SBarry Smith 1498a79cdbSBarry Smith Input Parameters: 1598a79cdbSBarry Smith + X - the matrix 1698a79cdbSBarry Smith . v - the vector 1798a79cdbSBarry Smith - c - the column requested 1898a79cdbSBarry Smith 19*15091d37SBarry Smith Level: advanced 20*15091d37SBarry Smith 2182bf6240SBarry Smith Contributed by: Denis Vanderstraeten 220925cdddSBarry Smith 2382bf6240SBarry Smith .keywords: matrix, column, get 240925cdddSBarry Smith 25*15091d37SBarry Smith .seealso: MatGetRow(), MatGetDiagonal() 26*15091d37SBarry Smith 270925cdddSBarry Smith @*/ 2882bf6240SBarry Smith int MatGetColumnVector(Mat A, Vec yy, int col) 290925cdddSBarry Smith { 3082bf6240SBarry Smith Scalar *y,*v,zero = 0.0; 3182bf6240SBarry Smith int ierr,i,j,nz,*idx,M,N,Mv,Rs,Re,rs,re; 320925cdddSBarry Smith 330925cdddSBarry Smith PetscFunctionBegin; 3482bf6240SBarry Smith PetscValidHeaderSpecific(A,MAT_COOKIE); 3582bf6240SBarry Smith PetscValidHeaderSpecific(yy,VEC_COOKIE); 360925cdddSBarry Smith 37596552b5SBarry Smith if (col < 0) SETERRQ1(1,1,"Requested negative column: %d",col); 3882bf6240SBarry Smith ierr = MatGetSize(A,&M,&N); CHKERRQ(ierr); 39596552b5SBarry Smith if (col >= N) SETERRQ2(1,1,"Requested column %d larger than number columns in matrix %d",col,N); 400925cdddSBarry Smith 4182bf6240SBarry Smith ierr = VecGetSize(yy,&Mv); CHKERRQ(ierr); 42596552b5SBarry Smith if (M != Mv) SETERRQ2(1,1,"Matrix does not have same number of columns %d as vector %d",M,Mv); 4382bf6240SBarry Smith 4482bf6240SBarry Smith ierr = MatGetOwnershipRange(A,&Rs,&Re);CHKERRQ(ierr); 4582bf6240SBarry Smith ierr = VecGetOwnershipRange(yy,&rs,&re);CHKERRQ(ierr); 46596552b5SBarry Smith 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); 4782bf6240SBarry Smith 4882bf6240SBarry Smith ierr = VecSet(&zero,yy);CHKERRQ(ierr); 4982bf6240SBarry Smith ierr = VecGetArray(yy,&y);CHKERRQ(ierr); 5082bf6240SBarry Smith 5182bf6240SBarry Smith for ( i=Rs; i<Re; i++ ) { 5282bf6240SBarry Smith ierr = MatGetRow(A,i,&nz,&idx,&v);CHKERRQ(ierr); 5382bf6240SBarry Smith if (nz && idx[0] <= col) { 5482bf6240SBarry Smith /* 5582bf6240SBarry Smith Should use faster search here 5682bf6240SBarry Smith */ 5782bf6240SBarry Smith for ( j=0; j<nz; j++ ) { 5882bf6240SBarry Smith if (idx[j] >= col) { 5982bf6240SBarry Smith if (idx[j] == col) y[i-rs] = v[j]; 6082bf6240SBarry Smith break; 610925cdddSBarry Smith } 620925cdddSBarry Smith } 630925cdddSBarry Smith } 6482bf6240SBarry Smith ierr = MatRestoreRow(A,i,&nz,&idx,&v);CHKERRQ(ierr); 650925cdddSBarry Smith } 660925cdddSBarry Smith 6782bf6240SBarry Smith ierr = VecRestoreArray(yy,&y);CHKERRQ(ierr); 680925cdddSBarry Smith PetscFunctionReturn(0); 690925cdddSBarry Smith } 70