How to do optimization for result calculated by using experimental data?
Show older comments
Hi All,
I have run a measurement for 28 fabricated samples based on few geometrical parameter (x(1), x(2), x(3) and x(4)). Each sample has a different combination of geometrical parameters. From the measurement result, i got list of 'alpha' and 'xi' (that varied with frequency) for each sample. Based on the 'alpha' and 'xi' values, I calculate the output, y (versus frequency) for each sample. With the data I have, I would like to find the optimized parameter for x(1), x(2), x(3) and x(4) when the output is maximized to 1 (y=1) at frequency 1600Hz.
I tried a few method but unfortunately it fails to achieve my goal.
I really appreaciate if someone can help me to solve this problem.
I have attached the files that I have for your kind reference. Kindly, please advice me the next step to achieve my goal.
Thank you.
Regards,
Nur Arafah
alpha = load ('alpha.txt'); %% alpha
xi = load ('xi.txt'); %% xi
c = 345.2;
n = 1.83e-5;
den = 1.189;
f = 500:1.5625:4000; %% frequency
%% for samples with different x(1)
x(1)=700e-6;
for i = 1:length(f)
w(i) = 2.*pi.*f(i);
a(i) = alpha(i,1);
E(i) = xi(i,1);
x(2) = 1000e-6
x(3) = 0.00615
x(4) = 0.0042
k(i) = x(1).*sqrt((w(i).*den)./(4.*n));
Rs(i) = (1./2).*sqrt(2.*n.*den.*w(i));
r1(i) = (((32.*n.*x(2))./(x(3).*den.*c.*(x(1).^2))).*sqrt(1+((k(i).^2)./32)))+((2.*a(i).*Rs(i))./(x(3).*den.*c));
m1(i) = ((((w(i).*x(2))./(x(3).*c)).*(1+(1./sqrt(9+((k(i).^2)./2)))))+(((w(i).^2).*(x(1).^2))./(8.*x(3).*(c.^2)))+((8.*x(1).*w(i))./(3.*pi.*x(3).*c)).*E(i))-cot((w(i).*x(4))./c);
z(i) = r1(i)+(j.*m1(i));
Ref(i) = (z(i)-1)./(z(i)+1);
Ref_re(i) = real(Ref(i));
Ref_im(i) = imag(Ref(i));
RF(i) = sqrt((Ref_re(i).^2)+(Ref_im(i).^2));
y(i) = 1 - abs(RF(i)).^2; %% output for sample #1
end
y1=y;
x(1)=750e-6;
for i = 1:length(f)
w(i) = 2.*pi.*f(i);
a(i) = alpha(i,2);
E(i) = xi(i,2);
x(2) = 1000e-6
x(3) = 0.00615
x(4) = 0.0042
k(i) = x(1).*sqrt((w(i).*den)./(4.*n));
Rs(i) = (1./2).*sqrt(2.*n.*den.*w(i));
r1(i) = (((32.*n.*x(2))./(x(3).*den.*c.*(x(1).^2))).*sqrt(1+((k(i).^2)./32)))+((2.*a(i).*Rs(i))./(x(3).*den.*c));
m1(i) = ((((w(i).*x(2))./(x(3).*c)).*(1+(1./sqrt(9+((k(i).^2)./2)))))+(((w(i).^2).*(x(1).^2))./(8.*x(3).*(c.^2)))+((8.*x(1).*w(i))./(3.*pi.*x(3).*c)).*E(i))-cot((w(i).*x(4))./c);
z(i) = r1(i)+(j.*m1(i));
Ref(i) = (z(i)-1)./(z(i)+1);
Ref_re(i) = real(Ref(i));
Ref_im(i) = imag(Ref(i));
RF(i) = sqrt((Ref_re(i).^2)+(Ref_im(i).^2));
y(i) = 1 - abs(RF(i)).^2; %% output for sample #2
end
y2=y;
x(1)=800e-6;
for i = 1:length(f)
w(i) = 2.*pi.*f(i);
a(i) = alpha(i,3);
E(i) = xi(i,3);
x(2) = 1000e-6
x(3) = 0.00615
x(4) = 0.0042
k(i) = x(1).*sqrt((w(i).*den)./(4.*n));
Rs(i) = (1./2).*sqrt(2.*n.*den.*w(i));
r1(i) = (((32.*n.*x(2))./(x(3).*den.*c.*(x(1).^2))).*sqrt(1+((k(i).^2)./32)))+((2.*a(i).*Rs(i))./(x(3).*den.*c));
m1(i) = ((((w(i).*x(2))./(x(3).*c)).*(1+(1./sqrt(9+((k(i).^2)./2)))))+(((w(i).^2).*(x(1).^2))./(8.*x(3).*(c.^2)))+((8.*x(1).*w(i))./(3.*pi.*x(3).*c)).*E(i))-cot((w(i).*x(4))./c);
z(i) = r1(i)+(j.*m1(i));
Ref(i) = (z(i)-1)./(z(i)+1);
Ref_re(i) = real(Ref(i));
Ref_im(i) = imag(Ref(i));
RF(i) = sqrt((Ref_re(i).^2)+(Ref_im(i).^2));
y(i) = 1 - abs(RF(i)).^2; %% output for sample #3
end
y3=y;
x(1)=850e-6;
for i = 1:length(f)
w(i) = 2.*pi.*f(i);
a(i) = alpha(i,4);
E(i) = xi(i,4);
x(2) = 1000e-6
x(3) = 0.00615
x(4) = 0.0042
k(i) = x(1).*sqrt((w(i).*den)./(4.*n));
Rs(i) = (1./2).*sqrt(2.*n.*den.*w(i));
r1(i) = (((32.*n.*x(2))./(x(3).*den.*c.*(x(1).^2))).*sqrt(1+((k(i).^2)./32)))+((2.*a(i).*Rs(i))./(x(3).*den.*c));
m1(i) = ((((w(i).*x(2))./(x(3).*c)).*(1+(1./sqrt(9+((k(i).^2)./2)))))+(((w(i).^2).*(x(1).^2))./(8.*x(3).*(c.^2)))+((8.*x(1).*w(i))./(3.*pi.*x(3).*c)).*E(i))-cot((w(i).*x(4))./c);
z(i) = r1(i)+(j.*m1(i));
Ref(i) = (z(i)-1)./(z(i)+1);
Ref_re(i) = real(Ref(i));
Ref_im(i) = imag(Ref(i));
RF(i) = sqrt((Ref_re(i).^2)+(Ref_im(i).^2));
y(i) = 1 - abs(RF(i)).^2; %% output for sample #4
end
y4=y;
x(1)=900e-6;
for i = 1:length(f)
w(i) = 2.*pi.*f(i);
a(i) = alpha(i,5);
E(i) = xi(i,5);
x(2) = 1000e-6
x(3) = 0.00615
x(4) = 0.0042
k(i) = x(1).*sqrt((w(i).*den)./(4.*n));
Rs(i) = (1./2).*sqrt(2.*n.*den.*w(i));
r1(i) = (((32.*n.*x(2))./(x(3).*den.*c.*(x(1).^2))).*sqrt(1+((k(i).^2)./32)))+((2.*a(i).*Rs(i))./(x(3).*den.*c));
m1(i) = ((((w(i).*x(2))./(x(3).*c)).*(1+(1./sqrt(9+((k(i).^2)./2)))))+(((w(i).^2).*(x(1).^2))./(8.*x(3).*(c.^2)))+((8.*x(1).*w(i))./(3.*pi.*x(3).*c)).*E(i))-cot((w(i).*x(4))./c);
z(i) = r1(i)+(j.*m1(i));
Ref(i) = (z(i)-1)./(z(i)+1);
Ref_re(i) = real(Ref(i));
Ref_im(i) = imag(Ref(i));
RF(i) = sqrt((Ref_re(i).^2)+(Ref_im(i).^2));
y(i) = 1 - abs(RF(i)).^2; %% output for sample #5
end
y5=y;
x(1)=950e-6;
for i = 1:length(f)
w(i) = 2.*pi.*f(i);
a(i) = alpha(i,6);
E(i) = xi(i,6);
x(2) = 1000e-6
x(3) = 0.00615
x(4) = 0.0042
k(i) = x(1).*sqrt((w(i).*den)./(4.*n));
Rs(i) = (1./2).*sqrt(2.*n.*den.*w(i));
r1(i) = (((32.*n.*x(2))./(x(3).*den.*c.*(x(1).^2))).*sqrt(1+((k(i).^2)./32)))+((2.*a(i).*Rs(i))./(x(3).*den.*c));
m1(i) = ((((w(i).*x(2))./(x(3).*c)).*(1+(1./sqrt(9+((k(i).^2)./2)))))+(((w(i).^2).*(x(1).^2))./(8.*x(3).*(c.^2)))+((8.*x(1).*w(i))./(3.*pi.*x(3).*c)).*E(i))-cot((w(i).*x(4))./c);
z(i) = r1(i)+(j.*m1(i));
Ref(i) = (z(i)-1)./(z(i)+1);
Ref_re(i) = real(Ref(i));
Ref_im(i) = imag(Ref(i));
RF(i) = sqrt((Ref_re(i).^2)+(Ref_im(i).^2));
y(i) = 1 - abs(RF(i)).^2; %% output for sample #6
end
y6=y;
x(1)=1000e-6;
for i = 1:length(f)
w(i) = 2.*pi.*f(i);
a(i) = alpha(i,7);
E(i) = xi(i,7);
x(2) = 1000e-6
x(3) = 0.00615
x(4) = 0.0042
k(i) = x(1).*sqrt((w(i).*den)./(4.*n));
Rs(i) = (1./2).*sqrt(2.*n.*den.*w(i));
r1(i) = (((32.*n.*x(2))./(x(3).*den.*c.*(x(1).^2))).*sqrt(1+((k(i).^2)./32)))+((2.*a(i).*Rs(i))./(x(3).*den.*c));
m1(i) = ((((w(i).*x(2))./(x(3).*c)).*(1+(1./sqrt(9+((k(i).^2)./2)))))+(((w(i).^2).*(x(1).^2))./(8.*x(3).*(c.^2)))+((8.*x(1).*w(i))./(3.*pi.*x(3).*c)).*E(i))-cot((w(i).*x(4))./c);
z(i) = r1(i)+(j.*m1(i));
Ref(i) = (z(i)-1)./(z(i)+1);
Ref_re(i) = real(Ref(i));
Ref_im(i) = imag(Ref(i));
RF(i) = sqrt((Ref_re(i).^2)+(Ref_im(i).^2));
y(i) = 1 - abs(RF(i)).^2; %% output for sample #7
end
y7=y;
%% for samples with different x(2)
x(2)=700e-6;
for i = 1:length(f)
w(i) = 2.*pi.*f(i);
a(i) = alpha(i,8);
E(i) = xi(i,8);
x(1) = 700e-6
x(3) = 0.00615
x(4) = 0.0042
k(i) = x(1).*sqrt((w(i).*den)./(4.*n));
Rs(i) = (1./2).*sqrt(2.*n.*den.*w(i));
r1(i) = (((32.*n.*x(2))./(x(3).*den.*c.*(x(1).^2))).*sqrt(1+((k(i).^2)./32)))+((2.*a(i).*Rs(i))./(x(3).*den.*c));
m1(i) = ((((w(i).*x(2))./(x(3).*c)).*(1+(1./sqrt(9+((k(i).^2)./2)))))+(((w(i).^2).*(x(1).^2))./(8.*x(3).*(c.^2)))+((8.*x(1).*w(i))./(3.*pi.*x(3).*c)).*E(i))-cot((w(i).*x(4))./c);
z(i) = r1(i)+(j.*m1(i));
Ref(i) = (z(i)-1)./(z(i)+1);
Ref_re(i) = real(Ref(i));
Ref_im(i) = imag(Ref(i));
RF(i) = sqrt((Ref_re(i).^2)+(Ref_im(i).^2));
y(i) = 1 - abs(RF(i)).^2; %% output for sample #8
end
y8=y;
x(2)=800e-6;
for i = 1:length(f)
w(i) = 2.*pi.*f(i);
a(i) = alpha(i,9);
E(i) = xi(i,9);
x(1) = 700e-6
x(3) = 0.00615
x(4) = 0.0042
k(i) = x(1).*sqrt((w(i).*den)./(4.*n));
Rs(i) = (1./2).*sqrt(2.*n.*den.*w(i));
r1(i) = (((32.*n.*x(2))./(x(3).*den.*c.*(x(1).^2))).*sqrt(1+((k(i).^2)./32)))+((2.*a(i).*Rs(i))./(x(3).*den.*c));
m1(i) = ((((w(i).*x(2))./(x(3).*c)).*(1+(1./sqrt(9+((k(i).^2)./2)))))+(((w(i).^2).*(x(1).^2))./(8.*x(3).*(c.^2)))+((8.*x(1).*w(i))./(3.*pi.*x(3).*c)).*E(i))-cot((w(i).*x(4))./c);
z(i) = r1(i)+(j.*m1(i));
Ref(i) = (z(i)-1)./(z(i)+1);
Ref_re(i) = real(Ref(i));
Ref_im(i) = imag(Ref(i));
RF(i) = sqrt((Ref_re(i).^2)+(Ref_im(i).^2));
y(i) = 1 - abs(RF(i)).^2; %% output for sample #9
end
y9=y;
x(2)=900e-6;
for i = 1:length(f)
w(i) = 2.*pi.*f(i);
a(i) = alpha(i,10);
E(i) = xi(i,10);
x(1) = 700e-6
x(3) = 0.00615
x(4) = 0.0042
k(i) = x(1).*sqrt((w(i).*den)./(4.*n));
Rs(i) = (1./2).*sqrt(2.*n.*den.*w(i));
r1(i) = (((32.*n.*x(2))./(x(3).*den.*c.*(x(1).^2))).*sqrt(1+((k(i).^2)./32)))+((2.*a(i).*Rs(i))./(x(3).*den.*c));
m1(i) = ((((w(i).*x(2))./(x(3).*c)).*(1+(1./sqrt(9+((k(i).^2)./2)))))+(((w(i).^2).*(x(1).^2))./(8.*x(3).*(c.^2)))+((8.*x(1).*w(i))./(3.*pi.*x(3).*c)).*E(i))-cot((w(i).*x(4))./c);
z(i) = r1(i)+(j.*m1(i));
Ref(i) = (z(i)-1)./(z(i)+1);
Ref_re(i) = real(Ref(i));
Ref_im(i) = imag(Ref(i));
RF(i) = sqrt((Ref_re(i).^2)+(Ref_im(i).^2));
y(i) = 1 - abs(RF(i)).^2; %% output for sample #10
end
y10=y;
x(2)=1000e-6;
for i = 1:length(f)
w(i) = 2.*pi.*f(i);
a(i) = alpha(i,11);
E(i) = xi(i,11);
x(1) = 700e-6
x(3) = 0.00615
x(4) = 0.0042
k(i) = x(1).*sqrt((w(i).*den)./(4.*n));
Rs(i) = (1./2).*sqrt(2.*n.*den.*w(i));
r1(i) = (((32.*n.*x(2))./(x(3).*den.*c.*(x(1).^2))).*sqrt(1+((k(i).^2)./32)))+((2.*a(i).*Rs(i))./(x(3).*den.*c));
m1(i) = ((((w(i).*x(2))./(x(3).*c)).*(1+(1./sqrt(9+((k(i).^2)./2)))))+(((w(i).^2).*(x(1).^2))./(8.*x(3).*(c.^2)))+((8.*x(1).*w(i))./(3.*pi.*x(3).*c)).*E(i))-cot((w(i).*x(4))./c);
z(i) = r1(i)+(j.*m1(i));
Ref(i) = (z(i)-1)./(z(i)+1);
Ref_re(i) = real(Ref(i));
Ref_im(i) = imag(Ref(i));
RF(i) = sqrt((Ref_re(i).^2)+(Ref_im(i).^2));
y(i) = 1 - abs(RF(i)).^2; %% output for sample #10
end
y11=y;
x(2)=1100e-6;
for i = 1:length(f)
w(i) = 2.*pi.*f(i);
a(i) = alpha(i,12);
E(i) = xi(i,12);
x(1) = 700e-6
x(3) = 0.00615
x(4) = 0.0042
k(i) = x(1).*sqrt((w(i).*den)./(4.*n));
Rs(i) = (1./2).*sqrt(2.*n.*den.*w(i));
r1(i) = (((32.*n.*x(2))./(x(3).*den.*c.*(x(1).^2))).*sqrt(1+((k(i).^2)./32)))+((2.*a(i).*Rs(i))./(x(3).*den.*c));
m1(i) = ((((w(i).*x(2))./(x(3).*c)).*(1+(1./sqrt(9+((k(i).^2)./2)))))+(((w(i).^2).*(x(1).^2))./(8.*x(3).*(c.^2)))+((8.*x(1).*w(i))./(3.*pi.*x(3).*c)).*E(i))-cot((w(i).*x(4))./c);
z(i) = r1(i)+(j.*m1(i));
Ref(i) = (z(i)-1)./(z(i)+1);
Ref_re(i) = real(Ref(i));
Ref_im(i) = imag(Ref(i));
RF(i) = sqrt((Ref_re(i).^2)+(Ref_im(i).^2));
y(i) = 1 - abs(RF(i)).^2; %% output for sample #12
end
y12=y;
x(2)=1200e-6;
for i = 1:length(f)
w(i) = 2.*pi.*f(i);
a(i) = alpha(i,13);
E(i) = xi(i,13);
x(1) = 700e-6
x(3) = 0.00615
x(4) = 0.0042
k(i) = x(1).*sqrt((w(i).*den)./(4.*n));
Rs(i) = (1./2).*sqrt(2.*n.*den.*w(i));
r1(i) = (((32.*n.*x(2))./(x(3).*den.*c.*(x(1).^2))).*sqrt(1+((k(i).^2)./32)))+((2.*a(i).*Rs(i))./(x(3).*den.*c));
m1(i) = ((((w(i).*x(2))./(x(3).*c)).*(1+(1./sqrt(9+((k(i).^2)./2)))))+(((w(i).^2).*(x(1).^2))./(8.*x(3).*(c.^2)))+((8.*x(1).*w(i))./(3.*pi.*x(3).*c)).*E(i))-cot((w(i).*x(4))./c);
z(i) = r1(i)+(j.*m1(i));
Ref(i) = (z(i)-1)./(z(i)+1);
Ref_re(i) = real(Ref(i));
Ref_im(i) = imag(Ref(i));
RF(i) = sqrt((Ref_re(i).^2)+(Ref_im(i).^2));
y(i) = 1 - abs(RF(i)).^2; %% output for sample #13
end
y13=y;
x(2)=1300e-6;
for i = 1:length(f)
w(i) = 2.*pi.*f(i);
a(i) = alpha(i,14);
E(i) = xi(i,14);
x(1) = 700e-6
x(3) = 0.00615
x(4) = 0.0042
k(i) = x(1).*sqrt((w(i).*den)./(4.*n));
Rs(i) = (1./2).*sqrt(2.*n.*den.*w(i));
r1(i) = (((32.*n.*x(2))./(x(3).*den.*c.*(x(1).^2))).*sqrt(1+((k(i).^2)./32)))+((2.*a(i).*Rs(i))./(x(3).*den.*c));
m1(i) = ((((w(i).*x(2))./(x(3).*c)).*(1+(1./sqrt(9+((k(i).^2)./2)))))+(((w(i).^2).*(x(1).^2))./(8.*x(3).*(c.^2)))+((8.*x(1).*w(i))./(3.*pi.*x(3).*c)).*E(i))-cot((w(i).*x(4))./c);
z(i) = r1(i)+(j.*m1(i));
Ref(i) = (z(i)-1)./(z(i)+1);
Ref_re(i) = real(Ref(i));
Ref_im(i) = imag(Ref(i));
RF(i) = sqrt((Ref_re(i).^2)+(Ref_im(i).^2));
y(i) = 1 - abs(RF(i)).^2; %% output for sample #14
end
y14=y;
%% for samples with different x(3)
x(3)=0.004;
for i = 1:length(f)
w(i) = 2.*pi.*f(i);
a(i) = alpha(i,15);
E(i) = xi(i,15);
x(1) = 700e-6
x(2) = 1000e-6
x(4) = 0.0042
k(i) = x(1).*sqrt((w(i).*den)./(4.*n));
Rs(i) = (1./2).*sqrt(2.*n.*den.*w(i));
r1(i) = (((32.*n.*x(2))./(x(3).*den.*c.*(x(1).^2))).*sqrt(1+((k(i).^2)./32)))+((2.*a(i).*Rs(i))./(x(3).*den.*c));
m1(i) = ((((w(i).*x(2))./(x(3).*c)).*(1+(1./sqrt(9+((k(i).^2)./2)))))+(((w(i).^2).*(x(1).^2))./(8.*x(3).*(c.^2)))+((8.*x(1).*w(i))./(3.*pi.*x(3).*c)).*E(i))-cot((w(i).*x(4))./c);
z(i) = r1(i)+(j.*m1(i));
Ref(i) = (z(i)-1)./(z(i)+1);
Ref_re(i) = real(Ref(i));
Ref_im(i) = imag(Ref(i));
RF(i) = sqrt((Ref_re(i).^2)+(Ref_im(i).^2));
y(i) = 1 - abs(RF(i)).^2; %% output for sample #15
end
y15=y;
x(3)=0.005;
for i = 1:length(f)
w(i) = 2.*pi.*f(i);
a(i) = alpha(i,16);
E(i) = xi(i,16);
x(1) = 700e-6
x(2) = 1000e-6
x(4) = 0.0042
k(i) = x(1).*sqrt((w(i).*den)./(4.*n));
Rs(i) = (1./2).*sqrt(2.*n.*den.*w(i));
r1(i) = (((32.*n.*x(2))./(x(3).*den.*c.*(x(1).^2))).*sqrt(1+((k(i).^2)./32)))+((2.*a(i).*Rs(i))./(x(3).*den.*c));
m1(i) = ((((w(i).*x(2))./(x(3).*c)).*(1+(1./sqrt(9+((k(i).^2)./2)))))+(((w(i).^2).*(x(1).^2))./(8.*x(3).*(c.^2)))+((8.*x(1).*w(i))./(3.*pi.*x(3).*c)).*E(i))-cot((w(i).*x(4))./c);
z(i) = r1(i)+(j.*m1(i));
Ref(i) = (z(i)-1)./(z(i)+1);
Ref_re(i) = real(Ref(i));
Ref_im(i) = imag(Ref(i));
RF(i) = sqrt((Ref_re(i).^2)+(Ref_im(i).^2));
y(i) = 1 - abs(RF(i)).^2; %% output for sample #16
end
y16=y;
x(3)=0.006;
for i = 1:length(f)
w(i) = 2.*pi.*f(i);
a(i) = alpha(i,17);
E(i) = xi(i,17);
x(1) = 700e-6
x(2) = 1000e-6
x(4) = 0.0042
k(i) = x(1).*sqrt((w(i).*den)./(4.*n));
Rs(i) = (1./2).*sqrt(2.*n.*den.*w(i));
r1(i) = (((32.*n.*x(2))./(x(3).*den.*c.*(x(1).^2))).*sqrt(1+((k(i).^2)./32)))+((2.*a(i).*Rs(i))./(x(3).*den.*c));
m1(i) = ((((w(i).*x(2))./(x(3).*c)).*(1+(1./sqrt(9+((k(i).^2)./2)))))+(((w(i).^2).*(x(1).^2))./(8.*x(3).*(c.^2)))+((8.*x(1).*w(i))./(3.*pi.*x(3).*c)).*E(i))-cot((w(i).*x(4))./c);
z(i) = r1(i)+(j.*m1(i));
Ref(i) = (z(i)-1)./(z(i)+1);
Ref_re(i) = real(Ref(i));
Ref_im(i) = imag(Ref(i));
RF(i) = sqrt((Ref_re(i).^2)+(Ref_im(i).^2));
y(i) = 1 - abs(RF(i)).^2; %% output for sample #17
end
y17=y;
x(3)=0.007;
for i = 1:length(f)
w(i) = 2.*pi.*f(i);
a(i) = alpha(i,18);
E(i) = xi(i,18);
x(1) = 700e-6
x(2) = 1000e-6
x(4) = 0.0042
k(i) = x(1).*sqrt((w(i).*den)./(4.*n));
Rs(i) = (1./2).*sqrt(2.*n.*den.*w(i));
r1(i) = (((32.*n.*x(2))./(x(3).*den.*c.*(x(1).^2))).*sqrt(1+((k(i).^2)./32)))+((2.*a(i).*Rs(i))./(x(3).*den.*c));
m1(i) = ((((w(i).*x(2))./(x(3).*c)).*(1+(1./sqrt(9+((k(i).^2)./2)))))+(((w(i).^2).*(x(1).^2))./(8.*x(3).*(c.^2)))+((8.*x(1).*w(i))./(3.*pi.*x(3).*c)).*E(i))-cot((w(i).*x(4))./c);
z(i) = r1(i)+(j.*m1(i));
Ref(i) = (z(i)-1)./(z(i)+1);
Ref_re(i) = real(Ref(i));
Ref_im(i) = imag(Ref(i));
RF(i) = sqrt((Ref_re(i).^2)+(Ref_im(i).^2));
y(i) = 1 - abs(RF(i)).^2; %% output for sample #18
end
y18=y;
x(3)=0.008;
for i = 1:length(f)
w(i) = 2.*pi.*f(i);
a(i) = alpha(i,19);
E(i) = xi(i,19);
x(1) = 700e-6
x(2) = 1000e-6
x(4) = 0.0042
k(i) = x(1).*sqrt((w(i).*den)./(4.*n));
Rs(i) = (1./2).*sqrt(2.*n.*den.*w(i));
r1(i) = (((32.*n.*x(2))./(x(3).*den.*c.*(x(1).^2))).*sqrt(1+((k(i).^2)./32)))+((2.*a(i).*Rs(i))./(x(3).*den.*c));
m1(i) = ((((w(i).*x(2))./(x(3).*c)).*(1+(1./sqrt(9+((k(i).^2)./2)))))+(((w(i).^2).*(x(1).^2))./(8.*x(3).*(c.^2)))+((8.*x(1).*w(i))./(3.*pi.*x(3).*c)).*E(i))-cot((w(i).*x(4))./c);
z(i) = r1(i)+(j.*m1(i));
Ref(i) = (z(i)-1)./(z(i)+1);
Ref_re(i) = real(Ref(i));
Ref_im(i) = imag(Ref(i));
RF(i) = sqrt((Ref_re(i).^2)+(Ref_im(i).^2));
y(i) = 1 - abs(RF(i)).^2; %% output for sample #19
end
y19=y;
x(3)=0.009;
for i = 1:length(f)
w(i) = 2.*pi.*f(i);
a(i) = alpha(i,20);
E(i) = xi(i,20);
x(1) = 700e-6
x(2) = 1000e-6
x(4) = 0.0042
k(i) = x(1).*sqrt((w(i).*den)./(4.*n));
Rs(i) = (1./2).*sqrt(2.*n.*den.*w(i));
r1(i) = (((32.*n.*x(2))./(x(3).*den.*c.*(x(1).^2))).*sqrt(1+((k(i).^2)./32)))+((2.*a(i).*Rs(i))./(x(3).*den.*c));
m1(i) = ((((w(i).*x(2))./(x(3).*c)).*(1+(1./sqrt(9+((k(i).^2)./2)))))+(((w(i).^2).*(x(1).^2))./(8.*x(3).*(c.^2)))+((8.*x(1).*w(i))./(3.*pi.*x(3).*c)).*E(i))-cot((w(i).*x(4))./c);
z(i) = r1(i)+(j.*m1(i));
Ref(i) = (z(i)-1)./(z(i)+1);
Ref_re(i) = real(Ref(i));
Ref_im(i) = imag(Ref(i));
RF(i) = sqrt((Ref_re(i).^2)+(Ref_im(i).^2));
y(i) = 1 - abs(RF(i)).^2; %% output for sample #20
end
y20=y;
x(3)=0.01;
for i = 1:length(f)
w(i) = 2.*pi.*f(i);
a(i) = alpha(i,21);
E(i) = xi(i,21);
x(1) = 700e-6
x(2) = 1000e-6
x(4) = 0.0042
k(i) = x(1).*sqrt((w(i).*den)./(4.*n));
Rs(i) = (1./2).*sqrt(2.*n.*den.*w(i));
r1(i) = (((32.*n.*x(2))./(x(3).*den.*c.*(x(1).^2))).*sqrt(1+((k(i).^2)./32)))+((2.*a(i).*Rs(i))./(x(3).*den.*c));
m1(i) = ((((w(i).*x(2))./(x(3).*c)).*(1+(1./sqrt(9+((k(i).^2)./2)))))+(((w(i).^2).*(x(1).^2))./(8.*x(3).*(c.^2)))+((8.*x(1).*w(i))./(3.*pi.*x(3).*c)).*E(i))-cot((w(i).*x(4))./c);
z(i) = r1(i)+(j.*m1(i));
Ref(i) = (z(i)-1)./(z(i)+1);
Ref_re(i) = real(Ref(i));
Ref_im(i) = imag(Ref(i));
RF(i) = sqrt((Ref_re(i).^2)+(Ref_im(i).^2));
y(i) = 1 - abs(RF(i)).^2; %% output for sample #21
end
y21=y;
%% for samples with different x(4)
x(4)=0.001;
for i = 1:length(f)
w(i) = 2.*pi.*f(i);
a(i) = alpha(i,22);
E(i) = xi(i,22);
x(1) = 700e-6
x(2) = 1000e-6
x(3) = 0.00615
k(i) = x(1).*sqrt((w(i).*den)./(4.*n));
Rs(i) = (1./2).*sqrt(2.*n.*den.*w(i));
r1(i) = (((32.*n.*x(2))./(x(3).*den.*c.*(x(1).^2))).*sqrt(1+((k(i).^2)./32)))+((2.*a(i).*Rs(i))./(x(3).*den.*c));
m1(i) = ((((w(i).*x(2))./(x(3).*c)).*(1+(1./sqrt(9+((k(i).^2)./2)))))+(((w(i).^2).*(x(1).^2))./(8.*x(3).*(c.^2)))+((8.*x(1).*w(i))./(3.*pi.*x(3).*c)).*E(i))-cot((w(i).*x(4))./c);
z(i) = r1(i)+(j.*m1(i));
Ref(i) = (z(i)-1)./(z(i)+1);
Ref_re(i) = real(Ref(i));
Ref_im(i) = imag(Ref(i));
RF(i) = sqrt((Ref_re(i).^2)+(Ref_im(i).^2));
y(i) = 1 - abs(RF(i)).^2; %% output for sample #22
end
y22=y;
x(4)=0.002;
for i = 1:length(f)
w(i) = 2.*pi.*f(i);
a(i) = alpha(i,23);
E(i) = xi(i,23);
x(1) = 700e-6
x(2) = 1000e-6
x(3) = 0.00615
k(i) = x(1).*sqrt((w(i).*den)./(4.*n));
Rs(i) = (1./2).*sqrt(2.*n.*den.*w(i));
r1(i) = (((32.*n.*x(2))./(x(3).*den.*c.*(x(1).^2))).*sqrt(1+((k(i).^2)./32)))+((2.*a(i).*Rs(i))./(x(3).*den.*c));
m1(i) = ((((w(i).*x(2))./(x(3).*c)).*(1+(1./sqrt(9+((k(i).^2)./2)))))+(((w(i).^2).*(x(1).^2))./(8.*x(3).*(c.^2)))+((8.*x(1).*w(i))./(3.*pi.*x(3).*c)).*E(i))-cot((w(i).*x(4))./c);
z(i) = r1(i)+(j.*m1(i));
Ref(i) = (z(i)-1)./(z(i)+1);
Ref_re(i) = real(Ref(i));
Ref_im(i) = imag(Ref(i));
RF(i) = sqrt((Ref_re(i).^2)+(Ref_im(i).^2));
y(i) = 1 - abs(RF(i)).^2; %% output for sample #23
end
y23=y;
x(4)=0.003;
for i = 1:length(f)
w(i) = 2.*pi.*f(i);
a(i) = alpha(i,24);
E(i) = xi(i,24);
x(1) = 700e-6
x(2) = 1000e-6
x(3) = 0.00615
k(i) = x(1).*sqrt((w(i).*den)./(4.*n));
Rs(i) = (1./2).*sqrt(2.*n.*den.*w(i));
r1(i) = (((32.*n.*x(2))./(x(3).*den.*c.*(x(1).^2))).*sqrt(1+((k(i).^2)./32)))+((2.*a(i).*Rs(i))./(x(3).*den.*c));
m1(i) = ((((w(i).*x(2))./(x(3).*c)).*(1+(1./sqrt(9+((k(i).^2)./2)))))+(((w(i).^2).*(x(1).^2))./(8.*x(3).*(c.^2)))+((8.*x(1).*w(i))./(3.*pi.*x(3).*c)).*E(i))-cot((w(i).*x(4))./c);
z(i) = r1(i)+(j.*m1(i));
Ref(i) = (z(i)-1)./(z(i)+1);
Ref_re(i) = real(Ref(i));
Ref_im(i) = imag(Ref(i));
RF(i) = sqrt((Ref_re(i).^2)+(Ref_im(i).^2));
y(i) = 1 - abs(RF(i)).^2; %% output for sample #24
end
y24=y;
x(4)=0.004;
for i = 1:length(f)
w(i) = 2.*pi.*f(i);
a(i) = alpha(i,25);
E(i) = xi(i,25);
x(1) = 700e-6
x(2) = 1000e-6
x(3) = 0.00615
k(i) = x(1).*sqrt((w(i).*den)./(4.*n));
Rs(i) = (1./2).*sqrt(2.*n.*den.*w(i));
r1(i) = (((32.*n.*x(2))./(x(3).*den.*c.*(x(1).^2))).*sqrt(1+((k(i).^2)./32)))+((2.*a(i).*Rs(i))./(x(3).*den.*c));
m1(i) = ((((w(i).*x(2))./(x(3).*c)).*(1+(1./sqrt(9+((k(i).^2)./2)))))+(((w(i).^2).*(x(1).^2))./(8.*x(3).*(c.^2)))+((8.*x(1).*w(i))./(3.*pi.*x(3).*c)).*E(i))-cot((w(i).*x(4))./c);
z(i) = r1(i)+(j.*m1(i));
Ref(i) = (z(i)-1)./(z(i)+1);
Ref_re(i) = real(Ref(i));
Ref_im(i) = imag(Ref(i));
RF(i) = sqrt((Ref_re(i).^2)+(Ref_im(i).^2));
y(i) = 1 - abs(RF(i)).^2; %% output for sample #25
end
y25=y;
x(4)=0.005;
for i = 1:length(f)
w(i) = 2.*pi.*f(i);
a(i) = alpha(i,26);
E(i) = xi(i,26);
x(1) = 700e-6
x(2) = 1000e-6
x(3) = 0.00615
k(i) = x(1).*sqrt((w(i).*den)./(4.*n));
Rs(i) = (1./2).*sqrt(2.*n.*den.*w(i));
r1(i) = (((32.*n.*x(2))./(x(3).*den.*c.*(x(1).^2))).*sqrt(1+((k(i).^2)./32)))+((2.*a(i).*Rs(i))./(x(3).*den.*c));
m1(i) = ((((w(i).*x(2))./(x(3).*c)).*(1+(1./sqrt(9+((k(i).^2)./2)))))+(((w(i).^2).*(x(1).^2))./(8.*x(3).*(c.^2)))+((8.*x(1).*w(i))./(3.*pi.*x(3).*c)).*E(i))-cot((w(i).*x(4))./c);
z(i) = r1(i)+(j.*m1(i));
Ref(i) = (z(i)-1)./(z(i)+1);
Ref_re(i) = real(Ref(i));
Ref_im(i) = imag(Ref(i));
RF(i) = sqrt((Ref_re(i).^2)+(Ref_im(i).^2));
y(i) = 1 - abs(RF(i)).^2; %% output for sample #26
end
y26=y;
x(4)=0.006;
for i = 1:length(f)
w(i) = 2.*pi.*f(i);
a(i) = alpha(i,27);
E(i) = xi(i,27);
x(1) = 700e-6
x(2) = 1000e-6
x(3) = 0.00615
k(i) = x(1).*sqrt((w(i).*den)./(4.*n));
Rs(i) = (1./2).*sqrt(2.*n.*den.*w(i));
r1(i) = (((32.*n.*x(2))./(x(3).*den.*c.*(x(1).^2))).*sqrt(1+((k(i).^2)./32)))+((2.*a(i).*Rs(i))./(x(3).*den.*c));
m1(i) = ((((w(i).*x(2))./(x(3).*c)).*(1+(1./sqrt(9+((k(i).^2)./2)))))+(((w(i).^2).*(x(1).^2))./(8.*x(3).*(c.^2)))+((8.*x(1).*w(i))./(3.*pi.*x(3).*c)).*E(i))-cot((w(i).*x(4))./c);
z(i) = r1(i)+(j.*m1(i));
Ref(i) = (z(i)-1)./(z(i)+1);
Ref_re(i) = real(Ref(i));
Ref_im(i) = imag(Ref(i));
RF(i) = sqrt((Ref_re(i).^2)+(Ref_im(i).^2));
y(i) = 1 - abs(RF(i)).^2; %% output for sample #27
end
y27=y;
x(4)=0.007;
for i = 1:length(f)
w(i) = 2.*pi.*f(i);
a(i) = alpha(i,28);
E(i) = xi(i,28);
x(1) = 700e-6
x(2) = 1000e-6
x(3) = 0.00615
k(i) = x(1).*sqrt((w(i).*den)./(4.*n));
Rs(i) = (1./2).*sqrt(2.*n.*den.*w(i));
r1(i) = (((32.*n.*x(2))./(x(3).*den.*c.*(x(1).^2))).*sqrt(1+((k(i).^2)./32)))+((2.*a(i).*Rs(i))./(x(3).*den.*c));
m1(i) = ((((w(i).*x(2))./(x(3).*c)).*(1+(1./sqrt(9+((k(i).^2)./2)))))+(((w(i).^2).*(x(1).^2))./(8.*x(3).*(c.^2)))+((8.*x(1).*w(i))./(3.*pi.*x(3).*c)).*E(i))-cot((w(i).*x(4))./c);
z(i) = r1(i)+(j.*m1(i));
Ref(i) = (z(i)-1)./(z(i)+1);
Ref_re(i) = real(Ref(i));
Ref_im(i) = imag(Ref(i));
RF(i) = sqrt((Ref_re(i).^2)+(Ref_im(i).^2));
y(i) = 1 - abs(RF(i)).^2; %% output for sample #28
end
y28=y;
%% all the output result from sample #1 to #28
f=[f;f;f;f;f;f;f;f;f;f;f;f;f;f;f;f;f;f;f;f;f;f;f;f;f;f;f;f]
y=[y1;y2;y3;y4;y5;y6;y7;y8;y9;y10;y11;y12;y13;y14;y15;y16;y17;y18;y19;y20;y21;y22;y23;y24;y25;y26;y27;y28]
plot(f,y);
10 Comments
Mario Malic
on 19 Jul 2020
If I understood your code correctly, you have vector for each sample that contains output values y for the frequencies in range [500, 4000].
You just need to find the value of i_freq that corresponds to frequency of 1600 Hz and find the maximum value of y(i_freq) for all the vector samples.
Since the code is quite repetitive, the only varibles that change are x(1) to x(4), write a function that you call for every sample
% You send a vector f and x variables, and receive y as a output vector and x as a vector of values
[y, x] = outputcalc(f, x(1), x(2), x(3), x(4))
Then you can save it in a array that each sample represents a row of a matrix. Then your problem would be just finding max value of a specific column in y, and with the index of it, you call the particular row of x array.
Jafar Nur Arafah
on 19 Jul 2020
Mario Malic
on 19 Jul 2020
Edited: Mario Malic
on 19 Jul 2020
- Please see the attached files. There's no need for this simple script to be this long.
- I hope the code represents what you want to do, in the other words, check if you get the same outputs
- For your and our sanity, work with arrays.
x = [700e-6, 1000e-6, 0.00615, 0.0042];
y = Sample_Calculation (x, f, alpha, xi, c, n, den, Sample_num);
EDIT 1: (Question is: what is your output? y is a (ixi) matrix and you mentioned how you want to maximize y on 1600 Hz. Most definitely is that you have expected a y as a vector so it means that you didn't write your code well.) I didn't clear variables, error on my side. Will update
EDIT 2: The second question is, do you need to use the optimization algorithm? Are you looking for minimum of y for given x(1), x(2), x(3), x(4) variables? Didn't read well either.
Mario Malic
on 19 Jul 2020
Edited: Mario Malic
on 19 Jul 2020
Here should be something good for you, at least to move further with your code.
I assume x variables need to be positively defined, so you need to include that in your constraints.
Can't upload more files, as I don't know how to delete the ones which are uploaded before, that's why it's written like this.
clc;
clear vars;
alpha = load ('alpha.txt');
xi = load ('xi.txt');
c = 345.2;
n = 1.83e-5;
den = 1.189;
f = 500:1.5625:4000;
Sample_Num = 1; % This doesn't change, unnecessary probably
Desired_Freq = find (f == 1600);
%% Calculating objective value for each sample
options = optimoptions('ga','TolFun', 1e-3, 'PlotFcn', @gaplotbestf);
x = [x(1), x(2), x(3), x(4)];
[X_Opt, fval, exitflag,output] = ga(@(x)Sample_Calculation(x, f, alpha, xi, c, n, den, Sample_Num, Desired_Freq), 4, options)
% y = Sample_Calculation (X_Opt, f, alpha, xi, c, n, den, Sample_Num, Desired_Freq)
Sample_Calculation function
function y = Sample_Calculation (x, f, alpha, xi, c, n, den, Sample_Num, Desired_Freq)
w = 2.*pi.*f;
a = alpha(:,Sample_Num)';
E = xi(:,Sample_Num)';
k = x(1).*sqrt((w.*den)./(4.*n));
Rs = (1./2).*sqrt(2.*n.*den.*w);
r1 = (((32.*n.*x(2))./(x(3).*den.*c.*(x(1).^2))).*sqrt(1+((k.^2)./32)))+((2.*a.*Rs)./(x(3).*den.*c));
m1 = ((((w.*x(2))./(x(3).*c)).*(1+(1./sqrt(9+((k.^2)./2)))))+(((w.^2).*(x(1).^2))./(8.*x(3).*(c.^2)))+((8.*x(1).*w)./(3.*pi.*x(3).*c)).*E)-cot((w.*x(4))./c);
z = r1+(j.*m1);
Ref = (z-1)./(z+1);
Ref_re = real(Ref);
Ref_im = imag(Ref);
RF = sqrt((Ref_re.^2)+(Ref_im.^2));
yn= 1 - abs(RF).^2; %% output for sample #1
y = 1 + yn(1,Desired_Freq); % Objective function, please CHECK this
% plot (f,yn); % only displays last iteration
end
Jafar Nur Arafah
on 19 Jul 2020
Mario Malic
on 19 Jul 2020
Edited: Mario Malic
on 19 Jul 2020
Here's code which is not working because y is not equal to the same variable from your code. I can't spend time troubleshooting where did things go wrong so, you have to figure it out yourself.
Also, reading your post history, I doubt this is what you actually wanted. Since this only gives you a maximum of your experiments, not the actual maximum.
format long
clc;
clearvars;
alpha = (load ('alpha.txt'))';
xi = (load ('xi.txt'))';
c = 345.2;
n = 1.83e-5;
den = 1.189;
f = 500:1.5625:4000;
Sample_Num = 28;
Desired_Freq = find (f == 1600);
% Defining X matrix
x = zeros(28,4);
x(1:7, 1) = [700e-6:50e-6:1000e-6];
x(8:28,1) = 700e-6;
x(1:7,2) = 1000e-6;
x(8:14, 2) = [700e-6:100e-6:1300e-6];
x(15:28,2) = 1000e-6;
x(1:14, 3) = 0.00615;
x(15:21, 3) = [0.004:0.001:0.010];
x(22:28,3) = 0.00615;
x(1:21, 4) = 0.0042;
x(22:28, 4) = [0.001:0.001:0.007];
y = zeros (Sample_Num, length(f));
for ii = 1 : 1 : Sample_Num
y(ii,:) = Sample_Calculation (x(ii,:), ii, f, alpha, xi, c, n, den);
end
[Max_Y, Index] = max(y(:,Desired_Freq))
X_Opt = x(Index,:)
Function
function y = Sample_Calculation (x, ii, f, alpha, xi, c, n, den)
w = 2.*pi.*f;
a = alpha(ii,:);
E = xi(ii,:);
k = x(1).*sqrt((w.*den)./(4.*n));
Rs = (1./2).*sqrt(2.*n.*den.*w);
r1 = (((32.*n.*x(2))./(x(3).*den.*c.*(x(1).^2))).*sqrt(1+((k.^2)./32)))+((2.*a.*Rs)./(x(3).*den.*c));
m1 = ((((w.*x(2))./(x(3).*c)).*(1+(1./sqrt(9+((k.^2)./2)))))+(((w.^2).*(x(1).^2))./(8.*x(3).*(c.^2)))+((8.*x(1).*w)./(3.*pi.*x(3).*c)).*E)-cot((w.*x(4))./c);
z = r1+(j.*m1);
Ref = (z-1)./(z+1);
Ref_re = real(Ref);
Ref_im = imag(Ref);
RF = sqrt((Ref_re.^2)+(Ref_im.^2));
y = 1 - abs(RF).^2; %% output for sample #1
end
Jafar Nur Arafah
on 20 Jul 2020
Mario Malic
on 20 Jul 2020
Edited: Mario Malic
on 20 Jul 2020
I am not so good (if at all) at optimization stuff. I would suggest rather spend some time reading about optimization functions and think how your problem may relate to these. Maybe work with simpler problems from the optimization documentation and make your way up.
It kind of feels that with every message the problem expands in some sort of way or changes the objective, so try to be clear what do you want to do and explain it in the best way possible.
"I got error "Your fitness function must return a scalar value.". The problem is i need to use all the value (scalar value) for optimization purpose.", for GA, by reffering to documentation you can see that objective function must return a scalar value. About the other optimization algorithms that accept vectors see here LINK.
Also, I would recommend you to write the code for optimization rather than use it through toolbox. You'll understand more while doing so.
Here's a sample code.
ii = 1;
options = optimoptions('ga','TolFun', 1e-3, 'PlotFcn', @gaplotbestf);
lb = [700e-6,700e-6,0.004,0.001];
ub = [1000e-6,1300e-6,0.010,0.007];
[X_Opt, fval, exitflag,output] = ga(@(x)Sample_Calculation(x, ii, f, alpha, xi, c, n, den), 4,[], [], [], [], lb, ub, [], options)
And if you change the output in function to return the scalar value like this or in any other way the code will work.
yn = 1 - abs(RF).^2; %% output for sample #1
y = 1 + yn(1,705); % obj fun
Sorry that I can not proceed further with the help, since I am not able to understand what you actually want to do.
Jafar Nur Arafah
on 20 Jul 2020
Mario Malic
on 20 Jul 2020
You're welcome. Good luck and have a great day as well!
Answers (0)
Categories
Find more on Multiobjective Optimization in Help Center and File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!