MPI-AMRVAC  2.0
The MPI - Adaptive Mesh Refinement - Versatile Advection Code
mod_nonlinear_phys.t
Go to the documentation of this file.
1 !> Module containing the physics routines for scalar nonlinear equation
3 
4  implicit none
5  private
6 
7  !> index of the single scalar unknown
8  integer, protected, public :: rho_ = 1
9 
10  !> switch between burgers (i.e. rho**2)
11  !> or nonconvex flux (i.e. rho**3)
12  integer, protected, public :: nonlinear_flux_type = 1
13 
14  !> whether the KdV source term is added
15  logical, protected, public :: kdv_source_term = .false.
16 
17  ! Public methods
18  public :: nonlinear_phys_init
19  public :: nonlinear_get_v
20 
21 contains
22 
23  !> Read this module's parameters from a file
24  subroutine nonlinear_params_read(files)
26  character(len=*), intent(in) :: files(:)
27  integer :: n
28 
29  namelist /nonlinear_list/ nonlinear_flux_type, kdv_source_term
30 
31  do n = 1, size(files)
32  open(unitpar, file=trim(files(n)), status='old')
33  read(unitpar, nonlinear_list, end=111)
34 111 close(unitpar)
35  end do
36 
37  end subroutine nonlinear_params_read
38 
39  !> Write this module's parameters to a snapshot
40  subroutine nonlinear_write_info(fh)
41  ! for nonlinear scalar equation, nothing to write
42  ! note: this is info only stored at end of dat files,
43  ! is never read/used for restarts, only expects
44  ! an integer (number of parameters) and
45  ! corresponding double values and character names
46  ! and is meant for use in the python tools
48  integer, intent(in) :: fh
49  integer, dimension(MPI_STATUS_SIZE) :: st
50  integer :: er
51 
52  ! Write zero parameters
53  call mpi_file_write(fh, 0, 1, mpi_integer, st, er)
54 
55  end subroutine nonlinear_write_info
56 
57  subroutine nonlinear_phys_init()
59  use mod_physics
60  use mod_kdv, only: kdv_init
61 
62 
63  call nonlinear_params_read(par_files)
64 
65  physics_type = "nonlinear"
66  phys_energy = .false.
67  ! Whether diagonal ghost cells are required for the physics
68  phys_req_diagonal = .false.
69 
70  rho_ = var_set_rho()
71 
72  ! Check whether custom flux types have been defined
73  if (.not. allocated(flux_type)) then
74  allocate(flux_type(ndir, nw))
76  else if (any(shape(flux_type) /= [ndir, nw])) then
77  call mpistop("phys_check error: flux_type has wrong shape")
78  end if
79 
89 
90  if (kdv_source_term) call kdv_init()
91 
92  end subroutine nonlinear_phys_init
93 
94  subroutine nonlinear_to_conserved(ixI^L, ixO^L, w, x)
96  integer, intent(in) :: ixI^L, ixO^L
97  double precision, intent(inout) :: w(ixi^s, nw)
98  double precision, intent(in) :: x(ixi^s, 1:^nd)
99 
100  ! Do nothing (primitive and conservative are equal for nonlinear module)
101  end subroutine nonlinear_to_conserved
102 
103  subroutine nonlinear_to_primitive(ixI^L, ixO^L, w, x)
105  integer, intent(in) :: ixI^L, ixO^L
106  double precision, intent(inout) :: w(ixi^s, nw)
107  double precision, intent(in) :: x(ixi^s, 1:^nd)
108 
109  ! Do nothing (primitive and conservative are equal for nonlinear module)
110  end subroutine nonlinear_to_primitive
111 
112  subroutine nonlinear_get_v(w, x, ixI^L, ixO^L, idim, v)
114  integer, intent(in) :: ixI^L, ixO^L, idim
115  double precision, intent(in) :: w(ixi^s, nw), x(ixi^s, 1:^nd)
116  double precision, intent(out) :: v(ixi^s)
117 
118  select case(nonlinear_flux_type)
119  case(1)
120  v(ixo^s)=w(ixo^s,rho_)
121  case(2)
122  v(ixo^s)=3.0d0*w(ixo^s,rho_)**2
123  case default
124  call mpistop('Undefined fluxtype: set nonlinear_flux_type to 1 or 2')
125  end select
126 
127  end subroutine nonlinear_get_v
128 
129  subroutine nonlinear_get_cmax(w, x, ixI^L, ixO^L, idim, cmax)
131  integer, intent(in) :: ixI^L, ixO^L, idim
132  double precision, intent(in) :: w(ixi^s, nw), x(ixi^s, 1:^nd)
133  double precision, intent(inout) :: cmax(ixi^s)
134 
135  call nonlinear_get_v(w, x, ixi^l, ixo^l, idim, cmax)
136 
137  cmax(ixo^s) = abs(cmax(ixo^s))
138 
139  end subroutine nonlinear_get_cmax
140 
141  subroutine nonlinear_get_cbounds(wLC, wRC, wLp, wRp, x, ixI^L, ixO^L, idim, cmax, cmin)
143  integer, intent(in) :: ixI^L, ixO^L, idim
144  double precision, intent(in) :: wLC(ixi^s, nw), wRC(ixi^s,nw)
145  double precision, intent(in) :: wLp(ixi^s, nw), wRp(ixi^s,nw)
146  double precision, intent(in) :: x(ixi^s, 1:^nd)
147  double precision, intent(inout) :: cmax(ixi^s)
148  double precision, intent(inout), optional :: cmin(ixi^s)
149 
150  double precision :: wmean(ixi^s,nw)
151 
152  ! since get_v depends on w, the first argument should be some average over the
153  ! left and right state
154  wmean(ixo^s,1:nwflux)=0.5d0*(wlc(ixo^s,1:nwflux)+wrc(ixo^s,1:nwflux))
155  call nonlinear_get_v(wmean, x, ixi^l, ixo^l, idim, cmax)
156 
157  if (present(cmin)) then
158  cmin(ixo^s) = min(cmax(ixo^s), zero)
159  cmax(ixo^s) = max(cmax(ixo^s), zero)
160  else
161  cmax(ixo^s) = maxval(abs(cmax(ixo^s)))
162  end if
163 
164  end subroutine nonlinear_get_cbounds
165 
166  subroutine nonlinear_get_dt(w, ixI^L, ixO^L, dtnew, dx^D, x)
168  use mod_kdv, only: kdv_get_dt
169 
170  integer, intent(in) :: ixI^L, ixO^L
171  double precision, intent(in) :: dx^D, x(ixi^s, 1:^nd)
172  double precision, intent(in) :: w(ixi^s, 1:nw)
173  double precision, intent(inout) :: dtnew
174 
175  dtnew = bigdouble
176 
177  if(kdv_source_term) then
178  call kdv_get_dt(w, ixi^l, ixo^l, dtnew, dx^d, x)
179  endif
180  end subroutine nonlinear_get_dt
181 
182  ! here we select the flux according to the nonlinear_flux_type parameter
183  subroutine nonlinear_get_flux(wC, w, x, ixI^L, ixO^L, idim, f)
185  integer, intent(in) :: ixI^L, ixO^L, idim
186  double precision, intent(in) :: wC(ixi^s, 1:nw)
187  double precision, intent(in) :: w(ixi^s, 1:nw)
188  double precision, intent(in) :: x(ixi^s, 1:^nd)
189  double precision, intent(out) :: f(ixi^s, nwflux)
190 
191  select case(nonlinear_flux_type)
192  case(1)
193  f(ixo^s,rho_)=half*w(ixo^s,rho_)**2
194  case(2)
195  f(ixo^s,rho_)=w(ixo^s,rho_)**3
196  case default
197  call mpistop('Undefined fluxtype: set nonlinear_flux_type to 1 or 2')
198  end select
199 
200  end subroutine nonlinear_get_flux
201 
202  subroutine nonlinear_add_source_geom(qdt, ixI^L, ixO^L, wCT, w, x)
204  ! Add geometrical source terms to w
205  ! There are no geometrical source terms
206 
208 
209  integer, intent(in) :: ixI^L, ixO^L
210  double precision, intent(in) :: qdt, x(ixi^s, 1:^nd)
211  double precision, intent(inout) :: wCT(ixi^s, 1:nw), w(ixi^s, 1:nw)
212 
213  end subroutine nonlinear_add_source_geom
214 
215  subroutine nonlinear_add_source(qdt,ixI^L,ixO^L,wCT,w,x,qsourcesplit,active)
217  use mod_kdv, only: kdv_add_source
218 
219  integer, intent(in) :: ixI^L, ixO^L
220  double precision, intent(in) :: qdt
221  double precision, intent(in) :: wCT(ixi^s, 1:nw), x(ixi^s, 1:ndim)
222  double precision, intent(inout) :: w(ixi^s, 1:nw)
223  logical, intent(in) :: qsourcesplit
224  logical, intent(inout) :: active
225 
226  if(kdv_source_term) then
227  call kdv_add_source(qdt,ixi^l,ixo^l,wct,w,x,qsourcesplit,active)
228  end if
229 
230  end subroutine nonlinear_add_source
231 
232 end module mod_nonlinear_phys
This module contains definitions of global parameters and variables and some generic functions/subrou...
integer, public, protected rho_
index of the single scalar unknown
character(len=name_len) physics_type
String describing the physics type of the simulation.
Definition: mod_physics.t:14
subroutine, public nonlinear_get_v(w, x, ixIL, ixOL, idim, v)
integer ndir
Number of spatial dimensions (components) for vector variables.
logical, public, protected kdv_source_term
whether the KdV source term is added
subroutine nonlinear_get_cbounds(wLC, wRC, wLp, wRp, x, ixIL, ixOL, idim, cmax, cmin)
subroutine kdv_get_dt(w, ixIL, ixOL, dtnew, dxD, x)
Definition: mod_kdv.t:95
logical phys_req_diagonal
Whether the physics routines require diagonal ghost cells, for example for computing a curl...
Definition: mod_physics.t:22
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
integer, public, protected nonlinear_flux_type
switch between burgers (i.e. rho**2) or nonconvex flux (i.e. rho**3)
subroutine nonlinear_add_source(qdt, ixIL, ixOL, wCT, w, x, qsourcesplit, active)
subroutine nonlinear_get_dt(w, ixIL, ixOL, dtnew, dxD, x)
subroutine nonlinear_to_conserved(ixIL, ixOL, w, x)
subroutine kdv_init()
Initialize the module.
Definition: mod_kdv.t:33
subroutine nonlinear_add_source_geom(qdt, ixIL, ixOL, wCT, w, x)
procedure(sub_convert), pointer phys_to_conserved
Definition: mod_physics.t:36
integer, parameter unitpar
integer, parameter flux_default
Indicates a normal flux.
Definition: mod_physics.t:29
procedure(sub_get_flux), pointer phys_get_flux
Definition: mod_physics.t:41
character(len=std_len), dimension(:), allocatable par_files
Which par files are used as input.
procedure(sub_get_dt), pointer phys_get_dt
Definition: mod_physics.t:42
procedure(sub_add_source_geom), pointer phys_add_source_geom
Definition: mod_physics.t:44
subroutine mpistop(message)
Exit MPI-AMRVAC with an error message.
Definition: comm_lib.t:160
subroutine, public nonlinear_phys_init()
subroutine nonlinear_to_primitive(ixIL, ixOL, w, x)
integer, parameter ndim
Number of spatial dimensions for grid variables.
subroutine nonlinear_write_info(fh)
Write this module's parameters to a snapshot.
Module for including kdv source term in simulations adds over dimensions i.
Definition: mod_kdv.t:3
procedure(sub_get_cmax), pointer phys_get_cmax
Definition: mod_physics.t:39
This module defines the procedures of a physics module. It contains function pointers for the various...
Definition: mod_physics.t:4
subroutine nonlinear_get_flux(wC, w, x, ixIL, ixOL, idim, f)
procedure(sub_add_source), pointer phys_add_source
Definition: mod_physics.t:45
subroutine nonlinear_get_cmax(w, x, ixIL, ixOL, idim, cmax)
procedure(sub_convert), pointer phys_to_primitive
Definition: mod_physics.t:37
Module containing the physics routines for scalar nonlinear equation.
logical phys_energy
Solve energy equation or not.
integer, dimension(:, :), allocatable flux_type
Array per direction per variable, which can be used to specify that certain fluxes have to be treated...
Definition: mod_physics.t:26
procedure(sub_get_cbounds), pointer phys_get_cbounds
Definition: mod_physics.t:40
procedure(sub_write_info), pointer phys_write_info
Definition: mod_physics.t:51