From 177cc7bbf481e4fdefc8c3833fa87d5454aa791c Mon Sep 17 00:00:00 2001 From: Eric Morway Date: Thu, 14 Dec 2023 09:21:00 -0800 Subject: [PATCH] chore(tdis.f90): cleanup docstrings in tdis class (#1492) * chore(tdis.f90): cleanup docstrings in tdis class * forgot to run fprettify * remove unused variables * specified wrong constant value --- src/Timing/tdis.f90 | 286 ++++++++++++++---------------------- src/Utilities/Constants.f90 | 7 + 2 files changed, 120 insertions(+), 173 deletions(-) diff --git a/src/Timing/tdis.f90 b/src/Timing/tdis.f90 index 44bf93a4b97..2b245c5345b 100644 --- a/src/Timing/tdis.f90 +++ b/src/Timing/tdis.f90 @@ -45,13 +45,9 @@ module TdisModule contains + !> @brief Create temporal discretization + !< subroutine tdis_cr(fname) -! ****************************************************************************** -! tdis_cr -- create temporal discretization. -! ****************************************************************************** -! -! SPECIFICATIONS: -! ------------------------------------------------------------------------------ ! -- modules use InputOutputModule, only: getunit, openfile use ConstantsModule, only: LINELENGTH, DZERO @@ -64,7 +60,6 @@ subroutine tdis_cr(fname) character(len=*), parameter :: fmtheader = & "(1X,/1X,'TDIS -- TEMPORAL DISCRETIZATION PACKAGE,', / & &' VERSION 1 : 11/13/2014 - INPUT READ FROM UNIT ',I4)" -! ------------------------------------------------------------------------------ ! ! -- Allocate the scalar variables call tdis_allocate_scalars() @@ -99,17 +94,13 @@ subroutine tdis_cr(fname) call ats_cr(inats, nper) end if ! - ! -- return + ! -- Return return end subroutine tdis_cr + !> @brief Set kstp and kper + !< subroutine tdis_set_counters() -! ****************************************************************************** -! tdis_set_counters -- Set kstp and kper -! ****************************************************************************** -! -! SPECIFICATIONS: -! ------------------------------------------------------------------------------ ! -- modules use ConstantsModule, only: DONE, DZERO, MNORMAL, MVALIDATE, DNODATA use SimVariablesModule, only: isim_mode @@ -131,7 +122,6 @@ subroutine tdis_set_counters() character(len=*), parameter :: fmtspits = & "(28X,'NUMBER OF TIME STEPS = ',I0,/ & &28X,'MULTIPLIER FOR DELT =',F10.3)" -! ------------------------------------------------------------------------------ ! ! -- Initialize variables for this step if (inats > 0) dtstable = DNODATA @@ -169,17 +159,13 @@ subroutine tdis_set_counters() end if end if ! - ! -- return + ! -- Return return end subroutine tdis_set_counters + !> @brief Set time step length + !< subroutine tdis_set_timestep() -! ****************************************************************************** -! tdis_set_timestep -- Set time step length -! ****************************************************************************** -! -! SPECIFICATIONS: -! ------------------------------------------------------------------------------ ! -- modules use ConstantsModule, only: DONE, DZERO use AdaptiveTimeStepModule, only: isAdaptivePeriod, & @@ -190,7 +176,6 @@ subroutine tdis_set_timestep() ! -- format character(len=*), parameter :: fmttsi = & "(28X,'INITIAL TIME STEP SIZE =',G15.7)" -! ------------------------------------------------------------------------------ ! ! -- Initialize adaptivePeriod = isAdaptivePeriod(kper) @@ -234,19 +219,16 @@ subroutine tdis_set_timestep() endofsimulation = .true. end if ! - ! -- return + ! -- Return return end subroutine tdis_set_timestep + !> @brief Reset delt and update timing variables and indicators + !! + !! This routine is called when a timestep fails to converge, and so it is + !! retried using a smaller time step (deltnew). + !< subroutine tdis_delt_reset(deltnew) -! ****************************************************************************** -! tdis_delt_reset -- reset delt and update timing variables and indicators. -! This routine is called when a timestep fails to converge, and so it is -! retried using a smaller time step (deltnew). -! ****************************************************************************** -! -! SPECIFICATIONS: -! ------------------------------------------------------------------------------ ! -- modules use ConstantsModule, only: DONE, DZERO use AdaptiveTimeStepModule, only: isAdaptivePeriod, & @@ -256,7 +238,6 @@ subroutine tdis_delt_reset(deltnew) real(DP), intent(in) :: deltnew ! -- local logical(LGP) :: adaptivePeriod -! ------------------------------------------------------------------------------ ! ! -- Set values adaptivePeriod = isAdaptivePeriod(kper) @@ -280,21 +261,15 @@ subroutine tdis_delt_reset(deltnew) totim = totalsimtime end if ! - ! -- return + ! -- Return return end subroutine tdis_delt_reset + !> @brief Set time step length + !< subroutine tdis_set_delt() -! ****************************************************************************** -! tdis_set_delt -- Set time step length -! ****************************************************************************** -! -! SPECIFICATIONS: -! ------------------------------------------------------------------------------ ! -- modules use ConstantsModule, only: DONE - ! -- local -! ------------------------------------------------------------------------------ ! if (kstp == 1) then ! -- Calculate the first value of delt for this stress period @@ -314,7 +289,7 @@ subroutine tdis_set_delt() delt = tsmult(kper) * delt end if ! - ! -- return + ! -- Return return end subroutine tdis_set_delt @@ -374,96 +349,91 @@ end subroutine tdis_set_delt ! totim = totalsimtime ! end if ! ! -! ! -- return +! ! -- Return ! return ! end subroutine tdis_set_delt_std + !> @brief Print simulation time + !< subroutine tdis_ot(iout) -! ****************************************************************************** -! PRINT SIMULATION TIME -! ****************************************************************************** -! -! SPECIFICATIONS: -! ------------------------------------------------------------------------------ + ! -- modules + use ConstantsModule, only: DZERO, DONE, DSIXTY, DSECPERHR, DHRPERDAY, & + DDYPERYR, DSECPERDY, DSECPERYR ! -- dummy integer(I4B), intent(in) :: iout ! -- local - real(DP) :: zero, cnv, delsec, totsec, persec, sixty, hrday, dayyr, & - delmn, delhr, totmn, tothr, totdy, totyr, permn, perhr, & - perdy, peryr, deldy, delyr -! ------------------------------------------------------------------------------ - WRITE (IOUT, 199) KSTP, KPER -199 FORMAT(1X, ///9X, 'TIME SUMMARY AT END OF TIME STEP', I5, & - & ' IN STRESS PERIOD ', I4) -!C -!C1------USE TIME UNIT INDICATOR TO GET FACTOR TO CONVERT TO SECONDS. - ZERO = 0.d0 - CNV = ZERO - IF (ITMUNI .EQ. 1) CNV = 1. - IF (ITMUNI .EQ. 2) CNV = 60. - IF (ITMUNI .EQ. 3) CNV = 3600. - IF (ITMUNI .EQ. 4) CNV = 86400. - IF (ITMUNI .EQ. 5) CNV = 31557600. -!C -!C2------IF FACTOR=0 THEN TIME UNITS ARE NON-STANDARD. - IF (CNV .NE. ZERO) GO TO 100 -!C -!C2A-----PRINT TIMES IN NON-STANDARD TIME UNITS. - WRITE (IOUT, 301) DELT, PERTIM, TOTIM -301 FORMAT(21X, ' TIME STEP LENGTH =', G15.6 / & - & 21X, ' STRESS PERIOD TIME =', G15.6 / & - & 21X, 'TOTAL SIMULATION TIME =', G15.6) -!C -!C2B-----RETURN - RETURN -!C -!C3------CALCULATE LENGTH OF TIME STEP & ELAPSED TIMES IN SECONDS. -100 DELSEC = CNV * DELT - TOTSEC = CNV * TOTIM - PERSEC = CNV * PERTIM -!C -!C4------CALCULATE TIMES IN MINUTES,HOURS,DAYS AND YEARS. - SIXTY = 60. - HRDAY = 24. - DAYYR = 365.25 - DELMN = DELSEC / SIXTY - DELHR = DELMN / SIXTY - DELDY = DELHR / HRDAY - DELYR = DELDY / DAYYR - TOTMN = TOTSEC / SIXTY - TOTHR = TOTMN / SIXTY - TOTDY = TOTHR / HRDAY - TOTYR = TOTDY / DAYYR - PERMN = PERSEC / SIXTY - PERHR = PERMN / SIXTY - PERDY = PERHR / HRDAY - PERYR = PERDY / DAYYR -!C -!C5------PRINT TIME STEP LENGTH AND ELAPSED TIMES IN ALL TIME UNITS. - WRITE (IOUT, 200) -200 FORMAT(19X, ' SECONDS MINUTES HOURS', 7X, & - & 'DAYS YEARS'/20X, 59('-')) - write (IOUT, 201) DELSEC, DELMN, DELHR, DELDY, DELYR -201 FORMAT(1X, ' TIME STEP LENGTH', 1P, 5G12.5) - WRITE (IOUT, 202) PERSEC, PERMN, PERHR, PERDY, PERYR -202 FORMAT(1X, 'STRESS PERIOD TIME', 1P, 5G12.5) - WRITE (IOUT, 203) TOTSEC, TOTMN, TOTHR, TOTDY, TOTYR -203 FORMAT(1X, ' TOTAL TIME', 1P, 5G12.5,/) -!C -!C6------RETURN - RETURN - END subroutine tdis_ot + real(DP) :: cnv, delsec, totsec, persec, delmn, delhr, totmn, tothr, & + totdy, totyr, permn, perhr, perdy, peryr, deldy, delyr + ! -- format + character(len=*), parameter :: fmttmsmry = "(1X, ///9X, & + &'TIME SUMMARY AT END OF TIME STEP', I5,' IN STRESS PERIOD ', I4)" + character(len=*), parameter :: fmttmstpmsg = & + &"(21X, ' TIME STEP LENGTH =', G15.6 / & + & 21X, ' STRESS PERIOD TIME =', G15.6 / & + & 21X, 'TOTAL SIMULATION TIME =', G15.6)" + character(len=*), parameter :: fmttottmmsg = & + &"(19X, ' SECONDS MINUTES HOURS', 7X, & + &'DAYS YEARS'/20X, 59('-'))" + character(len=*), parameter :: fmtdelttm = & + &"(1X, ' TIME STEP LENGTH', 1P, 5G12.5)" + character(len=*), parameter :: fmtpertm = & + &"(1X, 'STRESS PERIOD TIME', 1P, 5G12.5)" + character(len=*), parameter :: fmttottm = & + &"(1X, ' TOTAL TIME', 1P, 5G12.5,/)" + ! + ! -- Write header message for the information that follows + write (iout, fmttmsmry) kstp, kper + ! + ! -- Use time unit indicator to get factor to convert to seconds + cnv = DZERO + if (itmuni == 1) cnv = DONE + if (itmuni == 2) cnv = DSIXTY + if (itmuni == 3) cnv = DSECPERHR + if (itmuni == 4) cnv = DSECPERDY + if (itmuni == 5) cnv = DSECPERYR + ! + ! -- If FACTOR=0 then time units are non-standard + if (cnv == DZERO) then + ! -- Print times in non-standard time units + write (iout, fmttmstpmsg) delt, pertim, totim + else + ! -- Calculate length of time step & elapsed time in seconds + delsec = cnv * delt + totsec = cnv * totim + persec = cnv * pertim + ! + ! -- Calculate times in minutes, hours, days, and years + delmn = delsec / DSIXTY + delhr = delmn / DSIXTY + deldy = delhr / DHRPERDAY + delyr = deldy / DDYPERYR + totmn = totsec / DSIXTY + tothr = totmn / DSIXTY + totdy = tothr / DHRPERDAY + totyr = totdy / DDYPERYR + permn = persec / DSIXTY + perhr = permn / DSIXTY + perdy = perhr / DHRPERDAY + peryr = perdy / DDYPERYR + ! + ! -- Print time step length and elapsed times in all time units + write (iout, fmttottmmsg) + write (iout, fmtdelttm) delsec, delmn, delhr, deldy, delyr + write (iout, fmtpertm) persec, permn, perhr, perdy, peryr + write (iout, fmttottm) totsec, totmn, tothr, totdy, totyr + end if + ! + ! -- Return + return + end subroutine tdis_ot + !> @brief Deallocate memory + !< subroutine tdis_da() -! ****************************************************************************** -! tdis_da -- deallocate -! ****************************************************************************** -! -! SPECIFICATIONS: -! ------------------------------------------------------------------------------ + ! -- modules use MemoryManagerModule, only: mem_deallocate use AdaptiveTimeStepModule, only: ats_da -! ------------------------------------------------------------------------------ + ! ! -- ats if (inats > 0) call ats_da() ! @@ -498,17 +468,13 @@ subroutine tdis_da() return end subroutine tdis_da + !> @brief Read the timing discretization options + !< subroutine tdis_read_options() -! ****************************************************************************** -! tdis_read_options -- Read the options -! ****************************************************************************** -! -! SPECIFICATIONS: -! ------------------------------------------------------------------------------ + ! -- modules use ConstantsModule, only: LINELENGTH use SimModule, only: store_error use InputOutputModule, only: GetUnit, openfile - ! -- dummy ! -- local character(len=LINELENGTH) :: errmsg, keyword, fname integer(I4B) :: ierr @@ -520,7 +486,6 @@ subroutine tdis_read_options() character(len=*), parameter :: fmtdatetime0 = & &"(4x,'SIMULATION STARTING DATE AND TIME IS ',A)" !data -! ------------------------------------------------------------------------------ ! ! -- set variables itmuni = 0 @@ -600,17 +565,12 @@ subroutine tdis_read_options() return end subroutine tdis_read_options + !> @brief Read dimension NPER + !< subroutine tdis_allocate_scalars() -! ****************************************************************************** -! tdis_read_dimensions -- Read dimension NPER -! ****************************************************************************** -! -! SPECIFICATIONS: -! ------------------------------------------------------------------------------ ! -- modules use MemoryManagerModule, only: mem_allocate use ConstantsModule, only: DZERO -! ------------------------------------------------------------------------------ ! ! -- memory manager variables call mem_allocate(nper, 'NPER', 'TDIS') @@ -654,39 +614,30 @@ subroutine tdis_allocate_scalars() totalsimtime = DZERO datetime0 = '' ! - ! -- return + ! -- Return return end subroutine tdis_allocate_scalars + !> @brief Allocate tdis arrays + !< subroutine tdis_allocate_arrays() -! ****************************************************************************** -! tdis_allocate_arrays -- Allocate tdis arrays -! ****************************************************************************** -! -! SPECIFICATIONS: -! ------------------------------------------------------------------------------ ! -- modules use MemoryManagerModule, only: mem_allocate -! ------------------------------------------------------------------------------ ! call mem_allocate(perlen, nper, 'PERLEN', 'TDIS') call mem_allocate(nstp, nper, 'NSTP', 'TDIS') call mem_allocate(tsmult, nper, 'TSMULT', 'TDIS') ! - ! -- return + ! -- Return return end subroutine tdis_allocate_arrays + !> @brief Read dimension NPER + !< subroutine tdis_read_dimensions() -! ****************************************************************************** -! tdis_read_dimensions -- Read dimension NPER -! ****************************************************************************** -! -! SPECIFICATIONS: -! ------------------------------------------------------------------------------ + ! -- modules use ConstantsModule, only: LINELENGTH use SimModule, only: store_error - ! -- dummy ! -- local character(len=LINELENGTH) :: errmsg, keyword integer(I4B) :: ierr @@ -694,8 +645,6 @@ subroutine tdis_read_dimensions() ! -- formats character(len=*), parameter :: fmtnper = & "(1X,I4,' STRESS PERIOD(S) IN SIMULATION')" - !data -! ------------------------------------------------------------------------------ ! ! -- get DIMENSIONS block call parser%GetBlock('DIMENSIONS', isfound, ierr, & @@ -730,16 +679,12 @@ subroutine tdis_read_dimensions() return end subroutine tdis_read_dimensions + !> @brief Read timing information + !< subroutine tdis_read_timing() -! ****************************************************************************** -! tdis_read_timing -- Read timing information -! ****************************************************************************** -! -! SPECIFICATIONS: -! ------------------------------------------------------------------------------ + ! -- modules use ConstantsModule, only: LINELENGTH, DZERO use SimModule, only: store_error, count_errors - ! -- dummy ! -- local character(len=LINELENGTH) :: errmsg integer(I4B) :: ierr @@ -751,7 +696,6 @@ subroutine tdis_read_timing() &' MULTIPLIER FOR DELT',/1X,76('-'))" character(len=*), parameter :: fmtrow = & "(1X,I8,1PG21.7,I7,0PF25.3)" -! ------------------------------------------------------------------------------ ! ! -- get PERIODDATA block call parser%GetBlock('PERIODDATA', isfound, ierr, & @@ -789,16 +733,13 @@ subroutine tdis_read_timing() return end subroutine tdis_read_timing + !> @brief Check the tdis timing information + !! + !! Return back to tdis_read_timing if an error condition is found and let the + !! ustop routine be called there instead so the StoreErrorUnit routine can be + !! called to assign the correct file name. + !< subroutine check_tdis_timing(nper, perlen, nstp, tsmult) -! ****************************************************************************** -! check_tdis_timing -- Check the tdis timing information. Return back to -! tdis_read_timing if an error condition is found and let the ustop -! routine be called there instead so the StoreErrorUnit routine can be -! called to assign the correct file name. -! ****************************************************************************** -! -! SPECIFICATIONS: -! ------------------------------------------------------------------------------ ! -- modules use ConstantsModule, only: LINELENGTH, DZERO, DONE use SimModule, only: store_error, count_errors @@ -820,7 +761,6 @@ subroutine check_tdis_timing(nper, perlen, nstp, tsmult) character(len=*), parameter :: fmtdterror = & "('Time step length of ', G0, ' is too small in period ', I0, & &' and time step ', I0)" -! ------------------------------------------------------------------------------ ! ! -- Initialize tstart = DZERO diff --git a/src/Utilities/Constants.f90 b/src/Utilities/Constants.f90 index 784476c424f..718a9ae392c 100644 --- a/src/Utilities/Constants.f90 +++ b/src/Utilities/Constants.f90 @@ -81,6 +81,7 @@ module ConstantsModule real(DP), parameter :: DSIX = 6.0_DP !< real constant 6 real(DP), parameter :: DEIGHT = 8.0_DP !< real constant 8 real(DP), parameter :: DTEN = 1.0e1_DP !< real constant 10 + real(DP), parameter :: DSIXTY = 6.0e1_DP !< real constant 60 real(DP), parameter :: DHUNDRED = 1.0e2_DP !< real constant 100 real(DP), parameter :: DEP3 = 1.0e3_DP !< real constant 1000 @@ -92,6 +93,12 @@ module ConstantsModule real(DP), parameter :: DHDRY = -1.e30_DP !< real dry cell constant real(DP), parameter :: DNODATA = 3.0e30_DP !< real no data constant + real(DP), parameter :: DSECPERHR = 3.6e3_DP !< real constant representing number of seconds per hour (used in tdis) + real(DP), parameter :: DHRPERDAY = 2.4e1_DP !< real constant representing number of hours per day (used in tdis) + real(DP), parameter :: DDYPERYR = 3.6525e2_DP !< real constant representing the average number of days per year (used in tdis) + real(DP), parameter :: DSECPERDY = 8.64e4_DP !< real constant representing the number of seconds per day (used in tdis) + real(DP), parameter :: DSECPERYR = 3.1557600e7_DP !< real constant representing the average number of seconds per year (used in tdis) + real(DP), parameter :: DEM1 = 1.0e-1_DP !< real constant 1e-1 real(DP), parameter :: D5EM2 = 5.0e-2_DP !< real constant 5e-2 real(DP), parameter :: DEM2 = 1.0e-2_DP !< real constant 1e-2