Actual source code: ex1f90.F90

      1: ! test phase space (Maxwellian) mesh construction (serial)
      2: !
      3: ! Contributed by Mark Adams
      4: #include <petsc/finclude/petscts.h>
      5: #include <petsc/finclude/petscdmplex.h>
      6: program DMPlexTestLandauInterface
      7:   use petscts
      8:   use petscdmplex
      9:   implicit none

     11:   external DMPlexLandauIFunction
     12:   external DMPlexLandauIJacobian
     13:   DM dm
     14:   PetscInt, parameter :: dim = 2
     15:   PetscErrorCode ierr
     16:   TS ts
     17:   Vec X, X_0
     18:   Mat J
     19:   SNES snes
     20:   KSP ksp
     21:   PC pc
     22:   SNESLineSearch linesearch
     23:   PetscReal :: time
     24:   PetscScalar, parameter :: scalar = -1.0

     26:   PetscCallA(PetscInitialize(ierr))

     28:   ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
     29:   !  Create mesh (DM), read in parameters, create and add f_0 (X)
     30:   ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
     31:   PetscCallA(DMPlexLandauCreateVelocitySpace(PETSC_COMM_WORLD, dim, '', X, J, dm, ierr))
     32:   PetscCallA(DMSetUp(dm, ierr))
     33:   PetscCallA(VecDuplicate(X, X_0, ierr))
     34:   PetscCallA(VecCopy(X, X_0, ierr))
     35:   ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
     36:   !  View
     37:   ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
     38:   PetscCallA(DMPlexLandauPrintNorms(X, 0_PETSC_INT_KIND, ierr))
     39:   PetscCallA(DMSetOutputSequenceNumber(dm, 0_PETSC_INT_KIND, 0.0_PETSC_REAL_KIND, ierr))
     40:   ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
     41:   !    Create timestepping solver context
     42:   ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
     43:   PetscCallA(TSCreate(PETSC_COMM_SELF, ts, ierr))
     44:   PetscCallA(TSSetOptionsPrefix(ts, 'ex1_', ierr)) ! should get this from the dm or give it to the dm
     45:   PetscCallA(TSSetDM(ts, dm, ierr))
     46:   PetscCallA(TSGetSNES(ts, snes, ierr))
     47:   PetscCallA(SNESSetOptionsPrefix(snes, 'ex1_', ierr)) ! should get this from the dm or give it to the dm
     48:   PetscCallA(SNESGetLineSearch(snes, linesearch, ierr))
     49:   PetscCallA(SNESLineSearchSetType(linesearch, SNESLINESEARCHBASIC, ierr))
     50:   PetscCallA(TSSetIFunction(ts, PETSC_NULL_VEC, DMPlexLandauIFunction, PETSC_NULL_VEC, ierr))
     51:   PetscCallA(TSSetIJacobian(ts, J, J, DMPlexLandauIJacobian, PETSC_NULL_VEC, ierr))
     52:   PetscCallA(TSSetExactFinalTime(ts, TS_EXACTFINALTIME_STEPOVER, ierr))

     54:   PetscCallA(SNESGetKSP(snes, ksp, ierr))
     55:   PetscCallA(KSPSetOptionsPrefix(ksp, 'ex1_', ierr)) ! should get this from the dm or give it to the dm
     56:   PetscCallA(KSPGetPC(ksp, pc, ierr))
     57:   PetscCallA(PCSetOptionsPrefix(pc, 'ex1_', ierr)) ! should get this from the dm or give it to the dm

     59:   PetscCallA(TSSetFromOptions(ts, ierr))
     60:   PetscCallA(TSSetSolution(ts, X, ierr))
     61:   ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
     62:   !  Solve nonlinear system
     63:   ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
     64:   PetscCallA(TSSolve(ts, X, ierr))
     65:   PetscCallA(DMPlexLandauPrintNorms(X, 1_PETSC_INT_KIND, ierr))
     66:   PetscCallA(TSGetTime(ts, time, ierr))
     67:   PetscCallA(DMSetOutputSequenceNumber(dm, 1_PETSC_INT_KIND, time, ierr))
     68:   ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
     69:   !  remove f_0
     70:   ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
     71:   PetscCallA(VecAXPY(X, scalar, X_0, ierr))
     72:   PetscCallA(DMPlexLandauDestroyVelocitySpace(dm, ierr))
     73:   PetscCallA(TSDestroy(ts, ierr))
     74:   PetscCallA(VecDestroy(X, ierr))
     75:   PetscCallA(VecDestroy(X_0, ierr))
     76:   PetscCallA(PetscFinalize(ierr))
     77: end program DMPlexTestLandauInterface

     79: !/*TEST
     80: !
     81: !  build:
     82: !    requires: defined(PETSC_USE_DMLANDAU_2D)
     83: !
     84: !  test:
     85: !    suffix: 0
     86: !    requires: p4est !complex !kokkos_kernels !cuda
     87: !    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_time_step 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
     88: !
     89: !TEST*/