by using dataset for heat exchanger with model equation i written code but its not generating correct estimation values

clc
clear
close all

%% ===============================
% READ HEAT EXCHANGER DATASET
%% ===============================

data = readmatrix(‘exchanger.xlsx’);

% Column 1 = Input
g = data(:,1);

% Column 2 = Output
x = data(:,2);

% Number of estimation samples
M = 3863;

% Select first M samples
g = g(1:M);
x = x(1:M);

% Convert to column vectors
g = g(:);
x = x(:);

%% ===============================
% PSO PARAMETERS
%% ===============================

dim = 8;

nParticles = 100;
maxIter = 800;

lb = [-2.0 0.3 0.0 -0.4 -40 30 -30 55];
ub = [-1.0 0.9 0.4 -0.1 -25 50 -20 75];

%% ===============================
% INITIALIZATION
%% ===============================

particle = struct;

for i = 1:nParticles

particle(i).position = lb + rand(1,dim).*(ub-lb);

particle(i).velocity = zeros(1,dim);

particle(i).cost = fitness_function( …
particle(i).position,x,g,M);

particle(i).best.position = particle(i).position;

particle(i).best.cost = particle(i).cost;

end

costs = [particle.cost];

[~,idx] = min(costs);

global_best = particle(idx).best;

%% ===============================
% PSO PARAMETERS
%% ===============================

w = 0.6;
c1 = 2.5;
c2 = 2.5;

BestCost = zeros(maxIter,1);
theta_history = zeros(maxIter,dim);
mse_history = zeros(maxIter,1);

%% ===============================
% PSO MAIN LOOP
%% ===============================

for iter = 1:maxIter

for i = 1:nParticles

particle(i).velocity = …
w*particle(i).velocity …
+ c1*rand*(particle(i).best.position-particle(i).position) …
+ c2*rand*(global_best.position-particle(i).position);

particle(i).position = …
particle(i).position + particle(i).velocity;

particle(i).position = …
max(particle(i).position,lb);

particle(i).position = …
min(particle(i).position,ub);

particle(i).cost = fitness_function( …
particle(i).position,x,g,M);

if particle(i).cost < particle(i).best.cost

particle(i).best.position = particle(i).position;

particle(i).best.cost = particle(i).cost;

end

if particle(i).best.cost < global_best.cost

global_best = particle(i).best;

end

end

theta_history(iter,:) = global_best.position;

mse_history(iter) = global_best.cost;

BestCost(iter) = global_best.cost;

fprintf(‘Iteration %d Best Error = %en’, …
iter,global_best.cost);

end
%% ===============================
% Estimated Parameters
%% ===============================

xext = global_best.position;

disp(‘Estimated Parameters’)
disp(xext)

true_theta = [-1.6080 0.6385 0.2070 -0.1764 -31.55 41.73 -24.20 68.03];

disp(‘True Parameters’)
disp(true_theta)
%% ===============================
% Parameter Extraction
%% ===============================

sigma1 = xext(1);
sigma2 = xext(2);

delta1 = xext(3);
delta2 = xext(4);

l1 = xext(5);
l2 = xext(6);
l3 = xext(7);
l4 = xext(8);
%% ===============================
% Nonlinear Block
%% ===============================

khat = l1*g + l2*g.^2 + l3*g.^3 + l4*g.^4;
%% ===============================
% Estimated Output
%% ===============================

xhat = zeros(M,1);

xhat(1:5) = x(1:5);

for t = 6:M

xhat(t) = …
sigma1*(xhat(t-1)-0.3*xhat(t-2)-0.105*xhat(t-3)-0.0595*xhat(t-4)) …
+ sigma2*(xhat(t-2)-0.3*xhat(t-3)-0.105*xhat(t-4)-0.0595*xhat(t-5)) …
+ delta1*(khat(t-1)-0.3*khat(t-2)-0.105*khat(t-3)-0.0595*khat(t-4)) …
+ delta2*(khat(t-2)-0.3*khat(t-3)-0.105*khat(t-4)-0.0595*khat(t-5));

end

%% ===============================
% Performance
%% ===============================
disp(‘Best Cost from PSO’)
disp(global_best.cost)

err = x – xhat;

valid = isfinite(err) & (abs(err) < 1e6);

if any(valid)
mse_value = mean(err(valid).^2);
else
mse_value = Inf;
end

fprintf(‘Final MSE = %en’, mse_value);

%% ===============================
% Actual vs Estimated Output
%% ===============================

figure

plot(x,’b’,’LineWidth’,1.5)
hold on
plot(xhat,’r–‘,’LineWidth’,1.5)

xlabel(‘Samples’)
ylabel(‘Output’)
title(‘Actual vs Estimated Output’)

legend(‘Actual Output’,’Estimated Output’)
grid on

function j = fitness_function(xext,x,g,M)

%% ==========================================
% Column vectors
%% ==========================================

x = x(:);
g = g(:);

%% ==========================================
% Parameters
%% ==========================================

sigma1 = xext(1);
sigma2 = xext(2);

delta1 = xext(3);
delta2 = xext(4);

l1 = xext(5);
l2 = xext(6);
l3 = xext(7);
l4 = xext(8);

%% ==========================================
% Nonlinear Block
%% ==========================================

khat = l1*g + l2*g.^2 + l3*g.^3 + l4*g.^4;

%% ==========================================
% Free Run Simulation
%% ==========================================

xhat = zeros(M,1);

% Initial conditions
xhat(1:5) = x(1:5);

for t = 6:M

xhat(t) = …
sigma1*(xhat(t-1) …
-0.3*xhat(t-2) …
-0.105*xhat(t-3) …
-0.0595*xhat(t-4)) …
+ sigma2*(xhat(t-2) …
-0.3*xhat(t-3) …
-0.105*xhat(t-4) …
-0.0595*xhat(t-5)) …
+ delta1*(khat(t-1) …
-0.3*khat(t-2) …
-0.105*khat(t-3) …
-0.0595*khat(t-4)) …
+ delta2*(khat(t-2) …
-0.3*khat(t-3) …
-0.105*khat(t-4) …
-0.0595*khat(t-5));

end

%% ==========================================
% Mean Squared Error
%% ==========================================

err = x – xhat;

% Keep only reasonable errors
valid = isfinite(err) & (abs(err) < 1e6);

if sum(valid) < 0.80*M

j = 1e20;

else

err = err(valid);

j = mean(err.^2);

end

if isnan(j) || isinf(j)

j = 1e20;

end

endclc
clear
close all

%% ===============================
% READ HEAT EXCHANGER DATASET
%% ===============================

data = readmatrix(‘exchanger.xlsx’);

% Column 1 = Input
g = data(:,1);

% Column 2 = Output
x = data(:,2);

% Number of estimation samples
M = 3863;

% Select first M samples
g = g(1:M);
x = x(1:M);

% Convert to column vectors
g = g(:);
x = x(:);

%% ===============================
% PSO PARAMETERS
%% ===============================

dim = 8;

nParticles = 100;
maxIter = 800;

lb = [-2.0 0.3 0.0 -0.4 -40 30 -30 55];
ub = [-1.0 0.9 0.4 -0.1 -25 50 -20 75];

%% ===============================
% INITIALIZATION
%% ===============================

particle = struct;

for i = 1:nParticles

particle(i).position = lb + rand(1,dim).*(ub-lb);

particle(i).velocity = zeros(1,dim);

particle(i).cost = fitness_function( …
particle(i).position,x,g,M);

particle(i).best.position = particle(i).position;

particle(i).best.cost = particle(i).cost;

end

costs = [particle.cost];

[~,idx] = min(costs);

global_best = particle(idx).best;

%% ===============================
% PSO PARAMETERS
%% ===============================

w = 0.6;
c1 = 2.5;
c2 = 2.5;

BestCost = zeros(maxIter,1);
theta_history = zeros(maxIter,dim);
mse_history = zeros(maxIter,1);

%% ===============================
% PSO MAIN LOOP
%% ===============================

for iter = 1:maxIter

for i = 1:nParticles

particle(i).velocity = …
w*particle(i).velocity …
+ c1*rand*(particle(i).best.position-particle(i).position) …
+ c2*rand*(global_best.position-particle(i).position);

particle(i).position = …
particle(i).position + particle(i).velocity;

particle(i).position = …
max(particle(i).position,lb);

particle(i).position = …
min(particle(i).position,ub);

particle(i).cost = fitness_function( …
particle(i).position,x,g,M);

if particle(i).cost < particle(i).best.cost

particle(i).best.position = particle(i).position;

particle(i).best.cost = particle(i).cost;

end

if particle(i).best.cost < global_best.cost

global_best = particle(i).best;

end

end

theta_history(iter,:) = global_best.position;

mse_history(iter) = global_best.cost;

BestCost(iter) = global_best.cost;

fprintf(‘Iteration %d Best Error = %en’, …
iter,global_best.cost);

end
%% ===============================
% Estimated Parameters
%% ===============================

xext = global_best.position;

disp(‘Estimated Parameters’)
disp(xext)

true_theta = [-1.6080 0.6385 0.2070 -0.1764 -31.55 41.73 -24.20 68.03];

disp(‘True Parameters’)
disp(true_theta)
%% ===============================
% Parameter Extraction
%% ===============================

sigma1 = xext(1);
sigma2 = xext(2);

delta1 = xext(3);
delta2 = xext(4);

l1 = xext(5);
l2 = xext(6);
l3 = xext(7);
l4 = xext(8);
%% ===============================
% Nonlinear Block
%% ===============================

khat = l1*g + l2*g.^2 + l3*g.^3 + l4*g.^4;
%% ===============================
% Estimated Output
%% ===============================

xhat = zeros(M,1);

xhat(1:5) = x(1:5);

for t = 6:M

xhat(t) = …
sigma1*(xhat(t-1)-0.3*xhat(t-2)-0.105*xhat(t-3)-0.0595*xhat(t-4)) …
+ sigma2*(xhat(t-2)-0.3*xhat(t-3)-0.105*xhat(t-4)-0.0595*xhat(t-5)) …
+ delta1*(khat(t-1)-0.3*khat(t-2)-0.105*khat(t-3)-0.0595*khat(t-4)) …
+ delta2*(khat(t-2)-0.3*khat(t-3)-0.105*khat(t-4)-0.0595*khat(t-5));

end

%% ===============================
% Performance
%% ===============================
disp(‘Best Cost from PSO’)
disp(global_best.cost)

err = x – xhat;

valid = isfinite(err) & (abs(err) < 1e6);

if any(valid)
mse_value = mean(err(valid).^2);
else
mse_value = Inf;
end

fprintf(‘Final MSE = %en’, mse_value);

%% ===============================
% Actual vs Estimated Output
%% ===============================

figure

plot(x,’b’,’LineWidth’,1.5)
hold on
plot(xhat,’r–‘,’LineWidth’,1.5)

xlabel(‘Samples’)
ylabel(‘Output’)
title(‘Actual vs Estimated Output’)

legend(‘Actual Output’,’Estimated Output’)
grid on

function j = fitness_function(xext,x,g,M)

%% ==========================================
% Column vectors
%% ==========================================

x = x(:);
g = g(:);

%% ==========================================
% Parameters
%% ==========================================

sigma1 = xext(1);
sigma2 = xext(2);

delta1 = xext(3);
delta2 = xext(4);

l1 = xext(5);
l2 = xext(6);
l3 = xext(7);
l4 = xext(8);

%% ==========================================
% Nonlinear Block
%% ==========================================

khat = l1*g + l2*g.^2 + l3*g.^3 + l4*g.^4;

%% ==========================================
% Free Run Simulation
%% ==========================================

xhat = zeros(M,1);

% Initial conditions
xhat(1:5) = x(1:5);

for t = 6:M

xhat(t) = …
sigma1*(xhat(t-1) …
-0.3*xhat(t-2) …
-0.105*xhat(t-3) …
-0.0595*xhat(t-4)) …
+ sigma2*(xhat(t-2) …
-0.3*xhat(t-3) …
-0.105*xhat(t-4) …
-0.0595*xhat(t-5)) …
+ delta1*(khat(t-1) …
-0.3*khat(t-2) …
-0.105*khat(t-3) …
-0.0595*khat(t-4)) …
+ delta2*(khat(t-2) …
-0.3*khat(t-3) …
-0.105*khat(t-4) …
-0.0595*khat(t-5));

end

%% ==========================================
% Mean Squared Error
%% ==========================================

err = x – xhat;

% Keep only reasonable errors
valid = isfinite(err) & (abs(err) < 1e6);

if sum(valid) < 0.80*M

j = 1e20;

else

err = err(valid);

j = mean(err.^2);

end

if isnan(j) || isinf(j)

j = 1e20;

end

end clc
clear
close all

%% ===============================
% READ HEAT EXCHANGER DATASET
%% ===============================

data = readmatrix(‘exchanger.xlsx’);

% Column 1 = Input
g = data(:,1);

% Column 2 = Output
x = data(:,2);

% Number of estimation samples
M = 3863;

% Select first M samples
g = g(1:M);
x = x(1:M);

% Convert to column vectors
g = g(:);
x = x(:);

%% ===============================
% PSO PARAMETERS
%% ===============================

dim = 8;

nParticles = 100;
maxIter = 800;

lb = [-2.0 0.3 0.0 -0.4 -40 30 -30 55];
ub = [-1.0 0.9 0.4 -0.1 -25 50 -20 75];

%% ===============================
% INITIALIZATION
%% ===============================

particle = struct;

for i = 1:nParticles

particle(i).position = lb + rand(1,dim).*(ub-lb);

particle(i).velocity = zeros(1,dim);

particle(i).cost = fitness_function( …
particle(i).position,x,g,M);

particle(i).best.position = particle(i).position;

particle(i).best.cost = particle(i).cost;

end

costs = [particle.cost];

[~,idx] = min(costs);

global_best = particle(idx).best;

%% ===============================
% PSO PARAMETERS
%% ===============================

w = 0.6;
c1 = 2.5;
c2 = 2.5;

BestCost = zeros(maxIter,1);
theta_history = zeros(maxIter,dim);
mse_history = zeros(maxIter,1);

%% ===============================
% PSO MAIN LOOP
%% ===============================

for iter = 1:maxIter

for i = 1:nParticles

particle(i).velocity = …
w*particle(i).velocity …
+ c1*rand*(particle(i).best.position-particle(i).position) …
+ c2*rand*(global_best.position-particle(i).position);

particle(i).position = …
particle(i).position + particle(i).velocity;

particle(i).position = …
max(particle(i).position,lb);

particle(i).position = …
min(particle(i).position,ub);

particle(i).cost = fitness_function( …
particle(i).position,x,g,M);

if particle(i).cost < particle(i).best.cost

particle(i).best.position = particle(i).position;

particle(i).best.cost = particle(i).cost;

end

if particle(i).best.cost < global_best.cost

global_best = particle(i).best;

end

end

theta_history(iter,:) = global_best.position;

mse_history(iter) = global_best.cost;

BestCost(iter) = global_best.cost;

fprintf(‘Iteration %d Best Error = %en’, …
iter,global_best.cost);

end
%% ===============================
% Estimated Parameters
%% ===============================

xext = global_best.position;

disp(‘Estimated Parameters’)
disp(xext)

true_theta = [-1.6080 0.6385 0.2070 -0.1764 -31.55 41.73 -24.20 68.03];

disp(‘True Parameters’)
disp(true_theta)
%% ===============================
% Parameter Extraction
%% ===============================

sigma1 = xext(1);
sigma2 = xext(2);

delta1 = xext(3);
delta2 = xext(4);

l1 = xext(5);
l2 = xext(6);
l3 = xext(7);
l4 = xext(8);
%% ===============================
% Nonlinear Block
%% ===============================

khat = l1*g + l2*g.^2 + l3*g.^3 + l4*g.^4;
%% ===============================
% Estimated Output
%% ===============================

xhat = zeros(M,1);

xhat(1:5) = x(1:5);

for t = 6:M

xhat(t) = …
sigma1*(xhat(t-1)-0.3*xhat(t-2)-0.105*xhat(t-3)-0.0595*xhat(t-4)) …
+ sigma2*(xhat(t-2)-0.3*xhat(t-3)-0.105*xhat(t-4)-0.0595*xhat(t-5)) …
+ delta1*(khat(t-1)-0.3*khat(t-2)-0.105*khat(t-3)-0.0595*khat(t-4)) …
+ delta2*(khat(t-2)-0.3*khat(t-3)-0.105*khat(t-4)-0.0595*khat(t-5));

end

%% ===============================
% Performance
%% ===============================
disp(‘Best Cost from PSO’)
disp(global_best.cost)

err = x – xhat;

valid = isfinite(err) & (abs(err) < 1e6);

if any(valid)
mse_value = mean(err(valid).^2);
else
mse_value = Inf;
end

fprintf(‘Final MSE = %en’, mse_value);

%% ===============================
% Actual vs Estimated Output
%% ===============================

figure

plot(x,’b’,’LineWidth’,1.5)
hold on
plot(xhat,’r–‘,’LineWidth’,1.5)

xlabel(‘Samples’)
ylabel(‘Output’)
title(‘Actual vs Estimated Output’)

legend(‘Actual Output’,’Estimated Output’)
grid on

function j = fitness_function(xext,x,g,M)

%% ==========================================
% Column vectors
%% ==========================================

x = x(:);
g = g(:);

%% ==========================================
% Parameters
%% ==========================================

sigma1 = xext(1);
sigma2 = xext(2);

delta1 = xext(3);
delta2 = xext(4);

l1 = xext(5);
l2 = xext(6);
l3 = xext(7);
l4 = xext(8);

%% ==========================================
% Nonlinear Block
%% ==========================================

khat = l1*g + l2*g.^2 + l3*g.^3 + l4*g.^4;

%% ==========================================
% Free Run Simulation
%% ==========================================

xhat = zeros(M,1);

% Initial conditions
xhat(1:5) = x(1:5);

for t = 6:M

xhat(t) = …
sigma1*(xhat(t-1) …
-0.3*xhat(t-2) …
-0.105*xhat(t-3) …
-0.0595*xhat(t-4)) …
+ sigma2*(xhat(t-2) …
-0.3*xhat(t-3) …
-0.105*xhat(t-4) …
-0.0595*xhat(t-5)) …
+ delta1*(khat(t-1) …
-0.3*khat(t-2) …
-0.105*khat(t-3) …
-0.0595*khat(t-4)) …
+ delta2*(khat(t-2) …
-0.3*khat(t-3) …
-0.105*khat(t-4) …
-0.0595*khat(t-5));

end

%% ==========================================
% Mean Squared Error
%% ==========================================

err = x – xhat;

% Keep only reasonable errors
valid = isfinite(err) & (abs(err) < 1e6);

if sum(valid) < 0.80*M

j = 1e20;

else

err = err(valid);

j = mean(err.^2);

end

if isnan(j) || isinf(j)

j = 1e20;

end

end is my code correct can anyone help me MATLAB Answers — New Questions

​

Leave a Reply

Your email address will not be published. Required fields are marked *