MODFLOW 6  version 6.7.0.dev0
USGS Modular Hydrologic Model
virtualmodelmodule Module Reference

Data Types

interface  get_virtual_model
 
type  virtualmodeltype
 

Functions/Subroutines

subroutine vm_create (this, name, id, model)
 
subroutine init_virtual_data (this)
 
subroutine vm_prepare_stage (this, stage)
 
integer(i4b) function dis_get_nodeuser (this, node_reduced)
 Get user node number from reduced number. More...
 
subroutine dis_noder_to_string (this, node_reduced, node_str)
 
subroutine vm_destroy (this)
 
subroutine allocate_data (this)
 
subroutine deallocate_data (this)
 
class(virtualmodeltype) function, pointer, public get_virtual_model_from_list (model_list, idx)
 
class(virtualmodeltype) function, pointer, public cast_as_virtual_model (obj_ptr)
 
logical(lgp) function eq_virtual_model (this, v_model)
 
logical(lgp) function eq_numerical_model (this, num_model)
 
class(virtualmodeltype) function, pointer get_virtual_model_by_id (model_id)
 Returns a virtual model with the specified id. More...
 
class(virtualmodeltype) function, pointer get_virtual_model_by_name (model_name)
 Returns a virtual model with the specified name. More...
 

Function/Subroutine Documentation

◆ allocate_data()

subroutine virtualmodelmodule::allocate_data ( class(virtualmodeltype this)
private

Definition at line 233 of file VirtualModel.f90.

234  class(VirtualModelType) :: this
235 
236  allocate (this%con_ianglex)
237  allocate (this%con_ia)
238  allocate (this%con_ja)
239  allocate (this%con_jas)
240  allocate (this%con_ihc)
241  allocate (this%con_hwva)
242  allocate (this%con_cl1)
243  allocate (this%con_cl2)
244  allocate (this%con_anglex)
245  allocate (this%dis_ndim)
246  allocate (this%dis_nodes)
247  allocate (this%dis_nodesuser)
248  allocate (this%dis_nodeuser)
249  allocate (this%dis_nja)
250  allocate (this%dis_njas)
251  allocate (this%dis_xorigin)
252  allocate (this%dis_yorigin)
253  allocate (this%dis_angrot)
254  allocate (this%dis_xc)
255  allocate (this%dis_yc)
256  allocate (this%dis_top)
257  allocate (this%dis_bot)
258  allocate (this%dis_area)
259  allocate (this%moffset)
260  allocate (this%x)
261  allocate (this%x_old)
262  allocate (this%ibound)
263  allocate (this%idsoln)
264 

◆ cast_as_virtual_model()

class(virtualmodeltype) function, pointer, public virtualmodelmodule::cast_as_virtual_model ( class(*), pointer  obj_ptr)

Definition at line 316 of file VirtualModel.f90.

317  class(*), pointer :: obj_ptr
318  class(VirtualModelType), pointer :: v_model
319 
320  v_model => null()
321  select type (obj_ptr)
322  class is (virtualmodeltype)
323  v_model => obj_ptr
324  end select
325 
Here is the caller graph for this function:

◆ deallocate_data()

subroutine virtualmodelmodule::deallocate_data ( class(virtualmodeltype this)
private

Definition at line 267 of file VirtualModel.f90.

268  class(VirtualModelType) :: this
269 
270  ! CON
271  deallocate (this%con_ianglex)
272  deallocate (this%con_ia)
273  deallocate (this%con_ja)
274  deallocate (this%con_jas)
275  deallocate (this%con_ihc)
276  deallocate (this%con_hwva)
277  deallocate (this%con_cl1)
278  deallocate (this%con_cl2)
279  deallocate (this%con_anglex)
280  ! DIS
281  deallocate (this%dis_ndim)
282  deallocate (this%dis_nodes)
283  deallocate (this%dis_nodesuser)
284  deallocate (this%dis_nodeuser)
285  deallocate (this%dis_nja)
286  deallocate (this%dis_njas)
287  deallocate (this%dis_xorigin)
288  deallocate (this%dis_yorigin)
289  deallocate (this%dis_angrot)
290  deallocate (this%dis_xc)
291  deallocate (this%dis_yc)
292  deallocate (this%dis_top)
293  deallocate (this%dis_bot)
294  deallocate (this%dis_area)
295  ! Numerical model
296  deallocate (this%moffset)
297  deallocate (this%x)
298  deallocate (this%x_old)
299  deallocate (this%ibound)
300  ! Base model
301  deallocate (this%idsoln)
302 

◆ dis_get_nodeuser()

integer(i4b) function virtualmodelmodule::dis_get_nodeuser ( class(virtualmodeltype this,
integer(i4b), intent(in)  node_reduced 
)
private
Parameters
thisthis virtual model
[in]node_reducedthe reduced node number
Returns
the returned user node number

Definition at line 195 of file VirtualModel.f90.

196  class(VirtualModelType) :: this !< this virtual model
197  integer(I4B), intent(in) :: node_reduced !< the reduced node number
198  integer(I4B) :: node_user !< the returned user node number
199 
200  if (this%dis_nodes%get() < this%dis_nodesuser%get()) then
201  node_user = this%dis_nodeuser%get(node_reduced)
202  else
203  node_user = node_reduced
204  end if
205 

◆ dis_noder_to_string()

subroutine virtualmodelmodule::dis_noder_to_string ( class(virtualmodeltype this,
integer(i4b), intent(in)  node_reduced,
character(len=*), intent(inout)  node_str 
)
private
Parameters
thisthis virtual model
[in]node_reducedreduced node number
[in,out]node_strthe string representative of the user node number

Definition at line 208 of file VirtualModel.f90.

209  class(VirtualModelType) :: this !< this virtual model
210  integer(I4B), intent(in) :: node_reduced !< reduced node number
211  character(len=*), intent(inout) :: node_str !< the string representative of the user node number
212  ! local
213  character(len=11) :: nr_str
214 
215  if (this%is_local) then
216  call this%local_model%dis%noder_to_string(node_reduced, node_str)
217  else
218  ! for now this will look like: (102r)
219  write (nr_str, '(i0)') node_reduced
220  node_str = '('//trim(adjustl(nr_str))//'r)'
221  end if
222 

◆ eq_numerical_model()

logical(lgp) function virtualmodelmodule::eq_numerical_model ( class(virtualmodeltype), intent(in)  this,
class(numericalmodeltype), intent(in)  num_model 
)
private

Definition at line 337 of file VirtualModel.f90.

338  class(VirtualModelType), intent(in) :: this
339  class(NumericalModelType), intent(in) :: num_model
340  logical(LGP) :: is_equal
341 
342  is_equal = (this%id == num_model%id)
343 
Here is the caller graph for this function:

◆ eq_virtual_model()

logical(lgp) function virtualmodelmodule::eq_virtual_model ( class(virtualmodeltype), intent(in)  this,
class(virtualmodeltype), intent(in)  v_model 
)
private

Definition at line 328 of file VirtualModel.f90.

329  class(VirtualModelType), intent(in) :: this
330  class(VirtualModelType), intent(in) :: v_model
331  logical(LGP) :: is_equal
332 
333  is_equal = (this%id == v_model%id)
334 
Here is the caller graph for this function:

◆ get_virtual_model_by_id()

class(virtualmodeltype) function, pointer virtualmodelmodule::get_virtual_model_by_id ( integer(i4b)  model_id)
private

Definition at line 348 of file VirtualModel.f90.

350  integer(I4B) :: model_id
351  class(VirtualModelType), pointer :: virtual_model
352  ! local
353  integer(I4B) :: i
354  class(*), pointer :: vm
355 
356  virtual_model => null()
357  do i = 1, virtual_model_list%Count()
358  vm => virtual_model_list%GetItem(i)
359  select type (vm)
360  class is (virtualmodeltype)
361  if (vm%id == model_id) then
362  virtual_model => vm
363  return
364  end if
365  end select
366  end do
367 
type(listtype), public virtual_model_list

◆ get_virtual_model_by_name()

class(virtualmodeltype) function, pointer virtualmodelmodule::get_virtual_model_by_name ( character(len=*)  model_name)

Definition at line 372 of file VirtualModel.f90.

374  character(len=*) :: model_name
375  class(VirtualModelType), pointer :: virtual_model
376  ! local
377  integer(I4B) :: i
378  class(*), pointer :: vm
379 
380  virtual_model => null()
381  do i = 1, virtual_model_list%Count()
382  vm => virtual_model_list%GetItem(i)
383  select type (vm)
384  class is (virtualmodeltype)
385  if (vm%name == model_name) then
386  virtual_model => vm
387  return
388  end if
389  end select
390  end do
391 

◆ get_virtual_model_from_list()

class(virtualmodeltype) function, pointer, public virtualmodelmodule::get_virtual_model_from_list ( type(listtype model_list,
integer(i4b)  idx 
)

Definition at line 305 of file VirtualModel.f90.

306  type(ListType) :: model_list
307  integer(I4B) :: idx
308  class(VirtualModelType), pointer :: v_model
309  ! local
310  class(*), pointer :: obj_ptr
311 
312  obj_ptr => model_list%GetItem(idx)
313  v_model => cast_as_virtual_model(obj_ptr)
Here is the call graph for this function:
Here is the caller graph for this function:

◆ init_virtual_data()

subroutine virtualmodelmodule::init_virtual_data ( class(virtualmodeltype this)
private

Definition at line 93 of file VirtualModel.f90.

94  class(VirtualModelType) :: this
95 
96  ! CON
97  call this%set(this%con_ianglex%base(), 'IANGLEX', 'CON', map_all_type)
98  call this%set(this%con_ia%base(), 'IA', 'CON', map_all_type)
99  call this%set(this%con_ja%base(), 'JA', 'CON', map_all_type)
100  call this%set(this%con_jas%base(), 'JAS', 'CON', map_all_type)
101  call this%set(this%con_ihc%base(), 'IHC', 'CON', map_all_type)
102  call this%set(this%con_hwva%base(), 'HWVA', 'CON', map_all_type)
103  call this%set(this%con_cl1%base(), 'CL1', 'CON', map_all_type)
104  call this%set(this%con_cl2%base(), 'CL2', 'CON', map_all_type)
105  call this%set(this%con_anglex%base(), 'ANGLEX', 'CON', map_all_type)
106  ! DIS
107  call this%set(this%dis_ndim%base(), 'NDIM', 'DIS', map_all_type)
108  call this%set(this%dis_nodes%base(), 'NODES', 'DIS', map_all_type)
109  call this%set(this%dis_nodesuser%base(), 'NODESUSER', 'DIS', map_all_type)
110  call this%set(this%dis_nodeuser%base(), 'NODEUSER', 'DIS', map_all_type)
111  call this%set(this%dis_nja%base(), 'NJA', 'DIS', map_all_type)
112  call this%set(this%dis_njas%base(), 'NJAS', 'DIS', map_all_type)
113  call this%set(this%dis_xorigin%base(), 'XORIGIN', 'DIS', map_all_type)
114  call this%set(this%dis_yorigin%base(), 'YORIGIN', 'DIS', map_all_type)
115  call this%set(this%dis_angrot%base(), 'ANGROT', 'DIS', map_all_type)
116  call this%set(this%dis_xc%base(), 'XC', 'DIS', map_all_type)
117  call this%set(this%dis_yc%base(), 'YC', 'DIS', map_all_type)
118  call this%set(this%dis_top%base(), 'TOP', 'DIS', map_all_type)
119  call this%set(this%dis_bot%base(), 'BOT', 'DIS', map_all_type)
120  call this%set(this%dis_area%base(), 'AREA', 'DIS', map_all_type)
121  ! Numerical model
122  call this%set(this%moffset%base(), 'MOFFSET', '', map_all_type)
123  call this%set(this%x%base(), 'X', '', map_node_type)
124  call this%set(this%x_old%base(), 'XOLD', '', map_node_type)
125  call this%set(this%ibound%base(), 'IBOUND', '', map_node_type)
126  ! Base model
127  call this%set(this%idsoln%base(), 'IDSOLN', '', map_all_type)
128 

◆ vm_create()

subroutine virtualmodelmodule::vm_create ( class(virtualmodeltype this,
character(len=*)  name,
integer(i4b)  id,
class(numericalmodeltype), pointer  model 
)

Definition at line 75 of file VirtualModel.f90.

76  class(VirtualModelType) :: this
77  character(len=*) :: name
78  integer(I4B) :: id
79  class(NumericalModelType), pointer :: model
80  ! local
81  logical(LGP) :: is_local
82 
83  is_local = associated(model)
84  call this%VirtualDataContainerType%vdc_create(name, id, is_local)
85 
86  this%local_model => model
87 
88  call this%allocate_data()
89  call this%init_virtual_data()
90 

◆ vm_destroy()

subroutine virtualmodelmodule::vm_destroy ( class(virtualmodeltype this)
private

Definition at line 225 of file VirtualModel.f90.

226  class(VirtualModelType) :: this
227 
228  call this%VirtualDataContainerType%destroy()
229  call this%deallocate_data()
230 

◆ vm_prepare_stage()

subroutine virtualmodelmodule::vm_prepare_stage ( class(virtualmodeltype this,
integer(i4b)  stage 
)
private

Definition at line 131 of file VirtualModel.f90.

132  class(VirtualModelType) :: this
133  integer(I4B) :: stage
134  ! local
135  integer(I4B) :: nodes, nodesuser, nja, njas
136  logical(LGP) :: is_reduced
137 
138  if (stage == stg_aft_mdl_df) then
139 
140  call this%map(this%idsoln%base(), (/stg_aft_mdl_df/))
141  call this%map(this%con_ianglex%base(), (/stg_aft_mdl_df/))
142  call this%map(this%dis_ndim%base(), (/stg_aft_mdl_df/))
143  call this%map(this%dis_nodes%base(), (/stg_aft_mdl_df/))
144  call this%map(this%dis_nodesuser%base(), (/stg_aft_mdl_df/))
145  call this%map(this%dis_nja%base(), (/stg_aft_mdl_df/))
146  call this%map(this%dis_njas%base(), (/stg_aft_mdl_df/))
147 
148  else if (stage == stg_bfr_exg_ac) then
149 
150  nodes = this%dis_nodes%get()
151  nodesuser = this%dis_nodesuser%get()
152  is_reduced = (nodes /= nodesuser)
153  call this%map(this%moffset%base(), (/stg_bfr_exg_ac/))
154  if (is_reduced) then
155  call this%map(this%dis_nodeuser%base(), nodes, (/stg_bfr_exg_ac/))
156  else
157  ! no reduction, zero sized array, never synchronize
158  call this%map(this%dis_nodeuser%base(), 0, (/stg_never/))
159  end if
160 
161  else if (stage == stg_bfr_con_df) then
162 
163  nodes = this%dis_nodes%get()
164  nja = this%dis_nja%get()
165  njas = this%dis_njas%get()
166  ! DIS
167  call this%map(this%dis_xorigin%base(), (/stg_bfr_con_df/))
168  call this%map(this%dis_yorigin%base(), (/stg_bfr_con_df/))
169  call this%map(this%dis_angrot%base(), (/stg_bfr_con_df/))
170  call this%map(this%dis_xc%base(), nodes, (/stg_bfr_con_df/))
171  call this%map(this%dis_yc%base(), nodes, (/stg_bfr_con_df/))
172  call this%map(this%dis_top%base(), nodes, (/stg_bfr_con_df/))
173  call this%map(this%dis_bot%base(), nodes, (/stg_bfr_con_df/))
174  call this%map(this%dis_area%base(), nodes, (/stg_bfr_con_df/))
175  ! CON
176  call this%map(this%con_ia%base(), nodes + 1, (/stg_bfr_con_df/))
177  call this%map(this%con_ja%base(), nja, (/stg_bfr_con_df/))
178  call this%map(this%con_jas%base(), nja, (/stg_bfr_con_df/))
179  call this%map(this%con_ihc%base(), njas, (/stg_bfr_con_df/))
180  call this%map(this%con_hwva%base(), njas, (/stg_bfr_con_df/))
181  call this%map(this%con_cl1%base(), njas, (/stg_bfr_con_df/))
182  call this%map(this%con_cl2%base(), njas, (/stg_bfr_con_df/))
183  if (this%con_ianglex%get() > 0) then
184  call this%map(this%con_anglex%base(), njas, (/stg_bfr_con_df/))
185  else
186  call this%map(this%con_anglex%base(), 0, (/stg_never/))
187  end if
188 
189  end if
190