how can i avoid Nan in matlab expression and return 0
25 views (last 30 days)
Show older comments
The below three expression gives three different values,
and can someone advice me how can i avoid Nan values for 0 input and return 0 value for it
(sqrt((-220).^2)./-220+1)/2.*(10-5) + 5
(sqrt((220).^2)./220+1)/2.*(10-5) + 5
(sqrt((0).^2)./0+1)/2.*(10-5) + 5
>> (sqrt((-220).^2)./-220+1)/2.*(10-5) + 5
(sqrt((220).^2)./220+1)/2.*(10-5) + 5
(sqrt((0).^2)./0+1)/2.*(10-5) + 5
ans =
5
ans =
10
ans =
NaN
Answers (5)
Matt J
on 5 May 2021
Edited: Matt J
on 5 May 2021
Replace
sqrt(x.^2)./x
with
sign(x)
5 Comments
Walter Roberson
on 5 May 2021
Optimization variables cannot do logical multiplication. They can only use comparisons in the context of constraints.
Image Analyst
on 5 May 2021
Why can't you just simply assign it to a variable and check if that variable is nan, and if it is, assign it to zero?
result = (sqrt((0).^2)./0+1)/2.*(10-5) + 5;
if isnan(result)
result = 0;
end
1 Comment
Walter Roberson
on 6 May 2021
Because the context is Problem Based Optimization, which does not support if or isnan() and only permits comparisons as part of constraints.
Walter Roberson
on 6 May 2021
The following is not exactly right:
M = optimvar('M', 1, 5)
costa = 5
costb = 10
delta = eps(realmin)
part1a = (M - sqrt(M.^2))/2
part1b = part1a./(part1a - delta)
part2a = (M + sqrt(M.^2))/2
part2b = part2a./(part2a + delta)
cost = part1b * costa + part2b * costb
Mathematically it is wrong at exactly two points, M = -eps(realmin) and M = +eps(realmin) . For those two points, the output should be costa and costb respectively, but instead the formula mathematically gives costa/2 and costb/2 at those two points instead.
In practice, though, for values sufficiently close to +/- realmin, numeric evaluation might return costa+costb and for values sufficiently clost to eps(realmin) numeric evaluation might return 0 instead of costa or costb .
The exact result in a range close to +/- realmin is going to depend on the exact order of evaluation, which is not something that we have control over; optimization could potentially re-arrange the evaluation.
This code has been constructed so that it should never return NaN.
How important is it for your purposes that the values must be correct near +/- realmin, given that it is designed to return 0 for 0 exactly?
4 Comments
Walter Roberson
on 6 May 2021
Example. This needs to be upgraded to have Costb passed in as well, but your scripts and the Simulink model currently only pass in a single cost vector.
function [Pgrid,Pbatt,Ebatt] = battSolarOptimize(N,dt,Ppv,Pload,Einit,Cost,FinalWeight,batteryMinMax)
% Minimize the cost of power from the grid while meeting load with power
% from PV, battery and grid
prob = optimproblem;
% Decision variables
PgridV = optimvar('PgridV',N);
PbattV = optimvar('PbattV',N,'LowerBound',batteryMinMax.Pmin,'UpperBound',batteryMinMax.Pmax);
EbattV = optimvar('EbattV',N,'LowerBound',batteryMinMax.Emin,'UpperBound',batteryMinMax.Emax);
% Minimize cost of electricity from the grid
prob.ObjectiveSense = 'minimize';
%{
prob.Objective = dt*Cost'*PgridV - FinalWeight*EbattV(N);
%}
Costa = Cost(:,1);
Costb = Cost(:,end) + randn(size(Costa))/20;
r = 200;
%Cost = (PbattV<=r).*Costa+ (~PbattV>=r).*Costb;
Cost = fcn2optimexpr(@(PV) (PV < r) .* Costa + (PV > r) .* Costb, PbattV);
%P1 = dt*Cost'*PbattV;
%P2 = dt*Cost'*PbattV;
P = dt*Cost'*PbattV;
prob.Objective = dt*Cost'*PgridV - FinalWeight*EbattV(N);
% Power input/output to battery
prob.Constraints.energyBalance = optimconstr(N);
prob.Constraints.energyBalance(1) = EbattV(1) == Einit;
prob.Constraints.energyBalance(2:N) = EbattV(2:N) == EbattV(1:N-1) - PbattV(1:N-1)*dt;
% Satisfy power load with power from PV, grid and battery
prob.Constraints.loadBalance = Ppv + PgridV + PbattV == Pload;
x0.PgridV = ones(1,N);
x0.PbattV = batteryMinMax.Pmin * ones(1,N);
x0.EbattV = batteryMinMax.Emin * ones(1,N);
% Solve the linear program
options = optimoptions(prob.optimoptions, 'Display', 'iter', 'MaxFunctionEvaluations', 1e7, 'MaxIterations', 1024 );
[values,~,exitflag] = solve(prob, x0, 'Options', options);
% Parse optmization results
if exitflag <= 0
fprintf('warning: ran into iteration or function evaluation limit!\n');
end
Pgrid = values.PgridV;
Pbatt = values.PbattV;
Ebatt = values.EbattV;
Matt J
on 6 May 2021
Edited: Matt J
on 6 May 2021
Here's a way you can do it by adding some additional binary variables and linear constraints. It requires that x be bounded to the interval [-1,1]. You can easily introduce a variable z=A*x+B if you need a variable that spans a different range.
x=optimvar('x','LowerBound',-1,'UpperBound',+1);
b1=optimvar('b1','LowerBound',0,'UpperBound',1,'type','integer'); %Binary variables
b2=optimvar('b2','LowerBound',0,'UpperBound',1,'type','integer');
%%EDITED
con1(1)= b1>=x+eps(0); % b1 "ON" when x>=0
con2(1)= b1<=x+1; % b1 "OFF when x<0
con1(2)= b2>=x; % b2 "ON" when x>0
con2(2)= b2<=x+1-eps(1) % b2 "OFF when x<=0
y = 5*(1 + (b1+b2)/2 -3*(b1-b2)/2);
% x>0 ==> b1=b2=1 ==> y=10
% x<0 ==> b1=b2=0 ==> y=5
% x=0 ==> b1=1, b2=0 ==> y=0
3 Comments
Walter Roberson
on 7 May 2021
Each constraint is structed by constraint type (same for all elements in the constraint), and two pieces of data indicating what the constraint is operating on.
I do not know why they choose to implement that way..
See Also
Categories
Find more on Get Started with Optimization Toolbox in Help Center and File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!