init_common.F90 coverage: 30.00 %func 10.94 %block
1) module Init_Common_module
2)
3) use PFLOTRAN_Constants_module
4)
5) implicit none
6)
7) private
8)
9) #include "petsc/finclude/petscsys.h"
10)
11) #include "petsc/finclude/petscvec.h"
12) #include "petsc/finclude/petscvec.h90"
13) #include "petsc/finclude/petscmat.h"
14) #include "petsc/finclude/petscmat.h90"
15) #include "petsc/finclude/petscsnes.h"
16) #include "petsc/finclude/petscpc.h"
17) #include "petsc/finclude/petscts.h"
18)
19) public :: &
20) ! Init, &
21) InitCommonReadRegionFiles, &
22) #if defined(PETSC_HAVE_HDF5)
23) InitCommonReadVelocityField, &
24) #endif
25) InitCommonVerifyAllCouplers, &
26) setSurfaceFlowMode, &
27) InitCommonAddOutputWaypoints
28)
29) #if defined(SCORPIO)
30) public :: InitCommonCreateIOGroups
31) #endif
32)
33) contains
34)
35) ! ************************************************************************** !
36)
37) subroutine InitReadInputFilenames(option,filenames)
38) !
39) ! Reads filenames for multi-simulation runs
40) !
41) ! Author: Glenn Hammond
42) ! Date: 08/11/09
43) !
44)
45) use Option_module
46) use Input_Aux_module
47)
48) type(option_type) :: option
49) character(len=MAXSTRINGLENGTH), pointer :: filenames(:)
50)
51) character(len=MAXSTRINGLENGTH) :: string
52) character(len=MAXSTRINGLENGTH) :: filename
53) PetscInt :: filename_count
54) type(input_type), pointer :: input
55) PetscBool :: card_found
56)
57) input => InputCreate(IN_UNIT,option%input_filename,option)
58)
59) string = "FILENAMES"
60) call InputFindStringInFile(input,option,string)
61)
62) card_found = PETSC_FALSE
63) if (InputError(input)) then
64) ! if the FILENAMES card is not included, we will assume that only
65) ! filenames exist in the file.
66) rewind(input%fid)
67) else
68) card_found = PETSC_TRUE
69) endif
70)
71) filename_count = 0
72) do
73) call InputReadPflotranString(input,option)
74) if (InputError(input)) exit
75) if (InputCheckExit(input,option)) exit
76) call InputReadNChars(input,option,filename,MAXSTRINGLENGTH,PETSC_FALSE)
77) filename_count = filename_count + 1
78) enddo
79)
80) allocate(filenames(filename_count))
81) filenames = ''
82) rewind(input%fid)
83)
84) if (card_found) then
85) string = "FILENAMES"
86) call InputFindStringInFile(input,option,string)
87) endif
88)
89) filename_count = 0
90) do
91) call InputReadPflotranString(input,option)
92) if (InputError(input)) exit
93) if (InputCheckExit(input,option)) exit
94) call InputReadNChars(input,option,filename,MAXSTRINGLENGTH,PETSC_FALSE)
95) filename_count = filename_count + 1
96) filenames(filename_count) = filename
97) enddo
98)
99) call InputDestroy(input)
100)
101) end subroutine InitReadInputFilenames
102)
103) ! ************************************************************************** !
104)
105) subroutine setSurfaceFlowMode(option)
106) !
107) ! Sets the flow mode for surface (richards, th, etc.)
108) !
109) ! Author: Gautam Bisht
110) ! Date: 07/30/14
111) !
112)
113) use Option_module
114) use String_module
115)
116) implicit none
117)
118) type(option_type) :: option
119)
120) select case(option%iflowmode)
121) case(RICHARDS_MODE)
122) option%nsurfflowdof = ONE_INTEGER
123) case(TH_MODE)
124) option%nsurfflowdof = TWO_INTEGER
125) case default
126) write(option%io_buffer,*) option%iflowmode
127) option%io_buffer = 'Flow Mode ' // &
128) trim(option%io_buffer) // ' not recognized in setSurfaceFlowMode().'
129) call printErrMsg(option)
130) end select
131)
132) end subroutine setSurfaceFlowMode
133)
134) ! ************************************************************************** !
135)
136) subroutine InitCommonVerifyAllCouplers(realization)
137) !
138) ! Verifies the connectivity of a coupler
139) !
140) ! Author: Glenn Hammond
141) ! Date: 1/8/08
142) !
143)
144) use Realization_Subsurface_class
145) use Patch_module
146) use Coupler_module
147)
148) implicit none
149)
150) class(realization_subsurface_type) :: realization
151)
152) type(patch_type), pointer :: cur_patch
153)
154) cur_patch => realization%patch_list%first
155) do
156) if (.not.associated(cur_patch)) exit
157)
158) call InitCommonVerifyCoupler(realization,cur_patch, &
159) cur_patch%initial_condition_list)
160) call InitCommonVerifyCoupler(realization,cur_patch, &
161) cur_patch%boundary_condition_list)
162) call InitCommonVerifyCoupler(realization,cur_patch, &
163) cur_patch%source_sink_list)
164)
165) cur_patch => cur_patch%next
166) enddo
167)
168) end subroutine InitCommonVerifyAllCouplers
169)
170) ! ************************************************************************** !
171)
172) subroutine InitCommonVerifyCoupler(realization,patch,coupler_list)
173) !
174) ! Verifies the connectivity of a coupler
175) !
176) ! Author: Glenn Hammond
177) ! Date: 1/8/08
178) !
179)
180) use Realization_Subsurface_class
181) use Discretization_module
182) use Option_module
183) use Coupler_module
184) use Condition_module
185) use Grid_module
186) use Output_module
187) use Output_Tecplot_module, only : OutputVectorTecplot
188) use Patch_module
189)
190) implicit none
191)
192) class(realization_subsurface_type) :: realization
193) type(coupler_list_type), pointer :: coupler_list
194)
195) type(option_type), pointer :: option
196) type(grid_type), pointer :: grid
197) type(patch_type), pointer :: patch
198) type(coupler_type), pointer :: coupler
199) character(len=MAXWORDLENGTH) :: word
200) character(len=MAXSTRINGLENGTH) :: filename
201) character(len=MAXSTRINGLENGTH) :: dataset_name
202) PetscInt :: iconn, icell, local_id
203) Vec :: global_vec
204) PetscReal, pointer :: vec_ptr(:)
205) PetscErrorCode :: ierr
206)
207) patch => realization%patch
208) grid => patch%grid
209) option => realization%option
210)
211) if (.not.associated(coupler_list)) return
212)
213) call DiscretizationCreateVector(realization%discretization,ONEDOF, &
214) global_vec,GLOBAL,option)
215)
216) coupler => coupler_list%first
217)
218) do
219) if (.not.associated(coupler)) exit
220)
221) call VecZeroEntries(global_vec,ierr);CHKERRQ(ierr)
222) call VecGetArrayF90(global_vec,vec_ptr,ierr);CHKERRQ(ierr)
223) if (associated(coupler%connection_set)) then
224) do iconn = 1, coupler%connection_set%num_connections
225) local_id = coupler%connection_set%id_dn(iconn)
226) ! vec_ptr(local_id) = coupler%id
227) !geh: let's sum the # of connections
228) vec_ptr(local_id) = vec_ptr(local_id) + 1
229) enddo
230) else
231) if (associated(coupler%region)) then
232) do icell = 1, coupler%region%num_cells
233) local_id = coupler%region%cell_ids(icell)
234) ! vec_ptr(local_id) = coupler%id
235) vec_ptr(local_id) = vec_ptr(local_id) + 1
236) enddo
237) endif
238) endif
239) call VecRestoreArrayF90(global_vec,vec_ptr,ierr);CHKERRQ(ierr)
240) if (len_trim(coupler%flow_condition_name) > 0) then
241) dataset_name = coupler%flow_condition_name
242) else if (len_trim(coupler%tran_condition_name) > 0) then
243) dataset_name = coupler%tran_condition_name
244) endif
245) write(word,*) patch%id
246) dataset_name = trim(dataset_name) // '_' // &
247) trim(coupler%region%name) // '_' // &
248) trim(adjustl(word))
249) dataset_name = dataset_name(1:28)
250) filename = trim(dataset_name) // '.tec'
251) call OutputVectorTecplot(filename,dataset_name,realization,global_vec)
252)
253) coupler => coupler%next
254) enddo
255)
256) call VecDestroy(global_vec,ierr);CHKERRQ(ierr)
257)
258) end subroutine InitCommonVerifyCoupler
259)
260) ! ************************************************************************** !
261)
262) subroutine InitCommonReadRegionFiles(realization)
263) !
264) ! Reads in grid cell ids stored in files
265) !
266) ! Author: Glenn Hammond
267) ! Date: 1/03/08
268) !
269)
270) use Realization_Subsurface_class
271) use Region_module
272) use HDF5_module
273) use Option_module
274)
275) implicit none
276)
277) class(realization_subsurface_type) :: realization
278)
279) type(option_type), pointer :: option
280) type(region_type), pointer :: region
281) PetscBool :: cell_ids_exists
282) PetscBool :: face_ids_exists
283) PetscBool :: vert_ids_exists
284)
285) option => realization%option
286) region => realization%region_list%first
287) do
288) if (.not.associated(region)) exit
289) if (len_trim(region%filename) > 1) then
290) if (index(region%filename,'.h5') > 0) then
291) if (.not.region%hdf5_ugrid_kludge) then
292)
293) call HDF5QueryRegionDefinition(region, region%filename, &
294) realization%option, &
295) cell_ids_exists, face_ids_exists, vert_ids_exists)
296)
297) if ( (.not. cell_ids_exists) .and. &
298) (.not. face_ids_exists) .and. &
299) (.not. vert_ids_exists)) then
300)
301) option%io_buffer = '"Regions/' // trim(region%name) // &
302) ' is not defined by "Cell Ids" or "Face Ids" or "Vertex Ids".'
303) call printErrMsg(option)
304) end if
305)
306) if (cell_ids_exists .or. face_ids_exists) then
307) call HDF5ReadRegionFromFile(realization%patch%grid,region, &
308) region%filename,option)
309) else
310) call HDF5ReadRegionDefinedByVertex(realization%option, &
311) region, region%filename)
312) end if
313) else
314) !geh: Do not skip this subroutine if PETSC_HAVE_HDF5 is not
315) ! defined. The subroutine prints an error message if not defined
316) ! informing the user of the error. If you skip the subroutine,
317) ! no error message is printed and the user is unaware of the
318) ! region not being read.
319) call HDF5ReadUnstructuredGridRegionFromFile(realization%option, &
320) region, &
321) region%filename)
322) endif
323) else if (index(region%filename,'.ss') > 0) then
324) region%def_type = DEFINED_BY_SIDESET_UGRID
325) region%sideset => RegionCreateSideset()
326) call RegionReadFromFile(region%sideset,region%filename, &
327) realization%option)
328) else if (index(region%filename,'.ex') > 0) then
329) region%def_type = DEFINED_BY_FACE_UGRID_EXP
330) call RegionReadFromFile(region%explicit_faceset,region%cell_ids, &
331) region%filename,realization%option)
332) region%num_cells = size(region%cell_ids)
333) else
334) call RegionReadFromFile(region,realization%option, &
335) region%filename)
336) endif
337) endif
338) region => region%next
339) enddo
340)
341) end subroutine InitCommonReadRegionFiles
342)
343) ! ************************************************************************** !
344)
345) subroutine readVectorFromFile(realization,vector,filename,vector_type)
346) !
347) ! Reads data from a file into an associated vector
348) !
349) ! Author: Glenn Hammond
350) ! Date: 03/18/08
351) !
352)
353) use Realization_Subsurface_class
354) use Discretization_module
355) use Field_module
356) use Grid_module
357) use Option_module
358) use Patch_module
359) use Logging_module
360)
361) use HDF5_module
362)
363) implicit none
364)
365) class(realization_subsurface_type) :: realization
366) Vec :: vector
367) character(len=MAXWORDLENGTH) :: filename
368) PetscInt :: vector_type
369)
370) type(discretization_type), pointer :: discretization
371) type(field_type), pointer :: field
372) type(grid_type), pointer :: grid
373) type(option_type), pointer :: option
374) type(patch_type), pointer :: patch
375) PetscInt :: ghosted_id, natural_id, material_id
376) PetscInt :: fid = 86
377) PetscInt :: status
378) PetscErrorCode :: ierr
379) PetscInt :: count, read_count, i
380) PetscInt :: flag
381) PetscInt, pointer :: indices(:)
382) PetscReal, pointer :: values(:)
383) PetscInt, parameter :: block_size = 10000
384) Vec :: natural_vec, global_vec
385)
386) discretization => realization%discretization
387) field => realization%field
388) patch => realization%patch
389) grid => patch%grid
390) option => realization%option
391)
392) if (index(filename,'.h5') > 0) then
393) ! to be taken care of later in readPermeabilitiesFromFile()
394) else
395) open(unit=fid,file=filename,status="old",iostat=status)
396) if (status /= 0) then
397) option%io_buffer = 'File: ' // trim(filename) // ' not found.'
398) call printErrMsg(option)
399) endif
400) allocate(values(block_size))
401) allocate(indices(block_size))
402) call DiscretizationCreateVector(discretization,ONEDOF,natural_vec, &
403) NATURAL,option)
404) count = 0
405) do
406) if (count >= grid%nmax) exit
407) read_count = min(block_size,grid%nmax-count)
408) do i=1,read_count
409) indices(i) = count+i-1 ! zero-based indexing
410) enddo
411) ierr = 0
412) if (option%myrank == option%io_rank) &
413) read(fid,*,iostat=ierr) values(1:read_count)
414) flag = ierr
415) call MPI_Bcast(flag,ONE_INTEGER_MPI,MPIU_INTEGER,option%io_rank, &
416) option%mycomm,ierr)
417) if (flag /= 0) then
418) option%io_buffer = 'Insufficent data in file: ' // filename
419) call printErrMsg(option)
420) endif
421) if (option%myrank == option%io_rank) then
422) call VecSetValues(natural_vec,read_count,indices,values,INSERT_VALUES, &
423) ierr);CHKERRQ(ierr)
424) endif
425) count = count + read_count
426) enddo
427) call MPI_Bcast(count,ONE_INTEGER_MPI,MPIU_INTEGER,option%io_rank, &
428) option%mycomm,ierr)
429) if (count /= grid%nmax) then
430) write(option%io_buffer,'("Number of data in file (",i8, &
431) & ") does not match size of vector (",i8,")")') count, grid%nlmax
432) call printErrMsg(option)
433) endif
434) close(fid)
435) deallocate(values)
436) nullify(values)
437) deallocate(indices)
438) nullify(indices)
439) call VecAssemblyBegin(natural_vec,ierr);CHKERRQ(ierr)
440) call VecAssemblyEnd(natural_vec,ierr);CHKERRQ(ierr)
441) select case(vector_type)
442) case(LOCAL)
443) call DiscretizationCreateVector(discretization,ONEDOF,global_vec, &
444) GLOBAL,option)
445) call DiscretizationNaturalToGlobal(discretization,natural_vec, &
446) global_vec,ONEDOF)
447) call DiscretizationGlobalToLocal(discretization,global_vec, &
448) vector,ONEDOF)
449) call VecDestroy(global_vec,ierr);CHKERRQ(ierr)
450) case(GLOBAL)
451) call DiscretizationNaturalToGlobal(discretization,natural_vec, &
452) vector,ONEDOF)
453) end select
454) call VecDestroy(natural_vec,ierr);CHKERRQ(ierr)
455) endif
456)
457) end subroutine readVectorFromFile
458)
459) ! ************************************************************************** !
460)
461) subroutine InitCommonCreateIOGroups(option)
462) !
463) ! Create sub-communicators that are used in initialization
464) ! and output HDF5 routines.
465) !
466) ! Author: Vamsi Sripathi
467) ! Date: 07/14/09
468) !
469)
470) use Option_module
471) use Logging_module
472)
473) #if defined(SCORPIO)
474) use hdf5
475) #endif
476)
477) implicit none
478)
479) type(option_type) :: option
480) PetscErrorCode :: ierr
481)
482) #if defined(SCORPIO)
483)
484) PetscMPIInt :: numiogroups
485)
486) call PetscLogEventBegin(logging%event_create_iogroups,ierr);CHKERRQ(ierr)
487)
488) ! Initialize HDF interface to define global constants
489) call h5open_f(ierr)
490)
491) if (option%hdf5_read_group_size <= 0) then
492) write(option%io_buffer,&
493) '("The keyword HDF5_READ_GROUP_SIZE &
494) & in the input file (pflotran.in) is either not set or &
495) & its value is less than or equal to ZERO. &
496) & HDF5_READ_GROUP_SIZE = ",i6)') &
497) option%hdf5_read_group_size
498) !call printErrMsg(option)
499) call printMsg(option)
500) ! default is to let one process read and broadcast to everyone
501) option%hdf5_read_group_size = option%mycommsize
502) endif
503)
504) if (option%hdf5_write_group_size <= 0) then
505) write(option%io_buffer,&
506) '("The keyword HDF5_WRITE_GROUP_SIZE &
507) &in the input file (pflotran.in) is either not set or &
508) &its value is less than or equal to ZERO. &
509) &HDF5_WRITE_GROUP_SIZE = ",i6)') &
510) option%hdf5_write_group_size
511) !call printErrMsg(option)
512) call printMsg(option)
513) ! default is to let everyone write separately
514) option%hdf5_write_group_size = 1
515) endif
516)
517) ! create read IO groups
518) numiogroups = option%mycommsize/option%hdf5_read_group_size
519) call fscorpio_iogroup_init(numiogroups, option%mycomm, &
520) option%ioread_group_id, ierr)
521)
522) if ( option%hdf5_read_group_size == option%hdf5_write_group_size ) then
523) ! reuse read_group to use for writing too as both groups are same size
524) option%iowrite_group_id = option%ioread_group_id
525) else
526) ! create write IO groups
527) numiogroups = option%mycommsize/option%hdf5_write_group_size
528) call fscorpio_iogroup_init(numiogroups, option%mycomm, &
529) option%iowrite_group_id, ierr)
530) end if
531)
532) write(option%io_buffer, '(" Read group id : ", i6)') option%ioread_group_id
533) call printMsg(option)
534) write(option%io_buffer, '(" Write group id : ", i6)') &
535) option%iowrite_group_id
536) call printMsg(option)
537) call PetscLogEventEnd(logging%event_create_iogroups,ierr);CHKERRQ(ierr)
538) #endif
539) ! SCORPIO
540)
541) end subroutine InitCommonCreateIOGroups
542)
543) ! ************************************************************************** !
544)
545) subroutine InitCommonPrintPFLOTRANHeader(option,fid)
546) !
547) ! Initializes pflotran
548) !
549) ! Author: Glenn Hammond
550) ! Date: 10/23/07
551) !
552)
553) use Option_module
554)
555) implicit none
556)
557) PetscInt :: fid
558)
559) type(option_type) :: option
560)
561) write(fid,'(" PFLOTRAN Header")')
562)
563) end subroutine InitCommonPrintPFLOTRANHeader
564)
565) ! ************************************************************************** !
566)
567) #if defined(PETSC_HAVE_HDF5)
568)
569) subroutine InitCommonReadVelocityField(realization)
570) !
571) ! Reads fluxes in for transport with no flow.
572) !
573) ! Author: Glenn Hammond
574) ! Date: 02/05/13
575) !
576)
577) use Realization_Subsurface_class
578) use Patch_module
579) use Field_module
580) use Grid_module
581) use Option_module
582) use Coupler_module
583) use Connection_module
584) use Discretization_module
585) use HDF5_module
586) use HDF5_Aux_module
587)
588) implicit none
589)
590) class(realization_subsurface_type) :: realization
591) character(len=MAXSTRINGLENGTH) :: filename
592)
593) type(field_type), pointer :: field
594) type(patch_type), pointer :: patch
595) type(grid_type), pointer :: grid
596) type(discretization_type), pointer :: discretization
597) type(option_type), pointer :: option
598) character(len=MAXSTRINGLENGTH) :: group_name
599) character(len=MAXSTRINGLENGTH) :: dataset_name
600) PetscInt :: idir, iconn, sum_connection
601) PetscInt :: ghosted_id_up, local_id
602) PetscErrorCode :: ierr
603)
604) PetscReal, pointer :: vec_loc_p(:)
605) PetscReal, pointer :: vec_p(:)
606) type(coupler_type), pointer :: boundary_condition
607) type(connection_set_list_type), pointer :: connection_set_list
608) type(connection_set_type), pointer :: cur_connection_set
609)
610) field => realization%field
611) patch => realization%patch
612) grid => patch%grid
613) option => realization%option
614) discretization => realization%discretization
615)
616) filename = realization%nonuniform_velocity_filename
617)
618) group_name = ''
619) do idir = 1, 3
620) select case(idir)
621) case(1)
622) dataset_name = 'Internal Velocity X'
623) case(2)
624) dataset_name = 'Internal Velocity Y'
625) case(3)
626) dataset_name = 'Internal Velocity Z'
627) end select
628) if (.not.HDF5DatasetExists(filename,group_name,dataset_name,option)) then
629) option%io_buffer = 'Dataset "' // trim(group_name) // '/' // &
630) trim(dataset_name) // &
631) '" not found in HDF5 file "' // trim(filename) // '".'
632) call printErrMsg(option)
633) endif
634) call HDF5ReadCellIndexedRealArray(realization,field%work,filename, &
635) group_name,dataset_name,PETSC_FALSE)
636) call DiscretizationGlobalToLocal(discretization,field%work,field%work_loc, &
637) ONEDOF)
638) call VecGetArrayF90(field%work_loc,vec_loc_p,ierr);CHKERRQ(ierr)
639) connection_set_list => grid%internal_connection_set_list
640) cur_connection_set => connection_set_list%first
641) sum_connection = 0
642) do
643) if (.not.associated(cur_connection_set)) exit
644) do iconn = 1, cur_connection_set%num_connections
645) sum_connection = sum_connection + 1
646) ghosted_id_up = cur_connection_set%id_up(iconn)
647) if (cur_connection_set%dist(idir,iconn) > 0.9d0) then
648) patch%internal_velocities(1,sum_connection) = vec_loc_p(ghosted_id_up)
649) endif
650) enddo
651) cur_connection_set => cur_connection_set%next
652) enddo
653) call VecRestoreArrayF90(field%work_loc,vec_loc_p,ierr);CHKERRQ(ierr)
654) enddo
655)
656) boundary_condition => patch%boundary_condition_list%first
657) sum_connection = 0
658) do
659) if (.not.associated(boundary_condition)) exit
660) dataset_name = boundary_condition%name
661) if (.not.HDF5DatasetExists(filename,group_name,dataset_name,option)) then
662) option%io_buffer = 'Dataset "' // trim(group_name) // '/' // &
663) trim(dataset_name) // &
664) '" not found in HDF5 file "' // trim(filename) // '".'
665) call printErrMsg(option)
666) endif
667) call HDF5ReadCellIndexedRealArray(realization,field%work,filename, &
668) group_name,dataset_name,PETSC_FALSE)
669) call VecGetArrayF90(field%work,vec_p,ierr);CHKERRQ(ierr)
670) cur_connection_set => boundary_condition%connection_set
671) do iconn = 1, cur_connection_set%num_connections
672) sum_connection = sum_connection + 1
673) local_id = cur_connection_set%id_dn(iconn)
674) patch%boundary_velocities(1,sum_connection) = vec_p(local_id)
675) enddo
676) call VecRestoreArrayF90(field%work,vec_p,ierr);CHKERRQ(ierr)
677) boundary_condition => boundary_condition%next
678) enddo
679)
680) end subroutine InitCommonReadVelocityField
681)
682) #endif
683)
684) ! ************************************************************************** !
685)
686) subroutine InitCommonAddOutputWaypoints(option,output_option,waypoint_list)
687) !
688) ! Adds waypoints associated with output options to waypoint list
689) !
690) ! Author: Glenn Hammond
691) ! Date: 02/04/16
692) !
693) use Output_Aux_module
694) use Waypoint_module
695) use Option_module
696) use Utility_module
697)
698) implicit none
699)
700) type(option_type) :: option
701) type(output_option_type) :: output_option
702) type(waypoint_list_type) :: waypoint_list
703)
704) type(waypoint_type), pointer :: waypoint
705) character(len=MAXWORDLENGTH) :: word
706) PetscReal :: temp_real
707) PetscReal :: final_time
708) PetscReal :: num_waypoints, warning_num_waypoints
709) PetscInt :: k
710)
711) final_time = WaypointListGetFinalTime(waypoint_list)
712) warning_num_waypoints = 15000.0
713)
714) ! Add waypoints for periodic snapshot output
715) if (output_option%periodic_snap_output_time_incr > 0.d0) then
716) temp_real = 0.d0
717) num_waypoints = final_time / output_option%periodic_snap_output_time_incr
718) if ((num_waypoints > warning_num_waypoints) .and. &
719) OptionPrintToScreen(option)) then
720) write(word,*) floor(num_waypoints)
721) write(*,*) 'WARNING: Large number (' // trim(adjustl(word)) // &
722) ') of periodic snapshot output requested.'
723) write(*,'(a64)',advance='no') ' Creating periodic output &
724) &waypoints . . . Progress: 0%-'
725) endif
726) k = 0
727) do
728) k = k + 1
729) temp_real = temp_real + output_option%periodic_snap_output_time_incr
730) if (temp_real > final_time) exit
731) waypoint => WaypointCreate()
732) waypoint%time = temp_real
733) waypoint%print_snap_output = PETSC_TRUE
734) call WaypointInsertInList(waypoint,waypoint_list)
735) if ((num_waypoints > warning_num_waypoints) .and. &
736) OptionPrintToScreen(option)) then
737) call PrintProgressBarInt(floor(num_waypoints),10,k)
738) endif
739) enddo
740) endif
741)
742) ! Add waypoints for periodic observation output
743) if (output_option%periodic_obs_output_time_incr > 0.d0) then
744) temp_real = 0.d0
745) num_waypoints = final_time / output_option%periodic_obs_output_time_incr
746) if ((num_waypoints > warning_num_waypoints) .and. &
747) OptionPrintToScreen(option)) then
748) write(word,*) floor(num_waypoints)
749) write(*,*) 'WARNING: Large number (' // trim(adjustl(word)) // &
750) ') of periodic observation output requested.'
751) write(*,'(a64)',advance='no') ' Creating periodic output &
752) &waypoints . . . Progress: 0%-'
753) endif
754) k = 0
755) do
756) k = k + 1
757) temp_real = temp_real + output_option%periodic_obs_output_time_incr
758) if (temp_real > final_time) exit
759) waypoint => WaypointCreate()
760) waypoint%time = temp_real
761) waypoint%print_obs_output = PETSC_TRUE
762) call WaypointInsertInList(waypoint,waypoint_list)
763) if ((num_waypoints > warning_num_waypoints) .and. &
764) OptionPrintToScreen(option)) then
765) call PrintProgressBarInt(floor(num_waypoints),10,k)
766) endif
767) enddo
768) endif
769)
770) ! Add waypoints for periodic mass balance output
771) if (output_option%periodic_msbl_output_time_incr > 0.d0) then
772) temp_real = 0.d0
773) num_waypoints = final_time / output_option%periodic_msbl_output_time_incr
774) if ((num_waypoints > warning_num_waypoints) .and. &
775) OptionPrintToScreen(option)) then
776) write(word,*) floor(num_waypoints)
777) write(*,*) 'WARNING: Large number (' // trim(adjustl(word)) // &
778) ') of periodic mass balance output requested.'
779) write(*,'(a64)',advance='no') ' Creating periodic output &
780) &waypoints . . . Progress: 0%-'
781) endif
782) k = 0
783) do
784) k = k + 1
785) temp_real = temp_real + output_option%periodic_msbl_output_time_incr
786) if (temp_real > final_time) exit
787) waypoint => WaypointCreate()
788) waypoint%time = temp_real
789) waypoint%print_msbl_output = PETSC_TRUE
790) call WaypointInsertInList(waypoint,waypoint_list)
791) if ((num_waypoints > warning_num_waypoints) .and. &
792) OptionPrintToScreen(option)) then
793) call PrintProgressBarInt(floor(num_waypoints),10,k)
794) endif
795) enddo
796) endif
797)
798) end subroutine InitCommonAddOutputWaypoints
799)
800) end module Init_Common_module