MODFLOW 6  version 6.7.0.dev1
USGS Modular Hydrologic Model
mf6core.f90
Go to the documentation of this file.
1 !> @brief Core MODFLOW 6 module
2 !!
3 !! This module contains the core components for MODFLOW 6. This module
4 !! is used by the stand-alone executable and the share object versions
5 !! of MODFLOW 6.
6 !!
7 !<
9  use kindmodule, only: i4b, lgp
20  use simstagesmodule
21  use profilermodule
22  implicit none
23 
24  class(runcontroltype), pointer :: run_ctrl => null() !< the run controller for this simulation
25 
26 contains
27 
28  !> @brief Main controller
29  !!
30  !! This subroutine is the main controller for MODFLOW 6.
31  !!
32  !<
33  subroutine mf6run
34  ! -- modules
36  use tdismodule, only: endofsimulation
37  ! -- local
38  logical(LGP) :: hasConverged
39  !
40  ! -- parse any command line arguments
42  !
43  ! initialize simulation
44  call mf6initialize()
45  !
46  ! -- time loop
47  do while (.not. endofsimulation)
48 
49  ! perform a time step
50  hasconverged = mf6update()
51 
52  ! if not converged, break
53  if (.not. hasconverged) exit
54 
55  end do
56  !
57  ! -- finalize simulation
58  call mf6finalize()
59 
60  end subroutine mf6run
61 
62  !> @brief Initialize a simulation
63  !!
64  !! This subroutine initializes a MODFLOW 6 simulation. The subroutine:
65  !! - creates the simulation
66  !! - defines
67  !! - allocates and reads static data
68  !!
69  !<
70  subroutine mf6initialize()
71  ! -- modules
74  use sourceloadmodule, only: export_cr
75 
76  ! init timer and start
77  call g_prof%initialize()
78  call g_prof%start("Run", g_prof%tmr_run)
79  call g_prof%start("Initialize", g_prof%tmr_init)
80 
81  ! -- get the run controller for sequential or parallel builds
83  call run_ctrl%start()
84 
85  ! -- print info and start timer
86  call print_info()
87 
88  ! -- create mfsim.lst
89  call create_lstfile()
90 
91  ! -- load input context
92  call static_input_load()
93 
94  ! -- create
95  call simulation_cr()
96 
97  ! -- define
98  call simulation_df()
99 
100  ! -- allocate and read
101  call simulation_ar()
102 
103  ! -- create model exports
104  call export_cr()
105 
106  ! -- stop the timer
107  call g_prof%stop(g_prof%tmr_init)
108 
109  end subroutine mf6initialize
110 
111  !> @brief Run a time step
112  !!
113  !! This function runs a single time step to completion.
114  !!
115  !! @return hasConverged boolean indicating if convergence was achieved for the time step
116  !!
117  !<
118  function mf6update() result(hasConverged)
119  logical(LGP) :: hasconverged
120  ! start timer
121  call g_prof%start("Update", g_prof%tmr_update)
122  !
123  ! -- prepare timestep
124  call mf6preparetimestep()
125  !
126  ! -- do timestep
127  call mf6dotimestep()
128  !
129  ! -- after timestep
130  hasconverged = mf6finalizetimestep()
131 
132  ! stop timer
133  call g_prof%stop(g_prof%tmr_update)
134 
135  end function mf6update
136 
137  !> @brief Finalize the simulation
138  !!
139  !! This subroutine finalizes a simulation. Steps include:
140  !! - final processing
141  !! - deallocate memory
142  !!
143  !<
144  subroutine mf6finalize()
145  ! -- modules
146  use, intrinsic :: iso_fortran_env, only: output_unit
147  use listsmodule, only: lists_da
149  use tdismodule, only: tdis_da
150  use idmloadmodule, only: idm_da
151  use sourceloadmodule, only: export_da
152  use simvariablesmodule, only: iout
153  ! -- local variables
154  integer(I4B) :: im
155  integer(I4B) :: ic
156  integer(I4B) :: is
157  integer(I4B) :: isg
158  class(solutiongrouptype), pointer :: sgp => null()
159  class(basesolutiontype), pointer :: sp => null()
160  class(basemodeltype), pointer :: mp => null()
161  class(baseexchangetype), pointer :: ep => null()
162  class(spatialmodelconnectiontype), pointer :: mc => null()
163  integer(I4B) :: tmr_dealloc
164 
165  ! start timer
166  call g_prof%start("Finalize", g_prof%tmr_finalize)
167 
168  !
169  ! -- FINAL PROCESSING (FP)
170  ! -- Final processing for each model
171  do im = 1, basemodellist%Count()
172  mp => getbasemodelfromlist(basemodellist, im)
173  call mp%model_fp()
174  end do
175  !
176  ! -- Final processing for each exchange
177  do ic = 1, baseexchangelist%Count()
178  ep => getbaseexchangefromlist(baseexchangelist, ic)
179  call ep%exg_fp()
180  end do
181  !
182  ! -- Final processing for each solution
183  do is = 1, basesolutionlist%Count()
184  sp => getbasesolutionfromlist(basesolutionlist, is)
185  call sp%sln_fp()
186  end do
187 
188  ! start timer for deallocation
189  tmr_dealloc = -1
190  call g_prof%start("Deallocate", tmr_dealloc)
191 
192  !
193  ! -- DEALLOCATE (DA)
194  ! -- Deallocate tdis
195  call tdis_da()
196  !
197  ! -- Deallocate for each model
198  do im = 1, basemodellist%Count()
199  mp => getbasemodelfromlist(basemodellist, im)
200  call mp%model_da()
201  deallocate (mp)
202  end do
203  !
204  ! -- Deallocate for each exchange
205  do ic = 1, baseexchangelist%Count()
206  ep => getbaseexchangefromlist(baseexchangelist, ic)
207  call ep%exg_da()
208  deallocate (ep)
209  end do
210  !
211  ! -- Deallocate for each connection
212  do ic = 1, baseconnectionlist%Count()
213  mc => get_smc_from_list(baseconnectionlist, ic)
214  call mc%exg_da()
215  deallocate (mc)
216  end do
217  !
218  ! -- Deallocate for each solution
219  do is = 1, basesolutionlist%Count()
220  sp => getbasesolutionfromlist(basesolutionlist, is)
221  call sp%sln_da()
222  deallocate (sp)
223  end do
224  !
225  ! -- Deallocate solution group and simulation variables
226  do isg = 1, solutiongrouplist%Count()
227  sgp => getsolutiongroupfromlist(solutiongrouplist, isg)
228  call sgp%sgp_da()
229  deallocate (sgp)
230  end do
231  !
232  call idm_da(iout)
233  call export_da()
234  call simulation_da()
235  call lists_da()
236 
237  ! stop timer
238  call g_prof%stop(tmr_dealloc)
239 
240  ! finish gently (No calls after this)
241  ! timer is stopped inside because this call does not return
242  call run_ctrl%finish()
243 
244  end subroutine mf6finalize
245 
246  !> @brief print initial message
247  !<
248  subroutine print_info()
249  use simmodule, only: initial_message
250  use timermodule, only: print_start_time
251 
252  ! print initial message
253  call initial_message()
254 
255  ! get start time
256  call print_start_time()
257 
258  end subroutine print_info
259 
260  !> @brief Set up mfsim list file output logging
261  !!
262  !! This subroutine creates the mfsim list file
263  !! and writes the header.
264  !!
265  !<
266  subroutine create_lstfile()
267  use constantsmodule, only: linelength
270  use messagemodule, only: write_message
272  character(len=LINELENGTH) :: line
273  !
274  ! -- Open simulation list file
275  iout = getunit()
276  !
277  if (nr_procs > 1) then
279  end if
280  !
281  call openfile(iout, 0, simlstfile, 'LIST', filstat_opt='REPLACE')
282  !
283  ! -- write simlstfile to stdout
284  write (line, '(2(1x,A))') 'Writing simulation list file:', &
285  trim(adjustl(simlstfile))
286  !
287  call write_message(line)
289  end subroutine create_lstfile
290 
291  !> @brief Create simulation input context
292  !!
293  !! This subroutine creates the simulation input context
294  !!
295  !<
296  subroutine static_input_load()
297  ! -- modules
298  use constantsmodule, only: lenmempath
299  use simvariablesmodule, only: iout
300  use idmloadmodule, only: simnam_load, simtdis_load, &
304  use simvariablesmodule, only: iparamlog
305  !
306  ! -- load simnam input context
307  call simnam_load(iparamlog)
308  !
309  ! -- load tdis to input context
310  call simtdis_load()
311  !
312  ! -- load in scope models
313  call load_models(iout)
314  !
315  ! -- load in scope exchanges
316  call load_exchanges(iout)
317  end subroutine static_input_load
318 
319  !> @brief Define the simulation
320  !!
321  !! This subroutine defined the simulation. Steps include:
322  !! - define each model
323  !! - define each solution
324  !!
325  !<
326  subroutine simulation_df()
327  ! -- modules
328  use idmloadmodule, only: idm_df
329  ! -- local variables
330  integer(I4B) :: im
331  integer(I4B) :: ic
332  integer(I4B) :: is
333  class(basesolutiontype), pointer :: sp => null()
334  class(basemodeltype), pointer :: mp => null()
335  class(baseexchangetype), pointer :: ep => null()
336  class(spatialmodelconnectiontype), pointer :: mc => null()
337 
338  ! -- init virtual data environment
339  call run_ctrl%at_stage(stg_bfr_mdl_df)
340 
341  ! -- Define each model
342  do im = 1, basemodellist%Count()
344  call mp%model_df()
345  end do
346  !
347  ! -- synchronize
348  call run_ctrl%at_stage(stg_aft_mdl_df)
349  !
350  ! -- Define each exchange
351  do ic = 1, baseexchangelist%Count()
353  call ep%exg_df()
354  end do
355  !
356  ! -- synchronize
357  call run_ctrl%at_stage(stg_aft_exg_df)
358  !
359  ! -- when needed, this is where the interface models are
360  ! created and added to the numerical solutions
361  call connections_cr()
362  !
363  ! -- synchronize
364  call run_ctrl%at_stage(stg_aft_con_cr)
365  !
366  ! -- synchronize
367  call run_ctrl%at_stage(stg_bfr_con_df)
368  !
369  ! -- Define each connection
370  do ic = 1, baseconnectionlist%Count()
372  call mc%exg_df()
373  end do
374  !
375  ! -- synchronize
376  call run_ctrl%at_stage(stg_aft_con_df)
377  !
378  ! -- Define each solution
379  do is = 1, basesolutionlist%Count()
381  call sp%sln_df()
382  end do
383 
384  ! idm df
385  call idm_df()
386 
387  end subroutine simulation_df
388 
389  !> @brief Simulation allocate and read
390  !!
391  !! This subroutine allocates and reads static data for the simulation.
392  !! Steps include:
393  !! - allocate and read for each model
394  !! - allocate and read for each exchange
395  !! - allocate and read for each solution
396  !!
397  !<
398  subroutine simulation_ar()
400  ! -- local variables
401  integer(I4B) :: im
402  integer(I4B) :: ic
403  integer(I4B) :: is
404  class(basesolutiontype), pointer :: sp => null()
405  class(basemodeltype), pointer :: mp => null()
406  class(baseexchangetype), pointer :: ep => null()
407  class(spatialmodelconnectiontype), pointer :: mc => null()
408 
409  ! -- Allocate and read each model
410  do im = 1, basemodellist%Count()
412  call mp%model_ar()
413  end do
414  !
415  ! -- Allocate and read each exchange
416  do ic = 1, baseexchangelist%Count()
418  call ep%exg_ar()
419  end do
420  !
421  ! -- Synchronize
422  call run_ctrl%at_stage(stg_bfr_con_ar)
423  !
424  ! -- Allocate and read all model connections
425  do ic = 1, baseconnectionlist%Count()
427  call mc%exg_ar()
428  end do
429  !
430  ! -- Synchronize
431  call run_ctrl%at_stage(stg_aft_con_ar)
432  !
433  ! -- Allocate and read each solution
434  do is = 1, basesolutionlist%Count()
436  call sp%sln_ar()
437  end do
438  !
439  end subroutine simulation_ar
440 
441  !> @brief Create the model connections from the exchanges
442  !!
443  !! This will upgrade the numerical exchanges in the solution,
444  !! whenever the configuration requires this, to Connection
445  !! objects. Currently we anticipate:
446  !!
447  !! GWF-GWF => GwfGwfConnection
448  !! GWT-GWT => GwtGwtConecction
449  !<
450  subroutine connections_cr()
452  use simvariablesmodule, only: iout
453  use versionmodule, only: idevelopmode
454  integer(I4B) :: isol
455  type(connectionbuildertype) :: connectionBuilder
456  class(basesolutiontype), pointer :: sol => null()
457  integer(I4B) :: status
458  character(len=16) :: envvar
459 
460  write (iout, '(/a)') 'PROCESSING MODEL CONNECTIONS'
461 
462  if (baseexchangelist%Count() == 0) then
463  ! if this is not a coupled simulation in any way,
464  ! then we will not need model connections
465  return
466  end if
467 
468  if (idevelopmode == 1) then
469  call get_environment_variable('DEV_ALWAYS_USE_IFMOD', &
470  value=envvar, status=status)
471  if (status == 0 .and. envvar == '1') then
472  connectionbuilder%dev_always_ifmod = .true.
473  write (iout, '(/a)') "Development option: forcing interface model"
474  end if
475  end if
476 
477  do isol = 1, basesolutionlist%Count()
479  call connectionbuilder%processSolution(sol)
480  end do
481 
482  write (iout, '(a)') 'END OF MODEL CONNECTIONS'
483  end subroutine connections_cr
484 
485  !> @brief Read and prepare time step
486  !!
487  !! This subroutine reads and prepares period data for the simulation.
488  !! Steps include:
489  !! - read and prepare for each model
490  !! - read and prepare for each exchange
491  !! - reset convergence flag
492  !! - calculate maximum time step for each model
493  !! - calculate maximum time step for each exchange
494  !! - calculate maximum time step for each solution
495  !! - set time discretization timestep using smallest maximum timestep
496  !!
497  !<
498  subroutine mf6preparetimestep()
499  ! -- modules
500  use kindmodule, only: i4b
503  kstp, kper
508  use simmodule, only: converge_reset
509  use simvariablesmodule, only: isim_mode
510  use idmloadmodule, only: idm_rp
512  ! -- local variables
513  class(basemodeltype), pointer :: mp => null()
514  class(baseexchangetype), pointer :: ep => null()
515  class(spatialmodelconnectiontype), pointer :: mc => null()
516  class(basesolutiontype), pointer :: sp => null()
517  character(len=LINELENGTH) :: line
518  character(len=LINELENGTH) :: fmt
519  integer(I4B) :: im
520  integer(I4B) :: ie
521  integer(I4B) :: ic
522  integer(I4B) :: is
523 
524  ! start timer
525  call g_prof%start("Prepare time step", g_prof%tmr_prep_tstp)
526 
527  !
528  ! -- initialize fmt
529  fmt = "(/,a,/)"
530  !
531  ! -- period update
532  call tdis_set_counters()
533  !
534  ! -- set base line
535  write (line, '(a,i0,a,i0,a)') &
536  'start timestep kper="', kper, '" kstp="', kstp, '" mode="'
537  !
538  ! -- evaluate simulation mode
539  select case (isim_mode)
540  case (mvalidate)
541  line = trim(line)//'validate"'
542  case (mnormal)
543  line = trim(line)//'normal"'
544  end select
545 
546  ! -- load dynamic input
547  call idm_rp()
548 
549  ! -- Read and prepare each model
550  do im = 1, basemodellist%Count()
552  call mp%model_message(line, fmt=fmt)
553  call mp%model_rp()
554  end do
555  !
556  ! -- Synchronize
557  call run_ctrl%at_stage(stg_bfr_exg_rp)
558  !
559  ! -- Read and prepare each exchange
560  do ie = 1, baseexchangelist%Count()
562  call ep%exg_rp()
563  end do
564  !
565  ! -- Read and prepare each connection
566  do ic = 1, baseconnectionlist%Count()
567  mc => get_smc_from_list(baseconnectionlist, ic)
568  call mc%exg_rp()
569  end do
570  !
571  ! -- Synchronize
572  call run_ctrl%at_stage(stg_aft_con_rp)
573  !
574  ! -- reset simulation convergence flag
575  call converge_reset()
576  !
577  ! -- time update for each model
578  do im = 1, basemodellist%Count()
580  call mp%model_dt()
581  end do
582  !
583  ! -- time update for each exchange
584  do ie = 1, baseexchangelist%Count()
586  call ep%exg_dt()
587  end do
588  !
589  ! -- time update for each connection
590  do ic = 1, baseconnectionlist%Count()
591  mc => get_smc_from_list(baseconnectionlist, ic)
592  call mc%exg_dt()
593  end do
594  !
595  ! -- time update for each solution
596  do is = 1, basesolutionlist%Count()
597  sp => getbasesolutionfromlist(basesolutionlist, is)
598  call sp%sln_dt()
599  end do
600  !
601  ! -- update exports
602  call export_post_prepare()
603  !
604  ! -- set time step
605  call tdis_set_timestep()
606 
607  ! stop timer
608  call g_prof%stop(g_prof%tmr_prep_tstp)
609 
610  end subroutine mf6preparetimestep
611 
612  !> @brief Run time step
613  !!
614  !! This subroutine runs a single time step for the simulation.
615  !! Steps include:
616  !! - formulate the system of equations for each model and exchange
617  !! - solve each solution
618  !!
619  !<
620  subroutine mf6dotimestep()
621  ! -- modules
622  use kindmodule, only: i4b
623  use listsmodule, only: solutiongrouplist
626  use idmloadmodule, only: idm_ad
627  ! -- local variables
628  class(solutiongrouptype), pointer :: sgp => null()
629  integer(I4B) :: isg
630  logical :: finishedTrying
631 
632  ! start timer
633  call g_prof%start("Do time step", g_prof%tmr_do_tstp)
634 
635  ! -- By default, the solution groups will be solved once, and
636  ! may fail. But if adaptive stepping is active, then
637  ! the solution groups may be solved over and over with
638  ! progressively smaller time steps to see if convergence
639  ! can be obtained.
640  ifailedstepretry = 0
641  retryloop: do
642 
643  ! -- idm advance
644  call idm_ad()
645 
646  do isg = 1, solutiongrouplist%Count()
648  call sgp%sgp_ca()
649  end do
650 
651  call sim_step_retry(finishedtrying)
652  if (finishedtrying) exit retryloop
654 
655  end do retryloop
656 
657  ! stop timer
658  call g_prof%stop(g_prof%tmr_do_tstp)
659 
660  end subroutine mf6dotimestep
661 
662  !> @brief Rerun time step
663  !!
664  !! This subroutine reruns a single time step for the simulation when
665  !! the adaptive time step option is used.
666  !!
667  !<
668  subroutine sim_step_retry(finishedTrying)
669  ! -- modules
670  use kindmodule, only: dp
672  use simmodule, only: converge_reset
673  use tdismodule, only: kstp, kper, delt, tdis_delt_reset
675  ! -- dummy variables
676  logical, intent(out) :: finishedTrying !< boolean that indicates if no
677  ! additional reruns of the time step are required
678  !
679  ! -- Check with ats to reset delt and keep trying
680  finishedtrying = .true.
681  call ats_reset_delt(kstp, kper, laststepfailed, delt, finishedtrying)
682  !
683  if (.not. finishedtrying) then
684  !
685  ! -- Reset delt, which requires updating pertim, totim
686  ! and end of period and simulation indicators
687  call tdis_delt_reset(delt)
688  !
689  ! -- Reset state of the simulation convergence flag
690  call converge_reset()
691 
692  end if
693  end subroutine sim_step_retry
694 
695  !> @brief Finalize time step
696  !!
697  !! This function finalizes a single time step for the simulation
698  !! and writes output for the time step. Steps include:
699  !! - write output for each model
700  !! - write output for each exchange
701  !! - write output for each solutions
702  !! - perform a final convergence check and whether the simulation
703  !! can continue if convergence was not achieved
704  !!
705  !! @return hasConverged boolean indicating if convergence was achieved for the time step
706  !!
707  !<
708  function mf6finalizetimestep() result(hasConverged)
709  ! -- modules
710  use kindmodule, only: i4b
716  use simmodule, only: converge_check
717  use simvariablesmodule, only: isim_mode
719  ! -- return variable
720  logical(LGP) :: hasconverged
721  ! -- local variables
722  class(basesolutiontype), pointer :: sp => null()
723  class(basemodeltype), pointer :: mp => null()
724  class(baseexchangetype), pointer :: ep => null()
725  class(spatialmodelconnectiontype), pointer :: mc => null()
726  character(len=LINELENGTH) :: line
727  character(len=LINELENGTH) :: fmt
728  integer(I4B) :: im
729  integer(I4B) :: ix
730  integer(I4B) :: ic
731  integer(I4B) :: is
732  !
733  ! -- initialize format and line
734  fmt = "(/,a,/)"
735  line = 'end timestep'
736 
737  ! start timer
738  call g_prof%start("Finalize time step", g_prof%tmr_final_tstp)
739 
740  !
741  ! -- evaluate simulation mode
742  select case (isim_mode)
743  case (mvalidate)
744  !
745  ! -- Write final message for timestep for each model
746  do im = 1, basemodellist%Count()
748  call mp%model_message(line, fmt=fmt)
749  end do
750  case (mnormal)
751 
752  call g_prof%start("Write output", g_prof%tmr_output)
753  !
754  ! -- Write output and final message for timestep for each model
755  do im = 1, basemodellist%Count()
757  call mp%model_ot()
758  call mp%model_message(line, fmt=fmt)
759  end do
760  !
761  ! -- Write output for each exchange
762  do ix = 1, baseexchangelist%Count()
764  call ep%exg_ot()
765  end do
766  !
767  ! -- Write output for each connection
768  do ic = 1, baseconnectionlist%Count()
769  mc => get_smc_from_list(baseconnectionlist, ic)
770  call mc%exg_ot()
771  end do
772  !
773  ! -- Write output for each solution
774  do is = 1, basesolutionlist%Count()
776  call sp%sln_ot()
777  end do
778  !
779  ! -- update exports
780  call g_prof%start("NetCDF export", g_prof%tmr_nc_export)
781  call export_post_step()
782  call g_prof%stop(g_prof%tmr_nc_export)
783 
784  call g_prof%stop(g_prof%tmr_output)
785  end select
786  !
787  ! -- Check if we're done
788  call converge_check(hasconverged)
789 
790  ! stop timer
791  call g_prof%stop(g_prof%tmr_final_tstp)
792 
793  end function mf6finalizetimestep
794 
795 end module mf6coremodule
subroutine, public ats_reset_delt(kstp, kper, lastStepFailed, delt, finishedTrying)
@ brief Reset time step because failure has occurred
Definition: ats.f90:606
class(baseexchangetype) function, pointer, public getbaseexchangefromlist(list, idx)
Retrieve a specific BaseExchangeType object from a list.
class(basemodeltype) function, pointer, public getbasemodelfromlist(list, idx)
Definition: BaseModel.f90:172
class(basesolutiontype) function, pointer, public getbasesolutionfromlist(list, idx)
subroutine, public getcommandlinearguments()
Get command line arguments.
Definition: comarg.f90:29
This module contains simulation constants.
Definition: Constants.f90:9
integer(i4b), parameter linelength
maximum length of a standard line
Definition: Constants.f90:45
@ mvalidate
validation mode - do not run time steps
Definition: Constants.f90:205
@ mnormal
normal output mode
Definition: Constants.f90:206
integer(i4b), parameter lenmempath
maximum length of the memory path
Definition: Constants.f90:27
This module contains the IdmLoadModule.
Definition: IdmLoad.f90:7
subroutine, public simnam_load(paramlog)
MODFLOW 6 mfsim.nam input load routine.
Definition: IdmLoad.f90:441
subroutine, public idm_da(iout)
idm deallocate routine
Definition: IdmLoad.f90:74
subroutine, public idm_df()
advance package dynamic data for period steps
Definition: IdmLoad.f90:38
subroutine, public simtdis_load()
MODFLOW 6 tdis input load routine.
Definition: IdmLoad.f90:456
subroutine, public idm_rp()
load package dynamic data for period
Definition: IdmLoad.f90:50
subroutine, public idm_ad()
advance package dynamic data for period steps
Definition: IdmLoad.f90:62
subroutine, public load_models(iout)
load model namfiles and model package files
Definition: IdmLoad.f90:217
subroutine, public load_exchanges(iout)
load exchange files
Definition: IdmLoad.f90:284
integer(i4b) function, public getunit()
Get a free unit number.
subroutine, public append_processor_id(name, proc_id)
Append processor id to a string.
subroutine, public openfile(iu, iout, fname, ftype, fmtarg_opt, accarg_opt, filstat_opt, mode_opt)
Open a file.
Definition: InputOutput.f90:30
This module defines variable data types.
Definition: kind.f90:8
type(listtype), public basemodellist
Definition: mf6lists.f90:16
type(listtype), public baseexchangelist
Definition: mf6lists.f90:25
type(listtype), public solutiongrouplist
Definition: mf6lists.f90:22
type(listtype), public baseconnectionlist
Definition: mf6lists.f90:28
type(listtype), public basesolutionlist
Definition: mf6lists.f90:19
subroutine, public lists_da()
Definition: mf6lists.f90:34
character(len=lenmempath) function create_mem_path(component, subcomponent, context)
returns the path to the memory object
Store and issue logging messages to output units.
Definition: Message.f90:2
subroutine, public write_message(text, iunit, fmt, skipbefore, skipafter, advance)
Write a message to an output unit.
Definition: Message.f90:210
Core MODFLOW 6 module.
Definition: mf6core.f90:8
subroutine mf6dotimestep()
Run time step.
Definition: mf6core.f90:621
subroutine mf6run
Main controller.
Definition: mf6core.f90:34
subroutine static_input_load()
Create simulation input context.
Definition: mf6core.f90:297
subroutine mf6preparetimestep()
Read and prepare time step.
Definition: mf6core.f90:499
subroutine simulation_df()
Define the simulation.
Definition: mf6core.f90:327
subroutine simulation_ar()
Simulation allocate and read.
Definition: mf6core.f90:399
class(runcontroltype), pointer run_ctrl
the run controller for this simulation
Definition: mf6core.f90:24
logical(lgp) function mf6update()
Run a time step.
Definition: mf6core.f90:119
subroutine sim_step_retry(finishedTrying)
Rerun time step.
Definition: mf6core.f90:669
subroutine mf6initialize()
Initialize a simulation.
Definition: mf6core.f90:71
logical(lgp) function mf6finalizetimestep()
Finalize time step.
Definition: mf6core.f90:709
subroutine create_lstfile()
Set up mfsim list file output logging.
Definition: mf6core.f90:267
subroutine connections_cr()
Create the model connections from the exchanges.
Definition: mf6core.f90:451
subroutine mf6finalize()
Finalize the simulation.
Definition: mf6core.f90:145
subroutine print_info()
print initial message
Definition: mf6core.f90:249
type(profilertype), public g_prof
the global timer object (to reduce trivial lines of code)
Definition: Profiler.f90:65
class(runcontroltype) function, pointer, public create_run_control()
This module contains simulation methods.
Definition: Sim.f90:10
subroutine, public converge_reset()
Reset the simulation convergence flag.
Definition: Sim.f90:389
subroutine, public initial_message()
Print the header and initializes messaging.
Definition: Sim.f90:442
subroutine, public converge_check(hasConverged)
Simulation convergence check.
Definition: Sim.f90:402
integer(i4b), parameter, public stg_aft_con_ar
afterr connection allocate read
Definition: SimStages.f90:18
integer(i4b), parameter, public stg_aft_con_df
after connection define
Definition: SimStages.f90:15
integer(i4b), parameter, public stg_aft_mdl_df
after model define
Definition: SimStages.f90:11
integer(i4b), parameter, public stg_bfr_mdl_df
before model define
Definition: SimStages.f90:10
integer(i4b), parameter, public stg_aft_exg_df
after exchange define
Definition: SimStages.f90:12
integer(i4b), parameter, public stg_aft_con_cr
after connection create
Definition: SimStages.f90:13
integer(i4b), parameter, public stg_bfr_exg_rp
before exchange read prepare
Definition: SimStages.f90:19
integer(i4b), parameter, public stg_bfr_con_ar
before connection allocate read
Definition: SimStages.f90:17
integer(i4b), parameter, public stg_aft_con_rp
after connection read prepare
Definition: SimStages.f90:20
integer(i4b), parameter, public stg_bfr_con_df
before connection define
Definition: SimStages.f90:14
subroutine, public simulation_da()
Deallocate simulation variables.
subroutine, public simulation_cr()
Read the simulation name file and initialize the models, exchanges.
This module contains simulation variables.
Definition: SimVariables.f90:9
integer(i4b) ifailedstepretry
current retry for this time step
integer(i4b) nr_procs
integer(i4b) laststepfailed
flag indicating if the last step failed (1) if last step failed; (0) otherwise (set in converge_check...
integer(i4b) iout
file unit number for simulation output
integer(i4b) iparamlog
input (idm) parameter logging to simulation listing file
character(len=linelength) simlstfile
simulation listing file name
integer(i4b) proc_id
integer(i4b) isim_mode
simulation mode
class(solutiongrouptype) function, pointer, public getsolutiongroupfromlist(list, idx)
This module contains the SourceLoadModule.
Definition: SourceLoad.F90:8
subroutine, public export_da()
deallocate model export objects and list
Definition: SourceLoad.F90:336
subroutine, public export_cr()
create model exports list
Definition: SourceLoad.F90:301
subroutine, public export_post_prepare()
model exports post prepare step actions
Definition: SourceLoad.F90:322
subroutine, public export_post_step()
model exports post step actions
Definition: SourceLoad.F90:329
class(spatialmodelconnectiontype) function, pointer, public get_smc_from_list(list, idx)
Get the connection from a list.
logical(lgp), pointer, public endofsimulation
flag indicating end of simulation
Definition: tdis.f90:28
subroutine, public tdis_set_timestep()
Set time step length.
Definition: tdis.f90:153
subroutine, public tdis_set_counters()
Set kstp and kper.
Definition: tdis.f90:91
integer(i4b), pointer, public kstp
current time step number
Definition: tdis.f90:24
subroutine, public tdis_da()
Deallocate memory.
Definition: tdis.f90:345
subroutine, public tdis_delt_reset(deltnew)
Reset delt and update timing variables and indicators.
Definition: tdis.f90:213
integer(i4b), pointer, public kper
current stress period number
Definition: tdis.f90:23
real(dp), pointer, public delt
length of the current time step
Definition: tdis.f90:29
subroutine, public print_start_time()
Start simulation timer.
Definition: Timer.f90:19
This module contains version information.
Definition: version.f90:7
subroutine write_listfile_header(iout, cmodel_type, write_sys_command, write_kind_info)
@ brief Write program header
Definition: version.f90:98
integer(i4b), parameter idevelopmode
Definition: version.f90:19
Highest level model type. All models extend this parent type.
Definition: BaseModel.f90:13
Class to manage spatial connection of a model to one or more models of the same type....