Examples for `CVodeResizeHistory` when `y` changes size

[martinwang2002]

CVodeResizeHistory is a very nice function and I am glad it was released in v7.3.0. I'm working on a project that requires dynamically resizing the ODE system during integration using CVodeResizeHistory(). While the documentation explains the API, I’m unclear on how to handle some scenarios.

1. Resizing when the solution vector y changes size

• When the ODE system increases or decreases in size, how should y_hist and f_hist be constructed to accommodate the new variables?

• For example, consider this case:

  Time	t_n	t_{n-1}	t_{n-2}	t_{n-3}
  Number of equations	6	5	4	5
  Notes	y6 added	shrinks below original	intermediate (y1) or last (y5) removed	back to original size

• What should the sizes of y_hist and f_hist be in this case?
  At t_n, should the sizes of the vectors in history remain as [6, 5, 4, 5], or should all history vectors be resized to [6, 6, 6, 6] to match the new system size?

• If resizing history vectors, how should the new components be filled? Is zero-initialization acceptable, or should they be extrapolated or estimated in a specific way?

---

2. Clarification on test behavior

In this test case, I noticed two different behaviors:

sundials/test/unit_tests/cvode/CXX_serial/cv_test_resize_history.cpp

Lines 127 to 129 in 967f4d9

	// resize = 0 -- do not resize
	// resize = 1 -- call resize but with the same problem size
	// resize = 2 -- grow the problem by one element each time step

• Resizing while keeping the same problem size
• Growing the problem by one element each time step

Could you clarify the purpose of each approach? How do they differ in terms of results, stability, or accuracy?

---

3. Request for a complete working example

It would be incredibly helpful to see a full example that demonstrates a dynamic resize where the number of unknowns (NEQ) both expands and shrinks

---

Thank you for your work on this excellent integration library!

[gardner48]

@martinwang2002 could you describe the problem your solving? Depending on the application, using a one step method from ARKODE maybe a better fit as you only need to supply the current state at the new size instead of a history of values (see ARKodeResize).

1. Resizing when the solution vector y changes size

All the vectors supplied to CVodeResizeHistory must have the same length and, when adding entries, the new values need to be a "sufficiently accurate" approximation of the solution at the corresponding time in the history. How the values should be computed will depend on the application. For example, using interpolation if the unknowns are on a spatial mesh (however, one-step methods are typically applied in these cases). Other applications (e.g., solving the Kadanoff-Baym equations) can use symmetry properties of the system to obtain the needed history values.

2. Clarification on test behavior

The purpose of this test it to check that the computations within CVodeResizeHistory are performed correctly and it uses a problem with an analytic solution so history values can be easily computed for testing. So this is not necessarily a good example of how this function would be used in practice. The case that "resizes" the problem to the same system size should produce approximately the same solution as solving without resizing. These will not be exactly the same as resizing requires an additional evaluation the problem right-hand size function instead of approximating it using the correction vector. The other case adds a new copy of the problem each time step so each entry should have the same solution over time regardless of when it was added. That is, entry integrated from the very beginning and one added later should be essentially the same.

3. Request for a complete working example

Indeed this would be helpful and we'll try to add one in the near future (e.g., adapting the ARKODE example with resizing).

[martinwang2002]

@gardner48 Many thanks for your reply.

I'm working on atomic transitions of electrons, where the transition rates are governed by differential equations (df/dt). As temperature or potential changes, certain atomic states become inaccessible, reducing the number of equations (NEQ). Conversely, if the system is driven by a laser or other external sources, some atomic states may become populated, increasing the NEQ.

In such scenarios, the removed or added equations are not necessarily at the end of the system, so I need to understand how to correctly provide the N_Vector *y_hist and N_Vector *f_hist inputs to CVodeResizeHistory(). Specifically, do I need to fill zeros for the newly added equations, and how should I map the history vectors when the ordering of the system changes?

Currently, we handle this problem changes by fully reinitializing the CVODE memory, which unfortunately discards previous information about variable correlations and significantly slows down convergence.

---

Indeed, this would be helpful and we’ll try to add one in the near future (e.g., adapting the ARKODE example with resizing).

Having such an example would be incredibly helpful, particularly in understanding how to supply the inputs to CVodeResizeHistory() and maintain history consistency across a changing system size.

[gardner48]

Hi @martinwang2002, if inaccessible states become accessible again later in the simulation or new states are added during the simulation it might not be necessary to remove/add unknowns from the system but instead use a fixed state size (assuming you know the maximum number of possible states) and formulate the right-hand side function to account for the change in behavior. Depending on how the transition from accessible to inaccessible (or vice versa) is formulated, this could introduce is discontinuity. See this page and this page for more on handling discontinuities in the right-hand side.

If it is necessary to change the system size, there are no restrictions on how the unknowns are ordered when calling CVodeResizeHistory. The only requirement is the N_Vector *y_hist and N_Vector *f_hist arrays need to use a self-consistent ordering. So if you are inserting entries in an arbitrary location, the corresponding locations in the history vectors will need to be updated accordingly. Depending on the order of the method and the behavior of the system, filling values with zero may or may not be appropriate. For example if four history vectors are required when resizing, and a state exists in step n-4, is removed in n-3, does not exist in n-2, and is added back in n-1, then the state at n-4 for that unknown should be used in the last history vector.

[martinwang2002]

Hi @gardner48 , I quickly put together a multidimensional Gaussian integration (without correlation) to test the ideas.

The N_Vector *y_hist and N_Vector *f_hist arrays are provided self-consistently (via CVodeGetDky) and initialized with zeros when newly added:
https://github.com/martinwang2002/CVODE-ResizeHistory/blob/a4558cff146a9ed7fc9049d2d1e7f32df6c8a0e9/src/main.cpp#L183-L210

Does this represent the correct workflow, in the case of CV_NORMAL integration?

Also, I completely agree with your point about the potential discontinuity in both the solution y and the right-hand side function y'.