Revisiting Boundary Value Problems

[avik-pal]

With the recent developments around BoundaryValueDiffEq.jl and adjoint sensitivity computation for BVPs, there seem to be some clear shortcomings in how these problems and functions are defined.

Proposals

1. TwoPointBVPFunction: Should take 2 residuals, residuala to store the residual dependent on ua and similarly for ub

   • One thing that isn't very clear to me rn, is how do we initialize the separate residuals. Do we ask the user for len(residuala) & len(residualb)? If so, then we can simply slice the residual at a specific point and still maintain the current API.
   • Why do this? TwoPointBVPs have a specific residual structure if we know the exact bc dependence on ua and ub. So we won't have to do sparsity detection and can manually reorder the Jacobian Matrix to provide a fast path compared to standard BVPs
   • To smoothen the transition, should we check that if the bc has no 5 argument function, we construct a regular BVProblem?

2. Introduce a BVPFunction. This allows us to provide things like jac_prototype and even jac and similar stuff.

   • This is of particular interest if we want to support overconstrained and underconstrained BVPs. We can always use ODEFunction but it might be better to have a specific function for this.

3. Allow bc to be out of place (if we want to use AD through it peacefully). I see two options here:

   • If f is IIP then bc has to be IIP and similar for OOP f.
   • Allow f and bc to be independent, which means we will have to do more isinplace checks but also provides a smoother end-user experience.

Implemented Features

• Updated the TwoPointBVPFunction

cc @ChrisRackauckas @ErikQQY @EthanDecleyn

[ai-maintainer]

AI-Maintainer Review for PR - Revisiting Boundary Value Problems

Title and Description 👍

The Title and description are clear, concise and helpful

The title and description of the pull request are clear and concise. They effectively communicate the purpose of the changes, which is to revisit and improve the handling of boundary value problems. The description provides a detailed explanation of the shortcomings of the current implementation and proposes specific improvements to address them. It also mentions the implemented features related to the TwoPointBVPFunction.

Scope of Changes 👍

The changes are narrowly focused and address a specific set of issues

The changes proposed in this pull request are focused on addressing specific issues related to boundary value problems. The proposed improvements and implemented features are all related to enhancing the definition and functionality of boundary value problems. The author is not trying to tackle multiple unrelated problems simultaneously, which is a good practice.

Documentation 👎

Some newly added functions and methods lack docstrings

The following newly added functions and methods do not have docstrings:

• TwoPointBVPFunction (constructor)
• (f::TwoPointBVPFunction)(residuala, residualb, ua, ub, p) method
• (f::TwoPointBVPFunction)(residual::Tuple, u, p) method

These additions would benefit from having docstrings to describe their behavior, arguments, and return values. Please add these docstrings to improve code readability and maintainability.

Testing 👎

The description does not mention how the changes were tested

The description does not provide information about how the changes were tested. It would be beneficial to include details about the testing approach and any test cases that were used to validate the changes. This would provide transparency and ensure the changes have been adequately tested.

Suggested Changes

• Please add docstrings to the newly added functions and methods to describe their behavior, arguments, and return values. This will improve code readability and maintainability.
• Please provide information about how the changes were tested. This could include the testing approach and any test cases that were used to validate the changes.

Reviewed with AI Maintainer

[codecov]

Codecov Report

Merging #477 (593e7f0) into master (eb92995) will decrease coverage by 57.66%.
The diff coverage is 0.00%.

@@            Coverage Diff             @@
##           master    #477       +/-   ##
==========================================
- Coverage   57.65%   0.00%   -57.66%     
==========================================
  Files          50      50               
  Lines        3710    3635       -75     
==========================================
- Hits         2139       0     -2139     
- Misses       1571    3635     +2064     

Files Changed	Coverage Δ
src/problems/bvp_problems.jl	0.00% <0.00%> (-44.00%)	⬇️
src/scimlfunctions.jl	0.00% <0.00%> (-61.34%)	⬇️

... and 45 files with indirect coverage changes

📣 We’re building smart automated test selection to slash your CI/CD build times. Learn more

[ChrisRackauckas]

If f is IIP then bc has to be IIP and similar for OOP f.
Allow f and bc to be independent, which means we will have to do more isinplace checks but also provides a smoother end-user experience.

The first. We require matching in all other instances (for example, SDEs).

This is of particular interest if we want to support overconstrained and underconstrained BVPs. We can always use ODEFunction but it might be better to have a specific function for this.

It's just been a waiting for a BVPFunction. Yes we should do this. The reason why it was waiting though was to confirm whether we wanted two different jac and jac_prototypes, i.e. one also for bcjac and bcjac_prototype separate from the other one. I think we want jac to be the jac of f and bcjac as the jac of bc, and then the solver will make do with that information. W and W_prototype could be reserved for potentially the "global Jac" i.e. the jacobian of the nonlinear solve, which we would probably not expect most people to touch.

TwoPointBVPFunction: Should take 2 residuals, residuala to store the residual dependent on ua and similarly for ub

Well justified.

One thing that isn't very clear to me rn, is how do we initialize the separate residuals. Do we ask the user for len(residuala) & len(residualb)? If so, then we can simply slice the residual at a specific point and still maintain the current API.

Length assumes Array. We should ask for prototypes.

To smoothen the transition, should we check that if the bc has no 5 argument function, we construct a regular BVProblem?

Yeah having a transition path like that is good. That will keep old codes working, if unoptimized (but it wasn't well-optimized before anyways so whatever).

[ErikQQY]

I want to add a few more things here:

1. I implemented a draft of BVPFunction before: Add BVPFunction, I can add more features if we need.

2. Recently I am looking at how can we supply initial guess to the construction of BVProblem and TwoPointBVProblem(e.g. we pass some solution from ODE solving as the initial guess to our BVP problem, and use these initial guess as our initial state to avoid extensive computing from scratch), and I found the current ways of passing initial guess to BVP problem constructing is difficult to use. To be specific, let's see here:

SciMLBase.jl/src/problems/bvp_problems.jl

Lines 117 to 122 in cce1ce0

	# convenience interfaces:
	# Allow any previous timeseries solution
	function BVProblem(f::AbstractODEFunction, bc, sol::T, tspan::Tuple, p = NullParameters();
	kwargs...) where {T <: AbstractTimeseriesSolution}
	BVProblem(f, bc, sol.u, tspan, p)
	end

SciMLBase.jl/src/problems/bvp_problems.jl

Lines 153 to 160 in cce1ce0

	# Allow previous timeseries solution
	function TwoPointBVProblem(f::AbstractODEFunction,
	bc,
	sol::T,
	tspan::Tuple,
	p = NullParameters()) where {T <: AbstractTimeseriesSolution}
	TwoPointBVProblem(f, bc, sol.u, tspan, p)
	end

I think maybe we also need another property in BVProblem and TwoPointBVProblem to store the initial guess?

[ChrisRackauckas]

I think this looks good. Failing test though? I want this to be merged at the same time as #497

[avik-pal]

Give me a few hours, I am finishing up the ODEInterface integration to confirm we have everything we need. The tests here also need updating due to the breaking change

[avik-pal]

If this tests, pass we can bump it to 2.0 and then merge. Testing with 2.0 doesn't work due to version conflicts while installing StochasticDiffEq and others