ComputationalScienceLaboratory / ODE-Test-Problems

A MATLAB suite of initial value problems
https://computationalsciencelaboratory.github.io/ODE-Test-Problems/
MIT License
9 stars 2 forks source link

Add trivial adaptive DAE method #37

Open AndreyAPopov opened 2 years ago

AndreyAPopov commented 2 years ago

Adding an adaptive half-explicit method (one file, not over-engineered) for solving the pendulum with invariants would be very helpful.

reid-g commented 2 years ago

What order method do we want? For both Differential and Algebraic?

reid-g commented 2 years ago

Also HERK is good for problems which have a differential that can be treated explicitly. This may not be the case for all Index-2 DAEs and there may be a need for a Index-2 IRK method.

AndreyAPopov commented 2 years ago

I was literally thinking trapezoidal with Euler embedded. Call it dae12. if you have ready-made coefficients of higher order we can use them, but up to you.

I theoretically implemented enough derivatives in the pendulum problem to support a fully implicit method.

Steven-Roberts commented 2 years ago

I'm reluctant to implement integrators in OTP. To me it is slightly out of scope. Also it would need to be advanced enough to support the odeset options as that's what solve provides. I say for these DAEs of high-index, we should always aim to provide reductions to index-1 as the default RHS so built-in solvers can be used.

Steven-Roberts commented 2 years ago

But if we do want to provide an integrator, I think it should be an SDIRK method instead of a HERK. HERK requires knowledge of which variables are diff and alg, and the default RHS does not provide that.

AndreyAPopov commented 2 years ago

Ideally we would integrate with the new MATLODE, but this is a good stopgap

reid-g commented 2 years ago

Do we want to couple OTP with MATLODE?

reid-g commented 2 years ago

Also we would need to be carful with index reduced problems. Some index reductions are not stable.

reid-g commented 2 years ago

If we do provide an integrator for index-2 problems I agree it should just be an SDIRK method.

How in depth does it need to be ie how many options should it have?

Steven-Roberts commented 2 years ago

Do we want to couple OTP with MATLODE?

Couple: no. We can make them mutually compatible by having matlode use the same function signature as built in integrators. E.g. problem.solve('Solver', matlode.rk....)

Steven-Roberts commented 2 years ago

How in depth does it need to be ie how many options should it have?

It needs to accept odeset options, although it may ignore many. At the minimum is needs to support Jacobian and Mass option. Ideally it should support time and state dependant mass matrices.

AndreyAPopov commented 2 years ago

I have added a method that seems to work for the pendulum problem. There is a lot of cleanup that needs to be done, but it actually seems to work, which I am excited about.