xref: /petsc/src/mat/impls/aij/seq/aij.h (revision f22f69f01ec1ca0bb29797ff3752b65f1e1c47db)
1 
2 #if !defined(__AIJ_H)
3 #define __AIJ_H
4 
5 #include "private/matimpl.h"
6 
7 /*
8     Struct header shared by SeqAIJ, SeqBAIJ and SeqSBAIJ matrix formats
9 */
10 #define SEQAIJHEADER(datatype)	\
11   PetscTruth        roworiented;      /* if true, row-oriented input, default */\
12   PetscInt          nonew;            /* 1 don't add new nonzeros, -1 generate error on new */\
13   PetscInt          nounused;         /* -1 generate error on unused space */\
14   PetscTruth        singlemalloc;     /* if true a, i, and j have been obtained with one big malloc */\
15   PetscInt          maxnz;            /* allocated nonzeros */\
16   PetscInt          *imax;            /* maximum space allocated for each row */\
17   PetscInt          *ilen;            /* actual length of each row */\
18   PetscTruth        free_imax_ilen;  \
19   PetscInt          reallocs;         /* number of mallocs done during MatSetValues() \
20                                         as more values are set than were prealloced */\
21   PetscInt          rmax;             /* max nonzeros in any row */\
22   PetscTruth        keepnonzeropattern;   /* keeps matrix structure same in calls to MatZeroRows()*/\
23   PetscTruth        ignorezeroentries; \
24   PetscInt          *xtoy,*xtoyB;     /* map nonzero pattern of X into Y's, used by MatAXPY() */\
25   Mat               XtoY;             /* used by MatAXPY() */\
26   PetscTruth        free_ij;          /* free the column indices j and row offsets i when the matrix is destroyed */ \
27   PetscTruth        free_a;           /* free the numerical values when matrix is destroy */ \
28   Mat_CompressedRow compressedrow;    /* use compressed row format */                      \
29   PetscInt          nz;               /* nonzeros */                                       \
30   PetscInt          *i;               /* pointer to beginning of each row */               \
31   PetscInt          *j;               /* column values: j + i[k] - 1 is start of row k */  \
32   PetscInt          *diag;            /* pointers to diagonal elements */                  \
33   PetscTruth        free_diag;         \
34   datatype          *a;               /* nonzero elements */                               \
35   PetscScalar       *solve_work;      /* work space used in MatSolve */                    \
36   IS                row, col, icol;   /* index sets, used for reorderings */ \
37   PetscTruth        pivotinblocks;    /* pivot inside factorization of each diagonal block */ \
38   Mat               parent             /* set if this matrix was formed with MatDuplicate(...,MAT_SHARE_NONZERO_PATTERN,....);
39                                          means that this shares some data structures with the parent including diag, ilen, imax, i, j */
40 
41 /*
42   MATSEQAIJ format - Compressed row storage (also called Yale sparse matrix
43   format) or compressed sparse row (CSR).  The i[] and j[] arrays start at 0. For example,
44   j[i[k]+p] is the pth column in row k.  Note that the diagonal
45   matrix elements are stored with the rest of the nonzeros (not separately).
46 */
47 
48 /* Info about i-nodes (identical nodes) helper class for SeqAIJ */
49 typedef struct {
50   MatScalar   *bdiag,*ibdiag,*ssor_work;      /* diagonal blocks of matrix used for MatSOR_SeqAIJ_Inode() */
51   PetscInt    bdiagsize;                       /* length of bdiag and ibdiag */
52   PetscTruth  ibdiagvalid;                     /* do ibdiag[] and bdiag[] contain the most recent values */
53 
54   PetscTruth use;
55   PetscInt   node_count;                    /* number of inodes */
56   PetscInt   *size;                         /* size of each inode */
57   PetscInt   limit;                         /* inode limit */
58   PetscInt   max_limit;                     /* maximum supported inode limit */
59   PetscTruth checked;                       /* if inodes have been checked for */
60 } Mat_SeqAIJ_Inode;
61 
62 EXTERN PetscErrorCode MatView_SeqAIJ_Inode(Mat,PetscViewer);
63 EXTERN PetscErrorCode MatAssemblyEnd_SeqAIJ_Inode(Mat,MatAssemblyType);
64 EXTERN PetscErrorCode MatDestroy_SeqAIJ_Inode(Mat);
65 EXTERN PetscErrorCode MatCreate_SeqAIJ_Inode(Mat);
66 EXTERN PetscErrorCode MatSetOption_SeqAIJ_Inode(Mat,MatOption,PetscTruth);
67 EXTERN PetscErrorCode MatDuplicate_SeqAIJ_Inode(Mat,MatDuplicateOption,Mat*);
68 EXTERN PetscErrorCode MatLUFactorNumeric_SeqAIJ_Inode_inplace(Mat,Mat,const MatFactorInfo*);
69 EXTERN PetscErrorCode MatLUFactorNumeric_SeqAIJ_Inode(Mat,Mat,const MatFactorInfo*);
70 
71 typedef struct {
72   SEQAIJHEADER(MatScalar);
73   Mat_SeqAIJ_Inode inode;
74   MatScalar        *saved_values;             /* location for stashing nonzero values of matrix */
75 
76   PetscScalar      *idiag,*mdiag,*ssor_work;  /* inverse of diagonal entries, diagonal values and workspace for Eisenstat trick */
77   PetscTruth       idiagvalid;                     /* current idiag[] and mdiag[] are valid */
78   PetscScalar      fshift,omega;                   /* last used omega and fshift */
79 
80   ISColoring       coloring;                  /* set with MatADSetColoring() used by MatADSetValues() */
81 } Mat_SeqAIJ;
82 
83 /*
84     Frees the a, i, and j arrays from the XAIJ (AIJ, BAIJ, and SBAIJ) matrix types
85 */
86 #undef __FUNCT__
87 #define __FUNCT__ "MatSeqXAIJFreeAIJ"
88 PETSC_STATIC_INLINE PetscErrorCode MatSeqXAIJFreeAIJ(Mat AA,MatScalar **a,PetscInt **j,PetscInt **i)
89 {
90                                      PetscErrorCode ierr;
91                                      Mat_SeqAIJ     *A = (Mat_SeqAIJ*) AA->data;
92                                      if (A->singlemalloc) {
93                                        ierr = PetscFree3(*a,*j,*i);CHKERRQ(ierr);
94                                      } else {
95                                        if (A->free_a  && *a) {ierr = PetscFree(*a);CHKERRQ(ierr);}
96                                        if (A->free_ij && *j) {ierr = PetscFree(*j);CHKERRQ(ierr);}
97                                        if (A->free_ij && *i) {ierr = PetscFree(*i);CHKERRQ(ierr);}
98                                      }
99                                      *a = 0; *j = 0; *i = 0;
100                                      return 0;
101                                    }
102 
103 /*
104     Allocates larger a, i, and j arrays for the XAIJ (AIJ, BAIJ, and SBAIJ) matrix types
105 */
106 #define MatSeqXAIJReallocateAIJ(Amat,AM,BS2,NROW,ROW,COL,RMAX,AA,AI,AJ,RP,AP,AIMAX,NONEW,datatype) \
107   if (NROW >= RMAX) {\
108 	Mat_SeqAIJ *Ain = (Mat_SeqAIJ*)Amat->data;\
109         /* there is no extra room in row, therefore enlarge */ \
110         PetscInt   new_nz = AI[AM] + CHUNKSIZE,len,*new_i=0,*new_j=0; \
111         datatype   *new_a; \
112  \
113         if (NONEW == -2) SETERRQ2(PETSC_ERR_ARG_OUTOFRANGE,"New nonzero at (%D,%D) caused a malloc",ROW,COL); \
114         /* malloc new storage space */ \
115         ierr = PetscMalloc3(BS2*new_nz,datatype,&new_a,new_nz,PetscInt,&new_j,AM+1,PetscInt,&new_i);CHKERRQ(ierr);\
116  \
117         /* copy over old data into new slots */ \
118         for (ii=0; ii<ROW+1; ii++) {new_i[ii] = AI[ii];} \
119         for (ii=ROW+1; ii<AM+1; ii++) {new_i[ii] = AI[ii]+CHUNKSIZE;} \
120         ierr = PetscMemcpy(new_j,AJ,(AI[ROW]+NROW)*sizeof(PetscInt));CHKERRQ(ierr); \
121         len = (new_nz - CHUNKSIZE - AI[ROW] - NROW); \
122         ierr = PetscMemcpy(new_j+AI[ROW]+NROW+CHUNKSIZE,AJ+AI[ROW]+NROW,len*sizeof(PetscInt));CHKERRQ(ierr); \
123         ierr = PetscMemcpy(new_a,AA,BS2*(AI[ROW]+NROW)*sizeof(datatype));CHKERRQ(ierr); \
124         ierr = PetscMemzero(new_a+BS2*(AI[ROW]+NROW),BS2*CHUNKSIZE*sizeof(datatype));CHKERRQ(ierr);\
125         ierr = PetscMemcpy(new_a+BS2*(AI[ROW]+NROW+CHUNKSIZE),AA+BS2*(AI[ROW]+NROW),BS2*len*sizeof(datatype));CHKERRQ(ierr);  \
126         /* free up old matrix storage */ \
127         ierr = MatSeqXAIJFreeAIJ(A,&Ain->a,&Ain->j,&Ain->i);CHKERRQ(ierr);\
128         AA = new_a; \
129         Ain->a = (MatScalar*) new_a;		   \
130         AI = Ain->i = new_i; AJ = Ain->j = new_j;  \
131         Ain->singlemalloc = PETSC_TRUE; \
132  \
133         RP          = AJ + AI[ROW]; AP = AA + BS2*AI[ROW]; \
134         RMAX        = AIMAX[ROW] = AIMAX[ROW] + CHUNKSIZE; \
135         Ain->maxnz += CHUNKSIZE; \
136         Ain->reallocs++; \
137       } \
138 
139 EXTERN_C_BEGIN
140 EXTERN PetscErrorCode MatSeqAIJSetPreallocation_SeqAIJ(Mat,PetscInt,PetscInt*);
141 EXTERN_C_END
142 EXTERN PetscErrorCode MatILUFactorSymbolic_SeqAIJ_inplace(Mat,Mat,IS,IS,const MatFactorInfo*);
143 EXTERN PetscErrorCode MatILUFactorSymbolic_SeqAIJ(Mat,Mat,IS,IS,const MatFactorInfo*);
144 EXTERN PetscErrorCode MatILUFactorSymbolic_SeqAIJ_ilu0(Mat,Mat,IS,IS,const MatFactorInfo*);
145 
146 EXTERN PetscErrorCode MatICCFactorSymbolic_SeqAIJ_inplace(Mat,Mat,IS,const MatFactorInfo*);
147 EXTERN PetscErrorCode MatICCFactorSymbolic_SeqAIJ(Mat,Mat,IS,const MatFactorInfo*);
148 EXTERN PetscErrorCode MatCholeskyFactorSymbolic_SeqAIJ_inplace(Mat,Mat,IS,const MatFactorInfo*);
149 EXTERN PetscErrorCode MatCholeskyFactorSymbolic_SeqAIJ(Mat,Mat,IS,const MatFactorInfo*);
150 EXTERN PetscErrorCode MatCholeskyFactorNumeric_SeqAIJ_inplace(Mat,Mat,const MatFactorInfo*);
151 EXTERN PetscErrorCode MatCholeskyFactorNumeric_SeqAIJ(Mat,Mat,const MatFactorInfo*);
152 EXTERN PetscErrorCode MatDuplicate_SeqAIJ(Mat,MatDuplicateOption,Mat*);
153 EXTERN PetscErrorCode MatCopy_SeqAIJ(Mat,Mat,MatStructure);
154 EXTERN PetscErrorCode MatMissingDiagonal_SeqAIJ(Mat,PetscTruth*,PetscInt*);
155 EXTERN PetscErrorCode MatMarkDiagonal_SeqAIJ(Mat);
156 
157 EXTERN PetscErrorCode MatMult_SeqAIJ(Mat A,Vec,Vec);
158 EXTERN PetscErrorCode MatMultAdd_SeqAIJ(Mat A,Vec,Vec,Vec);
159 EXTERN PetscErrorCode MatMultTranspose_SeqAIJ(Mat A,Vec,Vec);
160 EXTERN PetscErrorCode MatMultTransposeAdd_SeqAIJ(Mat A,Vec,Vec,Vec);
161 EXTERN PetscErrorCode MatSOR_SeqAIJ(Mat,Vec,PetscReal,MatSORType,PetscReal,PetscInt,PetscInt,Vec);
162 
163 EXTERN PetscErrorCode MatSetColoring_SeqAIJ(Mat,ISColoring);
164 EXTERN PetscErrorCode MatSetValuesAdic_SeqAIJ(Mat,void*);
165 EXTERN PetscErrorCode MatSetValuesAdifor_SeqAIJ(Mat,PetscInt,void*);
166 
167 EXTERN PetscErrorCode MatGetSymbolicTranspose_SeqAIJ(Mat,PetscInt *[],PetscInt *[]);
168 EXTERN PetscErrorCode MatGetSymbolicTransposeReduced_SeqAIJ(Mat,PetscInt,PetscInt,PetscInt *[],PetscInt *[]);
169 EXTERN PetscErrorCode MatRestoreSymbolicTranspose_SeqAIJ(Mat,PetscInt *[],PetscInt *[]);
170 EXTERN PetscErrorCode MatToSymmetricIJ_SeqAIJ(PetscInt,PetscInt*,PetscInt*,PetscInt,PetscInt,PetscInt**,PetscInt**);
171 EXTERN PetscErrorCode MatLUFactorSymbolic_SeqAIJ_inplace(Mat,Mat,IS,IS,const MatFactorInfo*);
172 EXTERN PetscErrorCode MatLUFactorSymbolic_SeqAIJ(Mat,Mat,IS,IS,const MatFactorInfo*);
173 EXTERN PetscErrorCode MatLUFactorNumeric_SeqAIJ_inplace(Mat,Mat,const MatFactorInfo*);
174 EXTERN PetscErrorCode MatLUFactorNumeric_SeqAIJ(Mat,Mat,const MatFactorInfo*);
175 EXTERN PetscErrorCode MatLUFactorNumeric_SeqAIJ_InplaceWithPerm(Mat,Mat,const MatFactorInfo*);
176 EXTERN PetscErrorCode MatLUFactor_SeqAIJ(Mat,IS,IS,const MatFactorInfo*);
177 EXTERN PetscErrorCode MatSolve_SeqAIJ_inplace(Mat,Vec,Vec);
178 EXTERN PetscErrorCode MatSolve_SeqAIJ(Mat,Vec,Vec);
179 EXTERN PetscErrorCode MatSolve_SeqAIJ_Inode_inplace(Mat,Vec,Vec);
180 EXTERN PetscErrorCode MatSolve_SeqAIJ_Inode(Mat,Vec,Vec);
181 EXTERN PetscErrorCode MatSolve_SeqAIJ_NaturalOrdering_inplace(Mat,Vec,Vec);
182 EXTERN PetscErrorCode MatSolve_SeqAIJ_NaturalOrdering(Mat,Vec,Vec);
183 EXTERN PetscErrorCode MatSolve_SeqAIJ_InplaceWithPerm(Mat,Vec,Vec);
184 EXTERN PetscErrorCode MatSolveAdd_SeqAIJ_inplace(Mat,Vec,Vec,Vec);
185 EXTERN PetscErrorCode MatSolveAdd_SeqAIJ(Mat,Vec,Vec,Vec);
186 EXTERN PetscErrorCode MatSolveTranspose_SeqAIJ_inplace(Mat,Vec,Vec);
187 EXTERN PetscErrorCode MatSolveTranspose_SeqAIJ(Mat,Vec,Vec);
188 EXTERN PetscErrorCode MatSolveTransposeAdd_SeqAIJ_inplace(Mat,Vec,Vec,Vec);
189 EXTERN PetscErrorCode MatSolveTransposeAdd_SeqAIJ(Mat,Vec,Vec,Vec);
190 EXTERN PetscErrorCode MatMatSolve_SeqAIJ_inplace(Mat,Mat,Mat);
191 EXTERN PetscErrorCode MatMatSolve_SeqAIJ(Mat,Mat,Mat);
192 EXTERN PetscErrorCode MatEqual_SeqAIJ(Mat A,Mat B,PetscTruth* flg);
193 EXTERN PetscErrorCode MatFDColoringCreate_SeqAIJ(Mat,ISColoring,MatFDColoring);
194 EXTERN PetscErrorCode MatLoad_SeqAIJ(PetscViewer, const MatType,Mat*);
195 EXTERN PetscErrorCode RegisterApplyPtAPRoutines_Private(Mat);
196 EXTERN PetscErrorCode MatMatMultSymbolic_SeqAIJ_SeqAIJ(Mat,Mat,PetscReal,Mat*);
197 EXTERN PetscErrorCode MatMatMultNumeric_SeqAIJ_SeqAIJ(Mat,Mat,Mat);
198 EXTERN PetscErrorCode MatPtAPSymbolic_SeqAIJ(Mat,Mat,PetscReal,Mat*);
199 EXTERN PetscErrorCode MatPtAPNumeric_SeqAIJ(Mat,Mat,Mat);
200 EXTERN PetscErrorCode MatPtAPSymbolic_SeqAIJ_SeqAIJ(Mat,Mat,PetscReal,Mat*);
201 EXTERN PetscErrorCode MatPtAPNumeric_SeqAIJ_SeqAIJ(Mat,Mat,Mat);
202 EXTERN PetscErrorCode MatMatMultTranspose_SeqAIJ_SeqAIJ(Mat,Mat,MatReuse,PetscReal,Mat*);
203 EXTERN PetscErrorCode MatMatMultTransposeSymbolic_SeqAIJ_SeqAIJ(Mat,Mat,PetscReal,Mat*);
204 EXTERN PetscErrorCode MatMatMultTransposeNumeric_SeqAIJ_SeqAIJ(Mat,Mat,Mat);
205 EXTERN PetscErrorCode MatSetValues_SeqAIJ(Mat,PetscInt,const PetscInt[],PetscInt,const PetscInt[],const PetscScalar[],InsertMode);
206 EXTERN PetscErrorCode MatGetRow_SeqAIJ(Mat,PetscInt,PetscInt*,PetscInt**,PetscScalar**);
207 EXTERN PetscErrorCode MatRestoreRow_SeqAIJ(Mat,PetscInt,PetscInt*,PetscInt**,PetscScalar**);
208 EXTERN PetscErrorCode MatAXPY_SeqAIJ(Mat,PetscScalar,Mat,MatStructure);
209 EXTERN PetscErrorCode MatGetRowIJ_SeqAIJ(Mat,PetscInt,PetscTruth,PetscTruth,PetscInt*,PetscInt *[],PetscInt *[],PetscTruth *);
210 EXTERN PetscErrorCode MatRestoreRowIJ_SeqAIJ(Mat,PetscInt,PetscTruth,PetscTruth,PetscInt *,PetscInt *[],PetscInt *[],PetscTruth *);
211 EXTERN PetscErrorCode MatGetColumnIJ_SeqAIJ(Mat,PetscInt,PetscTruth,PetscTruth,PetscInt*,PetscInt *[],PetscInt *[],PetscTruth *);
212 EXTERN PetscErrorCode MatRestoreColumnIJ_SeqAIJ(Mat,PetscInt,PetscTruth,PetscTruth,PetscInt *,PetscInt *[],PetscInt *[],PetscTruth *);
213 EXTERN PetscErrorCode MatDestroy_SeqAIJ(Mat);
214 EXTERN PetscErrorCode MatView_SeqAIJ(Mat,PetscViewer);
215 
216 EXTERN PetscErrorCode Mat_CheckInode(Mat,PetscTruth);
217 EXTERN PetscErrorCode Mat_CheckInode_FactorLU(Mat,PetscTruth);
218 
219 EXTERN_C_BEGIN
220 EXTERN PetscErrorCode PETSCMAT_DLLEXPORT MatConvert_SeqAIJ_SeqSBAIJ(Mat,const MatType,MatReuse,Mat*);
221 EXTERN PetscErrorCode PETSCMAT_DLLEXPORT MatConvert_SeqAIJ_SeqBAIJ(Mat,const MatType,MatReuse,Mat*);
222 EXTERN PetscErrorCode PETSCMAT_DLLEXPORT MatConvert_SeqAIJ_SeqCSRPERM(Mat,const MatType,MatReuse,Mat*);
223 EXTERN PetscErrorCode PETSCMAT_DLLEXPORT MatReorderForNonzeroDiagonal_SeqAIJ(Mat,PetscReal,IS,IS);
224 EXTERN PetscErrorCode PETSCMAT_DLLEXPORT MatMatMult_SeqAIJ_SeqAIJ(Mat,Mat,MatReuse,PetscReal,Mat*);
225 EXTERN_C_END
226 
227 /*
228     PetscSparseDenseMinusDot - The inner kernel of triangular solves and Gauss-Siedel smoothing. \sum_i xv[i] * r[xi[i]] for CSR storage
229 
230   Input Parameters:
231 +  nnz - the number of entries
232 .  r - the array of vector values
233 .  xv - the matrix values for the row
234 -  xi - the column indices of the nonzeros in the row
235 
236   Output Parameter:
237 .  sum - negative the sum of results
238 
239   PETSc compile flags:
240 +   PETSC_KERNEL_USE_UNROLL_4 -   don't use this; it changes nnz and hence is WRONG
241 -   PETSC_KERNEL_USE_UNROLL_2 -
242 
243 .seealso: PetscSparseDensePlusDot()
244 
245 */
246 #ifdef PETSC_KERNEL_USE_UNROLL_4
247 #define PetscSparseDenseMinusDot(sum,r,xv,xi,nnz) {\
248 if (nnz > 0) {\
249 switch (nnz & 0x3) {\
250 case 3: sum -= *xv++ * r[*xi++];\
251 case 2: sum -= *xv++ * r[*xi++];\
252 case 1: sum -= *xv++ * r[*xi++];\
253 nnz -= 4;}\
254 while (nnz > 0) {\
255 sum -=  xv[0] * r[xi[0]] - xv[1] * r[xi[1]] -\
256 	xv[2] * r[xi[2]] - xv[3] * r[xi[3]];\
257 xv  += 4; xi += 4; nnz -= 4; }}}
258 
259 #elif defined(PETSC_KERNEL_USE_UNROLL_2)
260 #define PetscSparseDenseMinusDot(sum,r,xv,xi,nnz) {\
261 PetscInt __i,__i1,__i2;\
262 for(__i=0;__i<nnz-1;__i+=2) {__i1 = xi[__i]; __i2=xi[__i+1];\
263 sum -= (xv[__i]*r[__i1] + xv[__i+1]*r[__i2]);}\
264 if (nnz & 0x1) sum -= xv[__i] * r[xi[__i]];}
265 
266 #else
267 #define PetscSparseDenseMinusDot(sum,r,xv,xi,nnz) {\
268 PetscInt __i;\
269 for(__i=0;__i<nnz;__i++) sum -= xv[__i] * r[xi[__i]];}
270 #endif
271 
272 
273 
274 /*
275     PetscSparseDensePlusDot - The inner kernel of matrix-vector product \sum_i xv[i] * r[xi[i]] for CSR storage
276 
277   Input Parameters:
278 +  nnz - the number of entries
279 .  r - the array of vector values
280 .  xv - the matrix values for the row
281 -  xi - the column indices of the nonzeros in the row
282 
283   Output Parameter:
284 .  sum - the sum of results
285 
286   PETSc compile flags:
287 +   PETSC_KERNEL_USE_UNROLL_4 -  don't use this; it changes nnz and hence is WRONG
288 -   PETSC_KERNEL_USE_UNROLL_2 -
289 
290 .seealso: PetscSparseDenseMinusDot()
291 
292 */
293 #ifdef PETSC_KERNEL_USE_UNROLL_4
294 #define PetscSparseDensePlusDot(sum,r,xv,xi,nnz) {\
295 if (nnz > 0) {\
296 switch (nnz & 0x3) {\
297 case 3: sum += *xv++ * r[*xi++];\
298 case 2: sum += *xv++ * r[*xi++];\
299 case 1: sum += *xv++ * r[*xi++];\
300 nnz -= 4;}\
301 while (nnz > 0) {\
302 sum +=  xv[0] * r[xi[0]] + xv[1] * r[xi[1]] +\
303 	xv[2] * r[xi[2]] + xv[3] * r[xi[3]];\
304 xv  += 4; xi += 4; nnz -= 4; }}}
305 
306 #elif defined(PETSC_KERNEL_USE_UNROLL_2)
307 #define PetscSparseDensePlusDot(sum,r,xv,xi,nnz) {\
308 PetscInt __i,__i1,__i2;\
309 for(__i=0;__i<nnz-1;__i+=2) {__i1 = xi[__i]; __i2=xi[__i+1];\
310 sum += (xv[__i]*r[__i1] + xv[__i+1]*r[__i2]);}\
311 if (nnz & 0x1) sum += xv[__i] * r[xi[__i]];}
312 
313 #else
314 #define PetscSparseDensePlusDot(sum,r,xv,xi,nnz) {\
315  PetscInt __i;\
316 for(__i=0;__i<nnz;__i++) sum += xv[__i] * r[xi[__i]];}
317 #endif
318 
319 #endif
320