Actual source code: ex12.c
petsc-3.4.5 2014-06-29
2: /*
3: Simple example to show how PETSc programs can be run from MATLAB.
4: See ex12.m.
5: */
7: static char help[] = "Solves the one dimensional heat equation.\n\n";
9: #include <petscdmda.h>
13: int main(int argc,char **argv)
14: {
15: PetscMPIInt rank,size;
16: PetscInt M = 14,time_steps = 20,w=1,s=1,localsize,j,i,mybase,myend;
18: DM da;
19: Vec local,global,copy;
20: PetscScalar *localptr,*copyptr;
21: PetscReal h,k;
23: PetscInitialize(&argc,&argv,(char*)0,help);
25: PetscOptionsGetInt(NULL,"-M",&M,NULL);
26: PetscOptionsGetInt(NULL,"-time",&time_steps,NULL);
28: /* Set up the array */
29: DMDACreate1d(PETSC_COMM_WORLD,DMDA_BOUNDARY_NONE,M,w,s,NULL,&da);
30: DMCreateGlobalVector(da,&global);
31: DMCreateLocalVector(da,&local);
32: MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
33: MPI_Comm_size(PETSC_COMM_WORLD,&size);
35: /* Make copy of local array for doing updates */
36: VecDuplicate(local,©);
37: VecGetArray (copy,©ptr);
40: /* determine starting point of each processor */
41: VecGetOwnershipRange(global,&mybase,&myend);
43: /* Initialize the Array */
44: VecGetLocalSize (local,&localsize);
45: VecGetArray (local,&localptr);
47: localptr[0] = copyptr[0] = 0.0;
49: localptr[localsize-1] = copyptr[localsize-1] = 1.0;
50: for (i=1; i<localsize-1; i++) {
51: j = (i-1) + mybase;
53: localptr[i] = sin((PETSC_PI*j*6)/((PetscReal)M) + 1.2 * sin((PETSC_PI*j*2)/((PetscReal)M))) * 4+4;
54: }
56: VecRestoreArray (copy,©ptr);
57: VecRestoreArray(local,&localptr);
58: DMLocalToGlobalBegin(da,local,INSERT_VALUES,global);
59: DMLocalToGlobalEnd(da,local,INSERT_VALUES,global);
61: /* Assign Parameters */
62: h= 1.0/M;
63: k= h*h/2.2;
65: for (j=0; j<time_steps; j++) {
67: /* Global to Local */
68: DMGlobalToLocalBegin(da,global,INSERT_VALUES,local);
69: DMGlobalToLocalEnd(da,global,INSERT_VALUES,local);
71: /*Extract local array */
72: VecGetArray(local,&localptr);
73: VecGetArray (copy,©ptr);
75: /* Update Locally - Make array of new values */
76: /* Note: I don't do anything for the first and last entry */
77: for (i=1; i< localsize-1; i++) {
78: copyptr[i] = localptr[i] + (k/(h*h)) * (localptr[i+1]-2.0*localptr[i]+localptr[i-1]);
79: }
81: VecRestoreArray (copy,©ptr);
82: VecRestoreArray(local,&localptr);
84: /* Local to Global */
85: DMLocalToGlobalBegin(da,copy,INSERT_VALUES,global);
86: DMLocalToGlobalEnd(da,copy,INSERT_VALUES,global);
88: /* View Wave */
89: /* Set Up Display to Show Heat Graph */
90: #if defined(PETSC_USE_SOCKET_VIEWER)
91: VecView(global,PETSC_VIEWER_SOCKET_WORLD);
92: #endif
93: }
95: VecDestroy(©);
96: VecDestroy(&local);
97: VecDestroy(&global);
98: DMDestroy(&da);
99: PetscFinalize();
100: return 0;
101: }