182 class(GwfHfbType) :: this
183 integer(I4B) :: kiter
184 class(MatrixBaseType),
pointer :: matrix_sln
185 integer(I4B),
intent(in),
dimension(:) :: idxglo
186 real(DP),
intent(inout),
dimension(:) :: rhs
187 real(DP),
intent(inout),
dimension(:) :: hnew
189 integer(I4B) :: nodes, nja
190 integer(I4B) :: ihfb, n, m
192 integer(I4B) :: idiag, isymcon
193 integer(I4B) :: ixt3d
194 real(DP) :: cond, condhfb, aterm
195 real(DP) :: fawidth, faheight
196 real(DP) :: topn, topm, botn, botm
197 real(DP) :: viscratio
201 nodes = this%dis%nodes
202 nja = this%dis%con%nja
203 if (
associated(this%xt3d%ixt3d))
then
204 ixt3d = this%xt3d%ixt3d
211 do ihfb = 1, this%nhfb
212 n = min(this%noden(ihfb), this%nodem(ihfb))
213 m = max(this%noden(ihfb), this%nodem(ihfb))
215 if (this%ibound(n) == 0 .or. this%ibound(m) == 0) cycle
217 if (this%ivsc /= 0)
then
218 call this%vsc%get_visc_ratio(n, m, hnew(n), hnew(m), viscratio)
221 if (this%hydchr(ihfb) >
dzero)
then
222 if (this%inewton == 0)
then
223 ipos = this%idxloc(ihfb)
228 if (this%icelltype(n) == 1)
then
229 if (hnew(n) < topn) topn = hnew(n)
231 if (this%icelltype(m) == 1)
then
232 if (hnew(m) < topm) topm = hnew(m)
234 if (this%ihc(this%jas(ipos)) == 2)
then
235 faheight = min(topn, topm) - max(botn, botm)
237 faheight =
dhalf * ((topn - botn) + (topm - botm))
239 fawidth = this%hwva(this%jas(ipos))
240 condhfb = this%hydchr(ihfb) * viscratio * &
243 condhfb = this%hydchr(ihfb) * viscratio
246 condhfb = this%hydchr(ihfb)
249 call this%xt3d%xt3d_fhfb(kiter, nodes, nja, matrix_sln, idxglo, &
250 rhs, hnew, n, m, condhfb)
256 if (this%inewton == 0)
then
257 do ihfb = 1, this%nhfb
258 ipos = this%idxloc(ihfb)
259 aterm = matrix_sln%get_value_pos(idxglo(ipos))
262 if (this%ibound(n) == 0 .or. this%ibound(m) == 0) cycle
264 if (this%ivsc /= 0)
then
265 call this%vsc%get_visc_ratio(n, m, hnew(n), hnew(m), viscratio)
268 if (this%icelltype(n) == 1 .or. this%icelltype(m) == 1 .or. &
276 if (this%icelltype(n) == 1)
then
277 if (hnew(n) < topn) topn = hnew(n)
279 if (this%icelltype(m) == 1)
then
280 if (hnew(m) < topm) topm = hnew(m)
282 if (this%ihc(this%jas(ipos)) == 2)
then
283 faheight = min(topn, topm) - max(botn, botm)
285 faheight =
dhalf * ((topn - botn) + (topm - botm))
287 if (this%hydchr(ihfb) >
dzero)
then
288 fawidth = this%hwva(this%jas(ipos))
289 condhfb = this%hydchr(ihfb) * viscratio * &
291 cond = aterm * condhfb / (aterm + condhfb)
293 cond = -aterm * this%hydchr(ihfb)
297 this%condsav(ihfb) = cond
301 call matrix_sln%add_value_pos(idxglo(idiag), aterm - cond)
302 call matrix_sln%set_value_pos(idxglo(ipos), cond)
305 isymcon = this%isym(ipos)
307 call matrix_sln%add_value_pos(idxglo(idiag), aterm - cond)
308 call matrix_sln%set_value_pos(idxglo(isymcon), cond)