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