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 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 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(0); 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[1-3]/SNES iterations = 13/g" |sed -e "s/Linear iterations = 2[8-9]/Linear iterations = 29/g" |sed -e "s/Linear iterations = 3[0-1]/Linear iterations = 29/g" 121 122 TEST*/ 123