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