Update with BifurcationKit

[TorkelE]

This updates the bifurcation diagram tutorial to work with the latest BifurcationKit version. It also uses a slightly less powerful (and thus simpler) approach.

[TorkelE]

(Need to add an updated figure before we can merge)

[isaacsas]

Could you convert this to be a dynamic tutorial that uses @example blocks? (See some of the others.)

That would cause it to run when docs are built, and let us see immediately if it breaks.

(Also, maybe switch to using the preferred u0 / p style with mappings?)

[TorkelE]

Could you convert this to be a dynamic tutorial that uses @example blocks? (See some of the others.

Sounds good, will do.

I could use the mapping notation, however, I would then have to do e.g. first.() to extract a vector which made it all feel a bit silly. Ideally, at some point, it would be good if we could make the two packages more compatible, including permitting BifurcationKit use mapping input.

[isaacsas]

Right, but the point is even using first.(pmap) is implicitly then assuming you've defined the parameters with the correct ordering, and we no longer provide any guarantees about that. Calling varmap_to_vars on a symbolic mapping will give you the correct order. Alternatively you could just create an ODEProblem and extract p, u0 and the ODEFunction from that (instead of creating an ODEFunction directly).

[TorkelE]

Going via an ODEProblem seems like a good way to do it, agreed.

[isaacsas]

Yeah, that should probably be the way we suggest people get all the "stuff" out they need to interface with non-MT packages.

[isaacsas]

Looks like there are some issues in running the examples -- the output doesn't look all there:

https://catalyst.sciml.ai/previews/PR539/tutorials/bifurcation_diagram/

Have you tried building the docs locally to double check them?

[TorkelE]

I was unsure how to build them locally (I just ran the code), so I just used the preview function https://catalyst.sciml.ai/previews/PR539/tutorials/bifurcation_diagram/.

To me it seems that it is just the display of the actual diagram that fails:

The 3 before that doesn't display anything, but given that it would just be an incomprehensible mess, that is probably very much desired...

[isaacsas]

The final four example boxes have no output, so I think this is an issue.

You can build them locally by going into the docs folder, activating the environment there, resolving it, and then including the make.jl file.

[TorkelE]

Right, BifurcationKit wasn't available for the doc to build. I added it to the docs/Project.toml file.

[TorkelE]

Ok, some minor thing still needs fixing, will do that when I get home tonight.

[isaacsas]

Thanks for taking the time to make this dynamic. It will make it much easier to see when package updates, or Catalyst updates, break things going forward!

[TorkelE]

Think this is done now, unless something happens I'll merge by tomorrow.

[isaacsas]

If you make this a Dict then won't your p_bstart calculation become trivial? ModelingToolkit should support Dicts for the mappings too.

[isaacsas]

Just the one comment, maybe make the parameter mapping a dictionary. I think it would simplify the tutorial since it would avoid having to use setindex! or findfirst (if I've understood what you are doing).

[TorkelE]

Generally I think that would make sense, however, when I do:

bprob = BifurcationProblem(F, oprob.u0, oprob.p, (@lens _[findfirst(first.(p).==bif_par)]); recordFromSolution = (x, p) -> x[findfirst(first.(u0).==plot_var)], J=J)

I must supply the index of the bifurcation parameter and plotting variable in the resulting array, and I am not sure how to do this robustly, nor if it is even possible.

[isaacsas]

You can use https://mtk.sciml.ai/dev/basics/FAQ/#Getting-the-index-for-a-symbol to see the recommended way to get a symbolic variable's index. For the species replace parameters(sys) with states(sys).

[TorkelE]

How to I access this function? I got the latest version of MTK, but trying indexof, ModelingToolkit.indexof, or similar all fails. (I also tried importing Symbolics and Catalyst and using these).

[isaacsas]

You have to define that function as in the FAQ. It isn't pre-defined by ModelingToolkit unfortunately.

[TorkelE]

Ahh, I understand.

[TorkelE]

Those seems to only work well if you actually have declared the vars/pars as symbolic. After some tinkering the best I got was:

indexof(sym,syms) = findfirst(isequal(sym),syms)
par_idx = indexof(:S,Symbol.(ModelingToolkit.parameters(rn)))
var_idx = indexof(:X,Symbol.(getfield.(ModelingToolkit.states(rn),:f)))

but I think this will just make things messier than what we already have.

[isaacsas]

Just @unpack the symbols from the system to load them as symbolic variables.

[isaacsas]

If you want I can update it to handle this part tonight?

[TorkelE]

That might be better, yeah, I am still unsure exactly how we want this. Loads of thanks!