function?[?V??D?]?=??joint_diag（Ajthresh）
%?Joint?approximate?diagonalization
%?聯合近似對角化
%
%?Joint?approximate?of?n?（complex）?matrices?of?size?m*m?stored?in?the
%?m*mn?matrix?A?by?minimization?of?a?joint?diagonality?criterion
%?對儲存在A中的n個mxm矩陣做聯合對角化
%
%?Usage:??[?V??D?]?=??joint_diag（Ajthresh）
%
%?Input?:
%?*?the?m*nm?matrix?A?is?the?concatenation?of?n?matrices?with?size?m
%???by?m.?We?denote?A?=?[?A1?A2?....?An?]
%?*?threshold?is?an?optional?small?number?（typically?=?1.0e-8?see?the?M-file）.
%
%?Output?:
%?*?V?is?an?m*m?unitary?matrix.
%?*?D?=?V‘*A1*V??...??V‘*An*V?has?the?same?size?as?A?and?is?a
%???collection?of?diagonal?matrices?if?A1?...?An?are?exactly?jointly
%???unitarily?diagonalizable.
%
%?The?algorithm?finds?a?unitary?matrix?V?such?that?the?matrices
%?V‘*A1*V??...??V‘*An*V?are?as?diagonal?as?possible?providing?a
%?kind?of?‘average?eigen-structure‘?shared?by?the?matrices?A1?...An.
%?If?the?matrices?A1...An?do?have?an?exact?common?eigen-structure?ie
%?a?common?orthonormal?set?eigenvectors?then?the?algorithm?finds?it.
%?The?eigenvectors?THEN?are?the?column?vectors?of?V?and?D1?...Dn?are
%?diagonal?matrices.
%?聯合對角化的作用。
%
%?The?algorithm?implements?a?properly?extended?Jacobi?algorithm.??The
%?algorithm?stops?when?all?the?Givens?rotations?in?a?sweep?have?sines
%?smaller?than?‘threshold‘.?
%?算法在所有旋轉角度都小于給定閾值時停止
%
%?In?many?applications?the?notion?of?approximate?joint
%?diagonalization?is?ad?hoc?and?very?small?values?of?threshold?do?not
%?make?sense?because?the?diagonality?criterion?itself?is?ad?hoc.
%?Hence?it?is?often?not?necessary?in?applications?to?push?the
%?accuracy?of?the?rotation?matrix?V?to?the?machine?precision.
%?閾值的設置不用太經心。
%?PS:?If?a?numrical?analyst?knows?‘the?right?way‘?to?determine?jthresh
%?????in?terms?of?1）?machine?precision?and?2）?size?of?the?problem
%?????I?will?be?glad?to?hear?about?it.
%
%?This?version?of?the?code?is?for?complex?matrices?but?it?also?works
%?with?real?matrices.??However?simpler?implementations?are?possible
%?in?the?real?case.?
%?這個代碼適合復矩陣，當然也適合實矩陣；但實矩陣有更簡潔的程序實現。
%
%?See?more?info?references?and?version?history?at?the?bottom?of?this
%?m-file
%----------------------------------------------------------------
%?Version?1.2
%
%?Copyright?	:?Jean-Francois?Cardoso.
%?Author?	:?Jean-Francois?Cardoso.?cardoso@sig.enst.fr
%?Comments?bug?reports?etc?are?welcome.
%----------------------------------------------------------------

[mnm]?=?size（A）;
B?=?[?1?0?0?;?0?1?1?;?0?-1i?1i?]?;
Bt?=?B‘?;
V?=?eye（m）;
encore	=?1;
while?encore
????encore=0;
????for?p=1:m-1
????????Ip?=?p:m:nm;
????????for?q=p+1:m
????????????Iq?=?q:m:nm;
????????????
????????????%?Computing?the?Givens?angles?計算旋轉角度
????????????g?=?[A（pIp）-A（qIq）;?A（pIq）;?A（qIp）];
????????????[vcpD]?=?eig（real（B*（g*g‘）*Bt））;
????????????[~?K]?=?sort（diag（D））;
????????????angles??=?vcp（:K（3））;
????????????if?angles（1）<0
????????????????angles=?-angles;
????????????end
????????????c?=?sqrt（0.5+angles（1）/2）;
????????????s?=?0.5*（angles（2）-1j*

?屬性????????????大小?????日期????時間???名稱
-----------?---------??----------?-----??----
?????文件????????4903??2013-06-06?11:03??joint_diag.m
?????文件?????????165??2013-06-06?14:49??joint_diag_Note.txt
?????文件????????3813??2013-06-06?11:05??joint_diag_r.m
?????文件????????1021??2013-06-06?11:00??test_joint_diag.m
?????文件?????????970??2013-06-06?11:08??test_joint_diag_r.m