Skip to content

Feature: Higher-order Bautin normal form and unfolding#254

Open
mmbosschaert wants to merge 6 commits into
bifurcationkit:masterfrom
mmbosschaert:master
Open

Feature: Higher-order Bautin normal form and unfolding#254
mmbosschaert wants to merge 6 commits into
bifurcationkit:masterfrom
mmbosschaert:master

Conversation

@mmbosschaert

Copy link
Copy Markdown

Higher-Order Bautin Normal Form and Unfolding

This PR implements higher-order computation of the Bautin (Generalized Hopf) normal form and its associated predictor for fold of periodic orbits continuation.

Summary

The Bautin bifurcation occurs when the first Lyapunov coefficient vanishes at a Hopf point, requiring higher-order analysis. Previously, the code computed up to the second Lyapunov coefficient. This PR extends the computation to the third Lyapunov coefficient and implements a higher-order predictor that uses additional higher-order parameter-dependent terms for more accurate approximation to the fold of period orbits curve emanating from the Bautin bifurcation point.

This test validates the higher-order Bautin normal form computation and predictor:
- Uses Symbolics.jl to compute jet derivatives up to 7th order
- Tests Stuart-Landau equation with quintic and septic nonlinearities
- Verifies l2 and l3 coefficients match theoretical formulas (4·c5 and 8·c7)
- Confirms predictor actually uses c₃ (l3) by comparing predictions with/without it
- Shows parameter predictions differ by ~1e-7 when c₃ is included vs zeroed

The test proves the higher-order predictor is working correctly.
@rveltz

rveltz commented Oct 14, 2025

Copy link
Copy Markdown
Member

Thanks a lot for this massive PR! I will review this asap.

Best regards

Comment thread src/codim2/NormalForms.jl
~cv && @debug "[Bautin H3100] Linear solver for J did not converge. it = $it"

h2200 = D40(q0, q0, cq0, cq0) .+
4 .* C(q0, cq0, H1100) .+ C(cq0, cq0, H2000) .+ C(q0, q0, conj.(H2000)) .+

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

align with = sign (4 spaces)

Comment thread src/codim2/NormalForms.jl Outdated

G32 = dot(p0, D50(q0, q0, q0, cq0, cq0))
G32 += (dot(p0, D40(q0, q0, q0, conj.(H2000))) +
3dot(p0, D40(q0, cq0, cq0, H2000)) +

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

idem as above

@rveltz

rveltz commented Feb 4, 2026

Copy link
Copy Markdown
Member

You should modify the documentation as well https://github.com/bifurcationkit/BifurcationKitDocs.jl/blob/main/docs/src/bautin.md

You can probably add your name and a link to this PR there to underline you contribution

@rveltz

rveltz commented Feb 4, 2026

Copy link
Copy Markdown
Member

It seems your code only works for 2d systems as I have tried locally. I suggest minor changes in test/test_bautin_higherorder.jl to verify this claim.

@rveltz

rveltz commented Feb 5, 2026

Copy link
Copy Markdown
Member

If you want some help, I need to be able to push to your PR. This is possible as explained here

(4) On user-owned forks, if you want to allow anyone with push access to the upstream repository to make changes to your pull request, select Allow edits from maintainers.

@mmbosschaert

Copy link
Copy Markdown
Author

Thanks. It looks like ‘Allow edits and access to secrets by maintainers’ is already checked. Am I missing something?

@rveltz

rveltz commented Feb 26, 2026

Copy link
Copy Markdown
Member

I am not sure I am a maintener on your branch, am I?

Maikel Bosschaert added 2 commits May 21, 2026 11:59
Resolve conflicts in NormalForms.jl, periodicorbit/NormalForms.jl, Problems.jl,
and codim2/NormalForms.jl.

Widen _field_jet/_type_jet to include R20..R70 and R41/R51. The higher-order
Bautin code reads these jet fields directly and they were outside upstream's
constraint.

Route a user-supplied R01 function into jet.R01. Upstream made R01 an explicit
BifurcationProblem kwarg defaulting to AutoDiff(), which captured user
functions before they reached the jet. Only the dispatch when R01 is neither
AutoDiff() nor FiniteDifferences() changes.

Track upstream renames in the higher-order branch: bare dot to LA.dot, I to
LA.I, and alias prob_ma = 𝐌𝐚.
Conflict in src/periodicorbit/PeriodicOrbitCollocation.jl: upstream split
po_analytical_jacobian! into a thin wrapper plus inner
_collocation_analytical_jacobian! taking the time-slice matrix and period
directly, and added an optional _compute_borders. The higher-order Bautin
branch never touched this function, so the resolution takes upstream's
refactor verbatim.
@mmbosschaert

Copy link
Copy Markdown
Author

Is there something I can do for the PR to be approved at this stage?

@rveltz

rveltz commented Jun 6, 2026

Copy link
Copy Markdown
Member

Hi,

Have you adressed the comments/review in the code?
The main concern is that it only works for 2d system at the moment

Rotate the padded Stuart-Landau system by an orthogonal Q in the (1,3) plane so
the third coordinate is used while l2 and l3 keep their 2D values 4*c5 and 8*c7.
The test checks the eigenvector has a non-zero third component to make sure the
computation is not restricted to planar systems.
@mmbosschaert

Copy link
Copy Markdown
Author

I've added an example demonstrating that the normalization method also works in a higher dimensional system.

@mmbosschaert

Copy link
Copy Markdown
Author

Hi, have you been able to verify that the predictor works in higher-dimensional systems? If so, is there still work for me to do on the PR before it can be merged?

@rveltz

rveltz commented Jun 30, 2026

Copy link
Copy Markdown
Member

Hi,

I just pushed a few refactor but there is a conflict. Very strange...

Adapt the higher-order Bautin normal form and its test to the unified
codim-2 bifurcation point struct.
@mmbosschaert

Copy link
Copy Markdown
Author

Hi,

It should be fixed now. Please let me know if I can do anything else to help with the merge.

@rveltz

rveltz commented Jun 30, 2026

Copy link
Copy Markdown
Member

That is strange, I am not sure my changes have been pushed. Actually yes, it is here https://github.com/bifurcationkit/BifurcationKit.jl/blob/bautin-local/test/normal_forms/test_bautin_higherorder.jl

But it does not appear in the commits of this PR.

@rveltz

rveltz commented Jun 30, 2026

Copy link
Copy Markdown
Member

Please let me know if I can do anything else to help with the merge.

I need to read your code in detail.

@rveltz

rveltz commented Jun 30, 2026

Copy link
Copy Markdown
Member

I think there is a big issue:

julia> 1 +
       1
2

BUT

julia> 1 
       +1
1

this is used everywhere in your implementation

@mmbosschaert

Copy link
Copy Markdown
Author

Note that every multi-line sum is wrapped in parentheses, e.g., h4300 = ( ... ), and inside parentheses, Julia does continue across newlines, so no terms are dropped:

julia> ( 1
+ 1 )
2

The higher-order test relies on this: l3 is computed from h4300 (the largest of these formulas) and
comes out exact, l3 = 9.6 (and l2 = 1.2). Dropping a term there would fail the test.

The leading-operator style is fragile since it leans on the parentheses, so happy to move the operators
to the end of each line if you'd prefer.

@rveltz

rveltz commented Jul 1, 2026

Copy link
Copy Markdown
Member

Note that every multi-line sum is wrapped in parentheses

good catch

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

None yet

Projects

None yet

Development

Successfully merging this pull request may close these issues.

2 participants