Is there a julia package that uses Newton-Krylov method for systems of nonlinear equations?

[ChrisRackauckas]

@rveltz

The closest I can think of are methods for the computation of periodic orbits as BVP. I do not use the solvers of DE for that, except for the shooting algorithm where I call DE. These BVPs are quite specific (you have to add a constraint otherwise they are not well posed) so that’s why I dont call DE. I also use a lot of tricks to make the computation faster and the solvers are optimized for large scale problems. When the periodic orbit has a large period (close to homoclinic point), the methods show their weaknesses.

What are the interesting features of this BVP? It would probably be a good thing to have BVP solvers that specialize on this if it comes up a lot.

[rveltz]

There are like constrained BVP. You have to add a constraint to fix the phase of the periodic orbit (PO). This is explained here; here for the FD method, for the shooting here.

I tend to think that this is quite specific to the problem of computing PO but I may be wrong.

The standard shooting method (one of the methods here) is straightforward to implement using DE but useless if you cant pass a preconditionner.

In fact, the "only’’ difficulty for computing PO in large dimensions (with FD methods or Shooting ones) is to be able to invert cyclic-by blocks linear systems like those and they are notoriously difficult to solve in general. It is the main reason why computing PO in large dimensions is difficult.

[ChrisRackauckas]

@rveltz

There are like constrained BVP. You have to add a constraint to fix the phase of the periodic orbit (PO).

I tend to think that this is quite specific to the problem of computing PO but I may be wrong.

Ahh yeah, it’s a Boundary Value Differential-Algebraic Equation. I plan to get around to those, but indeed doing the general form in a very good way is fairly involved. The generalization of that finite difference type of method is the MIRK tableau form, but we need to work on adaptivity and stuff there. One of the things that was holding back this work though was that I didn’t get enough requests/interest from users about BVPs, but now that there’s a clear use case it’ll get a priority bump. Still probably a year or two away from something you’d really want to see though, but maybe we can collaborate a bit here. Though indeed, this is the kind of thing where it’ll need to specialize on a bunch of matrix forms and such in order to be good, which is exactly why it would take a bunch of work.

[rveltz]

@ChrisRackauckas

The generalization of that finite difference type of method is the MIRK tableau form,

The reason for the use of FD method is because the linear system is not too filled in and we can precondition it. Otherwise, I have no clue if it is possible.

One thing FD lacks is mesh adaptation, that’s a big thing if the PO has large period or a strange structure. I’ll have a look a the bvp mirk4 method.

[ChrisRackauckas]

@rveltz

The reason for the use of FD method is because the linear system is not too filled in and we can precondition it. Otherwise, I have no clue if it is possible.

Yes, we need to add a way to do preconditioning on it. Though IIRC a two-point BVP gets you a banded matrix? So if it’s a two-point BVP then you want to do a banded direct solver for sure, and then if you have different conditions you fallback to generic linear solver methods with preconditioners and all of that, and that’s what we haven’t exposed well in our BVP solvers.

[rveltz]

Not for PO. The boundary condition x(0) = x(1) makes the jacobian cyclic and not banded. That’s where the difficulty lies.