Actual source code: ex1f90.F90
1: ! test phase space (Maxwellian) mesh construction (serial)
2: !
3: ! Contributed by Mark Adams
4: program DMPlexTestLandauInterface
5: #include <petsc/finclude/petscts.h>
6: #include <petsc/finclude/petscdmplex.h>
7: use petscts
8: use petscdmplex
9: implicit none
11: external DMPlexLandauIFunction
12: external DMPlexLandauIJacobian
13: DM dm
14: PetscInt dim
15: PetscInt ii
16: PetscErrorCode ierr
17: TS ts
18: Vec X,X_0
19: Mat J
20: SNES snes
21: KSP ksp
22: PC pc
23: SNESLineSearch linesearch
24: PetscReal mone
25: PetscScalar scalar
27: PetscCallA(PetscInitialize(ierr))
29: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
30: ! Create mesh (DM), read in parameters, create and add f_0 (X)
31: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
32: dim = 2
33: PetscCallA(DMPlexLandauCreateVelocitySpace(PETSC_COMM_WORLD, dim, '', X, J, dm, ierr))
34: PetscCallA(DMSetUp(dm,ierr))
35: PetscCallA(VecDuplicate(X,X_0,ierr))
36: PetscCallA(VecCopy(X,X_0,ierr))
37: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
38: ! View
39: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
40: ii = 0
41: PetscCallA(DMPlexLandauPrintNorms(X,ii,ierr))
42: mone = 0;
43: PetscCallA(DMSetOutputSequenceNumber(dm, ii, mone, ierr))
44: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
45: ! Create timestepping solver context
46: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
47: PetscCallA(TSCreate(PETSC_COMM_SELF,ts,ierr))
48: PetscCallA(TSSetOptionsPrefix(ts, 'ex1_', ierr)) ! should get this from the dm or give it to the dm
49: PetscCallA(TSSetDM(ts,dm,ierr))
50: PetscCallA(TSGetSNES(ts,snes,ierr))
51: PetscCallA(SNESSetOptionsPrefix(snes, 'ex1_', ierr)) ! should get this from the dm or give it to the dm
52: PetscCallA(SNESGetLineSearch(snes,linesearch,ierr))
53: PetscCallA(SNESLineSearchSetType(linesearch,SNESLINESEARCHBASIC,ierr))
54: PetscCallA(TSSetIFunction(ts,PETSC_NULL_VEC,DMPlexLandauIFunction,PETSC_NULL_VEC,ierr))
55: PetscCallA(TSSetIJacobian(ts,J,J,DMPlexLandauIJacobian,PETSC_NULL_VEC,ierr))
56: PetscCallA(TSSetExactFinalTime(ts,TS_EXACTFINALTIME_STEPOVER,ierr))
58: PetscCallA(SNESGetKSP(snes,ksp,ierr))
59: PetscCallA(KSPSetOptionsPrefix(ksp, 'ex1_', ierr)) ! should get this from the dm or give it to the dm
60: PetscCallA(KSPGetPC(ksp,pc,ierr))
61: PetscCallA(PCSetOptionsPrefix(pc, 'ex1_', ierr)) ! should get this from the dm or give it to the dm
63: PetscCallA(TSSetFromOptions(ts,ierr))
64: PetscCallA(TSSetSolution(ts,X,ierr))
65: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
66: ! Solve nonlinear system
67: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
68: PetscCallA(TSSolve(ts,X,ierr))
69: ii = 1
70: PetscCallA(DMPlexLandauPrintNorms(X,ii,ierr))
71: PetscCallA(TSGetTime(ts, mone, ierr))
72: PetscCallA(DMSetOutputSequenceNumber(dm, ii, mone, ierr))
73: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
74: ! remove f_0
75: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
76: scalar = -1.
77: PetscCallA(VecAXPY(X,scalar,X_0,ierr))
78: PetscCallA(DMPlexLandauDestroyVelocitySpace(dm, ierr))
79: PetscCallA(TSDestroy(ts, ierr))
80: PetscCallA(VecDestroy(X, ierr))
81: PetscCallA(VecDestroy(X_0, ierr))
82: PetscCallA(PetscFinalize(ierr))
83: end program DMPlexTestLandauInterface
85: !/*TEST
86: ! build:
87: ! requires: defined(PETSC_USING_F90FREEFORM) defined(PETSC_USE_DMLANDAU_2D)
88: !
89: ! test:
90: ! suffix: 0
91: ! requires: p4est !complex !kokkos_kernels !cuda
92: ! 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 unlimited -ex1_ts_rtol 1e-1 -ex1_ts_dt 1.e-3 -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 -dm_landau_device_type cpu -dm_landau_verbose 4 -dm_landau_normalization_grid 1
93: !
94: !TEST*/