%?l1qc_logbarrier.m
%
%?Solve?quadratically?constrained?l1?minimization:
%?min?||x||_1???s.t.??||Ax?-?b||_2?<=?\epsilon
%
%?Reformulate?as?the?second-order?cone?program
%?min_{xu}??sum（u）???s.t.????x?-?u?<=?0
%????????????????????????????-x?-?u?<=?0
%??????1/2（||Ax-b||^2?-?\epsilon^2）?<=?0
%?and?use?a?log?barrier?algorithm.
%
%?Usage:??xp?=?l1qc_logbarrier（x0?A?At?b?epsilon?lbtol?mu?cgtol?cgmaxiter）
%
%?x0?-?Nx1?vector?initial?point.
%
%?A?-?Either?a?handle?to?a?function?that?takes?a?N?vector?and?returns?a?K?
%?????vector??or?a?KxN?matrix.??If?A?is?a?function?handle?the?algorithm
%?????operates?in?“largescale“?mode?solving?the?Newton?systems?via?the
%?????Conjugate?Gradients?algorithm.
%
%?At?-?Handle?to?a?function?that?takes?a?K?vector?and?returns?an?N?vector.
%??????If?A?is?a?KxN?matrix?At?is?ignored.
%
%?b?-?Kx1?vector?of?observations.
%
%?epsilon?-?scalar?constraint?relaxation?parameter
%
%?lbtol?-?The?log?barrier?algorithm?terminates?when?the?duality?gap?<=?lbtol.
%?????????Also?the?number?of?log?barrier?iterations?is?completely
%?????????determined?by?lbtol.
%?????????Default?=?1e-3.
%
%?mu?-?Factor?by?which?to?increase?the?barrier?constant?at?each?iteration.
%??????Default?=?10.
%
%?cgtol?-?Tolerance?for?Conjugate?Gradients;?ignored?if?A?is?a?matrix.
%?????Default?=?1e-8.
%
%?cgmaxiter?-?Maximum?number?of?iterations?for?Conjugate?Gradients;?ignored
%?????if?A?is?a?matrix.
%?????Default?=?200.
%
%?Written?by:?Justin?Romberg?Caltech
%?Email:?jrom@acm.caltech.edu
%?Created:?October?2005
%

function?xp?=?l1qc_logbarrier（x0?A?At?b?epsilon?lbtol?mu?cgtol?cgmaxiter）??

largescale?=?isa（A‘function_handle‘）;

if?（nargin?if?（nargin?if?（nargin?if?（nargin?
newtontol?=?lbtol;
newtonmaxiter?=?50;

N?=?length（x0）;

%?starting?point?---?make?sure?that?it?is?feasible
if?（largescale）
??if?（norm（A（x0）-b）?>?epsilon）
????disp（‘Starting?point?infeasible;?using?x0?=?At*inv（AAt）*y.‘）;
????AAt?=?@（z）?A（At（z））;
????[w?cgres]?=?cgsolve（AAt?b?cgtol?cgmaxiter?0）;
????if?（cgres?>?1/2）
??????disp（‘A*At?is?ill-conditioned:?cannot?find?starting?point‘）;
??????xp?=?x0;
??????return;
????end
????x0?=?At（w）;
??end
else
??if?（norm（A*x0-b）?>?epsilon）
????disp（‘Starting?point?infeasible;?using?x0?=?At*inv（AAt）*y.‘）;
????opts.POSDEF?=?true;?opts.SYM?=?true;
????[w?hcond]?=?linsolve（A*A‘?b?opts）;
????if?（hcond???????disp（‘A*At?is?ill-conditioned:?cannot?find?starting?point‘）;
??????xp?=?x0;
??????return;
????end
????x0?=?A‘*w;
??end??
end
x?=?x0;
u?=?（0.95）*abs（x0）?+?（0.10）*max（abs（x0））;

disp（sprintf（‘Original?l1?norm?=?%.3f?original?functional?=?%.3f‘?sum（abs（x0））?sum（u）））;

%?choose?initial?value?of?tau?so?that?the?duality?gap?after?the?first
%?step?will?be?about?the?origial?norm
tau?=?max（（2*N+1）/sum（abs（x0））?1）;
?????????????????????????????

?屬性????????????大小?????日期????時間???名稱
-----------?---------??----------?-----??----
?????目錄???????????0??2020-01-17?14:30??src\
?????文件????????3633??2020-01-17?14:18??src\l1qc_logbarrier.m
?????文件????????4565??2020-01-17?14:16??src\l1qc_newton.m
?????文件????????1453??2020-01-17?14:11??src\Main.m