17 character(len=LENPACKAGENAME) ::
text =
' PRTFMI'
21 integer(I4B),
public :: max_faces
22 real(dp),
allocatable,
public :: sourceflows(:)
23 real(dp),
allocatable,
public :: sinkflows(:)
24 real(dp),
allocatable,
public :: storageflows(:)
25 real(dp),
allocatable,
public :: boundaryflows(:)
26 integer(I4B),
allocatable,
public :: boundaryfaces(:)
43 subroutine fmi_cr(fmiobj, name_model, input_mempath, inunit, iout)
46 character(len=*),
intent(in) :: name_model
47 character(len=*),
intent(in) :: input_mempath
48 integer(I4B),
intent(inout) :: inunit
49 integer(I4B),
intent(in) :: iout
55 call fmiobj%set_names(1, name_model,
'FMI',
'FMI', input_mempath)
59 call fmiobj%allocate_scalars()
62 fmiobj%inunit = inunit
66 fmiobj%depvartype =
'TRACKS '
78 character(len=15) :: nodestr
79 character(len=*),
parameter :: fmtdry = &
80 &
"(/1X,'WARNING: DRY CELL ENCOUNTERED AT ',a,'; RESET AS INACTIVE')"
81 character(len=*),
parameter :: fmtrewet = &
82 &
"(/1X,'DRY CELL REACTIVATED AT ', a)"
87 this%iflowsupdated = 1
90 if (this%iubud /= 0)
call this%advance_bfr()
93 if (this%iuhds /= 0)
call this%advance_hfr()
96 if (this%iumvr /= 0) &
97 call this%mvrbudobj%bfr_advance(this%dis, this%iout)
100 call this%accumulate_flows()
103 do n = 1, this%dis%nodes
106 if (this%gwfsat(n) > dzero)
then
107 this%ibdgwfsat0(n) = 1
109 this%ibdgwfsat0(n) = 0
113 if (this%ibound(n) > 0)
then
114 if (this%gwfhead(n) ==
dhdry)
then
117 call this%dis%noder_to_string(n, nodestr)
118 write (this%iout, fmtdry) trim(nodestr)
123 if (this%ibound(n) == 0)
then
124 if (this%gwfhead(n) /=
dhdry)
then
127 call this%dis%noder_to_string(n, nodestr)
128 write (this%iout, fmtrewet) trim(nodestr)
139 integer(I4B),
intent(in) :: idryinactive
141 call this%FlowModelInterfaceType%fmi_df(dis, idryinactive)
143 this%max_faces = this%dis%get_max_npolyverts() + 2
144 allocate (this%StorageFlows(this%dis%nodes))
145 allocate (this%SourceFlows(this%dis%nodes))
146 allocate (this%SinkFlows(this%dis%nodes))
147 allocate (this%BoundaryFlows(this%dis%nodes * this%max_faces))
148 allocate (this%BoundaryFaces(this%dis%nodes))
157 integer(I4B) :: j, i, ip, ib
158 integer(I4B) :: ioffset, iflowface, iauxiflowface, icellface
160 character(len=LENAUXNAME) :: auxname
163 this%StorageFlows =
dzero
164 if (this%igwfstrgss /= 0) &
165 this%StorageFlows = this%StorageFlows + this%gwfstrgss
166 if (this%igwfstrgsy /= 0) &
167 this%StorageFlows = this%StorageFlows + this%gwfstrgsy
169 this%SourceFlows =
dzero
170 this%SinkFlows =
dzero
171 this%BoundaryFlows =
dzero
172 this%BoundaryFaces = 0
173 do ip = 1, this%nflowpack
175 naux = this%gwfpackages(ip)%naux
178 auxname = this%gwfpackages(ip)%auxname(j)
179 if (trim(adjustl(auxname)) ==
"IFLOWFACE")
then
185 do ib = 1, this%gwfpackages(ip)%nbound
186 i = this%gwfpackages(ip)%nodelist(ib)
188 if (this%ibound(i) <= 0) cycle
189 qbnd = this%gwfpackages(ip)%get_flow(ib)
193 if (iauxiflowface > 0)
then
194 iflowface = nint(this%gwfpackages(ip)%auxvar(iauxiflowface, ib))
195 icellface = this%iflowface_to_icellface(iflowface)
197 if (icellface > 0)
then
198 call this%mark_boundary_face(i, icellface)
199 ioffset = (i - 1) * this%max_faces
200 this%BoundaryFlows(ioffset + icellface) = &
201 this%BoundaryFlows(ioffset + icellface) + qbnd
202 else if (qbnd .gt.
dzero)
then
203 this%SourceFlows(i) = this%SourceFlows(i) + qbnd
204 else if (qbnd .lt.
dzero)
then
205 this%SinkFlows(i) = this%SinkFlows(i) + qbnd
215 integer(I4B),
intent(in) :: ic
216 integer(I4B),
intent(in) :: iface
218 integer(I4B) :: bit_pos
220 if (ic <= 0 .or. ic > this%dis%nodes)
then
221 print *,
'Invalid cell number: ', ic
222 print *,
'Expected a value in range [1, ', this%dis%nodes,
']'
226 print *,
'Invalid face number: ', iface
227 print *,
'Expected a value in range [1, ', this%max_faces,
']'
231 if (bit_pos < 0 .or. bit_pos > 31)
then
232 print *,
'Invalid bitmask position: ', iface
233 print *,
'Expected a value in range [0, 31]'
236 this%BoundaryFaces(ic) = ibset(this%BoundaryFaces(ic), bit_pos)
242 integer(I4B),
intent(in) :: ic
243 integer(I4B),
intent(in) :: iface
244 logical(LGP) :: is_boundary
246 integer(I4B) :: bit_pos
248 is_boundary = .false.
249 if (ic <= 0 .or. ic > this%dis%nodes)
then
250 print *,
'Invalid cell number: ', ic
251 print *,
'Expected a value in range [1, ', this%dis%nodes,
']'
255 print *,
'Invalid face number: ', iface
256 print *,
'Expected a value in range [1, ', this%max_faces,
']'
260 if (bit_pos < 0 .or. bit_pos > 31)
then
261 print *,
'Invalid bitmask position: ', iface
262 print *,
'Expected a value in range [0, 31]'
265 is_boundary = btest(this%BoundaryFaces(ic), bit_pos)
271 integer(I4B),
intent(in) :: ic
272 integer(I4B),
intent(in) :: iface
273 logical(LGP) :: is_net_out_boundary
275 integer(I4B) :: ioffset
277 is_net_out_boundary = .false.
278 if (.not. this%is_boundary_face(ic, iface))
return
279 ioffset = (ic - 1) * this%max_faces
280 if (this%BoundaryFlows(ioffset + iface) <
dzero) is_net_out_boundary = .true.
287 integer(I4B),
intent(in) :: iflowface
288 integer(I4B) :: iface
291 if (iface < 0) iface = iface + this%max_faces + 1
This module contains simulation constants.
real(dp), parameter dhdry
real dry cell constant
integer(i4b), parameter lenpackagename
maximum length of the package name
integer(i4b), parameter lenvarname
maximum length of a variable name
integer(i4b), parameter lenauxname
maximum length of a aux variable
real(dp), parameter dzero
real constant zero
subroutine pstop(status, message)
Stop the program, optionally specifying an error status code.
This module defines variable data types.
subroutine accumulate_flows(this)
Accumulate flows.
integer(i4b) function iflowface_to_icellface(this, iflowface)
Convert an iflowface number to a cell face number. Maps bottom (-2) -> max_faces - 1,...
subroutine fmi_ad(this)
Time step advance.
subroutine mark_boundary_face(this, ic, iface)
Mark a face as a boundary face.
subroutine prtfmi_df(this, dis, idryinactive)
Define the flow model interface.
logical(lgp) function is_net_out_boundary_face(this, ic, iface)
Check if a face is an assigned boundary with net outflow.
character(len=lenpackagename) text
subroutine, public fmi_cr(fmiobj, name_model, input_mempath, inunit, iout)
Create a new PrtFmi object.
logical(lgp) function is_boundary_face(this, ic, iface)
Check if a face is assigned to a boundary package.
This module contains simulation methods.
subroutine, public store_error(msg, terminate)
Store an error message.
This module contains simulation variables.
character(len=maxcharlen) errmsg
error message string