News | Application Package Repository Telkom University

Category: News

cooperative logistics – profit maximization

PuTI / 2025-06-20

I need help about this code dealing with a problem of cooperative logistics. I have to maximize profits, that is, demonstrate that with the three carriers together you get a profit greater than the sum of the profits given by the three carriers taken individually. MatLab takes over 30 minutes and gives infeasible. Thank you so much!
clc
clear all

%%%%%%%%%%%%%%%% PARAMETERS %%%%%%%%%%%%%%%%%%%%%

Nr=20; % number of trips
Nc=3; % number of carriers
Vr=22; % number of vehicles

alpha=1; % extra charge on revenue;
cE= 0.2;
beta=0.8;
mu=80;

% dummy quantity
M=1000000;

% data about trips 20
d= [85, 80, 123, 111, 153, 112, 150, 107, 181, 135, 109, 106, 113, 148, 164, 127, 148, 141, 108, 118];
delta=[200, 200, 200, 200, 200, 200, 200, 200, 300, 300, 300, 300, 300, 300, 300, 400, 400, 400, 400, 400];
p=[25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25, 25];

% data about carriers, INVENTATI
I=[3, 3, 3];
cc=[0.8, 1.0, 0.9];

% data about vehicles, 22
T= [8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8];
CAP= [30, 30, 30, 30, 30, 30, 30, 30, 30, 30, 30, 30, 30, 30, 30, 30, 30, 30, 30, 30, 30, 30];
Eempty= [0.3, 0.2, 0.3, 0.3, 0.2, 0.3, 0.3, 0.2, 0.3, 0.2, 0.2, 0.3, 0.2, 0.2, 0.2, 0.2, 0.3, 0.2, 0.2, 0.3, 0.3, 0.2];
Efull= [0.3, 0.4, 0.3, 0.3, 0.3, 0.4, 0.3, 0.4, 0.3, 0.4, 0.4, 0.4, 0.4, 0.4, 0.3, 0.3, 0.3, 0.4, 0.4, 0.3, 0.4, 0.4]; % it is already Efull-Eempty
IsE=[0, 0, 0, 1, 0, 0, 1, 1, 1, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 1, 0, 0]; %1 if the vehicle is electric
Owner= [1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 2, 2, 2, 2, 2, 2, 2, 3, 3, 3, 3, 3];

% repositioing distances between trips 20×20
eps= [14 16 6 19 26 27 21 10 16 18 29 7 30 27 28 10 8 11 17 10;
27 9 20 11 9 27 22 18 27 14 29 22 5 25 9 21 17 18 11 16;
25 11 10 30 19 14 7 9 19 26 27 22 29 7 12 29 24 7 9 30;
13 16 6 19 7 7 24 5 27 28 6 18 27 9 17 30 24 14 18 10;
22 23 21 13 16 16 25 21 14 8 14 13 18 17 11 7 16 29 6 6;
5 21 9 17 26 30 27 11 24 28 16 6 25 22 27 11 18 30 13 7;
13 11 15 24 20 30 26 11 7 21 7 17 29 5 11 12 22 22 8 5;
28 7 14 20 7 25 20 12 15 5 27 10 5 27 20 26 22 25 14 28;
16 19 11 12 8 25 24 9 26 26 17 15 8 5 20 19 7 21 28 24;
29 8 13 18 17 13 15 12 21 11 10 25 17 29 22 11 28 27 14 13;
10 22 26 16 22 19 7 28 10 14 12 28 27 10 12 6 23 30 28 26;
18 11 15 25 20 7 21 17 21 9 16 5 11 5 18 16 17 26 21 19;
14 15 25 19 11 16 8 25 21 9 11 8 13 5 10 7 5 28 30 8;
26 15 11 19 21 23 9 24 9 26 16 17 9 27 25 27 20 21 15 13;
27 27 14 19 28 9 22 24 18 11 22 15 15 26 27 27 27 19 17 19;
25 11 17 24 6 30 15 23 21 30 16 5 18 25 12 20 9 11 25 27;
29 9 17 11 24 8 29 28 21 16 20 15 6 28 8 22 7 13 26 22;
23 30 28 15 11 20 13 10 25 18 29 13 12 11 16 20 24 11 6 30;
28 23 12 29 16 29 26 22 28 17 20 8 30 22 8 25 12 7 13 9;
7 7 22 11 23 9 17 25 27 26 21 29 23 14 8 6 6 26 24 30];

% 1 if the pair of trips can be combined INVENTATO 20×20
O=[0 0 0 0 0 0 0 0 0 0 1 1 1 1 1 1 1 1 1 1;
0 0 0 0 0 0 0 0 0 0 1 1 1 1 1 1 1 1 1 1;
0 0 0 0 0 0 0 0 0 0 1 1 1 1 1 1 1 1 1 1;
0 0 0 0 0 0 0 0 0 0 1 1 1 1 1 1 1 1 1 1;
0 0 0 0 0 0 0 0 0 0 1 1 1 1 1 1 1 1 1 1;
0 0 0 0 0 0 0 0 0 0 1 1 1 1 1 1 1 1 1 1;
0 0 0 0 0 0 0 0 0 0 1 1 1 1 1 1 1 1 1 1;
0 0 0 0 0 0 0 0 0 0 1 1 1 1 1 1 1 1 1 1;
0 0 0 0 0 0 0 0 0 0 1 1 1 1 1 1 1 1 1 1;
0 0 0 0 0 0 0 0 0 0 1 1 1 1 1 1 1 1 1 1;
1 1 1 1 1 1 1 1 1 1 0 0 0 0 0 0 0 0 0 0;
1 1 1 1 1 1 1 1 1 1 0 0 0 0 0 0 0 0 0 0;
1 1 1 1 1 1 1 1 1 1 0 0 0 0 0 0 0 0 0 0;
1 1 1 1 1 1 1 1 1 1 0 0 0 0 0 0 0 0 0 0;
1 1 1 1 1 1 1 1 1 1 0 0 0 0 0 0 0 0 0 0;
1 1 1 1 1 1 1 1 1 1 0 0 0 0 0 0 0 0 0 0;
1 1 1 1 1 1 1 1 1 1 0 0 0 0 0 0 0 0 0 0;
1 1 1 1 1 1 1 1 1 1 0 0 0 0 0 0 0 0 0 0;
1 1 1 1 1 1 1 1 1 1 0 0 0 0 0 0 0 0 0 0;
1 1 1 1 1 1 1 1 1 1 0 0 0 0 0 0 0 0 0 0];
%O=[0 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1;
%1 0 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1;
%1 1 0 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1;
%1 1 1 0 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1;
%1 1 1 1 0 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1;
%1 1 1 1 1 0 1 1 1 1 1 1 1 1 1 1 1 1 1 1;
%1 1 1 1 1 1 0 1 1 1 1 1 1 1 1 1 1 1 1 1;
%1 1 1 1 1 1 1 0 1 1 1 1 1 1 1 1 1 1 1 1;
%1 1 1 1 1 1 1 1 0 1 1 1 1 1 1 1 1 1 1 1;
%1 1 1 1 1 1 1 1 1 0 1 1 1 1 1 1 1 1 1 1;
%1 1 1 1 1 1 1 1 1 1 0 1 1 1 1 1 1 1 1 1;
%1 1 1 1 1 1 1 1 1 1 1 0 1 1 1 1 1 1 1 1;
%1 1 1 1 1 1 1 1 1 1 1 1 0 1 1 1 1 1 1 1;
%1 1 1 1 1 1 1 1 1 1 1 1 1 0 1 1 1 1 1 1;
%1 1 1 1 1 1 1 1 1 1 1 1 1 1 0 1 1 1 1 1;
%1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 0 1 1 1 1;
%1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 0 1 1 1;
%1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 0 1 1;
%1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 0 1;
%1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 0];

% distances between the depot of vehicle v (column) and
% the origin of trip n (row) 20×22 dati mancanti metto 1000 da 1 a 10
dO= [3 8 4 2 2 4 3 9 3 4 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000;
9 3 5 7 8 2 5 9 1 4 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000;
8 9 5 1 3 9 5 6 8 5 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000;
5 8 2 2 9 2 8 2 7 7 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000;
1 6 9 5 6 5 6 3 1 9 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000;
1 5 7 7 4 8 1 5 2 3 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000;
1 3 1 5 8 2 6 1 5 2 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000;
5 8 7 6 3 7 7 9 1 6 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000;
1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 3 8 8 2 7 7 7 1000 1000 1000 1000 1000;
1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 7 7 8 2 8 3 9 1000 1000 1000 1000 1000;
1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 5 6 7 2 3 3 6 1000 1000 1000 1000 1000;
1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 2 6 8 5 4 2 8 1000 1000 1000 1000 1000;
1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1 2 6 5 6 7 7 1000 1000 1000 1000 1000;
1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 9 1 2 7 8 8 5 1000 1000 1000 1000 1000;
1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 8 9 8 9 8 9 5 1000 1000 1000 1000 1000;
1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 8 6 7 9 8;
1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1 5 5 7 9;
1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 5 1 8 6 4;
1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 2 4 1 3 8;
1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 2 8 9 7 3];

% distances between the destination of trip n (row)
% and the depot of vehicle v (column) 20×22 dati mancanti metto 1000 da 1 a 10

dD=[7 5 1 3 9 1 5 2 1 9 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000;
7 2 7 3 3 1 5 8 1 1 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000;
1 7 2 5 3 5 8 4 7 8 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000;
7 7 4 8 3 7 3 1 7 1 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000;
5 1 1 6 7 2 4 9 9 9 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000;
6 7 1 1 5 2 4 3 9 5 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000;
8 7 2 9 5 6 6 5 8 9 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000;
6 8 9 7 6 8 2 5 8 8 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000;
1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 2 7 7 6 3 9 1 1000 1000 1000 1000 1000;
1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 2 9 7 3 6 9 6 1000 1000 1000 1000 1000;
1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 3 4 6 6 5 9 9 1000 1000 1000 1000 1000;
1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 7 7 1 7 2 1 7 1000 1000 1000 1000 1000;
1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 5 6 7 9 8 6 2 1000 1000 1000 1000 1000;
1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 4 1 5 4 1 3 3 1000 1000 1000 1000 1000;
1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 6 6 6 4 6 9 5 1000 1000 1000 1000 1000;
1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1 5 6 7 2;
1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 9 9 6 8 2;
1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 2 9 4 6 2;
1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 5 5 4 9 1;
1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 4 8 6 9 9];

% operative cost 3×22
c=[0.7, 0.8, 0.8, 0.6, 0.9, 0.8, 0.6, 0.9, 0.6, 0.8, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100;
100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 0.9, 0.8, 0.8, 1, 0.7, 0.7, 0.8, 100, 100, 100, 100, 100;
100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 0.9, 0.9, 0.6, 0.7, 0.7];

% 1 if the trip (row) is owned by the carrier (column) INVENTATO 20×3
z=[1 0 0;
1 0 0;
1 0 0;
1 0 0;
1 0 0;
1 0 0;
1 0 0;
1 0 0;
0 1 0;
0 1 0;
0 1 0;
0 1 0;
0 1 0
0 1 0;
0 1 0;
0 0 1;
0 0 1;
0 0 1;
0 0 1;
0 0 1];

% initial profit of carriers INVENTATO dati ottenuti dal prob1 per ciascun
% carrier
S0= [3380 2382 1561];
%S0 = [1000 1000 1000];

TOA=zeros(Nr,Vr);

% computation of TOA(n,v) for single trips;
for n=1:Nr
for v=1:Vr
TOA(n,v) = dO(n,v)+d(n)+dD(n,v);
end
end

E=zeros(Nr,Vr);
% computation of cost E;
for n=1:Nr
for v=1:Vr
E(n,v)= TOA(n,v)*Eempty(v) + (Efull(v)*d(n)*p(n)/CAP(v));
end
end

TOAC=zeros(Nr,Nr,Vr);
% computation of TOAC(n,k,v) for pairs of trips;
for n=1:Nr
for k=1:Nr
for v=1:Vr
TOAC(n,k,v) = dO(n,v)+d(n)+d(k)+eps(n,k)+dD(k,v);
end
end
end

EC=zeros(Nr,Nr,Vr);
% computation of cost E for combined trips;
for n=1:Nr
for k=1:Nr
for v=1:Vr
EC(n,k,v)= TOAC(n,k,v)*Eempty(v) + (Efull(v)*(d(n)*p(n)+d(k)*p(k))/CAP(v));
end
end
end

% % %%%%%%%%%%%%%%%% MAX %%%%%%%%%%%%%%%%%%%%%
%
%TO BE NOTED: I have not included the carrier in decision variables
% because here we are not includeing vehicles in the data structrue related
% to carriers but we are considering the whole list of vehicles
% than we associate vehicles with carriers by using the parameters Owner(v)

% decision variables with lower and upper bounds
x = optimvar(‘x’,Nr,Vr,’Type’,’integer’,’LowerBound’,0, ‘UpperBound’,1);
y = optimvar(‘y’,Nr,Nr,Vr,’Type’,’integer’,’LowerBound’,0, ‘UpperBound’,1);

% problema definition
prob = optimproblem(‘ObjectiveSense’,’max’);

% objective function computation

obj=0;
% %profit S1 of carriers;

for r=1:Nc
for n=1:Nr
for v=1:Vr
if (Owner(v) == r)
obj=obj+(alpha* d(n)*I(r) – TOA(n,v)*c(r,v)-(1-IsE(v))*beta*cE*E(n,v))*x(n,v);
end
end
end
end
%
%profit S2 of carriers;

for r=1:Nc
for n=1:Nr
for k=1:Nr
if (k ~= n)
for v=1:Vr
if (Owner(v) == r)
obj= obj+ (alpha*(d(n)+d(k))*I(r)-TOAC(n,k,v)*c(r,v)-(1-IsE(v))*beta*cE*EC(n,k,v))*y(n,k,v);
end
end
end
end
end
end

%profit S3 of carriers;

for r=1:Nc
%compensation profits and costs for single trips
for n=1:Nr
for rr=1:Nc
if (rr ~= r)
for v=1:Vr
if (Owner(v) == rr)
obj=obj+ delta(n)*z(n,r)*x(n,v);
else
if (Owner(v) == r)
obj=obj- delta(n)*z(n,rr)*x(n,v);
end
end
end
end
end
end
%compensation profits and costs for combined trips
for n=1:Nr
for k=1:Nr
for rr=1:Nc
if (rr ~= r)
for v=1:Vr
if (Owner(v) == rr)
obj = obj+ y(n,k,v)*(delta(n)*z(n,r)+delta(k)*z(k,r));
else
if (Owner(v) == r)
obj = obj- y(n,k,v)*(delta(n)*z(n,rr)+delta(k)*z(k,rr));
end
end
end
end
end
end
end
end

prob.Objective = obj;

%%% CONSTRAINTS

% profits higher than initial profits

constr_profit = optimineq(Nc,1);
%
for r=1:Nc
S=0;
%S1 computation
for n=1:Nr
for v=1:Vr
if (Owner(v) == r)
S=S+(alpha* d(n)*I(r) – TOA(n,v)*c(r,v)-(1-IsE(v))*beta*cE*E(n,v))*x(n,v);
end
end
end

%S2 computation
for n=1:Nr
for k=1:Nr
if (k ~= n)
for v=1:Vr
if (Owner(v) == r)
S= S + (alpha*(d(n)+d(k))*I(r)-TOAC(n,k,v)*c(r,v)-(1-IsE(v))*beta*cE*EC(n,k,v))*y(n,k,v);
end
end
end
end
end

% S3 computation;

%compensation profits and costs for single trips
for n=1:Nr
for rr=1:Nc
if (rr ~= r)
for v=1:Vr
if (Owner(v) == rr)
S=S+ delta(n)*z(n,r)*x(n,v);
else
if (Owner(v) == r)
S=S- delta(n)*z(n,rr)*x(n,v);
end
end
end
end
end
end
%compensation profits and costs for combined trips
for n=1:Nr
for k=1:Nr
for rr=1:Nc
if (rr ~= r)
for v=1:Vr
if (Owner(v) == rr)
S = S+ y(n,k,v)*(delta(n)*z(n,r)+delta(k)*z(k,r));
else
if (Owner(v) == r)
S = S- y(n,k,v)*(delta(n)*z(n,rr)+delta(k)*z(k,rr));
end
end
end
end
end
end
end
constr_profit(r)=S>=S0(r);
end

prob.Constraints.constr_profit = constr_profit;

% all demand is executed (14);

constr_demand_sat=optimeq(1,1);
%
for n=1:Nr
temp=0;
for v=1:Vr
temp=temp+x(n,v);
for k=1:Nr
if (k ~= n)
temp=temp+y(n,k,v)+y(k,n,v);
end
end
end
constr_demand_sat(n)=temp==1;
end

prob.Constraints.constr_demand_sat = constr_demand_sat;

% Every vehicle must be assigned to at most one trip;

constr_assign=optimineq(Vr,1);

for v=1:Vr
temp=0;
for n=1:Nr
temp=temp+x(n,v);
for k=1:Nr
if (k ~= n)
temp=temp+y(n,k,v);
end
end

end
constr_assign(v)=temp<=1;
end

prob.Constraints.constr_assign = constr_assign;

% duration of the batteries of electric trucks;

constr_batteries_single=optimineq(Nr,Vr);

for n=1:Nr
for v=1:Vr
if (IsE(v) == 1)
constr_batteries_single(n,v)=x(n,v)*TOA(n,v)<=T(v)*mu;
else
constr_batteries_single(n,v)=x(n,v)*TOA(n,v)<=M;
end
end
end

prob.Constraints.constr_batteries_single = constr_batteries_single;

constr_batteries_combined=optimineq(Nr,Nr,Vr);

for n=1:Nr
for k=1:Nr
if (k ~= n)
for v=1:Vr
if (IsE(v) == 1)
constr_batteries_combined(n,k,v)=y(n,k,v)*TOAC(n,k,v)<=T(v)*mu;
else
constr_batteries_combined(n,k,v)=y(n,k,v)*TOAC(n,k,v)<=M;
end
end
end
end
end

prob.Constraints.constr_batteries_combined = constr_batteries_combined;

% constraints on trips combination;

constr_combined=optimineq(Nr,Nr,Vr,1);

for n=1:Nr
for k=1:Nr
for v=1:Vr
constr_combined(n,k,v)=y(n,k,v)<=O(n,k);
end
end
end

prob.Constraints.constr_combined = constr_combined;

%
% % problem solution
show(prob)
[sol,val] = solve(prob);

% retrieve of optimal values
opt_variables_x=sol.x;
opt_variables_y=sol.y;

opt_obj=val;

% profits higher than initial profits

% %SS=zeros(Nc,1);
% %for r=1:Nc
% SS(r)=0;
% %S1 computation
% for n=1:Nr
% for v=1:Vr
% if (Owner(v) == r)
% SS(r)=SS(r)+(alpha* d(n)*I(r) – TOA(n,v)*c(r,v)-(1-IsE(v))*beta*cE*E(n,v))*sol.x(n,v);
% end
% end
% end
%
%
% %S2 computation
% for n=1:Nr
% for k=1:Nr
% if (k ~= n)
% for v=1:Vr
% if (Owner(v) == r)
% SS(r)= SS(r) + (alpha*(d(n)+d(k))*I(r)-TOAC(n,k,v)*c(r,v)-(1-IsE(v))*beta*cE*EC(n,k,v))*sol.y(n,k,v);
% end
% end
% end
% end
% end
%
%
% % S3 computation;
%
% %compensation profits and costs for single trips
% for n=1:Nr
% for rr=1:Nc
% if (rr ~= r)
% for v=1:Vr
% if (Owner(v) == rr)
% SS(r)=SS(r)+ delta(n)*z(n,r)*sol.x(n,v);
% else
% if (Owner(v) == r)
% SS(r)=SS(r)- delta(n)*z(n,rr)*sol.x(n,v);
% end
% end
% end
% end
% end
% end
% %compensation profits and costs for combined trips
% for n=1:Nr
% for k=1:Nr
% for rr=1:Nc
% if (rr ~= r)
% for v=1:Vr
% if (Owner(v) == rr)
% SS(r) = SS(r)+ sol.y(n,k,v)*(delta(n)*z(n,r)+delta(k)*z(k,r));
% else
% if (Owner(v) == r)
% SS(r) = SS(r)- sol.y(n,k,v)*(delta(n)*z(n,rr)+delta(k)*z(k,rr));
% end
% end
% end
% end
% end
% end
% end
%
% y2=0;
% for i=1:Nr
% for v=1:Vr
% y2=y2+opt_variables_y(i,2,v)+opt_variables_y(2,i,v);
% end;
% end;
%
% end
—-I need help about this code dealing with a problem of cooperative logistics. I have to maximize profits, that is, demonstrate that with the three carriers together you get a profit greater than the sum of the profits given by the three carriers taken individually. MatLab takes over 30 minutes and gives infeasible. Thank you so much!
clc
clear all

%%%%%%%%%%%%%%%% PARAMETERS %%%%%%%%%%%%%%%%%%%%%

Nr=20; % number of trips
Nc=3; % number of carriers
Vr=22; % number of vehicles

alpha=1; % extra charge on revenue;
cE= 0.2;
beta=0.8;
mu=80;

% dummy quantity
M=1000000;

% data about carriers, INVENTATI
I=[3, 3, 3];
cc=[0.8, 1.0, 0.9];

% distances between the destination of trip n (row)
% and the depot of vehicle v (column) 20×22 dati mancanti metto 1000 da 1 a 10

% initial profit of carriers INVENTATO dati ottenuti dal prob1 per ciascun
% carrier
S0= [3380 2382 1561];
%S0 = [1000 1000 1000];

TOA=zeros(Nr,Vr);

% computation of TOA(n,v) for single trips;
for n=1:Nr
for v=1:Vr
TOA(n,v) = dO(n,v)+d(n)+dD(n,v);
end
end

E=zeros(Nr,Vr);
% computation of cost E;
for n=1:Nr
for v=1:Vr
E(n,v)= TOA(n,v)*Eempty(v) + (Efull(v)*d(n)*p(n)/CAP(v));
end
end

TOAC=zeros(Nr,Nr,Vr);
% computation of TOAC(n,k,v) for pairs of trips;
for n=1:Nr
for k=1:Nr
for v=1:Vr
TOAC(n,k,v) = dO(n,v)+d(n)+d(k)+eps(n,k)+dD(k,v);
end
end
end

EC=zeros(Nr,Nr,Vr);
% computation of cost E for combined trips;
for n=1:Nr
for k=1:Nr
for v=1:Vr
EC(n,k,v)= TOAC(n,k,v)*Eempty(v) + (Efull(v)*(d(n)*p(n)+d(k)*p(k))/CAP(v));
end
end
end

% problema definition
prob = optimproblem(‘ObjectiveSense’,’max’);

% objective function computation

obj=0;
% %profit S1 of carriers;

for r=1:Nc
for n=1:Nr
for v=1:Vr
if (Owner(v) == r)
obj=obj+(alpha* d(n)*I(r) – TOA(n,v)*c(r,v)-(1-IsE(v))*beta*cE*E(n,v))*x(n,v);
end
end
end
end
%
%profit S2 of carriers;

for r=1:Nc
for n=1:Nr
for k=1:Nr
if (k ~= n)
for v=1:Vr
if (Owner(v) == r)
obj= obj+ (alpha*(d(n)+d(k))*I(r)-TOAC(n,k,v)*c(r,v)-(1-IsE(v))*beta*cE*EC(n,k,v))*y(n,k,v);
end
end
end
end
end
end

%profit S3 of carriers;

prob.Objective = obj;

%%% CONSTRAINTS

% profits higher than initial profits

constr_profit = optimineq(Nc,1);
%
for r=1:Nc
S=0;
%S1 computation
for n=1:Nr
for v=1:Vr
if (Owner(v) == r)
S=S+(alpha* d(n)*I(r) – TOA(n,v)*c(r,v)-(1-IsE(v))*beta*cE*E(n,v))*x(n,v);
end
end
end

%S2 computation
for n=1:Nr
for k=1:Nr
if (k ~= n)
for v=1:Vr
if (Owner(v) == r)
S= S + (alpha*(d(n)+d(k))*I(r)-TOAC(n,k,v)*c(r,v)-(1-IsE(v))*beta*cE*EC(n,k,v))*y(n,k,v);
end
end
end
end
end

% S3 computation;

prob.Constraints.constr_profit = constr_profit;

% all demand is executed (14);

constr_demand_sat=optimeq(1,1);
%
for n=1:Nr
temp=0;
for v=1:Vr
temp=temp+x(n,v);
for k=1:Nr
if (k ~= n)
temp=temp+y(n,k,v)+y(k,n,v);
end
end
end
constr_demand_sat(n)=temp==1;
end

prob.Constraints.constr_demand_sat = constr_demand_sat;

% Every vehicle must be assigned to at most one trip;

constr_assign=optimineq(Vr,1);

for v=1:Vr
temp=0;
for n=1:Nr
temp=temp+x(n,v);
for k=1:Nr
if (k ~= n)
temp=temp+y(n,k,v);
end
end

end
constr_assign(v)=temp<=1;
end

prob.Constraints.constr_assign = constr_assign;

% duration of the batteries of electric trucks;

constr_batteries_single=optimineq(Nr,Vr);

for n=1:Nr
for v=1:Vr
if (IsE(v) == 1)
constr_batteries_single(n,v)=x(n,v)*TOA(n,v)<=T(v)*mu;
else
constr_batteries_single(n,v)=x(n,v)*TOA(n,v)<=M;
end
end
end

prob.Constraints.constr_batteries_single = constr_batteries_single;

constr_batteries_combined=optimineq(Nr,Nr,Vr);

prob.Constraints.constr_batteries_combined = constr_batteries_combined;

% constraints on trips combination;

constr_combined=optimineq(Nr,Nr,Vr,1);

for n=1:Nr
for k=1:Nr
for v=1:Vr
constr_combined(n,k,v)=y(n,k,v)<=O(n,k);
end
end
end

prob.Constraints.constr_combined = constr_combined;

%
% % problem solution
show(prob)
[sol,val] = solve(prob);

% retrieve of optimal values
opt_variables_x=sol.x;
opt_variables_y=sol.y;

opt_obj=val;

% profits higher than initial profits

% %SS=zeros(Nc,1);
% %for r=1:Nc
% SS(r)=0;
% %S1 computation
% for n=1:Nr
% for v=1:Vr
% if (Owner(v) == r)
% SS(r)=SS(r)+(alpha* d(n)*I(r) – TOA(n,v)*c(r,v)-(1-IsE(v))*beta*cE*E(n,v))*sol.x(n,v);
% end
% end
% end
%
%
% %S2 computation
% for n=1:Nr
% for k=1:Nr
% if (k ~= n)
% for v=1:Vr
% if (Owner(v) == r)
% SS(r)= SS(r) + (alpha*(d(n)+d(k))*I(r)-TOAC(n,k,v)*c(r,v)-(1-IsE(v))*beta*cE*EC(n,k,v))*sol.y(n,k,v);
% end
% end
% end
% end
% end
%
%
% % S3 computation;
%
% %compensation profits and costs for single trips
% for n=1:Nr
% for rr=1:Nc
% if (rr ~= r)
% for v=1:Vr
% if (Owner(v) == rr)
% SS(r)=SS(r)+ delta(n)*z(n,r)*sol.x(n,v);
% else
% if (Owner(v) == r)
% SS(r)=SS(r)- delta(n)*z(n,rr)*sol.x(n,v);
% end
% end
% end
% end
% end
% end
% %compensation profits and costs for combined trips
% for n=1:Nr
% for k=1:Nr
% for rr=1:Nc
% if (rr ~= r)
% for v=1:Vr
% if (Owner(v) == rr)
% SS(r) = SS(r)+ sol.y(n,k,v)*(delta(n)*z(n,r)+delta(k)*z(k,r));
% else
% if (Owner(v) == r)
% SS(r) = SS(r)- sol.y(n,k,v)*(delta(n)*z(n,rr)+delta(k)*z(k,rr));
% end
% end
% end
% end
% end
% end
% end
%
% y2=0;
% for i=1:Nr
% for v=1:Vr
% y2=y2+opt_variables_y(i,2,v)+opt_variables_y(2,i,v);
% end;
% end;
%
% end
—- I need help about this code dealing with a problem of cooperative logistics. I have to maximize profits, that is, demonstrate that with the three carriers together you get a profit greater than the sum of the profits given by the three carriers taken individually. MatLab takes over 30 minutes and gives infeasible. Thank you so much!
clc
clear all

%%%%%%%%%%%%%%%% PARAMETERS %%%%%%%%%%%%%%%%%%%%%

Nr=20; % number of trips
Nc=3; % number of carriers
Vr=22; % number of vehicles

alpha=1; % extra charge on revenue;
cE= 0.2;
beta=0.8;
mu=80;

% dummy quantity
M=1000000;

% data about carriers, INVENTATI
I=[3, 3, 3];
cc=[0.8, 1.0, 0.9];

% distances between the destination of trip n (row)
% and the depot of vehicle v (column) 20×22 dati mancanti metto 1000 da 1 a 10

% initial profit of carriers INVENTATO dati ottenuti dal prob1 per ciascun
% carrier
S0= [3380 2382 1561];
%S0 = [1000 1000 1000];

TOA=zeros(Nr,Vr);

% computation of TOA(n,v) for single trips;
for n=1:Nr
for v=1:Vr
TOA(n,v) = dO(n,v)+d(n)+dD(n,v);
end
end

E=zeros(Nr,Vr);
% computation of cost E;
for n=1:Nr
for v=1:Vr
E(n,v)= TOA(n,v)*Eempty(v) + (Efull(v)*d(n)*p(n)/CAP(v));
end
end

TOAC=zeros(Nr,Nr,Vr);
% computation of TOAC(n,k,v) for pairs of trips;
for n=1:Nr
for k=1:Nr
for v=1:Vr
TOAC(n,k,v) = dO(n,v)+d(n)+d(k)+eps(n,k)+dD(k,v);
end
end
end

EC=zeros(Nr,Nr,Vr);
% computation of cost E for combined trips;
for n=1:Nr
for k=1:Nr
for v=1:Vr
EC(n,k,v)= TOAC(n,k,v)*Eempty(v) + (Efull(v)*(d(n)*p(n)+d(k)*p(k))/CAP(v));
end
end
end

% problema definition
prob = optimproblem(‘ObjectiveSense’,’max’);

% objective function computation

obj=0;
% %profit S1 of carriers;

for r=1:Nc
for n=1:Nr
for v=1:Vr
if (Owner(v) == r)
obj=obj+(alpha* d(n)*I(r) – TOA(n,v)*c(r,v)-(1-IsE(v))*beta*cE*E(n,v))*x(n,v);
end
end
end
end
%
%profit S2 of carriers;

for r=1:Nc
for n=1:Nr
for k=1:Nr
if (k ~= n)
for v=1:Vr
if (Owner(v) == r)
obj= obj+ (alpha*(d(n)+d(k))*I(r)-TOAC(n,k,v)*c(r,v)-(1-IsE(v))*beta*cE*EC(n,k,v))*y(n,k,v);
end
end
end
end
end
end

%profit S3 of carriers;

prob.Objective = obj;

%%% CONSTRAINTS

% profits higher than initial profits

constr_profit = optimineq(Nc,1);
%
for r=1:Nc
S=0;
%S1 computation
for n=1:Nr
for v=1:Vr
if (Owner(v) == r)
S=S+(alpha* d(n)*I(r) – TOA(n,v)*c(r,v)-(1-IsE(v))*beta*cE*E(n,v))*x(n,v);
end
end
end

%S2 computation
for n=1:Nr
for k=1:Nr
if (k ~= n)
for v=1:Vr
if (Owner(v) == r)
S= S + (alpha*(d(n)+d(k))*I(r)-TOAC(n,k,v)*c(r,v)-(1-IsE(v))*beta*cE*EC(n,k,v))*y(n,k,v);
end
end
end
end
end

% S3 computation;

prob.Constraints.constr_profit = constr_profit;

% all demand is executed (14);

constr_demand_sat=optimeq(1,1);
%
for n=1:Nr
temp=0;
for v=1:Vr
temp=temp+x(n,v);
for k=1:Nr
if (k ~= n)
temp=temp+y(n,k,v)+y(k,n,v);
end
end
end
constr_demand_sat(n)=temp==1;
end

prob.Constraints.constr_demand_sat = constr_demand_sat;

% Every vehicle must be assigned to at most one trip;

constr_assign=optimineq(Vr,1);

for v=1:Vr
temp=0;
for n=1:Nr
temp=temp+x(n,v);
for k=1:Nr
if (k ~= n)
temp=temp+y(n,k,v);
end
end

end
constr_assign(v)=temp<=1;
end

prob.Constraints.constr_assign = constr_assign;

% duration of the batteries of electric trucks;

constr_batteries_single=optimineq(Nr,Vr);

for n=1:Nr
for v=1:Vr
if (IsE(v) == 1)
constr_batteries_single(n,v)=x(n,v)*TOA(n,v)<=T(v)*mu;
else
constr_batteries_single(n,v)=x(n,v)*TOA(n,v)<=M;
end
end
end

prob.Constraints.constr_batteries_single = constr_batteries_single;

constr_batteries_combined=optimineq(Nr,Nr,Vr);

prob.Constraints.constr_batteries_combined = constr_batteries_combined;

% constraints on trips combination;

constr_combined=optimineq(Nr,Nr,Vr,1);

for n=1:Nr
for k=1:Nr
for v=1:Vr
constr_combined(n,k,v)=y(n,k,v)<=O(n,k);
end
end
end

prob.Constraints.constr_combined = constr_combined;

%
% % problem solution
show(prob)
[sol,val] = solve(prob);

% retrieve of optimal values
opt_variables_x=sol.x;
opt_variables_y=sol.y;

opt_obj=val;

% profits higher than initial profits

% %SS=zeros(Nc,1);
% %for r=1:Nc
% SS(r)=0;
% %S1 computation
% for n=1:Nr
% for v=1:Vr
% if (Owner(v) == r)
% SS(r)=SS(r)+(alpha* d(n)*I(r) – TOA(n,v)*c(r,v)-(1-IsE(v))*beta*cE*E(n,v))*sol.x(n,v);
% end
% end
% end
%
%
% %S2 computation
% for n=1:Nr
% for k=1:Nr
% if (k ~= n)
% for v=1:Vr
% if (Owner(v) == r)
% SS(r)= SS(r) + (alpha*(d(n)+d(k))*I(r)-TOAC(n,k,v)*c(r,v)-(1-IsE(v))*beta*cE*EC(n,k,v))*sol.y(n,k,v);
% end
% end
% end
% end
% end
%
%
% % S3 computation;
%
% %compensation profits and costs for single trips
% for n=1:Nr
% for rr=1:Nc
% if (rr ~= r)
% for v=1:Vr
% if (Owner(v) == rr)
% SS(r)=SS(r)+ delta(n)*z(n,r)*sol.x(n,v);
% else
% if (Owner(v) == r)
% SS(r)=SS(r)- delta(n)*z(n,rr)*sol.x(n,v);
% end
% end
% end
% end
% end
% end
% %compensation profits and costs for combined trips
% for n=1:Nr
% for k=1:Nr
% for rr=1:Nc
% if (rr ~= r)
% for v=1:Vr
% if (Owner(v) == rr)
% SS(r) = SS(r)+ sol.y(n,k,v)*(delta(n)*z(n,r)+delta(k)*z(k,r));
% else
% if (Owner(v) == r)
% SS(r) = SS(r)- sol.y(n,k,v)*(delta(n)*z(n,rr)+delta(k)*z(k,rr));
% end
% end
% end
% end
% end
% end
% end
%
% y2=0;
% for i=1:Nr
% for v=1:Vr
% y2=y2+opt_variables_y(i,2,v)+opt_variables_y(2,i,v);
% end;
% end;
%
% end
—- #logistics #cooperation MATLAB Answers — New Questions

Matlab News

I am trying to use the PCBwriter function to output a Gerber file for PCB stack object with greater than 2 layers.

PuTI / 2025-06-20

I am trying to use the PCBwriter function to output a gerber file for pcb stack object with greater than 2 metal layers and I am getting the error:
[A, g] = gerberWrite(PcbWriteObj); Error using PCBWriter More than 2 layers of metal is not supported for generating gerber files.
Is there a way around this for PCBs with inner copper layers?
Does MATLAB support export of mulitlayer PCB designs (>2 Layers) as Gerber files?I am trying to use the PCBwriter function to output a gerber file for pcb stack object with greater than 2 metal layers and I am getting the error:
[A, g] = gerberWrite(PcbWriteObj); Error using PCBWriter More than 2 layers of metal is not supported for generating gerber files.
Is there a way around this for PCBs with inner copper layers?
Does MATLAB support export of mulitlayer PCB designs (>2 Layers) as Gerber files? I am trying to use the PCBwriter function to output a gerber file for pcb stack object with greater than 2 metal layers and I am getting the error:
[A, g] = gerberWrite(PcbWriteObj); Error using PCBWriter More than 2 layers of metal is not supported for generating gerber files.
Is there a way around this for PCBs with inner copper layers?
Does MATLAB support export of mulitlayer PCB designs (>2 Layers) as Gerber files? pcbwriter, gerber file MATLAB Answers — New Questions

Matlab News

Compare same frequency components of two signals

PuTI / 2025-06-20

Hello,
I have two signals and I want to split each one of them in same number of components. After that, I want to find the components with the (almost) same frequency and compare their plots.
Any idea how could I find same frequency components of the signals?

Thanks.Hello,
I have two signals and I want to split each one of them in same number of components. After that, I want to find the components with the (almost) same frequency and compare their plots.
Any idea how could I find same frequency components of the signals?

Thanks. Hello,
I have two signals and I want to split each one of them in same number of components. After that, I want to find the components with the (almost) same frequency and compare their plots.
Any idea how could I find same frequency components of the signals?

Thanks. signal processing, frequency MATLAB Answers — New Questions

Matlab News

Interp1 is not working after applying unique because of rounding off

PuTI / 2025-06-20

I have a double array with many values. The first value is 3.569990000000001, the next multiple values is 3.570000000000001, then comes 3.570010000000001 and it works like this. I can see only these long values which i click on that cell, otherwise every value is being displayed as 3.57. Now I applied unique, so the multiple 3.570000000000001 is removed and there is now only one 3.570000000000001 and again and so on, but still all values are seen as 3.57. So if compared the first two values to see if it equal it is showing that they are not equal even though me as the user sees both values as 3.57. Now here comes the problem. I need to use interp1 on these values, but here interp1 is not reading them as 3.569990000000001 and 3.570000000000001. Interp1 is reading both of them as 3.57 and is giving me an error saying values need to be unique. I dont know what to do. I tried googling like converting to format long or something like that, i couldnt find anything. What can I do. Thank you.I have a double array with many values. The first value is 3.569990000000001, the next multiple values is 3.570000000000001, then comes 3.570010000000001 and it works like this. I can see only these long values which i click on that cell, otherwise every value is being displayed as 3.57. Now I applied unique, so the multiple 3.570000000000001 is removed and there is now only one 3.570000000000001 and again and so on, but still all values are seen as 3.57. So if compared the first two values to see if it equal it is showing that they are not equal even though me as the user sees both values as 3.57. Now here comes the problem. I need to use interp1 on these values, but here interp1 is not reading them as 3.569990000000001 and 3.570000000000001. Interp1 is reading both of them as 3.57 and is giving me an error saying values need to be unique. I dont know what to do. I tried googling like converting to format long or something like that, i couldnt find anything. What can I do. Thank you. I have a double array with many values. The first value is 3.569990000000001, the next multiple values is 3.570000000000001, then comes 3.570010000000001 and it works like this. I can see only these long values which i click on that cell, otherwise every value is being displayed as 3.57. Now I applied unique, so the multiple 3.570000000000001 is removed and there is now only one 3.570000000000001 and again and so on, but still all values are seen as 3.57. So if compared the first two values to see if it equal it is showing that they are not equal even though me as the user sees both values as 3.57. Now here comes the problem. I need to use interp1 on these values, but here interp1 is not reading them as 3.569990000000001 and 3.570000000000001. Interp1 is reading both of them as 3.57 and is giving me an error saying values need to be unique. I dont know what to do. I tried googling like converting to format long or something like that, i couldnt find anything. What can I do. Thank you. matlab MATLAB Answers — New Questions

Matlab News

Readtable and Readmatrix Ignore Specified Range and Produce Extra Variables

PuTI / 2025-06-20

I have this very simple code:
stupid_matlab = readtable(‘PATH’, ‘ReadVariableNames’, false, ‘Range’, strcat(‘C2:D’, string(ri_length)));
However, it completely ignores the column portion of the range I give it. Changing the number part of the specified range changes the number of rows, but it always produces four extra columns. If I specify C2:C3, I get five columns. If I specify C2:D3, I get six columns. I have not idea where it is getting the extra columns from. I have also tried "readmatrix" and it does the exact same thing. I’ve attatched the .csv, but I’ve opened it directly and in excel and don’t see anything wrong.I have this very simple code:
stupid_matlab = readtable(‘PATH’, ‘ReadVariableNames’, false, ‘Range’, strcat(‘C2:D’, string(ri_length)));
However, it completely ignores the column portion of the range I give it. Changing the number part of the specified range changes the number of rows, but it always produces four extra columns. If I specify C2:C3, I get five columns. If I specify C2:D3, I get six columns. I have not idea where it is getting the extra columns from. I have also tried "readmatrix" and it does the exact same thing. I’ve attatched the .csv, but I’ve opened it directly and in excel and don’t see anything wrong. I have this very simple code:
stupid_matlab = readtable(‘PATH’, ‘ReadVariableNames’, false, ‘Range’, strcat(‘C2:D’, string(ri_length)));
However, it completely ignores the column portion of the range I give it. Changing the number part of the specified range changes the number of rows, but it always produces four extra columns. If I specify C2:C3, I get five columns. If I specify C2:D3, I get six columns. I have not idea where it is getting the extra columns from. I have also tried "readmatrix" and it does the exact same thing. I’ve attatched the .csv, but I’ve opened it directly and in excel and don’t see anything wrong. readtable, readmatrix, range MATLAB Answers — New Questions

Matlab News

Clean and automate indents in scripts an live scripts

PuTI / 2025-06-19

Hi everyone,
Please, do you have any tips to automate indentation in scripts?
Example by the picture below:Hi everyone,
Please, do you have any tips to automate indentation in scripts?
Example by the picture below: Hi everyone,
Please, do you have any tips to automate indentation in scripts?
Example by the picture below: automatic, indentation, alignment, text MATLAB Answers — New Questions

Matlab News

Varying device latency for ASIO Device

PuTI / 2025-06-19

Hi,
my goal is to measure impulse responses using Matlab and an ASIO device. I encountered the problem, that for different runs of my code, the peak of my impulse response is at a different location on the time axis. I was able to reproduce this behavior outside of my code in the "Impulse Respones Measurer". Looping back the output channel to an input channel and capturing a fiew imulse reponses while leaving everything else untouched gives me shifted impulse responses.

I also followed this example: Measure Audio Latency, which gave me the same results.
When measuring the same impulse response (i.e. same device, same configuration of channels etc.) with another software (REW), the delay remained constant between different runs.
I was also able to reproduce this behavior vor R2020a and R2024b
Any suggestions how to keep a constant latency between several runs of audioPlayerRecorder are appreciated.
Kind Regards,
SebastianHi,
my goal is to measure impulse responses using Matlab and an ASIO device. I encountered the problem, that for different runs of my code, the peak of my impulse response is at a different location on the time axis. I was able to reproduce this behavior outside of my code in the "Impulse Respones Measurer". Looping back the output channel to an input channel and capturing a fiew imulse reponses while leaving everything else untouched gives me shifted impulse responses.

I also followed this example: Measure Audio Latency, which gave me the same results.
When measuring the same impulse response (i.e. same device, same configuration of channels etc.) with another software (REW), the delay remained constant between different runs.
I was also able to reproduce this behavior vor R2020a and R2024b
Any suggestions how to keep a constant latency between several runs of audioPlayerRecorder are appreciated.
Kind Regards,
Sebastian Hi,
my goal is to measure impulse responses using Matlab and an ASIO device. I encountered the problem, that for different runs of my code, the peak of my impulse response is at a different location on the time axis. I was able to reproduce this behavior outside of my code in the "Impulse Respones Measurer". Looping back the output channel to an input channel and capturing a fiew imulse reponses while leaving everything else untouched gives me shifted impulse responses.

Matlab News

COnnecting points on opposite end of mask

PuTI / 2025-06-19

Given an mask image i have certain points as seen in the image, I want to connect the point to the point at the opposite side of it given the list of points/array of point.Well the output would look something like.Given an mask image i have certain points as seen in the image, I want to connect the point to the point at the opposite side of it given the list of points/array of point.Well the output would look something like. Given an mask image i have certain points as seen in the image, I want to connect the point to the point at the opposite side of it given the list of points/array of point.Well the output would look something like. image analysis, image processing, digital image processing, computer vision, contour MATLAB Answers — New Questions

Matlab News

Mass importing .txt files with sequential names to plot multiple graphs on same figure for easy comparison

PuTI / 2025-06-19

I have 30 .txt files numbered like so:
FreeField1.txt
FreeField2.txt
…
FreeField30.txt

I need to find a way to import the data from all the files at once so I can then plot a graph with all of them overlaid to allow for quick comparison.

The current method I have found which inputs all the data is like so
numfiles = 30;
S = dir(‘*.txt’)
N = length(S)
for i = 1:numel(S)
S(i).data = readtable(S(i).name)
end

But this gives a structure which I do not how to progress from. I am looking for any possible methods which could be used for this mass importation to then build off for the plotting of graphs.

Any help would be appreciated as I am very new to MATLab and am wishing to understand it more. I have attached a sample text file for convenience.I have 30 .txt files numbered like so:
FreeField1.txt
FreeField2.txt
…
FreeField30.txt

I need to find a way to import the data from all the files at once so I can then plot a graph with all of them overlaid to allow for quick comparison.

The current method I have found which inputs all the data is like so
numfiles = 30;
S = dir(‘*.txt’)
N = length(S)
for i = 1:numel(S)
S(i).data = readtable(S(i).name)
end

But this gives a structure which I do not how to progress from. I am looking for any possible methods which could be used for this mass importation to then build off for the plotting of graphs.

Any help would be appreciated as I am very new to MATLab and am wishing to understand it more. I have attached a sample text file for convenience. I have 30 .txt files numbered like so:
FreeField1.txt
FreeField2.txt
…
FreeField30.txt

I need to find a way to import the data from all the files at once so I can then plot a graph with all of them overlaid to allow for quick comparison.

The current method I have found which inputs all the data is like so
numfiles = 30;
S = dir(‘*.txt’)
N = length(S)
for i = 1:numel(S)
S(i).data = readtable(S(i).name)
end

But this gives a structure which I do not how to progress from. I am looking for any possible methods which could be used for this mass importation to then build off for the plotting of graphs.

Any help would be appreciated as I am very new to MATLab and am wishing to understand it more. I have attached a sample text file for convenience. data import, graph, text file MATLAB Answers — New Questions

Matlab News

3D surf plot for more than two quantities

PuTI / 2025-06-19

I want to plot the value in the same 3D graph .

syms x t r b %alpha
% Parameter values
a=(pi)/3;
g=9.8;
U=2.5;
O=7.29*10^(-5);
f=2*O*sin(a);
H=-(f/g)*U;
alpha=0.75; % fractional order
%%%%%%%%%initalization of variable
u_l=sym(zeros(1));
v_l=zeros(1,’sym’);
h_l=zeros(1,’sym’);
A_l=zeros(1,2,’sym’);
B_l=zeros(1,2,’sym’);
C_l=zeros(1,2,’sym’);
D_l=zeros(1,2,’sym’);
series1_l(x,t)=sym(zeros(1,1));
series2_l(x,t)=sym(zeros(1,1));
series3_l(x,t)=sym(zeros(1,1));
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
u_up=sym(zeros(1));
v_up=zeros(1,’sym’);
h_up=zeros(1,’sym’);
A_up=zeros(1,2,’sym’);
B_up=zeros(1,2,’sym’);
C_up=zeros(1,2,’sym’);
D_up=zeros(1,2,’sym’);
series1_up(x,t)=sym(zeros(1,1));
series2_up(x,t)=sym(zeros(1,1));
series3_up(x,t)=sym(zeros(1,1));
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
u_cr=sym(zeros(1));
v_cr=zeros(1,’sym’);
h_cr=zeros(1,’sym’);
A_cr=zeros(1,2,’sym’);
B_cr=zeros(1,2,’sym’);
C_cr=zeros(1,2,’sym’);
D_cr=zeros(1,2,’sym’);
series1_cr(x,t)=sym(zeros(1,1));
series2_cr(x,t)=sym(zeros(1,1));
series3_cr(x,t)=sym(zeros(1,1));
%%%%%%%% Initial condition fuzzy condition
R=0.5;
b_l=0 % lower bound
b_cr=0.5 % middle value
b_u=1 % upper value
%%%%%%%%%%%%%%%%%% lOWER VALUE%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
u_l(1)=(2*b_l*(1-R)+R)*exp(x)*(sech(x))^2;
v_l(1)=(2*b_l*(1-R)+R+1)*2*x*(sech(2*x))^2;
h_l(1)=(2*b_l*(1-R)+R)*x^2*(sech(2*x))^2;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%UPPER VALUE
u_up(1)=(2*b_u*(1-R)+R)*exp(x)*(sech(x))^2;
v_up(1)=(2*b_u*(1-R)+R+1)*2*x*(sech(2*x))^2;
h_up(1)=(2*b_u*(1-R)+R)*x^2*(sech(2*x))^2;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%Middle valur %%%%%%%%%%%%%%%%%%%%
u_cr(1)=(2*b_cr*(1-R)+R)*exp(x)*(sech(x))^2;
v_cr(1)=(2*b_cr*(1-R)+R+1)*2*x*(sech(2*x))^2;
h_cr(1)=(2*b_cr*(1-R)+R)*x^2*(sech(2*x))^2;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%lOWER VALUR %%%%%%%%%%%%%%%%%%%%%%%%%
for k=1:2
A_l=0;
B_l=0;
C_l=0;
D_l=0;
if k==1
E=1;
else
E=0;
end
for r=1:k
A_l=A_l+u_l(r)*diff(u_l(k-r+1),x,1);
B_l=B_l+u_l(r)*diff(v_l(k-r+1),x,1);
C_l=C_l+u_l(r)*diff(h_l(k-r+1),x,1);
D_l=D_l+h_l(r)*diff(u_l(k-r+1),x,1);
end
u_l(k+1)=simplify(gamma(((k-1)*alpha)+1)/gamma((alpha*(k+1-1))+1)*(-A_l+f*v_l(k)-g*diff(h_l(k),x,1)));
v_l(k+1)=gamma(((k-1)*alpha)+1)/gamma((alpha*(k+1-1))+1)*(-B_l-f*u_l(k)-g*H*E);
h_l(k+1)=gamma(((k-1)*alpha)+1)/gamma((alpha*(k+1-1))+1)*(-C_l-H*v_l(k)-D_l);
end
% t1=simplify(u(2))
% var2 = vpa(t1)
% t2=simplify(v(2))
% var = vpa(t2)
% t3=simplify(h(2))
for k=1:3
series1_l(x,t)=simplify(series1_l(x,t)+u_l(k)*(power(t,k-1)));
series2_l(x,t)=simplify(series2_l(x,t)+v_l(k)*(power(t,k-1)));
series3_l(x,t)=simplify(series3_l(x,t)+h_l(k)*(power(t,k-1)));
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% UPPER VALUE %%%%%%%%%%%%%%%%
for k=1:2
A_up=0;
B_up=0;
C_up=0;
D_up=0;
if k==1
E=1;
else
E=0;
end
for r=1:k
A_up=A_up+u_up(r)*diff(u_up(k-r+1),x,1);
B_up=B_up+u_up(r)*diff(v_up(k-r+1),x,1);
C_up=C_up+u_up(r)*diff(h_up(k-r+1),x,1);
D_up=D_up+h_up(r)*diff(u_up(k-r+1),x,1);
end
u_up(k+1)=simplify(gamma(((k-1)*alpha)+1)/gamma((alpha*(k+1-1))+1)*(-A_up+f*v_up(k)-g*diff(h_up(k),x,1)));
v_up(k+1)=gamma(((k-1)*alpha)+1)/gamma((alpha*(k+1-1))+1)*(-B_up-f*u_up(k)-g*H*E);
h_up(k+1)=gamma(((k-1)*alpha)+1)/gamma((alpha*(k+1-1))+1)*(-C_up-H*v_up(k)-D_up);
end
% t1=simplify(u(2))
% var2 = vpa(t1)_u
% t2=simplify(v(2))
% var = vpa(t2)
% t3=simplify(h(2))
for k=1:3
series1_up(x,t)=simplify(series1_up(x,t)+u_up(k)*(power(t,k-1)));
series2_up(x,t)=simplify(series2_up(x,t)+v_up(k)*(power(t,k-1)));
series3_up(x,t)=simplify(series3_up(x,t)+h_up(k)*(power(t,k-1)));
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%middle value %%%%%%%%%%%%%%
for k=1:2
A_cr=0;
B_cr=0;
C_cr=0;
D_cr=0;
if k==1
E=1;
else
E=0;
end
for r=1:k
A_cr=A_cr+u_cr(r)*diff(u_up(k-r+1),x,1);
B_cr=B_cr+u_cr(r)*diff(v_up(k-r+1),x,1);
C_cr=C_cr+u_cr(r)*diff(h_up(k-r+1),x,1);
D_cr=D_up+h_cr(r)*diff(u_up(k-r+1),x,1);
end
u_cr(k+1)=simplify(gamma(((k-1)*alpha)+1)/gamma((alpha*(k+1-1))+1)*(-A_cr+f*v_up(k)-g*diff(h_cr(k),x,1)));
v_cr(k+1)=gamma(((k-1)*alpha)+1)/gamma((alpha*(k+1-1))+1)*(-B_cr-f*u_cr(k)-g*H*E);
h_cr(k+1)=gamma(((k-1)*alpha)+1)/gamma((alpha*(k+1-1))+1)*(-C_cr-H*v_cr(k)-D_cr);
end
% t1=simplify(u(2))
% var2 = vpa(t1)_u
% t2=simplify(v(2))
% var = vpa(t2)
% t3=simplify(h(2))
for k=1:3
series1_cr(x,t)=simplify(series1_cr(x,t)+u_cr(k)*(power(t,k-1)));
series2_cr(x,t)=simplify(series2_cr(x,t)+v_cr(k)*(power(t,k-1)));
series3_cr(x,t)=simplify(series3_cr(x,t)+h_cr(k)*(power(t,k-1)));
end

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
C1_l=zeros(1)
C2_l=zeros(1)
C3_l=zeros(1)
C1_up=zeros(1)
C2_up=zeros(1)
C3_up=zeros(1)
C1_cr=zeros(1)
C2_cr=zeros(1)
C3_cr=zeros(1)
row=0;
x=0:0.04:2
t=0:0.002:0.1
for i=1:length(x)
row=row+1;
col=0;
for j=1:length(t)
col=col+1;
C1_l(row,col)=series1_l(x(i),t(j));
C2_l(row,col)=series2_l(x(i),t(j));
C3_l(row,col)=series3_l(x(i),t(j));
%———————————-
C1_up(row,col)=series1_up(x(i),t(j));
C2_up(row,col)=series2_up(x(i),t(j));
C3_up(row,col)=series3_up(x(i),t(j));
%————————————
C1_cr(row,col)=series1_cr(x(i),t(j));
C2_cr(row,col)=series2_cr(x(i),t(j));
C3_cr(row,col)=series3_cr(x(i),t(j));

end
end
%——————————————————–
surf(x,t,C1_l,C1_up,C1_cr)
surf(x,t,C2_l,C2_up,C2_cr)
surf(x,t,C3_l,C3_up,C3_cr)
I want to plot all the graph in the same 3D plot. I dont want different graph using subplot. and I am using MATLAB version 2023aI want to plot the value in the same 3D graph .

Matlab News

using spice models in simlink

PuTI / 2025-06-19

Hello,

I have simple spice model defined as described in the attached image.
I want to use this script in the simulink modelling environment.
Is this possible?
Can anyone help with this?

ThanksHello,

I have simple spice model defined as described in the attached image.
I want to use this script in the simulink modelling environment.
Is this possible?
Can anyone help with this?

Thanks Hello,

I have simple spice model defined as described in the attached image.
I want to use this script in the simulink modelling environment.
Is this possible?
Can anyone help with this?

Thanks simulink MATLAB Answers — New Questions

News

Microsoft to Block Third-Party App Access to User Sites and Files

Tony Redmond / 2025-06-19

New Microsoft-Managed App Consent Policy to Control User Consent for Apps

Message center notification MC1097272 (17 June 2025) announces Microsoft’s intention to restrict access to some legacy protocols and introduce a new managed app consent policy to the ability of users to grant consent to third-party apps that want access to files and sites.

Microsoft says that they are updating default settings to help Microsoft 365 tenants “meet the minimum security benchmark and harden your tenant’s security posture.” As far as I can tell, this appears to be a reference to section IM-2 of the Microsoft cloud security benchmark. For good measure, Microsoft throws in the Secure Future Initiative and Secure by Default principle to provide further justification for the change.

No Problem with Blocking Obsolete and Insecure Protocols

I don’t think anyone will complain about blocking browser access to SharePoint and OneDrive via the Relying Party Suite (RPS – another relatively unknown component for most Microsoft 365 tenants). Legacy protocols are blocked in the SharePoint tenant configuration, and this change reinforces the block.

Get-SPOTenant | Select-Object LegacyBrowserAuthProtocolsEnabled

LegacyBrowserAuthProtocolsEnabled
---------------------------------
                             True

Likewise, I don’t think anyone will complain about blocking the FrontPage Remote Procedure Call (FPRPC) protocol for Office file opens. It’s an outdated protocol that attackers have leveraged (here’s an example).

App Consent Policy to Prevent Third-Party Access to Files and Sites

My interest was drawn to the third block, which will introduce a Microsoft-managed app consent policy to require administrator consent for third-party apps that access files and sites. There are a bunch of app consent policies already present in tenants that you can see by running the Get-MgPolicyPermissionGrantPolicy cmdlet from the Microsoft Graph PowerShell SDK (any policy prefixed by “microsoft” is a Microsoft-managed app consent policy):

Get-MgPolicyPermissionGrantPolicy | Format-Table Id, DisplayName, Description -AutoSize

Like many other Microsoft 365 policies, the policy is a container, and the real settings (“condition sets”) are found by running the Get-MgPolicyPermissionGrantPolicyInclude cmdlet. For example, this app consent policy allows administrators to manage all aspects of all apps in a tenant:

Get-MgPolicyPermissionGrantPolicyInclude -PermissionGrantPolicyId "microsoft-application-admin" | Format-List

ClientApplicationIds                        : {all}
ClientApplicationPublisherIds               : {all}
ClientApplicationTenantIds                  : {all}
ClientApplicationsFromVerifiedPublisherOnly : False
Id                                          : 811d2da7-443c-43da-96e7-28d285b234e9
PermissionClassification                    : all
PermissionType                              : application
Permissions                                 : {all}
ResourceApplication                         : any
AdditionalProperties                        : {}

ClientApplicationIds                        : {all}
ClientApplicationPublisherIds               : {all}
ClientApplicationTenantIds                  : {all}
ClientApplicationsFromVerifiedPublisherOnly : False
Id                                          : 60461179-740e-4d8b-9e00-1456a338c44b
PermissionClassification                    : all
PermissionType                              : delegated
Permissions                                 : {all}
ResourceApplication                         : any
AdditionalProperties                        : {}

For more details, see the Graph documentation for permission grant policies. There’s no UX in the Entra admin center to manage app consent policies. This article throws more light onto how to build your own app consent policies.

I don’t believe that users should be able to grant consent to use any app within a tenant. Disabling the ability for users to register apps in Entra user settings is recommended (Figure 1). If some users need to register apps, do what the documentation says and assign the Application Developer role to their accounts.

Blocking the ability for users to register applications.App consent policy. — Figure 1: Blocking the ability for users to register applications

If you do allow users to register apps, it’s likely that they will need to grant consent for delegated permissions to allow apps to access the user’s data. This is an area to monitor to make sure that apps are not asking for unexplainable permissions. It’s easy to check permission consents through audit records.

If you follow my advice and don’t allow users to register applications, you’ll need to make sure that the admin consent request workflow is operational. Obviously, you should monitor and report consent approvals (the article links to a PowerShell script to do the job).

The ChatGPT Conundrum

What’s interesting about Microsoft’s move is that it neatly blocks the ability of users to grant consent for the ChatGPT app that allows them to upload files from SharePoint Online and OneDrive for Business for processing by one of the ChatGPT models. A cynic might say that Microsoft is taking this step to make sure that Microsoft 365 Copilot has sole access to files stored in SharePoint Online and OneDrive for Business. A more benign reading is that Microsoft is simply making sure that users can’t inadvertently grant access to third-party apps to access and read their Microsoft 365 files.

In any case, I don’t think people should upload files to ChatGPT because this activity creates all sorts of security concerns. Fortunately, it’s easy to find and block the ChatGPT app if it’s already in a tenant. In addition, ChatGPT cannot process encrypted files protected by sensitivity labels because it doesn’t have the access right needed to open protected files.

Don’t Drop Your Guard

No can argue that we do need to do better to secure tenants, so the changes proposed by Microsoft are welcome. The changes will begin rolling out in mid-July and are due to be in all tenants sometime in August 2025.

There are still too many tenants that don’t protect user accounts with multifactor authentication, which is why bad actors keep using password spray attackers in an attempt to compromise accounts. A recent report describes a password spray attack by a group called SneakyStrike against Entra ID accounts. The report is a little overhyped, but it’s a good reminder that attackers still patiently look for weak tenants to penetrate.

Insight like this doesn’t come easily. You’ve got to know the technology and understand how to look behind the scenes. Benefit from the knowledge and experience of the Office 365 for IT Pros team by subscribing to the best eBook covering Office 365 and the wider Microsoft 365 ecosystem.

News

Updating the Entra ID Custom Banned Password List with PowerShell

Tony Redmond / 2025-06-19

Use Microsoft Graph PowerShell SDK Cmdlets to Maintain the Entra ID Custom Banned Password List

Vasil Michev is busy these days. Apart from his day job, he’s doing the technical reviews for the Office 365 for IT Pros (2026 edition) and Automating Microsoft 365 with PowerShell (2^nd edition) eBooks, both due for release on July 1, 2025. Technical editing is an important part of our publication process because it’s an annual end-to-end review of all content to help authors refine their chapters, eliminate old and unnecessary text, and consider what they should be covering.

And still Vasil finds time for his own writing, such as a recent article about using the Microsoft Graph PowerShell SDK to update the banned password list for Entra ID accounts. Given that the Graph PowerShell SDK is a major topic for Automating Microsoft with PowerShell, my attention was immediately drawn to the article to understand what it described and consider it for inclusion in the book. It is now, along with 350+ pages of other PowerShell content about automating different aspects of Microsoft 365 activities.

Global Banned Password List

The Entra ID password protection feature maintains a global list of banned passwords. Microsoft maintains the list and updates it on an ongoing basis from telemetry for Entra ID authentication. All attempts to change account passwords are checked against the global banned list to make sure that the new password is reasonably strong. In other words, it’s not something like “Mypassword” or “Cats.” Tenant administrators cannot affect how Entra ID uses the global list of banned passwords, nor can they add or remove values from the list. It’s just part of how Entra ID works, and this part of password protection is included in the version of Entra ID included with all Microsoft 365 tenants.

Custom Banned Password List

If a tenant has Entra P1 or P2 licenses, they can implement a custom banned password list. The custom list supplements the global banned password list. The custom list is limited to 1,000 entries, but those entries are “key base terms” of between 4 and 16 characters. In other words, Entra ID blocks variations and combinations of the terms in the custom banned password list.

When a custom banned password list is available, Entra ID combines its entries with the global banned password list. The idea is that tenants might want to stop people using organization-specific terms like the names of locations or buildings in passwords because these terms might be easy for attackers to guess in a spray attack. Of course, you shouldn’t be depending on passwords and should deploy multifactor authentication to protect accounts, but it’s worthwhile protecting passwords as much as possible.

Blocking Passwords

Figure 1 shows some of the entries in the custom banned password list as viewed through the Entra admin center. You can see that the last entry is for “VictorMeldrew.” This is a key base term for password checking.

The custom banned password list in the Entra admin center. — Figure 1: The custom banned password list in the Entra admin center

In Figure 2, an administrator has attempted to change an account password through the Microsoft 365 admin center. The password looks strong, but Entra ID rejects it because it includes a key base term. Telling the administrator that the password is easily guessable is just the way Microsoft chose to say: “can’t use that password.”

The custom banned password list stops a password based on a key base term. — Figure 2: : The custom banned password list stops a password based on a key base term

Updating the Custom Banned Password List with a Script

Vasil’s article covers the basics of creating a directory settings object to hold password protection settings, including the custom banned password list. I used that information to create a script that’s more like something you might use as production code, which you can download from GitHub.

The code:

Checks if the correct permission (Directory.ReadWrite.All) is available to read, create, and update directory settings. This is a very high-level permission that should be restricted as tightly as possible. You should also monitor the apps that hold this permission to make sure that they are used correctly.
Import a list of key base terms from a CSV file and checks that each term is at least 4 and no more than 16 characters long.
Uses the Get-MgBetaDirectorySetting cmdlet to check if a directory setting object for password protection is defined in the tenant. If not, the script runs the New-MgBetaDirectorySetting cmdlet to create and populate a new directory setting object with the list of key base terms (and other default values). The directory setting object is derived from the directory settings template for password rules. The template always has an identifier of 5cf42378-d67d-4f36-ba46-e8b86229381d.
If a directory setting object for password protection is available, fetch the list of current key base terms and combine it with the new list to generate a combined list. The Update-MgBetaDirectorySetting cmdlet then updates the directory setting object with the combined list.
Export the newly-updated list to a CSV file.

If you prefer to use the input CSV file as the definitive set of key base terms and not combine the input set with the current set, it’s easy to comment out the two lines that create a combined list.

The only semi-weird thing about the list of key base terms is that it uses tabs for delimitation (which is why the code splits the list using [char]9).

Hopefully the script is of some use. If not, I won’t be offended. Check out the 320-plus scripts in the Office365Itpros GitHub repository. You might find something more useful there!

Need some assistance to write and manage PowerShell scripts for Microsoft 365? Get a copy of the Automating Microsoft 365 with PowerShell eBook, available standalone or as part of the Office 365 for IT Pros eBook bundle.

Matlab News

Heatmap inside custom 2d geometry shape

PuTI / 2025-06-18

Hi,
I have a specific 2d geometry of an object with multiple fields. I also know the data values for each field. How can I plot a heatmap with my values inside the fileds using the shape of my 2d object. I could get the CAD of the object for example, create a smiple drawing myself, or imitate my 2d shape with blocks or somehow else ….
Thanks.Hi,
I have a specific 2d geometry of an object with multiple fields. I also know the data values for each field. How can I plot a heatmap with my values inside the fileds using the shape of my 2d object. I could get the CAD of the object for example, create a smiple drawing myself, or imitate my 2d shape with blocks or somehow else ….
Thanks. Hi,
I have a specific 2d geometry of an object with multiple fields. I also know the data values for each field. How can I plot a heatmap with my values inside the fileds using the shape of my 2d object. I could get the CAD of the object for example, create a smiple drawing myself, or imitate my 2d shape with blocks or somehow else ….
Thanks. heatmap, plot, 2d, shape, individual, geometry, fields, labels MATLAB Answers — New Questions

Matlab News

I’m sending 170 buffered samples from an ESP32 with MPU6050 to ThingSpeak, but the total entries increase inconsistently with a few samples missing in each interval.

PuTI / 2025-06-18

I’m using an ESP32 to collect data from an MPU6050 sensor. The ESP32 buffers 170 samples at a time and uploads them to ThingSpeak using the bulk update API.
However, the data entries received on ThingSpeak are inconsistent. Here’s the pattern of total entries I see on the channel after each upload:
After 1st upload: 170 entries (correct)
After 2nd upload: 338 entries (should be 340)
After 3rd upload: 499 entries (should be 510)
After 4th upload: 657 entries (should be 680)
After 5th upload: 818 entries (should be 850)
As you can see, the expected total should increase by 170 each time, but it’s losing a few samples in each batch, and the number of missing entries increases over time. The losses are not consistent — it seems like 1 or 2 samples get dropped randomly in each interval.
Some things I’ve considered:
The ESP32 is collecting and sending the full buffer (verified via serial debug).
Could it be a payload size limit, HTTP timeout, or formatting issue?
Maybe ThingSpeak is silently dropping malformed lines from the bulk update?
Has anyone else run into this with buffered uploads to ThingSpeak? Any tips for debugging the exact cause of dropped samples?
Thanks for any help!I’m using an ESP32 to collect data from an MPU6050 sensor. The ESP32 buffers 170 samples at a time and uploads them to ThingSpeak using the bulk update API.
However, the data entries received on ThingSpeak are inconsistent. Here’s the pattern of total entries I see on the channel after each upload:
After 1st upload: 170 entries (correct)
After 2nd upload: 338 entries (should be 340)
After 3rd upload: 499 entries (should be 510)
After 4th upload: 657 entries (should be 680)
After 5th upload: 818 entries (should be 850)
As you can see, the expected total should increase by 170 each time, but it’s losing a few samples in each batch, and the number of missing entries increases over time. The losses are not consistent — it seems like 1 or 2 samples get dropped randomly in each interval.
Some things I’ve considered:
The ESP32 is collecting and sending the full buffer (verified via serial debug).
Could it be a payload size limit, HTTP timeout, or formatting issue?
Maybe ThingSpeak is silently dropping malformed lines from the bulk update?
Has anyone else run into this with buffered uploads to ThingSpeak? Any tips for debugging the exact cause of dropped samples?
Thanks for any help! I’m using an ESP32 to collect data from an MPU6050 sensor. The ESP32 buffers 170 samples at a time and uploads them to ThingSpeak using the bulk update API.
However, the data entries received on ThingSpeak are inconsistent. Here’s the pattern of total entries I see on the channel after each upload:
After 1st upload: 170 entries (correct)
After 2nd upload: 338 entries (should be 340)
After 3rd upload: 499 entries (should be 510)
After 4th upload: 657 entries (should be 680)
After 5th upload: 818 entries (should be 850)
As you can see, the expected total should increase by 170 each time, but it’s losing a few samples in each batch, and the number of missing entries increases over time. The losses are not consistent — it seems like 1 or 2 samples get dropped randomly in each interval.
Some things I’ve considered:
The ESP32 is collecting and sending the full buffer (verified via serial debug).
Could it be a payload size limit, HTTP timeout, or formatting issue?
Maybe ThingSpeak is silently dropping malformed lines from the bulk update?
Has anyone else run into this with buffered uploads to ThingSpeak? Any tips for debugging the exact cause of dropped samples?
Thanks for any help! thingspeak, buffer, burst MATLAB Answers — New Questions

Matlab News

sprintf vs. compose performance for large arrays

PuTI / 2025-06-18

My understanding is that compose is based on sprintf, yet I’m noticing a huge discrepancy in performance. Consider the following example.
M = 32;
N = 32;
A = randn(M,N);
[m,n] = ndgrid(1:M,1:N);
Nrep = 100;

%% compose and no for loops
t1 = zeros(Nrep,1);
for kr=1:Nrep
tic;
s1 = compose("%d,%d,%.3g",m(:),n(:),A(:));
s1 = reshape(s1,M,N);
t1(kr) = toc;
end

%% sprintf and two nested for loops
t2 = zeros(Nrep,1);
for kr=1:Nrep
s2 = repelem("",M,N);
tic;
for j=1:N
for i=1:M
s2(i,j) = sprintf("%d,%d,%.3g",i,j,A(i,j));
end
end
t2(kr) = toc;
end

%% compose and two nested for loops
t3 = zeros(Nrep,1);
for kr=1:Nrep
tic;
s3 = repelem("",M,N);
for j=1:N
for i=1:M
s3(i,j) = compose("%d,%d,%.3g",i,j,A(i,j));
end
end
t3(kr) = toc;
end

fprintf(" min | mean | maxn");
fprintf("compose, no loops: %.6f | %.6f | %.6fn",min(t1),mean(t1),max(t1));
fprintf("sprintf, 2 loops: %.6f | %.6f | %.6fn",min(t2),mean(t2),max(t2));
fprintf("compose, 2 loops: %.6f | %.6f | %.6fn",min(t3),mean(t3),max(t3));
This code produces the following output on my machine.
min | mean | max
compose, no loops: 1.008151 | 1.035783 | 1.076671
sprintf, 2 loops: 0.004265 | 0.004384 | 0.008990
compose, 2 loops: 1.018900 | 1.050923 | 1.464713
It seems that sprintf is 60x – 70x faster in this example. Any idea why that is?
Here is the output of ver on my machine:
—————————————————————————————–
MATLAB Version: 25.1.0.2943329 (R2025a)
MATLAB License Number:
Operating System: macOS Version: 15.5 Build: 24F74
Java Version: Java 11.0.27+6-LTS with Amazon.com Inc. OpenJDK 64-Bit Server VM mixed mode
—————————————————————————————–My understanding is that compose is based on sprintf, yet I’m noticing a huge discrepancy in performance. Consider the following example.
M = 32;
N = 32;
A = randn(M,N);
[m,n] = ndgrid(1:M,1:N);
Nrep = 100;

%% compose and no for loops
t1 = zeros(Nrep,1);
for kr=1:Nrep
tic;
s1 = compose("%d,%d,%.3g",m(:),n(:),A(:));
s1 = reshape(s1,M,N);
t1(kr) = toc;
end

%% sprintf and two nested for loops
t2 = zeros(Nrep,1);
for kr=1:Nrep
s2 = repelem("",M,N);
tic;
for j=1:N
for i=1:M
s2(i,j) = sprintf("%d,%d,%.3g",i,j,A(i,j));
end
end
t2(kr) = toc;
end

%% compose and two nested for loops
t3 = zeros(Nrep,1);
for kr=1:Nrep
tic;
s3 = repelem("",M,N);
for j=1:N
for i=1:M
s3(i,j) = compose("%d,%d,%.3g",i,j,A(i,j));
end
end
t3(kr) = toc;
end

fprintf(" min | mean | maxn");
fprintf("compose, no loops: %.6f | %.6f | %.6fn",min(t1),mean(t1),max(t1));
fprintf("sprintf, 2 loops: %.6f | %.6f | %.6fn",min(t2),mean(t2),max(t2));
fprintf("compose, 2 loops: %.6f | %.6f | %.6fn",min(t3),mean(t3),max(t3));
This code produces the following output on my machine.
min | mean | max
compose, no loops: 1.008151 | 1.035783 | 1.076671
sprintf, 2 loops: 0.004265 | 0.004384 | 0.008990
compose, 2 loops: 1.018900 | 1.050923 | 1.464713
It seems that sprintf is 60x – 70x faster in this example. Any idea why that is?
Here is the output of ver on my machine:
—————————————————————————————–
MATLAB Version: 25.1.0.2943329 (R2025a)
MATLAB License Number:
Operating System: macOS Version: 15.5 Build: 24F74
Java Version: Java 11.0.27+6-LTS with Amazon.com Inc. OpenJDK 64-Bit Server VM mixed mode
—————————————————————————————– My understanding is that compose is based on sprintf, yet I’m noticing a huge discrepancy in performance. Consider the following example.
M = 32;
N = 32;
A = randn(M,N);
[m,n] = ndgrid(1:M,1:N);
Nrep = 100;

%% compose and no for loops
t1 = zeros(Nrep,1);
for kr=1:Nrep
tic;
s1 = compose("%d,%d,%.3g",m(:),n(:),A(:));
s1 = reshape(s1,M,N);
t1(kr) = toc;
end

%% sprintf and two nested for loops
t2 = zeros(Nrep,1);
for kr=1:Nrep
s2 = repelem("",M,N);
tic;
for j=1:N
for i=1:M
s2(i,j) = sprintf("%d,%d,%.3g",i,j,A(i,j));
end
end
t2(kr) = toc;
end

%% compose and two nested for loops
t3 = zeros(Nrep,1);
for kr=1:Nrep
tic;
s3 = repelem("",M,N);
for j=1:N
for i=1:M
s3(i,j) = compose("%d,%d,%.3g",i,j,A(i,j));
end
end
t3(kr) = toc;
end

fprintf(" min | mean | maxn");
fprintf("compose, no loops: %.6f | %.6f | %.6fn",min(t1),mean(t1),max(t1));
fprintf("sprintf, 2 loops: %.6f | %.6f | %.6fn",min(t2),mean(t2),max(t2));
fprintf("compose, 2 loops: %.6f | %.6f | %.6fn",min(t3),mean(t3),max(t3));
This code produces the following output on my machine.
min | mean | max
compose, no loops: 1.008151 | 1.035783 | 1.076671
sprintf, 2 loops: 0.004265 | 0.004384 | 0.008990
compose, 2 loops: 1.018900 | 1.050923 | 1.464713
It seems that sprintf is 60x – 70x faster in this example. Any idea why that is?
Here is the output of ver on my machine:
—————————————————————————————–
MATLAB Version: 25.1.0.2943329 (R2025a)
MATLAB License Number:
Operating System: macOS Version: 15.5 Build: 24F74
Java Version: Java 11.0.27+6-LTS with Amazon.com Inc. OpenJDK 64-Bit Server VM mixed mode
—————————————————————————————– string, strings, sprintf, compose, text MATLAB Answers — New Questions

Matlab News

defining function which generates 3d array within a class. Different behavior in class versus command window

PuTI / 2025-06-18

This code:
lam = @(t) 3*(1 + 0.8*cos(2*pi*t));
A1 = @(t) reshape([lam(t(:).’); lam(t(:).’); zeros(1, numel(t)); lam(t(:).’)], [2, 2, numel(t)]);
tt=0:1/400:1-1/400;
A1stack=A1(tt);
generates a 2x2x400 array. This is what I want to happen. When I embed the code within a class function however, it returns a 2×800 array. How do I fix this?This code:
lam = @(t) 3*(1 + 0.8*cos(2*pi*t));
A1 = @(t) reshape([lam(t(:).’); lam(t(:).’); zeros(1, numel(t)); lam(t(:).’)], [2, 2, numel(t)]);
tt=0:1/400:1-1/400;
A1stack=A1(tt);
generates a 2x2x400 array. This is what I want to happen. When I embed the code within a class function however, it returns a 2×800 array. How do I fix this? This code:
lam = @(t) 3*(1 + 0.8*cos(2*pi*t));
A1 = @(t) reshape([lam(t(:).’); lam(t(:).’); zeros(1, numel(t)); lam(t(:).’)], [2, 2, numel(t)]);
tt=0:1/400:1-1/400;
A1stack=A1(tt);
generates a 2x2x400 array. This is what I want to happen. When I embed the code within a class function however, it returns a 2×800 array. How do I fix this? array function, 3d array, oop MATLAB Answers — New Questions

Matlab News

Word document “Saveas2()” method no longer working

PuTI / 2025-06-18

As of the past few days, I seem to be no longer able to save word documents created through matlab. Several months ago I ran the following script which opens an existing word document and saves it as another document.
reportPath = ‘Table_out.docx’;
templateName = append(pwd, ‘Table_Inject_Test.docx’);
actx_word = actxserver(‘Word.Application’);
actx_word.Visible = true;
trace(actx_word.Visible);
wordTemplate = actx_word.Documents.Open(templateName); % Open template
wordTemplate.SaveAs2(reportPath); % Save copy as report
Once I started getting this error in another project, I reverted back to this script to check it wasn’t something I had done within the word template to cause the error.
However when I try to run the same script now, I get the following error:
% Unrecognized method, property, or field ‘SaveAs2’ for class
% ‘Interface.0002096B_0000_0000_C000_000000000046’.

% Error in tableinjecttest2 (line 7)
% wordTemplate.SaveAs2(reportPath); % Save copy as report

Has something changed within the WordInterface object? I see there is a ‘saveobj’ method but that doesn’t seem to save the document either. That just gives me a different error:

% Unable to resolve the name ‘wordParent.saveobj’.
%
% Error in tableInjectTest (line 22)
% wordTemplate.saveobj(reportPath)As of the past few days, I seem to be no longer able to save word documents created through matlab. Several months ago I ran the following script which opens an existing word document and saves it as another document.
reportPath = ‘Table_out.docx’;
templateName = append(pwd, ‘Table_Inject_Test.docx’);
actx_word = actxserver(‘Word.Application’);
actx_word.Visible = true;
trace(actx_word.Visible);
wordTemplate = actx_word.Documents.Open(templateName); % Open template
wordTemplate.SaveAs2(reportPath); % Save copy as report
Once I started getting this error in another project, I reverted back to this script to check it wasn’t something I had done within the word template to cause the error.
However when I try to run the same script now, I get the following error:
% Unrecognized method, property, or field ‘SaveAs2’ for class
% ‘Interface.0002096B_0000_0000_C000_000000000046’.

% Error in tableinjecttest2 (line 7)
% wordTemplate.SaveAs2(reportPath); % Save copy as report

Has something changed within the WordInterface object? I see there is a ‘saveobj’ method but that doesn’t seem to save the document either. That just gives me a different error:

% Unable to resolve the name ‘wordParent.saveobj’.
%
% Error in tableInjectTest (line 22)
% wordTemplate.saveobj(reportPath) As of the past few days, I seem to be no longer able to save word documents created through matlab. Several months ago I ran the following script which opens an existing word document and saves it as another document.
reportPath = ‘Table_out.docx’;
templateName = append(pwd, ‘Table_Inject_Test.docx’);
actx_word = actxserver(‘Word.Application’);
actx_word.Visible = true;
trace(actx_word.Visible);
wordTemplate = actx_word.Documents.Open(templateName); % Open template
wordTemplate.SaveAs2(reportPath); % Save copy as report
Once I started getting this error in another project, I reverted back to this script to check it wasn’t something I had done within the word template to cause the error.
However when I try to run the same script now, I get the following error:
% Unrecognized method, property, or field ‘SaveAs2’ for class
% ‘Interface.0002096B_0000_0000_C000_000000000046’.

% Error in tableinjecttest2 (line 7)
% wordTemplate.SaveAs2(reportPath); % Save copy as report

Has something changed within the WordInterface object? I see there is a ‘saveobj’ method but that doesn’t seem to save the document either. That just gives me a different error:

% Unable to resolve the name ‘wordParent.saveobj’.
%
% Error in tableInjectTest (line 22)
% wordTemplate.saveobj(reportPath) word application, actxserver MATLAB Answers — New Questions

News

Microsoft Pushes European Sovereign Solutions

Tony Redmond / 2025-06-18

Marked Lack of Detail around Microsoft 365 Local

Microsoft’s June 16 announcement about “sovereign solutions empowering European organizations” (Figure 1) is obviously an attempt by Microsoft to reassure European customers that continuing to use Microsoft (U.S.-based) technology is a safe choice at a time when many question the policies of the current U.S. administration.

Microsoft sovereign clouds, including Microsoft 365 Local. — Figure 1: Microsoft sovereign clouds, including Microsoft 365 Local (source: Microsoft)

To be fair to Microsoft, they’ve been on the path to respect data sovereignty for many years, starting with the original “Black Forest” implementation of Office 365 for German customers to a point where multiple national-level datacenter regions are available within Europe. Microsoft’s continued efforts to provide comfort to customers who want their data stored in-country and under the control of European law is commendable.

However, the announcement of Microsoft 365 Local confused everyone. According to the announcement, “Microsoft 365 Local provides customers with additional choice by bringing together Microsoft’s productivity server software into an Azure Local environment that can run entirely in a customer’s own datacenter.”

Apart from the Name, No Trace of Microsoft 365

Applying the Microsoft 365 branding to the offering implies some form of connection to Microsoft 365. But apart from a need to connect to Azure., this solution seems to have nothing much to do with Microsoft 365 cloud services. Instead, it appears to be the on-premises versions of Exchange Server, SharePoint Server, and Skype for Business Server running on an Azure Local instance, defined as “a machine or a cluster of machines running the Azure Stack HCI operating system and connected to Azure.”

At this point, Microsoft hasn’t shared details of how the services connect together, but I assume that Active Directory is in the mix too. We also don’t know if the Azure-based local infrastructure operates as a separate deployment, can be integrated into an existing on-premises organization, or operate as part of a hybrid organization.

In other words, Microsoft 365 Local is a modernized example of a packaged Azure-based installation of Exchange, SharePoint, and Skype for Business built according to a reference architecture and accessed via the same kind of clients that people use today to connect to on-premises servers. Unsurprisingly, Microsoft 365 Local doesn’t include Teams because Teams relies so heavily on services from Exchange, SharePoint, OneDrive, Planner, and a bunch of Azure microservices.

The packaging might be innovative, and Microsoft marketing will certainly call the announcement a triumph for branding, but it has nothing to do with Microsoft 365. Anyone who steps back from using Exchange Online with its close integration with SharePoint Online will quickly discover how different things are.

Some Organizations Will Love Microsoft 365 Local

Although I hate the name, a place exists for a solution like Microsoft 365 Local. Some companies want to control their own destiny, which is why they continue running on-premises software; others don’t have sufficient external network capacity to be dependent on cloud services.

Other companies simply want to not have to deal with the blizzard of changes that Microsoft 365 customers have to cope with, or the constant nagging from Microsoft to adopt and use its AI-based tools like Microsoft 365 Copilot. European customers have a strong track record of respecting user privacy, and solutions like the recently-launched AI-powered People Skills are unlikely to be popular with unions or works councils.

Being able to purchase a packaged solution that is hopefully better integrated out-of-the-box is a nicer option than having to convince Exchange Server and SharePoint Server (for instance) to work together, an exercise that is usually guaranteed to frustrate. Presumably the solution leverages the subscription version of the three on-premises servers and will be paid for via an Azure subscription in the same manner as Azure Local.

Lack of Detail is Frustrating

The trouble is the total lack of detail currently available about Microsoft 365 Local. The above is inspired guesswork based on reading between the lines of Microsoft’s announcement. Many questions remain unanswered. Customers will need pricing and availability details from the various hardware vendors listed in the announcement are before they can decide if Microsoft 365 Local is for them. Migration from current on-premises deployments is another issue to resolve as is deployment alongside existing deployments.

The lack of detail is frustrating, but this is a classic marketing playbook: announce a product to gauge interest and follow up if the interest is there. It will be interesting to see what Microsoft 365 Local can deliver and at what cost.

Matlab News

Finding roots of a complex function

PuTI / 2025-06-17

Hi,
I’m trying to find the roots of a complex function z=x+iyz in MATLAB. However, it seems that the root-finding routine is missing some roots, which leads to inaccurate results. Could you advise on how I can improve the code to ensure more reliable root detection? I have the two matlab files.
Thanks,
clearvars
close all

% set parameter values
pars.gamma1=0.1093;
pars.alpha3=-0.1104e-2;
pars.K1=6e-12;
pars.d=0.2e-3;
pars.eta1=0.240e-1;
pars.chia=1.219e-6;
pars.alpha=1-pars.alpha3^2/(pars.gamma1*pars.eta1);
pars.Ha=pi*sqrt(pars.K1/pars.chia)/pars.d;

% set lists of u (field) and xi (activity) values
uvals=0:0.05:3;
xivals=-0.3:0.01:0.3;

nu=length(uvals);
nxi=length(xivals);

% initiate arrays for output
taumin=ones(nu,nxi);
wavenummin=ones(nu,nxi);

% start timer
tic
disp(‘Starting u and xi loops’);

% start loop around u values
for i=1:nu

pars.H=uvals(i)*pars.Ha;

% start loop around xi value
for j=1:nxi

xi=xivals(j);

% set initial tau values for root finding, tau is a complex
% variable
tauRvals=-50:0.1:50;
tauIvals=0.1*ones(size(tauRvals));
ntau=length(tauRvals);

% plot equation (projected onto real line) to solve for these values of u and xi
fig1 = figure(1);
y=zeros(size(tauRvals));
for ii=1:ntau
tau=tauRvals(ii);
y(ii) = (pars.H ^ 2 * sin(pars.K1 ^ (-0.1e1 / 0.2e1) * pars.eta1 ^ (-0.1e1 / 0.2e1) * tau ^ (-0.1e1 / 0.2e1) * sqrt(pars.H ^ 2 * pars.chia * pars.eta1 * tau + pars.alpha3 * tau * xi – pars.alpha3 ^ 2 + pars.eta1 * pars.gamma1) * pars.d) * pars.chia * pars.d * tau ^ (0.3e1 / 0.2e1) * sqrt(pars.eta1) * sqrt(pars.H ^ 2 * pars.chia * pars.eta1 * tau + pars.alpha3 * tau * xi – pars.alpha3 ^ 2 + pars.eta1 * pars.gamma1) + sin(pars.K1 ^ (-0.1e1 / 0.2e1) * pars.eta1 ^ (-0.1e1 / 0.2e1) * tau ^ (-0.1e1 / 0.2e1) * sqrt(pars.H ^ 2 * pars.chia * pars.eta1 * tau + pars.alpha3 * tau * xi – pars.alpha3 ^ 2 + pars.eta1 * pars.gamma1) * pars.d) * pars.d * pars.gamma1 * sqrt(pars.eta1) * sqrt(tau) * sqrt(pars.H ^ 2 * pars.chia * pars.eta1 * tau + pars.alpha3 * tau * xi – pars.alpha3 ^ 2 + pars.eta1 * pars.gamma1) – 0.2e1 * sqrt(pars.K1) * cos(pars.K1 ^ (-0.1e1 / 0.2e1) * pars.eta1 ^ (-0.1e1 / 0.2e1) * tau ^ (-0.1e1 / 0.2e1) * sqrt(pars.H ^ 2 * pars.chia * pars.eta1 * tau + pars.alpha3 * tau * xi – pars.alpha3 ^ 2 + pars.eta1 * pars.gamma1) * pars.d) * pars.alpha3 * tau ^ 2 * xi + 0.2e1 * sqrt(pars.K1) * cos(pars.K1 ^ (-0.1e1 / 0.2e1) * pars.eta1 ^ (-0.1e1 / 0.2e1) * tau ^ (-0.1e1 / 0.2e1) * sqrt(pars.H ^ 2 * pars.chia * pars.eta1 * tau + pars.alpha3 * tau * xi – pars.alpha3 ^ 2 + pars.eta1 * pars.gamma1) * pars.d) * pars.alpha3 ^ 2 * tau + 0.2e1 * sqrt(pars.K1) * pars.alpha3 * tau ^ 2 * xi – 0.2e1 * sqrt(pars.K1) * pars.alpha3 ^ 2 * tau) * pars.eta1 ^ (-0.1e1 / 0.2e1) * tau ^ (-0.3e1 / 0.2e1) * (pars.H ^ 2 * pars.chia * pars.eta1 * tau + pars.alpha3 * tau * xi – pars.alpha3 ^ 2 + pars.eta1 * pars.gamma1) ^ (-0.1e1 / 0.2e1) / pars.alpha3 / sin(pars.K1 ^ (-0.1e1 / 0.2e1) * pars.eta1 ^ (-0.1e1 / 0.2e1) * tau ^ (-0.1e1 / 0.2e1) * sqrt(pars.H ^ 2 * pars.chia * pars.eta1 * tau + pars.alpha3 * tau * xi – pars.alpha3 ^ 2 + pars.eta1 * pars.gamma1) * pars.d);
end
plot(tauRvals,y);
hold on
plot(tauRvals,zeros(size(tauRvals)),’-k’);
xline(0);
yline(0);
xlabel(‘tau’);
ylabel(‘tau equation’);
axis([min(tauRvals) max(tauRvals) -1e-2 1e-2]);
drawnow
hold off

% loop around initial tau values for root finding
tausol=zeros(1,ntau);
flag=zeros(1,ntau);
wavenum=zeros(1,ntau);
for k=1:ntau

% set function, options and inital tau value
fun = @(x)rootsolver_complex(x,xi,pars);
options = optimset(‘TolFun’,1e-15,’MaxFunEvals’,1e5,’Maxiter’,1e5,’Display’,’none’);
tauinit=[tauRvals(k),tauIvals(k)];

% find root of tau equation
[x,fval,exitflag,output] = fsolve(fun,tauinit,options);
% save complex tau solution
tausol(k)=complex(x(1),x(2));
% set solve flag (if exitflag>0 the root finder has solved)
flag(k)=(exitflag>0);
% calulate wavenumber (real part of) using equation from Maple
% file
wavenum(k)=real(sqrt((pars.H ^ 2 * pars.chia * pars.eta1 * tau + pars.alpha3 * tau * xi – pars.alpha3 ^ 2 + pars.eta1 * pars.gamma1) / pars.K1 / pars.eta1 / tau) * pars.d / pi / 0.2e1);
end

tauflag=[tausol’,flag’];
tausol_found=tauflag(flag==1);
tausolR=real(tausol_found);
tauRvals=tauRvals(flag==1);
wavenum=wavenum(flag==1);

% plot real part of tau solution versus initial tau
fig2 = figure(2);
plot(tauRvals,tausolR,’g-‘)
hold on
plot(tauRvals,tauRvals,’k-‘)
xlabel(‘initial $tau$’, ‘Interpreter’,’latex’);
ylabel(‘$tau$ solution’, ‘Interpreter’,’latex’);

hold off

% calculate min value of real part of tau and the wavenumber at
% that min value of tau
taumin(i,j)=min(tausolR);
wavenummins=wavenum(tausolR==min(tausolR));
wavenummin(i,j)=wavenummins(1);

% filled contour plot of minimum tau value (negative tau means instability)
fig3 = figure(3);
[Xi,U] = meshgrid(xivals,uvals);
N=[0:0.1:1];
map = [0.95*(1-N’) 0.95*(1-N’) N’];
contourf(U,Xi,taumin,[-100:10:100])
colormap(map)
colorbar
xlabel(‘Orienting field, $u$’, ‘Interpreter’,’latex’);
ylabel(‘Activity, $xi$ [Pa]’, ‘Interpreter’,’latex’);
%title(‘minimum tau’);
drawnow

% filled contour plot of minimum tau value (negative tau means instability)
fig4 =figure(4);
contourf(U,Xi,wavenummin,10)
colormap(map)
colorbar
xlabel(‘Orienting field, $u$’, ‘Interpreter’,’latex’);
ylabel(‘Activity, $xi$ [Pa]’, ‘Interpreter’,’latex’);
%title(‘wavenumber at minimum tau’);
drawnow

% stability domain in (u,xi) plane
fig5 = figure(5);
S = 25; % size of symbols in pixels
% normalize colouring vector to go from zero to 1
normtau = (taumin>0);
normtau=reshape(normtau,nu*nxi,1);
C = [0.95*(1-normtau) 0.95*(1-normtau) normtau];
scatter(reshape(U,nu*nxi,1),reshape(Xi,nu*nxi,1),S,C,’filled’,’Marker’,’o’)
xlabel(‘Orienting field, $u$’, ‘Interpreter’,’latex’);
ylabel(‘Activity, $xi$ [Pa]’, ‘Interpreter’,’latex’);
title(‘Blue = stable, Yellow = unstable’, ‘Interpreter’,’latex’);
drawnow

end

% display time taken and percentage complete
toc
disp([‘Progress: ‘ num2str(round(100*(i*(j-1)+j)/(nu*nxi))) ‘ % completed’]);
end

function F = rootsolver_complex(x,xi,pars)
% function to provide right-hand-side of the equation for tau

gamma1=pars.gamma1;
alpha3=pars.alpha3;
K1=pars.K1;
d=pars.d;
eta1=pars.eta1;
chia=pars.chia;
alpha=pars.alpha;
H=pars.H;

tau=complex(x(1),x(2));

% equation taken directly from Maple file eq.mw
y = (H ^ 2 * sin(K1 ^ (-0.1e1 / 0.2e1) * eta1 ^ (-0.1e1 / 0.2e1) * tau ^ (-0.1e1 / 0.2e1) * sqrt(H ^ 2 * chia * eta1 * tau + alpha3 * tau * xi – alpha3 ^ 2 + eta1 * gamma1) * d) * chia * d * tau ^ (0.3e1 / 0.2e1) * sqrt(eta1) * sqrt(H ^ 2 * chia * eta1 * tau + alpha3 * tau * xi – alpha3 ^ 2 + eta1 * gamma1) + sin(K1 ^ (-0.1e1 / 0.2e1) * eta1 ^ (-0.1e1 / 0.2e1) * tau ^ (-0.1e1 / 0.2e1) * sqrt(H ^ 2 * chia * eta1 * tau + alpha3 * tau * xi – alpha3 ^ 2 + eta1 * gamma1) * d) * d * gamma1 * sqrt(eta1) * sqrt(tau) * sqrt(H ^ 2 * chia * eta1 * tau + alpha3 * tau * xi – alpha3 ^ 2 + eta1 * gamma1) – 0.2e1 * sqrt(K1) * cos(K1 ^ (-0.1e1 / 0.2e1) * eta1 ^ (-0.1e1 / 0.2e1) * tau ^ (-0.1e1 / 0.2e1) * sqrt(H ^ 2 * chia * eta1 * tau + alpha3 * tau * xi – alpha3 ^ 2 + eta1 * gamma1) * d) * alpha3 * tau ^ 2 * xi + 0.2e1 * sqrt(K1) * cos(K1 ^ (-0.1e1 / 0.2e1) * eta1 ^ (-0.1e1 / 0.2e1) * tau ^ (-0.1e1 / 0.2e1) * sqrt(H ^ 2 * chia * eta1 * tau + alpha3 * tau * xi – alpha3 ^ 2 + eta1 * gamma1) * d) * alpha3 ^ 2 * tau + 0.2e1 * sqrt(K1) * alpha3 * tau ^ 2 * xi – 0.2e1 * sqrt(K1) * alpha3 ^ 2 * tau) * eta1 ^ (-0.1e1 / 0.2e1) * tau ^ (-0.3e1 / 0.2e1) * (H ^ 2 * chia * eta1 * tau + alpha3 * tau * xi – alpha3 ^ 2 + eta1 * gamma1) ^ (-0.1e1 / 0.2e1) / alpha3 / sin(K1 ^ (-0.1e1 / 0.2e1) * eta1 ^ (-0.1e1 / 0.2e1) * tau ^ (-0.1e1 / 0.2e1) * sqrt(H ^ 2 * chia * eta1 * tau + alpha3 * tau * xi – alpha3 ^ 2 + eta1 * gamma1) * d);
F=[real(y),imag(y)];

endHi,
I’m trying to find the roots of a complex function z=x+iyz in MATLAB. However, it seems that the root-finding routine is missing some roots, which leads to inaccurate results. Could you advise on how I can improve the code to ensure more reliable root detection? I have the two matlab files.
Thanks,
clearvars
close all

% set lists of u (field) and xi (activity) values
uvals=0:0.05:3;
xivals=-0.3:0.01:0.3;

nu=length(uvals);
nxi=length(xivals);

% initiate arrays for output
taumin=ones(nu,nxi);
wavenummin=ones(nu,nxi);

% start timer
tic
disp(‘Starting u and xi loops’);

% start loop around u values
for i=1:nu

pars.H=uvals(i)*pars.Ha;

% start loop around xi value
for j=1:nxi

xi=xivals(j);

% set initial tau values for root finding, tau is a complex
% variable
tauRvals=-50:0.1:50;
tauIvals=0.1*ones(size(tauRvals));
ntau=length(tauRvals);

% loop around initial tau values for root finding
tausol=zeros(1,ntau);
flag=zeros(1,ntau);
wavenum=zeros(1,ntau);
for k=1:ntau

tauflag=[tausol’,flag’];
tausol_found=tauflag(flag==1);
tausolR=real(tausol_found);
tauRvals=tauRvals(flag==1);
wavenum=wavenum(flag==1);

hold off

% calculate min value of real part of tau and the wavenumber at
% that min value of tau
taumin(i,j)=min(tausolR);
wavenummins=wavenum(tausolR==min(tausolR));
wavenummin(i,j)=wavenummins(1);

end

% display time taken and percentage complete
toc
disp([‘Progress: ‘ num2str(round(100*(i*(j-1)+j)/(nu*nxi))) ‘ % completed’]);
end

function F = rootsolver_complex(x,xi,pars)
% function to provide right-hand-side of the equation for tau

gamma1=pars.gamma1;
alpha3=pars.alpha3;
K1=pars.K1;
d=pars.d;
eta1=pars.eta1;
chia=pars.chia;
alpha=pars.alpha;
H=pars.H;

tau=complex(x(1),x(2));

end Hi,
I’m trying to find the roots of a complex function z=x+iyz in MATLAB. However, it seems that the root-finding routine is missing some roots, which leads to inaccurate results. Could you advise on how I can improve the code to ensure more reliable root detection? I have the two matlab files.
Thanks,
clearvars
close all

% set lists of u (field) and xi (activity) values
uvals=0:0.05:3;
xivals=-0.3:0.01:0.3;

nu=length(uvals);
nxi=length(xivals);

% initiate arrays for output
taumin=ones(nu,nxi);
wavenummin=ones(nu,nxi);

% start timer
tic
disp(‘Starting u and xi loops’);

% start loop around u values
for i=1:nu

pars.H=uvals(i)*pars.Ha;

% start loop around xi value
for j=1:nxi

xi=xivals(j);

% set initial tau values for root finding, tau is a complex
% variable
tauRvals=-50:0.1:50;
tauIvals=0.1*ones(size(tauRvals));
ntau=length(tauRvals);

% loop around initial tau values for root finding
tausol=zeros(1,ntau);
flag=zeros(1,ntau);
wavenum=zeros(1,ntau);
for k=1:ntau

tauflag=[tausol’,flag’];
tausol_found=tauflag(flag==1);
tausolR=real(tausol_found);
tauRvals=tauRvals(flag==1);
wavenum=wavenum(flag==1);

hold off

% calculate min value of real part of tau and the wavenumber at
% that min value of tau
taumin(i,j)=min(tausolR);
wavenummins=wavenum(tausolR==min(tausolR));
wavenummin(i,j)=wavenummins(1);

end

% display time taken and percentage complete
toc
disp([‘Progress: ‘ num2str(round(100*(i*(j-1)+j)/(nu*nxi))) ‘ % completed’]);
end

function F = rootsolver_complex(x,xi,pars)
% function to provide right-hand-side of the equation for tau

gamma1=pars.gamma1;
alpha3=pars.alpha3;
K1=pars.K1;
d=pars.d;
eta1=pars.eta1;
chia=pars.chia;
alpha=pars.alpha;
H=pars.H;

tau=complex(x(1),x(2));

end roots finding, complex finding MATLAB Answers — New Questions

Previous 1 … 11 12 13 14 15 … 56 Next

Cart

Cart

Category: News

cooperative logistics – profit maximization

I am trying to use the PCBwriter function to output a Gerber file for PCB stack object with greater than 2 layers.

Compare same frequency components of two signals

Interp1 is not working after applying unique because of rounding off

Readtable and Readmatrix Ignore Specified Range and Produce Extra Variables

Clean and automate indents in scripts an live scripts

Varying device latency for ASIO Device

COnnecting points on opposite end of mask

Mass importing .txt files with sequential names to plot multiple graphs on same figure for easy comparison

3D surf plot for more than two quantities

using spice models in simlink

Microsoft to Block Third-Party App Access to User Sites and Files

New Microsoft-Managed App Consent Policy to Control User Consent for Apps

No Problem with Blocking Obsolete and Insecure Protocols

App Consent Policy to Prevent Third-Party Access to Files and Sites

The ChatGPT Conundrum

Don’t Drop Your Guard

Updating the Entra ID Custom Banned Password List with PowerShell

Use Microsoft Graph PowerShell SDK Cmdlets to Maintain the Entra ID Custom Banned Password List

Global Banned Password List

Custom Banned Password List

Blocking Passwords

Updating the Custom Banned Password List with a Script

Heatmap inside custom 2d geometry shape

I’m sending 170 buffered samples from an ESP32 with MPU6050 to ThingSpeak, but the total entries increase inconsistently with a few samples missing in each interval.

sprintf vs. compose performance for large arrays

defining function which generates 3d array within a class. Different behavior in class versus command window

Word document “Saveas2()” method no longer working

Microsoft Pushes European Sovereign Solutions

Marked Lack of Detail around Microsoft 365 Local

Apart from the Name, No Trace of Microsoft 365

Some Organizations Will Love Microsoft 365 Local

Lack of Detail is Frustrating

Finding roots of a complex function

Information

Contact Us

All Categories

Search

Cart

All Categories

Search

Cart

New Microsoft-Managed App Consent Policy to Control User Consent for Apps

No Problem with Blocking Obsolete and Insecure Protocols

App Consent Policy to Prevent Third-Party Access to Files and Sites

The ChatGPT Conundrum

Don’t Drop Your Guard

Use Microsoft Graph PowerShell SDK Cmdlets to Maintain the Entra ID Custom Banned Password List

Global Banned Password List

Custom Banned Password List

Blocking Passwords

Updating the Custom Banned Password List with a Script

Marked Lack of Detail around Microsoft 365 Local

Apart from the Name, No Trace of Microsoft 365

Some Organizations Will Love Microsoft 365 Local

Lack of Detail is Frustrating

Sign Up For Newsletters

Information

Contact Us