xref: /petsc/src/ts/utils/dmplexlandau/tutorials/ex1.c (revision e0eea49528ddcce32a775a862994a11c74fcac87)
1*e0eea495SMark static char help[] = "Landau collision operator driver\n\n";
2*e0eea495SMark 
3*e0eea495SMark #include <petscts.h>
4*e0eea495SMark #include <petsclandau.h>
5*e0eea495SMark 
6*e0eea495SMark int main(int argc, char **argv)
7*e0eea495SMark {
8*e0eea495SMark   DM             dm;
9*e0eea495SMark   Vec            X,X_0;
10*e0eea495SMark   PetscErrorCode ierr;
11*e0eea495SMark   PetscInt       dim=2;
12*e0eea495SMark   TS             ts;
13*e0eea495SMark   Mat            J;
14*e0eea495SMark   SNES           snes;
15*e0eea495SMark   KSP            ksp;
16*e0eea495SMark   PC             pc;
17*e0eea495SMark   SNESLineSearch linesearch;
18*e0eea495SMark   PetscReal      time;
19*e0eea495SMark 
20*e0eea495SMark   ierr = PetscInitialize(&argc, &argv, NULL,help);if (ierr) return ierr;
21*e0eea495SMark   ierr = PetscOptionsGetInt(NULL,NULL, "-dim", &dim, NULL);CHKERRQ(ierr);
22*e0eea495SMark   /* Create a mesh */
23*e0eea495SMark   ierr = LandauCreateVelocitySpace(PETSC_COMM_SELF, dim, "", &X, &J, &dm);CHKERRQ(ierr);
24*e0eea495SMark   ierr = LandauCreateMassMatrix(dm, NULL);CHKERRQ(ierr);
25*e0eea495SMark   ierr = DMSetUp(dm);CHKERRQ(ierr);
26*e0eea495SMark   ierr = VecDuplicate(X,&X_0);CHKERRQ(ierr);
27*e0eea495SMark   ierr = VecCopy(X,X_0);CHKERRQ(ierr);
28*e0eea495SMark   ierr = LandauPrintNorms(X,0);CHKERRQ(ierr);
29*e0eea495SMark   ierr = DMSetOutputSequenceNumber(dm, 0, 0.0);CHKERRQ(ierr);
30*e0eea495SMark   ierr = DMViewFromOptions(dm,NULL,"-dm_view");CHKERRQ(ierr);
31*e0eea495SMark   ierr = VecViewFromOptions(X,NULL,"-vec_view");CHKERRQ(ierr);
32*e0eea495SMark   /* Create timestepping solver context */
33*e0eea495SMark   ierr = TSCreate(PETSC_COMM_SELF,&ts);CHKERRQ(ierr);
34*e0eea495SMark   ierr = TSSetOptionsPrefix(ts, "ex1_");CHKERRQ(ierr);  /* should get this from the dm or give it to the dm */
35*e0eea495SMark   ierr = TSSetDM(ts,dm);CHKERRQ(ierr);
36*e0eea495SMark   ierr = TSGetSNES(ts,&snes);CHKERRQ(ierr);
37*e0eea495SMark   ierr = SNESSetOptionsPrefix(snes, "ex1_");CHKERRQ(ierr);  /* should get this from the dm or give it to the dm */
38*e0eea495SMark   ierr = SNESGetLineSearch(snes,&linesearch);CHKERRQ(ierr);
39*e0eea495SMark   ierr = SNESLineSearchSetType(linesearch,SNESLINESEARCHBASIC);CHKERRQ(ierr);
40*e0eea495SMark   ierr = TSSetIFunction(ts,NULL,LandauIFunction,NULL);CHKERRQ(ierr);
41*e0eea495SMark   ierr = TSSetIJacobian(ts,J,J,LandauIJacobian,NULL);CHKERRQ(ierr);
42*e0eea495SMark   ierr = TSSetExactFinalTime(ts,TS_EXACTFINALTIME_STEPOVER);CHKERRQ(ierr);
43*e0eea495SMark   ierr = SNESGetKSP(snes,&ksp);CHKERRQ(ierr);
44*e0eea495SMark   ierr = KSPSetOptionsPrefix(ksp, "ex1_");CHKERRQ(ierr);  /* should get this from the dm or give it to the dm */
45*e0eea495SMark   ierr = KSPGetPC(ksp,&pc);CHKERRQ(ierr);
46*e0eea495SMark   ierr = PCSetOptionsPrefix(pc, "ex1_");CHKERRQ(ierr);  /* should get this from the dm or give it to the dm */
47*e0eea495SMark   ierr = TSSetFromOptions(ts);CHKERRQ(ierr);
48*e0eea495SMark   ierr = TSSetSolution(ts,X);CHKERRQ(ierr);
49*e0eea495SMark   ierr = TSSolve(ts,X);CHKERRQ(ierr);
50*e0eea495SMark   ierr = LandauPrintNorms(X,1);CHKERRQ(ierr);
51*e0eea495SMark   ierr = TSGetTime(ts, &time);CHKERRQ(ierr);
52*e0eea495SMark   ierr = DMSetOutputSequenceNumber(dm, 1, time);CHKERRQ(ierr);
53*e0eea495SMark   ierr = VecViewFromOptions(X,NULL,"-vec_view");CHKERRQ(ierr);
54*e0eea495SMark   ierr = VecAXPY(X,-1,X_0);CHKERRQ(ierr);
55*e0eea495SMark   /* clean up */
56*e0eea495SMark   ierr = LandauDestroyVelocitySpace(&dm);CHKERRQ(ierr);
57*e0eea495SMark   ierr = TSDestroy(&ts);CHKERRQ(ierr);
58*e0eea495SMark   ierr = VecDestroy(&X);CHKERRQ(ierr);
59*e0eea495SMark   ierr = VecDestroy(&X_0);CHKERRQ(ierr);
60*e0eea495SMark   ierr = PetscFinalize();
61*e0eea495SMark   return ierr;
62*e0eea495SMark }
63*e0eea495SMark 
64*e0eea495SMark 
65*e0eea495SMark /*TEST
66*e0eea495SMark 
67*e0eea495SMark   test:
68*e0eea495SMark     suffix: 0
69*e0eea495SMark     requires: p4est !complex
70*e0eea495SMark     args: -petscspace_degree 4 -petscspace_poly_tensor 1 -dm_landau_type p4est -info :dm,tsadapt -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-6 -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-4 -ex1_ts_dt 1.e-6 -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 7 -dm_landau_domain_radius 5 -dm_landau_amr_re_levels 0 -dm_landau_re_radius 1 -dm_landau_amr_z_refine1 1 -dm_landau_amr_z_refine2 0 -dm_landau_amr_post_refine 0 -dm_landau_z_radius1 .1 -dm_landau_z_radius2 .1 -dm_refine 1
71*e0eea495SMark 
72*e0eea495SMark TEST*/
73