# Can someone help me implement my Heaviside function into my Kuramoto model please?

2 views (last 30 days)
Newbie on 21 Sep 2021
Edited: Paul on 21 Sep 2021
Hi there,
I have a Heaviside function that is 1 only when 5.76<t_pert<5.96 and I want to implement that in my Kuramoto model. However, this does not work as it tell me as my Heaviside is symbolic variable but my dde23 inputs are numerical. Please can anyone help?
Here is my code so far:
t1=5.76;
t2=5.96;
t_pert=linspace(0,20*pi,1000);
% Define the Heaviside function
syms t_pert
H=1-(heaviside(t_pert-t1)*heaviside(t2-t_pert));
fplot(H,[5.6,6])
title('Heaviside function')
N=2;
omega=2*pi;
tau=11*pi/36;
h=0.1; % step size
iter=100; % number of iterations
t=1:h:h*iter; % time
K_coeff=[10 10];
K_guess=[0 K_coeff(1,1);K_coeff(1,2) 0]; % full coupling matrix
tspan=[0 20*pi]; % interval of integration
s=ones(2,1); % history of DDEs
sol=dde23(@(t,theta,Z)(my_kuramoto_pert(t,theta,Z,omega,H,K_coeff)),tau,s,tspan)
thetalag1=Z(:,1);
omega+K_coeff(1,2)*sin(thetalag1(1)-theta(2))];
end
##### 2 CommentsShowHide 1 older comment
Newbie on 21 Sep 2021
Thank you for the link, I tried using this but it does not work for me.

Paul on 21 Sep 2021
Edited: Paul on 21 Sep 2021
The function H seems to do the opposite of what you want it do, i.e., it's zero over the desired range and one outside that range:
t1=5.76;
t2=5.96;
% Define the Heaviside function
syms t_pert
H=1-(heaviside(t_pert-t1)*heaviside(t2-t_pert));
fplot(H,[5.6,6],'-o') Once you have H defined the way you want, you can do this (the code executes, don't know if the result is what you expect):
Hfunc = matlabFunction(H); % convert H to a function that can be evaluated with numeric input
You could also define Hfunc directly without going through the symbolic stuff:
Hfunc1 = @(t)(double(t<5.76 | t > 5.96));
fplot(Hfunc1,[5.6 6],'-o') N=2;
omega=2*pi;
tau=11*pi/36;
h=0.1; % step size
iter=100; % number of iterations
t=1:h:h*iter; % time
K_coeff=[10 10];
K_guess=[0 K_coeff(1,1);K_coeff(1,2) 0]; % full coupling matrix
tspan=[0 20*pi]; % interval of integration
s=ones(2,1); % history of DDEs
sol=dde23(@(t,theta,Z)(my_kuramoto_pert(t,theta,Z,omega,Hfunc,K_coeff)),tau,s,tspan); % use Hfunc instead of H
thetalag1=Z(:,1);
dthetadt_pert=[omega+H(t)+K_coeff(1,1)*sin(theta(2)-theta(1)); % t is the input argument to H, which is really Hfunc
omega+K_coeff(1,2)*sin(thetalag1(1)-theta(2))];
end