Livescript results in incorrect Math output

2 次查看(过去 30 天)
When using a LiveScript for symbolic math, I noticed a presumed bug in which the script outputs the incorrect value for an evaluation. When running the same equation(s) in the command window, the correct value is calculated. I tried restarting MATLAB but no change. It is concerning to me that the LiveScript is somehow providing erroneous calculations. What causes this?
Offending line of code is at the bottom; the rest is just in case you want to run it as a script yourself.
clear all
sympref('AbbreviateOutput',false); sympref('MatrixWithSquareBrackets',true);
syms x1 x2
f = -x1*x2 %objective function
g1 = (1/9)*x1^2 + (1/4)*x2^2 - 1 %given inequality constraint
g2 = -x1 %bounds set as inequality constraints
g3 = x1 - 3
g4 = -x2
g5 = x2 - 2
lambda_j = 1 %initialize inequality lagrange parameter
rp = 1 %initialize penalty parameter
b = 20 %initialize penalty adjustment parameter 'beta'
k = 1; %first ALM iteration
x = [2.5; 1.6] %initial design point
psi_j = max(subs(g1,[x1;x2],x), -lambda_j/(2*rp))
if subs(g1,[x1;x2],x) > 0 %if inequality is violated
A = f + lambda_j*psi_j + rp*psi_j^2
else %if no constraints are violated
A = f
end
c = subs(gradient(A),[x1;x2],x) %gradient of objective function at current point
d = -c
% Analytically solve for astar
syms alpha_star
xnew = x + alpha_star*d %update design point with alpha_star
% Solve for alpha using necessary condition df/da = 0
fa = subs(A,[x1;x2],xnew) %create alpha function
df_da = diff(fa)
astar = solve(df_da,alpha_star);
displaySymSolution(astar);
s = astar*d %calclulate the change in design
xnew = x + astar*d %update design point numerically
x = xnew;
%lambda_j + 2*rp*max(subs(g1,[x1;x2],x), -lambda_j/(2*rp)) some issue with
%the live script here makes the lambda_j update incorrect. Must manually
%update.
lambda_j = 1.6689
Lambda_j should == 1.6689; however, the LiveScript outputs '0' while the command window will output '1.6689' for the same code.
  2 个评论
Cris LaPierre
Cris LaPierre 2021-11-8
编辑:Cris LaPierre 2021-11-8
What code are you running in the command window vs the live script? I get the same result running the following code in both places:
clear all
sympref('AbbreviateOutput',false); sympref('MatrixWithSquareBrackets',true);
syms x1 x2
f = -x1*x2 %objective function
f = 
g1 = (1/9)*x1^2 + (1/4)*x2^2 - 1 %given inequality constraint
g1 = 
g2 = -x1 %bounds set as inequality constraints
g2 = 
g3 = x1 - 3
g3 = 
g4 = -x2
g4 = 
g5 = x2 - 2
g5 = 
lambda_j = 1 %initialize inequality lagrange parameter
lambda_j = 1
rp = 1 %initialize penalty parameter
rp = 1
b = 20 %initialize penalty adjustment parameter 'beta'
b = 20
k = 1; %first ALM iteration
x = [2.5; 1.6] %initial design point
x = 2×1
2.5000 1.6000
psi_j = max(subs(g1,[x1;x2],x), -lambda_j/(2*rp))
psi_j = 
if subs(g1,[x1;x2],x) > 0 %if inequality is violated
A = f + lambda_j*psi_j + rp*psi_j^2
else %if no constraints are violated
A = f
end
A = 
c = subs(gradient(A),[x1;x2],x) %gradient of objective function at current point
c = 
d = -c
d = 
% Analytically solve for astar
syms alpha_star
xnew = x + alpha_star*d %update design point with alpha_star
xnew = 
% Solve for alpha using necessary condition df/da = 0
fa = subs(A,[x1;x2],xnew) %create alpha function
fa = 
df_da = diff(fa)
df_da = 
astar = solve(df_da,alpha_star);
displaySymSolution(astar);
astar = 
s = astar*d %calclulate the change in design
s = 
xnew = x + astar*d %update design point numerically
xnew = 
x = xnew;
lambda_j + 2*rp*max(subs(g1,[x1;x2],x), -lambda_j/(2*rp))
ans = 
0
Colton Campbell
Colton Campbell 2021-11-8
Ah yes I see... I was not updating the x-value when executing in the command window. Thank you!

请先登录,再进行评论。

回答(1 个)

Chris
Chris 2021-11-8
编辑:Chris 2021-11-8
The commented-out line is:
lambda_j + 2*rp*max(subs(g1,[x1;x2],x), -lambda_j/(2*rp))
lambda_j is 1. The second term (which uses the second argument to max(), containing lambda_j) is -1. The sum of the two is zero.
  1 个评论
Colton Campbell
Colton Campbell 2021-11-8
Ah yes I see... I was not updating the x-value when executing in the command window. Thank you!

请先登录,再进行评论。

类别

Help CenterFile Exchange 中查找有关 Symbolic Math Toolbox 的更多信息

产品


版本

R2021b

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by