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