MATLAB Answers

How to implement optimal control using Pontryagin's Minimum Principle for HESS in EV application as in attached paper?

87 views (last 30 days)
Praveen Kumar
Praveen Kumar on 18 Jan 2019
Commented: Praveen Kumar on 28 Jan 2019
I am trying to implement a paper based on optimal control of energy management of HESS in EV. It's based on Pontryagin's Minimum Principle using Hamiltonian, state and costate equations. How to simulate this? Is it better to use matlab m file or simulink?
Please find the paper attached.

  0 Comments

Sign in to comment.

Answers (1)

Stephan
Stephan on 18 Jan 2019
Edited: Stephan on 18 Jan 2019
Hi,
there are some examples here on answers, where people used Matlab with Symbolic Math Toolbox to solve this kind of problems. I do not remember Simulink modesl in this context - but maybe there are some. A newer example using Matlab with a working code you can find here:
Best regards
Stephan

  2 Comments

Praveen Kumar
Praveen Kumar on 28 Jan 2019
clear;
%Initializations and values of constant variables
C = 28;
Rs = 0.01;
Rp = 10e3;
Vsmax = 305;
Preq = 50e3;
k = 1;
Vb = 360;
Rb = 0.3;
p1 = -3;
tf = 5;
% State equations
syms x1 Pb p1; %Pb control variable & x1 = SOC_s is state variable
Dx1 = -1/C*((1/(2*Rs)+1/Rp)*x1-1/(2*Rs)*sqrt(x1^2-4*Rs/Vsmax^2*(Preq-Pb)));
% Cost function inside integral
syms g;
g = k*((Vb - sqrt(Vb^2 - 4*Rb*Pb))/(2*Rb))^2;
% Hamiltonian
syms p1 H;
H = g + p1*Dx1;
% Costate equations
% Dx1 = diff(H,p1);
Dp1 = -diff(H,x1);
% Solve for control Pb
dPb = diff(H,Pb);
sol_Pb = solve(dPb,Pb);
% Substitute Pb to state equations
Dx1 = subs(Dx1,Pb,sol_Pb);
% Convert symbolic objects to strings for using 'dsolve'
eq1 = strcat('Dx1=',char(Dx1));
eq2 = strcat('Dp1=',char(Dp1));
sol_h = dsolve(eq1,eq2);
% Use boundary conditions to determine the coefficients
% Case a: (a)Dx1(0)=90 Dx1(tf)=50
conA1 = 'x1(0) = 0';
conA2 = 'x1(10) = 1';
sol_a = dsolve(eq1,eq2,conA1,conA2);
% plot solutions
figure(1)
ezplot(sol_a.x1,[0 1]); hold on;
ezplot(sol_a.p1,[0 1]); % plot the control u=-p1
axis([0 20 0 1]);
xlabel('time');
ylabel('states');
title('Solution of PMP')
I have tried to implement my example using the framework given in the above link. I am getting an error for dsolve.
Error using mupadengine/feval (line 187)
Invalid number of arguments.
Error in dsolve>mupadDsolve (line 340)
T = feval(symengine,'symobj::dsolve',sys,x,options);
Error in dsolve (line 194)
sol = mupadDsolve(args, options);
Error in PMP_18012019 (line 33)
sol_h = dsolve(eq1,eq2);.

Sign in to comment.

Sign in to answer this question.