MPI-AMRVAC  2.0
The MPI - Adaptive Mesh Refinement - Versatile Advection Code
boundary_conditions.t
Go to the documentation of this file.
1 !=============================================================================
2 subroutine bc_phys(iside,idims,time,qdt,w,x,ixG^L,ixB^L)
4 use mod_bc_data, only: bc_data_set
6 
7 integer, intent(in) :: iside, idims, ixG^L,ixB^L
8 double precision, intent(in) :: time,qdt
9 double precision, intent(inout) :: w(ixg^s,1:nw)
10 double precision, intent(in) :: x(ixg^s,1:ndim)
11 double precision :: wtmp(ixg^s,1:nwflux)
12 
13 integer :: iw, iB, ix^D, ixI^L, ixM^L, nghostcellsi,iib^D
14 logical :: isphysbound
15 !-----------------------------------------------------------------------------
16 select case (idims)
17 {case (^d)
18  if (iside==2) then
19  ! maximal boundary
20  ib=ismax^d
21  iximin^dd=ixbmax^d+1-nghostcells^d%ixImin^dd=ixbmin^dd;
22  iximax^dd=ixbmax^dd;
23  ! cont/symm/asymm types
24  do iw=1,nwflux+nwaux
25  select case (typeboundary(iw,ib))
26  case ("symm")
27  w(ixi^s,iw) = w(iximin^d-1:iximin^d-nghostcells:-1^d%ixI^s,iw)
28  case ("asymm")
29  w(ixi^s,iw) =-w(iximin^d-1:iximin^d-nghostcells:-1^d%ixI^s,iw)
30  case ("cont")
31  do ix^d=iximin^d,iximax^d
32  w(ix^d^d%ixI^s,iw) = w(iximin^d-1^d%ixI^s,iw)
33  end do
34  case("noinflow")
35  if (iw==1+^d)then
36  do ix^d=iximin^d,iximax^d
37  w(ix^d^d%ixI^s,iw) = max(w(iximin^d-1^d%ixI^s,iw),zero)
38  end do
39  else
40  do ix^d=iximin^d,iximax^d
41  w(ix^d^d%ixI^s,iw) = w(iximin^d-1^d%ixI^s,iw)
42  end do
43  end if
44  case ("special", "bc_data")
45  ! skip it here, do AFTER all normal type boundaries are set
46  case ("character")
47  ! skip it here, do AFTER all normal type boundaries are set
48  case ("aperiodic")
49  !this just multiplies the variables with (-), they have been set from neighbors just like periodic.
50  w(ixi^s,iw) = - w(ixi^s,iw)
51  case ("periodic")
52 ! call mpistop("periodic bc info should come from neighbors")
53  case default
54  write (unitterm,*) "Undefined boundarytype ",typeboundary(iw,ib), &
55  "for variable iw=",iw," and side iB=",ib
56  end select
57  end do
58  else
59  ! minimal boundary
60  ib=ismin^d
61  iximin^dd=ixbmin^dd;
62  iximax^dd=ixbmin^d-1+nghostcells^d%ixImax^dd=ixbmax^dd;
63  ! cont/symm/asymm types
64  do iw=1,nwflux+nwaux
65  select case (typeboundary(iw,ib))
66  case ("symm")
67  w(ixi^s,iw) = w(iximax^d+nghostcells:iximax^d+1:-1^d%ixI^s,iw)
68  case ("asymm")
69  w(ixi^s,iw) =-w(iximax^d+nghostcells:iximax^d+1:-1^d%ixI^s,iw)
70  case ("cont")
71  do ix^d=iximin^d,iximax^d
72  w(ix^d^d%ixI^s,iw) = w(iximax^d+1^d%ixI^s,iw)
73  end do
74  case("noinflow")
75  if (iw==1+^d)then
76  do ix^d=iximin^d,iximax^d
77  w(ix^d^d%ixI^s,iw) = min(w(iximax^d+1^d%ixI^s,iw),zero)
78  end do
79  else
80  do ix^d=iximin^d,iximax^d
81  w(ix^d^d%ixI^s,iw) = w(iximax^d+1^d%ixI^s,iw)
82  end do
83  end if
84  case ("special", "bc_data")
85  ! skip it here, do AFTER all normal type boundaries are set
86  case ("character")
87  ! skip it here, do AFTER all normal type boundaries are set
88  case ("aperiodic")
89  !this just multiplies the variables with (-), they have been set from neighbors just like periodic.
90  w(ixi^s,iw) = - w(ixi^s,iw)
91  case ("periodic")
92 ! call mpistop("periodic bc info should come from neighbors")
93  case default
94  write (unitterm,*) "Undefined boundarytype ",typeboundary(iw,ib), &
95  "for variable iw=",iw," and side iB=",ib
96  end select
97  end do
98  end if \}
99 end select
100 
101 ! do special case AFTER all normal cases are set
102 !do iw=1,nwflux+nwaux
103 ! opedit: iw==0 since this breaks fewest of setups.
104 if (any(typeboundary(1:nwflux+nwaux,ib)=="special")) then
105  if (.not. associated(usr_special_bc)) &
106  call mpistop("usr_special_bc not defined")
107  call usr_special_bc(time,ixg^l,ixi^l,ib,w,x)
108 end if
109 
110 if (any(typeboundary(1:nwflux+nwaux,ib)=="bc_data")) then
111  call bc_data_set(time,ixg^l,ixi^l,ib,w,x)
112 end if
113 
114 {#IFDEF EVOLVINGBOUNDARY
115 if (any(typeboundary(1:nwflux,ib)=="character")) then
116  ixm^l=ixm^ll;
117  if(ixgmax1==ixghi1) then
118  nghostcellsi=nghostcells
119  else
120  nghostcellsi=ceiling(nghostcells*0.5d0)
121  end if
122  select case (idims)
123  {case (^d)
124  if (iside==2) then
125  ! maximal boundary
126  iximin^dd=ixgmax^d+1-nghostcellsi^d%ixImin^dd=ixbmin^dd;
127  iximax^dd=ixbmax^dd;
128  if(all(w(ixi^s,1:nwflux)==0.d0)) then
129  do ix^d=iximin^d,iximax^d
130  w(ix^d^d%ixI^s,1:nwflux) = w(iximin^d-1^d%ixI^s,1:nwflux)
131  end do
132  end if
133  if(qdt>0.d0.and.ixgmax^d==ixghi^d) then
134  iximin^dd=iximin^d^d%ixImin^dd=ixmmin^dd;
135  iximax^dd=iximax^d^d%ixImax^dd=ixmmax^dd;
136  wtmp(ixg^s,1:nw)=pso(saveigrid)%w(ixg^s,1:nw)
137  call characteristic_project(idims,iside,ixg^l,ixi^l,wtmp,x,dxlevel,qdt)
138  w(ixi^s,1:nwflux)=wtmp(ixi^s,1:nwflux)
139  end if
140  else
141  ! minimal boundary
142  iximin^dd=ixbmin^dd;
143  iximax^dd=ixgmin^d-1+nghostcellsi^d%ixImax^dd=ixbmax^dd;
144  if(all(w(ixi^s,1:nwflux)==0.d0)) then
145  do ix^d=iximin^d,iximax^d
146  w(ix^d^d%ixI^s,1:nwflux) = w(iximax^d+1^d%ixI^s,1:nwflux)
147  end do
148  end if
149  if(qdt>0.d0.and.ixgmax^d==ixghi^d) then
150  iximin^dd=iximin^d^d%ixImin^dd=ixmmin^dd;
151  iximax^dd=iximax^d^d%ixImax^dd=ixmmax^dd;
152  wtmp(ixg^s,1:nw)=pso(saveigrid)%w(ixg^s,1:nw)
153  call characteristic_project(idims,iside,ixg^l,ixi^l,wtmp,x,dxlevel,qdt)
154  w(ixi^s,1:nwflux)=wtmp(ixi^s,1:nwflux)
155  end if
156  end if \}
157  end select
158  if(ixgmax1==ixghi1) then
159  call identifyphysbound(saveigrid,isphysbound,iib^d)
160  if(iib1==-1.and.iib2==-1) then
161  do ix2=nghostcells,1,-1
162  do ix1=nghostcells,1,-1
163  w(ix^d,1:nwflux)=(w(ix1+1,ix2+1,1:nwflux)+w(ix1+1,ix2,1:nwflux)+w(ix1,ix2+1,1:nwflux))/3.d0
164  end do
165  end do
166  end if
167  if(iib1== 1.and.iib2==-1) then
168  do ix2=nghostcells,1,-1
169  do ix1=ixmmax1+1,ixgmax1
170  w(ix^d,1:nwflux)=(w(ix1-1,ix2+1,1:nwflux)+w(ix1-1,ix2,1:nwflux)+w(ix1,ix2+1,1:nwflux))/3.d0
171  end do
172  end do
173  end if
174  if(iib1==-1.and.iib2== 1) then
175  do ix2=ixmmax2+1,ixgmax2
176  do ix1=nghostcells,1,-1
177  w(ix^d,1:nwflux)=(w(ix1+1,ix2-1,1:nwflux)+w(ix1+1,ix2,1:nwflux)+w(ix1,ix2-1,1:nwflux))/3.d0
178  end do
179  end do
180  end if
181  if(iib1== 1.and.iib2== 1) then
182  do ix2=ixmmax2+1,ixgmax2
183  do ix1=ixmmax1+1,ixgmax1
184  w(ix^d,1:nwflux)=(w(ix1-1,ix2-1,1:nwflux)+w(ix1-1,ix2,1:nwflux)+w(ix1,ix2-1,1:nwflux))/3.d0
185  end do
186  end do
187  end if
188  end if
189 end if
190 }
191 !end do
192 
193 end subroutine bc_phys
194 !=============================================================================
195 subroutine getintbc(time,ixG^L)
198 
199 double precision, intent(in) :: time
200 integer, intent(in) :: ixG^L
201 
202 ! .. local ..
203 integer :: iigrid, igrid, ixO^L,level
204 !----------------------------------------------------------------------------
205 ixo^l=ixg^l^lsubnghostcells;
206 
207 do iigrid=1,igridstail_active; igrid=igrids_active(iigrid);
208 !do iigrid=1,igridstail; igrid=igrids(iigrid);
209  ^d&dxlevel(^d)=rnode(rpdx^d_,igrid);
210  block=>ps(igrid)
213  level=node(plevel_,igrid)
214  saveigrid=igrid
215 
216  if (associated(usr_internal_bc)) then
217  call usr_internal_bc(level,time,ixg^l,ixo^l,ps(igrid)%w,ps(igrid)%x)
218  end if
219 end do
220 
221 end subroutine getintbc
222 !=============================================================================
This module contains definitions of global parameters and variables and some generic functions/subrou...
procedure(internal_bc), pointer usr_internal_bc
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.
subroutine, public bc_data_set(qt, ixIL, ixOL, iB, w, x)
Set boundary conditions according to user data.
Definition: mod_bc_data.t:125
Module to set boundary conditions from user data.
Definition: mod_bc_data.t:2
integer nghostcells
Number of ghost cells surrounding a grid.
Module with all the methods that users can customize in AMRVAC.
integer ixghi
Upper index of grid block arrays.
character(len=std_len), dimension(:, :), allocatable typeboundary
Array indicating the type of boundary condition per variable and per physical boundary.
integer, parameter unitterm
Unit for standard output.
integer, dimension(:), allocatable, parameter d
subroutine characteristic_project(idims, iside, ixIL, ixOL, w, x, dxndim, qdt)
Definition: eigensystem.t:182
subroutine mpistop(message)
Exit MPI-AMRVAC with an error message.
Definition: comm_lib.t:160
double precision, dimension(ndim) dxlevel
integer, parameter ndim
Number of spatial dimensions for grid variables.
integer, dimension(:,:), allocatable node
double precision, dimension(:,:), allocatable rnode
Corner coordinates.
subroutine getintbc(time, ixGL)
procedure(special_bc), pointer usr_special_bc
subroutine bc_phys(iside, idims, time, qdt, w, x, ixGL, ixBL)