MATLAB Answers

0

Find the minimum of a function and use it to plot a graph.

Asked by JACINTA ONWUKA on 9 Jul 2019
Latest activity Answered by Robert U
on 10 Jul 2019
My mission is to find Mycp such that Jcp is minimum and then plot a graph with minimum of Mycp as a function of rho.
I just started using matlab so dont understand why the graph is not showing? Can someone please help me?
%close all
L=1;
T=100;
r=0.03;
I1=0.1;
p=0.05;
epsilon=0.3;
beta=0.1;
rho=1000;
Myrho=0:100:1000
MycpRho=zeros(numel(Myrho),1);
for j = 1:numel(Myrho)
MyrhoCurrent=Myrho(j);
Mycp = 0:1:100
n = zeros(numel(Mycp),1 );
n2 = zeros(numel(Mycp),1 );
n3 = zeros(numel(Mycp),1 );
Jcp = zeros(numel(Mycp),1 );
for i = 1:numel(Mycp )
MycpCurrent=Mycp(i);
delta = 1-MycpCurrent/100 ;
tau = (1/(beta*(L+delta*p)))*log((L*(I1+delta*p))/(delta*p*(L-I1 )));
t05 =(1/(beta*(L+delta*p)))*log((L*(0.05*L+delta*p))/(delta*p*(L-0.05*L )));
I2= @(t)(L*delta*p*(exp (beta*(L+delta*p)*t)-1)) ./ (L + delta*p* exp(beta*(L+delta*p)*t ));
I3= @(t)(L*(I1+delta*p)*exp((epsilon*beta)*(L+delta*p)*(t-tau))-...
delta*p*(L -I1))./(L-I1+(I1+delta*p)*exp(epsilon*beta*(L+delta*p)*(t-tau)));
fun = @(t,MycpCurrent) MycpCurrent*L*exp(-r*t);
fun2=@(t)rho*I2(t).*exp(-r*t);
fun3=@(t)rho*I3(t).*exp(-r*t);
n(i) = integral(@(t)fun(t,MycpCurrent),0,100, 'ArrayValued',1);
n2(i)= integral(fun2,t05,tau);
n3(i)= integral(fun3,tau,100);
JCp(i)= n(i)+n2(i)+n3(i);
end
MycpBest=Mycp(JCp==min(JCp));
JCpBest=min(JCp);
plot(Myrho,MycpBest)
hold on
scatter(MycpBest,JCpBest)
hold off
end

  0 Comments

Sign in to comment.

1 Answer

Answer by Robert U
on 10 Jul 2019

Hi JACINTA ONWUKA,
I do not know what you expect as output. The graph you are trying to plot is produced with every iteration. Since you have the hold off command in every loop iteration, every new call of plot() is overwriting the actual active axes. Depending on what you want, you can either include
  • the command figure that would produce a new figure window every loop iteration or
  • you delete the hold off command in order to continue drawing into your initially opened axis.
Remark:
The command plot(Myrho, MycpBest) is trying to plot the vector Myrho against a scalar value. That gives you an empty plot. You can plot the current Myrho values against the current MycpBest values as dots in the graph by changing the last lines in your shared code:
L=1;
T=100;
r=0.03;
I1=0.1;
p=0.05;
epsilon=0.3;
beta=0.1;
rho=1000;
Myrho=0:100:1000;
MycpRho=zeros(numel(Myrho),1);
Mycp = 0:1:100;
n = zeros( numel(Mycp),1 );
n2 = zeros( numel(Mycp),1 );
n3 = zeros( numel(Mycp),1 );
Jcp = zeros( numel(Mycp),1 );
for j = 1:numel(Myrho)
MyrhoCurrent=Myrho(j);
for i = 1:numel(Mycp )
MycpCurrent=Mycp(i);
delta = 1-MycpCurrent/100 ;
tau = (1/(beta*(L+delta*p)))*log((L*(I1+delta*p))/(delta*p*(L-I1 )));
t05 =(1/(beta*(L+delta*p)))*log((L*(0.05*L+delta*p))/(delta*p*(L-0.05*L )));
I2 = @(t)(L*delta*p*(exp (beta*(L+delta*p)*t)-1)) ./ (L + delta*p*exp(beta*(L+delta*p)*t ));
I3 = @(t)(L*(I1+delta*p)*exp((epsilon*beta)*(L+delta*p)*(t-tau))-...
delta*p*(L -I1))./(L-I1+(I1+delta*p)*exp(epsilon*beta*(L+delta*p)*(t-tau)));
fun = @(t,MycpCurrent) MycpCurrent*L*exp(-r*t);
fun2 = @(t)rho*I2(t).*exp(-r*t);
fun3 = @(t)rho*I3(t).*exp(-r*t);
n(i) = integral(@(t)fun(t,MycpCurrent),0,100, 'ArrayValued',1);
n2(i)= integral(fun2,t05,tau);
n3(i)= integral(fun3,tau,100);
JCp(i)= n(i)+n2(i)+n3(i);
end
MycpBest = Mycp( JCp==min(JCp) );
JCpBest = min(JCp);
plot(Myrho(j),MycpBest,'.blue','MarkerSize',20)
hold on
scatter(MycpBest,JCpBest)
end
Kind regards,
Robert

  0 Comments

Sign in to comment.