Actual source code: ex1.c
1: static char help[] = "Landau collision operator driver\n\n";
3: #include <petscts.h>
4: #include <petsclandau.h>
6: int main(int argc, char **argv)
7: {
8: DM dm;
9: Vec X,X_0;
11: PetscInt dim=2;
12: TS ts;
13: Mat J;
14: SNES snes;
15: KSP ksp;
16: PC pc;
17: SNESLineSearch linesearch;
18: PetscReal time;
20: PetscInitialize(&argc, &argv, NULL,help);if (ierr) return ierr;
21: PetscOptionsGetInt(NULL,NULL, "-dim", &dim, NULL);
22: /* Create a mesh */
23: LandauCreateVelocitySpace(PETSC_COMM_SELF, dim, "", &X, &J, &dm);
24: LandauCreateMassMatrix(dm, NULL);
25: DMSetUp(dm);
26: VecDuplicate(X,&X_0);
27: VecCopy(X,X_0);
28: LandauPrintNorms(X,0);
29: DMSetOutputSequenceNumber(dm, 0, 0.0);
30: DMViewFromOptions(dm,NULL,"-dm_view");
31: VecViewFromOptions(X,NULL,"-vec_view");
32: /* Create timestepping solver context */
33: TSCreate(PETSC_COMM_SELF,&ts);
34: TSSetOptionsPrefix(ts, "ex1_"); /* should get this from the dm or give it to the dm */
35: TSSetDM(ts,dm);
36: TSGetSNES(ts,&snes);
37: SNESSetOptionsPrefix(snes, "ex1_"); /* should get this from the dm or give it to the dm */
38: SNESGetLineSearch(snes,&linesearch);
39: SNESLineSearchSetType(linesearch,SNESLINESEARCHBASIC);
40: TSSetIFunction(ts,NULL,LandauIFunction,NULL);
41: TSSetIJacobian(ts,J,J,LandauIJacobian,NULL);
42: TSSetExactFinalTime(ts,TS_EXACTFINALTIME_STEPOVER);
43: SNESGetKSP(snes,&ksp);
44: KSPSetOptionsPrefix(ksp, "ex1_"); /* should get this from the dm or give it to the dm */
45: KSPGetPC(ksp,&pc);
46: PCSetOptionsPrefix(pc, "ex1_"); /* should get this from the dm or give it to the dm */
47: TSSetFromOptions(ts);
48: TSSetSolution(ts,X);
49: TSSolve(ts,X);
50: LandauPrintNorms(X,1);
51: TSGetTime(ts, &time);
52: DMSetOutputSequenceNumber(dm, 1, time);
53: VecViewFromOptions(X,NULL,"-vec_view");
54: VecAXPY(X,-1,X_0);
55: /* clean up */
56: LandauDestroyVelocitySpace(&dm);
57: TSDestroy(&ts);
58: VecDestroy(&X);
59: VecDestroy(&X_0);
60: PetscFinalize();
61: return ierr;
62: }
64: /*TEST
66: test:
67: suffix: 0
68: requires: p4est !complex
69: args: -dm_landau_num_species_grid 1,2 -petscspace_degree 3 -petscspace_poly_tensor 1 -dm_landau_type p4est -dm_landau_ion_masses 2,4 -dm_landau_ion_charges 1,18 -dm_landau_thermal_temps 5,5,.5 -dm_landau_n 1.00018,1,1e-5 -dm_landau_n_0 1e20 -ex1_ts_monitor -ex1_snes_rtol 1.e-14 -ex1_snes_stol 1.e-14 -ex1_snes_monitor -ex1_snes_converged_reason -ex1_ts_type arkimex -ex1_ts_arkimex_type 1bee -ex1_ts_max_snes_failures -1 -ex1_ts_rtol 1e-1 -ex1_ts_dt 1.e-1 -ex1_ts_max_time 1 -ex1_ts_adapt_clip .5,1.25 -ex1_ts_adapt_scale_solve_failed 0.75 -ex1_ts_adapt_time_step_increase_delay 5 -ex1_ts_max_steps 1 -ex1_pc_type lu -ex1_ksp_type preonly -dm_landau_amr_levels_max 2,1
71: TEST*/