1 static const char help[] = "Minimum surface problem in 2D.\n\
2 Uses 2-dimensional distributed arrays.\n\
3 \n\
4 Solves the linear systems via multilevel methods \n\
5 \n\n";
6
7 /*
8
9 This example models the partial differential equation
10
11 - Div((1 + ||GRAD T||^2)^(1/2) (GRAD T)) = 0.
12
13 in the unit square, which is uniformly discretized in each of x and
14 y in this simple encoding. The degrees of freedom are vertex centered
15
16 A finite difference approximation with the usual 5-point stencil
17 is used to discretize the boundary value problem to obtain a
18 nonlinear system of equations.
19
20 */
21
22 #include <petscsnes.h>
23 #include <petscdm.h>
24 #include <petscdmda.h>
25
26 extern PetscErrorCode FormFunctionLocal(DMDALocalInfo *, PetscScalar **, PetscScalar **, void *);
27
main(int argc,char ** argv)28 int main(int argc, char **argv)
29 {
30 SNES snes;
31 PetscInt its, lits;
32 PetscReal litspit;
33 DM da;
34
35 PetscFunctionBeginUser;
36 PetscCall(PetscInitialize(&argc, &argv, NULL, help));
37 /*
38 Set the DMDA (grid structure) for the grids.
39 */
40 PetscCall(DMDACreate2d(PETSC_COMM_WORLD, DM_BOUNDARY_NONE, DM_BOUNDARY_NONE, DMDA_STENCIL_STAR, 5, 5, PETSC_DECIDE, PETSC_DECIDE, 1, 1, 0, 0, &da));
41 PetscCall(DMSetFromOptions(da));
42 PetscCall(DMSetUp(da));
43 PetscCall(DMDASNESSetFunctionLocal(da, INSERT_VALUES, (PetscErrorCode (*)(DMDALocalInfo *, void *, void *, void *))FormFunctionLocal, NULL));
44 PetscCall(SNESCreate(PETSC_COMM_WORLD, &snes));
45 PetscCall(SNESSetDM(snes, da));
46 PetscCall(DMDestroy(&da));
47
48 PetscCall(SNESSetFromOptions(snes));
49
50 PetscCall(SNESSolve(snes, 0, 0));
51 PetscCall(SNESGetIterationNumber(snes, &its));
52 PetscCall(SNESGetLinearSolveIterations(snes, &lits));
53 litspit = ((PetscReal)lits) / ((PetscReal)its);
54 PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Number of SNES iterations = %" PetscInt_FMT "\n", its));
55 PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Number of Linear iterations = %" PetscInt_FMT "\n", lits));
56 PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Average Linear its / SNES = %e\n", (double)litspit));
57
58 PetscCall(SNESDestroy(&snes));
59 PetscCall(PetscFinalize());
60 return 0;
61 }
62
FormFunctionLocal(DMDALocalInfo * info,PetscScalar ** t,PetscScalar ** f,void * ptr)63 PetscErrorCode FormFunctionLocal(DMDALocalInfo *info, PetscScalar **t, PetscScalar **f, void *ptr)
64 {
65 PetscInt i, j;
66 PetscScalar hx, hy;
67 PetscScalar gradup, graddown, gradleft, gradright, gradx, grady;
68 PetscScalar coeffup, coeffdown, coeffleft, coeffright;
69
70 PetscFunctionBeginUser;
71 hx = 1.0 / (PetscReal)(info->mx - 1);
72 hy = 1.0 / (PetscReal)(info->my - 1);
73
74 /* Evaluate function */
75 for (j = info->ys; j < info->ys + info->ym; j++) {
76 for (i = info->xs; i < info->xs + info->xm; i++) {
77 if (i == 0 || i == info->mx - 1 || j == 0 || j == info->my - 1) {
78 f[j][i] = t[j][i] - (1.0 - (2.0 * hx * (PetscReal)i - 1.0) * (2.0 * hx * (PetscReal)i - 1.0));
79 } else {
80 gradup = (t[j + 1][i] - t[j][i]) / hy;
81 graddown = (t[j][i] - t[j - 1][i]) / hy;
82 gradright = (t[j][i + 1] - t[j][i]) / hx;
83 gradleft = (t[j][i] - t[j][i - 1]) / hx;
84
85 gradx = .5 * (t[j][i + 1] - t[j][i - 1]) / hx;
86 grady = .5 * (t[j + 1][i] - t[j - 1][i]) / hy;
87
88 coeffup = 1.0 / PetscSqrtScalar(1.0 + gradup * gradup + gradx * gradx);
89 coeffdown = 1.0 / PetscSqrtScalar(1.0 + graddown * graddown + gradx * gradx);
90
91 coeffleft = 1.0 / PetscSqrtScalar(1.0 + gradleft * gradleft + grady * grady);
92 coeffright = 1.0 / PetscSqrtScalar(1.0 + gradright * gradright + grady * grady);
93
94 f[j][i] = (coeffup * gradup - coeffdown * graddown) * hx + (coeffright * gradright - coeffleft * gradleft) * hy;
95 }
96 }
97 }
98 PetscFunctionReturn(PETSC_SUCCESS);
99 }
100
101 /*TEST
102
103 test:
104 args: -pc_type mg -da_refine 1 -ksp_type fgmres
105
106 test:
107 suffix: 2
108 nsize: 2
109 args: -pc_type mg -da_refine 1 -ksp_type fgmres
110
111 test:
112 suffix: 3
113 nsize: 2
114 args: -pc_type mg -da_refine 1 -ksp_type fgmres -snes_type newtontrdc -snes_trdc_use_cauchy false
115
116 test:
117 suffix: 4
118 nsize: 2
119 args: -pc_type mg -da_refine 1 -ksp_type fgmres -snes_type newtontrdc
120 filter: sed -e "s/SNES iterations = 1[0-3]/SNES iterations = 13/g" |sed -e "s/Linear iterations = 2[7-9]/Linear iterations = 29/g" |sed -e "s/Linear iterations = 3[0-1]/Linear iterations = 29/g"
121
122 TEST*/
123