%第一步：寫一個m文件 brownvv.m 計算目標函數值、梯度以及海色矩陣值。
%brownvv.m 是計算目標函數文件，它需傳遞給優化函數fmincon。
%type brownvv
%因為我們使用brownfgh計算目標函數值，梯度和海色矩陣，你需要在使用optimset命
%令時使用GradObj 和海色矩陣參數。

%第二步：定義一個函數hmfleq1，它利用在ibrownvv計算的矩陣Hinfo和V，計算海
%色矩陣乘積W。
W = H*Y = (Hinfo - V*V')*Y.
%此函數的形式為：
W = hmfleq1(Hinfo,Y,p1,p2...)
%第一個參數即為函數返回的第三個參數，第二個參數為矩陣Y。
%因為函數fmincon使用參數Y去計算海色矩陣的乘積，Y是一個n行矩陣，n即為
%此問題的維數。Y的列數是可變的。最后，任何附加的參數都會傳遞給函數hmfleq1，
%于是hmfleq1接受了相同的參數例如矩陣V。

function W = hmfleq1(Hinfo,Y,V);
W = Hinfo*Y - V*(V'*Y);

%第三步：使用命令fmincon 解決此問題。
%從文件中取得V, Aeq, beq的值。
%使用命令 optimset 設置GradObj 和Hessian參數為 'on' ，設定HessMult
％參數為函數句柄，指向函數 hmfleq1。令目標函數和矩陣V為優化函數的參數，計算最優結
%果。

% 得到 V, Aeq, beq值
load fleq1 
% 問題維數
n = 1000; 
mtxmpy = @hmfleq1; % Function handle to function hmfleq1
xstart = -ones(n,1); 
xstart(2:2:n,1) = ones(length(2:2:n),1);
options = optimset('GradObj','on','Hessian','on',...
'HessMult',mtxmpy,'Display','iter');
[x,fval,exitflag,output] = fmincon(@brownvv,xstart,[],[],Aeq,beq,[],[],[], ... options,V);