xref: /petsc/src/ksp/pc/impls/gamg/classical.c (revision bfde193ff4e761bc853abc118b6965d7dd4b968a) 
1 #include <../src/ksp/pc/impls/gamg/gamg.h>        /*I "petscpc.h" I*/
2 #include <petsc-private/kspimpl.h>
3 
4 typedef struct {
5   PetscReal dummy; /* empty struct; save for later */
6 } PC_GAMG_Classical;
7 
8 
9 #undef __FUNCT__
10 #define __FUNCT__ "PCGAMGClassicalCreateGhostVector_Private"
11 PetscErrorCode PCGAMGClassicalCreateGhostVector_Private(Mat G,Vec *gvec,PetscInt **global)
12 {
13   Mat_MPIAIJ     *aij = (Mat_MPIAIJ*)G->data;
14   PetscErrorCode ierr;
15   PetscBool      isMPIAIJ;
16 
17   PetscFunctionBegin;
18   ierr = PetscObjectTypeCompare((PetscObject)G, MATMPIAIJ, &isMPIAIJ ); CHKERRQ(ierr);
19   if (isMPIAIJ) {
20     if (gvec)ierr = VecDuplicate(aij->lvec,gvec);CHKERRQ(ierr);
21     if (global)*global = aij->garray;
22   } else {
23     /* no off-processor nodes */
24     if (gvec)*gvec = NULL;
25     if (global)*global = NULL;
26   }
27   PetscFunctionReturn(0);
28 }
29 
30 #undef __FUNCT__
31 #define __FUNCT__ "PCGAMGClassicalGraphSplitting_Private"
32 /*
33  Split the relevant graph into diagonal and off-diagonal parts in local numbering; for now this
34  a roundabout private interface to the mats' internal diag and offdiag mats.
35  */
36 PetscErrorCode PCGAMGClassicalGraphSplitting_Private(Mat G,Mat *Gd, Mat *Go)
37 {
38   Mat_MPIAIJ     *aij = (Mat_MPIAIJ*)G->data;
39   PetscErrorCode ierr;
40   PetscBool      isMPIAIJ;
41   PetscFunctionBegin;
42   ierr = PetscObjectTypeCompare((PetscObject)G, MATMPIAIJ, &isMPIAIJ ); CHKERRQ(ierr);
43   if (isMPIAIJ) {
44     *Gd = aij->A;
45     *Go = aij->B;
46   } else {
47     *Gd = G;
48     *Go = NULL;
49   }
50   PetscFunctionReturn(0);
51 }
52 
53 #undef __FUNCT__
54 #define __FUNCT__ "PCGAMGGraph_Classical"
55 PetscErrorCode PCGAMGGraph_Classical(PC pc,Mat A,Mat *G)
56 {
57   PetscInt          s,f,idx;
58   PetscInt          r,c,ncols,*gcol;
59   const PetscInt    *rcol;
60   const PetscScalar *rval;
61   PetscScalar       *gval;
62   PetscReal         rmax,Amax;
63   PetscInt          cmax;
64   PC_MG             *mg;
65   PC_GAMG           *gamg;
66   PetscErrorCode    ierr;
67   PetscInt          *gsparse,*lsparse;
68   Mat               lA,gA;
69   MatType           mtype;
70 
71   PetscFunctionBegin;
72   mg   = (PC_MG *)pc->data;
73   gamg = (PC_GAMG *)mg->innerctx;
74 
75   ierr = MatGetOwnershipRange(A,&s,&f);CHKERRQ(ierr);
76 
77   ierr = PetscMalloc(sizeof(PetscInt)*(f - s),&lsparse);CHKERRQ(ierr);
78   ierr = PetscMalloc(sizeof(PetscInt)*(f - s),&gsparse);CHKERRQ(ierr);
79 
80   ierr = PCGAMGClassicalGraphSplitting_Private(A,&lA,&gA);CHKERRQ(ierr);
81 
82   /* find the maximum off-diagonal entry in the matrix */
83   rmax = 0.;
84   cmax = 0;
85   for (r = s;r < f;r++) {
86     ierr = MatGetRow(A,r,&ncols,&rcol,&rval);CHKERRQ(ierr);
87     for (c = 0; c < ncols; c++) {
88       if (rcol[c] != r)
89         if (PetscAbsScalar(rval[c]) > rmax) rmax = PetscAbsScalar(rval[c]);
90     }
91     if (ncols > cmax) cmax = ncols;
92     ierr = MatRestoreRow(A,r,&ncols,&rcol,&rval);CHKERRQ(ierr);
93   }
94   ierr = PetscMalloc(sizeof(PetscScalar)*cmax,&gval);CHKERRQ(ierr);
95   ierr = PetscMalloc(sizeof(PetscInt)*cmax,&gcol);CHKERRQ(ierr);
96   ierr = MPI_Reduce(&rmax,&Amax,1,MPI_DOUBLE,MPI_MAX,0,PetscObjectComm((PetscObject)pc));
97 
98   ierr = PetscInfo1(pc,"Maximum off-diagonal value in classical AMG graph: %f",rmax);CHKERRQ(ierr);
99 
100   for (r = 0;r < f-s;r++) {
101     lsparse[r] = 0;
102     gsparse[r] = 0;
103   }
104 
105   /* for now this recreates the entire matrix due to a bug in MatCoarsen */
106   for (r = 0;r < f-s;r++) {
107     ierr = MatGetRow(lA,r,&ncols,&rcol,&rval);CHKERRQ(ierr);
108     idx = 0;
109     for (c = 0; c < ncols; c++) {
110       if (PetscAbsScalar(rval[c]) > gamg->threshold*Amax) {
111         idx++;
112       } else {
113         idx++;
114       }
115     }
116     ierr = MatRestoreRow(lA,r,&ncols,&rcol,&rval);CHKERRQ(ierr);
117     lsparse[r] = idx;
118   }
119   if (gA) {
120     for (r = 0;r < f-s;r++) {
121       ierr = MatGetRow(gA,r,&ncols,&rcol,&rval);CHKERRQ(ierr);
122       idx = 0;
123       for (c = 0; c < ncols; c++) {
124         if (PetscAbsScalar(rval[c]) > gamg->threshold*Amax) {
125           idx++;
126         } else {
127           idx++;
128         }
129       }
130       ierr = MatRestoreRow(gA,r,&ncols,&rcol,&rval);CHKERRQ(ierr);
131       gsparse[r] = idx;
132     }
133   }
134   ierr = MatCreate(PetscObjectComm((PetscObject)A),G); CHKERRQ(ierr);
135   ierr = MatGetType(A,&mtype);CHKERRQ(ierr);
136   ierr = MatSetType(*G,mtype);CHKERRQ(ierr);
137   ierr = MatSetSizes(*G,f-s,f-s,PETSC_DETERMINE,PETSC_DETERMINE);CHKERRQ(ierr);
138   ierr = MatMPIAIJSetPreallocation(*G,0,lsparse,0,gsparse);CHKERRQ(ierr);
139   ierr = MatSeqAIJSetPreallocation(*G,0,lsparse);CHKERRQ(ierr);
140   for (r = s;r < f;r++) {
141     ierr = MatGetRow(A,r,&ncols,&rcol,&rval);CHKERRQ(ierr);
142     idx = 0;
143     for (c = 0; c < ncols; c++) {
144       /* classical strength of connection */
145       if (PetscAbsScalar(rval[c]) > gamg->threshold*Amax) {
146         gcol[idx] = rcol[c];
147         gval[idx] = rval[c];
148         idx++;
149       } else {
150         gcol[idx] = rcol[c];
151         gval[idx] = 0.;
152         idx++;
153       }
154     }
155     ierr = MatSetValues(*G,1,&r,idx,gcol,gval,INSERT_VALUES);CHKERRQ(ierr);
156     ierr = MatRestoreRow(A,r,&ncols,&rcol,&rval);CHKERRQ(ierr);
157   }
158   ierr = MatAssemblyBegin(*G, MAT_FINAL_ASSEMBLY); CHKERRQ(ierr);
159   ierr = MatAssemblyEnd(*G, MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);
160 
161 
162 
163   ierr = PetscFree(gval);CHKERRQ(ierr);
164   ierr = PetscFree(gcol);CHKERRQ(ierr);
165   ierr = PetscFree(lsparse);CHKERRQ(ierr);
166   ierr = PetscFree(gsparse);CHKERRQ(ierr);
167   PetscFunctionReturn(0);
168 }
169 
170 
171 #undef __FUNCT__
172 #define __FUNCT__ "PCGAMGCoarsen_Classical"
173 PetscErrorCode PCGAMGCoarsen_Classical(PC pc,Mat *G,PetscCoarsenData **agg_lists)
174 {
175   PetscErrorCode   ierr;
176   MatCoarsen       crs;
177   MPI_Comm         fcomm = ((PetscObject)pc)->comm;
178 
179   PetscFunctionBegin;
180 
181 
182   /* construct the graph if necessary */
183   if (!G) {
184     SETERRQ(fcomm,PETSC_ERR_ARG_WRONGSTATE,"Must set Graph in PC in PCGAMG before coarsening");
185   }
186 
187   ierr = MatCoarsenCreate(fcomm,&crs);CHKERRQ(ierr);
188   ierr = MatCoarsenSetFromOptions(crs);CHKERRQ(ierr);
189   ierr = MatCoarsenSetAdjacency(crs,*G);CHKERRQ(ierr);
190   ierr = MatCoarsenSetStrictAggs(crs,PETSC_TRUE);CHKERRQ(ierr);
191   ierr = MatCoarsenApply(crs);CHKERRQ(ierr);
192   ierr = MatCoarsenGetData(crs,agg_lists);CHKERRQ(ierr);
193   ierr = MatCoarsenDestroy(&crs);CHKERRQ(ierr);
194 
195   PetscFunctionReturn(0);
196 }
197 
198 #undef __FUNCT__
199 #define __FUNCT__ "PCGAMGClassicalGhost_Private"
200 /*
201  Find all ghost nodes that are coarse and output the fine/coarse splitting for those as well
202 
203  Input:
204  G - graph;
205  gvec - Global Vector
206  avec - Local part of the scattered vec
207  bvec - Global part of the scattered vec
208 
209  Output:
210  findx - indirection t
211 
212  */
213 PetscErrorCode PCGAMGClassicalGhost_Private(Mat G,Vec v,Vec gv)
214 {
215   PetscErrorCode ierr;
216   Mat_MPIAIJ     *aij = (Mat_MPIAIJ*)G->data;
217   PetscBool      isMPIAIJ;
218 
219   PetscFunctionBegin;
220   ierr = PetscObjectTypeCompare((PetscObject)G, MATMPIAIJ, &isMPIAIJ ); CHKERRQ(ierr);
221   if (isMPIAIJ) {
222     ierr = VecScatterBegin(aij->Mvctx,v,gv,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr);
223     ierr = VecScatterEnd(aij->Mvctx,v,gv,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr);
224   }
225   PetscFunctionReturn(0);
226 }
227 
228 #undef __FUNCT__
229 #define __FUNCT__ "PCGAMGProlongator_Classical"
230 PetscErrorCode PCGAMGProlongator_Classical(PC pc, Mat A, Mat G, PetscCoarsenData *agg_lists,Mat *P)
231 {
232   PetscErrorCode    ierr;
233   MPI_Comm          comm;
234   Mat               lG,gG,lA,gA;     /* on and off diagonal matrices */
235   PetscInt          fn;                        /* fine local blocked sizes */
236   PetscInt          cn;                        /* coarse local blocked sizes */
237   PetscInt          gn;                        /* size of the off-diagonal fine vector */
238   PetscInt          fs,fe;                     /* fine (row) ownership range*/
239   PetscInt          cs,ce;                     /* coarse (column) ownership range */
240   PetscInt          i,j,k;                     /* indices! */
241   PetscBool         iscoarse;                  /* flag for determining if a node is coarse */
242   PetscInt          *lcid,*gcid;               /* on and off-processor coarse unknown IDs */
243   PetscInt          *lsparse,*gsparse;         /* on and off-processor sparsity patterns for prolongator */
244   PetscScalar       pij;
245   const PetscScalar *rval;
246   const PetscInt    *rcol;
247   PetscScalar       g_pos,g_neg,a_pos,a_neg,diag,invdiag,alpha,beta;
248   Vec               F;   /* vec of coarse size */
249   Vec               C;   /* vec of fine size */
250   Vec               gF;  /* vec of off-diagonal fine size */
251   MatType           mtype;
252   PetscInt          c_indx;
253   const PetscScalar *vcols;
254   const PetscInt    *icols;
255   PetscScalar       c_scalar;
256   PetscInt          ncols,col;
257   PetscInt          row_f,row_c;
258   PetscInt          cmax=0,ncolstotal,idx;
259   PetscScalar       *pvals;
260   PetscInt          *pcols;
261 
262   PetscFunctionBegin;
263   comm = ((PetscObject)pc)->comm;
264   ierr = MatGetOwnershipRange(A,&fs,&fe); CHKERRQ(ierr);
265   fn = (fe - fs);
266 
267   ierr = MatGetVecs(A,&F,NULL);CHKERRQ(ierr);
268 
269   /* get the number of local unknowns and the indices of the local unknowns */
270 
271   ierr = PetscMalloc(sizeof(PetscInt)*fn,&lsparse);CHKERRQ(ierr);
272   ierr = PetscMalloc(sizeof(PetscInt)*fn,&gsparse);CHKERRQ(ierr);
273   ierr = PetscMalloc(sizeof(PetscInt)*fn,&lcid);CHKERRQ(ierr);
274 
275   /* count the number of coarse unknowns */
276   cn = 0;
277   for (i=0;i<fn;i++) {
278     /* filter out singletons */
279     ierr = PetscCDEmptyAt(agg_lists,i,&iscoarse); CHKERRQ(ierr);
280     lcid[i] = -1;
281     if (!iscoarse) {
282       cn++;
283     }
284   }
285 
286    /* create the coarse vector */
287   ierr = VecCreateMPI(comm,cn,PETSC_DECIDE,&C);CHKERRQ(ierr);
288   ierr = VecGetOwnershipRange(C,&cs,&ce);CHKERRQ(ierr);
289 
290   /* construct a global vector indicating the global indices of the coarse unknowns */
291   cn = 0;
292   for (i=0;i<fn;i++) {
293     ierr = PetscCDEmptyAt(agg_lists,i,&iscoarse); CHKERRQ(ierr);
294     if (!iscoarse) {
295       lcid[i] = cs+cn;
296       cn++;
297     } else {
298       lcid[i] = -1;
299     }
300     c_scalar = (PetscScalar)lcid[i];
301     c_indx = fs+i;
302     ierr = VecSetValues(F,1,&c_indx,&c_scalar,INSERT_VALUES);CHKERRQ(ierr);
303   }
304 
305   ierr = VecAssemblyBegin(F);CHKERRQ(ierr);
306   ierr = VecAssemblyEnd(F);CHKERRQ(ierr);
307 
308   /* split the graph into two */
309   ierr = PCGAMGClassicalGraphSplitting_Private(G,&lG,&gG);CHKERRQ(ierr);
310   ierr = PCGAMGClassicalGraphSplitting_Private(A,&lA,&gA);CHKERRQ(ierr);
311 
312   /* scatter to the ghost vector */
313   ierr = PCGAMGClassicalCreateGhostVector_Private(G,&gF,NULL);CHKERRQ(ierr);
314   ierr = PCGAMGClassicalGhost_Private(G,F,gF);CHKERRQ(ierr);
315 
316   if (gG) {
317     ierr = VecGetSize(gF,&gn);CHKERRQ(ierr);
318     ierr = PetscMalloc(sizeof(PetscInt)*gn,&gcid);CHKERRQ(ierr);
319     for (i=0;i<gn;i++) {
320       ierr = VecGetValues(gF,1,&i,&c_scalar);CHKERRQ(ierr);
321       gcid[i] = (PetscInt)PetscRealPart(c_scalar);
322     }
323   }
324 
325   ierr = VecDestroy(&F);CHKERRQ(ierr);
326   ierr = VecDestroy(&gF);CHKERRQ(ierr);
327   ierr = VecDestroy(&C);CHKERRQ(ierr);
328 
329   /* count the on and off processor sparsity patterns for the prolongator */
330 
331   for (i=0;i<fn;i++) {
332     /* on */
333     ierr = MatGetRow(lG,i,&ncols,&icols,&vcols);CHKERRQ(ierr);
334     ncolstotal = ncols;
335     lsparse[i] = 0;
336     if (lcid[i] >= 0) {
337       lsparse[i] = 1;
338       gsparse[i] = 0;
339     } else {
340       for (j = 0;j < ncols;j++) {
341         col = icols[j];
342         if (lcid[col] >= 0 && vcols[j] != 0.) {
343           lsparse[i] += 1;
344         }
345       }
346 
347       ierr = MatRestoreRow(lG,i,&ncols,&icols,&vcols);CHKERRQ(ierr);
348       ncolstotal += ncols;
349       /* off */
350       gsparse[i] = 0;
351       if (gG) {
352         ierr = MatGetRow(gG,i,&ncols,&icols,&vcols);CHKERRQ(ierr);
353         for (j = 0; j < ncols; j++) {
354           col = icols[j];
355           if (gcid[col] >= 0 && vcols[j] != 0.) {
356             gsparse[i] += 1;
357           }
358         }
359         ierr = MatRestoreRow(gG,i,&ncols,NULL,NULL);CHKERRQ(ierr);
360       }
361       if (ncolstotal > cmax) cmax = ncolstotal;
362     }
363   }
364 
365   ierr = PetscMalloc(sizeof(PetscInt)*cmax,&pcols);CHKERRQ(ierr);
366   ierr = PetscMalloc(sizeof(PetscScalar)*cmax,&pvals);CHKERRQ(ierr);
367 
368   /* preallocate and create the prolongator */
369   ierr = MatCreate(comm,P); CHKERRQ(ierr);
370   ierr = MatGetType(G,&mtype);CHKERRQ(ierr);
371   ierr = MatSetType(*P,mtype);CHKERRQ(ierr);
372 
373   ierr = MatSetSizes(*P,fn,cn,PETSC_DETERMINE,PETSC_DETERMINE);CHKERRQ(ierr);
374   ierr = MatMPIAIJSetPreallocation(*P,0,lsparse,0,gsparse);CHKERRQ(ierr);
375   ierr = MatSeqAIJSetPreallocation(*P,0,lsparse);CHKERRQ(ierr);
376 
377   /* loop over local fine nodes -- get the diagonal, the sum of positive and negative strong and weak weights, and set up the row */
378   for (i = 0;i < fn;i++) {
379     /* determine on or off */
380     row_f = i + fs;
381     row_c = lcid[i];
382     if (row_c >= 0) {
383       pij = 1.;
384       ierr = MatSetValues(*P,1,&row_f,1,&row_c,&pij,INSERT_VALUES);CHKERRQ(ierr);
385     } else {
386       g_pos = 0.;
387       g_neg = 0.;
388       a_pos = 0.;
389       a_neg = 0.;
390       diag = 0.;
391 
392       /* local strong connections */
393       ierr = MatGetRow(lG,i,&ncols,&rcol,&rval);CHKERRQ(ierr);
394       for (k = 0; k < ncols; k++) {
395         if (PetscRealPart(rval[k]) > 0) {
396           g_pos += rval[k];
397         } else {
398           g_neg += rval[k];
399         }
400       }
401       ierr = MatRestoreRow(lG,i,&ncols,&rcol,&rval);CHKERRQ(ierr);
402 
403 
404       /* ghosted strong connections */
405       if (gG) {
406         ierr = MatGetRow(gG,i,&ncols,&rcol,&rval);CHKERRQ(ierr);
407         for (k = 0; k < ncols; k++) {
408           if (PetscRealPart(rval[k]) > 0.) {
409             g_pos += rval[k];
410           } else {
411             g_neg += rval[k];
412           }
413         }
414         ierr = MatRestoreRow(gG,i,&ncols,&rcol,&rval);CHKERRQ(ierr);
415       }
416 
417       /* local all connections */
418       ierr = MatGetRow(lA,i,&ncols,&rcol,&rval);CHKERRQ(ierr);
419       for (k = 0; k < ncols; k++) {
420         if (PetscRealPart(rval[k]) > 0) {
421           a_pos += rval[k];
422         } else {
423           a_neg += rval[k];
424         }
425         if (rcol[k] == i) {
426           diag = rval[k];
427         }
428       }
429       ierr = MatRestoreRow(lA,i,&ncols,&rcol,&rval);CHKERRQ(ierr);
430 
431       /* ghosted all connections */
432       if (gA) {
433         ierr = MatGetRow(gA,i,&ncols,&rcol,&rval);CHKERRQ(ierr);
434         for (k = 0; k < ncols; k++) {
435           if (PetscRealPart(rval[k]) > 0.) {
436             a_pos += PetscRealPart(rval[k]);
437           } else {
438             a_neg += PetscRealPart(rval[k]);
439           }
440         }
441         ierr = MatRestoreRow(gA,i,&ncols,&rcol,&rval);CHKERRQ(ierr);
442       }
443 
444       if (g_neg == 0.) {
445         alpha = 0.;
446       } else {
447         alpha = -a_neg/g_neg;
448       }
449 
450       if (g_pos == 0.) {
451         diag += a_pos;
452         beta = 0.;
453       } else {
454         beta = -a_pos/g_pos;
455       }
456 
457       invdiag = 1. / diag;
458       /* on */
459       ierr = MatGetRow(lG,i,&ncols,&icols,&vcols);CHKERRQ(ierr);
460       idx = 0;
461       for (j = 0;j < ncols;j++) {
462         col = icols[j];
463         if (lcid[col] >= 0 && vcols[j] != 0.) {
464           row_f = i + fs;
465           row_c = lcid[col];
466           /* set the values for on-processor ones */
467           if (PetscRealPart(vcols[j]) < 0.) {
468             pij = vcols[j]*alpha*invdiag;
469           } else {
470             pij = vcols[j]*beta*invdiag;
471           }
472           if (PetscAbsScalar(pij) != 0.) {
473             pvals[idx] = pij;
474             pcols[idx] = row_c;
475             idx++;
476           }
477         }
478       }
479       ierr = MatRestoreRow(lG,i,&ncols,&icols,&vcols);CHKERRQ(ierr);
480       /* off */
481       if (gG) {
482         ierr = MatGetRow(gG,i,&ncols,&icols,&vcols);CHKERRQ(ierr);
483         for (j = 0; j < ncols; j++) {
484           col = icols[j];
485           if (gcid[col] >= 0 && vcols[j] != 0.) {
486             row_f = i + fs;
487             row_c = gcid[col];
488             /* set the values for on-processor ones */
489             if (PetscRealPart(vcols[j]) < 0.) {
490               pij = vcols[j]*alpha*invdiag;
491             } else {
492               pij = vcols[j]*beta*invdiag;
493             }
494             if (PetscAbsScalar(pij) != 0.) {
495               pvals[idx] = pij;
496               pcols[idx] = row_c;
497               idx++;
498             }
499           }
500         }
501         ierr = MatSetValues(*P,1,&row_f,idx,pcols,pvals,INSERT_VALUES);CHKERRQ(ierr);
502         ierr = MatRestoreRow(gG,i,&ncols,NULL,NULL);CHKERRQ(ierr);
503       }
504     }
505   }
506   ierr = MatAssemblyBegin(*P, MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);
507   ierr = MatAssemblyEnd(*P, MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);
508 
509   ierr = PetscFree(lsparse);CHKERRQ(ierr);
510   ierr = PetscFree(gsparse);CHKERRQ(ierr);
511   ierr = PetscFree(pcols);CHKERRQ(ierr);
512   ierr = PetscFree(pvals);CHKERRQ(ierr);
513   ierr = PetscFree(lcid);CHKERRQ(ierr);
514   if (gG) {ierr = PetscFree(gcid);CHKERRQ(ierr);}
515 
516   PetscFunctionReturn(0);
517 }
518 
519 #undef __FUNCT__
520 #define __FUNCT__ "PCGAMGDestroy_Classical"
521 PetscErrorCode PCGAMGDestroy_Classical(PC pc)
522 {
523   PetscErrorCode ierr;
524   PC_MG          *mg          = (PC_MG*)pc->data;
525   PC_GAMG        *pc_gamg     = (PC_GAMG*)mg->innerctx;
526 
527   PetscFunctionBegin;
528   ierr = PetscFree(pc_gamg->subctx);CHKERRQ(ierr);
529   PetscFunctionReturn(0);
530 }
531 
532 #undef __FUNCT__
533 #define __FUNCT__ "PCGAMGSetFromOptions_Classical"
534 PetscErrorCode PCGAMGSetFromOptions_Classical(PC pc)
535 {
536   PetscFunctionBegin;
537   PetscFunctionReturn(0);
538 }
539 
540 #undef __FUNCT__
541 #define __FUNCT__ "PCGAMGSetData_Classical"
542 PetscErrorCode PCGAMGSetData_Classical(PC pc, Mat A)
543 {
544   PC_MG          *mg      = (PC_MG*)pc->data;
545   PC_GAMG        *pc_gamg = (PC_GAMG*)mg->innerctx;
546 
547   PetscFunctionBegin;
548   /* no data for classical AMG */
549   pc_gamg->data           = NULL;
550   pc_gamg->data_cell_cols = 1;
551   pc_gamg->data_cell_rows = 1;
552   pc_gamg->data_sz = 0;
553   PetscFunctionReturn(0);
554 }
555 
556 /* -------------------------------------------------------------------------- */
557 /*
558    PCCreateGAMG_Classical
559 
560 */
561 #undef __FUNCT__
562 #define __FUNCT__ "PCCreateGAMG_Classical"
563 PetscErrorCode  PCCreateGAMG_Classical(PC pc)
564 {
565   PetscErrorCode ierr;
566   PC_MG             *mg      = (PC_MG*)pc->data;
567   PC_GAMG           *pc_gamg = (PC_GAMG*)mg->innerctx;
568   PC_GAMG_Classical *pc_gamg_classical;
569 
570   PetscFunctionBegin;
571   if (pc_gamg->subctx) {
572     /* call base class */
573     ierr = PCDestroy_GAMG(pc);CHKERRQ(ierr);
574   }
575 
576   /* create sub context for SA */
577   ierr = PetscNewLog(pc, PC_GAMG_Classical, &pc_gamg_classical);CHKERRQ(ierr);
578   pc_gamg->subctx = pc_gamg_classical;
579   pc->ops->setfromoptions = PCGAMGSetFromOptions_Classical;
580   /* reset does not do anything; setup not virtual */
581 
582   /* set internal function pointers */
583   pc_gamg->ops->destroy     = PCGAMGDestroy_Classical;
584   pc_gamg->ops->graph       = PCGAMGGraph_Classical;
585   pc_gamg->ops->coarsen     = PCGAMGCoarsen_Classical;
586   pc_gamg->ops->prolongator = PCGAMGProlongator_Classical;
587   pc_gamg->ops->optprol     = NULL;
588   pc_gamg->ops->formkktprol = NULL;
589 
590   pc_gamg->ops->createdefaultdata = PCGAMGSetData_Classical;
591   PetscFunctionReturn(0);
592 }
593