# Delay differential equations where every variable can take on multiple delays

조회 수: 7 (최근 30일)
Yasir Çatal . 2023년 10월 3일
댓글: Yasir Çatal . 2023년 10월 3일
Hello,
I'm trying to code the Kuramoto model as a delay differential equation. I have a connectivity matrix and a delay matrix. So the equation looks like this: The problem is the delay matrix D has 0s on the diagonal and when formulating the DDE all the time delays has to be positive or I get the error message "The lags must all be positive". One way to circumvent this is to do an if / else loop:
if i ~= j
theta_j = Z(:, ???)
else
theta_j = theta(j);
end
Then the problem is how to formulate lags and Z. One can do the following:
lags = D(find(triu(D, 1)));
Which would only give the non-diagonal upper triangle of D. Then one would need to find a mapping between ij and lags to put in as a Z index. I am inclined to use sub2ind but they would assume that the matrix is symmetrical and square, not only the upper triangle. Is there a way to make this work?
For context, here is the full code for the Kuramoto equation:
function dtheta = krm_dde(t, theta, Z, p)
for i = 1:length(p.omega)
sum_coupling = 0;
% Do the sigma term
for j = 1:length(p.omega)
if i ~= j
theta_j = Z(:, ???);
else
theta_j = theta(j);
end
sum_coupling = sum_coupling + p.C(i, j) * sin(theta_j - theta(i));
end
% Do the RHS
dtheta(i) = p.omega(i) + K * sum_coupling;
end
end
Tangentially related: I would also appreciate any tips to make this code vectorized. It is horribly slow in this state.
Thanks,
Yasir

댓글을 달려면 로그인하십시오.

### 채택된 답변

Torsten 2023년 10월 3일
편집: Torsten 님. 2023년 10월 3일
I don't understand your mathematical formulation. Should the summation be over j instead of i ? And shouldn't there be a bracket around the sin arguments ?
If this is the case, I guess the following should work if you define the lags as tau(1)=D12,tau(2)=D13, tau(3)= D14,...,tau(n-1)=D1n,tau(n) = D21,tau(n+1) = D23,... etc.
delay_number = 0;
for i = 1:length(p.omega)
sum_coupling = 0;
% Do the sigma term
for j = 1:length(p.omega)
if j ~= i
delay_number = delay_number + 1;
theta_j = Z(j, delay_number);
else
theta_j = theta(j);
end
sum_coupling = sum_coupling + p.C(i, j) * sin(theta_j-theta(i));
end
% Do the RHS
dtheta(i) = p.omega(i) + K * sum_coupling;
end
##### 댓글 수: 2없음 표시없음 숨기기
Yasir Çatal 2023년 10월 3일
You're right about the mathematical formulation, sorry; sleep deprivation is not good for math. Editing the formula now. I will check your recommendation today.
Yasir Çatal 2023년 10월 3일
Thanks! I implemented the following and it worked:
tau = zeros(length(D) * (length(D)-1), 1);
c = 1;
for i = 1:length(D)
for j = 1:length(D)
if i~=j
tau(c) = D(i, j);
c = c + 1;
end
end
end

댓글을 달려면 로그인하십시오.

### 카테고리

Help CenterFile Exchange에서 Correlation and Convolution에 대해 자세히 알아보기

R2023a

### Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!