Actual source code: p2.c
2: static char help[] = "Model single-physics solver. Modified from ex19.c \n\\n";
4: /* ------------------------------------------------------------------------
6: See ex19.c for discussion of the problem
8: ------------------------------------------------------------------------- */
9: #include mp.h
17: int main(int argc,char **argv)
18: {
19: DMMG *dmmg; /* multilevel grid structure */
20: AppCtx user; /* user-defined work context */
21: PetscInt mx,my,its;
23: MPI_Comm comm;
24: SNES snes;
25: DA da2;
27: PetscInitialize(&argc,&argv,(char *)0,help);
28: PetscLogEventRegister("FormFunc2", 0,&EVENT_FORMFUNCTIONLOCAL2);
29: comm = PETSC_COMM_WORLD;
31: /* Problem parameters (velocity of lid, prandtl, and grashof numbers) */
32: PetscOptionsGetReal(PETSC_NULL,"-lidvelocity",&user.lidvelocity,PETSC_NULL);
33: PetscOptionsGetReal(PETSC_NULL,"-prandtl",&user.prandtl,PETSC_NULL);
34: PetscOptionsGetReal(PETSC_NULL,"-grashof",&user.grashof,PETSC_NULL);
36: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
37: Create user context, set problem data, create vector data structures.
38: Also, compute the initial guess.
39: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
41: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
42: Setup Physics 2:
43: - Lap(T) + PR*Div([U*T,V*T]) = 0
44: where U and V are given by the given x.u and x.v
45: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
46: DACreate2d(comm,DA_NONPERIODIC,DA_STENCIL_STAR,-4,-4,PETSC_DECIDE,PETSC_DECIDE,1,1,0,0,&da2);
47: DASetFieldName(da2,0,"temperature");
49: /* Create the solver object and attach the grid/physics info */
50: DMMGCreate(comm,1,&user,&dmmg);
51: DMMGSetDM(dmmg,(DM)da2);
52: DMMGSetISColoringType(dmmg,IS_COLORING_GLOBAL);
54: DMMGSetInitialGuess(dmmg,FormInitialGuess);
55: DMMGSetSNES(dmmg,FormFunction,0);
56: DMMGSetFromOptions(dmmg);
58: DAGetInfo(da2,PETSC_NULL,&mx,&my,0,0,0,0,0,0,0,0);
59: user.lidvelocity = 1.0/(mx*my);
60: user.prandtl = 1.0;
61: user.grashof = 1.0;
63: /* Solve the nonlinear system */
64: DMMGSolve(dmmg);
65: snes = DMMGGetSNES(dmmg);
66: SNESGetIterationNumber(snes,&its);
67: PetscPrintf(comm,"Physics 2: Number of Newton iterations = %D\n\n", its);
69: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
70: Free spaces
71: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
72: DADestroy(da2);
73: DMMGDestroy(dmmg);
74: PetscFinalize();
75: return 0;
76: }
80: PetscErrorCode FormInitialGuess(DMMG dmmg,Vec X)
81: {
83: AppCtx *user = (AppCtx*)dmmg->user;
84: DA da2 = (DA)dmmg->dm;
85: Field2 **x2;
86: DALocalInfo info2;
89: /* Access the arrays inside of X */
90: DAVecGetArray(da2,X,(void**)&x2);
92: DAGetLocalInfo(da2,&info2);
94: /* Evaluate local user provided function */
95: FormInitialGuessLocal2(&info2,x2,user);
97: DAVecRestoreArray(da2,X,(void**)&x2);
98: return(0);
99: }
103: PetscErrorCode FormFunction(SNES snes,Vec X,Vec F,void *ctx)
104: {
106: DMMG dmmg = (DMMG)ctx;
107: AppCtx *user = (AppCtx*)dmmg->user;
108: DA da2 = (DA)dmmg->dm;
109: DALocalInfo info2;
110: Field2 **x2,**f2;
111: Vec X2;
114: DAGetLocalInfo(da2,&info2);
116: /* Get local vectors to hold ghosted parts of X */
117: DAGetLocalVector(da2,&X2);
118: DAGlobalToLocalBegin(da2,X,INSERT_VALUES,X2);
119: DAGlobalToLocalEnd(da2,X,INSERT_VALUES,X2);
121: /* Access the array inside of X1 */
122: DAVecGetArray(da2,X2,(void**)&x2);
124: /* Access the subvectors in F.
125: These are not ghosted so directly access the memory locations in F */
126: DAVecGetArray(da2,F,(void**)&f2);
128: /* Evaluate local user provided function */
129: FormFunctionLocal2(&info2,0,x2,f2,(void**)user);
131: DAVecRestoreArray(da2,F,(void**)&f2);
132: DAVecRestoreArray(da2,X2,(void**)&x2);
133: DARestoreLocalVector(da2,&X2);
134: return(0);
135: }