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: int DMMGComputeJacobian_Multigrid(SNES snes,Vec X,Mat *J,Mat *B,MatStructure *flag,void *ptr)
12: {
13: DMMG *dmmg = (DMMG*)ptr;
14: int ierr,i,nlevels = dmmg[0]->nlevels;
15: SLES sles,lsles;
16: PC pc;
17: PetscTruth ismg;
18: Vec W;
21: if (!dmmg) SETERRQ(1,"Passing null as user context which should contain DMMG");
23: /* compute Jacobian on finest grid */
24: (*DMMGGetFine(dmmg)->computejacobian)(snes,X,J,B,flag,DMMGGetFine(dmmg));
25: MatSNESMFSetBase(DMMGGetFine(dmmg)->J,X);
27: /* create coarser grid Jacobians for preconditioner if multigrid is the preconditioner */
28: SNESGetSLES(snes,&sles);
29: SLESGetPC(sles,&pc);
30: PetscTypeCompare((PetscObject)pc,PCMG,&ismg);
31: if (ismg) {
33: MGGetSmoother(pc,nlevels-1,&lsles);
34: SLESSetOperators(lsles,DMMGGetFine(dmmg)->J,DMMGGetFine(dmmg)->B,*flag);
36: for (i=nlevels-1; i>0; i--) {
38: if (!dmmg[i-1]->w) {
39: VecDuplicate(dmmg[i-1]->x,&dmmg[i-1]->w);
40: }
42: W = dmmg[i-1]->w;
43: /* restrict X to coarser grid */
44: MatRestrict(dmmg[i]->R,X,W);
45: X = W;
47: /* scale to "natural" scaling for that grid */
48: VecPointwiseMult(dmmg[i]->Rscale,X,X);
50: /* tell the base vector for matrix free multiplies */
51: MatSNESMFSetBase(dmmg[i-1]->J,X);
53: /* compute Jacobian on coarse grid */
54: (*dmmg[i-1]->computejacobian)(snes,X,&dmmg[i-1]->J,&dmmg[i-1]->B,flag,dmmg[i-1]);
56: MGGetSmoother(pc,i-1,&lsles);
57: SLESSetOperators(lsles,dmmg[i-1]->J,dmmg[i-1]->B,*flag);
58: }
59: }
60: return(0);
61: }
63: /* ---------------------------------------------------------------------------*/
65: /*
66: DMMGFormFunction - This is a universal global FormFunction used by the DMMG code
67: when the user provides a local function.
69: Input Parameters:
70: + snes - the SNES context
71: . X - input vector
72: - ptr - optional user-defined context, as set by SNESSetFunction()
74: Output Parameter:
75: . F - function vector
77: */
78: int DMMGFormFunction(SNES snes,Vec X,Vec F,void *ptr)
79: {
80: DMMG dmmg = (DMMG)ptr;
81: int ierr;
82: Vec localX;
83: DA da = (DA)dmmg->dm;
86: DAGetLocalVector(da,&localX);
87: /*
88: Scatter ghost points to local vector, using the 2-step process
89: DAGlobalToLocalBegin(), DAGlobalToLocalEnd().
90: */
91: DAGlobalToLocalBegin(da,X,INSERT_VALUES,localX);
92: DAGlobalToLocalEnd(da,X,INSERT_VALUES,localX);
93: DAFormFunction1(da,localX,F,dmmg->user);
94: DARestoreLocalVector(da,&localX);
95: return(0);
96: }
98: /*@C
99: SNESDAFormFunction - This is a universal function evaluation routine that
100: may be used with SNESSetFunction() as long as the user context has a DA
101: as its first record and the user has called DASetLocalFunction().
103: Collective on SNES
105: Input Parameters:
106: + snes - the SNES context
107: . X - input vector
108: . F - function vector
109: - ptr - pointer to a structure that must have a DA as its first entry. For example this
110: could be a DMMG
112: Level: intermediate
114: .seealso: DASetLocalFunction(), DASetLocalJacobian(), DASetLocalAdicFunction(), DASetLocalAdicMFFunction(),
115: SNESSetFunction(), SNESSetJacobian()
117: @*/
118: int SNESDAFormFunction(SNES snes,Vec X,Vec F,void *ptr)
119: {
120: int ierr;
121: Vec localX;
122: DA da = *(DA*)ptr;
125: DAGetLocalVector(da,&localX);
126: /*
127: Scatter ghost points to local vector, using the 2-step process
128: DAGlobalToLocalBegin(), DAGlobalToLocalEnd().
129: */
130: DAGlobalToLocalBegin(da,X,INSERT_VALUES,localX);
131: DAGlobalToLocalEnd(da,X,INSERT_VALUES,localX);
132: DAFormFunction1(da,localX,F,ptr);
133: DARestoreLocalVector(da,&localX);
134: return(0);
135: }
137: /* ---------------------------------------------------------------------------------------------------------------------------*/
139: int DMMGComputeJacobianWithFD(SNES snes,Vec x1,Mat *J,Mat *B,MatStructure *flag,void *ctx)
140: {
141: int ierr;
142: DMMG dmmg = (DMMG)ctx;
143:
145: SNESDefaultComputeJacobianColor(snes,x1,J,B,flag,dmmg->fdcoloring);
146: return(0);
147: }
149: int DMMGComputeJacobianWithMF(SNES snes,Vec x1,Mat *J,Mat *B,MatStructure *flag,void *ctx)
150: {
151: int ierr;
152:
154: MatAssemblyBegin(*J,MAT_FINAL_ASSEMBLY);
155: MatAssemblyEnd(*J,MAT_FINAL_ASSEMBLY);
156: return(0);
157: }
159: /*
160: DMMGComputeJacobianWithAdic - Evaluates the Jacobian via Adic when the user has provided
161: a local function evaluation routine.
162: */
163: int DMMGComputeJacobianWithAdic(SNES snes,Vec X,Mat *J,Mat *B,MatStructure *flag,void *ptr)
164: {
165: DMMG dmmg = (DMMG) ptr;
166: int ierr;
167: Vec localX;
168: DA da = (DA) dmmg->dm;
171: DAGetLocalVector(da,&localX);
172: DAGlobalToLocalBegin(da,X,INSERT_VALUES,localX);
173: DAGlobalToLocalEnd(da,X,INSERT_VALUES,localX);
174: DAComputeJacobian1WithAdic(da,localX,*B,dmmg->user);
175: DARestoreLocalVector(da,&localX);
176: /* Assemble true Jacobian; if it is different */
177: if (*J != *B) {
178: ierr = MatAssemblyBegin(*J,MAT_FINAL_ASSEMBLY);
179: ierr = MatAssemblyEnd(*J,MAT_FINAL_ASSEMBLY);
180: }
181: ierr = MatSetOption(*B,MAT_NEW_NONZERO_LOCATION_ERR);
182: *flag = SAME_NONZERO_PATTERN;
183: return(0);
184: }
186: /*@
187: SNESDAComputeJacobianWithAdic - This is a universal Jacobian evaluation routine
188: that may be used with SNESSetJacobian() as long as the user context has a DA as
189: its first record and DASetLocalAdicFunction() has been called.
191: Collective on SNES
193: Input Parameters:
194: + snes - the SNES context
195: . X - input vector
196: . J - Jacobian
197: . B - Jacobian used in preconditioner (usally same as J)
198: . flag - indicates if the matrix changed its structure
199: - ptr - optional user-defined context, as set by SNESSetFunction()
201: Level: intermediate
203: .seealso: DASetLocalFunction(), DASetLocalAdicFunction(), SNESSetFunction(), SNESSetJacobian()
205: @*/
206: int SNESDAComputeJacobianWithAdic(SNES snes,Vec X,Mat *J,Mat *B,MatStructure *flag,void *ptr)
207: {
208: DA da = *(DA*) ptr;
209: int ierr;
210: Vec localX;
213: DAGetLocalVector(da,&localX);
214: DAGlobalToLocalBegin(da,X,INSERT_VALUES,localX);
215: DAGlobalToLocalEnd(da,X,INSERT_VALUES,localX);
216: DAComputeJacobian1WithAdic(da,localX,*B,ptr);
217: DARestoreLocalVector(da,&localX);
218: /* Assemble true Jacobian; if it is different */
219: if (*J != *B) {
220: ierr = MatAssemblyBegin(*J,MAT_FINAL_ASSEMBLY);
221: ierr = MatAssemblyEnd(*J,MAT_FINAL_ASSEMBLY);
222: }
223: ierr = MatSetOption(*B,MAT_NEW_NONZERO_LOCATION_ERR);
224: *flag = SAME_NONZERO_PATTERN;
225: return(0);
226: }
228: /*
229: SNESDAComputeJacobianWithAdifor - This is a universal Jacobian evaluation routine
230: that may be used with SNESSetJacobian() from Fortran as long as the user context has
231: a DA as its first record and DASetLocalAdiforFunction() has been called.
233: Collective on SNES
235: Input Parameters:
236: + snes - the SNES context
237: . X - input vector
238: . J - Jacobian
239: . B - Jacobian used in preconditioner (usally same as J)
240: . flag - indicates if the matrix changed its structure
241: - ptr - optional user-defined context, as set by SNESSetFunction()
243: Level: intermediate
245: .seealso: DASetLocalFunction(), DASetLocalAdicFunction(), SNESSetFunction(), SNESSetJacobian()
247: */
248: int SNESDAComputeJacobianWithAdifor(SNES snes,Vec X,Mat *J,Mat *B,MatStructure *flag,void *ptr)
249: {
250: DA da = *(DA*) ptr;
251: int ierr;
252: Vec localX;
255: DAGetLocalVector(da,&localX);
256: DAGlobalToLocalBegin(da,X,INSERT_VALUES,localX);
257: DAGlobalToLocalEnd(da,X,INSERT_VALUES,localX);
258: DAComputeJacobian1WithAdifor(da,localX,*B,ptr);
259: DARestoreLocalVector(da,&localX);
260: /* Assemble true Jacobian; if it is different */
261: if (*J != *B) {
262: ierr = MatAssemblyBegin(*J,MAT_FINAL_ASSEMBLY);
263: ierr = MatAssemblyEnd(*J,MAT_FINAL_ASSEMBLY);
264: }
265: ierr = MatSetOption(*B,MAT_NEW_NONZERO_LOCATION_ERR);
266: *flag = SAME_NONZERO_PATTERN;
267: return(0);
268: }
270: /*
271: SNESDAComputeJacobian - This is a universal Jacobian evaluation routine for a
272: locally provided Jacobian.
274: Collective on SNES
276: Input Parameters:
277: + snes - the SNES context
278: . X - input vector
279: . J - Jacobian
280: . B - Jacobian used in preconditioner (usally same as J)
281: . flag - indicates if the matrix changed its structure
282: - ptr - optional user-defined context, as set by SNESSetFunction()
284: Level: intermediate
286: .seealso: DASetLocalFunction(), DASetLocalJacobian(), SNESSetFunction(), SNESSetJacobian()
288: */
289: int SNESDAComputeJacobian(SNES snes,Vec X,Mat *J,Mat *B,MatStructure *flag,void *ptr)
290: {
291: DA da = *(DA*) ptr;
292: int ierr;
293: Vec localX;
296: DAGetLocalVector(da,&localX);
297: DAGlobalToLocalBegin(da,X,INSERT_VALUES,localX);
298: DAGlobalToLocalEnd(da,X,INSERT_VALUES,localX);
299: DAComputeJacobian1(da,localX,*B,ptr);
300: DARestoreLocalVector(da,&localX);
301: /* Assemble true Jacobian; if it is different */
302: if (*J != *B) {
303: ierr = MatAssemblyBegin(*J,MAT_FINAL_ASSEMBLY);
304: ierr = MatAssemblyEnd(*J,MAT_FINAL_ASSEMBLY);
305: }
306: ierr = MatSetOption(*B,MAT_NEW_NONZERO_LOCATION_ERR);
307: *flag = SAME_NONZERO_PATTERN;
308: return(0);
309: }
311: int DMMGSolveSNES(DMMG *dmmg,int level)
312: {
313: int ierr,nlevels = dmmg[0]->nlevels,its;
316: dmmg[0]->nlevels = level+1;
317: ierr = SNESSolve(dmmg[level]->snes,dmmg[level]->x,&its);
318: dmmg[0]->nlevels = nlevels;
319: return(0);
320: }
322: /* ===========================================================================================================*/
324: /*@C
325: DMMGSetSNES - Sets the nonlinear function that defines the nonlinear set of equations
326: to be solved using the grid hierarchy.
328: Collective on DMMG
330: Input Parameter:
331: + dmmg - the context
332: . function - the function that defines the nonlinear system
333: - jacobian - optional function to compute Jacobian
335: Options Database:
336: + -dmmg_snes_monitor
337: . -dmmg_jacobian_fd
338: . -dmmg_jacobian_ad
339: . -dmmg_jacobian_mf_fd_operator
340: . -dmmg_jacobian_mf_fd
341: . -dmmg_jacobian_mf_ad_operator
342: - -dmmg_jacobian_mf_ad
344: Level: advanced
346: .seealso DMMGCreate(), DMMGDestroy, DMMGSetSLES(), DMMGSetSNESLocal()
348: @*/
349: int DMMGSetSNES(DMMG *dmmg,int (*function)(SNES,Vec,Vec,void*),int (*jacobian)(SNES,Vec,Mat*,Mat*,MatStructure*,void*))
350: {
351: int ierr,i,nlevels = dmmg[0]->nlevels;
352: PetscTruth snesmonitor,mffdoperator,mffd,fdjacobian;
353: #if defined(PETSC_HAVE_ADIC) && !defined(PETSC_USE_COMPLEX) && !defined(PETSC_USE_SINGLE)
354: PetscTruth mfadoperator,mfad,adjacobian;
355: #endif
356: SLES sles;
357: PetscViewer ascii;
358: MPI_Comm comm;
361: if (!dmmg) SETERRQ(1,"Passing null as DMMG");
362: if (!jacobian) jacobian = DMMGComputeJacobianWithFD;
364: PetscOptionsBegin(dmmg[0]->comm,PETSC_NULL,"DMMG Options","SNES");
365: PetscOptionsName("-dmmg_snes_monitor","Monitor nonlinear convergence","SNESSetMonitor",&snesmonitor);
368: PetscOptionsName("-dmmg_jacobian_fd","Compute sparse Jacobian explicitly with finite differencing","DMMGSetSNES",&fdjacobian);
369: if (fdjacobian) jacobian = DMMGComputeJacobianWithFD;
370: #if defined(PETSC_HAVE_ADIC) && !defined(PETSC_USE_COMPLEX) && !defined(PETSC_USE_SINGLE)
371: PetscOptionsName("-dmmg_jacobian_ad","Compute sparse Jacobian explicitly with ADIC (automatic differentiation)","DMMGSetSNES",&adjacobian);
372: if (adjacobian) jacobian = DMMGComputeJacobianWithAdic;
373: #endif
375: PetscOptionsLogicalGroupBegin("-dmmg_jacobian_mf_fd_operator","Apply Jacobian via matrix free finite differencing","DMMGSetSNES",&mffdoperator);
376: PetscOptionsLogicalGroupEnd("-dmmg_jacobian_mf_fd","Apply Jacobian via matrix free finite differencing even in computing preconditioner","DMMGSetSNES",&mffd);
377: if (mffd) mffdoperator = PETSC_TRUE;
378: #if defined(PETSC_HAVE_ADIC) && !defined(PETSC_USE_COMPLEX) && !defined(PETSC_USE_SINGLE)
379: PetscOptionsLogicalGroupBegin("-dmmg_jacobian_mf_ad_operator","Apply Jacobian via matrix free ADIC (automatic differentiation)","DMMGSetSNES",&mfadoperator);
380: PetscOptionsLogicalGroupEnd("-dmmg_jacobian_mf_ad","Apply Jacobian via matrix free ADIC (automatic differentiation) even in computing preconditioner","DMMGSetSNES",&mfad);
381: if (mfad) mfadoperator = PETSC_TRUE;
382: #endif
383: PetscOptionsEnd();
385: /* create solvers for each level */
386: for (i=0; i<nlevels; i++) {
387: SNESCreate(dmmg[i]->comm,&dmmg[i]->snes);
388: if (snesmonitor) {
389: PetscObjectGetComm((PetscObject)dmmg[i]->snes,&comm);
390: PetscViewerASCIIOpen(comm,"stdout",&ascii);
391: PetscViewerASCIISetTab(ascii,nlevels-i);
392: SNESSetMonitor(dmmg[i]->snes,SNESDefaultMonitor,ascii,(int(*)(void*))PetscViewerDestroy);
393: }
395: if (mffdoperator) {
396: MatCreateSNESMF(dmmg[i]->snes,dmmg[i]->x,&dmmg[i]->J);
397: VecDuplicate(dmmg[i]->x,&dmmg[i]->work1);
398: VecDuplicate(dmmg[i]->x,&dmmg[i]->work2);
399: MatSNESMFSetFunction(dmmg[i]->J,dmmg[i]->work1,function,dmmg[i]);
400: if (mffd) {
401: dmmg[i]->B = dmmg[i]->J;
402: jacobian = DMMGComputeJacobianWithMF;
403: }
404: #if defined(PETSC_HAVE_ADIC) && !defined(PETSC_USE_COMPLEX) && !defined(PETSC_USE_SINGLE)
405: } else if (mfadoperator) {
406: MatRegisterDAAD();
407: MatCreateDAAD((DA)dmmg[i]->dm,&dmmg[i]->J);
408: MatDAADSetCtx(dmmg[i]->J,dmmg[i]->user);
409: if (mfad) {
410: dmmg[i]->B = dmmg[i]->J;
411: jacobian = DMMGComputeJacobianWithMF;
412: }
413: #endif
414: }
415:
416: if (!dmmg[i]->B) {
417: DMGetMatrix(dmmg[i]->dm,MATMPIAIJ,&dmmg[i]->B);
418: }
419: if (!dmmg[i]->J) {
420: dmmg[i]->J = dmmg[i]->B;
421: }
423: SNESGetSLES(dmmg[i]->snes,&sles);
424: DMMGSetUpLevel(dmmg,sles,i+1);
425:
426: /*
427: if the number of levels is > 1 then we want the coarse solve in the grid sequencing to use LU
428: when possible
429: */
430: if (nlevels > 1 && i == 0) {
431: PC pc;
432: SLES csles;
433: PetscTruth flg1,flg2,flg3;
435: SLESGetPC(sles,&pc);
436: MGGetCoarseSolve(pc,&csles);
437: SLESGetPC(csles,&pc);
438: PetscTypeCompare((PetscObject)pc,PCILU,&flg1);
439: PetscTypeCompare((PetscObject)pc,PCSOR,&flg2);
440: PetscTypeCompare((PetscObject)pc,PETSC_NULL,&flg3);
441: if (flg1 || flg2 || flg3) {
442: PCSetType(pc,PCLU);
443: }
444: }
446: SNESSetFromOptions(dmmg[i]->snes);
447: dmmg[i]->solve = DMMGSolveSNES;
448: dmmg[i]->computejacobian = jacobian;
449: dmmg[i]->computefunction = function;
450: }
453: if (jacobian == DMMGComputeJacobianWithFD) {
454: ISColoring iscoloring;
455: for (i=0; i<nlevels; i++) {
456: DMGetColoring(dmmg[i]->dm,IS_COLORING_LOCAL,&iscoloring);
457: MatFDColoringCreate(dmmg[i]->B,iscoloring,&dmmg[i]->fdcoloring);
458: ISColoringDestroy(iscoloring);
459: MatFDColoringSetFunction(dmmg[i]->fdcoloring,(int(*)(void))function,dmmg[i]);
460: MatFDColoringSetFromOptions(dmmg[i]->fdcoloring);
461: }
462: #if defined(PETSC_HAVE_ADIC) && !defined(PETSC_USE_COMPLEX) && !defined(PETSC_USE_SINGLE)
463: } else if (jacobian == DMMGComputeJacobianWithAdic) {
464: for (i=0; i<nlevels; i++) {
465: ISColoring iscoloring;
466: DMGetColoring(dmmg[i]->dm,IS_COLORING_GHOSTED,&iscoloring);
467: MatSetColoring(dmmg[i]->B,iscoloring);
468: ISColoringDestroy(iscoloring);
469: }
470: #endif
471: }
473: for (i=0; i<nlevels; i++) {
474: SNESSetJacobian(dmmg[i]->snes,dmmg[i]->J,dmmg[i]->B,DMMGComputeJacobian_Multigrid,dmmg);
475: SNESSetFunction(dmmg[i]->snes,dmmg[i]->b,function,dmmg[i]);
476: }
478: /* Create interpolation scaling */
479: for (i=1; i<nlevels; i++) {
480: DMGetInterpolationScale(dmmg[i-1]->dm,dmmg[i]->dm,dmmg[i]->R,&dmmg[i]->Rscale);
481: }
483: return(0);
484: }
486: /*@C
487: DMMGSetInitialGuess - Sets the function that computes an initial guess, if not given
488: uses 0.
490: Collective on DMMG and SNES
492: Input Parameter:
493: + dmmg - the context
494: - guess - the function
496: Level: advanced
498: .seealso DMMGCreate(), DMMGDestroy, DMMGSetSLES()
500: @*/
501: int DMMGSetInitialGuess(DMMG *dmmg,int (*guess)(SNES,Vec,void*))
502: {
503: int i,nlevels = dmmg[0]->nlevels;
506: for (i=0; i<nlevels; i++) {
507: dmmg[i]->initialguess = guess;
508: }
509: return(0);
510: }
513: /*M
514: DMMGSetSNESLocal - Sets the local user function that defines the nonlinear set of equations
515: that will use the grid hierarchy and (optionally) its derivative.
517: Collective on DMMG
519: Synopsis:
520: int DMMGSetSNESLocal(DMMG *dmmg,DALocalFunction1 function, DALocalFunction1 jacobian,
521: DALocalFunction1 ad_function, DALocalFunction1 admf_function);
523: Input Parameter:
524: + dmmg - the context
525: . function - the function that defines the nonlinear system
526: . jacobian - function defines the local part of the Jacobian (not currently supported)
527: . ad_function - the name of the function with an ad_ prefix. This is ignored if ADIC is
528: not installed
529: - admf_function - the name of the function with an ad_ prefix. This is ignored if ADIC is
530: not installed
532: Level: intermediate
534: Notes:
535: If ADIC or ADIFOR have been installed, this routine can use ADIC or ADIFOR to compute
536: the derivative; however, that function cannot call other functions except those in
537: standard C math libraries.
539: If ADIC/ADIFOR have not been installed and the Jacobian is not provided, this routine
540: uses finite differencing to approximate the Jacobian.
542: .seealso DMMGCreate(), DMMGDestroy, DMMGSetSLES(), DMMGSetSNES()
544: M*/
546: int DMMGSetSNESLocal_Private(DMMG *dmmg,DALocalFunction1 function,DALocalFunction1 jacobian,DALocalFunction1 ad_function,DALocalFunction1 admf_function)
547: {
548: int ierr,i,nlevels = dmmg[0]->nlevels;
549: int (*computejacobian)(SNES,Vec,Mat*,Mat*,MatStructure*,void*) = 0;
553: if (jacobian) computejacobian = SNESDAComputeJacobian;
554: #if defined(PETSC_HAVE_ADIC) && !defined(PETSC_USE_COMPLEX) && !defined(PETSC_USE_SINGLE)
555: else if (ad_function) computejacobian = DMMGComputeJacobianWithAdic;
556: #endif
558: DMMGSetSNES(dmmg,DMMGFormFunction,computejacobian);
559: for (i=0; i<nlevels; i++) {
560: DASetLocalFunction((DA)dmmg[i]->dm,function);
561: DASetLocalJacobian((DA)dmmg[i]->dm,jacobian);
562: DASetLocalAdicFunction((DA)dmmg[i]->dm,ad_function);
563: DASetLocalAdicMFFunction((DA)dmmg[i]->dm,admf_function);
564: }
565: return(0);
566: }
568: static int DMMGFunctioni(int i,Vec u,PetscScalar* r,void* ctx)
569: {
570: DMMG dmmg = (DMMG)ctx;
571: Vec U = dmmg->lwork1;
572: int ierr;
573: VecScatter gtol;
576: /* copy u into interior part of U */
577: DAGetScatter((DA)dmmg->dm,0,>ol,0);
578: VecScatterBegin(u,U,INSERT_VALUES,SCATTER_FORWARD_LOCAL,gtol);
579: VecScatterEnd(u,U,INSERT_VALUES,SCATTER_FORWARD_LOCAL,gtol);
581: DAFormFunctioni1((DA)dmmg->dm,i,U,r,dmmg->user);
582: return(0);
583: }
585: static int DMMGFunctioniBase(Vec u,void* ctx)
586: {
587: DMMG dmmg = (DMMG)ctx;
588: Vec U = dmmg->lwork1;
589: int ierr;
592: DAGlobalToLocalBegin((DA)dmmg->dm,u,INSERT_VALUES,U);
593: DAGlobalToLocalEnd((DA)dmmg->dm,u,INSERT_VALUES,U);
594: return(0);
595: }
597: 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*))
598: {
599: int ierr,i,nlevels = dmmg[0]->nlevels;
602: CHKMEMQ;
603: for (i=0; i<nlevels; i++) {
604: CHKMEMQ;
605: DASetLocalFunctioni((DA)dmmg[i]->dm,functioni);
606: CHKMEMQ;
607: DASetLocalAdicFunctioni((DA)dmmg[i]->dm,adi);
608: CHKMEMQ;
609: DASetLocalAdicMFFunctioni((DA)dmmg[i]->dm,adimf);
610: CHKMEMQ;
611: MatSNESMFSetFunctioni(dmmg[i]->J,DMMGFunctioni);
612: CHKMEMQ;
613: MatSNESMFSetFunctioniBase(dmmg[i]->J,DMMGFunctioniBase);
614: CHKMEMQ;
615: DACreateLocalVector((DA)dmmg[i]->dm,&dmmg[i]->lwork1);
616: }
617: return(0);
618: }
621: #if defined(PETSC_HAVE_ADIC) && !defined(PETSC_USE_COMPLEX) && !defined(PETSC_USE_SINGLE)
622: EXTERN_C_BEGIN
623: #include "adic/ad_utils.h"
624: EXTERN_C_END
626: int PetscADView(int N,int nc,double *ptr,PetscViewer viewer)
627: {
628: int i,j,nlen = PetscADGetDerivTypeSize();
629: char *cptr = (char*)ptr;
630: double *values;
633: for (i=0; i<N; i++) {
634: printf("Element %d value %g derivatives: ",i,*(double*)cptr);
635: values = PetscADGetGradArray(cptr);
636: for (j=0; j<nc; j++) {
637: printf("%g ",*values++);
638: }
639: printf("n");
640: cptr += nlen;
641: }
643: return(0);
644: }
646: #endif