2D AMRCLAW
step2x.f90
Go to the documentation of this file.
1 !
2 ! -------------------------------------------------------------
3 !
4 subroutine step2x(maxm,meqn,maux,mbc,mx,my,qold,aux,dx,dt,cflgrid,fm,fp,rpn2)
5 !
6 ! clawpack routine ... modified for AMRCLAW
7 !
8 ! Take one time step, updating q.
9 ! On entry, qold gives
10 ! initial data for this step
11 ! and is unchanged in this version.
12 !
13 ! fm, fp are fluxes to left and right of single cell edge
14 ! See the flux2 documentation for more information.
15 !
16 ! Converted to f90 2012-1-04 (KTM)
17 !
18 
19  use amr_module
20 
21  implicit none
22 
23  external rpn2
24 
25  ! Arguments
26  integer, intent(in) :: maxm,meqn,maux,mbc,mx,my
27  real(kind=8), intent(in) :: dx,dt
28  real(kind=8), intent(inout) :: cflgrid
29  real(kind=8), intent(inout) :: qold(meqn, 1-mbc:mx+mbc, 1-mbc:my+mbc)
30  real(kind=8), intent(inout) :: aux(maux,1-mbc:mx+mbc, 1-mbc:my+mbc)
31  real(kind=8), intent(inout) :: fm(meqn, 1-mbc:mx+mbc, 1-mbc:my+mbc)
32  real(kind=8), intent(inout) :: fp(meqn,1-mbc:mx+mbc, 1-mbc:my+mbc)
33 
34 
35  ! Local storage for flux accumulation
36  real(kind=8) :: faddm(meqn,1-mbc:maxm+mbc)
37  real(kind=8) :: faddp(meqn,1-mbc:maxm+mbc)
38 
39 
40  ! Scratch storage for Sweeps and Riemann problems
41  real(kind=8) :: q1d(meqn,1-mbc:maxm+mbc)
42  real(kind=8) :: aux2(maux,1-mbc:maxm+mbc)
43  real(kind=8) :: dtdx1d(1-mbc:maxm+mbc)
44 
45 
46  real(kind=8) :: wave(meqn, mwaves, 1-mbc:maxm+mbc)
47  real(kind=8) :: s(mwaves, 1-mbc:maxm + mbc)
48  real(kind=8) :: cqxx(meqn,1-mbc:maxm + mbc)
49 
50 
51  ! Looping scalar storage
52  integer :: j
53  real(kind=8) :: dtdx,cfl1d
54 
55 
56  cflgrid = 0.d0
57  dtdx = dt/dx
58 
59  fm = 0.d0
60  fp = 0.d0
61 
62  ! ============================================================================
63  ! Perform X-Sweeps (for Godunov split method, called from stepgrid_dimSplit)
64  ! loop indices assume x sweep goes first
65  do j = 1-mbc,my+mbc
66 
67  ! Copy old q into 1d slice
68  q1d(:,1-mbc:mx+mbc) = qold(:,1-mbc:mx+mbc,j)
69 
70  ! Set dtdx slice if a capacity array exists
71  if (mcapa > 0) then
72  dtdx1d(1-mbc:mx+mbc) = dtdx / aux(mcapa,1-mbc:mx+mbc,j)
73  else
74  dtdx1d = dtdx
75  endif
76 
77  ! Copy aux array into slices
78  if (maux > 0) then
79 
80  aux2(:,1-mbc:mx+mbc) = aux(:,1-mbc:mx+mbc,j )
81 
82  endif
83 
84  ! Compute modifications fadd and gadd to fluxes along this slice:
85 !!! ::: dont need aux1 and aux3 for dimensionally split method
86 !!! ::: but flux2 needs 3 aux arrays, so passing in aux2 3 times
87  call flux2_dimsplit(1,maxm,meqn,maux,mbc,mx,q1d,dtdx1d,aux2, &
88  faddm,faddp,cfl1d,wave,s,cqxx,rpn2)
89 
90  cflgrid = max(cflgrid,cfl1d)
91 
92  ! Update fluxes
93  fm(:,1:mx+1,j) = fm(:,1:mx+1,j) + faddm(:,1:mx+1)
94  fp(:,1:mx+1,j) = fp(:,1:mx+1,j) + faddp(:,1:mx+1)
95 
96  enddo
97 
98 
99 end subroutine step2x
subroutine flux2_dimsplit(ixy, maxm, meqn, maux, mbc, mx, q1d, dtdx1d, aux2, faddm, faddp, cfl1d, wave, s, cqxx, rpn2)
Definition: flux2_dimSplit.f:8
integer mcapa
Definition: amr_module.f90:253
subroutine step2x(maxm, meqn, maux, mbc, mx, my, qold, aux, dx, dt, cflgrid, fm, fp, rpn2)
Definition: step2x.f90:5
subroutine rpn2(ixy, maxm, meqn, mwaves, maux, mbc, mx, ql, qr, auxl, auxr, wave, s, amdq, apdq)
Definition: rpn2from1d.f:7
The module contains the definition of a "node descriptor" as well as other global variables used duri...
Definition: amr_module.f90:21
integer mwaves
Definition: amr_module.f90:253