xref: /petsc/src/mat/utils/getcolv.c (revision e6e75211d226c622f451867f53ce5d558649ff4f)
1 
2 #include <petsc/private/matimpl.h>  /*I   "petscmat.h"  I*/
3 
4 #undef __FUNCT__
5 #define __FUNCT__ "MatGetColumnVector"
6 /*@
7    MatGetColumnVector - Gets the values from a given column of a matrix.
8 
9    Not Collective
10 
11    Input Parameters:
12 +  A - the matrix
13 .  yy - the vector
14 -  c - the column requested (in global numbering)
15 
16    Level: advanced
17 
18    Notes:
19    Each processor for which this is called gets the values for its rows.
20 
21    Since PETSc matrices are usually stored in compressed row format, this routine
22    will generally be slow.
23 
24    The vector must have the same parallel row layout as the matrix.
25 
26    Contributed by: Denis Vanderstraeten
27 
28 .keywords: matrix, column, get
29 
30 .seealso: MatGetRow(), MatGetDiagonal()
31 
32 @*/
33 PetscErrorCode  MatGetColumnVector(Mat A,Vec yy,PetscInt col)
34 {
35   PetscScalar       *y;
36   const PetscScalar *v;
37   PetscErrorCode    ierr;
38   PetscInt          i,j,nz,N,Rs,Re,rs,re;
39   const PetscInt    *idx;
40 
41   PetscFunctionBegin;
42   PetscValidHeaderSpecific(A,MAT_CLASSID,1);
43   PetscValidHeaderSpecific(yy,VEC_CLASSID,2);
44   if (col < 0) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Requested negative column: %D",col);
45   ierr = MatGetSize(A,NULL,&N);CHKERRQ(ierr);
46   if (col >= N) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Requested column %D larger than number columns in matrix %D",col,N);
47   ierr = MatGetOwnershipRange(A,&Rs,&Re);CHKERRQ(ierr);
48   ierr = VecGetOwnershipRange(yy,&rs,&re);CHKERRQ(ierr);
49   if (Rs != rs || Re != re) SETERRQ4(PETSC_COMM_SELF,PETSC_ERR_ARG_INCOMP,"Matrix %D %D does not have same ownership range (size) as vector %D %D",Rs,Re,rs,re);
50 
51   if (A->ops->getcolumnvector) {
52     ierr = (*A->ops->getcolumnvector)(A,yy,col);CHKERRQ(ierr);
53   } else {
54     ierr = VecSet(yy,0.0);CHKERRQ(ierr);
55     ierr = VecGetArray(yy,&y);CHKERRQ(ierr);
56 
57     for (i=Rs; i<Re; i++) {
58       ierr = MatGetRow(A,i,&nz,&idx,&v);CHKERRQ(ierr);
59       if (nz && idx[0] <= col) {
60         /*
61           Should use faster search here
62         */
63         for (j=0; j<nz; j++) {
64           if (idx[j] >= col) {
65             if (idx[j] == col) y[i-rs] = v[j];
66             break;
67           }
68         }
69       }
70       ierr = MatRestoreRow(A,i,&nz,&idx,&v);CHKERRQ(ierr);
71     }
72     ierr = VecRestoreArray(yy,&y);CHKERRQ(ierr);
73   }
74   PetscFunctionReturn(0);
75 }
76 
77 
78 
79 
80 #undef __FUNCT__
81 #define __FUNCT__ "MatGetColumnNorms"
82 /*@
83     MatGetColumnNorms - Gets the norms of each column of a sparse or dense matrix.
84 
85   Input Parameter:
86 +  A - the matrix
87 -  type - NORM_2, NORM_1 or NORM_INFINITY
88 
89   Output Parameter:
90 .  norms - an array as large as the TOTAL number of columns in the matrix
91 
92    Level: intermediate
93 
94    Notes: Each process has ALL the column norms after the call. Because of the way this is computed each process gets all the values,
95     if each process wants only some of the values it should extract the ones it wants from the array.
96 
97 .seealso: MatGetColumns()
98 
99 @*/
100 PetscErrorCode MatGetColumnNorms(Mat A,NormType type,PetscReal norms[])
101 {
102   PetscErrorCode ierr;
103 
104   PetscFunctionBegin;
105   PetscValidHeaderSpecific(A,MAT_CLASSID,1);
106   if (A->ops->getcolumnnorms) {
107     ierr = (*A->ops->getcolumnnorms)(A,type,norms);CHKERRQ(ierr);
108   } else SETERRQ(PetscObjectComm((PetscObject)A),PETSC_ERR_SUP,"Not coded for this matrix type");
109   PetscFunctionReturn(0);
110 }
111 
112