Actual source code: ex29.c
1: static char help[] = "Tests TS time span \n\n";
3: #include <petscts.h>
5: static PetscErrorCode RHSFunction(TS ts, PetscReal t, Vec X, Vec F, void *ctx)
6: {
7: PetscInt i, n;
8: const PetscScalar *xx;
9: PetscScalar *ff;
12: VecGetLocalSize(X, &n);
13: VecGetArrayRead(X, &xx);
14: VecGetArray(F, &ff);
15: if (n >= 1) ff[0] = 1;
16: for (i = 1; i < n; i++) ff[i] = (i + 1) * (xx[i - 1] + PetscPowReal(t, i)) / 2;
17: VecRestoreArrayRead(X, &xx);
18: VecRestoreArray(F, &ff);
19: return 0;
20: }
22: int main(int argc, char *argv[])
23: {
24: TS ts;
25: Vec X, *Xs;
26: PetscInt i, n, N = 9;
27: PetscReal tspan[8] = {16.0, 16.1, 16.2, 16.3, 16.4, 16.5, 16.6, 16.7};
28: const PetscReal *tspan2;
31: PetscInitialize(&argc, &argv, NULL, help);
32: TSCreate(PETSC_COMM_SELF, &ts);
33: TSSetType(ts, TSRK);
34: TSSetRHSFunction(ts, NULL, RHSFunction, NULL);
35: VecCreateSeq(PETSC_COMM_SELF, N, &X);
36: VecZeroEntries(X);
37: TSSetTimeStep(ts, 0.001);
38: TSSetTimeSpan(ts, 8, tspan);
39: TSSetExactFinalTime(ts, TS_EXACTFINALTIME_MATCHSTEP);
40: TSSetFromOptions(ts);
41: TSSolve(ts, X);
42: TSGetTimeSpanSolutions(ts, &n, &Xs);
43: TSGetTimeSpan(ts, &n, &tspan2);
44: PetscPrintf(PETSC_COMM_WORLD, "Time Span: ");
45: for (i = 0; i < n; i++) PetscPrintf(PETSC_COMM_WORLD, " %g", (double)tspan2[i]);
46: PetscPrintf(PETSC_COMM_WORLD, "\n");
47: TSDestroy(&ts);
48: VecDestroy(&X);
49: PetscFinalize();
50: return 0;
51: }
53: /*TEST
55: testset:
56: test:
57: suffix: 1
58: args: -ts_monitor
59: test:
60: suffix: 2
61: requires: !single
62: args: -ts_monitor -ts_adapt_type none
63: TEST*/