MODFLOW 6  version 6.7.0.dev1
USGS Modular Hydrologic Model
simulationcreatemodule Module Reference

Functions/Subroutines

subroutine, public simulation_cr ()
 Read the simulation name file and initialize the models, exchanges. More...
 
subroutine, public simulation_da ()
 Deallocate simulation variables. More...
 
subroutine source_simulation_nam ()
 Source the simulation name file. More...
 
subroutine options_create ()
 Set the simulation options. More...
 
subroutine timing_create ()
 Set the timing module to be used for the simulation. More...
 
subroutine models_create ()
 Set the models to be used for the simulation. More...
 
subroutine exchanges_create ()
 Set the exchanges to be used for the simulation. More...
 
subroutine solution_group_check (sgp, sgid, isgpsoln)
 Check a solution_group to be used for the simulation. More...
 
subroutine solution_groups_create ()
 Set the solution_groups to be used for the simulation. More...
 
subroutine check_model_assignment ()
 Check for dangling models, and break with error when found. More...
 
subroutine assign_exchanges ()
 Assign exchanges to solutions. More...
 
subroutine check_model_name (mtype, mname)
 Check that the model name is valid. More...
 

Function/Subroutine Documentation

◆ assign_exchanges()

subroutine simulationcreatemodule::assign_exchanges
private

This assigns NumericalExchanges to NumericalSolutions, based on the link between the models in the solution and those exchanges. The BaseExchangeconnects_model() function should be overridden to indicate if such a link exists.

Definition at line 760 of file SimulationCreate.f90.

761  ! -- local
762  class(BaseSolutionType), pointer :: sp
763  class(BaseExchangeType), pointer :: ep
764  class(BaseModelType), pointer :: mp
765  type(ListType), pointer :: models_in_solution
766  integer(I4B) :: is, ie, im
767 
768  do is = 1, basesolutionlist%Count()
769  sp => getbasesolutionfromlist(basesolutionlist, is)
770  !
771  ! -- now loop over exchanges
772  do ie = 1, baseexchangelist%Count()
773  ep => getbaseexchangefromlist(baseexchangelist, ie)
774  !
775  ! -- and add when it affects (any model in) the solution matrix
776  models_in_solution => sp%get_models()
777  do im = 1, models_in_solution%Count()
778  mp => getbasemodelfromlist(models_in_solution, im)
779  if (ep%connects_model(mp)) then
780  !
781  ! -- add to solution (and only once)
782  call sp%add_exchange(ep)
783  exit
784  end if
785  end do
786  end do
787  end do
Here is the call graph for this function:
Here is the caller graph for this function:

◆ check_model_assignment()

subroutine simulationcreatemodule::check_model_assignment

Definition at line 734 of file SimulationCreate.f90.

735  character(len=LINELENGTH) :: errmsg
736  class(BaseModelType), pointer :: mp
737  integer(I4B) :: im
738 
739  do im = 1, basemodellist%Count()
740  mp => getbasemodelfromlist(basemodellist, im)
741  if (mp%idsoln == 0) then
742  write (errmsg, '(a,a)') &
743  'Model was not assigned to a solution: ', mp%name
744  call store_error(errmsg)
745  end if
746  end do
747  if (count_errors() > 0) then
748  call store_error_filename('mfsim.nam')
749  end if
750 
Here is the call graph for this function:
Here is the caller graph for this function:

◆ check_model_name()

subroutine simulationcreatemodule::check_model_name ( character(len=*), intent(in)  mtype,
character(len=*), intent(inout)  mname 
)
private

Definition at line 792 of file SimulationCreate.f90.

793  ! -- dummy
794  character(len=*), intent(in) :: mtype
795  character(len=*), intent(inout) :: mname
796  ! -- local
797  integer :: ilen
798  integer :: i
799  character(len=LINELENGTH) :: errmsg
800  logical :: terminate = .true.
801 
802  ilen = len_trim(mname)
803  if (ilen > lenmodelname) then
804  write (errmsg, '(a,a)') 'Invalid model name: ', trim(mname)
805  call store_error(errmsg)
806  write (errmsg, '(a,i0,a,i0)') &
807  'Name length of ', ilen, ' exceeds maximum length of ', &
808  lenmodelname
809  call store_error(errmsg, terminate)
810  end if
811  do i = 1, ilen
812  if (mname(i:i) == ' ') then
813  write (errmsg, '(a,a)') 'Invalid model name: ', trim(mname)
814  call store_error(errmsg)
815  write (errmsg, '(a)') &
816  'Model name cannot have spaces within it.'
817  call store_error(errmsg, terminate)
818  end if
819  end do
Here is the call graph for this function:
Here is the caller graph for this function:

◆ exchanges_create()

subroutine simulationcreatemodule::exchanges_create

Definition at line 353 of file SimulationCreate.f90.

354  ! -- modules
369  ! use VirtualPrtExchangeModule, only: add_virtual_prt_exchange
370  ! -- dummy
371  ! -- locals
372  character(len=LENMEMPATH) :: input_mempath
373  type(CharacterStringType), dimension(:), contiguous, &
374  pointer :: etypes !< exg types
375  type(CharacterStringType), dimension(:), contiguous, &
376  pointer :: efiles !< exg file names
377  type(CharacterStringType), dimension(:), contiguous, &
378  pointer :: emnames_a !< model a names
379  type(CharacterStringType), dimension(:), contiguous, &
380  pointer :: emnames_b !< model b names
381  type(CharacterStringType), dimension(:), contiguous, &
382  pointer :: emempaths
383  character(len=LINELENGTH) :: exgtype
384  integer(I4B) :: exg_id
385  integer(I4B) :: m1_id, m2_id
386  character(len=LINELENGTH) :: fname, name1, name2
387  character(len=LENEXCHANGENAME) :: exg_name
388  character(len=LENMEMPATH) :: exg_mempath
389  integer(I4B) :: n
390  character(len=LINELENGTH) :: errmsg
391  logical(LGP) :: terminate = .true.
392  logical(LGP) :: both_remote, both_local
393  ! -- formats
394  character(len=*), parameter :: fmtmerr = "('Error in simulation control ', &
395  &'file. Could not find model: ', a)"
396  !
397  ! -- set input memory path
398  input_mempath = create_mem_path('SIM', 'NAM', idm_context)
399  !
400  ! -- set pointers to input context exchange attribute arrays
401  call mem_setptr(etypes, 'EXGTYPE', input_mempath)
402  call mem_setptr(efiles, 'EXGFILE', input_mempath)
403  call mem_setptr(emnames_a, 'EXGMNAMEA', input_mempath)
404  call mem_setptr(emnames_b, 'EXGMNAMEB', input_mempath)
405  call mem_setptr(emempaths, 'EXGMEMPATHS', input_mempath)
406  !
407  ! -- open exchange logging block
408  write (iout, '(/1x,a)') 'READING SIMULATION EXCHANGES'
409  !
410  ! -- initialize
411  exg_id = 0
412  !
413  ! -- create exchanges
414  do n = 1, size(etypes)
415  !
416  ! -- attributes for this exchange
417  exgtype = etypes(n)
418  fname = efiles(n)
419  name1 = emnames_a(n)
420  name2 = emnames_b(n)
421  exg_mempath = emempaths(n)
422 
423  exg_id = exg_id + 1
424 
425  ! find model index in list
426  m1_id = ifind(model_names, name1)
427  if (m1_id < 0) then
428  write (errmsg, fmtmerr) trim(name1)
429  call store_error(errmsg, terminate)
430  end if
431  m2_id = ifind(model_names, name2)
432  if (m2_id < 0) then
433  write (errmsg, fmtmerr) trim(name2)
434  call store_error(errmsg, terminate)
435  end if
436 
437  ! both models on other process? then don't create it here...
438  both_remote = (model_loc_idx(m1_id) == -1 .and. &
439  model_loc_idx(m2_id) == -1)
440  both_local = (model_loc_idx(m1_id) > 0 .and. &
441  model_loc_idx(m2_id) > 0)
442  if (.not. both_remote) then
443  write (iout, '(4x,a,a,i0,a,i0,a,i0)') trim(exgtype), ' exchange ', &
444  exg_id, ' will be created to connect model ', m1_id, &
445  ' with model ', m2_id
446  end if
447 
448  select case (exgtype)
449  case ('CHF6-GWF6')
450  write (exg_name, '(a,i0)') 'CHF-GWF_', exg_id
451  if (both_local) then
452  call chfgwf_cr(fname, exg_name, exg_id, m1_id, m2_id, exg_mempath)
453  end if
454  case ('GWF6-GWF6')
455  write (exg_name, '(a,i0)') 'GWF-GWF_', exg_id
456  if (.not. both_remote) then
457  call gwfexchange_create(fname, exg_name, exg_id, m1_id, m2_id, &
458  exg_mempath)
459  end if
460  call add_virtual_gwf_exchange(exg_name, exg_id, m1_id, m2_id)
461  case ('GWF6-GWT6')
462  if (both_local) then
463  call gwfgwt_cr(fname, exg_id, m1_id, m2_id)
464  end if
465  case ('GWF6-GWE6')
466  if (both_local) then
467  call gwfgwe_cr(fname, exg_id, m1_id, m2_id)
468  end if
469  case ('GWF6-PRT6')
470  call gwfprt_cr(fname, exg_id, m1_id, m2_id)
471  case ('GWT6-GWT6')
472  write (exg_name, '(a,i0)') 'GWT-GWT_', exg_id
473  if (.not. both_remote) then
474  call gwtexchange_create(fname, exg_name, exg_id, m1_id, m2_id, &
475  exg_mempath)
476  end if
477  call add_virtual_gwt_exchange(exg_name, exg_id, m1_id, m2_id)
478  case ('GWE6-GWE6')
479  write (exg_name, '(a,i0)') 'GWE-GWE_', exg_id
480  if (.not. both_remote) then
481  call gweexchange_create(fname, exg_name, exg_id, m1_id, m2_id, &
482  exg_mempath)
483  end if
484  call add_virtual_gwe_exchange(exg_name, exg_id, m1_id, m2_id)
485  case ('OLF6-GWF6')
486  write (exg_name, '(a,i0)') 'OLF-GWF_', exg_id
487  if (both_local) then
488  call olfgwf_cr(fname, exg_name, exg_id, m1_id, m2_id, exg_mempath)
489  end if
490  case default
491  write (errmsg, '(a,a)') &
492  'Unknown simulation exchange type: ', trim(exgtype)
493  call store_error(errmsg, terminate)
494  end select
495  end do
496  !
497  ! -- close exchange logging block
498  write (iout, '(1x,a)') 'END OF SIMULATION EXCHANGES'
This module contains the ChfGwfExchangeModule Module.
Definition: exg-chfgwf.f90:6
subroutine, public chfgwf_cr(filename, name, id, m1_id, m2_id, input_mempath)
@ brief Create CHF GWF exchange
Definition: exg-chfgwf.f90:41
This module contains the GweGweExchangeModule Module.
Definition: exg-gwegwe.f90:10
subroutine, public gweexchange_create(filename, name, id, m1_id, m2_id, input_mempath)
@ brief Create GWT GWT exchange
Definition: exg-gwegwe.f90:111
subroutine, public gwfgwe_cr(filename, id, m1_id, m2_id)
Create a new GWF to GWE exchange object.
Definition: exg-gwfgwe.f90:47
This module contains the GwfGwfExchangeModule Module.
Definition: exg-gwfgwf.f90:10
subroutine, public gwfexchange_create(filename, name, id, m1_id, m2_id, input_mempath)
@ brief Create GWF GWF exchange
Definition: exg-gwfgwf.f90:122
subroutine, public gwfgwt_cr(filename, id, m1_id, m2_id)
Create a new GWF to GWT exchange object.
Definition: exg-gwfgwt.f90:49
subroutine, public gwfprt_cr(filename, id, m1id, m2id)
Create a new GWF to PRT exchange object.
Definition: exg-gwfprt.f90:40
This module contains the GwtGwtExchangeModule Module.
Definition: exg-gwtgwt.f90:10
subroutine, public gwtexchange_create(filename, name, id, m1_id, m2_id, input_mempath)
@ brief Create GWT GWT exchange
Definition: exg-gwtgwt.f90:110
character(len=lenmempath) function create_mem_path(component, subcomponent, context)
returns the path to the memory object
This module contains the OlfGwfExchangeModule Module.
Definition: exg-olfgwf.f90:6
subroutine, public olfgwf_cr(filename, name, id, m1_id, m2_id, input_mempath)
@ brief Create OLF GWF exchange
Definition: exg-olfgwf.f90:41
This module contains simulation variables.
Definition: SimVariables.f90:9
character(len=linelength) idm_context
subroutine, public add_virtual_gwe_exchange(name, exchange_id, model1_id, model2_id)
Add a virtual GWE-GWE exchange to the simulation.
subroutine, public add_virtual_gwf_exchange(name, exchange_id, model1_id, model2_id)
Add a virtual GWF-GWF exchange to the simulation.
subroutine, public add_virtual_gwt_exchange(name, exchange_id, model1_id, model2_id)
Add a virtual GWT-GWT exchange to the simulation.
Here is the call graph for this function:
Here is the caller graph for this function:

◆ models_create()

subroutine simulationcreatemodule::models_create

Definition at line 200 of file SimulationCreate.f90.

201  ! -- modules
206  use chfmodule, only: chf_cr
207  use gwfmodule, only: gwf_cr
208  use gwtmodule, only: gwt_cr
209  use gwemodule, only: gwe_cr
210  use olfmodule, only: olf_cr
211  use prtmodule, only: prt_cr
216  ! use VirtualPrtModelModule, only: add_virtual_prt_model
217  use constantsmodule, only: lenmodelname
218  ! -- dummy
219  ! -- locals
220  type(DistributedSimType), pointer :: ds
221  character(len=LENMEMPATH) :: input_mempath
222  type(CharacterStringType), dimension(:), contiguous, &
223  pointer :: mtypes !< model types
224  type(CharacterStringType), dimension(:), contiguous, &
225  pointer :: mfnames !< model file names
226  type(CharacterStringType), dimension(:), contiguous, &
227  pointer :: mnames !< model names
228  integer(I4B) :: im
229  class(NumericalModelType), pointer :: num_model
230  character(len=LINELENGTH) :: model_type
231  character(len=LINELENGTH) :: fname, model_name
232  integer(I4B) :: n, nr_models_glob
233  integer(I4B), dimension(:), pointer :: model_ranks => null()
234  logical :: terminate = .true.
235  !
236  ! -- set input memory path
237  input_mempath = create_mem_path('SIM', 'NAM', idm_context)
238  !
239  ! -- set pointers to input context model attribute arrays
240  call mem_setptr(mtypes, 'MTYPE', input_mempath)
241  call mem_setptr(mfnames, 'MFNAME', input_mempath)
242  call mem_setptr(mnames, 'MNAME', input_mempath)
243  !
244  ! -- allocate global arrays
245  nr_models_glob = size(mnames)
246  allocate (model_names(nr_models_glob))
247  allocate (model_loc_idx(nr_models_glob))
248  !
249  ! -- get model-to-cpu assignment (in serial all to rank 0)
250  ds => get_dsim()
251  model_ranks => ds%get_load_balance()
252  !
253  ! -- open model logging block
254  write (iout, '(/1x,a)') 'READING SIMULATION MODELS'
255  !
256  ! -- create models
257  im = 0
258  do n = 1, size(mtypes)
259  !
260  ! -- attributes for this model
261  model_type = mtypes(n)
262  fname = mfnames(n)
263  model_name = mnames(n)
264  !
265  call check_model_name(model_type, model_name)
266  !
267  ! increment global model id
268  model_names(n) = model_name(1:lenmodelname)
269  model_loc_idx(n) = -1
270  num_model => null()
271  !
272  ! -- add a new (local or global) model
273  select case (model_type)
274  case ('GWF6')
275  if (model_ranks(n) == proc_id) then
276  im = im + 1
277  write (iout, '(4x,2a,i0,a)') trim(model_type), ' model ', &
278  n, ' will be created'
279  call gwf_cr(fname, n, model_names(n))
280  num_model => getnumericalmodelfromlist(basemodellist, im)
281  model_loc_idx(n) = im
282  end if
283  call add_virtual_gwf_model(n, model_names(n), num_model)
284  case ('GWT6')
285  if (model_ranks(n) == proc_id) then
286  im = im + 1
287  write (iout, '(4x,2a,i0,a)') trim(model_type), ' model ', &
288  n, ' will be created'
289  call gwt_cr(fname, n, model_names(n))
290  num_model => getnumericalmodelfromlist(basemodellist, im)
291  model_loc_idx(n) = im
292  end if
293  call add_virtual_gwt_model(n, model_names(n), num_model)
294  case ('GWE6')
295  if (model_ranks(n) == proc_id) then
296  im = im + 1
297  write (iout, '(4x,2a,i0,a)') trim(model_type), ' model ', &
298  n, ' will be created'
299  call gwe_cr(fname, n, model_names(n))
300  num_model => getnumericalmodelfromlist(basemodellist, im)
301  model_loc_idx(n) = im
302  end if
303  call add_virtual_gwe_model(n, model_names(n), num_model)
304  case ('CHF6')
305  if (model_ranks(n) == proc_id) then
306  im = im + 1
307  write (iout, '(4x,2a,i0,a)') trim(model_type), " model ", &
308  n, " will be created"
309  call chf_cr(fname, n, model_names(n))
310  call dev_feature('CHF is still under development, install the &
311  &nightly build or compile from source with IDEVELOPMODE = 1.')
312  num_model => getnumericalmodelfromlist(basemodellist, im)
313  model_loc_idx(n) = im
314  end if
315  case ('OLF6')
316  if (model_ranks(n) == proc_id) then
317  im = im + 1
318  write (iout, '(4x,2a,i0,a)') trim(model_type), " model ", &
319  n, " will be created"
320  call olf_cr(fname, n, model_names(n))
321  call dev_feature('OLF is still under development, install the &
322  &nightly build or compile from source with IDEVELOPMODE = 1.')
323  num_model => getnumericalmodelfromlist(basemodellist, im)
324  model_loc_idx(n) = im
325  end if
326  case ('PRT6')
327  im = im + 1
328  write (iout, '(4x,2a,i0,a)') trim(model_type), ' model ', &
329  n, ' will be created'
330  call prt_cr(fname, n, model_names(n))
331  num_model => getnumericalmodelfromlist(basemodellist, im)
332  model_loc_idx(n) = im
333  case default
334  write (errmsg, '(a,a)') &
335  'Unknown simulation model type: ', trim(model_type)
336  call store_error(errmsg, terminate)
337  end select
338  end do
339  !
340  ! -- close model logging block
341  write (iout, '(1x,a)') 'END OF SIMULATION MODELS'
342  !
343  ! -- sanity check
344  if (simulation_mode == 'PARALLEL' .and. im == 0) then
345  write (errmsg, '(a, i0)') &
346  'No MODELS assigned to process ', proc_id
347  call store_error(errmsg, terminate)
348  end if
Channel Flow (CHF) Module.
Definition: chf.f90:3
subroutine, public chf_cr(filename, id, modelname)
Create a new surface water flow model object.
Definition: chf.f90:56
This module contains simulation constants.
Definition: Constants.f90:9
integer(i4b), parameter lenmodelname
maximum length of the model name
Definition: Constants.f90:22
class(distributedsimtype) function, pointer, public get_dsim()
Get pointer to the distributed simulation object.
Definition: gwe.f90:3
subroutine, public gwe_cr(filename, id, modelname)
Create a new groundwater energy transport model object.
Definition: gwe.f90:98
Definition: gwf.f90:1
subroutine, public gwf_cr(filename, id, modelname)
Create a new groundwater flow model object.
Definition: gwf.f90:138
Definition: gwt.f90:8
subroutine, public gwt_cr(filename, id, modelname)
Create a new groundwater transport model object.
Definition: gwt.f90:101
class(numericalmodeltype) function, pointer, public getnumericalmodelfromlist(list, idx)
Channel Flow (OLF) Module.
Definition: olf.f90:3
subroutine, public olf_cr(filename, id, modelname)
Create a new overland flow model object.
Definition: olf.f90:56
Definition: prt.f90:1
subroutine, public prt_cr(filename, id, modelname)
Create a new particle tracking model object.
Definition: prt.f90:120
character(len=maxcharlen) errmsg
error message string
subroutine, public add_virtual_gwe_model(model_id, model_name, model)
subroutine, public add_virtual_gwf_model(model_id, model_name, model)
Add virtual GWF model.
subroutine, public add_virtual_gwt_model(model_id, model_name, model)
Here is the call graph for this function:
Here is the caller graph for this function:

◆ options_create()

subroutine simulationcreatemodule::options_create
private

Definition at line 95 of file SimulationCreate.f90.

96  ! -- modules
101  use profilermodule, only: g_prof
103  ! -- dummy
104  ! -- locals
105  character(len=LENMEMPATH) :: input_mempath
106  integer(I4B), pointer :: simcontinue, nocheck, maxerror
107  character(len=:), pointer :: prmem, prprof
108  character(len=LINELENGTH) :: errmsg
109  !
110  ! -- set input memory path
111  input_mempath = create_mem_path('SIM', 'NAM', idm_context)
112  !
113  ! -- set pointers to input context option params
114  call mem_setptr(simcontinue, 'CONTINUE', input_mempath)
115  call mem_setptr(nocheck, 'NOCHECK', input_mempath)
116  call mem_setptr(prmem, 'PRMEM', input_mempath)
117  call mem_setptr(prprof, 'PRPROF', input_mempath)
118  call mem_setptr(maxerror, 'MAXERRORS', input_mempath)
119  !
120  ! -- update sim options
121  isimcontinue = simcontinue
122  isimcheck = nocheck
123  call maxerrors(maxerror)
124  !
125  if (prmem /= '') then
126  errmsg = ''
127  call mem_set_print_option(iout, prmem, errmsg)
128  if (errmsg /= '') then
129  call store_error(errmsg, .true.)
130  end if
131  end if
132 
133  ! set profiler print option
134  call g_prof%set_print_option(prprof)
135 
136  !
137  ! -- log values to list file
138  if (iout > 0) then
139  write (iout, '(/1x,a)') 'READING SIMULATION OPTIONS'
140  !
141  if (isimcontinue == 1) then
142  write (iout, '(4x, a)') &
143  'SIMULATION WILL CONTINUE EVEN IF THERE IS NONCONVERGENCE.'
144  end if
145  !
146  if (isimcheck == 0) then
147  write (iout, '(4x, a)') &
148  'MODEL DATA WILL NOT BE CHECKED FOR ERRORS.'
149  end if
150  !
151  write (iout, '(4x, a, i0)') &
152  'MAXIMUM NUMBER OF ERRORS THAT WILL BE STORED IS ', maxerror
153  !
154  if (prmem /= '') then
155  write (iout, '(4x, a, a, a)') &
156  'MEMORY_PRINT_OPTION SET TO "', trim(prmem), '".'
157  end if
158  !
159  write (iout, '(1x,a)') 'END OF SIMULATION OPTIONS'
160  end if
subroutine, public mem_set_print_option(iout, keyword, error_msg)
Set the memory print option.
type(profilertype), public g_prof
the global timer object (to reduce trivial lines of code)
Definition: Profiler.f90:65
integer(i4b) isimcheck
simulation input check flag (1) to check input, (0) to ignore checks
integer(i4b) isimcontinue
simulation continue flag (1) to continue if isimcnvg = 0, (0) to terminate
Here is the call graph for this function:
Here is the caller graph for this function:

◆ simulation_cr()

subroutine, public simulationcreatemodule::simulation_cr

Definition at line 34 of file SimulationCreate.f90.

35  ! -- modules
36  ! -- local
37  !
38  ! -- Source simulation nam input context and create objects
39  call source_simulation_nam()
Here is the call graph for this function:
Here is the caller graph for this function:

◆ simulation_da()

subroutine, public simulationcreatemodule::simulation_da

Definition at line 44 of file SimulationCreate.f90.

45  ! -- modules
48  ! -- local
49  type(DistributedSimType), pointer :: ds
50 
51  ! -- variables
52  ds => get_dsim()
53  call ds%destroy()
54  !
55  deallocate (model_names)
56  deallocate (model_loc_idx)
Here is the call graph for this function:
Here is the caller graph for this function:

◆ solution_group_check()

subroutine simulationcreatemodule::solution_group_check ( type(solutiongrouptype), intent(inout), pointer  sgp,
integer(i4b), intent(in)  sgid,
integer(i4b), intent(in)  isgpsoln 
)

Definition at line 503 of file SimulationCreate.f90.

504  ! -- modules
505  ! -- dummy
506  type(SolutionGroupType), pointer, intent(inout) :: sgp
507  integer(I4B), intent(in) :: sgid
508  integer(I4B), intent(in) :: isgpsoln
509  ! -- local
510  character(len=LINELENGTH) :: errmsg
511  logical :: terminate = .true.
512  ! -- formats
513  character(len=*), parameter :: fmterrmxiter = &
514  "('MXITER is set to ', i0, ' but there is only one solution', &
515  &' in SOLUTION GROUP ', i0, '. Set MXITER to 1 in simulation control', &
516  &' file.')"
517  !
518  ! -- error check completed group
519  if (sgid > 0) then
520  !
521  ! -- Make sure there is a solution in this solution group
522  if (isgpsoln == 0) then
523  write (errmsg, '(a,i0)') &
524  'There are no solutions for solution group ', sgid
525  call store_error(errmsg, terminate)
526  end if
527  !
528  ! -- If there is only one solution then mxiter should be 1.
529  if (isgpsoln == 1 .and. sgp%mxiter > 1) then
530  write (errmsg, fmterrmxiter) sgp%mxiter, isgpsoln
531  call store_error(errmsg, terminate)
532  end if
533  end if
Here is the call graph for this function:
Here is the caller graph for this function:

◆ solution_groups_create()

subroutine simulationcreatemodule::solution_groups_create
private

Definition at line 538 of file SimulationCreate.f90.

539  ! -- modules
547  use basemodelmodule, only: basemodeltype
550  ! -- dummy
551  ! -- local
552  character(len=LENMEMPATH) :: input_mempath
553  type(CharacterStringType), dimension(:), contiguous, &
554  pointer :: slntype
555  type(CharacterStringType), dimension(:), contiguous, &
556  pointer :: slnfname
557  type(CharacterStringType), dimension(:), contiguous, &
558  pointer :: slnmnames
559  integer(I4B), dimension(:), contiguous, pointer :: blocknum
560  character(len=LINELENGTH) :: stype, fname
561  character(len=:), allocatable :: mnames
562  type(SolutionGroupType), pointer :: sgp
563  class(BaseSolutionType), pointer :: sp
564  class(BaseModelType), pointer :: mp
565  integer(I4B) :: isoln
566  integer(I4B) :: isgpsoln
567  integer(I4B) :: sgid
568  integer(I4B) :: glo_mid
569  integer(I4B) :: loc_idx
570  integer(I4B) :: i, j, istat, mxiter
571  integer(I4B) :: nwords
572  character(len=LENMODELNAME), dimension(:), allocatable :: words
573  character(len=:), allocatable :: parse_str
574  character(len=LINELENGTH) :: errmsg
575  logical :: terminate = .true.
576  !
577  ! -- set memory path
578  input_mempath = create_mem_path('SIM', 'NAM', idm_context)
579  !
580  ! -- set pointers to input context solution attribute arrays
581  call mem_setptr(slntype, 'SLNTYPE', input_mempath)
582  call mem_setptr(slnfname, 'SLNFNAME', input_mempath)
583  call mem_setptr(slnmnames, 'SLNMNAMES', input_mempath)
584  call mem_setptr(blocknum, 'SOLUTIONGROUPNUM', input_mempath)
585  !
586  ! -- open solution group logging block
587  write (iout, '(/1x,a)') 'READING SOLUTIONGROUP'
588  !
589  ! -- initialize
590  sgid = 0 ! integer id of soln group, tracks with blocknum
591  isoln = 0 ! cumulative solution number
592  !
593  ! -- create solution groups
594  do i = 1, size(blocknum)
595  !
596  ! -- allocate slnmnames string
597  allocate (character(slnmnames(i)%strlen()) :: mnames)
598  !
599  ! -- attributes for this solution
600  stype = slntype(i)
601  fname = slnfname(i)
602  mnames = slnmnames(i)
603 
604  if (blocknum(i) /= sgid) then
605  !
606  ! -- check for new soln group
607  if (blocknum(i) == sgid + 1) then
608  !
609  ! -- error check completed group
610  call solution_group_check(sgp, sgid, isgpsoln)
611  !
612  ! -- reinitialize
613  nullify (sgp)
614  isgpsoln = 0 ! solution counter for this solution group
615  !
616  ! -- set sgid
617  sgid = blocknum(i)
618  !
619  ! -- create new soln group and add to global list
620  call solutiongroup_create(sgp, sgid)
621  call addsolutiongrouptolist(solutiongrouplist, sgp)
622  else
623  write (errmsg, '(a,i0,a,i0,a)') &
624  'Solution group blocks are not listed consecutively. Found ', &
625  blocknum(i), ' when looking for ', sgid + 1, '.'
626  call store_error(errmsg, terminate)
627  end if
628  end if
629  !
630  ! --
631  select case (stype)
632  !
633  case ('MXITER')
634  read (fname, *, iostat=istat) mxiter
635  if (istat == 0) then
636  sgp%mxiter = mxiter
637  end if
638  case ('IMS6')
639  !
640  ! -- increment solution counters
641  isoln = isoln + 1
642  isgpsoln = isgpsoln + 1
643  !
644  ! -- create soln and add to group
645  sp => create_ims_solution(simulation_mode, fname, isoln)
646  call sgp%add_solution(isoln, sp)
647  !
648  ! -- parse model names
649  parse_str = trim(mnames)//' '
650  call parseline(parse_str, nwords, words)
651  !
652  ! -- Find each model id and get model
653  do j = 1, nwords
654  call upcase(words(j))
655  glo_mid = ifind(model_names, words(j))
656  if (glo_mid == -1) then
657  write (errmsg, '(a,a)') 'Invalid model name: ', trim(words(j))
658  call store_error(errmsg, terminate)
659  end if
660  !
661  loc_idx = model_loc_idx(glo_mid)
662  if (loc_idx == -1) then
663  if (simulation_mode == 'PARALLEL') then
664  ! this is still ok
665  cycle
666  end if
667  end if
668  !
669  mp => getbasemodelfromlist(basemodellist, loc_idx)
670  !
671  ! -- Add the model to the solution
672  call sp%add_model(mp)
673  mp%idsoln = isoln
674  end do
675  case ('EMS6')
676  !
677  ! -- increment solution counters
678  isoln = isoln + 1
679  isgpsoln = isgpsoln + 1
680  !
681  ! -- create soln and add to group
682  sp => create_ems_solution(simulation_mode, fname, isoln)
683  call sgp%add_solution(isoln, sp)
684  !
685  ! -- parse model names
686  parse_str = trim(mnames)//' '
687  call parseline(parse_str, nwords, words)
688  !
689  ! -- Find each model id and get model
690  do j = 1, nwords
691  call upcase(words(j))
692  glo_mid = ifind(model_names, words(j))
693  if (glo_mid == -1) then
694  write (errmsg, '(a,a)') 'Invalid model name: ', trim(words(j))
695  call store_error(errmsg, terminate)
696  end if
697  !
698  loc_idx = model_loc_idx(glo_mid)
699  if (loc_idx == -1) then
700  if (simulation_mode == 'PARALLEL') then
701  ! this is still ok
702  cycle
703  end if
704  end if
705  !
706  mp => getbasemodelfromlist(basemodellist, loc_idx)
707  !
708  ! -- Add the model to the solution
709  call sp%add_model(mp)
710  mp%idsoln = isoln
711  end do
712  case default
713  end select
714  !
715  ! -- clean up
716  deallocate (mnames)
717  end do
718  !
719  ! -- error check final group
720  call solution_group_check(sgp, sgid, isgpsoln)
721  !
722  ! -- close exchange logging block
723  write (iout, '(1x,a)') 'END OF SOLUTIONGROUP'
724  !
725  ! -- Check and make sure at least one solution group was found
726  if (solutiongrouplist%Count() == 0) then
727  call store_error('There are no solution groups.', terminate)
728  end if
subroutine, public parseline(line, nwords, words, inunit, filename)
Parse a line into words.
subroutine, public upcase(word)
Convert to upper case.
character(len=linelength) simulation_mode
class(basesolutiontype) function, pointer, public create_ims_solution(sim_mode, filename, sol_id)
Create an IMS solution of type NumericalSolution for serial runs or its sub-type ParallelSolution for...
class(basesolutiontype) function, pointer, public create_ems_solution(sim_mode, filename, sol_id)
Create an EMS solution of type ExplicitSolution for serial runs or its sub-type ParallelSolution for.
subroutine, public solutiongroup_create(sgp, id)
Create a new solution group.
Highest level model type. All models extend this parent type.
Definition: BaseModel.f90:13
Here is the call graph for this function:
Here is the caller graph for this function:

◆ source_simulation_nam()

subroutine simulationcreatemodule::source_simulation_nam

Source from the simulation nam input context to initialize the models, exchanges, solutions, solutions groups. Then add the exchanges to the appropriate solutions.

Definition at line 66 of file SimulationCreate.f90.

67  ! -- dummy
68  ! -- local
69  !
70  ! -- Process OPTIONS block in namfile
71  call options_create()
72  !
73  ! -- Process TIMING block in namfile
74  call timing_create()
75  !
76  ! -- Process MODELS block in namfile
77  call models_create()
78  !
79  ! -- Process EXCHANGES block in namfile
80  call exchanges_create()
81  !
82  ! -- Process SOLUTION_GROUPS blocks in namfile
83  call solution_groups_create()
84  !
85  ! -- Go through each model and make sure that it has been assigned to
86  ! a solution.
87  call check_model_assignment()
88  !
89  ! -- Go through each solution and assign exchanges accordingly
90  call assign_exchanges()
Here is the call graph for this function:
Here is the caller graph for this function:

◆ timing_create()

subroutine simulationcreatemodule::timing_create

Definition at line 165 of file SimulationCreate.f90.

166  ! -- modules
170  use tdismodule, only: tdis_cr
171  ! -- dummy
172  ! -- locals
173  character(len=LENMEMPATH) :: input_mempath
174  character(len=LENMEMPATH) :: tdis_input_mempath
175  character(len=:), pointer :: tdis6
176  logical :: terminate = .true.
177  !
178  ! -- set input memory path
179  input_mempath = create_mem_path('SIM', 'NAM', idm_context)
180  tdis_input_mempath = create_mem_path('SIM', 'TDIS', idm_context)
181  !
182  write (iout, '(/1x,a)') 'READING SIMULATION TIMING'
183  !
184  ! -- set pointers to input context timing params
185  call mem_setptr(tdis6, 'TDIS6', input_mempath)
186  !
187  ! -- create timing
188  if (tdis6 /= '') then
189  call tdis_cr(tdis6, tdis_input_mempath)
190  else
191  call store_error('TIMING block variable TDIS6 is unset'// &
192  ' in simulation control input.', terminate)
193  end if
194  !
195  write (iout, '(1x,a)') 'END OF SIMULATION TIMING'
subroutine, public tdis_cr(fname, inmempath)
Create temporal discretization.
Definition: tdis.f90:50
Here is the call graph for this function:
Here is the caller graph for this function: