MATLAB: Fmincon linear contraint doens’t work for me

fminconlinear contraintsMATLABOptimization Toolbox

Hello,
I'm trying to use fmincon to minimize the error between a designed filter and physical data by manipulating the filter FSTOP, FPASS and order. The contraint is that all F are in increasing order. So I gave the following contraint:
x0=[Fstop1 Fpass1 Fpass2 Fstop2 flt_ord];
% Inequalities
A=zeros(3,5); b=0*ones(3,1);
A(1,1)=1;A(1,2)=-1; %x1<x2
A(2,2)=1;A(2,3)=-1; %x2<x3
A(3,3)=1;A(3,4)=-1; %x3<x4
then I run fmincon(@myfun,x0,A,b).
My problem is that for some reason, i get that the F's are not in ascending order, and the filter designer shows and error and stops the entire script.
What am I doing wrong ??
tnx, Sivan

Best Answer

You could try as below. In other words, solve instead for p1=x1, p2=x2-x1, p3=x3-x2, p4=x4-x3 and just put positivity bounds on the increments p(2:4). As Alan said, fmincon won't enforce all constraints at all iterations, but the simple lb,ub bounds are an exception - under default settings, they will be enforced.
Some other remarks. Do not include flt_ord as one of the unknowns, since fmincon can't handle integer constraints - just do a looping search for that. Finally, use pchip interpolation in interp1, to keep myfun differentiable.
flt_orders=5:30; 
lb=[-inf,0,0,0];
ub=+inf(1,4);
p0=diff([0,x0(:).']);
for i=1:numel(flt_orders)
   [p{i} ,fval(i)] = fmincon(@(p)myfun(p, flt_orders(i)), p0, [],[],[],[],lb,ub)
end 
[~,imin]=min(fval);
x=cumsum(p{imin});
and
function f = myfun(p, flt_ord)
%physical data
Fs=22050;
load SinForTom
f=sinPlayed1*2/Fs; % sin freq[Hz]*2/sample rate [Hz]
a=rmsMS1./rmsLS1; %attentuation in linear scale
%filter params
Fstop1 = p(1);         % First Stopband Frequency
Fpass1 = Fstop1 + p(2);         % First Passband Frequency
Fpass2 = Fpass1 + p(3);         % Second Passband Frequency
Fstop2 = Fpass2 + p(4);         % Second Stopband Frequency
bp = designfilt('bandpassfir', 'FilterOrder', flt_ord, 'StopbandFrequency1', ...
                Fstop1, 'PassbandFrequency1', Fpass1, ...
                'PassbandFrequency2', Fpass2, 'StopbandFrequency2', ...
                Fstop2);
[h,w] = freqz(bp,1024);
habs=abs(h);
h_bela = interp1(w,habs,f,'pchip'); %find h at f points
f=sum((abs(a)-h_bela).^2); %calc lms
Related Question