Skip to content
This repository was archived by the owner on Aug 28, 2024. It is now read-only.

Commit

Permalink
Added ability to stretch time during output steps
Browse files Browse the repository at this point in the history
  • Loading branch information
daminton committed Feb 25, 2023
1 parent 15e766f commit 65eaabf
Show file tree
Hide file tree
Showing 4 changed files with 67 additions and 8 deletions.
34 changes: 29 additions & 5 deletions python/swiftest/swiftest/simulation_class.py
Original file line number Diff line number Diff line change
Expand Up @@ -111,14 +111,20 @@ def __init__(self,read_param: bool = False,
istep_out : int, optional
The number of time steps between output saves to file. *Note*: only `istep_out` or `tstep_out` can be set.
Parameter input file equivalent: `ISTEP_OUT`
tstep_out : float, optional
The approximate time between when outputs are written to file. Passing this computes
`istep_out = floor(tstep_out/dt)`. *Note*: only `istep_out` or `tstep_out` can be set.
Parameter input file equivalent: None
nstep_out : int, optional
The total number of times that outputs are written to file. Passing this allows for a geometric progression of output steps:
`TSTART, f**0 * TSTEP_OUT, f**1 * TSTEP_OUT, f**2 * TSTEP_OUT, ..., f**(nstep_out-1) * TSTEP_OUT`,
where `f` is a factor that can stretch (or shrink) the time between outputs. Setting `nstep_out = int((tstart - tstop) / (tstep_out))` is
equivalent to the standard linear output (i.e. `f==1`) and is the same as not passing anything for this argument.
*Note*: Passing `nstep_out` requires passing either `istep_out` or `tstep_out` as well.
dump_cadence : int, optional
The number of output steps (given by `istep_out`) between when the saved data is dumped to a file. Setting it to 0
is equivalent to only dumping data to file at the end of the simulation. Default value is 10.
Parameter input file equivalent: `DUMP_CADENCE`
tstep_out : float, optional
The approximate time between when outputs are written to file. Passing this computes
`istep_out = floor(tstep_out/dt)`. *Note*: only `istep_out` or `tstep_out` can be set.
Parameter input file equivalent: None
init_cond_file_type : {"NETCDF_DOUBLE", "NETCDF_FLOAT", "ASCII"}, default "NETCDF_DOUBLE"
The file type containing initial conditions for the simulation:
* NETCDF_DOUBLE: A single initial conditions input file in NetCDF file format of type NETCDF_DOUBLE.
Expand Down Expand Up @@ -571,6 +577,7 @@ def set_simulation_time(self,
dt: float | None = None,
istep_out: int | None = None,
tstep_out: float | None = None,
nstep_out: int | None = None,
dump_cadence: int | None = None,
verbose: bool | None = None,
**kwargs: Any
Expand All @@ -588,10 +595,18 @@ def set_simulation_time(self,
dt : float, optional
The step size of the simulation. `dt` must be less than or equal to `tstop-dstart`.
istep_out : int, optional
The number of time steps between outputs to file. *Note*: only `istep_out` or `tstep_out` can be set.
The number of time steps between output saves to file. *Note*: only `istep_out` or `tstep_out` can be set.
Parameter input file equivalent: `ISTEP_OUT`
tstep_out : float, optional
The approximate time between when outputs are written to file. Passing this computes
`istep_out = floor(tstep_out/dt)`. *Note*: only `istep_out` or `tstep_out` can be set.
Parameter input file equivalent: None
nstep_out : int, optional
The total number of times that outputs are written to file. Passing this allows for a geometric progression of output steps:
`TSTART, f**0 * TSTEP_OUT, f**1 * TSTEP_OUT, f**2 * TSTEP_OUT, ..., f**(nstep_out-1) * TSTEP_OUT`,
where `f` is a factor that can stretch (or shrink) the time between outputs. Setting `nstep_out = int((tstart - tstop) / (tstep_out))` is
equivalent to the standard linear output (i.e. `f==1`) and is the same as not passing anything for this argument.
*Note*: Passing `nstep_out` requires passing either `istep_out` or `tstep_out` as well.
dump_cadence : int, optional
The number of output steps (given by `istep_out`) between when the saved data is dumped to a file. Setting it to 0
is equivalent to only dumping data to file at the end of the simulation. Default value is 10.
Expand Down Expand Up @@ -672,6 +687,12 @@ def set_simulation_time(self,

if istep_out is not None:
self.param['ISTEP_OUT'] = int(istep_out)

if nstep_out is not None:
if istep_out is None:
warnings.warn("nstep_out requires either istep_out or tstep_out to also be set", stacklevel=2)
else:
self.param['NSTEP_OUT'] = int(nstep_out)

if dump_cadence is None:
dump_cadence = self.param.pop("DUMP_CADENCE", 1)
Expand Down Expand Up @@ -714,6 +735,7 @@ def get_simulation_time(self, arg_list: str | List[str] | None = None, verbose:
"tstop": "TSTOP",
"dt": "DT",
"istep_out": "ISTEP_OUT",
"nstep_out": "NSTEP_OUT",
"dump_cadence": "DUMP_CADENCE",
}

Expand All @@ -723,6 +745,7 @@ def get_simulation_time(self, arg_list: str | List[str] | None = None, verbose:
"dt": self.TU_name,
"tstep_out": self.TU_name,
"istep_out": "",
"nstep_out": "",
"dump_cadence": ""}

tstep_out = None
Expand Down Expand Up @@ -772,6 +795,7 @@ def set_parameter(self, verbose: bool = True, **kwargs):
"dt": None,
"istep_out": 1,
"tstep_out": None,
"nstep_out": None,
"dump_cadence": 10,
"init_cond_file_type": "NETCDF_DOUBLE",
"init_cond_file_name": None,
Expand Down
2 changes: 2 additions & 0 deletions src/base/base_module.f90
Original file line number Diff line number Diff line change
Expand Up @@ -35,6 +35,8 @@ module base
character(STRMAX) :: in_type = "NETCDF_DOUBLE" !! Data representation type of input data files
character(STRMAX) :: in_form = "XV" !! Format of input data files ("EL" or ["XV"])
integer(I4B) :: istep_out = -1 !! Number of time steps between saved outputs
integer(I4B) :: nstep_out = -1 !! Total number of saved outputs
real(DP) :: fstep_out = 1.0_DP !! The output step time stretching factor
character(STRMAX) :: outfile = BIN_OUTFILE !! Name of output binary file
character(STRMAX) :: out_type = "NETCDF_DOUBLE" !! Binary format of output file
character(STRMAX) :: out_form = "XVEL" !! Data to write to output file
Expand Down
6 changes: 3 additions & 3 deletions src/misc/solver_module.f90
Original file line number Diff line number Diff line change
Expand Up @@ -290,7 +290,7 @@ end function f
integer(I4B) :: i,j, ev,br, niter, error
integer(I4B),parameter :: NTRY = 50
integer(I4B),parameter :: NBRACKET =20
real(DP),parameter :: FIRSTFACTOR = 1.6_DP
real(DP),parameter :: FIRSTFACTOR = 1.1_DP
real(DP) :: f1, f2, fmin
integer(I4B) :: numev,numbr

Expand All @@ -313,7 +313,7 @@ end function f
bracket: do br = 1, NBRACKET
numbr = numbr + 1
x1 = startx1
x2 = startx1 + abs(startx1 * factor)
x2 = startx1 + abs(startx1 * (1.0_DP - factor))

! First bracket the root
f1 = f(x1)
Expand All @@ -339,7 +339,7 @@ end function f
end if
end do
x1 = x2
factor = 0.5_DP * (factor + 1._DP)
factor = factor + 0.5_DP * (factor - 1._DP)
end do bracket

! Now do a Brent's method to find the root
Expand Down
33 changes: 33 additions & 0 deletions src/swiftest/swiftest_io.f90
Original file line number Diff line number Diff line change
Expand Up @@ -1883,6 +1883,7 @@ module subroutine swiftest_io_param_reader(self, unit, iotype, v_list, iostat, i
integer(I4B) :: nseeds, nseeds_from_file
logical :: seed_set = .false. !! Is the random seed set in the input file?
character(len=:), allocatable :: integrator
real(DP) :: tratio, y


! Parse the file line by line, extracting tokens then matching them up with known parameters if possible
Expand Down Expand Up @@ -1931,6 +1932,8 @@ module subroutine swiftest_io_param_reader(self, unit, iotype, v_list, iostat, i
param%in_form = param_value
case ("ISTEP_OUT")
read(param_value, *) param%istep_out
case ("NSTEP_OUT")
read(param_value, *) param%nstep_out
case ("BIN_OUT")
param%outfile = param_value
case ("OUT_TYPE")
Expand Down Expand Up @@ -2126,6 +2129,15 @@ module subroutine swiftest_io_param_reader(self, unit, iotype, v_list, iostat, i
iostat = -1
return
end if
if (param%nstep_out <= 0) then
param%nstep_out = int((param%tstop - param%t0) / (param%istep_out * param%dt))
param%fstep_out = 1.0_DP ! Linear output time
else
param%fstep_out = 1._DP
tratio = (param%TSTOP - param%T0) / (param%istep_out * param%dt)
y = time_stretcher(param%fstep_out)
call solve_roots(time_stretcher,param%fstep_out)
end if
if (param%dump_cadence < 0) then
write(iomsg,*) 'Invalid DUMP_CADENCE. Must be a positive integer or 0.'
iostat = -1
Expand Down Expand Up @@ -2332,6 +2344,27 @@ module subroutine swiftest_io_param_reader(self, unit, iotype, v_list, iostat, i
return
667 continue
write(*,*) "Error reading param file: ", trim(adjustl(iomsg))

contains
function time_stretcher(fstep_out) result(ans)
!! author: David A. Minton
!!
!! Equation for the time stretchinf function. Solving the roots of this equation gives the time stretching factor for non-linear file output cadence.
implicit none
! Arguments
real(DP), intent(in) :: fstep_out
! Result
real(DP) :: ans

if (abs(fstep_out-1.0_DP) < epsilon(1.0_DP)) then
ans = self%nstep_out - tratio
else
ans = (1.0_DP - fstep_out**(self%nstep_out))/ (1.0_DP - fstep_out) - tratio
end if

return
end function time_stretcher

end subroutine swiftest_io_param_reader


Expand Down

0 comments on commit 65eaabf

Please sign in to comment.