Livescript results in incorrect Math output

4 vues (au cours des 30 derniers jours)
Colton Campbell
Colton Campbell le 8 Nov 2021
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 commentaires
Cris LaPierre
Cris LaPierre le 8 Nov 2021
Modifié(e) : Cris LaPierre le 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 le 8 Nov 2021
Ah yes I see... I was not updating the x-value when executing in the command window. Thank you!

Connectez-vous pour commenter.

Réponses (1)

Chris
Chris le 8 Nov 2021
Modifié(e) : Chris le 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 commentaire
Colton Campbell
Colton Campbell le 8 Nov 2021
Ah yes I see... I was not updating the x-value when executing in the command window. Thank you!

Connectez-vous pour commenter.

Produits


Version

R2021b

Community Treasure Hunt

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

Start Hunting!

Translated by