xref: /petsc/src/ksp/pc/impls/gamg/classical.c (revision 3c9ab2c330197cde143110315a271226c1b750df) 
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 = PCGAMGClassicalGraphSplitting_Private(A,&lA,&gA);CHKERRQ(ierr);
78 
79   ierr = PetscMalloc(sizeof(PetscInt)*(f - s),&lsparse);CHKERRQ(ierr);
80   if (gA) {ierr = PetscMalloc(sizeof(PetscInt)*(f - s),&gsparse);CHKERRQ(ierr);}
81   else {
82     gsparse = NULL;
83   }
84 
85   /* find the maximum off-diagonal entry in the matrix */
86   rmax = 0.;
87   cmax = 0;
88   for (r = s;r < f;r++) {
89     ierr = MatGetRow(A,r,&ncols,&rcol,&rval);CHKERRQ(ierr);
90     for (c = 0; c < ncols; c++) {
91       if (rcol[c] != r)
92         if (PetscAbsScalar(rval[c]) > rmax) rmax = PetscAbsScalar(rval[c]);
93     }
94     if (ncols > cmax) cmax = ncols;
95     ierr = MatRestoreRow(A,r,&ncols,&rcol,&rval);CHKERRQ(ierr);
96   }
97   ierr = PetscMalloc(sizeof(PetscScalar)*cmax,&gval);CHKERRQ(ierr);
98   ierr = PetscMalloc(sizeof(PetscInt)*cmax,&gcol);CHKERRQ(ierr);
99   ierr = MPI_Reduce(&rmax,&Amax,1,MPI_DOUBLE,MPI_MAX,0,PetscObjectComm((PetscObject)pc));
100 
101   ierr = PetscInfo1(pc,"Maximum off-diagonal value in classical AMG graph: %f",rmax);CHKERRQ(ierr);
102 
103   for (r = 0;r < f-s;r++) {
104     lsparse[r] = 0;
105     if (gsparse) gsparse[r] = 0;
106   }
107 
108   /* for now this recreates the entire matrix due to a bug in MatCoarsen */
109   for (r = 0;r < f-s;r++) {
110     ierr = MatGetRow(lA,r,&ncols,&rcol,&rval);CHKERRQ(ierr);
111     idx = 0;
112     for (c = 0; c < ncols; c++) {
113       if (PetscAbsScalar(rval[c]) > gamg->threshold*Amax) {
114         idx++;
115       } else {
116         idx++;
117       }
118     }
119     ierr = MatRestoreRow(lA,r,&ncols,&rcol,&rval);CHKERRQ(ierr);
120     lsparse[r] = idx;
121   }
122   if (gA) {
123     for (r = 0;r < f-s;r++) {
124       ierr = MatGetRow(gA,r,&ncols,&rcol,&rval);CHKERRQ(ierr);
125       idx = 0;
126       for (c = 0; c < ncols; c++) {
127         if (PetscAbsScalar(rval[c]) > gamg->threshold*Amax) {
128           idx++;
129         } else {
130           idx++;
131         }
132       }
133       ierr = MatRestoreRow(gA,r,&ncols,&rcol,&rval);CHKERRQ(ierr);
134       gsparse[r] = idx;
135     }
136   }
137   ierr = MatCreate(PetscObjectComm((PetscObject)A),G); CHKERRQ(ierr);
138   ierr = MatGetType(A,&mtype);CHKERRQ(ierr);
139   ierr = MatSetType(*G,mtype);CHKERRQ(ierr);
140   ierr = MatSetSizes(*G,f-s,f-s,PETSC_DETERMINE,PETSC_DETERMINE);CHKERRQ(ierr);
141   ierr = MatMPIAIJSetPreallocation(*G,0,lsparse,0,gsparse);CHKERRQ(ierr);
142   ierr = MatSeqAIJSetPreallocation(*G,0,lsparse);CHKERRQ(ierr);
143   for (r = s;r < f;r++) {
144     ierr = MatGetRow(A,r,&ncols,&rcol,&rval);CHKERRQ(ierr);
145     idx = 0;
146     for (c = 0; c < ncols; c++) {
147       /* classical strength of connection */
148       if (PetscAbsScalar(rval[c]) > gamg->threshold*Amax) {
149         gcol[idx] = rcol[c];
150         gval[idx] = rval[c];
151         idx++;
152       } else {
153         gcol[idx] = rcol[c];
154         gval[idx] = 0.;
155         idx++;
156       }
157     }
158     ierr = MatSetValues(*G,1,&r,idx,gcol,gval,INSERT_VALUES);CHKERRQ(ierr);
159     ierr = MatRestoreRow(A,r,&ncols,&rcol,&rval);CHKERRQ(ierr);
160   }
161   ierr = MatAssemblyBegin(*G, MAT_FINAL_ASSEMBLY); CHKERRQ(ierr);
162   ierr = MatAssemblyEnd(*G, MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);
163 
164   ierr = PetscFree(gval);CHKERRQ(ierr);
165   ierr = PetscFree(gcol);CHKERRQ(ierr);
166   ierr = PetscFree(lsparse);CHKERRQ(ierr);
167   ierr = PetscFree(gsparse);CHKERRQ(ierr);
168   PetscFunctionReturn(0);
169 }
170 
171 
172 #undef __FUNCT__
173 #define __FUNCT__ "PCGAMGCoarsen_Classical"
174 PetscErrorCode PCGAMGCoarsen_Classical(PC pc,Mat *G,PetscCoarsenData **agg_lists)
175 {
176   PetscErrorCode   ierr;
177   MatCoarsen       crs;
178   MPI_Comm         fcomm = ((PetscObject)pc)->comm;
179 
180   PetscFunctionBegin;
181 
182 
183   /* construct the graph if necessary */
184   if (!G) {
185     SETERRQ(fcomm,PETSC_ERR_ARG_WRONGSTATE,"Must set Graph in PC in PCGAMG before coarsening");
186   }
187 
188   ierr = MatCoarsenCreate(fcomm,&crs);CHKERRQ(ierr);
189   ierr = MatCoarsenSetFromOptions(crs);CHKERRQ(ierr);
190   ierr = MatCoarsenSetAdjacency(crs,*G);CHKERRQ(ierr);
191   ierr = MatCoarsenSetStrictAggs(crs,PETSC_TRUE);CHKERRQ(ierr);
192   ierr = MatCoarsenApply(crs);CHKERRQ(ierr);
193   ierr = MatCoarsenGetData(crs,agg_lists);CHKERRQ(ierr);
194   ierr = MatCoarsenDestroy(&crs);CHKERRQ(ierr);
195 
196   PetscFunctionReturn(0);
197 }
198 
199 #undef __FUNCT__
200 #define __FUNCT__ "PCGAMGClassicalGhost_Private"
201 /*
202  Find all ghost nodes that are coarse and output the fine/coarse splitting for those as well
203 
204  Input:
205  G - graph;
206  gvec - Global Vector
207  avec - Local part of the scattered vec
208  bvec - Global part of the scattered vec
209 
210  Output:
211  findx - indirection t
212 
213  */
214 PetscErrorCode PCGAMGClassicalGhost_Private(Mat G,Vec v,Vec gv)
215 {
216   PetscErrorCode ierr;
217   Mat_MPIAIJ     *aij = (Mat_MPIAIJ*)G->data;
218   PetscBool      isMPIAIJ;
219 
220   PetscFunctionBegin;
221   ierr = PetscObjectTypeCompare((PetscObject)G, MATMPIAIJ, &isMPIAIJ ); CHKERRQ(ierr);
222   if (isMPIAIJ) {
223     ierr = VecScatterBegin(aij->Mvctx,v,gv,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr);
224     ierr = VecScatterEnd(aij->Mvctx,v,gv,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr);
225   }
226   PetscFunctionReturn(0);
227 }
228 
229 #undef __FUNCT__
230 #define __FUNCT__ "PCGAMGProlongator_Classical"
231 PetscErrorCode PCGAMGProlongator_Classical(PC pc, Mat A, Mat G, PetscCoarsenData *agg_lists,Mat *P)
232 {
233   PetscErrorCode    ierr;
234   MPI_Comm          comm;
235   Mat               lG,gG,lA,gA;     /* on and off diagonal matrices */
236   PetscInt          fn;                        /* fine local blocked sizes */
237   PetscInt          cn;                        /* coarse local blocked sizes */
238   PetscInt          gn;                        /* size of the off-diagonal fine vector */
239   PetscInt          fs,fe;                     /* fine (row) ownership range*/
240   PetscInt          cs,ce;                     /* coarse (column) ownership range */
241   PetscInt          i,j,k;                     /* indices! */
242   PetscBool         iscoarse;                  /* flag for determining if a node is coarse */
243   PetscInt          *lcid,*gcid;               /* on and off-processor coarse unknown IDs */
244   PetscInt          *lsparse,*gsparse;         /* on and off-processor sparsity patterns for prolongator */
245   PetscScalar       pij;
246   const PetscScalar *rval;
247   const PetscInt    *rcol;
248   PetscScalar       g_pos,g_neg,a_pos,a_neg,diag,invdiag,alpha,beta;
249   Vec               F;   /* vec of coarse size */
250   Vec               C;   /* vec of fine size */
251   Vec               gF;  /* vec of off-diagonal fine size */
252   MatType           mtype;
253   PetscInt          c_indx;
254   const PetscScalar *vcols;
255   const PetscInt    *icols;
256   PetscScalar       c_scalar;
257   PetscInt          ncols,col;
258   PetscInt          row_f,row_c;
259   PetscInt          cmax=0,ncolstotal,idx;
260   PetscScalar       *pvals;
261   PetscInt          *pcols;
262 
263   PetscFunctionBegin;
264   comm = ((PetscObject)pc)->comm;
265   ierr = MatGetOwnershipRange(A,&fs,&fe); CHKERRQ(ierr);
266   fn = (fe - fs);
267 
268   ierr = MatGetVecs(A,&F,NULL);CHKERRQ(ierr);
269 
270   /* get the number of local unknowns and the indices of the local unknowns */
271 
272   ierr = PetscMalloc(sizeof(PetscInt)*fn,&lsparse);CHKERRQ(ierr);
273   ierr = PetscMalloc(sizeof(PetscInt)*fn,&gsparse);CHKERRQ(ierr);
274   ierr = PetscMalloc(sizeof(PetscInt)*fn,&lcid);CHKERRQ(ierr);
275 
276   /* count the number of coarse unknowns */
277   cn = 0;
278   for (i=0;i<fn;i++) {
279     /* filter out singletons */
280     ierr = PetscCDEmptyAt(agg_lists,i,&iscoarse); CHKERRQ(ierr);
281     lcid[i] = -1;
282     if (!iscoarse) {
283       cn++;
284     }
285   }
286 
287    /* create the coarse vector */
288   ierr = VecCreateMPI(comm,cn,PETSC_DECIDE,&C);CHKERRQ(ierr);
289   ierr = VecGetOwnershipRange(C,&cs,&ce);CHKERRQ(ierr);
290 
291   /* construct a global vector indicating the global indices of the coarse unknowns */
292   cn = 0;
293   for (i=0;i<fn;i++) {
294     ierr = PetscCDEmptyAt(agg_lists,i,&iscoarse); CHKERRQ(ierr);
295     if (!iscoarse) {
296       lcid[i] = cs+cn;
297       cn++;
298     } else {
299       lcid[i] = -1;
300     }
301     c_scalar = (PetscScalar)lcid[i];
302     c_indx = fs+i;
303     ierr = VecSetValues(F,1,&c_indx,&c_scalar,INSERT_VALUES);CHKERRQ(ierr);
304   }
305 
306   ierr = VecAssemblyBegin(F);CHKERRQ(ierr);
307   ierr = VecAssemblyEnd(F);CHKERRQ(ierr);
308 
309   /* split the graph into two */
310   ierr = PCGAMGClassicalGraphSplitting_Private(G,&lG,&gG);CHKERRQ(ierr);
311   ierr = PCGAMGClassicalGraphSplitting_Private(A,&lA,&gA);CHKERRQ(ierr);
312 
313   /* scatter to the ghost vector */
314   ierr = PCGAMGClassicalCreateGhostVector_Private(G,&gF,NULL);CHKERRQ(ierr);
315   ierr = PCGAMGClassicalGhost_Private(G,F,gF);CHKERRQ(ierr);
316 
317   if (gG) {
318     ierr = VecGetSize(gF,&gn);CHKERRQ(ierr);
319     ierr = PetscMalloc(sizeof(PetscInt)*gn,&gcid);CHKERRQ(ierr);
320     for (i=0;i<gn;i++) {
321       ierr = VecGetValues(gF,1,&i,&c_scalar);CHKERRQ(ierr);
322       gcid[i] = (PetscInt)PetscRealPart(c_scalar);
323     }
324   }
325 
326   ierr = VecDestroy(&F);CHKERRQ(ierr);
327   ierr = VecDestroy(&gF);CHKERRQ(ierr);
328   ierr = VecDestroy(&C);CHKERRQ(ierr);
329 
330   /* count the on and off processor sparsity patterns for the prolongator */
331 
332   for (i=0;i<fn;i++) {
333     /* on */
334     ierr = MatGetRow(lG,i,&ncols,&icols,&vcols);CHKERRQ(ierr);
335     ncolstotal = ncols;
336     lsparse[i] = 0;
337     if (lcid[i] >= 0) {
338       lsparse[i] = 1;
339       gsparse[i] = 0;
340     } else {
341       for (j = 0;j < ncols;j++) {
342         col = icols[j];
343         if (lcid[col] >= 0 && vcols[j] != 0.) {
344           lsparse[i] += 1;
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       /* ghosted strong connections */
404       if (gG) {
405         ierr = MatGetRow(gG,i,&ncols,&rcol,&rval);CHKERRQ(ierr);
406         for (k = 0; k < ncols; k++) {
407           if (PetscRealPart(rval[k]) > 0.) {
408             g_pos += rval[k];
409           } else {
410             g_neg += rval[k];
411           }
412         }
413         ierr = MatRestoreRow(gG,i,&ncols,&rcol,&rval);CHKERRQ(ierr);
414       }
415 
416       /* local all connections */
417       ierr = MatGetRow(lA,i,&ncols,&rcol,&rval);CHKERRQ(ierr);
418       for (k = 0; k < ncols; k++) {
419         if (PetscRealPart(rval[k]) > 0) {
420           a_pos += rval[k];
421         } else {
422           a_neg += rval[k];
423         }
424         if (rcol[k] == i) {
425           diag = rval[k];
426         }
427       }
428       ierr = MatRestoreRow(lA,i,&ncols,&rcol,&rval);CHKERRQ(ierr);
429 
430       /* ghosted all connections */
431       if (gA) {
432         ierr = MatGetRow(gA,i,&ncols,&rcol,&rval);CHKERRQ(ierr);
433         for (k = 0; k < ncols; k++) {
434           if (PetscRealPart(rval[k]) > 0.) {
435             a_pos += PetscRealPart(rval[k]);
436           } else {
437             a_neg += PetscRealPart(rval[k]);
438           }
439         }
440         ierr = MatRestoreRow(gA,i,&ncols,&rcol,&rval);CHKERRQ(ierr);
441       }
442 
443       if (g_neg == 0.) {
444         alpha = 0.;
445       } else {
446         alpha = -a_neg/g_neg;
447       }
448 
449       if (g_pos == 0.) {
450         diag += a_pos;
451         beta = 0.;
452       } else {
453         beta = -a_pos/g_pos;
454       }
455 
456       invdiag = 1. / diag;
457       /* on */
458       ierr = MatGetRow(lG,i,&ncols,&icols,&vcols);CHKERRQ(ierr);
459       idx = 0;
460       for (j = 0;j < ncols;j++) {
461         col = icols[j];
462         if (lcid[col] >= 0 && vcols[j] != 0.) {
463           row_f = i + fs;
464           row_c = lcid[col];
465           /* set the values for on-processor ones */
466           if (PetscRealPart(vcols[j]) < 0.) {
467             pij = vcols[j]*alpha*invdiag;
468           } else {
469             pij = vcols[j]*beta*invdiag;
470           }
471           if (PetscAbsScalar(pij) != 0.) {
472             pvals[idx] = pij;
473             pcols[idx] = row_c;
474             idx++;
475           }
476         }
477       }
478       ierr = MatRestoreRow(lG,i,&ncols,&icols,&vcols);CHKERRQ(ierr);
479       /* off */
480       if (gG) {
481         ierr = MatGetRow(gG,i,&ncols,&icols,&vcols);CHKERRQ(ierr);
482         for (j = 0; j < ncols; j++) {
483           col = icols[j];
484           if (gcid[col] >= 0 && vcols[j] != 0.) {
485             row_f = i + fs;
486             row_c = gcid[col];
487             /* set the values for on-processor ones */
488             if (PetscRealPart(vcols[j]) < 0.) {
489               pij = vcols[j]*alpha*invdiag;
490             } else {
491               pij = vcols[j]*beta*invdiag;
492             }
493             if (PetscAbsScalar(pij) != 0.) {
494               pvals[idx] = pij;
495               pcols[idx] = row_c;
496               idx++;
497             }
498           }
499         }
500         ierr = MatRestoreRow(gG,i,&ncols,&icols,&vcols);CHKERRQ(ierr);
501       }
502       ierr = MatSetValues(*P,1,&row_f,idx,pcols,pvals,INSERT_VALUES);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