Actual source code: mg.c

  1: /*$Id: mg.c,v 1.122 2001/03/23 23:23:10 balay Exp $*/
  2: /*
  3:     Defines the multigrid preconditioner interface.
  4: */
  5: #include "src/sles/pc/impls/mg/mgimpl.h"                    /*I "petscmg.h" I*/


  8: /*
  9:        MGMCycle_Private - Given an MG structure created with MGCreate() runs 
 10:                   one multiplicative cycle down through the levels and
 11:                   back up.

 13:     Input Parameter:
 14: .   mg - structure created with  MGCreate().
 15: */
 16: int MGMCycle_Private(MG *mglevels)
 17: {
 18:   MG     mg = *mglevels,mgc = *(mglevels - 1);
 19:   int    cycles = mg->cycles,ierr,its;
 20:   Scalar zero = 0.0;

 23:   if (!mg->level) {  /* coarse grid */
 24:     SLESSolve(mg->smoothd,mg->b,mg->x,&its);
 25:   } else {
 26:     while (cycles--) {
 27:       SLESSolve(mg->smoothd,mg->b,mg->x,&its);
 28:       (*mg->residual)(mg->A,mg->b,mg->x,mg->r);
 29:       MatRestrict(mg->restrct,mg->r,mgc->b);
 30:       VecSet(&zero,mgc->x);
 31:       MGMCycle_Private(mglevels-1);
 32:       MatInterpolateAdd(mg->interpolate,mgc->x,mg->x,mg->x);
 33:       SLESSolve(mg->smoothu,mg->b,mg->x,&its);
 34:     }
 35:   }
 36:   return(0);
 37: }

 39: /*
 40:        MGCreate_Private - Creates a MG structure for use with the
 41:                multigrid code. Level 0 is the coarsest. (But the 
 42:                finest level is stored first in the array).

 44: */
 45: static int MGCreate_Private(MPI_Comm comm,int levels,PC pc,MPI_Comm *comms,MG **result)
 46: {
 47:   MG   *mg;
 48:   int  i,ierr,size;
 49:   char *prefix;
 50:   KSP  ksp;
 51:   PC   ipc;

 54:   PetscMalloc(levels*sizeof(MG),&mg);
 55:   PetscLogObjectMemory(pc,levels*(sizeof(MG)+sizeof(struct _MG)));

 57:   PCGetOptionsPrefix(pc,&prefix);

 59:   for (i=0; i<levels; i++) {
 60:     PetscNew(struct _MG,&mg[i]);
 61:     PetscMemzero(mg[i],sizeof(struct _MG));
 62:     mg[i]->level  = i;
 63:     mg[i]->levels = levels;
 64:     mg[i]->cycles = 1;

 66:     if (comms) comm = comms[i];
 67:     SLESCreate(comm,&mg[i]->smoothd);
 68:     SLESGetKSP(mg[i]->smoothd,&ksp);
 69:     KSPSetTolerances(ksp,PETSC_DEFAULT,PETSC_DEFAULT,PETSC_DEFAULT,1);
 70:     SLESSetOptionsPrefix(mg[i]->smoothd,prefix);

 72:     /* do special stuff for coarse grid */
 73:     if (!i && levels > 1) {
 74:       SLESAppendOptionsPrefix(mg[0]->smoothd,"mg_coarse_");

 76:       /* coarse solve is (redundant) LU by default */
 77:       SLESGetKSP(mg[0]->smoothd,&ksp);
 78:       KSPSetType(ksp,KSPPREONLY);
 79:       SLESGetPC(mg[0]->smoothd,&ipc);
 80:       MPI_Comm_size(comm,&size);
 81:       if (size > 1) {
 82:         PCSetType(ipc,PCREDUNDANT);
 83:         PCRedundantGetPC(ipc,&ipc);
 84:       }
 85:       PCSetType(ipc,PCLU);

 87:     } else {
 88:       SLESAppendOptionsPrefix(mg[i]->smoothd,"mg_levels_");
 89:     }
 90:     PetscLogObjectParent(pc,mg[i]->smoothd);
 91:     mg[i]->smoothu         = mg[i]->smoothd;
 92:     mg[i]->default_smoothu = 10000;
 93:     mg[i]->default_smoothd = 10000;
 94:   }
 95:   *result = mg;
 96:   return(0);
 97: }

 99: static int PCDestroy_MG(PC pc)
100: {
101:   MG  *mg = (MG*)pc->data;
102:   int i,n = mg[0]->levels,ierr;

105:   for (i=0; i<n; i++) {
106:     if (mg[i]->smoothd != mg[i]->smoothu) {
107:       SLESDestroy(mg[i]->smoothd);
108:     }
109:     SLESDestroy(mg[i]->smoothu);
110:     PetscFree(mg[i]);
111:   }
112:   PetscFree(mg);
113:   return(0);
114: }



118: EXTERN int MGACycle_Private(MG*);
119: EXTERN int MGFCycle_Private(MG*);
120: EXTERN int MGKCycle_Private(MG*);

122: /*
123:    PCApply_MG - Runs either an additive, multiplicative, Kaskadic
124:              or full cycle of multigrid. 

126:   Note: 
127:   A simple wrapper which calls MGMCycle(),MGACycle(), or MGFCycle(). 
128: */
129: static int PCApply_MG(PC pc,Vec b,Vec x)
130: {
131:   MG     *mg = (MG*)pc->data;
132:   Scalar zero = 0.0;
133:   int    levels = mg[0]->levels,ierr;

136:   mg[levels-1]->b = b;
137:   mg[levels-1]->x = x;
138:   if (mg[0]->am == MGMULTIPLICATIVE) {
139:     VecSet(&zero,x);
140:     MGMCycle_Private(mg+levels-1);
141:   }
142:   else if (mg[0]->am == MGADDITIVE) {
143:     MGACycle_Private(mg);
144:   }
145:   else if (mg[0]->am == MGKASKADE) {
146:     MGKCycle_Private(mg);
147:   }
148:   else {
149:     MGFCycle_Private(mg);
150:   }
151:   return(0);
152: }

154: static int PCApplyRichardson_MG(PC pc,Vec b,Vec x,Vec w,int its)
155: {
156:   MG  *mg = (MG*)pc->data;
157:   int ierr,levels = mg[0]->levels;

160:   mg[levels-1]->b = b;
161:   mg[levels-1]->x = x;
162:   while (its--) {
163:     MGMCycle_Private(mg+levels-1);
164:   }
165:   return(0);
166: }

168: static int PCSetFromOptions_MG(PC pc)
169: {
170:   int        ierr,m,levels = 1;
171:   PetscTruth flg;
172:   char       buff[16],*type[] = {"additive","multiplicative","full","cascade"};


176:   PetscOptionsHead("Multigrid options");
177:     if (!pc->data) {
178:       PetscOptionsInt("-pc_mg_levels","Number of Levels","MGSetLevels",levels,&levels,&flg);
179:       MGSetLevels(pc,levels,PETSC_NULL);
180:     }
181:     PetscOptionsInt("-pc_mg_cycles","1 for V cycle, 2 for W-cycle","MGSetCycles",1,&m,&flg);
182:     if (flg) {
183:       MGSetCycles(pc,m);
184:     }
185:     PetscOptionsInt("-pc_mg_smoothup","Number of post-smoothing steps","MGSetNumberSmoothUp",1,&m,&flg);
186:     if (flg) {
187:       MGSetNumberSmoothUp(pc,m);
188:     }
189:     PetscOptionsInt("-pc_mg_smoothdown","Number of pre-smoothing steps","MGSetNumberSmoothDown",1,&m,&flg);
190:     if (flg) {
191:       MGSetNumberSmoothDown(pc,m);
192:     }
193:     PetscOptionsEList("-pc_mg_type","Multigrid type","MGSetType",type,4,"multiplicative",buff,15,&flg);
194:     if (flg) {
195:       MGType     mg;
196:       PetscTruth isadd,ismult,isfull,iskask,iscasc;

198:       PetscStrcmp(buff,type[0],&isadd);
199:       PetscStrcmp(buff,type[1],&ismult);
200:       PetscStrcmp(buff,type[2],&isfull);
201:       PetscStrcmp(buff,type[3],&iscasc);
202:       PetscStrcmp(buff,"kaskade",&iskask);

204:       if      (isadd)  mg = MGADDITIVE;
205:       else if (ismult) mg = MGMULTIPLICATIVE;
206:       else if (isfull) mg = MGFULL;
207:       else if (iskask) mg = MGKASKADE;
208:       else if (iscasc) mg = MGKASKADE;
209:       else SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,"Unknown type: %s",buff);
210:       MGSetType(pc,mg);
211:     }
212:   PetscOptionsTail();
213:   return(0);
214: }

216: static int PCView_MG(PC pc,PetscViewer viewer)
217: {
218:   MG         *mg = (MG*)pc->data;
219:   KSP        kspu,kspd;
220:   int        itu,itd,ierr,levels = mg[0]->levels,i;
221:   PetscReal  dtol,atol,rtol;
222:   char       *cstring;
223:   PetscTruth isascii;

226:   PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&isascii);
227:   if (isascii) {
228:     SLESGetKSP(mg[0]->smoothu,&kspu);
229:     SLESGetKSP(mg[0]->smoothd,&kspd);
230:     KSPGetTolerances(kspu,&dtol,&atol,&rtol,&itu);
231:     KSPGetTolerances(kspd,&dtol,&atol,&rtol,&itd);
232:     if (mg[0]->am == MGMULTIPLICATIVE) cstring = "multiplicative";
233:     else if (mg[0]->am == MGADDITIVE)  cstring = "additive";
234:     else if (mg[0]->am == MGFULL)      cstring = "full";
235:     else if (mg[0]->am == MGKASKADE)   cstring = "Kaskade";
236:     else cstring = "unknown";
237:     PetscViewerASCIIPrintf(viewer,"  MG: type is %s, levels=%d cycles=%d, pre-smooths=%d, post-smooths=%dn",
238:                       cstring,levels,mg[0]->cycles,mg[0]->default_smoothu,mg[0]->default_smoothd);
239:     for (i=0; i<levels; i++) {
240:       PetscViewerASCIIPrintf(viewer,"Down solver (pre-smoother) on level %d -------------------------------n",i);
241:       PetscViewerASCIIPushTab(viewer);
242:       SLESView(mg[i]->smoothd,viewer);
243:       PetscViewerASCIIPopTab(viewer);
244:       if (mg[i]->smoothd == mg[i]->smoothu) {
245:         PetscViewerASCIIPrintf(viewer,"Up solver same as down solvern");
246:       } else {
247:         PetscViewerASCIIPrintf(viewer,"Up solver (post-smoother) on level %d -------------------------------n",i);
248:         PetscViewerASCIIPushTab(viewer);
249:         SLESView(mg[i]->smoothu,viewer);
250:         PetscViewerASCIIPopTab(viewer);
251:       }
252:     }
253:   } else {
254:     SETERRQ1(1,"Viewer type %s not supported for PCMG",((PetscObject)viewer)->type_name);
255:   }
256:   return(0);
257: }

259: /*
260:     Calls setup for the SLES on each level
261: */
262: static int PCSetUp_MG(PC pc)
263: {
264:   MG          *mg = (MG*)pc->data;
265:   int         ierr,i,n = mg[0]->levels;
266:   KSP         ksp;
267:   PC          cpc;
268:   PetscTruth  preonly,lu,redundant,monitor = PETSC_FALSE,dump;
269:   PetscViewer ascii;
270:   MPI_Comm    comm;

273:   /*
274:      temporarily stick pc->vec into mg[0]->b and x so that 
275:    SLESSetUp is happy. Since currently those slots are empty.
276:   */
277:   mg[n-1]->x = pc->vec;
278:   mg[n-1]->b = pc->vec;

280:   if (pc->setupcalled == 0) {
281:     PetscOptionsHasName(0,"-pc_mg_monitor",&monitor);
282: 
283:     for (i=1; i<n; i++) {
284:       if (mg[i]->smoothd) {
285:         if (monitor) {
286:           SLESGetKSP(mg[i]->smoothd,&ksp);
287:           PetscObjectGetComm((PetscObject)ksp,&comm);
288:           PetscViewerASCIIOpen(comm,"stdout",&ascii);
289:           PetscViewerASCIISetTab(ascii,n-i);
290:           KSPSetMonitor(ksp,KSPDefaultMonitor,ascii,(int(*)(void*))PetscViewerDestroy);
291:         }
292:         SLESSetFromOptions(mg[i]->smoothd);
293:       }
294:     }
295:     for (i=0; i<n; i++) {
296:       if (mg[i]->smoothu && mg[i]->smoothu != mg[i]->smoothd) {
297:         if (monitor) {
298:           SLESGetKSP(mg[i]->smoothu,&ksp);
299:           PetscObjectGetComm((PetscObject)ksp,&comm);
300:           PetscViewerASCIIOpen(comm,"stdout",&ascii);
301:           PetscViewerASCIISetTab(ascii,n-i);
302:           KSPSetMonitor(ksp,KSPDefaultMonitor,ascii,(int(*)(void*))PetscViewerDestroy);
303:         }
304:         SLESSetFromOptions(mg[i]->smoothu);
305:       }
306:     }
307:   }

309:   for (i=1; i<n; i++) {
310:     if (mg[i]->smoothd) {
311:       SLESGetKSP(mg[i]->smoothd,&ksp);
312:       KSPSetInitialGuessNonzero(ksp);
313:       SLESSetUp(mg[i]->smoothd,mg[i]->b,mg[i]->x);
314:     }
315:   }
316:   for (i=0; i<n; i++) {
317:     if (mg[i]->smoothu && mg[i]->smoothu != mg[i]->smoothd) {
318:       SLESGetKSP(mg[i]->smoothu,&ksp);
319:       KSPSetInitialGuessNonzero(ksp);
320:       SLESSetUp(mg[i]->smoothu,mg[i]->b,mg[i]->x);
321:     }
322:   }

324:   /*
325:       If coarse solver is not direct method then DO NOT USE preonly 
326:   */
327:   SLESGetKSP(mg[0]->smoothd,&ksp);
328:   PetscTypeCompare((PetscObject)ksp,KSPPREONLY,&preonly);
329:   if (preonly) {
330:     SLESGetPC(mg[0]->smoothd,&cpc);
331:     PetscTypeCompare((PetscObject)cpc,PCLU,&lu);
332:     PetscTypeCompare((PetscObject)cpc,PCREDUNDANT,&redundant);
333:     if (!lu && !redundant) {
334:       KSPSetType(ksp,KSPGMRES);
335:     }
336:   }

338:   if (pc->setupcalled == 0) {
339:     if (monitor) {
340:       SLESGetKSP(mg[0]->smoothd,&ksp);
341:       PetscObjectGetComm((PetscObject)ksp,&comm);
342:       PetscViewerASCIIOpen(comm,"stdout",&ascii);
343:       PetscViewerASCIISetTab(ascii,n);
344:       KSPSetMonitor(ksp,KSPDefaultMonitor,ascii,(int(*)(void*))PetscViewerDestroy);
345:     }
346:     SLESSetFromOptions(mg[0]->smoothd);
347:   }

349:   SLESSetUp(mg[0]->smoothd,mg[0]->b,mg[0]->x);

351:   /*
352:      Dump the interpolation/restriction matrices to matlab plus the 
353:    Jacobian/stiffness on each level. This allows Matlab users to 
354:    easily check if the Galerkin condition A_c = R A_f R^T is satisfied */
355:   PetscOptionsHasName(pc->prefix,"-pc_mg_dump_matlab",&dump);
356:   if (dump) {
357:     for (i=1; i<n; i++) {
358:       MatView(mg[i]->restrct,PETSC_VIEWER_SOCKET_(pc->comm));
359:     }
360:     for (i=0; i<n; i++) {
361:       SLESGetPC(mg[i]->smoothd,&pc);
362:       MatView(pc->mat,PETSC_VIEWER_SOCKET_(pc->comm));
363:     }
364:   }

366:   return(0);
367: }

369: /* -------------------------------------------------------------------------------------*/

371: /*@C
372:    MGSetLevels - Sets the number of levels to use with MG.
373:    Must be called before any other MG routine.

375:    Collective on PC

377:    Input Parameters:
378: +  pc - the preconditioner context
379: .  levels - the number of levels
380: -  comms - optional communicators for each level; this is to allow solving the coarser problems
381:            on smaller sets of processors. Use PETSC_NULL_OBJECT for default in Fortran

383:    Level: intermediate

385:    Notes:
386:      If the number of levels is one then the multigrid uses the -mg_levels prefix
387:   for setting the level options rather than the -mg_coarse prefix.

389: .keywords: MG, set, levels, multigrid

391: .seealso: MGSetType(), MGGetLevels()
392: @*/
393: int MGSetLevels(PC pc,int levels,MPI_Comm *comms)
394: {
396:   MG  *mg;


401:   if (pc->data) {
402:     SETERRQ(1,"Number levels already set for MGn
403:     make sure that you call MGSetLevels() before SLESSetFromOptions()");
404:   }
405:   ierr                     = MGCreate_Private(pc->comm,levels,pc,comms,&mg);
406:   mg[0]->am                = MGMULTIPLICATIVE;
407:   pc->data                 = (void*)mg;
408:   pc->ops->applyrichardson = PCApplyRichardson_MG;
409:   return(0);
410: }

412: /*@
413:    MGGetLevels - Gets the number of levels to use with MG.

415:    Not Collective

417:    Input Parameter:
418: .  pc - the preconditioner context

420:    Output parameter:
421: .  levels - the number of levels

423:    Level: advanced

425: .keywords: MG, get, levels, multigrid

427: .seealso: MGSetLevels()
428: @*/
429: int MGGetLevels(PC pc,int *levels)
430: {
431:   MG  *mg;


436:   mg      = (MG*)pc->data;
437:   *levels = mg[0]->levels;
438:   return(0);
439: }

441: /*@
442:    MGSetType - Determines the form of multigrid to use:
443:    multiplicative, additive, full, or the Kaskade algorithm.

445:    Collective on PC

447:    Input Parameters:
448: +  pc - the preconditioner context
449: -  form - multigrid form, one of MGMULTIPLICATIVE, MGADDITIVE,
450:    MGFULL, MGKASKADE

452:    Options Database Key:
453: .  -pc_mg_type <form> - Sets <form>, one of multiplicative,
454:    additive, full, kaskade   

456:    Level: advanced

458: .keywords: MG, set, method, multiplicative, additive, full, Kaskade, multigrid

460: .seealso: MGSetLevels()
461: @*/
462: int MGSetType(PC pc,MGType form)
463: {
464:   MG *mg;

468:   mg = (MG*)pc->data;

470:   if (!mg) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling");
471:   mg[0]->am = form;
472:   if (form == MGMULTIPLICATIVE) pc->ops->applyrichardson = PCApplyRichardson_MG;
473:   else pc->ops->applyrichardson = 0;
474:   return(0);
475: }

477: /*@
478:    MGSetCycles - Sets the type cycles to use.  Use MGSetCyclesOnLevel() for more 
479:    complicated cycling.

481:    Collective on PC

483:    Input Parameters:
484: +  mg - the multigrid context 
485: -  n - the number of cycles

487:    Options Database Key:
488: $  -pc_mg_cycles n - 1 denotes a V-cycle; 2 denotes a W-cycle.

490:    Level: advanced

492: .keywords: MG, set, cycles, V-cycle, W-cycle, multigrid

494: .seealso: MGSetCyclesOnLevel()
495: @*/
496: int MGSetCycles(PC pc,int n)
497: {
498:   MG  *mg;
499:   int i,levels;

503:   mg     = (MG*)pc->data;
504:   if (!mg) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling");
505:   levels = mg[0]->levels;

507:   for (i=0; i<levels; i++) {
508:     mg[i]->cycles  = n;
509:   }
510:   return(0);
511: }

513: /*@
514:    MGCheck - Checks that all components of the MG structure have 
515:    been set.

517:    Collective on PC

519:    Input Parameters:
520: .  mg - the MG structure

522:    Level: advanced

524: .keywords: MG, check, set, multigrid
525: @*/
526: int MGCheck(PC pc)
527: {
528:   MG  *mg;
529:   int i,n,count = 0;

533:   mg = (MG*)pc->data;

535:   if (!mg) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling");

537:   n = mg[0]->levels;

539:   for (i=1; i<n; i++) {
540:     if (!mg[i]->restrct) {
541:       (*PetscErrorPrintf)("No restrict set level %d n",n-i); count++;
542:     }
543:     if (!mg[i]->interpolate) {
544:       (*PetscErrorPrintf)("No interpolate set level %d n",n-i); count++;
545:     }
546:     if (!mg[i]->residual) {
547:       (*PetscErrorPrintf)("No residual set level %d n",n-i); count++;
548:     }
549:     if (!mg[i]->smoothu) {
550:       (*PetscErrorPrintf)("No smoothup set level %d n",n-i); count++;
551:     }
552:     if (!mg[i]->smoothd) {
553:       (*PetscErrorPrintf)("No smoothdown set level %d n",n-i); count++;
554:     }
555:     if (!mg[i]->r) {
556:       (*PetscErrorPrintf)("No r set level %d n",n-i); count++;
557:     }
558:     if (!mg[i-1]->x) {
559:       (*PetscErrorPrintf)("No x set level %d n",n-i); count++;
560:     }
561:     if (!mg[i-1]->b) {
562:       (*PetscErrorPrintf)("No b set level %d n",n-i); count++;
563:     }
564:   }
565:   PetscFunctionReturn(count);
566: }


569: /*@
570:    MGSetNumberSmoothDown - Sets the number of pre-smoothing steps to
571:    use on all levels. Use MGGetSmootherDown() to set different 
572:    pre-smoothing steps on different levels.

574:    Collective on PC

576:    Input Parameters:
577: +  mg - the multigrid context 
578: -  n - the number of smoothing steps

580:    Options Database Key:
581: .  -pc_mg_smoothdown <n> - Sets number of pre-smoothing steps

583:    Level: advanced

585: .keywords: MG, smooth, down, pre-smoothing, steps, multigrid

587: .seealso: MGSetNumberSmoothUp()
588: @*/
589: int MGSetNumberSmoothDown(PC pc,int n)
590: {
591:   MG  *mg;
592:   int i,levels,ierr;
593:   KSP ksp;

597:   mg     = (MG*)pc->data;
598:   if (!mg) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling");
599:   levels = mg[0]->levels;

601:   for (i=0; i<levels; i++) {
602:     SLESGetKSP(mg[i]->smoothd,&ksp);
603:     KSPSetTolerances(ksp,PETSC_DEFAULT,PETSC_DEFAULT,PETSC_DEFAULT,n);
604:     mg[i]->default_smoothd = n;
605:   }
606:   return(0);
607: }

609: /*@
610:    MGSetNumberSmoothUp - Sets the number of post-smoothing steps to use 
611:    on all levels. Use MGGetSmootherUp() to set different numbers of 
612:    post-smoothing steps on different levels.

614:    Collective on PC

616:    Input Parameters:
617: +  mg - the multigrid context 
618: -  n - the number of smoothing steps

620:    Options Database Key:
621: .  -pc_mg_smoothup <n> - Sets number of post-smoothing steps

623:    Level: advanced

625: .keywords: MG, smooth, up, post-smoothing, steps, multigrid

627: .seealso: MGSetNumberSmoothDown()
628: @*/
629: int  MGSetNumberSmoothUp(PC pc,int n)
630: {
631:   MG  *mg;
632:   int i,levels,ierr;
633:   KSP ksp;

637:   mg     = (MG*)pc->data;
638:   if (!mg) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must set MG levels before calling");
639:   levels = mg[0]->levels;

641:   for (i=0; i<levels; i++) {
642:     SLESGetKSP(mg[i]->smoothu,&ksp);
643:     KSPSetTolerances(ksp,PETSC_DEFAULT,PETSC_DEFAULT,PETSC_DEFAULT,n);
644:     mg[i]->default_smoothu = n;
645:   }
646:   return(0);
647: }

649: /* ----------------------------------------------------------------------------------------*/

651: EXTERN_C_BEGIN
652: int PCCreate_MG(PC pc)
653: {
655:   pc->ops->apply          = PCApply_MG;
656:   pc->ops->setup          = PCSetUp_MG;
657:   pc->ops->destroy        = PCDestroy_MG;
658:   pc->ops->setfromoptions = PCSetFromOptions_MG;
659:   pc->ops->view           = PCView_MG;

661:   pc->data                = (void*)0;
662:   return(0);
663: }
664: EXTERN_C_END