Actual source code: ex17.c
petsc-3.3-p7 2013-05-11
1: static const char help[] = "Newton's method to solve a two-variable system, sequentially.\n"
2: "The same problem is solved twice - i) fully assembled system + ii) block system\n\n";
4: /*T
5: Concepts: SNES^basic uniprocessor example, block objects
6: Processors: 1
7: T*/
9: /*
10: Include "petscsnes.h" so that we can use SNES solvers. Note that this
11: file automatically includes:
12: petscsys.h - base PETSc routines petscvec.h - vectors
13: petscsys.h - system routines petscmat.h - matrices
14: petscis.h - index sets petscksp.h - Krylov subspace methods
15: petscviewer.h - viewers petscpc.h - preconditioners
16: petscksp.h - linear solvers
17: */
18: #include <petscsnes.h>
20: /*
21: This example is block version of the test found at
22: ${PETSC_DIR}/src/snes/examples/tutorials/ex1.c
23: In this test we replace the Jacobian systems
24: [J]{x} = {F}
25: where
27: [J] = ( j_00, j_01 ), {x} = ( x_0, x_1 )^T, {F} = ( f_0, f_1 )^T
28: ( j_10, j_11 )
29: where [J] \in \mathbb^{2 \times 2}, {x},{F} \in \mathbb^{2 \times 1},
31: with a block system in which each block is of length 1.
32: i.e. The block system is thus
34: [J] = ( [j00], [j01] ), {x} = ( {x0}, {x1} )^T, {F} = ( {f0}, {f1} )^T
35: ( [j10], [j11] )
36: where
37: [j00], [j01], [j10], [j11] \in \mathbb^{1 \times 1}
38: {x0}, {x1}, {f0}, {f1} \in \mathbb^{1 \times 1}
40: In practice we would not bother defing blocks of size one, and would instead assemble the
41: full system. This is just a simple test to illustrate how to manipulate the blocks and
42: to confirm the implementation is correct.
43: */
45: /*
46: User-defined routines
47: */
48: static PetscErrorCode FormJacobian1(SNES,Vec,Mat*,Mat*,MatStructure*,void*);
49: static PetscErrorCode FormFunction1(SNES,Vec,Vec,void*);
50: static PetscErrorCode FormJacobian2(SNES,Vec,Mat*,Mat*,MatStructure*,void*);
51: static PetscErrorCode FormFunction2(SNES,Vec,Vec,void*);
52: static PetscErrorCode FormJacobian1_block(SNES,Vec,Mat*,Mat*,MatStructure*,void*);
53: static PetscErrorCode FormFunction1_block(SNES,Vec,Vec,void*);
54: static PetscErrorCode FormJacobian2_block(SNES,Vec,Mat*,Mat*,MatStructure*,void*);
55: static PetscErrorCode FormFunction2_block(SNES,Vec,Vec,void*);
60: static PetscErrorCode assembled_system(void)
61: {
62: SNES snes; /* nonlinear solver context */
63: KSP ksp; /* linear solver context */
64: PC pc; /* preconditioner context */
65: Vec x,r; /* solution, residual vectors */
66: Mat J; /* Jacobian matrix */
68: PetscInt its;
69: PetscScalar pfive = .5,*xx;
70: PetscBool flg;
73: PetscPrintf( PETSC_COMM_WORLD, "\n\n========================= Assembled system =========================\n\n" );
75: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
76: Create nonlinear solver context
77: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
79: SNESCreate(PETSC_COMM_WORLD,&snes);
81: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
82: Create matrix and vector data structures; set corresponding routines
83: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
85: /*
86: Create vectors for solution and nonlinear function
87: */
88: VecCreateSeq(PETSC_COMM_SELF,2,&x);
89: VecDuplicate(x,&r);
91: /*
92: Create Jacobian matrix data structure
93: */
94: MatCreate(PETSC_COMM_SELF,&J);
95: MatSetSizes(J,PETSC_DECIDE,PETSC_DECIDE,2,2);
96: MatSetFromOptions(J);
97: MatSetUp(J);
99: PetscOptionsHasName(PETSC_NULL,"-hard",&flg);
100: if (!flg) {
101: /*
102: Set function evaluation routine and vector.
103: */
104: SNESSetFunction(snes,r,FormFunction1,PETSC_NULL);
106: /*
107: Set Jacobian matrix data structure and Jacobian evaluation routine
108: */
109: SNESSetJacobian(snes,J,J,FormJacobian1,PETSC_NULL);
110: } else {
111: SNESSetFunction(snes,r,FormFunction2,PETSC_NULL);
112: SNESSetJacobian(snes,J,J,FormJacobian2,PETSC_NULL);
113: }
115: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
116: Customize nonlinear solver; set runtime options
117: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
119: /*
120: Set linear solver defaults for this problem. By extracting the
121: KSP, KSP, and PC contexts from the SNES context, we can then
122: directly call any KSP, KSP, and PC routines to set various options.
123: */
124: SNESGetKSP(snes,&ksp);
125: KSPGetPC(ksp,&pc);
126: PCSetType(pc,PCNONE);
127: KSPSetTolerances(ksp,1.e-4,PETSC_DEFAULT,PETSC_DEFAULT,20);
129: /*
130: Set SNES/KSP/KSP/PC runtime options, e.g.,
131: -snes_view -snes_monitor -ksp_type <ksp> -pc_type <pc>
132: These options will override those specified above as long as
133: SNESSetFromOptions() is called _after_ any other customization
134: routines.
135: */
136: SNESSetFromOptions(snes);
138: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
139: Evaluate initial guess; then solve nonlinear system
140: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
141: if (!flg) {
142: VecSet(x,pfive);
143: } else {
144: VecGetArray(x,&xx);
145: xx[0] = 2.0; xx[1] = 3.0;
146: VecRestoreArray(x,&xx);
147: }
148: /*
149: Note: The user should initialize the vector, x, with the initial guess
150: for the nonlinear solver prior to calling SNESSolve(). In particular,
151: to employ an initial guess of zero, the user should explicitly set
152: this vector to zero by calling VecSet().
153: */
155: SNESSolve(snes,PETSC_NULL,x);
156: SNESGetIterationNumber(snes,&its);
157: if (flg) {
158: Vec f;
159: VecView(x,PETSC_VIEWER_STDOUT_WORLD);
160: SNESGetFunction(snes,&f,0,0);
161: VecView(r,PETSC_VIEWER_STDOUT_WORLD);
162: }
164: PetscPrintf(PETSC_COMM_SELF,"number of SNES iterations = %D\n\n",its);
166: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
167: Free work space. All PETSc objects should be destroyed when they
168: are no longer needed.
169: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
171: VecDestroy(&x); VecDestroy(&r);
172: MatDestroy(&J); SNESDestroy(&snes);
173: return(0);
174: }
175: /* ------------------------------------------------------------------- */
178: /*
179: FormFunction1 - Evaluates nonlinear function, F(x).
181: Input Parameters:
182: . snes - the SNES context
183: . x - input vector
184: . dummy - optional user-defined context (not used here)
186: Output Parameter:
187: . f - function vector
188: */
189: static PetscErrorCode FormFunction1(SNES snes,Vec x,Vec f,void *dummy)
190: {
192: PetscScalar *xx,*ff;
195: /*
196: Get pointers to vector data.
197: - For default PETSc vectors, VecGetArray() returns a pointer to
198: the data array. Otherwise, the routine is implementation dependent.
199: - You MUST call VecRestoreArray() when you no longer need access to
200: the array.
201: */
202: VecGetArray(x,&xx);
203: VecGetArray(f,&ff);
205: /*
206: Compute function
207: */
208: ff[0] = xx[0]*xx[0] + xx[0]*xx[1] - 3.0;
209: ff[1] = xx[0]*xx[1] + xx[1]*xx[1] - 6.0;
212: /*
213: Restore vectors
214: */
215: VecRestoreArray(x,&xx);
216: VecRestoreArray(f,&ff);
217: return(0);
218: }
219: /* ------------------------------------------------------------------- */
222: /*
223: FormJacobian1 - Evaluates Jacobian matrix.
225: Input Parameters:
226: . snes - the SNES context
227: . x - input vector
228: . dummy - optional user-defined context (not used here)
230: Output Parameters:
231: . jac - Jacobian matrix
232: . B - optionally different preconditioning matrix
233: . flag - flag indicating matrix structure
234: */
235: static PetscErrorCode FormJacobian1(SNES snes,Vec x,Mat *jac,Mat *B,MatStructure *flag,void *dummy)
236: {
237: PetscScalar *xx,A[4];
239: PetscInt idx[2] = {0,1};
242: /*
243: Get pointer to vector data
244: */
245: VecGetArray(x,&xx);
247: /*
248: Compute Jacobian entries and insert into matrix.
249: - Since this is such a small problem, we set all entries for
250: the matrix at once.
251: */
252: A[0] = 2.0*xx[0] + xx[1]; A[1] = xx[0];
253: A[2] = xx[1]; A[3] = xx[0] + 2.0*xx[1];
254: MatSetValues(*jac,2,idx,2,idx,A,INSERT_VALUES);
255: *flag = SAME_NONZERO_PATTERN;
257: /*
258: Restore vector
259: */
260: VecRestoreArray(x,&xx);
262: /*
263: Assemble matrix
264: */
265: MatAssemblyBegin(*jac,MAT_FINAL_ASSEMBLY);
266: MatAssemblyEnd(*jac,MAT_FINAL_ASSEMBLY);
267: return(0);
268: }
271: /* ------------------------------------------------------------------- */
274: static PetscErrorCode FormFunction2(SNES snes,Vec x,Vec f,void *dummy)
275: {
277: PetscScalar *xx,*ff;
280: /*
281: Get pointers to vector data.
282: - For default PETSc vectors, VecGetArray() returns a pointer to
283: the data array. Otherwise, the routine is implementation dependent.
284: - You MUST call VecRestoreArray() when you no longer need access to
285: the array.
286: */
287: VecGetArray(x,&xx);
288: VecGetArray(f,&ff);
290: /*
291: Compute function
292: */
293: ff[0] = PetscSinScalar(3.0*xx[0]) + xx[0];
294: ff[1] = xx[1];
296: /*
297: Restore vectors
298: */
299: VecRestoreArray(x,&xx);
300: VecRestoreArray(f,&ff);
301: return(0);
302: }
303: /* ------------------------------------------------------------------- */
306: static PetscErrorCode FormJacobian2(SNES snes,Vec x,Mat *jac,Mat *B,MatStructure *flag,void *dummy)
307: {
308: PetscScalar *xx,A[4];
310: PetscInt idx[2] = {0,1};
313: /*
314: Get pointer to vector data
315: */
316: VecGetArray(x,&xx);
318: /*
319: Compute Jacobian entries and insert into matrix.
320: - Since this is such a small problem, we set all entries for
321: the matrix at once.
322: */
323: A[0] = 3.0*PetscCosScalar(3.0*xx[0]) + 1.0; A[1] = 0.0;
324: A[2] = 0.0; A[3] = 1.0;
325: MatSetValues(*jac,2,idx,2,idx,A,INSERT_VALUES);
326: *flag = SAME_NONZERO_PATTERN;
328: /*
329: Restore vector
330: */
331: VecRestoreArray(x,&xx);
333: /*
334: Assemble matrix
335: */
336: MatAssemblyBegin(*jac,MAT_FINAL_ASSEMBLY);
337: MatAssemblyEnd(*jac,MAT_FINAL_ASSEMBLY);
338: return(0);
339: }
343: static int block_system(void)
344: {
345: SNES snes; /* nonlinear solver context */
346: KSP ksp; /* linear solver context */
347: PC pc; /* preconditioner context */
348: Vec x,r; /* solution, residual vectors */
349: Mat J; /* Jacobian matrix */
351: PetscInt its;
352: PetscScalar pfive = .5;
353: PetscBool flg;
355: Mat j11, j12, j21, j22;
356: Vec x1, x2, r1, r2;
357: Vec bv;
358: Vec bx[2];
359: Mat bA[2][2];
362: PetscPrintf( PETSC_COMM_WORLD, "\n\n========================= Block system =========================\n\n" );
364: SNESCreate(PETSC_COMM_WORLD,&snes);
366: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
367: Create matrix and vector data structures; set corresponding routines
368: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
370: /*
371: Create sub vectors for solution and nonlinear function
372: */
373: VecCreateSeq(PETSC_COMM_SELF,1,&x1);
374: VecDuplicate(x1,&r1);
376: VecCreateSeq(PETSC_COMM_SELF,1,&x2);
377: VecDuplicate(x2,&r2);
379: /*
380: Create the block vectors
381: */
382: bx[0] = x1;
383: bx[1] = x2;
384: VecCreateNest(PETSC_COMM_WORLD,2,PETSC_NULL,bx,&x);
385: VecAssemblyBegin(x);
386: VecAssemblyEnd(x);
387: VecDestroy(&x1);
388: VecDestroy(&x2);
390: bx[0] = r1;
391: bx[1] = r2;
392: VecCreateNest(PETSC_COMM_WORLD,2,PETSC_NULL,bx,&r);
393: VecDestroy(&r1);
394: VecDestroy(&r2);
395: VecAssemblyBegin(r);
396: VecAssemblyEnd(r);
398: /*
399: Create sub Jacobian matrix data structure
400: */
401: MatCreate( PETSC_COMM_WORLD, &j11 );
402: MatSetSizes( j11, 1, 1, 1, 1 );
403: MatSetType( j11, MATSEQAIJ );
404: MatSetUp(j11);
406: MatCreate( PETSC_COMM_WORLD, &j12 );
407: MatSetSizes( j12, 1, 1, 1, 1 );
408: MatSetType( j12, MATSEQAIJ );
409: MatSetUp(j12);
411: MatCreate( PETSC_COMM_WORLD, &j21 );
412: MatSetSizes( j21, 1, 1, 1, 1 );
413: MatSetType( j21, MATSEQAIJ );
414: MatSetUp(j21);
416: MatCreate( PETSC_COMM_WORLD, &j22 );
417: MatSetSizes( j22, PETSC_DECIDE, PETSC_DECIDE, 1, 1 );
418: MatSetType( j22, MATSEQAIJ );
419: MatSetUp(j22);
420: /*
421: Create block Jacobian matrix data structure
422: */
423: bA[0][0] = j11;
424: bA[0][1] = j12;
425: bA[1][0] = j21;
426: bA[1][1] = j22;
427: MatCreateNest(PETSC_COMM_WORLD,2,PETSC_NULL,2,PETSC_NULL,&bA[0][0],&J);
428: MatSetUp(J);
429: MatNestSetVecType(J,VECNEST);
430: MatDestroy(&j11);
431: MatDestroy(&j12);
432: MatDestroy(&j21);
433: MatDestroy(&j22);
435: MatAssemblyBegin(J,MAT_FINAL_ASSEMBLY);
436: MatAssemblyEnd(J,MAT_FINAL_ASSEMBLY);
438: PetscOptionsHasName(PETSC_NULL,"-hard",&flg);
439: if (!flg) {
440: /*
441: Set function evaluation routine and vector.
442: */
443: SNESSetFunction(snes,r,FormFunction1_block,PETSC_NULL);
445: /*
446: Set Jacobian matrix data structure and Jacobian evaluation routine
447: */
448: SNESSetJacobian(snes,J,J,FormJacobian1_block,PETSC_NULL);
449: } else {
450: SNESSetFunction(snes,r,FormFunction2_block,PETSC_NULL);
451: SNESSetJacobian(snes,J,J,FormJacobian2_block,PETSC_NULL);
452: }
454: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
455: Customize nonlinear solver; set runtime options
456: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
458: /*
459: Set linear solver defaults for this problem. By extracting the
460: KSP, KSP, and PC contexts from the SNES context, we can then
461: directly call any KSP, KSP, and PC routines to set various options.
462: */
463: SNESGetKSP(snes,&ksp);
464: KSPGetPC(ksp,&pc);
465: PCSetType(pc,PCNONE);
466: KSPSetTolerances(ksp,1.e-4,PETSC_DEFAULT,PETSC_DEFAULT,20);
468: /*
469: Set SNES/KSP/KSP/PC runtime options, e.g.,
470: -snes_view -snes_monitor -ksp_type <ksp> -pc_type <pc>
471: These options will override those specified above as long as
472: SNESSetFromOptions() is called _after_ any other customization
473: routines.
474: */
475: SNESSetFromOptions(snes);
477: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
478: Evaluate initial guess; then solve nonlinear system
479: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
480: if (!flg) {
481: VecSet(x,pfive);
482: } else {
483: Vec *vecs;
484: VecNestGetSubVecs( x, PETSC_NULL, &vecs );
485: bv = vecs[0];
486: // VecBlockGetSubVec( x, 0, &bv );
487: VecSetValue( bv, 0, 2.0, INSERT_VALUES ); /* xx[0] = 2.0; */
488: VecAssemblyBegin(bv);
489: VecAssemblyEnd(bv);
491: // VecBlockGetSubVec( x, 1, &bv );
492: bv = vecs[1];
493: VecSetValue( bv, 0, 3.0, INSERT_VALUES ); /* xx[1] = 3.0; */
494: VecAssemblyBegin(bv);
495: VecAssemblyEnd(bv);
496: }
497: /*
498: Note: The user should initialize the vector, x, with the initial guess
499: for the nonlinear solver prior to calling SNESSolve(). In particular,
500: to employ an initial guess of zero, the user should explicitly set
501: this vector to zero by calling VecSet().
502: */
503: SNESSolve(snes,PETSC_NULL,x);
504: SNESGetIterationNumber(snes,&its);
505: if (flg) {
506: Vec f;
507: VecView(x,PETSC_VIEWER_STDOUT_WORLD);
508: SNESGetFunction(snes,&f,0,0);
509: VecView(r,PETSC_VIEWER_STDOUT_WORLD);
510: }
512: PetscPrintf(PETSC_COMM_SELF,"number of SNES iterations = %D\n\n",its);
514: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
515: Free work space. All PETSc objects should be destroyed when they
516: are no longer needed.
517: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
518: VecDestroy(&x); VecDestroy(&r);
519: MatDestroy(&J); SNESDestroy(&snes);
521: return(0);
522: }
523: /* ------------------------------------------------------------------- */
526: static PetscErrorCode FormFunction1_block(SNES snes,Vec x,Vec f,void *dummy)
527: {
528: Vec *xx, *ff, x1,x2, f1,f2;
529: PetscScalar ff_0, ff_1;
530: PetscScalar xx_0, xx_1;
531: PetscInt index,nb;
535: /* get blocks for function */
536: VecNestGetSubVecs( f, &nb, &ff );
537: f1 = ff[0]; f2 = ff[1];
539: /* get blocks for solution */
540: VecNestGetSubVecs( x, &nb, &xx );
541: x1 = xx[0]; x2 = xx[1];
543: /* get solution values */
544: index = 0;
545: VecGetValues( x1,1, &index, &xx_0 );
546: VecGetValues( x2,1, &index, &xx_1 );
548: /* Compute function */
549: ff_0 = xx_0*xx_0 + xx_0*xx_1 - 3.0;
550: ff_1 = xx_0*xx_1 + xx_1*xx_1 - 6.0;
552: /* set function values */
553: VecSetValue( f1, index, ff_0, INSERT_VALUES );
555: VecSetValue( f2, index, ff_1, INSERT_VALUES );
557: VecAssemblyBegin(f);
558: VecAssemblyEnd(f);
560: return(0);
561: }
562: /* ------------------------------------------------------------------- */
565: static PetscErrorCode FormJacobian1_block(SNES snes,Vec x,Mat *jac,Mat *B,MatStructure *flag,void *dummy)
566: {
567: Vec *xx, x1,x2;
568: PetscScalar xx_0, xx_1;
569: PetscInt index,nb;
570: PetscScalar A_00, A_01, A_10, A_11;
571: Mat j11, j12, j21, j22;
572: Mat **mats;
576: /* get blocks for solution */
577: VecNestGetSubVecs( x, &nb, &xx );
578: x1 = xx[0]; x2 = xx[1];
580: /* get solution values */
581: index = 0;
582: VecGetValues( x1,1, &index, &xx_0 );
583: VecGetValues( x2,1, &index, &xx_1 );
585: /* get block matrices */
586: MatNestGetSubMats(*jac,PETSC_NULL,PETSC_NULL,&mats);
587: j11 = mats[0][0];
588: j12 = mats[0][1];
589: j21 = mats[1][0];
590: j22 = mats[1][1];
592: /* compute jacobian entries */
593: A_00 = 2.0*xx_0 + xx_1;
594: A_01 = xx_0;
595: A_10 = xx_1;
596: A_11 = xx_0 + 2.0*xx_1;
598: /* set jacobian values */
599: MatSetValue( j11, 0,0, A_00, INSERT_VALUES);
600: MatSetValue( j12, 0,0, A_01, INSERT_VALUES);
601: MatSetValue( j21, 0,0, A_10, INSERT_VALUES);
602: MatSetValue( j22, 0,0, A_11, INSERT_VALUES);
604: *flag = SAME_NONZERO_PATTERN;
606: /* Assemble sub matrix */
607: MatAssemblyBegin(*jac,MAT_FINAL_ASSEMBLY);
608: MatAssemblyEnd(*jac,MAT_FINAL_ASSEMBLY);
609: return(0);
610: }
613: /* ------------------------------------------------------------------- */
616: static PetscErrorCode FormFunction2_block(SNES snes,Vec x,Vec f,void *dummy)
617: {
619: PetscScalar *xx,*ff;
622: /*
623: Get pointers to vector data.
624: - For default PETSc vectors, VecGetArray() returns a pointer to
625: the data array. Otherwise, the routine is implementation dependent.
626: - You MUST call VecRestoreArray() when you no longer need access to
627: the array.
628: */
629: VecGetArray(x,&xx);
630: VecGetArray(f,&ff);
632: /*
633: Compute function
634: */
635: ff[0] = PetscSinScalar(3.0*xx[0]) + xx[0];
636: ff[1] = xx[1];
638: /*
639: Restore vectors
640: */
641: VecRestoreArray(x,&xx);
642: VecRestoreArray(f,&ff);
643: return(0);
644: }
645: /* ------------------------------------------------------------------- */
648: static PetscErrorCode FormJacobian2_block(SNES snes,Vec x,Mat *jac,Mat *B,MatStructure *flag,void *dummy)
649: {
650: PetscScalar *xx,A[4];
652: PetscInt idx[2] = {0,1};
655: /*
656: Get pointer to vector data
657: */
658: VecGetArray(x,&xx);
660: /*
661: Compute Jacobian entries and insert into matrix.
662: - Since this is such a small problem, we set all entries for
663: the matrix at once.
664: */
665: A[0] = 3.0*PetscCosScalar(3.0*xx[0]) + 1.0; A[1] = 0.0;
666: A[2] = 0.0; A[3] = 1.0;
667: MatSetValues(*jac,2,idx,2,idx,A,INSERT_VALUES);
668: *flag = SAME_NONZERO_PATTERN;
670: /*
671: Restore vector
672: */
673: VecRestoreArray(x,&xx);
675: /*
676: Assemble matrix
677: */
678: MatAssemblyBegin(*jac,MAT_FINAL_ASSEMBLY);
679: MatAssemblyEnd(*jac,MAT_FINAL_ASSEMBLY);
680: return(0);
681: }
685: int main(int argc,char **argv)
686: {
687: PetscMPIInt size;
690: PetscInitialize(&argc,&argv,(char *)0,help);
692: MPI_Comm_size(PETSC_COMM_WORLD,&size);
693: if (size != 1) SETERRQ(PETSC_COMM_WORLD, 1,"This is a uniprocessor example only!");
695: assembled_system();
697: block_system();
699: PetscFinalize();
700: return 0;
701: }