10925cdddSBarry Smith 2e090d566SSatish Balay #include "src/mat/matimpl.h" /*I "petscmat.h" I*/ 30925cdddSBarry Smith 44a2ae208SSatish Balay #undef __FUNCT__ 54a2ae208SSatish Balay #define __FUNCT__ "MatGetColumnVector" 60925cdddSBarry Smith /*@ 782bf6240SBarry Smith MatGetColumnVector - Gets the values from a given column of a matrix. 80925cdddSBarry Smith 972257631SBarry Smith Not Collective 10fee21e36SBarry Smith 1198a79cdbSBarry Smith Input Parameters: 125ed6d96aSBarry Smith + A - the matrix 135ed6d96aSBarry Smith . yy - the vector 145ed6d96aSBarry Smith - c - the column requested (in global numbering) 1598a79cdbSBarry Smith 1615091d37SBarry Smith Level: advanced 1715091d37SBarry Smith 189d006df2SBarry Smith Notes: 199d006df2SBarry Smith Each processor for which this is called gets the values for its rows. 209d006df2SBarry Smith 219d006df2SBarry Smith Since PETSc matrices are usually stored in compressed row format, this routine 229d006df2SBarry Smith will generally be slow. 239d006df2SBarry Smith 243d81755aSBarry Smith The vector must have the same parallel row layout as the matrix. 253d81755aSBarry Smith 2682bf6240SBarry Smith Contributed by: Denis Vanderstraeten 270925cdddSBarry Smith 2882bf6240SBarry Smith .keywords: matrix, column, get 290925cdddSBarry Smith 3015091d37SBarry Smith .seealso: MatGetRow(), MatGetDiagonal() 3115091d37SBarry Smith 320925cdddSBarry Smith @*/ 33*38baddfdSBarry Smith PetscErrorCode MatGetColumnVector(Mat A,Vec yy,PetscInt col) 340925cdddSBarry Smith { 35b3cc6726SBarry Smith PetscScalar *y,zero = 0.0; 36b3cc6726SBarry Smith const PetscScalar *v; 37dfbe8321SBarry Smith PetscErrorCode ierr; 38*38baddfdSBarry Smith PetscInt i,j,nz,N,Rs,Re,rs,re; 39*38baddfdSBarry Smith const PetscInt *idx; 40011b8408SBarry Smith MPI_Comm comm; 410925cdddSBarry Smith 420925cdddSBarry Smith PetscFunctionBegin; 434482741eSBarry Smith PetscValidHeaderSpecific(A,MAT_COOKIE,1); 444482741eSBarry Smith PetscValidHeaderSpecific(yy,VEC_COOKIE,2); 450925cdddSBarry Smith 46*38baddfdSBarry Smith if (col < 0) SETERRQ1(1,"Requested negative column: %d",(int)col); 47011b8408SBarry Smith ierr = MatGetSize(A,PETSC_NULL,&N);CHKERRQ(ierr); 48*38baddfdSBarry Smith if (col >= N) SETERRQ2(1,"Requested column %d larger than number columns in matrix %d",(int)col,(int)N); 490925cdddSBarry Smith 5082bf6240SBarry Smith ierr = MatGetOwnershipRange(A,&Rs,&Re);CHKERRQ(ierr); 51011b8408SBarry Smith 52011b8408SBarry Smith ierr = PetscObjectGetComm((PetscObject)yy,&comm);CHKERRQ(ierr); 5382bf6240SBarry Smith ierr = VecGetOwnershipRange(yy,&rs,&re);CHKERRQ(ierr); 54*38baddfdSBarry Smith if (Rs != rs || Re != re) SETERRQ4(1,"Matrix %d %d does not have same ownership range (size) as vector %d %d",(int)Rs,(int)Re,(int)rs,(int)re); 5582bf6240SBarry Smith 5682bf6240SBarry Smith ierr = VecSet(&zero,yy);CHKERRQ(ierr); 5782bf6240SBarry Smith ierr = VecGetArray(yy,&y);CHKERRQ(ierr); 5882bf6240SBarry Smith 5982bf6240SBarry Smith for (i=Rs; i<Re; i++) { 6082bf6240SBarry Smith ierr = MatGetRow(A,i,&nz,&idx,&v);CHKERRQ(ierr); 6182bf6240SBarry Smith if (nz && idx[0] <= col) { 6282bf6240SBarry Smith /* 6382bf6240SBarry Smith Should use faster search here 6482bf6240SBarry Smith */ 6582bf6240SBarry Smith for (j=0; j<nz; j++) { 6682bf6240SBarry Smith if (idx[j] >= col) { 6782bf6240SBarry Smith if (idx[j] == col) y[i-rs] = v[j]; 6882bf6240SBarry Smith break; 690925cdddSBarry Smith } 700925cdddSBarry Smith } 710925cdddSBarry Smith } 7282bf6240SBarry Smith ierr = MatRestoreRow(A,i,&nz,&idx,&v);CHKERRQ(ierr); 730925cdddSBarry Smith } 740925cdddSBarry Smith 7582bf6240SBarry Smith ierr = VecRestoreArray(yy,&y);CHKERRQ(ierr); 760925cdddSBarry Smith PetscFunctionReturn(0); 770925cdddSBarry Smith } 78