function r = nonlinconstr(prob)
global globstat
prob.cnl = [];
prob.cnl.lb = prob.lb;
prob.cnl.ub = prob.ub;
if ~isempty(prob.b) || ~isempty(prob.beq) || ...
(~isempty(prob.lb) && any(isfinite(prob.lb))) || ...
(~isempty(prob.ub) && any(isfinite(prob.ub)))
prob = lin2nonlin(prob);
end
%setting Lagrange multipliers - AFTER lin2nonlin()!
lambda = zeros(prob.nh,1);%beware! - zeros can yield unbounded solution
mu = zeros(prob.nc,1);
maxUpdateMultipliers = 100;
ExtractRoutineParamsFromProb;
if ~isfield(prob, 'lambda'); prob.lambda = lambda; end
if ~isfield(prob, 'mu'); prob.mu = mu; end
FUNCS = {'f' 'd2f' 'c' 'h' 'dc' 'dh' 'd2c' 'd2h'};
% NO df!! - see df handling below!!
prob.cnlf = prob.f;
prob.cnldf = prob.df;
if ~isempty(prob.df) && isempty(prob.primal.df)
prob.cnl.df = @(x, prob, varargin) prob.cnldf(x, prob, [], prob.cnlf);
prob.df = @(x, prob, varargin) cnl_df(x, prob, varargin, prob.cnlf);
else
FUNCS{end+1} = 'df';
end
% for fn = {'df' 'f' 'd2f' 'c' 'h' 'dc' 'dh' 'd2c' 'd2h'} % NO df!! - see df handling above!!
for fn = FUNCS
if ~isfield(prob, fn{1}) || isempty(getfield(prob, fn{1})); continue; end %#ok<GFLD>
F = getfield(prob, fn{1}); %#ok<GFLD>
prob.cnl = setfield(prob.cnl, fn{1}, F); %#ok<SFLD>
eval(['prob.' fn{1} '= @cnl_' fn{1} ';']); %#ok<SFLD>
%debug
if ~ismember(fn, {'f' 'c' 'h' 'df' 'dc' 'dh'}) && ~isempty(getfield(prob.primal, fn{1})) %#ok<GFLD>
prob.err(['nonlinconstr wrapper: ' fn{1} ' isn''t implemented yet'])
end
%debug end
end
%todo: remove (make auto before 1st call)
prob = updateLM(prob.x0, prob);
prob.cnl.iterfcn = prob.iterfcn;
prob.iterfcn = @nonLinConstrIterFcn;
% prob.cnl.iterfcn = prob.iterfcn;
% prob.iterfcn = @LMiter;
prob.nc = 0;
prob.nh = 0;
% todo: add prob.cnl.nc, prob.cnl.nh
nUpdateMultipliers = 0;
prob.f0 = prob.f(prob.x0, prob);
% prob.advanced.ignoreIStop(end+1)=5;
nn=10;
if isequal(lower(prob.solverName), 'hpso')
prob.MaxIter = prob.MaxIter/nn;
end
for itn=1:maxUpdateMultipliers
r = prob.solver(prob);
% if prob.debug; disp(['istop: ' num2str(r.istop)]); end
% disp(r.xf)
if r.istop == 13; return; end
if ~stopcase(r) % MaxTime, MaxCPUTime, MaxIter, MaxFunEvals or etc is exeeded
if ~isequal(lower(prob.solverName), 'hpso') || itn==nn
break
end
end
if ~isempty(prob.c)%todo: if r includes r.c info
c = prob.cnl.c(r.xf, prob);
else c = [];
end
if ~isempty(prob.h)
h = prob.cnl.h(r.xf, prob);
else h = [];
end
bViolationH = abs(h)>prob.TolCon;%matrix
bViolationC = c>prob.TolCon;%matrix
indViolationH = find(bViolationH);
indViolationC = find(bViolationC);%#ok<EFIND> %matrix
if isempty(indViolationH) && isempty(indViolationC) && r.istop ~= -8.8888
r.advanced.nUpdateMultipliers = nUpdateMultipliers;
% todo: hint for setting start lambda, mu
break
end
nUpdateMultipliers = nUpdateMultipliers + 1;
if nUpdateMultipliers>maxUpdateMultipliers
r.istop = -13;
return
end
% if prob.advanced.debug
if prob.iterPrint>0 && ~isequal(lower(prob.solverName), 'hpso')
prob.info(['nIter: ' int2str(globstat.nIter) ' nUpdateLagrangeMultipliers = ' int2str(nUpdateMultipliers)])
end
if prob.debug
disp(['istop: ' num2str(r.istop)])
if r.istop == -5; disp('No line-search minimum is found'); end
s = '';
if ~isempty(c); s = horzcat(s, ' ', int2str(length(indViolationC)), ' ineq'); end
if ~isempty(h); s = horzcat(s, ' ', int2str(length(indViolationH)), ' eq'); end
prob.info(['constraints violation: ' s])
end
% [prob.lambda prob.mu] = updateLM(prob);
% if ~isempty(bViolationH)
% j = find(bViolationH & ~prob.lambda);
% prob.lambda(j) = prob.lambda(j) + 1;
% end
% if ~isempty(bViolationC)
% j = find(bViolationC & ~prob.mu);
% prob.mu(j) = prob.mu(j) + 1;
% end
% if r.istop ~= -5
prob.x0 = r.xf;
% end
prob = updateLM(prob.x0, prob);
% globstat.xPrev=[];%checkme!
if r.istop ~= -5%no line-search found
prob.f0 = prob.f(prob.x0, prob);
end
if ismember(prob.solverName, {'ralg' 'hPSO'}) && isfield(r, 'advanced') && isfield(r.advanced, 'ralg') && isfield(r.advanced.ralg, 'hs')
% prob.ralg.b = r.advanced.ralg.b;
prob.ralg.hs = r.advanced.ralg.hs;
end
% prob.f0 = r.ff;
globstat.isFinished = 0;
if isfield(globstat, 'lastX')
globstat.lastX = inf*ones(size(globstat.lastX));
end
globstat.ff = inf;
globstat.bestOptimFVal = inf;
if isfield(globstat, 'lastOptimFVal')
globstat.lastOptimFVal = inf;
end
end
if itn==maxUpdateMultipliers && (~(isempty(indViolationH) && isempty(indViolationC)) || r.istop == -8.8888)
r.istop = -13;
end
if ~isempty(prob.lambda); r.lambda = prob.lambda; end
if ~isempty(prob.mu); r.mu = prob.mu; end
end
% function prob = updateLM(x, prob)
% %linear constraints must be alreay converted to non-lin
% % f = prob.f(x, prob);
% % ndf = norm(prob.df(x, prob),1);
% DF = prob.cnldf(x, prob);
% % ndf = norm(DF,1);
%
% if isfield(prob.cnl, 'dc')
% c = prob.cnl.c(x, prob);
% % Nc = norm(c, 1)+realmin;
% dc = prob.cnl.dc(x, prob);
% else
% dc = [];
% end
% if isfield(prob.cnl, 'dh')
% h = prob.cnl.h(x, prob);
% % Nh = norm(h, 1)+realmin;
% dh = prob.cnl.dh(x, prob);
% else
% dh = [];
% end
%
% p = [c; abs(h)];
% p(p<prob.TolCon) = 0;
% if ~any(p); return; end
% dh2 = dh .* repmat(sign(h), 1, size(dh,1)).';
% A = [dc dh2];
% b = DF;
% % A = sum(A);
% % b = sum(b);
% x0 = [prob.mu; prob.lambda];
% [x lam exitflag output] = linprog(p,A,b,[],[],zeros(length(x0),1),[],[]);
% if exitflag<=0
% disp(exitflag)
% disp(output)
% keyboard
% end
% prob.mu = x(1:length(prob.mu));
% prob.lambda = x(1+length(prob.mu):end);
% end