Actual source code: ex26.c
1: static char help[] = "Solves the trival ODE 2 du/dt = 1, u(0) = 0. \n\n";
3: #include <petscts.h>
4: #include <petscpc.h>
6: PetscErrorCode IFunction(TS,PetscReal,Vec,Vec,Vec,void*);
7: PetscErrorCode IJacobian(TS,PetscReal,Vec,Vec,PetscReal,Mat,Mat,void*);
9: int main(int argc,char **argv)
10: {
11: TS ts;
12: Vec x;
13: Vec f;
14: Mat A;
15: PetscErrorCode ierr;
17: PetscInitialize(&argc,&argv,(char*)0,help);
19: TSCreate(PETSC_COMM_WORLD,&ts);
20: TSSetEquationType(ts,TS_EQ_ODE_IMPLICIT);
21: VecCreate(PETSC_COMM_WORLD,&f);
22: VecSetSizes(f,1,PETSC_DECIDE);
23: VecSetFromOptions(f);
24: VecSetUp(f);
25: TSSetIFunction(ts,f,IFunction,NULL);
26: VecDestroy(&f);
28: MatCreate(PETSC_COMM_WORLD,&A);
29: MatSetSizes(A,1,1,PETSC_DECIDE,PETSC_DECIDE);
30: MatSetFromOptions(A);
31: MatSetUp(A);
32: MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);
33: MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);
34: /* ensure that the Jacobian matrix has diagonal entries since that is required by TS */
35: MatShift(A,(PetscReal)1);
36: MatShift(A,(PetscReal)-1);
37: TSSetIJacobian(ts,A,A,IJacobian,NULL);
38: MatDestroy(&A);
40: VecCreate(PETSC_COMM_WORLD,&x);
41: VecSetSizes(x,1,PETSC_DECIDE);
42: VecSetFromOptions(x);
43: VecSetUp(x);
44: TSSetSolution(ts,x);
45: VecDestroy(&x);
46: TSSetFromOptions(ts);
48: TSSetStepNumber(ts,0);
49: TSSetTimeStep(ts,1);
50: TSSetTime(ts,0);
51: TSSetMaxTime(ts,PETSC_MAX_REAL);
52: TSSetMaxSteps(ts,3);
54: /*
55: When an ARKIMEX scheme with an explicit stage is used this will error with a message informing the user it is not possible to use
56: a non-trivial mass matrix with ARKIMEX schemes with explicit stages.
57: */
58: TSSolve(ts,NULL);
59: if (ierr != PETSC_ERR_ARG_INCOMP) ierr;
61: TSDestroy(&ts);
62: PetscFinalize();
63: return 0;
64: }
66: PetscErrorCode IFunction(TS ts,PetscReal t,Vec x,Vec xdot,Vec f,void *ctx)
67: {
68: VecCopy(xdot,f);
69: VecScale(f,2.0);
70: VecShift(f,-1.0);
71: return 0;
72: }
74: PetscErrorCode IJacobian(TS ts,PetscReal t,Vec x,Vec xdot,PetscReal shift,Mat A,Mat B,void *ctx)
75: {
76: PetscScalar j;
78: j = shift*2.0;
79: MatSetValue(B,0,0,j,INSERT_VALUES);
80: MatAssemblyBegin(B,MAT_FINAL_ASSEMBLY);
81: MatAssemblyEnd(B,MAT_FINAL_ASSEMBLY);
82: return 0;
83: }
85: /*TEST
87: test:
88: suffix: arkimex_explicit_stage
89: requires: defined(PETSC_USE_DEBUG)
90: args: -ts_type arkimex -error_output_stdout
91: filter: egrep -v "(Petsc|on a| at |Configure)"
93: test:
94: suffix: arkimex_implicit_stage
95: args: -ts_type arkimex -ts_arkimex_type l2 -ts_monitor_solution -ts_monitor
97: TEST*/