MPI-AMRVAC 3.1
The MPI - Adaptive Mesh Refinement - Versatile Advection Code (development version)
Loading...
Searching...
No Matches
mod_nonlinear_phys.t
Go to the documentation of this file.
1!> Module containing the physics routines for scalar nonlinear equation
3 use mod_comm_lib, only: mpistop
4
5 implicit none
6 private
7
8 !> index of the single scalar unknown
9 integer, protected, public :: rho_ = 1
10
11 !> switch between burgers (i.e. rho**2)
12 !> or nonconvex flux (i.e. rho**3)
13 integer, protected, public :: nonlinear_flux_type = 1
14
15 !> whether the KdV source term is added
16 logical, protected, public :: kdv_source_term = .false.
17
18 !> Whether particles module is added
19 logical, public, protected :: nonlinear_particles = .false.
20
21 ! Public methods
22 public :: nonlinear_phys_init
23 public :: nonlinear_get_v
24
25contains
26
27 !> Read this module's parameters from a file
28 subroutine nonlinear_params_read(files)
30 character(len=*), intent(in) :: files(:)
31 integer :: n
32
34
35 do n = 1, size(files)
36 open(unitpar, file=trim(files(n)), status='old')
37 read(unitpar, nonlinear_list, end=111)
38111 close(unitpar)
39 end do
40
41 end subroutine nonlinear_params_read
42
43 !> Write this module's parameters to a snapshot
44 subroutine nonlinear_write_info(fh)
45 ! for nonlinear scalar equation, nothing to write
46 ! note: this is info only stored at end of dat files,
47 ! is never read/used for restarts, only expects
48 ! an integer (number of parameters) and
49 ! corresponding double values and character names
50 ! and is meant for use in the python tools
52 integer, intent(in) :: fh
53 integer, dimension(MPI_STATUS_SIZE) :: st
54 integer :: er
55
56 ! Write zero parameters
57 call mpi_file_write(fh, 0, 1, mpi_integer, st, er)
58
59 end subroutine nonlinear_write_info
60
63 use mod_physics
64 use mod_kdv, only: kdv_init
66
67 call nonlinear_params_read(par_files)
68
69 physics_type = "nonlinear"
70 phys_energy = .false.
72
73 allocate(start_indices(number_species),stop_indices(number_species))
74 ! set the index of the first flux variable for species 1
75 start_indices(1)=1
76 rho_ = var_set_rho()
77
78 ! set number of variables which need update ghostcells
79 nwgc=nwflux
80
81 ! set the index of the last flux variable for species 1
82 stop_indices(1)=nwflux
83
84 ! Check whether custom flux types have been defined
85 if (.not. allocated(flux_type)) then
86 allocate(flux_type(ndir, nw))
88 else if (any(shape(flux_type) /= [ndir, nw])) then
89 call mpistop("phys_check error: flux_type has wrong shape")
90 end if
91
92 phys_get_cmax => nonlinear_get_cmax
93 phys_get_cbounds => nonlinear_get_cbounds
94 phys_get_flux => nonlinear_get_flux
95 phys_add_source_geom => nonlinear_add_source_geom
96 phys_add_source => nonlinear_add_source
97 phys_to_conserved => nonlinear_to_conserved
98 phys_to_primitive => nonlinear_to_primitive
99 phys_get_dt => nonlinear_get_dt
100 phys_write_info => nonlinear_write_info
101
102 if (kdv_source_term) call kdv_init()
103
104 ! Initialize particles module
105 if (nonlinear_particles) then
106 call particles_init()
107 end if
108
109 end subroutine nonlinear_phys_init
110
111 subroutine nonlinear_to_conserved(ixI^L, ixO^L, w, x)
113 integer, intent(in) :: ixi^l, ixo^l
114 double precision, intent(inout) :: w(ixi^s, nw)
115 double precision, intent(in) :: x(ixi^s, 1:^nd)
116
117 ! Do nothing (primitive and conservative are equal for nonlinear module)
118 end subroutine nonlinear_to_conserved
119
120 subroutine nonlinear_to_primitive(ixI^L, ixO^L, w, x)
122 integer, intent(in) :: ixi^l, ixo^l
123 double precision, intent(inout) :: w(ixi^s, nw)
124 double precision, intent(in) :: x(ixi^s, 1:^nd)
125
126 ! Do nothing (primitive and conservative are equal for nonlinear module)
127 end subroutine nonlinear_to_primitive
128
129 subroutine nonlinear_get_v(w, x, ixI^L, ixO^L, idim, v)
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(out) :: v(ixi^s)
134
135 select case(nonlinear_flux_type)
136 case(1)
137 v(ixo^s)=w(ixo^s,rho_)
138 case(2)
139 v(ixo^s)=3.0d0*w(ixo^s,rho_)**2
140 case default
141 call mpistop('Undefined fluxtype: set nonlinear_flux_type to 1 or 2')
142 end select
143
144 end subroutine nonlinear_get_v
145
146 subroutine nonlinear_get_cmax(w, x, ixI^L, ixO^L, idim, cmax)
148 integer, intent(in) :: ixi^l, ixo^l, idim
149 double precision, intent(in) :: w(ixi^s, nw), x(ixi^s, 1:^nd)
150 double precision, intent(inout) :: cmax(ixi^s)
151
152 call nonlinear_get_v(w, x, ixi^l, ixo^l, idim, cmax)
153
154 cmax(ixo^s) = abs(cmax(ixo^s))
155
156 end subroutine nonlinear_get_cmax
157
158 subroutine nonlinear_get_cbounds(wLC, wRC, wLp, wRp, x, ixI^L, ixO^L, idim,Hspeed, cmax, cmin)
160 use mod_variables
161 integer, intent(in) :: ixi^l, ixo^l, idim
162 double precision, intent(in) :: wlc(ixi^s, nw), wrc(ixi^s,nw)
163 double precision, intent(in) :: wlp(ixi^s, nw), wrp(ixi^s,nw)
164 double precision, intent(in) :: x(ixi^s, 1:^nd)
165 double precision, intent(inout) :: cmax(ixi^s,1:number_species)
166 double precision, intent(inout), optional :: cmin(ixi^s,1:number_species)
167 double precision, intent(in) :: hspeed(ixi^s,1:number_species)
168 double precision :: wmean(ixi^s,nw)
169
170 ! since get_v depends on w, the first argument should be some average over the
171 ! left and right state
172 wmean(ixo^s,1:nwflux)=0.5d0*(wlc(ixo^s,1:nwflux)+wrc(ixo^s,1:nwflux))
173 call nonlinear_get_v(wmean, x, ixi^l, ixo^l, idim, cmax(ixi^s,1))
174
175 if (present(cmin)) then
176 cmin(ixo^s,1) = min(cmax(ixo^s,1), zero)
177 cmax(ixo^s,1) = max(cmax(ixo^s,1), zero)
178 else
179 cmax(ixo^s,1) = maxval(abs(cmax(ixo^s,1)))
180 end if
181
182 end subroutine nonlinear_get_cbounds
183
184 subroutine nonlinear_get_dt(w, ixI^L, ixO^L, dtnew, dx^D, x)
186 use mod_kdv, only: kdv_get_dt
187
188 integer, intent(in) :: ixi^l, ixo^l
189 double precision, intent(in) :: dx^d, x(ixi^s, 1:^nd)
190 double precision, intent(in) :: w(ixi^s, 1:nw)
191 double precision, intent(inout) :: dtnew
192
193 dtnew = bigdouble
194
195 if(kdv_source_term) then
196 call kdv_get_dt(w, ixi^l, ixo^l, dtnew, dx^d, x)
197 endif
198 end subroutine nonlinear_get_dt
199
200 ! here we select the flux according to the nonlinear_flux_type parameter
201 subroutine nonlinear_get_flux(wC, w, x, ixI^L, ixO^L, idim, f)
203 integer, intent(in) :: ixi^l, ixo^l, idim
204 double precision, intent(in) :: wc(ixi^s, 1:nw)
205 double precision, intent(in) :: w(ixi^s, 1:nw)
206 double precision, intent(in) :: x(ixi^s, 1:^nd)
207 double precision, intent(out) :: f(ixi^s, nwflux)
208
209 select case(nonlinear_flux_type)
210 case(1)
211 f(ixo^s,rho_)=half*w(ixo^s,rho_)**2
212 case(2)
213 f(ixo^s,rho_)=w(ixo^s,rho_)**3
214 case default
215 call mpistop('Undefined fluxtype: set nonlinear_flux_type to 1 or 2')
216 end select
217
218 end subroutine nonlinear_get_flux
219
220 subroutine nonlinear_add_source_geom(qdt, dtfactor, ixI^L, ixO^L, wCT,wprim, w, x)
221
222 ! Add geometrical source terms to w
223 ! There are no geometrical source terms
224
226
227 integer, intent(in) :: ixi^l, ixo^l
228 double precision, intent(in) :: qdt, dtfactor, x(ixi^s, 1:^nd)
229 double precision, intent(inout) :: wct(ixi^s, 1:nw), wprim(ixi^s,1:nw),w(ixi^s, 1:nw)
230
231 end subroutine nonlinear_add_source_geom
232
233 subroutine nonlinear_add_source(qdt, dtfactor, ixI^L,ixO^L,wCT,wCTprim,w,x,qsourcesplit,active)
235 use mod_kdv, only: kdv_add_source
236
237 integer, intent(in) :: ixi^l, ixo^l
238 double precision, intent(in) :: qdt,dtfactor
239 double precision, intent(in) :: wct(ixi^s, 1:nw),wctprim(ixi^s,1:nw), x(ixi^s, 1:ndim)
240 double precision, intent(inout) :: w(ixi^s, 1:nw)
241 logical, intent(in) :: qsourcesplit
242 logical, intent(inout) :: active
243
244 if(kdv_source_term) then
245 call kdv_add_source(qdt,ixi^l,ixo^l,wct,w,x,qsourcesplit,active)
246 end if
247
248 end subroutine nonlinear_add_source
249
250end module mod_nonlinear_phys
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...
integer, parameter unitpar
file handle for IO
integer, parameter ndim
Number of spatial dimensions for grid variables.
logical use_particles
Use particles module or not.
character(len=std_len), dimension(:), allocatable par_files
Which par files are used as input.
double precision, dimension(:), allocatable, parameter d
integer ndir
Number of spatial dimensions (components) for vector variables.
double precision, dimension(:,:), allocatable dx
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
subroutine kdv_get_dt(w, ixil, ixol, dtnew, dxd, x)
Definition mod_kdv.t:96
Module containing the physics routines for scalar nonlinear equation.
logical, public, protected nonlinear_particles
Whether particles module is added.
integer, public, protected nonlinear_flux_type
switch between burgers (i.e. rho**2) or nonconvex flux (i.e. rho**3)
logical, public, protected kdv_source_term
whether the KdV source term is added
integer, public, protected rho_
index of the single scalar unknown
subroutine, public nonlinear_phys_init()
subroutine, public nonlinear_get_v(w, x, ixil, ixol, idim, v)
Module containing all the particle routines.
subroutine particles_init()
Initialize particle data and parameters.
This module defines the procedures of a physics module. It contains function pointers for the various...
Definition mod_physics.t:4
procedure(sub_convert), pointer phys_to_primitive
Definition mod_physics.t:55
procedure(sub_write_info), pointer phys_write_info
Definition mod_physics.t:77
procedure(sub_get_flux), pointer phys_get_flux
Definition mod_physics.t:64
procedure(sub_get_cbounds), pointer phys_get_cbounds
Definition mod_physics.t:63
procedure(sub_get_dt), pointer phys_get_dt
Definition mod_physics.t:67
procedure(sub_add_source_geom), pointer phys_add_source_geom
Definition mod_physics.t:68
integer, parameter flux_default
Indicates a normal flux.
Definition mod_physics.t:24
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:21
procedure(sub_convert), pointer phys_to_conserved
Definition mod_physics.t:54
character(len=name_len) physics_type
String describing the physics type of the simulation.
Definition mod_physics.t:50
procedure(sub_add_source), pointer phys_add_source
Definition mod_physics.t:69
procedure(sub_get_cmax), pointer phys_get_cmax
Definition mod_physics.t:57
logical phys_energy
Solve energy equation or not.
Definition mod_physics.t:35
integer nw
Total number of variables.
integer number_species
number of species: each species has different characterictic speeds and should be used accordingly in...
integer nwflux
Number of flux variables.