MPI-AMRVAC  2.2
The MPI - Adaptive Mesh Refinement - Versatile Advection Code
mod_source.t
Go to the documentation of this file.
1 !> Module for handling split source terms (split from the fluxes)
2 module mod_source
3 
4  implicit none
5  private
6 
7  public :: add_split_source
8  public :: addsource2
9 
10 contains
11 
12  subroutine add_split_source(prior)
17 
18  logical, intent(in) :: prior
19 
20  double precision :: qdt, qt
21  integer :: iigrid, igrid, i^D
22  logical :: src_active
23 
24  ! add thermal conduction
26 
27  src_active = .false.
28 
29  if (prior .and. associated(phys_global_source_before)) then
30  call phys_global_source_before(dt, qt, src_active)
31  end if
32 
33  if ((.not.prior).and.&
34  (typesourcesplit=='sf' .or. typesourcesplit=='ssf')) return
35 
36  if (prior) then
37  qt=global_time
38  else
39  qt=global_time+dt
40  end if
41 
42  !$OMP PARALLEL DO PRIVATE(igrid,qdt,i^D)
43  do iigrid=1,igridstail_active; igrid=igrids_active(iigrid);
44  qdt=dt_grid(igrid)
45  block=>ps(igrid)
46  call addsource1_grid(igrid,qdt,qt,ps(igrid)%w,src_active)
47  end do
48  !$OMP END PARALLEL DO
49 
50  if (.not. prior .and. associated(phys_global_source_after)) then
51  call phys_global_source_after(dt, qt, src_active)
52  end if
53 
54  if (src_active) then
55  call getbc(qt,0.d0,ps,1,nwflux+nwaux,phys_req_diagonal)
56  end if
57 
58  end subroutine add_split_source
59 
60  subroutine addsource1_grid(igrid,qdt,qt,w,src_active)
61 
63 
64  integer, intent(in) :: igrid
65  double precision, intent(in) :: qdt, qt
66  double precision, intent(inout) :: w(ixg^t,nw)
67  logical, intent(inout) :: src_active
68 
69  double precision :: w1(ixg^t,nw)
70 
71  saveigrid=igrid
74 
75  ^d&dxlevel(^d)=rnode(rpdx^d_,igrid);
76 
77  w1(ixg^t,1:nw)=w(ixg^t,1:nw)
78 
79  select case (typesourcesplit)
80  case ('sf')
81  call addsource2(qdt ,ixg^ll,ixm^ll,1,nw,qt,w1,qt,w,&
82  ps(igrid)%x,.true.,src_active)
83  case ('sfs')
84  call addsource2(qdt/2,ixg^ll,ixm^ll,1,nw,qt,w1,qt,w,&
85  ps(igrid)%x,.true.,src_active)
86  case ('ssf')
87  call addsource2(qdt/2,ixg^ll,ixg^ll,1,nw,qt,w,qt,w1,&
88  ps(igrid)%x,.true.,src_active)
89  call addsource2(qdt ,ixg^ll,ixm^ll,1,nw,qt,w1,qt,w,&
90  ps(igrid)%x,.true.,src_active)
91  case ('ssfss')
92  call addsource2(qdt/4,ixg^ll,ixg^ll,1,nw,qt,w,qt,w1,&
93  ps(igrid)%x,.true.,src_active)
94  call addsource2(qdt/2,ixg^ll,ixm^ll,1,nw,qt,w1,qt,w,&
95  ps(igrid)%x,.true.,src_active)
96  case default
97  write(unitterm,*)'No such typesourcesplit=',typesourcesplit
98  call mpistop("Error: Unknown typesourcesplit!")
99  end select
100 
101  end subroutine addsource1_grid
102 
103  !> Add source within ixO for iws: w=w+qdt*S[wCT]
104  subroutine addsource2(qdt,ixI^L,ixO^L,iw^LIM,qtC,wCT,qt,&
105  w,x,qsourcesplit,src_active)
107  use mod_physics, only: phys_add_source
108  use mod_usr_methods, only: usr_source
109  ! differences with VAC is in iw^LIM and in declaration of ranges for wCT,w
110 
111  integer, intent(in) :: ixI^L, ixO^L, iw^LIM
112  double precision, intent(in) :: qdt, qtC, qt
113  double precision, intent(in) :: wCT(ixi^s,1:nw), x(ixi^s,1:ndim)
114  double precision, intent(inout) :: w(ixi^s,1:nw)
115  logical, intent(in) :: qsourcesplit
116  logical, intent(inout), optional :: src_active
117  logical :: tmp_active
118 
119  tmp_active = .false.
120 
121  ! physics defined sources, typically explicitly added,
122  ! along with geometrical source additions
123  call phys_add_source(qdt,ixi^l,ixo^l,wct,w,x,qsourcesplit,tmp_active)
124 
125  ! user defined sources, typically explicitly added
126  if ((qsourcesplit .eqv. source_split_usr) .and. associated(usr_source)) then
127  tmp_active = .true.
128  call usr_source(qdt,ixi^l,ixo^l,iw^lim,qtc,wct,qt,w,x)
129  end if
130 
131  if (present(src_active)) src_active = src_active .or. tmp_active
132 
133  end subroutine addsource2
134 
135 end module mod_source
This module contains definitions of global parameters and variables and some generic functions/subrou...
Module for handling split source terms (split from the fluxes)
Definition: mod_source.t:2
procedure(sub_global_source), pointer phys_global_source_after
Definition: mod_physics.t:76
update ghost cells of all blocks including physical boundaries
integer, dimension(:), allocatable type_gradient_limiter
Type of slope limiter used for computing gradients or divergences, when typegrad or typediv are set t...
integer, parameter plevel_
integer, dimension(:), allocatable type_limiter
Type of slope limiter used for reconstructing variables on cell edges.
logical phys_req_diagonal
Whether the physics routines require diagonal ghost cells, for example for computing a curl...
Definition: mod_physics.t:22
Module with all the methods that users can customize in AMRVAC.
character(len=std_len) typesourcesplit
How to apply dimensional splitting to the source terms, see Discretization.
procedure(source), pointer usr_source
subroutine, public add_split_source(prior)
Definition: mod_source.t:13
Thermal conduction for HD and MHD.
integer ixm
the mesh range (within a block with ghost cells)
double precision, dimension(:), allocatable dt_grid
subroutine getbc(time, qdt, psb, nwstart, nwbc, req_diag)
do update ghost cells of all blocks including physical boundaries
integer, parameter unitterm
Unit for standard output.
integer, dimension(:), allocatable, parameter d
double precision global_time
The global simulation time.
subroutine addsource1_grid(igrid, qdt, qt, w, src_active)
Definition: mod_source.t:61
procedure(sub_global_source), pointer phys_global_source_before
Definition: mod_physics.t:75
subroutine mpistop(message)
Exit MPI-AMRVAC with an error message.
Definition: comm_lib.t:198
double precision, dimension(ndim) dxlevel
integer, parameter ndim
Number of spatial dimensions for grid variables.
procedure(thermal_conduction), pointer phys_thermal_conduction
This module defines the procedures of a physics module. It contains function pointers for the various...
Definition: mod_physics.t:4
procedure(sub_add_source), pointer phys_add_source
Definition: mod_physics.t:74
logical source_split_usr
Use split or unsplit way to add user's source terms, default: unsplit.
integer, dimension(:,:), allocatable node
double precision, dimension(:,:), allocatable rnode
Corner coordinates.
subroutine, public addsource2(qdt, ixIL, ixOL, iwLIM, qtC, wCT, qt, w, x, qsourcesplit, src_active)
Add source within ixO for iws: w=w+qdt*S[wCT].
Definition: mod_source.t:106