Actual source code: damgsnes.c
1: /*$Id: damgsnes.c,v 1.51 2001/08/06 21:18:06 bsmith Exp $*/
2:
3: #include petscda.h
4: #include petscmg.h
7: /*
8: Evaluates the Jacobian on all of the grids. It is used by DMMG to provide the
9: ComputeJacobian() function that SNESSetJacobian() requires.
10: */
11: #undef __FUNCT__
13: int DMMGComputeJacobian_Multigrid(SNES snes,Vec X,Mat *J,Mat *B,MatStructure *flag,void *ptr)
14: {
15: DMMG *dmmg = (DMMG*)ptr;
16: int ierr,i,nlevels = dmmg[0]->nlevels;
17: SLES sles,lsles;
18: PC pc;
19: PetscTruth ismg;
20: Vec W;
23: if (!dmmg) SETERRQ(1,"Passing null as user context which should contain DMMG");
25: /* compute Jacobian on finest grid */
26: (*DMMGGetFine(dmmg)->computejacobian)(snes,X,J,B,flag,DMMGGetFine(dmmg));
27: MatSNESMFSetBase(DMMGGetFine(dmmg)->J,X);
29: /* create coarser grid Jacobians for preconditioner if multigrid is the preconditioner */
30: SNESGetSLES(snes,&sles);
31: SLESGetPC(sles,&pc);
32: PetscTypeCompare((PetscObject)pc,PCMG,&ismg);
33: if (ismg) {
35: MGGetSmoother(pc,nlevels-1,&lsles);
36: SLESSetOperators(lsles,DMMGGetFine(dmmg)->J,DMMGGetFine(dmmg)->B,*flag);
38: for (i=nlevels-1; i>0; i--) {
40: if (!dmmg[i-1]->w) {
41: VecDuplicate(dmmg[i-1]->x,&dmmg[i-1]->w);
42: }
44: W = dmmg[i-1]->w;
45: /* restrict X to coarser grid */
46: MatRestrict(dmmg[i]->R,X,W);
47: X = W;
49: /* scale to "natural" scaling for that grid */
50: VecPointwiseMult(dmmg[i]->Rscale,X,X);
52: /* tell the base vector for matrix free multiplies */
53: MatSNESMFSetBase(dmmg[i-1]->J,X);
55: /* compute Jacobian on coarse grid */
56: (*dmmg[i-1]->computejacobian)(snes,X,&dmmg[i-1]->J,&dmmg[i-1]->B,flag,dmmg[i-1]);
58: MGGetSmoother(pc,i-1,&lsles);
59: SLESSetOperators(lsles,dmmg[i-1]->J,dmmg[i-1]->B,*flag);
60: }
61: }
62: return(0);
63: }
65: /* ---------------------------------------------------------------------------*/
67: #undef __FUNCT__
69: /*
70: DMMGFormFunction - This is a universal global FormFunction used by the DMMG code
71: when the user provides a local function.
73: Input Parameters:
74: + snes - the SNES context
75: . X - input vector
76: - ptr - optional user-defined context, as set by SNESSetFunction()
78: Output Parameter:
79: . F - function vector
81: */
82: int DMMGFormFunction(SNES snes,Vec X,Vec F,void *ptr)
83: {
84: DMMG dmmg = (DMMG)ptr;
85: int ierr;
86: Vec localX;
87: DA da = (DA)dmmg->dm;
90: DAGetLocalVector(da,&localX);
91: /*
92: Scatter ghost points to local vector, using the 2-step process
93: DAGlobalToLocalBegin(), DAGlobalToLocalEnd().
94: */
95: DAGlobalToLocalBegin(da,X,INSERT_VALUES,localX);
96: DAGlobalToLocalEnd(da,X,INSERT_VALUES,localX);
97: DAFormFunction1(da,localX,F,dmmg->user);
98: DARestoreLocalVector(da,&localX);
99: return(0);
100: }
102: #undef __FUNCT__
104: /*@C
105: SNESDAFormFunction - This is a universal function evaluation routine that
106: may be used with SNESSetFunction() as long as the user context has a DA
107: as its first record and the user has called DASetLocalFunction().
109: Collective on SNES
111: Input Parameters:
112: + snes - the SNES context
113: . X - input vector
114: . F - function vector
115: - ptr - pointer to a structure that must have a DA as its first entry. For example this
116: could be a DMMG
118: Level: intermediate
120: .seealso: DASetLocalFunction(), DASetLocalJacobian(), DASetLocalAdicFunction(), DASetLocalAdicMFFunction(),
121: SNESSetFunction(), SNESSetJacobian()
123: @*/
124: int SNESDAFormFunction(SNES snes,Vec X,Vec F,void *ptr)
125: {
126: int ierr;
127: Vec localX;
128: DA da = *(DA*)ptr;
131: DAGetLocalVector(da,&localX);
132: /*
133: Scatter ghost points to local vector, using the 2-step process
134: DAGlobalToLocalBegin(), DAGlobalToLocalEnd().
135: */
136: DAGlobalToLocalBegin(da,X,INSERT_VALUES,localX);
137: DAGlobalToLocalEnd(da,X,INSERT_VALUES,localX);
138: DAFormFunction1(da,localX,F,ptr);
139: DARestoreLocalVector(da,&localX);
140: return(0);
141: }
143: /* ---------------------------------------------------------------------------------------------------------------------------*/
145: #undef __FUNCT__
147: int DMMGComputeJacobianWithFD(SNES snes,Vec x1,Mat *J,Mat *B,MatStructure *flag,void *ctx)
148: {
149: int ierr;
150: DMMG dmmg = (DMMG)ctx;
151:
153: SNESDefaultComputeJacobianColor(snes,x1,J,B,flag,dmmg->fdcoloring);
154: return(0);
155: }
157: #undef __FUNCT__
159: int DMMGComputeJacobianWithMF(SNES snes,Vec x1,Mat *J,Mat *B,MatStructure *flag,void *ctx)
160: {
161: int ierr;
162:
164: MatAssemblyBegin(*J,MAT_FINAL_ASSEMBLY);
165: MatAssemblyEnd(*J,MAT_FINAL_ASSEMBLY);
166: return(0);
167: }
169: #undef __FUNCT__
171: /*
172: DMMGComputeJacobianWithAdic - Evaluates the Jacobian via Adic when the user has provided
173: a local function evaluation routine.
174: */
175: int DMMGComputeJacobianWithAdic(SNES snes,Vec X,Mat *J,Mat *B,MatStructure *flag,void *ptr)
176: {
177: DMMG dmmg = (DMMG) ptr;
178: int ierr;
179: Vec localX;
180: DA da = (DA) dmmg->dm;
183: DAGetLocalVector(da,&localX);
184: DAGlobalToLocalBegin(da,X,INSERT_VALUES,localX);
185: DAGlobalToLocalEnd(da,X,INSERT_VALUES,localX);
186: DAComputeJacobian1WithAdic(da,localX,*B,dmmg->user);
187: DARestoreLocalVector(da,&localX);
188: /* Assemble true Jacobian; if it is different */
189: if (*J != *B) {
190: ierr = MatAssemblyBegin(*J,MAT_FINAL_ASSEMBLY);
191: ierr = MatAssemblyEnd(*J,MAT_FINAL_ASSEMBLY);
192: }
193: ierr = MatSetOption(*B,MAT_NEW_NONZERO_LOCATION_ERR);
194: *flag = SAME_NONZERO_PATTERN;
195: return(0);
196: }
198: #undef __FUNCT__
200: /*@
201: SNESDAComputeJacobianWithAdic - This is a universal Jacobian evaluation routine
202: that may be used with SNESSetJacobian() as long as the user context has a DA as
203: its first record and DASetLocalAdicFunction() has been called.
205: Collective on SNES
207: Input Parameters:
208: + snes - the SNES context
209: . X - input vector
210: . J - Jacobian
211: . B - Jacobian used in preconditioner (usally same as J)
212: . flag - indicates if the matrix changed its structure
213: - ptr - optional user-defined context, as set by SNESSetFunction()
215: Level: intermediate
217: .seealso: DASetLocalFunction(), DASetLocalAdicFunction(), SNESSetFunction(), SNESSetJacobian()
219: @*/
220: int SNESDAComputeJacobianWithAdic(SNES snes,Vec X,Mat *J,Mat *B,MatStructure *flag,void *ptr)
221: {
222: DA da = *(DA*) ptr;
223: int ierr;
224: Vec localX;
227: DAGetLocalVector(da,&localX);
228: DAGlobalToLocalBegin(da,X,INSERT_VALUES,localX);
229: DAGlobalToLocalEnd(da,X,INSERT_VALUES,localX);
230: DAComputeJacobian1WithAdic(da,localX,*B,ptr);
231: DARestoreLocalVector(da,&localX);
232: /* Assemble true Jacobian; if it is different */
233: if (*J != *B) {
234: ierr = MatAssemblyBegin(*J,MAT_FINAL_ASSEMBLY);
235: ierr = MatAssemblyEnd(*J,MAT_FINAL_ASSEMBLY);
236: }
237: ierr = MatSetOption(*B,MAT_NEW_NONZERO_LOCATION_ERR);
238: *flag = SAME_NONZERO_PATTERN;
239: return(0);
240: }
242: #undef __FUNCT__
244: /*
245: SNESDAComputeJacobianWithAdifor - This is a universal Jacobian evaluation routine
246: that may be used with SNESSetJacobian() from Fortran as long as the user context has
247: a DA as its first record and DASetLocalAdiforFunction() has been called.
249: Collective on SNES
251: Input Parameters:
252: + snes - the SNES context
253: . X - input vector
254: . J - Jacobian
255: . B - Jacobian used in preconditioner (usally same as J)
256: . flag - indicates if the matrix changed its structure
257: - ptr - optional user-defined context, as set by SNESSetFunction()
259: Level: intermediate
261: .seealso: DASetLocalFunction(), DASetLocalAdicFunction(), SNESSetFunction(), SNESSetJacobian()
263: */
264: int SNESDAComputeJacobianWithAdifor(SNES snes,Vec X,Mat *J,Mat *B,MatStructure *flag,void *ptr)
265: {
266: DA da = *(DA*) ptr;
267: int ierr;
268: Vec localX;
271: DAGetLocalVector(da,&localX);
272: DAGlobalToLocalBegin(da,X,INSERT_VALUES,localX);
273: DAGlobalToLocalEnd(da,X,INSERT_VALUES,localX);
274: DAComputeJacobian1WithAdifor(da,localX,*B,ptr);
275: DARestoreLocalVector(da,&localX);
276: /* Assemble true Jacobian; if it is different */
277: if (*J != *B) {
278: ierr = MatAssemblyBegin(*J,MAT_FINAL_ASSEMBLY);
279: ierr = MatAssemblyEnd(*J,MAT_FINAL_ASSEMBLY);
280: }
281: ierr = MatSetOption(*B,MAT_NEW_NONZERO_LOCATION_ERR);
282: *flag = SAME_NONZERO_PATTERN;
283: return(0);
284: }
286: #undef __FUNCT__
288: /*
289: SNESDAComputeJacobian - This is a universal Jacobian evaluation routine for a
290: locally provided Jacobian.
292: Collective on SNES
294: Input Parameters:
295: + snes - the SNES context
296: . X - input vector
297: . J - Jacobian
298: . B - Jacobian used in preconditioner (usally same as J)
299: . flag - indicates if the matrix changed its structure
300: - ptr - optional user-defined context, as set by SNESSetFunction()
302: Level: intermediate
304: .seealso: DASetLocalFunction(), DASetLocalJacobian(), SNESSetFunction(), SNESSetJacobian()
306: */
307: int SNESDAComputeJacobian(SNES snes,Vec X,Mat *J,Mat *B,MatStructure *flag,void *ptr)
308: {
309: DA da = *(DA*) ptr;
310: int ierr;
311: Vec localX;
314: DAGetLocalVector(da,&localX);
315: DAGlobalToLocalBegin(da,X,INSERT_VALUES,localX);
316: DAGlobalToLocalEnd(da,X,INSERT_VALUES,localX);
317: DAComputeJacobian1(da,localX,*B,ptr);
318: DARestoreLocalVector(da,&localX);
319: /* Assemble true Jacobian; if it is different */
320: if (*J != *B) {
321: ierr = MatAssemblyBegin(*J,MAT_FINAL_ASSEMBLY);
322: ierr = MatAssemblyEnd(*J,MAT_FINAL_ASSEMBLY);
323: }
324: ierr = MatSetOption(*B,MAT_NEW_NONZERO_LOCATION_ERR);
325: *flag = SAME_NONZERO_PATTERN;
326: return(0);
327: }
329: #undef __FUNCT__
331: int DMMGSolveSNES(DMMG *dmmg,int level)
332: {
333: int ierr,nlevels = dmmg[0]->nlevels,its;
336: dmmg[0]->nlevels = level+1;
337: ierr = SNESSolve(dmmg[level]->snes,dmmg[level]->x,&its);
338: dmmg[0]->nlevels = nlevels;
339: return(0);
340: }
342: /* ===========================================================================================================*/
344: #undef __FUNCT__
346: /*@C
347: DMMGSetSNES - Sets the nonlinear function that defines the nonlinear set of equations
348: to be solved using the grid hierarchy.
350: Collective on DMMG
352: Input Parameter:
353: + dmmg - the context
354: . function - the function that defines the nonlinear system
355: - jacobian - optional function to compute Jacobian
357: Options Database Keys:
358: + -dmmg_snes_monitor
359: . -dmmg_jacobian_fd
360: . -dmmg_jacobian_ad
361: . -dmmg_jacobian_mf_fd_operator
362: . -dmmg_jacobian_mf_fd
363: . -dmmg_jacobian_mf_ad_operator
364: - -dmmg_jacobian_mf_ad
366: Level: advanced
368: .seealso DMMGCreate(), DMMGDestroy, DMMGSetSLES(), DMMGSetSNESLocal()
370: @*/
371: int DMMGSetSNES(DMMG *dmmg,int (*function)(SNES,Vec,Vec,void*),int (*jacobian)(SNES,Vec,Mat*,Mat*,MatStructure*,void*))
372: {
373: int ierr,i,nlevels = dmmg[0]->nlevels;
374: PetscTruth snesmonitor,mffdoperator,mffd,fdjacobian;
375: #if defined(PETSC_HAVE_ADIC) && !defined(PETSC_USE_COMPLEX) && !defined(PETSC_USE_SINGLE)
376: PetscTruth mfadoperator,mfad,adjacobian;
377: #endif
378: SLES sles;
379: PetscViewer ascii;
380: MPI_Comm comm;
383: if (!dmmg) SETERRQ(1,"Passing null as DMMG");
384: if (!jacobian) jacobian = DMMGComputeJacobianWithFD;
386: PetscOptionsBegin(dmmg[0]->comm,PETSC_NULL,"DMMG Options","SNES");
387: PetscOptionsName("-dmmg_snes_monitor","Monitor nonlinear convergence","SNESSetMonitor",&snesmonitor);
390: PetscOptionsName("-dmmg_jacobian_fd","Compute sparse Jacobian explicitly with finite differencing","DMMGSetSNES",&fdjacobian);
391: if (fdjacobian) jacobian = DMMGComputeJacobianWithFD;
392: #if defined(PETSC_HAVE_ADIC) && !defined(PETSC_USE_COMPLEX) && !defined(PETSC_USE_SINGLE)
393: PetscOptionsName("-dmmg_jacobian_ad","Compute sparse Jacobian explicitly with ADIC (automatic differentiation)","DMMGSetSNES",&adjacobian);
394: if (adjacobian) jacobian = DMMGComputeJacobianWithAdic;
395: #endif
397: PetscOptionsLogicalGroupBegin("-dmmg_jacobian_mf_fd_operator","Apply Jacobian via matrix free finite differencing","DMMGSetSNES",&mffdoperator);
398: PetscOptionsLogicalGroupEnd("-dmmg_jacobian_mf_fd","Apply Jacobian via matrix free finite differencing even in computing preconditioner","DMMGSetSNES",&mffd);
399: if (mffd) mffdoperator = PETSC_TRUE;
400: #if defined(PETSC_HAVE_ADIC) && !defined(PETSC_USE_COMPLEX) && !defined(PETSC_USE_SINGLE)
401: PetscOptionsLogicalGroupBegin("-dmmg_jacobian_mf_ad_operator","Apply Jacobian via matrix free ADIC (automatic differentiation)","DMMGSetSNES",&mfadoperator);
402: PetscOptionsLogicalGroupEnd("-dmmg_jacobian_mf_ad","Apply Jacobian via matrix free ADIC (automatic differentiation) even in computing preconditioner","DMMGSetSNES",&mfad);
403: if (mfad) mfadoperator = PETSC_TRUE;
404: #endif
405: PetscOptionsEnd();
407: /* create solvers for each level */
408: for (i=0; i<nlevels; i++) {
409: SNESCreate(dmmg[i]->comm,&dmmg[i]->snes);
410: if (snesmonitor) {
411: PetscObjectGetComm((PetscObject)dmmg[i]->snes,&comm);
412: PetscViewerASCIIOpen(comm,"stdout",&ascii);
413: PetscViewerASCIISetTab(ascii,nlevels-i);
414: SNESSetMonitor(dmmg[i]->snes,SNESDefaultMonitor,ascii,(int(*)(void*))PetscViewerDestroy);
415: }
417: if (mffdoperator) {
418: MatCreateSNESMF(dmmg[i]->snes,dmmg[i]->x,&dmmg[i]->J);
419: VecDuplicate(dmmg[i]->x,&dmmg[i]->work1);
420: VecDuplicate(dmmg[i]->x,&dmmg[i]->work2);
421: MatSNESMFSetFunction(dmmg[i]->J,dmmg[i]->work1,function,dmmg[i]);
422: if (mffd) {
423: dmmg[i]->B = dmmg[i]->J;
424: jacobian = DMMGComputeJacobianWithMF;
425: }
426: #if defined(PETSC_HAVE_ADIC) && !defined(PETSC_USE_COMPLEX) && !defined(PETSC_USE_SINGLE)
427: } else if (mfadoperator) {
428: MatRegisterDAAD();
429: MatCreateDAAD((DA)dmmg[i]->dm,&dmmg[i]->J);
430: MatDAADSetCtx(dmmg[i]->J,dmmg[i]->user);
431: if (mfad) {
432: dmmg[i]->B = dmmg[i]->J;
433: jacobian = DMMGComputeJacobianWithMF;
434: }
435: #endif
436: }
437:
438: if (!dmmg[i]->B) {
439: DMGetMatrix(dmmg[i]->dm,MATMPIAIJ,&dmmg[i]->B);
440: }
441: if (!dmmg[i]->J) {
442: dmmg[i]->J = dmmg[i]->B;
443: }
445: SNESGetSLES(dmmg[i]->snes,&sles);
446: DMMGSetUpLevel(dmmg,sles,i+1);
447:
448: /*
449: if the number of levels is > 1 then we want the coarse solve in the grid sequencing to use LU
450: when possible
451: */
452: if (nlevels > 1 && i == 0) {
453: PC pc;
454: SLES csles;
455: PetscTruth flg1,flg2,flg3;
457: SLESGetPC(sles,&pc);
458: MGGetCoarseSolve(pc,&csles);
459: SLESGetPC(csles,&pc);
460: PetscTypeCompare((PetscObject)pc,PCILU,&flg1);
461: PetscTypeCompare((PetscObject)pc,PCSOR,&flg2);
462: PetscTypeCompare((PetscObject)pc,PETSC_NULL,&flg3);
463: if (flg1 || flg2 || flg3) {
464: PCSetType(pc,PCLU);
465: }
466: }
468: SNESSetFromOptions(dmmg[i]->snes);
469: dmmg[i]->solve = DMMGSolveSNES;
470: dmmg[i]->computejacobian = jacobian;
471: dmmg[i]->computefunction = function;
472: }
475: if (jacobian == DMMGComputeJacobianWithFD) {
476: ISColoring iscoloring;
477: for (i=0; i<nlevels; i++) {
478: DMGetColoring(dmmg[i]->dm,IS_COLORING_LOCAL,&iscoloring);
479: MatFDColoringCreate(dmmg[i]->B,iscoloring,&dmmg[i]->fdcoloring);
480: ISColoringDestroy(iscoloring);
481: MatFDColoringSetFunction(dmmg[i]->fdcoloring,(int(*)(void))function,dmmg[i]);
482: MatFDColoringSetFromOptions(dmmg[i]->fdcoloring);
483: }
484: #if defined(PETSC_HAVE_ADIC) && !defined(PETSC_USE_COMPLEX) && !defined(PETSC_USE_SINGLE)
485: } else if (jacobian == DMMGComputeJacobianWithAdic) {
486: for (i=0; i<nlevels; i++) {
487: ISColoring iscoloring;
488: DMGetColoring(dmmg[i]->dm,IS_COLORING_GHOSTED,&iscoloring);
489: MatSetColoring(dmmg[i]->B,iscoloring);
490: ISColoringDestroy(iscoloring);
491: }
492: #endif
493: }
495: for (i=0; i<nlevels; i++) {
496: SNESSetJacobian(dmmg[i]->snes,dmmg[i]->J,dmmg[i]->B,DMMGComputeJacobian_Multigrid,dmmg);
497: SNESSetFunction(dmmg[i]->snes,dmmg[i]->b,function,dmmg[i]);
498: }
500: /* Create interpolation scaling */
501: for (i=1; i<nlevels; i++) {
502: DMGetInterpolationScale(dmmg[i-1]->dm,dmmg[i]->dm,dmmg[i]->R,&dmmg[i]->Rscale);
503: }
505: return(0);
506: }
508: #undef __FUNCT__
510: /*@C
511: DMMGSetInitialGuess - Sets the function that computes an initial guess, if not given
512: uses 0.
514: Collective on DMMG and SNES
516: Input Parameter:
517: + dmmg - the context
518: - guess - the function
520: Level: advanced
522: .seealso DMMGCreate(), DMMGDestroy, DMMGSetSLES()
524: @*/
525: int DMMGSetInitialGuess(DMMG *dmmg,int (*guess)(SNES,Vec,void*))
526: {
527: int i,nlevels = dmmg[0]->nlevels;
530: for (i=0; i<nlevels; i++) {
531: dmmg[i]->initialguess = guess;
532: }
533: return(0);
534: }
537: /*M
538: DMMGSetSNESLocal - Sets the local user function that defines the nonlinear set of equations
539: that will use the grid hierarchy and (optionally) its derivative.
541: Collective on DMMG
543: Synopsis:
544: int DMMGSetSNESLocal(DMMG *dmmg,DALocalFunction1 function, DALocalFunction1 jacobian,
545: DALocalFunction1 ad_function, DALocalFunction1 admf_function);
547: Input Parameter:
548: + dmmg - the context
549: . function - the function that defines the nonlinear system
550: . jacobian - function defines the local part of the Jacobian (not currently supported)
551: . ad_function - the name of the function with an ad_ prefix. This is ignored if ADIC is
552: not installed
553: - admf_function - the name of the function with an ad_ prefix. This is ignored if ADIC is
554: not installed
556: Options Database Keys:
557: + -dmmg_snes_monitor
558: . -dmmg_jacobian_fd
559: . -dmmg_jacobian_ad
560: . -dmmg_jacobian_mf_fd_operator
561: . -dmmg_jacobian_mf_fd
562: . -dmmg_jacobian_mf_ad_operator
563: - -dmmg_jacobian_mf_ad
566: Level: intermediate
568: Notes:
569: If ADIC or ADIFOR have been installed, this routine can use ADIC or ADIFOR to compute
570: the derivative; however, that function cannot call other functions except those in
571: standard C math libraries.
573: If ADIC/ADIFOR have not been installed and the Jacobian is not provided, this routine
574: uses finite differencing to approximate the Jacobian.
576: .seealso DMMGCreate(), DMMGDestroy, DMMGSetSLES(), DMMGSetSNES()
578: M*/
580: #undef __FUNCT__
582: int DMMGSetSNESLocal_Private(DMMG *dmmg,DALocalFunction1 function,DALocalFunction1 jacobian,DALocalFunction1 ad_function,DALocalFunction1 admf_function)
583: {
584: int ierr,i,nlevels = dmmg[0]->nlevels;
585: int (*computejacobian)(SNES,Vec,Mat*,Mat*,MatStructure*,void*) = 0;
589: if (jacobian) computejacobian = SNESDAComputeJacobian;
590: #if defined(PETSC_HAVE_ADIC) && !defined(PETSC_USE_COMPLEX) && !defined(PETSC_USE_SINGLE)
591: else if (ad_function) computejacobian = DMMGComputeJacobianWithAdic;
592: #endif
594: DMMGSetSNES(dmmg,DMMGFormFunction,computejacobian);
595: for (i=0; i<nlevels; i++) {
596: DASetLocalFunction((DA)dmmg[i]->dm,function);
597: DASetLocalJacobian((DA)dmmg[i]->dm,jacobian);
598: DASetLocalAdicFunction((DA)dmmg[i]->dm,ad_function);
599: DASetLocalAdicMFFunction((DA)dmmg[i]->dm,admf_function);
600: }
601: return(0);
602: }
604: #undef __FUNCT__
606: static int DMMGFunctioni(int i,Vec u,PetscScalar* r,void* ctx)
607: {
608: DMMG dmmg = (DMMG)ctx;
609: Vec U = dmmg->lwork1;
610: int ierr;
611: VecScatter gtol;
614: /* copy u into interior part of U */
615: DAGetScatter((DA)dmmg->dm,0,>ol,0);
616: VecScatterBegin(u,U,INSERT_VALUES,SCATTER_FORWARD_LOCAL,gtol);
617: VecScatterEnd(u,U,INSERT_VALUES,SCATTER_FORWARD_LOCAL,gtol);
619: DAFormFunctioni1((DA)dmmg->dm,i,U,r,dmmg->user);
620: return(0);
621: }
623: #undef __FUNCT__
625: static int DMMGFunctioniBase(Vec u,void* ctx)
626: {
627: DMMG dmmg = (DMMG)ctx;
628: Vec U = dmmg->lwork1;
629: int ierr;
632: DAGlobalToLocalBegin((DA)dmmg->dm,u,INSERT_VALUES,U);
633: DAGlobalToLocalEnd((DA)dmmg->dm,u,INSERT_VALUES,U);
634: return(0);
635: }
637: #undef __FUNCT__
639: int DMMGSetSNESLocali_Private(DMMG *dmmg,int (*functioni)(DALocalInfo*,MatStencil*,void*,PetscScalar*,void*),int (*adi)(DALocalInfo*,MatStencil*,void*,void*,void*),int (*adimf)(DALocalInfo*,MatStencil*,void*,void*,void*))
640: {
641: int ierr,i,nlevels = dmmg[0]->nlevels;
644: CHKMEMQ;
645: for (i=0; i<nlevels; i++) {
646: CHKMEMQ;
647: DASetLocalFunctioni((DA)dmmg[i]->dm,functioni);
648: CHKMEMQ;
649: DASetLocalAdicFunctioni((DA)dmmg[i]->dm,adi);
650: CHKMEMQ;
651: DASetLocalAdicMFFunctioni((DA)dmmg[i]->dm,adimf);
652: CHKMEMQ;
653: MatSNESMFSetFunctioni(dmmg[i]->J,DMMGFunctioni);
654: CHKMEMQ;
655: MatSNESMFSetFunctioniBase(dmmg[i]->J,DMMGFunctioniBase);
656: CHKMEMQ;
657: DACreateLocalVector((DA)dmmg[i]->dm,&dmmg[i]->lwork1);
658: }
659: return(0);
660: }
663: #if defined(PETSC_HAVE_ADIC) && !defined(PETSC_USE_COMPLEX) && !defined(PETSC_USE_SINGLE)
664: EXTERN_C_BEGIN
665: #include "adic/ad_utils.h"
666: EXTERN_C_END
668: #undef __FUNCT__
670: int PetscADView(int N,int nc,double *ptr,PetscViewer viewer)
671: {
672: int i,j,nlen = PetscADGetDerivTypeSize();
673: char *cptr = (char*)ptr;
674: double *values;
677: for (i=0; i<N; i++) {
678: printf("Element %d value %g derivatives: ",i,*(double*)cptr);
679: values = PetscADGetGradArray(cptr);
680: for (j=0; j<nc; j++) {
681: printf("%g ",*values++);
682: }
683: printf("n");
684: cptr += nlen;
685: }
687: return(0);
688: }
690: #endif