Pendulum with algebraic differential equation

Contents

Consider a pendulum subject to external harmonic moment, which is governed by

$\ddot{\theta}+c\dot{\theta}+\sin\theta=\epsilon\cos\Omega t$

Define $x_1=\theta$, $x_2=\dot{\theta}$, $x_3=\sin\theta$ and $x_4=\cos\theta-1$, we have

$\dot{x}_1=x_2,\\\dot{x}_2=-cx_2-x_3+\epsilon\cos\Omega t,\\\dot{x}_3=(1+x_4)x_2,\\x_3^2+(1+x_4)^2=1.$

Lyapunov subcenter manifold

clear all
c=0.001;
B = diag([1 1 1 0]);
A = [0 1 0 0;0 -c -1 0;0 1 0 0;0 0 0 2];
subs2 = [3 2 4
    4 3 3
    4 4 4];
vals2 = [1 1 1]';
F2 = sptensor(subs2,vals2,[4 4 4]);

Setup model

order = 1;
DS = DynamicalSystem(order);
set(DS,'B',B,'A',A,'fnl',{F2});
set(DS.Options,'Emax',5,'Nmax',10,'notation','multiindex')

Linear Modal analysis

[V,D,W] = DS.linear_spectral_analysis();

1 nan/inf eigenvalues are removed
1 zero eigenvalues are removed

 The first 2 nonzero eigenvalues are given as 
  -0.0005 + 1.0000i
  -0.0005 - 1.0000i

Autonomous SSM

S = SSM(DS);
set(S.Options, 'reltol', 1,'notation','multiindex');
resonant_modes = [1 2]; % choose master spectral subspace
order = 35;                  % SSM expansion order
S.choose_E(resonant_modes)
[W0,R0] = S.compute_whisker(order);
% Reduced dynamics in symbolic form

lamdMaster = DS.spectrum.Lambda(resonant_modes);
options = struct();
options.isauto = true;
options.isdamped = true;
options.numDigits = 4;
y = reduced_dynamics_symbolic(lamdMaster,R0,options)

sigma_out = 0
sigma_in = 1
Manifold computation time at order 2 = 00:00:00
Estimated memory usage at order  2 = 2.98E-02 MB
Manifold computation time at order 3 = 00:00:00
Estimated memory usage at order  3 = 3.10E-02 MB
Manifold computation time at order 4 = 00:00:00
Estimated memory usage at order  4 = 3.35E-02 MB
Manifold computation time at order 5 = 00:00:00
Estimated memory usage at order  5 = 3.61E-02 MB
Manifold computation time at order 6 = 00:00:00
Estimated memory usage at order  6 = 4.01E-02 MB
Manifold computation time at order 7 = 00:00:00
Estimated memory usage at order  7 = 4.45E-02 MB
Manifold computation time at order 8 = 00:00:00
Estimated memory usage at order  8 = 5.07E-02 MB
Manifold computation time at order 9 = 00:00:00
Estimated memory usage at order  9 = 5.74E-02 MB
Manifold computation time at order 10 = 00:00:00
Estimated memory usage at order  10 = 6.62E-02 MB
Manifold computation time at order 11 = 00:00:00
Estimated memory usage at order  11 = 7.56E-02 MB
Manifold computation time at order 12 = 00:00:00
Estimated memory usage at order  12 = 8.76E-02 MB
Manifold computation time at order 13 = 00:00:00
Estimated memory usage at order  13 = 1.00E-01 MB
Manifold computation time at order 14 = 00:00:00
Estimated memory usage at order  14 = 1.16E-01 MB
Manifold computation time at order 15 = 00:00:00
Estimated memory usage at order  15 = 1.32E-01 MB
Manifold computation time at order 16 = 00:00:00
Estimated memory usage at order  16 = 1.52E-01 MB
Manifold computation time at order 17 = 00:00:00
Estimated memory usage at order  17 = 1.73E-01 MB
Manifold computation time at order 18 = 00:00:00
Estimated memory usage at order  18 = 1.97E-01 MB
Manifold computation time at order 19 = 00:00:00
Estimated memory usage at order  19 = 2.22E-01 MB
Manifold computation time at order 20 = 00:00:00
Estimated memory usage at order  20 = 2.52E-01 MB
Manifold computation time at order 21 = 00:00:00
Estimated memory usage at order  21 = 2.82E-01 MB
Manifold computation time at order 22 = 00:00:00
Estimated memory usage at order  22 = 3.17E-01 MB
Manifold computation time at order 23 = 00:00:00
Estimated memory usage at order  23 = 3.54E-01 MB
Manifold computation time at order 24 = 00:00:00
Estimated memory usage at order  24 = 3.95E-01 MB
Manifold computation time at order 25 = 00:00:00
Estimated memory usage at order  25 = 4.37E-01 MB
Manifold computation time at order 26 = 00:00:00
Estimated memory usage at order  26 = 4.85E-01 MB
Manifold computation time at order 27 = 00:00:00
Estimated memory usage at order  27 = 5.34E-01 MB
Manifold computation time at order 28 = 00:00:00
Estimated memory usage at order  28 = 5.88E-01 MB
Manifold computation time at order 29 = 00:00:00
Estimated memory usage at order  29 = 6.45E-01 MB
Manifold computation time at order 30 = 00:00:00
Estimated memory usage at order  30 = 7.07E-01 MB
Manifold computation time at order 31 = 00:00:00
Estimated memory usage at order  31 = 7.71E-01 MB
Manifold computation time at order 32 = 00:00:00
Estimated memory usage at order  32 = 8.41E-01 MB
Manifold computation time at order 33 = 00:00:00
Estimated memory usage at order  33 = 9.13E-01 MB
Manifold computation time at order 34 = 00:00:00
Estimated memory usage at order  34 = 9.92E-01 MB
Manifold computation time at order 35 = 00:00:00
Estimated memory usage at order  35 = 1.07E+00 MB
 
y =
 
2.045*rho_1^35 + 1.027*rho_1^33 + 0.5184*rho_1^31 + 0.2632*rho_1^29 
+ 0.1344*rho_1^27 + 0.06919*rho_1^25 + 0.03592*rho_1^23 + 0.01883*rho_1^21 
+ 0.009996*rho_1^19 + 0.005383*rho_1^17 + 0.002951*rho_1^15 
+ 0.001656*rho_1^13 + 0.0009584*rho_1^11 + 0.0005789*rho_1^9 
+ 0.0003731*rho_1^7 + 0.000269*rho_1^5 + 0.0002498*rho_1^3 
- 0.0005*rho_1 - 686.8*rho_1^34 - 347.6*rho_1^32 
- 177.0*rho_1^30 - 90.73*rho_1^28 - 46.88*rho_1^26 
- 24.44*rho_1^24 - 12.88*rho_1^22 - 6.873*rho_1^20 
- 3.725*rho_1^18 - 2.058*rho_1^16 - 1.165*rho_1^14 
- 0.681*rho_1^12 - 0.4159*rho_1^10 - 0.2709*rho_1^8 
- 0.1957*rho_1^6 - 0.1715*rho_1^4 - 0.2498*rho_1^2 + 1.0
 

Symbolic string to latex

% sympref('FloatingPointOutput',true);
% rho1dot = latex(y(1))
% theta1dot = latex(y(2))

Convergence of backbone curve

% 3.467*rho_1^6 - 8.963*rho_1^4 + 0.8168*rho_1^2 + 2.0
syms rho_1 positive
[coeffs,powers]=coeffs(y(2));
tmp = simplify(log(powers));
exp_idx = double(tmp./log(rho_1));

rhosamp = 0:0.01:0.95;
% orders = [3 5 7 9 11 13];
orders = [3 7 11 15 19 23 27 31 35];
plot_backbone_curves(double(coeffs),exp_idx,rhosamp,orders)
% Transient response validation
% We take an initial condition on SSM and perform forward simulation using both
% the reduced-order model and the full model

tf = 50;
nsteps = 1000;
q0 = 0.5*exp(1i*0);
q0 = [q0;conj(q0)];
z0 = reduced_to_full_traj(0,q0,W0);
traj = transient_traj_on_auto_ssm(DS, resonant_modes, W0, R0, tf, nsteps, 1:4, [], q0);

Reference solution from forward simulation

options = odeset('RelTol',1e-8,'AbsTol',1e-10,'Events',@zero_crossing_event);
[tInt1, zInt1,te,ye] = ode15s(@(t,x) pend_ode_auto(x,0), [0 tf], z0(1:2),options); % Transients

figure;
plot(traj.time,traj.phy(:,1),'r-'); hold on
plot(tInt1,zInt1(:,1),'b--');
legend('SSM-prediction','Full system')
xlabel('$t$','Interpreter',"latex",'FontSize',14);
ylabel('$z_1$','Interpreter',"latex",'FontSize',14);
% Backbone curve
% SSM prediction

set(S.FRCOptions,'outDOF',1);
freqRange = [0.1 1.1];
rhomax = 1.5;
BB = S.extract_backbone(resonant_modes, freqRange, [3 11 19 27 35],rhomax);
figbc = gcf;

sigma_out = 0
sigma_in = 1
Manifold computation time at order 2 = 00:00:00
Estimated memory usage at order  2 = 6.85E-03 MB
Manifold computation time at order 3 = 00:00:00
Estimated memory usage at order  3 = 8.10E-03 MB
gamma = 
   0.0002 - 0.2498i

Manifold computation time at order 2 = 00:00:00
Estimated memory usage at order  2 = 1.11E-02 MB
Manifold computation time at order 3 = 00:00:00
Estimated memory usage at order  3 = 1.24E-02 MB
Manifold computation time at order 4 = 00:00:00
Estimated memory usage at order  4 = 1.48E-02 MB
Manifold computation time at order 5 = 00:00:00
Estimated memory usage at order  5 = 1.74E-02 MB
Manifold computation time at order 6 = 00:00:00
Estimated memory usage at order  6 = 2.15E-02 MB
Manifold computation time at order 7 = 00:00:00
Estimated memory usage at order  7 = 2.58E-02 MB
Manifold computation time at order 8 = 00:00:00
Estimated memory usage at order  8 = 3.21E-02 MB
Manifold computation time at order 9 = 00:00:00
Estimated memory usage at order  9 = 3.87E-02 MB
Manifold computation time at order 10 = 00:00:00
Estimated memory usage at order  10 = 4.76E-02 MB
Manifold computation time at order 11 = 00:00:00
Estimated memory usage at order  11 = 5.70E-02 MB
gamma = 
   0.0002 - 0.2498i
   0.0003 - 0.1715i
   0.0004 - 0.1957i
   0.0006 - 0.2709i
   0.0010 - 0.4159i

Manifold computation time at order 2 = 00:00:00
Estimated memory usage at order  2 = 1.64E-02 MB
Manifold computation time at order 3 = 00:00:00
Estimated memory usage at order  3 = 1.76E-02 MB
Manifold computation time at order 4 = 00:00:00
Estimated memory usage at order  4 = 2.01E-02 MB
Manifold computation time at order 5 = 00:00:00
Estimated memory usage at order  5 = 2.26E-02 MB
Manifold computation time at order 6 = 00:00:00
Estimated memory usage at order  6 = 2.67E-02 MB
Manifold computation time at order 7 = 00:00:00
Estimated memory usage at order  7 = 3.11E-02 MB
Manifold computation time at order 8 = 00:00:00
Estimated memory usage at order  8 = 3.73E-02 MB
Manifold computation time at order 9 = 00:00:00
Estimated memory usage at order  9 = 4.39E-02 MB
Manifold computation time at order 10 = 00:00:00
Estimated memory usage at order  10 = 5.28E-02 MB
Manifold computation time at order 11 = 00:00:00
Estimated memory usage at order  11 = 6.22E-02 MB
Manifold computation time at order 12 = 00:00:00
Estimated memory usage at order  12 = 7.42E-02 MB
Manifold computation time at order 13 = 00:00:00
Estimated memory usage at order  13 = 8.69E-02 MB
Manifold computation time at order 14 = 00:00:00
Estimated memory usage at order  14 = 1.03E-01 MB
Manifold computation time at order 15 = 00:00:00
Estimated memory usage at order  15 = 1.19E-01 MB
Manifold computation time at order 16 = 00:00:00
Estimated memory usage at order  16 = 1.39E-01 MB
Manifold computation time at order 17 = 00:00:00
Estimated memory usage at order  17 = 1.59E-01 MB
Manifold computation time at order 18 = 00:00:00
Estimated memory usage at order  18 = 1.84E-01 MB
Manifold computation time at order 19 = 00:00:00
Estimated memory usage at order  19 = 2.09E-01 MB
gamma = 
   0.0002 - 0.2498i
   0.0003 - 0.1715i
   0.0004 - 0.1957i
   0.0006 - 0.2709i
   0.0010 - 0.4159i
   0.0017 - 0.6810i
   0.0030 - 1.1650i
   0.0054 - 2.0580i
   0.0100 - 3.7252i

Manifold computation time at order 2 = 00:00:00
Estimated memory usage at order  2 = 2.26E-02 MB
Manifold computation time at order 3 = 00:00:00
Estimated memory usage at order  3 = 2.38E-02 MB
Manifold computation time at order 4 = 00:00:00
Estimated memory usage at order  4 = 2.63E-02 MB
Manifold computation time at order 5 = 00:00:00
Estimated memory usage at order  5 = 2.89E-02 MB
Manifold computation time at order 6 = 00:00:00
Estimated memory usage at order  6 = 3.29E-02 MB
Manifold computation time at order 7 = 00:00:00
Estimated memory usage at order  7 = 3.73E-02 MB
Manifold computation time at order 8 = 00:00:00
Estimated memory usage at order  8 = 4.35E-02 MB
Manifold computation time at order 9 = 00:00:00
Estimated memory usage at order  9 = 5.02E-02 MB
Manifold computation time at order 10 = 00:00:00
Estimated memory usage at order  10 = 5.90E-02 MB
Manifold computation time at order 11 = 00:00:00
Estimated memory usage at order  11 = 6.84E-02 MB
Manifold computation time at order 12 = 00:00:00
Estimated memory usage at order  12 = 8.04E-02 MB
Manifold computation time at order 13 = 00:00:00
Estimated memory usage at order  13 = 9.31E-02 MB
Manifold computation time at order 14 = 00:00:00
Estimated memory usage at order  14 = 1.09E-01 MB
Manifold computation time at order 15 = 00:00:00
Estimated memory usage at order  15 = 1.25E-01 MB
Manifold computation time at order 16 = 00:00:00
Estimated memory usage at order  16 = 1.45E-01 MB
Manifold computation time at order 17 = 00:00:00
Estimated memory usage at order  17 = 1.65E-01 MB
Manifold computation time at order 18 = 00:00:00
Estimated memory usage at order  18 = 1.90E-01 MB
Manifold computation time at order 19 = 00:00:00
Estimated memory usage at order  19 = 2.15E-01 MB
Manifold computation time at order 20 = 00:00:00
Estimated memory usage at order  20 = 2.45E-01 MB
Manifold computation time at order 21 = 00:00:00
Estimated memory usage at order  21 = 2.75E-01 MB
Manifold computation time at order 22 = 00:00:00
Estimated memory usage at order  22 = 3.10E-01 MB
Manifold computation time at order 23 = 00:00:00
Estimated memory usage at order  23 = 3.46E-01 MB
Manifold computation time at order 24 = 00:00:00
Estimated memory usage at order  24 = 3.87E-01 MB
Manifold computation time at order 25 = 00:00:00
Estimated memory usage at order  25 = 4.30E-01 MB
Manifold computation time at order 26 = 00:00:00
Estimated memory usage at order  26 = 4.77E-01 MB
Manifold computation time at order 27 = 00:00:00
Estimated memory usage at order  27 = 5.27E-01 MB
gamma = 
   0.0002 - 0.2498i
   0.0003 - 0.1715i
   0.0004 - 0.1957i
   0.0006 - 0.2709i
   0.0010 - 0.4159i
   0.0017 - 0.6810i
   0.0030 - 1.1650i
   0.0054 - 2.0580i
   0.0100 - 3.7252i
   0.0188 - 6.8731i
   0.0359 -12.8785i
   0.0692 -24.4401i
   0.1344 -46.8780i

Manifold computation time at order 2 = 00:00:00
Estimated memory usage at order  2 = 2.98E-02 MB
Manifold computation time at order 3 = 00:00:00
Estimated memory usage at order  3 = 3.10E-02 MB
Manifold computation time at order 4 = 00:00:00
Estimated memory usage at order  4 = 3.35E-02 MB
Manifold computation time at order 5 = 00:00:00
Estimated memory usage at order  5 = 3.61E-02 MB
Manifold computation time at order 6 = 00:00:00
Estimated memory usage at order  6 = 4.01E-02 MB
Manifold computation time at order 7 = 00:00:00
Estimated memory usage at order  7 = 4.45E-02 MB
Manifold computation time at order 8 = 00:00:00
Estimated memory usage at order  8 = 5.07E-02 MB
Manifold computation time at order 9 = 00:00:00
Estimated memory usage at order  9 = 5.74E-02 MB
Manifold computation time at order 10 = 00:00:00
Estimated memory usage at order  10 = 6.62E-02 MB
Manifold computation time at order 11 = 00:00:00
Estimated memory usage at order  11 = 7.56E-02 MB
Manifold computation time at order 12 = 00:00:00
Estimated memory usage at order  12 = 8.76E-02 MB
Manifold computation time at order 13 = 00:00:00
Estimated memory usage at order  13 = 1.00E-01 MB
Manifold computation time at order 14 = 00:00:00
Estimated memory usage at order  14 = 1.16E-01 MB
Manifold computation time at order 15 = 00:00:00
Estimated memory usage at order  15 = 1.32E-01 MB
Manifold computation time at order 16 = 00:00:00
Estimated memory usage at order  16 = 1.52E-01 MB
Manifold computation time at order 17 = 00:00:00
Estimated memory usage at order  17 = 1.73E-01 MB
Manifold computation time at order 18 = 00:00:00
Estimated memory usage at order  18 = 1.97E-01 MB
Manifold computation time at order 19 = 00:00:00
Estimated memory usage at order  19 = 2.22E-01 MB
Manifold computation time at order 20 = 00:00:00
Estimated memory usage at order  20 = 2.52E-01 MB
Manifold computation time at order 21 = 00:00:00
Estimated memory usage at order  21 = 2.82E-01 MB
Manifold computation time at order 22 = 00:00:00
Estimated memory usage at order  22 = 3.17E-01 MB
Manifold computation time at order 23 = 00:00:00
Estimated memory usage at order  23 = 3.54E-01 MB
Manifold computation time at order 24 = 00:00:00
Estimated memory usage at order  24 = 3.95E-01 MB
Manifold computation time at order 25 = 00:00:00
Estimated memory usage at order  25 = 4.37E-01 MB
Manifold computation time at order 26 = 00:00:00
Estimated memory usage at order  26 = 4.85E-01 MB
Manifold computation time at order 27 = 00:00:00
Estimated memory usage at order  27 = 5.34E-01 MB
Manifold computation time at order 28 = 00:00:00
Estimated memory usage at order  28 = 5.88E-01 MB
Manifold computation time at order 29 = 00:00:00
Estimated memory usage at order  29 = 6.45E-01 MB
Manifold computation time at order 30 = 00:00:00
Estimated memory usage at order  30 = 7.07E-01 MB
Manifold computation time at order 31 = 00:00:00
Estimated memory usage at order  31 = 7.71E-01 MB
Manifold computation time at order 32 = 00:00:00
Estimated memory usage at order  32 = 8.41E-01 MB
Manifold computation time at order 33 = 00:00:00
Estimated memory usage at order  33 = 9.13E-01 MB
Manifold computation time at order 34 = 00:00:00
Estimated memory usage at order  34 = 9.92E-01 MB
Manifold computation time at order 35 = 00:00:00
Estimated memory usage at order  35 = 1.07E+00 MB
gamma = 
   1.0e+02 *

   0.0000 - 0.0025i
   0.0000 - 0.0017i
   0.0000 - 0.0020i
   0.0000 - 0.0027i
   0.0000 - 0.0042i
   0.0000 - 0.0068i
   0.0000 - 0.0117i
   0.0001 - 0.0206i
   0.0001 - 0.0373i
   0.0002 - 0.0687i
   0.0004 - 0.1288i
   0.0007 - 0.2444i
   0.0013 - 0.4688i
   0.0026 - 0.9073i
   0.0052 - 1.7700i
   0.0103 - 3.4762i
   0.0204 - 6.8681i

Total time spent on backbone curve computation = 00:00:06

[~,idx1] = min(abs(tInt1-te(1)));
[~,idx2] = min(abs(tInt1-te(2)));

prob = coco_prob();
prob = coco_set(prob, 'cont', 'NAdapt', 2, 'h_max', 0.5, 'PtMX', 250);
prob = coco_set(prob, 'coll', 'NTST', 20);
prob = coco_set(prob, 'po', 'bifus', 'off');
funcs = {@pend_ode_auto};
coll_args = [funcs, {tInt1(idx1:idx2)-tInt1(idx1), zInt1(idx1:idx2,:), {'c'}, 0}];
prob = ode_isol2po(prob, '', coll_args{:});
[data, uidx] = coco_get_func_data(prob, 'po.orb.coll', 'data', 'uidx');
maps = data.coll_seg.maps;
ampdata.dof  = 1;
ampdata.zdim = 2;
prob = coco_add_func(prob, 'amp1', @amplitude, ampdata, 'regular', 'x1',...
    'uidx', uidx(maps.xbp_idx), 'remesh', @amplitude_remesh);

Tmin = 2*pi/max([BB.Omega]);
Tmax = 2*pi/min([BB.Omega]);
cont_args = {1, {'po.period' 'x1' 'c'}, [0.95*Tmin,1.1*Tmax]};

fprintf('\n Run=''%s'': Continue primary family of periodic orbits.\n', ...
  'freq_resp');

bd1  = coco(prob, 'auto_freq_resp', [], cont_args{:});

 Run='freq_resp': Continue primary family of periodic orbits.

    STEP   DAMPING               NORMS              COMPUTATION TIMES
  IT SIT     GAMMA     ||d||     ||f||     ||U||   F(x)  DF(x)  SOLVE
   0                          1.01e-02  1.45e+01    0.0    0.0    0.0
   1   1  1.00e+00  3.28e-02  7.04e-06  1.45e+01    0.0    0.0    0.0
   2   1  1.00e+00  4.10e-05  3.93e-12  1.45e+01    0.0    0.0    0.0
   3   1  1.00e+00  2.05e-11  3.94e-15  1.45e+01    0.0    0.0    0.0

 STEP      TIME        ||U||  LABEL  TYPE     po.period           x1            c
    0  00:00:00   1.4493e+01      1  EP      6.8146e+00   1.1200e+00  -2.0138e-11
   10  00:00:01   1.8185e+01      2          7.4394e+00   1.5843e+00   1.6095e-12
   20  00:00:01   2.3932e+01      3          8.7858e+00   2.1442e+00  -6.9035e-12
   30  00:00:01   2.9654e+01      4          1.0785e+01   2.5759e+00  -2.7126e-13
   40  00:00:01   3.4819e+01      5          1.3260e+01   2.8461e+00   1.1024e-13
   50  00:00:02   4.0077e+01      6          1.6028e+01   2.9953e+00  -9.1383e-12
   60  00:00:02   4.3346e+01      7          1.8989e+01   3.0720e+00  -3.2853e-11
   70  00:00:02   4.8434e+01      8          2.2077e+01   3.1095e+00   1.7286e-12
   75  00:00:02   4.9876e+01      9  EP      2.3466e+01   3.1189e+00   4.2188e-12

 STEP      TIME        ||U||  LABEL  TYPE     po.period           x1            c
    0  00:00:02   1.4493e+01     10  EP      6.8146e+00   1.1200e+00  -2.0138e-11
   10  00:00:02   1.0957e+01     11          6.4607e+00   6.6379e-01   6.5051e-14
   20  00:00:02   8.8986e+00     12          6.2845e+00   5.7263e-02  -2.9720e-16
   22  00:00:02   8.8861e+00     13  BP      6.2832e+00   1.0154e-02  -4.4047e-16
   30  00:00:03   1.0340e+01     14          6.4461e+00   6.3655e-01   4.0767e-13
   40  00:00:03   1.4501e+01     15          6.9930e+00   1.2788e+00  -3.6412e-13
   50  00:00:03   2.0897e+01     16          8.0479e+00   1.8827e+00   7.0132e-12
   60  00:00:03   2.6786e+01     17          9.6723e+00   2.3721e+00   2.1728e-13
   70  00:00:03   3.2222e+01     18          1.1928e+01   2.7242e+00  -1.1087e-14
   80  00:00:03   3.7227e+01     19          1.4572e+01   2.9302e+00  -2.9327e-12
   90  00:00:03   4.2399e+01     20          1.7451e+01   3.0394e+00  -2.7846e-12
  100  00:00:04   4.5679e+01     21          2.0478e+01   3.0937e+00  -4.3443e-11
  110  00:00:04   5.0575e+01     22  EP      2.3466e+01   3.1189e+00  -3.1892e-12

bd1 = coco_bd_read('auto_freq_resp');
figure(figbc); hold on
% plot([BB.Omega],[BB.Aout],'b-'); hold on
% bd = coco_bd_read('backbone_coco');
amp_auto_coco = coco_bd_col(bd1,'x1');
om_auto_coco = 2*pi./coco_bd_col(bd1,'po.period');
plot(om_auto_coco(1:4:end),amp_auto_coco(1:4:end),'ko','LineWidth',1.5,'MarkerSize',6,'DisplayName','ODE (COCO)');
xlim([0.4 1]);
ylim([0,pi]);

Published with MATLAB® R2023a