MODFLOW 6  version 6.7.0.dev2
USGS Modular Hydrologic Model
UpstreamScheme.f90
Go to the documentation of this file.
2  use kindmodule, only: dp, i4b
3  use constantsmodule, only: dzero
6  use basedismodule, only: disbasetype
7  use tspfmimodule, only: tspfmitype
8 
9  implicit none
10  private
11 
12  public :: upstreamschemetype
13 
15  private
16  class(disbasetype), pointer :: dis
17  type(tspfmitype), pointer :: fmi
18  contains
19  procedure :: compute
20  end type upstreamschemetype
21 
22  interface upstreamschemetype
23  module procedure constructor
24  end interface upstreamschemetype
25 
26 contains
27  function constructor(dis, fmi) result(interpolation_scheme)
28  ! -- return
29  type(upstreamschemetype) :: interpolation_scheme
30  ! --dummy
31  class(disbasetype), pointer, intent(in) :: dis
32  type(tspfmitype), pointer, intent(in) :: fmi
33 
34  interpolation_scheme%dis => dis
35  interpolation_scheme%fmi => fmi
36 
37  end function constructor
38 
39  function compute(this, n, m, iposnm, phi) result(phi_face)
40  !-- return
41  type(coefficientstype) :: phi_face
42  ! -- dummy
43  class(upstreamschemetype), target :: this
44  integer(I4B), intent(in) :: n
45  integer(I4B), intent(in) :: m
46  integer(I4B), intent(in) :: iposnm
47  real(dp), intent(in), dimension(:) :: phi
48  ! -- local
49  real(dp) :: qnm
50 
51  ! -- Compute the coefficients for the upwind scheme
52  qnm = this%fmi%gwfflowja(iposnm)
53 
54  if (qnm < dzero) then
55  phi_face%c_n = 1.0_dp
56  else
57  phi_face%c_m = 1.0_dp
58  end if
59 
60  end function compute
61 end module upstreamschememodule
This module contains simulation constants.
Definition: Constants.f90:9
real(dp), parameter dzero
real constant zero
Definition: Constants.f90:65
This module defines variable data types.
Definition: kind.f90:8
type(coefficientstype) function compute(this, n, m, iposnm, phi)
type(upstreamschemetype) function constructor(dis, fmi)