Livescript results in incorrect Math output

2 views (last 30 days)
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 Comments
Cris LaPierre
Cris LaPierre on 8 Nov 2021
Edited: Cris LaPierre on 8 Nov 2021
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 on 8 Nov 2021
Ah yes I see... I was not updating the x-value when executing in the command window. Thank you!

Sign in to comment.

Answers (1)

Chris
Chris on 8 Nov 2021
Edited: Chris on 8 Nov 2021
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 Comment
Colton Campbell
Colton Campbell on 8 Nov 2021
Ah yes I see... I was not updating the x-value when executing in the command window. Thank you!

Sign in to comment.

Categories

Find more on Symbolic Math Toolbox in Help Center and File Exchange

Products


Release

R2021b

Community Treasure Hunt

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

Start Hunting!