MPI-AMRVAC 3.1
The MPI - Adaptive Mesh Refinement - Versatile Advection Code (development version)
Loading...
Searching...
No Matches
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
3module 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
14contains
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)
26111 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)
43 use mod_comm_lib, only: mpistop
44
45 integer, intent(in) :: ixI^L, ixO^L
46 double precision, intent(in) :: qdt, x(ixI^S,1:ndim)
47 double precision, intent(in) :: wCT(ixI^S,1:nw)
48 double precision, intent(inout) :: w(ixI^S,1:nw)
49 logical, intent(in) :: qsourcesplit
50 logical, intent(inout) :: active
51
52 integer :: idir, lx^L,kx^L,jx^L,hx^L,gx^L,fx^L
53 double precision :: skdv(ixI^S)
54
55 if(qsourcesplit .eqv. kdv_split) then
56 active = .true.
57 skdv(ixo^s)=zero
58 select case(kdv_order)
59 case(1)
60 do idir=1,ndim
61 ! The source is based on the time centered wCT
62 kx^l=ixo^l+2*kr(idir,^d);
63 jx^l=ixo^l+kr(idir,^d);
64 hx^l=ixo^l-kr(idir,^d);
65 gx^l=ixo^l-2*kr(idir,^d);
66 ! 2nd order centered difference for -\partial_xxx \rho
67 ! warning: needs 2 ghostcells, equidistant grid
68 skdv(ixo^s)=skdv(ixo^s)+(wct(kx^s,iw_rho)-2.0d0*wct(jx^s,iw_rho) &
69 +2.0d0*wct(hx^s,iw_rho)-wct(gx^s,iw_rho)) &
70 /(2.0d0 *dxlevel(idir)**3)
71 enddo
72 case(2)
73 do idir=1,ndim
74 ! The source is based on the time centered wCT
75 lx^l=ixo^l+3*kr(idir,^d);
76 kx^l=ixo^l+2*kr(idir,^d);
77 jx^l=ixo^l+kr(idir,^d);
78 hx^l=ixo^l-kr(idir,^d);
79 gx^l=ixo^l-2*kr(idir,^d);
80 fx^l=ixo^l-3*kr(idir,^d);
81 ! 4th order centered difference for -\partial_xxx \rho
82 ! warning: needs 3 ghostcells, equidistant grid
83 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) &
84 +13.0d0*wct(hx^s,iw_rho)-8.0d0*wct(gx^s,iw_rho)+wct(fx^s,iw_rho)) &
85 /(8.0d0 *dxlevel(idir)**3)
86 enddo
87 case default
88 call mpistop('undefined kdv_order parameter: see mod_kdv.t')
89 end select
90 w(ixo^s,iw_rho) = w(ixo^s,iw_rho) - qdt*kdv_delta**2*skdv(ixo^s)
91 end if
92
93 end subroutine kdv_add_source
94
95 subroutine kdv_get_dt(w,ixI^L,ixO^L,dtnew,dx^D,x)
98
99 integer, intent(in) :: ixI^L, ixO^L
100 double precision, intent(in) :: dx^D, x(ixI^S,1:ndim), w(ixI^S,1:nw)
101 double precision, intent(inout) :: dtnew
102
103 double precision :: dxarr(ndim), max_sinefactor
104
105 !> Time step constraint for leap-frog time stepping combined with central 2nd order FD
106 !> see e.g. Chun-Te Lee et al, Journal of Mathematics Research vol. 9, no.4, 2017
107 !> ISSN 1916-9795
108 ^d&dxarr(^d)=dx^d;
109 max_sinefactor=3.0d0*dsqrt(3.0d0)/2.0d0
110 dtnew=dtdiffpar*minval(dxarr(1:ndim))**3/(max_sinefactor*kdv_delta**2)
111
112 end subroutine kdv_get_dt
113
114end module mod_kdv
subroutine, public mpistop(message)
Exit MPI-AMRVAC with an error message.
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.
double precision, dimension(:), allocatable, parameter d
double precision, dimension(^nd) dxlevel
store unstretched cell size of current level
Module for including kdv source term in simulations adds over dimensions i.
Definition mod_kdv.t:3
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
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
subroutine kdv_get_dt(w, ixil, ixol, dtnew, dxd, x)
Definition mod_kdv.t:96
double precision kdv_delta
forefactor of term
Definition mod_kdv.t:9
logical kdv_split
source split or not
Definition mod_kdv.t:7
Module with all the methods that users can customize in AMRVAC.