MODFLOW 6  version 6.7.0.dev3
USGS Modular Hydrologic Model
centraldifferenceschememodule Module Reference

Data Types

interface  centraldifferenceschemetype
 

Functions/Subroutines

type(centraldifferenceschemetype) function constructor (dis, fmi)
 
subroutine set_field (this, phi)
 Set the scalar field for which interpolation will be computed. More...
 
type(coefficientstype) function compute (this, n, m, iposnm)
 

Function/Subroutine Documentation

◆ compute()

type(coefficientstype) function centraldifferenceschememodule::compute ( class(centraldifferenceschemetype), target  this,
integer(i4b), intent(in)  n,
integer(i4b), intent(in)  m,
integer(i4b), intent(in)  iposnm 
)
private

Definition at line 54 of file CentralDifferenceScheme.f90.

55  !-- return
56  type(CoefficientsType) :: phi_face
57  ! -- dummy
58  class(CentralDifferenceSchemeType), target :: this
59  integer(I4B), intent(in) :: n
60  integer(I4B), intent(in) :: m
61  integer(I4B), intent(in) :: iposnm
62  ! -- local
63  real(DP) :: lnm, lmn, omega
64 
65  ! -- calculate weight based on distances between nodes and the shared
66  ! face of the connection
67  if (this%dis%con%ihc(this%dis%con%jas(iposnm)) == 0) then
68  ! -- vertical connection; assume cell is fully saturated
69  lnm = dhalf * (this%dis%top(n) - this%dis%bot(n))
70  lmn = dhalf * (this%dis%top(m) - this%dis%bot(m))
71  else
72  ! -- horizontal connection
73  lnm = this%dis%con%cl1(this%dis%con%jas(iposnm))
74  lmn = this%dis%con%cl2(this%dis%con%jas(iposnm))
75  end if
76 
77  omega = lmn / (lnm + lmn)
78 
79  phi_face%c_n = omega
80  phi_face%c_m = done - omega
81 

◆ constructor()

type(centraldifferenceschemetype) function centraldifferenceschememodule::constructor ( class(disbasetype), intent(in), pointer  dis,
type(tspfmitype), intent(in), pointer  fmi 
)
private

Definition at line 29 of file CentralDifferenceScheme.f90.

30  ! -- return
31  type(CentralDifferenceSchemeType) :: interpolation_scheme
32  ! --dummy
33  class(DisBaseType), pointer, intent(in) :: dis
34  type(TspFmiType), pointer, intent(in) :: fmi
35 
36  interpolation_scheme%dis => dis
37  interpolation_scheme%fmi => fmi
38 
Here is the caller graph for this function:

◆ set_field()

subroutine centraldifferenceschememodule::set_field ( class(centraldifferenceschemetype), target  this,
real(dp), dimension(:), intent(in), pointer  phi 
)
private

This method establishes a pointer to the scalar field data for subsequent central difference interpolation computations.

Definition at line 46 of file CentralDifferenceScheme.f90.

47  ! -- dummy
48  class(CentralDifferenceSchemeType), target :: this
49  real(DP), intent(in), dimension(:), pointer :: phi
50 
51  this%phi => phi