# Continuation of folds of periodic orbits

(c) DDE-BIFTOOL v. 3.1.1(75), 31/12/2014

The extension ddebiftool_extra_psol is able to continue local bifurcations of periodic orbits in two parameters. This demo shows how one can continue folds (saddle-nodes) of periodic orbits i nthe neuron problem. (requires running of demo1_psol.html first).

## Contents

## Add extension folder to path

The extension is installed in a separate folder. In addition, a folder with utilities for user convenience is loaded.

addpath('../../ddebiftool_extra_psol/'); addpath('../../ddebiftool_utilities/'); %#ok<*ASGLU,*NOPTS,*NASGU>

## Speed up computations by vectorization

The functions `neuron_sys_rhs` and `neuron_sys_deri` are not vectorized. In order to speed up computations we re-define `neuron_sys_rhs`, and replace `neuron_sys_seri` with the default finite-difference approximation

neuron_sys_rhs=@(xx,par)[... -par(1)*xx(1,1,:)+par(2)*tanh(xx(1,4,:))+par(3)*tanh(xx(2,3,:));.... -par(1)*xx(2,1,:)+par(2)*tanh(xx(2,4,:))+par(4)*tanh(xx(1,2,:))]; vfuncs=set_funcs(... 'sys_rhs',neuron_sys_rhs,... 'sys_tau',@()[5,6,7],... 'x_vectorized',true);

## Find initial guess

For the branch computed in demo1_psol.html one fold occured at the maximal parameter value. So we extract its index `indmax`.

[dummy,indmax]=max(arrayfun(@(x)x.parameter(ind_a21),branch5.point));

## Initialize branch and set up extended system

Then we call `SetupPOfold` to initialize the branch and to create the functions for the extended system. For the core DDE-Biftool routines the fold of periodic orbits is of the same type as a standard periodic orbit. However, the user-provided right-hand side has been extended (eg, foldfuncs.sys_rhs is different from funcs.sys_rhs). SetupPOfold has three mandatory parameters:

`funcs`containing the user-problem functions (in the demo:`vfuncs`),`branch`, the branch on which the fold occurred (in the demo:`branch5`), and`ind`, the index in`branch.point`near which the fold occurs (in the demo`indmax`).

Output (to be fed into `br_contn`):

`foldfuncs`: functions for the extended system (derived from`funcs`)`FPObranch`: branch of folds of periodic orbits with two points

The other parameters are optional name-value pairs. Important are

`'contpar'`: typically two integers ,the two continuation parameters,`'dir'`: in which direction the first step is taken), and`'step'`: the length of the first step. Note that name-value pairs are also passed on to fields of the`branch`structure: the output`FPObranch`inherits all fields from the input`branch`. Additional name-value argument pairs canbe used to change selected fields.

[foldfuncs,branch7]=SetupPOfold(vfuncs,branch5,indmax,'contpar',[ind_a21,ind_taus],... 'dir',ind_taus,'print_residual_info',1,'step',0.01,'plot_measure',[],... 'min_bound',[ind_a21,0],'max_bound',[ind_a21,4; ind_taus,10],... 'max_step',[ind_a21,0.2; ind_taus,0.5]);

it=1, res=0.000574744 it=2, res=0.0083463 it=3, res=0.000139207 it=4, res=4.6582e-08 it=5, res=2.22616e-10 it=1, res=0.000267241 it=2, res=2.96116e-10 it=1, res=0.0699482 it=2, res=0.014741 it=3, res=6.29929e-06 it=4, res=2.28572e-10 it=1, res=0.000175272 it=2, res=8.75278e-10

## Branch continuation

The output of `SetupPOfold` can be fed directly into `br_contn` to perform a branch continuation. Note that the computation is typically slower than a standard continuation of periodic orbits because the system dimension is twice as large and additional delays have been introduced.

figure(6); branch7.method.point.print_residual_info=0; branch7=br_contn(foldfuncs,branch7,100); xlabel('a21');ylabel('tau_s'); title('Continuation of fold of periodic orbits');

BR_CONTN warning: boundary hit.

## Extracting solution components

Since the fold continuation solves an extended system (with additional components) it does not make sense to compute the linear stability of the resulting orbits directly. However, `foldfuncs` has an additional field `'get_comp'`. This is a function that extacts the components of the extended system that correspond to the solution of the original system. The output of `get_comp` is an array of periodic orbits that are located precisely on the fold.

```
pf_orbits=foldfuncs.get_comp(branch7.point,'solution');
```

## Stability of fold orbits

We compute Floquet multipliers for these orbits using the utility function `GetStability(psol_array,...)`. Its optional inputs instruct it to ignore the two eigenvalues closest to 1 when computing stability. It also returns `triv_defect` which measures the distance of the supposedly trivial eigenvalues (as computed) to their correct values. `GetStability` needs the optional input `'funcs'` if its first argument does not yet have a stability structure.

[nunst_pf,dom,triv_defect,pf_orbits]=GetStability(pf_orbits,'funcs',vfuncs,... 'exclude_trivial',true,'locate_trivial',@(x)[1,1]); fprintf('max number of unstable Floquet multipliers: %d\n',max(nunst_pf)); n_orbits=length(pf_orbits); figure(16); a21_pfold=arrayfun(@(x)x.parameter(4),branch7.point); taus_pfold=arrayfun(@(x)x.parameter(7),branch7.point); plot(1:n_orbits,triv_defect); xlabel('point number'); ylabel('defect of trivial Floquet multiplier');

max number of unstable Floquet multipliers: 0

## Save results (end of tutorial demo, but try also demo1_hcli.html)

```
save('demo1_POfold_results.mat')
```