Reproducibility when restarting generalized-alpha
Hi, I've recently been helping some co-workers with restarting PETSc time integrators from saved solution data. It looks like the only supported path for restarting the generalized-alpha integrator for 2nd-order-in-time systems (`TSALPHA2`) is to follow the same procedure as initialization, in which two first-order-accurate half-steps are used to estimate an acceleration from the given displacement and velocity. However, the resulting acceleration is not exactly equivalent to the intermediate one that would have been used by the integrator if the integration simply proceeded without restarting. This prevents exact reproducibility of computations from saved intermediate results. (An analogous issue would also affect `TSALPHA` for first-order-in-time problems, where velocity is estimated on initialization/restart.) Am I misunderstanding this, or missing a better method of restarting the 2nd-order generalized-alpha integrator? If not, would there be interest in adding an alternate initialization/restart option to the `TSALPHA2` integrator that takes a user-provided `Vec` for the initial/intermediate acceleration, and skips over the half-step estimation procedure? Thanks, David
Note that BDF has the same issue. I think the correct way to handle this is to support storing/loading these extra vectors via TSView()/TSLoad(). How are you currently restarting the simulation? Il giorno gio 2 gen 2025 alle ore 19:25 David Kamensky via petsc-users < [email protected]> ha scritto:
Hi,
I've recently been helping some co-workers with restarting PETSc time integrators from saved solution data.
It looks like the only supported path for restarting the generalized-alpha integrator for 2nd-order-in-time systems (`TSALPHA2`) is to follow the same procedure as initialization, in which two first-order-accurate half-steps are used to estimate an acceleration from the given displacement and velocity. However, the resulting acceleration is not exactly equivalent to the intermediate one that would have been used by the integrator if the integration simply proceeded without restarting. This prevents exact reproducibility of computations from saved intermediate results. (An analogous issue would also affect `TSALPHA` for first-order-in-time problems, where velocity is estimated on initialization/restart.)
Am I misunderstanding this, or missing a better method of restarting the 2nd-order generalized-alpha integrator? If not, would there be interest in adding an alternate initialization/restart option to the `TSALPHA2` integrator that takes a user-provided `Vec` for the initial/intermediate acceleration, and skips over the half-step estimation procedure?
Thanks, David
-- Stefano
How are you currently restarting the simulation?
I just reviewed the code, and we're not currently using the `TSView/Load` functions. We're just manually (de)serializing displacement, velocity, and acceleration data using a neutral format, populating PETSc `Vec`s with this data, and associating them with a new `TS` object via `TS2SetSolution` (and setting other relevant data, like time, time step size, etc.). However, `TS2SetSolution` only accepts displacement and velocity. I think the correct way to handle this is to support storing/loading these
extra vectors via TSView()/TSLoad().
I took a quick look at the implementations of `TSView/Load`, and it looks like the "base class" (to borrow some OOP terminology) implementation in `ts/interface/ts.c` only saves/loads the solution vector, while the subclass-specific logic from `TSView_Alpha` in `ts/impls/implicit/alpha/alpha2.c` only adds some additional output writing the generalized-alpha parameters to ASCII viewers (and similar for BDF). So, following the `TSView/Load` path, I don't see where it would even save/load the velocity vector for 2nd-order-in-time integrators. Is it the case that this functionality is known to be incomplete, and you're suggesting that the best path forward would be to update it? Thanks, David On Thu, Jan 2, 2025 at 10:41 AM Stefano Zampini <[email protected]> wrote:
Note that BDF has the same issue. I think the correct way to handle this is to support storing/loading these extra vectors via TSView()/TSLoad(). How are you currently restarting the simulation?
Il giorno gio 2 gen 2025 alle ore 19:25 David Kamensky via petsc-users < [email protected]> ha scritto:
Hi,
I've recently been helping some co-workers with restarting PETSc time integrators from saved solution data.
It looks like the only supported path for restarting the generalized-alpha integrator for 2nd-order-in-time systems (`TSALPHA2`) is to follow the same procedure as initialization, in which two first-order-accurate half-steps are used to estimate an acceleration from the given displacement and velocity. However, the resulting acceleration is not exactly equivalent to the intermediate one that would have been used by the integrator if the integration simply proceeded without restarting. This prevents exact reproducibility of computations from saved intermediate results. (An analogous issue would also affect `TSALPHA` for first-order-in-time problems, where velocity is estimated on initialization/restart.)
Am I misunderstanding this, or missing a better method of restarting the 2nd-order generalized-alpha integrator? If not, would there be interest in adding an alternate initialization/restart option to the `TSALPHA2` integrator that takes a user-provided `Vec` for the initial/intermediate acceleration, and skips over the half-step estimation procedure?
Thanks, David
-- Stefano
David, I think Stefano was saying the TSView/Load approach should be improved to save the additional vector(s) and use them in the restart. Are you up to trying this by adding this functionality to TSView_*/TSLoad_*, or should we try to fit in time to add this (needed) support? Barry
On Jan 2, 2025, at 2:34 PM, David Kamensky via petsc-users <[email protected]> wrote:
How are you currently restarting the simulation?
I just reviewed the code, and we're not currently using the `TSView/Load` functions. We're just manually (de)serializing displacement, velocity, and acceleration data using a neutral format, populating PETSc `Vec`s with this data, and associating them with a new `TS` object via `TS2SetSolution` (and setting other relevant data, like time, time step size, etc.). However, `TS2SetSolution` only accepts displacement and velocity.
I think the correct way to handle this is to support storing/loading these extra vectors via TSView()/TSLoad().
I took a quick look at the implementations of `TSView/Load`, and it looks like the "base class" (to borrow some OOP terminology) implementation in `ts/interface/ts.c` only saves/loads the solution vector, while the subclass-specific logic from `TSView_Alpha` in `ts/impls/implicit/alpha/alpha2.c` only adds some additional output writing the generalized-alpha parameters to ASCII viewers (and similar for BDF). So, following the `TSView/Load` path, I don't see where it would even save/load the velocity vector for 2nd-order-in-time integrators. Is it the case that this functionality is known to be incomplete, and you're suggesting that the best path forward would be to update it?
Thanks, David
On Thu, Jan 2, 2025 at 10:41 AM Stefano Zampini <[email protected] <mailto:[email protected]>> wrote:
Note that BDF has the same issue. I think the correct way to handle this is to support storing/loading these extra vectors via TSView()/TSLoad(). How are you currently restarting the simulation?
Il giorno gio 2 gen 2025 alle ore 19:25 David Kamensky via petsc-users <[email protected] <mailto:[email protected]>> ha scritto:
Hi,
I've recently been helping some co-workers with restarting PETSc time integrators from saved solution data.
It looks like the only supported path for restarting the generalized-alpha integrator for 2nd-order-in-time systems (`TSALPHA2`) is to follow the same procedure as initialization, in which two first-order-accurate half-steps are used to estimate an acceleration from the given displacement and velocity. However, the resulting acceleration is not exactly equivalent to the intermediate one that would have been used by the integrator if the integration simply proceeded without restarting. This prevents exact reproducibility of computations from saved intermediate results. (An analogous issue would also affect `TSALPHA` for first-order-in-time problems, where velocity is estimated on initialization/restart.)
Am I misunderstanding this, or missing a better method of restarting the 2nd-order generalized-alpha integrator? If not, would there be interest in adding an alternate initialization/restart option to the `TSALPHA2` integrator that takes a user-provided `Vec` for the initial/intermediate acceleration, and skips over the half-step estimation procedure?
Thanks, David
-- Stefano
Are you up to trying this by adding this functionality to TSView_*/TSLoad_*, or should we try to fit in time to add this (needed) support?
I'll have to consult with the team to decide what direction we want to go. We may want to keep our restart data in a neutral format to support other uses of it (e.g., solution postprocessing), or to maintain consistency and interoperability with some of our non-PETSc time steppers (so that we can, e.g., restart from an intermediate step of a PETSc `TS` using a non-PETSc integrator). If we do stick with our more manual re-initialization of the `TS`, it might be preferable for us to implement an initialization option that allows us to provide an acceleration. (This would be useful even for purposes other than restarting, if a user has a cheaper problem-specific method for computing the initial acceleration, e.g., inverting a mass matrix against the initial configuration's force vector.) In any case, getting consistent view/load behavior across all time integrators and testing it thoroughly would be significantly larger in scope than what we need. Best, David On Thu, Jan 2, 2025 at 2:16 PM Barry Smith <[email protected]> wrote:
David,
I think Stefano was saying the TSView/Load approach should be improved to save the additional vector(s) and use them in the restart.
Are you up to trying this by adding this functionality to TSView_*/TSLoad_*, or should we try to fit in time to add this (needed) support?
Barry
On Jan 2, 2025, at 2:34 PM, David Kamensky via petsc-users < [email protected]> wrote:
How are you currently restarting the simulation?
I just reviewed the code, and we're not currently using the `TSView/Load` functions. We're just manually (de)serializing displacement, velocity, and acceleration data using a neutral format, populating PETSc `Vec`s with this data, and associating them with a new `TS` object via `TS2SetSolution` (and setting other relevant data, like time, time step size, etc.). However, `TS2SetSolution` only accepts displacement and velocity.
I think the correct way to handle this is to support storing/loading these
extra vectors via TSView()/TSLoad().
I took a quick look at the implementations of `TSView/Load`, and it looks like the "base class" (to borrow some OOP terminology) implementation in `ts/interface/ts.c` only saves/loads the solution vector, while the subclass-specific logic from `TSView_Alpha` in `ts/impls/implicit/alpha/alpha2.c` only adds some additional output writing the generalized-alpha parameters to ASCII viewers (and similar for BDF). So, following the `TSView/Load` path, I don't see where it would even save/load the velocity vector for 2nd-order-in-time integrators. Is it the case that this functionality is known to be incomplete, and you're suggesting that the best path forward would be to update it?
Thanks, David
On Thu, Jan 2, 2025 at 10:41 AM Stefano Zampini <[email protected]> wrote:
Note that BDF has the same issue. I think the correct way to handle this is to support storing/loading these extra vectors via TSView()/TSLoad(). How are you currently restarting the simulation?
Il giorno gio 2 gen 2025 alle ore 19:25 David Kamensky via petsc-users < [email protected]> ha scritto:
Hi,
I've recently been helping some co-workers with restarting PETSc time integrators from saved solution data.
It looks like the only supported path for restarting the generalized-alpha integrator for 2nd-order-in-time systems (`TSALPHA2`) is to follow the same procedure as initialization, in which two first-order-accurate half-steps are used to estimate an acceleration from the given displacement and velocity. However, the resulting acceleration is not exactly equivalent to the intermediate one that would have been used by the integrator if the integration simply proceeded without restarting. This prevents exact reproducibility of computations from saved intermediate results. (An analogous issue would also affect `TSALPHA` for first-order-in-time problems, where velocity is estimated on initialization/restart.)
Am I misunderstanding this, or missing a better method of restarting the 2nd-order generalized-alpha integrator? If not, would there be interest in adding an alternate initialization/restart option to the `TSALPHA2` integrator that takes a user-provided `Vec` for the initial/intermediate acceleration, and skips over the half-step estimation procedure?
Thanks, David
-- Stefano
Sure, I was suggesting you might implement it for your one-needed TSView_*/TSLoad_*. With a single template done we can easily add the rest. I agree having an API to start with multiple solutions is a good idea and would be needed for any TSView_*/TSLoad_*. so that may be the simplest way for you to get started. Or we can look at providing the new API if it is out of scope for you. Barry
On Jan 2, 2025, at 5:47 PM, David Kamensky <[email protected]> wrote:
Are you up to trying this by adding this functionality to TSView_*/TSLoad_*, or should we try to fit in time to add this (needed) support?
I'll have to consult with the team to decide what direction we want to go. We may want to keep our restart data in a neutral format to support other uses of it (e.g., solution postprocessing), or to maintain consistency and interoperability with some of our non-PETSc time steppers (so that we can, e.g., restart from an intermediate step of a PETSc `TS` using a non-PETSc integrator). If we do stick with our more manual re-initialization of the `TS`, it might be preferable for us to implement an initialization option that allows us to provide an acceleration. (This would be useful even for purposes other than restarting, if a user has a cheaper problem-specific method for computing the initial acceleration, e.g., inverting a mass matrix against the initial configuration's force vector.)
In any case, getting consistent view/load behavior across all time integrators and testing it thoroughly would be significantly larger in scope than what we need.
Best, David
On Thu, Jan 2, 2025 at 2:16 PM Barry Smith <[email protected] <mailto:[email protected]>> wrote:
David,
I think Stefano was saying the TSView/Load approach should be improved to save the additional vector(s) and use them in the restart.
Are you up to trying this by adding this functionality to TSView_*/TSLoad_*, or should we try to fit in time to add this (needed) support?
Barry
On Jan 2, 2025, at 2:34 PM, David Kamensky via petsc-users <[email protected] <mailto:[email protected]>> wrote:
How are you currently restarting the simulation?
I just reviewed the code, and we're not currently using the `TSView/Load` functions. We're just manually (de)serializing displacement, velocity, and acceleration data using a neutral format, populating PETSc `Vec`s with this data, and associating them with a new `TS` object via `TS2SetSolution` (and setting other relevant data, like time, time step size, etc.). However, `TS2SetSolution` only accepts displacement and velocity.
I think the correct way to handle this is to support storing/loading these extra vectors via TSView()/TSLoad().
I took a quick look at the implementations of `TSView/Load`, and it looks like the "base class" (to borrow some OOP terminology) implementation in `ts/interface/ts.c` only saves/loads the solution vector, while the subclass-specific logic from `TSView_Alpha` in `ts/impls/implicit/alpha/alpha2.c` only adds some additional output writing the generalized-alpha parameters to ASCII viewers (and similar for BDF). So, following the `TSView/Load` path, I don't see where it would even save/load the velocity vector for 2nd-order-in-time integrators. Is it the case that this functionality is known to be incomplete, and you're suggesting that the best path forward would be to update it?
Thanks, David
On Thu, Jan 2, 2025 at 10:41 AM Stefano Zampini <[email protected] <mailto:[email protected]>> wrote:
Note that BDF has the same issue. I think the correct way to handle this is to support storing/loading these extra vectors via TSView()/TSLoad(). How are you currently restarting the simulation?
Il giorno gio 2 gen 2025 alle ore 19:25 David Kamensky via petsc-users <[email protected] <mailto:[email protected]>> ha scritto:
Hi,
I've recently been helping some co-workers with restarting PETSc time integrators from saved solution data.
It looks like the only supported path for restarting the generalized-alpha integrator for 2nd-order-in-time systems (`TSALPHA2`) is to follow the same procedure as initialization, in which two first-order-accurate half-steps are used to estimate an acceleration from the given displacement and velocity. However, the resulting acceleration is not exactly equivalent to the intermediate one that would have been used by the integrator if the integration simply proceeded without restarting. This prevents exact reproducibility of computations from saved intermediate results. (An analogous issue would also affect `TSALPHA` for first-order-in-time problems, where velocity is estimated on initialization/restart.)
Am I misunderstanding this, or missing a better method of restarting the 2nd-order generalized-alpha integrator? If not, would there be interest in adding an alternate initialization/restart option to the `TSALPHA2` integrator that takes a user-provided `Vec` for the initial/intermediate acceleration, and skips over the half-step estimation procedure?
Thanks, David
-- Stefano
Just to provide an update here, the conclusion of our internal discussion was to backlog this on our end for now, in favor of more urgent tasks. Perhaps we can open an issue for this on the PETSc issue tracker (or comment on an existing issue for BDF, if it exists), and I can check in there if/when we come back to it. Thanks, David On Thu, Jan 2, 2025 at 7:02 PM Barry Smith <[email protected]> wrote:
Sure, I was suggesting you might implement it for your one-needed TSView_*/TSLoad_*. With a single template done we can easily add the rest.
I agree having an API to start with multiple solutions is a good idea and would be needed for any TSView_*/TSLoad_*. so that may be the simplest way for you to get started. Or we can look at providing the new API if it is out of scope for you.
Barry
On Jan 2, 2025, at 5:47 PM, David Kamensky <[email protected]> wrote:
Are you up to trying this by adding this functionality to
TSView_*/TSLoad_*, or should we try to fit in time to add this (needed) support?
I'll have to consult with the team to decide what direction we want to go. We may want to keep our restart data in a neutral format to support other uses of it (e.g., solution postprocessing), or to maintain consistency and interoperability with some of our non-PETSc time steppers (so that we can, e.g., restart from an intermediate step of a PETSc `TS` using a non-PETSc integrator). If we do stick with our more manual re-initialization of the `TS`, it might be preferable for us to implement an initialization option that allows us to provide an acceleration. (This would be useful even for purposes other than restarting, if a user has a cheaper problem-specific method for computing the initial acceleration, e.g., inverting a mass matrix against the initial configuration's force vector.)
In any case, getting consistent view/load behavior across all time integrators and testing it thoroughly would be significantly larger in scope than what we need.
Best, David
On Thu, Jan 2, 2025 at 2:16 PM Barry Smith <[email protected]> wrote:
David,
I think Stefano was saying the TSView/Load approach should be improved to save the additional vector(s) and use them in the restart.
Are you up to trying this by adding this functionality to TSView_*/TSLoad_*, or should we try to fit in time to add this (needed) support?
Barry
On Jan 2, 2025, at 2:34 PM, David Kamensky via petsc-users < [email protected]> wrote:
How are you currently restarting the simulation?
I just reviewed the code, and we're not currently using the `TSView/Load` functions. We're just manually (de)serializing displacement, velocity, and acceleration data using a neutral format, populating PETSc `Vec`s with this data, and associating them with a new `TS` object via `TS2SetSolution` (and setting other relevant data, like time, time step size, etc.). However, `TS2SetSolution` only accepts displacement and velocity.
I think the correct way to handle this is to support storing/loading
these extra vectors via TSView()/TSLoad().
I took a quick look at the implementations of `TSView/Load`, and it looks like the "base class" (to borrow some OOP terminology) implementation in `ts/interface/ts.c` only saves/loads the solution vector, while the subclass-specific logic from `TSView_Alpha` in `ts/impls/implicit/alpha/alpha2.c` only adds some additional output writing the generalized-alpha parameters to ASCII viewers (and similar for BDF). So, following the `TSView/Load` path, I don't see where it would even save/load the velocity vector for 2nd-order-in-time integrators. Is it the case that this functionality is known to be incomplete, and you're suggesting that the best path forward would be to update it?
Thanks, David
On Thu, Jan 2, 2025 at 10:41 AM Stefano Zampini < [email protected]> wrote:
Note that BDF has the same issue. I think the correct way to handle this is to support storing/loading these extra vectors via TSView()/TSLoad(). How are you currently restarting the simulation?
Il giorno gio 2 gen 2025 alle ore 19:25 David Kamensky via petsc-users < [email protected]> ha scritto:
Hi,
I've recently been helping some co-workers with restarting PETSc time integrators from saved solution data.
It looks like the only supported path for restarting the generalized-alpha integrator for 2nd-order-in-time systems (`TSALPHA2`) is to follow the same procedure as initialization, in which two first-order-accurate half-steps are used to estimate an acceleration from the given displacement and velocity. However, the resulting acceleration is not exactly equivalent to the intermediate one that would have been used by the integrator if the integration simply proceeded without restarting. This prevents exact reproducibility of computations from saved intermediate results. (An analogous issue would also affect `TSALPHA` for first-order-in-time problems, where velocity is estimated on initialization/restart.)
Am I misunderstanding this, or missing a better method of restarting the 2nd-order generalized-alpha integrator? If not, would there be interest in adding an alternate initialization/restart option to the `TSALPHA2` integrator that takes a user-provided `Vec` for the initial/intermediate acceleration, and skips over the half-step estimation procedure?
Thanks, David
-- Stefano
participants (3)
-
Barry Smith -
David Kamensky -
Stefano Zampini