Lines Matching +full:test +full:- +full:basic
8 "y_dot = -0.2*x\n"
10 … "7 polynomials (dir=+1) with zeros: 1+2^i, i=-3,...3, on ranks=(i+3)%size\n"
11 … "7 polynomials (dir=-1) with zeros: 1+(8-2^i), i=-3,...3, on ranks=(i+3)%size\n"
12 "(t-5)^2 * sin(pi*t), with zeros = 1,2,...10, on rank-0\n"
13 " 0.5 * cos(pi*t), with zeros = 0.5,1.5,...9.5, on last-rank\n"
15 "-dir d : zero-crossing direction for events\n"
16 "-flg : additional output in Postevent\n"
17 … "-errtol e : error tolerance, for printing 'pass/fail' for located events (1e-5 by default)\n"
18 "-restart : flag for TSRestartStep() in PostEvent\n"
19 … "-dtpost x : if x > 0, then on even PostEvent calls 1st-post-event-step = x is set,\n"
20 … " on odd PostEvent calls 1st-post-event-step = PETSC_DECIDE is set,\n"
30 PetscReal evres[MAX_NEV]; // times of found zero-crossings
31 PetscReal ref[MAX_NEV]; // reference times of zero-crossings, for checking
35 …errtol; // error tolerance, for printing 'pass/fail' for located events (1e-5 by default)
37 PetscReal dtpost; // post-event step
40 PetscInt dir0; // desired zero-crossing direction
64 ctx.pi = PetscAcosReal(-1.0); in main()
68 ctx.errtol = 1e-5; in main()
79 vals[2] = -0.2; in main()
110 …PetscCall(PetscOptionsGetInt(NULL, NULL, "-dir", &ctx.dir0, NULL)); // desired zero-crossi… in main()
111 …PetscCall(PetscOptionsHasName(NULL, NULL, "-flg", &ctx.flg)); // flag for additional… in main()
112 …PetscCall(PetscOptionsGetReal(NULL, NULL, "-errtol", &ctx.errtol, NULL)); // error tolerance for… in main()
113 …PetscCall(PetscOptionsGetBool(NULL, NULL, "-restart", &ctx.restart, NULL)); // flag for TSRestartS… in main()
114 PetscCall(PetscOptionsGetReal(NULL, NULL, "-dtpost", &ctx.dtpost, NULL)); // post-event step in main()
118 for (PetscInt i = -3; i <= 3; i++) { // pos-polynomials in main()
127 for (PetscInt i = -3; i <= 3; i++) { // neg-polynomials in main()
131 if (ctx.dir0 <= 0) ctx.ref[ctx.cntref++] = 9.0 - aux; in main()
135 if (ctx.rank == 0) { // sin-event -- on rank-0 in main()
138 for (PetscInt i = 1; i < MAX_NEV / 2 - 10; i++) { in main()
143 if (ctx.rank == ctx.size - 1) { // cos-event -- on last rank in main()
146 for (PetscInt i = 1; i < MAX_NEV / 2 - 10; i++) { in main()
147 if (i % 2 == 1 && ctx.dir0 <= 0) ctx.ref[ctx.cntref++] = i - 0.5; in main()
148 if (i % 2 == 0 && ctx.dir0 >= 0) ctx.ref[ctx.cntref++] = i - 0.5; in main()
153 SetVtols(ctx.rank, ctx.size, 1e-8, 1e-8, ctx.vtol); in main()
162 if (j < ctx.cntref) err = PetscAbsReal(ctx.evres[j] - ctx.ref[j]); in main()
176 User callback for defining the event-functions
185 // for the test purposes, event-functions are defined based on t in EventFunction()
186 for (PetscInt i = -3; i <= 3; i++) { // pos-polynomials in EventFunction()
187 if (Ctx->rank == (i + 3) % Ctx->size) { in EventFunction()
189 if (t < 2 + P) gval[n++] = 1 - PetscPowReal(2 + P - t, i + 5); in EventFunction()
193 for (PetscInt i = -3; i <= 3; i++) { // neg-polynomials in EventFunction()
194 if (Ctx->rank == (i + 3) % Ctx->size) { in EventFunction()
196 if (t > 8 - P) gval[n++] = 1 - PetscPowReal(t - 8 + P, i + 5); in EventFunction()
200 …if (Ctx->rank == 0) gval[n++] = (t - 5) * (t - 5) * PetscSinReal(Ctx->pi * t); // sin-event -- on … in EventFunction()
201 …if (Ctx->rank == Ctx->size - 1) gval[n++] = 0.5 * PetscCosReal(Ctx->pi * t); // cos-event -- on … in EventFunction()
206 User callback for the post-event stuff
213 if (Ctx->flg) { in Postevent()
214 PetscCallBack("EventFunction", EventFunction(ts, t, U, Ctx->fvals, ctx)); in Postevent()
215 …, "[%d] At t = %20.16g : %" PetscInt_FMT " events triggered, fvalues =", Ctx->rank, (double)t, nev… in Postevent()
216 …j++) PetscCall(PetscSynchronizedPrintf(PETSC_COMM_WORLD, "\t%g", (double)Ctx->fvals[evs_zero[j]])); in Postevent()
221 if (Ctx->cnt + nev_zero < MAX_NEV) in Postevent()
222 …for (PetscInt i = 0; i < nev_zero; i++) Ctx->evres[Ctx->cnt++] = t; // save the repeating zeros se… in Postevent()
225 Ctx->postcnt++; // sync in Postevent()
226 if (Ctx->dtpost > 0) { in Postevent()
227 if (Ctx->postcnt % 2 == 0) PetscCall(TSSetPostEventStep(ts, Ctx->dtpost)); in Postevent()
232 …if ((Ctx->dir0 == 0 && PetscAbsReal(t - (PetscReal)4.0) < 0.01) || (Ctx->dir0 == -1 && PetscAbsRea… in Postevent()
233 …SetVtols(Ctx->rank, Ctx->size, 1e-8, 1e-26, Ctx->vtol); // for better resolution of sin-event at t… in Postevent()
234 PetscCall(TSSetEventTolerances(ts, PETSC_DECIDE, Ctx->vtol)); in Postevent()
236 if (PetscAbsReal(t - (PetscReal)5.0) < 0.01) { in Postevent()
237 SetVtols(Ctx->rank, Ctx->size, 1e-8, 1e-8, Ctx->vtol); // back to normal in Postevent()
238 PetscCall(TSSetEventTolerances(ts, PETSC_DECIDE, Ctx->vtol)); in Postevent()
241 if (Ctx->restart) PetscCall(TSRestartStep(ts)); in Postevent()
249 for (PetscInt i = -3; i <= 3; i++) in SetVtols()
250 if (rank == (i + 3) % size) vtol[n++] = tol0; // pos-polynomials in SetVtols()
251 for (PetscInt i = -3; i <= 3; i++) in SetVtols()
252 if (rank == (i + 3) % size) vtol[n++] = tol0; // neg-polynomials in SetVtols()
253 if (rank == 0) vtol[n++] = tolsin; // sin-event -- on rank-0 in SetVtols()
254 if (rank == size - 1) vtol[n++] = tol0; // cos-event -- on last rank in SetVtols()
256 /*---------------------------------------------------------------------------------------------*/
258 Note, in the tests below, -ts_event_post_event_step is occasionally set to -1,
260 explicitly specify -1 in this option. Rather, if PETSC_DECIDE behaviour is needed,
264 /*TEST
265 test:
268 args: -dir 1 -ts_event_dt_min 1e-6
269 args: -restart 1
270 args: -dtpost {{0 0.25}}
271 args: -ts_event_post_event_step 0.31
272 args: -ts_type rk
273 args: -ts_adapt_type {{none basic}}
276 test:
279 args: -dir 1 -ts_event_dt_min 1e-6 -ts_time_step 0.25
280 args: -restart 0
281 args: -dtpost 0
282 args: -ts_event_post_event_step -1
283 args: -ts_type {{beuler rk}}
284 args: -ts_adapt_type {{none basic}}
289 test:
292 args: -dir 0 -ts_event_dt_min 1e-6
293 args: -restart 0
294 args: -dtpost {{0 0.25}}
295 args: -ts_event_post_event_step -1
296 args: -ts_type rk
297 args: -ts_adapt_type {{none basic}}
300 test:
303 args: -dir 0 -ts_event_dt_min 1e-6 -ts_time_step 0.25
304 args: -dtpost 0
305 args: -ts_event_post_event_step {{-1 0.29}}
306 args: -ts_event_post_event_second_step {{-1 0.31}}
307 args: -ts_type rk
308 args: -ts_adapt_type {{none basic}}
313 test:
316 args: -dir -1 -ts_event_dt_min 1e-6
317 args: -restart 1
318 args: -dtpost {{0 0.25}}
319 args: -ts_event_post_event_step 0.31
320 args: -ts_type beuler
321 args: -ts_adapt_type {{none basic}}
326 test:
329 args: -dir -1 -ts_event_dt_min 1e-6 -ts_time_step 0.25
330 args: -restart 0
331 args: -dtpost 0
332 args: -ts_event_post_event_step -1
333 args: -ts_type {{beuler rk}}
334 args: -ts_adapt_type {{none basic}}
339 TEST*/