MPI-AMRVAC  3.1
The MPI - Adaptive Mesh Refinement - Versatile Advection Code (development version)
mod_kdv.t
Go to the documentation of this file.
1 !> Module for including kdv source term in simulations
2 !> adds \f$-\delta^2*\sum_i \partial_{iii} \rho \f$ over dimensions i
3 module mod_kdv
4  implicit none
5 
6  !> source split or not
7  logical :: kdv_split= .false.
8  !> forefactor \f$ \delta^2\f$ of \f$ \partial_{iii} \f$ term
9  double precision :: kdv_delta = 1.0d0
10  !> switch for second order [1] or fourth order [2] central FD for \f$ \partial_{iii}\f$
11  !> Note: fourth order needs 3 nghostcells, all assume equidistant grid
12  integer :: kdv_order = 1
13 
14 contains
15  !> Read this module's parameters from a file
16  subroutine kdv_params_read(files)
18  character(len=*), intent(in) :: files(:)
19  integer :: n
20 
21  namelist /kdv_list/ kdv_split, kdv_delta, kdv_order
22 
23  do n = 1, size(files)
24  open(unitpar, file=trim(files(n)), status="old")
25  read(unitpar, kdv_list, end=111)
26 111 close(unitpar)
27  end do
28 
29  end subroutine kdv_params_read
30 
31  !> Initialize the module
32  subroutine kdv_init()
34 
36 
37  end subroutine kdv_init
38 
39  !> w[iw]=w[iw]+qdt*S[wCT,qtC,x] where S is the source based on wCT within ixO
40  subroutine kdv_add_source(qdt,ixI^L,ixO^L,wCT,w,x,qsourcesplit,active)
42  use mod_usr_methods
43 
44  integer, intent(in) :: ixI^L, ixO^L
45  double precision, intent(in) :: qdt, x(ixI^S,1:ndim)
46  double precision, intent(in) :: wCT(ixI^S,1:nw)
47  double precision, intent(inout) :: w(ixI^S,1:nw)
48  logical, intent(in) :: qsourcesplit
49  logical, intent(inout) :: active
50 
51  integer :: idir, lx^L,kx^L,jx^L,hx^L,gx^L,fx^L
52  double precision :: skdv(ixI^S)
53 
54  if(qsourcesplit .eqv. kdv_split) then
55  active = .true.
56  skdv(ixo^s)=zero
57  select case(kdv_order)
58  case(1)
59  do idir=1,ndim
60  ! The source is based on the time centered wCT
61  kx^l=ixo^l+2*kr(idir,^d);
62  jx^l=ixo^l+kr(idir,^d);
63  hx^l=ixo^l-kr(idir,^d);
64  gx^l=ixo^l-2*kr(idir,^d);
65  ! 2nd order centered difference for -\partial_xxx \rho
66  ! warning: needs 2 ghostcells, equidistant grid
67  skdv(ixo^s)=skdv(ixo^s)+(wct(kx^s,iw_rho)-2.0d0*wct(jx^s,iw_rho) &
68  +2.0d0*wct(hx^s,iw_rho)-wct(gx^s,iw_rho)) &
69  /(2.0d0 *dxlevel(idir)**3)
70  enddo
71  case(2)
72  do idir=1,ndim
73  ! The source is based on the time centered wCT
74  lx^l=ixo^l+3*kr(idir,^d);
75  kx^l=ixo^l+2*kr(idir,^d);
76  jx^l=ixo^l+kr(idir,^d);
77  hx^l=ixo^l-kr(idir,^d);
78  gx^l=ixo^l-2*kr(idir,^d);
79  fx^l=ixo^l-3*kr(idir,^d);
80  ! 4th order centered difference for -\partial_xxx \rho
81  ! warning: needs 3 ghostcells, equidistant grid
82  skdv(ixo^s)=skdv(ixo^s)+(-wct(lx^s,iw_rho)+8.0d0*wct(kx^s,iw_rho)-13.0d0*wct(jx^s,iw_rho) &
83  +13.0d0*wct(hx^s,iw_rho)-8.0d0*wct(gx^s,iw_rho)+wct(fx^s,iw_rho)) &
84  /(8.0d0 *dxlevel(idir)**3)
85  enddo
86  case default
87  call mpistop('undefined kdv_order parameter: see mod_kdv.t')
88  end select
89  w(ixo^s,iw_rho) = w(ixo^s,iw_rho) - qdt*kdv_delta**2*skdv(ixo^s)
90  end if
91 
92  end subroutine kdv_add_source
93 
94  subroutine kdv_get_dt(w,ixI^L,ixO^L,dtnew,dx^D,x)
96  use mod_usr_methods
97 
98  integer, intent(in) :: ixI^L, ixO^L
99  double precision, intent(in) :: dx^D, x(ixI^S,1:ndim), w(ixI^S,1:nw)
100  double precision, intent(inout) :: dtnew
101 
102  double precision :: dxarr(ndim), max_sinefactor
103 
104  !> Time step constraint for leap-frog time stepping combined with central 2nd order FD
105  !> see e.g. Chun-Te Lee et al, Journal of Mathematics Research vol. 9, no.4, 2017
106  !> ISSN 1916-9795
107  ^d&dxarr(^d)=dx^d;
108  max_sinefactor=3.0d0*dsqrt(3.0d0)/2.0d0
109  dtnew=dtdiffpar*minval(dxarr(1:ndim))**3/(max_sinefactor*kdv_delta**2)
110 
111  end subroutine kdv_get_dt
112 
113 end module mod_kdv
This module contains definitions of global parameters and variables and some generic functions/subrou...
double precision dtdiffpar
For resistive MHD, the time step is also limited by the diffusion time: .
integer, parameter unitpar
file handle for IO
integer, dimension(3, 3) kr
Kronecker delta tensor.
character(len=std_len), dimension(:), allocatable par_files
Which par files are used as input.
integer, dimension(:), allocatable, parameter d
double precision, dimension(ndim) dxlevel
Module for including kdv source term in simulations adds over dimensions i.
Definition: mod_kdv.t:3
subroutine kdv_get_dt(w, ixIL, ixOL, dtnew, dxD, x)
Definition: mod_kdv.t:95
subroutine kdv_init()
Initialize the module.
Definition: mod_kdv.t:33
integer kdv_order
switch for second order [1] or fourth order [2] central FD for Note: fourth order needs 3 nghostcell...
Definition: mod_kdv.t:12
subroutine kdv_params_read(files)
Read this module's parameters from a file.
Definition: mod_kdv.t:17
double precision kdv_delta
forefactor of term
Definition: mod_kdv.t:9
logical kdv_split
source split or not
Definition: mod_kdv.t:7
subroutine kdv_add_source(qdt, ixIL, ixOL, wCT, w, x, qsourcesplit, active)
w[iw]=w[iw]+qdt*S[wCT,qtC,x] where S is the source based on wCT within ixO
Definition: mod_kdv.t:41
Module with all the methods that users can customize in AMRVAC.