Design a 1/f spectrum shaping (pink-noise) filter
n = 40;
m = 15*n;
wa = 0.01*pi; wb = pi;
wl = logspace(log10(wa),log10(wb),m)';
D = 1./sqrt(wl);
Al = [ones(m,1) 2*cos(kron(wl,[1:n-1]))];
cvx_begin
variable r(n,1)
variable R(m,1)
minimize( max( max( [R./(D.^2) (D.^2).*inv_pos(R)]' ) ) )
subject to
R == Al*r;
cvx_end
disp(['Problem is ' cvx_status])
if ~strfind(cvx_status,'Solved')
return
end
h = spectral_fact(r);
figure(1)
H = exp(-j*kron(wl,[0:n-1]))*h;
loglog(wl,abs(H),wl,D,'r--')
set(gca,'XLim',[wa pi])
xlabel('freq w')
ylabel('mag H(w) and D(w)')
legend('optimized','desired')
Calling sedumi: 3000 variables, 641 equality constraints
For improved efficiency, sedumi is solving the dual problem.
------------------------------------------------------------
SeDuMi 1.21 by AdvOL, 2005-2008 and Jos F. Sturm, 1998-2003.
Alg = 2: xz-corrector, Adaptive Step-Differentiation, theta = 0.250, beta = 0.500
eqs m = 641, order n = 1803, dim = 3602, blocks = 602
nnz(A) = 27149 + 0, nnz(ADA) = 51481, nnz(L) = 26061
Handling 1 + 1 dense columns.
it : b*y gap delta rate t/tP* t/tD* feas cg cg prec
0 : 5.19E+01 0.000
1 : -5.87E+00 2.04E+01 0.000 0.3939 0.9000 0.9000 -0.97 1 1 5.9E+02
2 : -3.67E+01 3.72E+00 0.000 0.1819 0.9000 0.9000 -0.86 1 1 4.1E+02
3 : -5.23E+01 1.81E-01 0.000 0.0488 0.9900 0.9900 0.23 1 1 2.5E+01
4 : -9.93E+00 8.50E-02 0.000 0.4687 0.9000 0.9000 8.65 1 1 2.3E+00
5 : -2.83E+00 3.12E-02 0.000 0.3665 0.9000 0.9000 3.70 1 1 3.4E-01
6 : -1.65E+00 1.14E-02 0.000 0.3645 0.9000 0.9000 1.34 1 1 1.1E-01
7 : -1.31E+00 4.78E-03 0.000 0.4205 0.9000 0.9000 0.96 1 1 5.0E-02
8 : -1.22E+00 2.30E-03 0.000 0.4809 0.6154 0.9000 0.95 1 1 2.4E-02
9 : -1.21E+00 1.28E-03 0.000 0.5593 0.9000 0.6822 1.00 1 1 1.4E-02
10 : -1.19E+00 3.98E-04 0.000 0.3102 0.9000 0.9000 1.00 1 1 4.2E-03
11 : -1.19E+00 1.31E-04 0.000 0.3291 0.9000 0.9000 1.00 1 1 1.4E-03
12 : -1.19E+00 4.30E-05 0.000 0.3281 0.9000 0.9000 1.00 1 1 4.6E-04
13 : -1.19E+00 1.17E-05 0.000 0.2713 0.9024 0.9000 1.00 1 1 1.2E-04
14 : -1.19E+00 2.32E-06 0.000 0.1987 0.9084 0.9000 1.00 1 1 1.9E-05
15 : -1.19E+00 5.26E-07 0.000 0.2270 0.9088 0.9000 1.00 1 1 3.6E-06
16 : -1.19E+00 4.17E-08 0.465 0.0792 0.9900 0.9901 1.00 1 1 3.1E-07
17 : -1.19E+00 1.02E-08 0.000 0.2451 0.9000 0.9097 1.00 1 1 7.9E-08
18 : -1.19E+00 1.41E-10 0.000 0.0138 0.9901 0.9900 1.00 2 2 1.0E-09
iter seconds digits c*x b*y
18 1.4 9.6 -1.1873310834e+00 -1.1873310837e+00
|Ax-b| = 5.0e-10, [Ay-c]_+ = 1.2E-11, |x|= 1.2e+01, |y|= 2.6e+02
Detailed timing (sec)
Pre IPM Post
5.000E-02 1.390E+00 1.000E-02
Max-norms: ||b||=3.183099e-01, ||c|| = 2,
Cholesky |add|=0, |skip| = 0, ||L.L|| = 6.28317.
------------------------------------------------------------
Status: Solved
Optimal value (cvx_optval): +1.18733
Problem is Solved