In rete: http://www.elegio.it/omnia/ht/ortogonalemia.html

Tabella della matrice ortogonale...

Documento dedicato ad una matrice ortogonale ( https://it.wikipedia.org/wiki/Matrice_ortogonale ) poco usata ma comunque interessante come lo sono tutte le matrici ortogonali come l'utilissima matrice di Alston Scott Householder ( https://it.wikipedia.org/wiki/Alston_Scott_Householder ideatore di https://it.wikipedia.org/wiki/Trasformazione_di_Householder ) Anche questa matrice può essere usata per fare la decomposizione QR di una qualsiasi matrice ... ( https://it.wikipedia.org/wiki/Decomposizione_QR ) analogamente alla trasformazione di https://en.wikipedia.org/wiki/Wallace_Givens che lavorò nell' Argonne National Laboratory ( https://it.wikipedia.org/wiki/Argonne_National_Laboratory ) frequetato da Enrico Fermi quando, nel 1942, realizzò il primo reattore nucleare...

http://www.elegio.it/images/francobolli/
Qui mi dilungo molto sperando di essere chiaro nel descrivere questa insolita matrice ortogonale...

Per calcolare i coefficienti q1 ... q8 e f1 ... f8   dove gli fi vengono usati solo invertiti e 1/fi si annulla se un Vi è nullo

 q12 = V12;  f1= −q1*q2/V2
 q22 = q12 + V22;  f2= −q2*q3/V3
 q32 = q22 + V32;  f3= −q3*q4/V4
 q42 = q32 + V42;  f4= −q4*q5/V5
 q52 = q42 + V52;  f5= −q5*q6/V6
 q62 = q52 + V62;  f6= −q6*q7/V7
 q72 = q62 + V72;  f7= −q7*q8/V8
 q82 = q72 + V82;           f8= q8     

La matrice M è questa ipotizzando che il V1 sia l'elemento di massimo modulo del vettore assegnato

M =
┌
│
│
│
│
│
│
│
│
│
│
│
│
│
└		 V6/f8V7/f8V8/f8V1/f8V2/f8V3/f8V4/f8V5/f8 ┐
│
│
│
│
│
│
│
│
│
│
│
│
│
┘
000V1/f1q1/q2000
000V1/f2V2/f2q2/q300
000V1/f3V2/f3V3/f3q3/q40
000V1/f4V2/f4V3/f4V4/f4q4/q5
q5/q600V1/f5V2/f5V3/f5V4/f5V5/f5
V6/f6q6/q70V1/f6V2/f6V3/f6V4/f6V5/f6
V6/f7V7/f7q7/q8V1/f7V2/f7V3/f7V4/f7V5/f7

Dunque bisogna fare tante radici quadrate che ovviamente hanno un costo computazionale

q1 = |V1|
q2 = sqrt(V12 + V22)
q3 = sqrt(V12 + V22 + V32)
q4 = sqrt(V12 + V22 + V32 + V42)
q5 = sqrt(V12 + V22 + V32 + V42 + V52)
q6 = sqrt(V12 + V22 + V32 + V42 + V52 + V62)
q7 = sqrt(V12 + V22 + V32 + V42 + V52 + V62 + V72)
q8 = sqrt(V12 + V22 + V32 + V42 + V52 + V62 + V72 + V82)

Questa funzione calcola sia la matrice ortogonale sia i dati necessari per fare il prodotto della matrice ortogonale per un vettore qualsiasi in modo veloce ossia non facendo tante moltiplicazioni quanto il quadrato dell'ordine della matrice ma solo il triplo dell'ordine.

var fuma=function(Vt){
    var j,k,kn=0,uma=[[0]],uni=[[0]];
    var Q=[0],Qb=[0],Vq=0,Yn=-1;
    var uso=[[0],[0],[0]],V=[],nord=Vt.length;
    for(j=0;nord>j;j++){V[0]=Math.abs(Vt[j]);
       if(V[0]>Yn){kn=j;Yn=V[0]}}
    for(j=kn;nord>j;j++)V[j-kn]=Vt[j];
    for(j=0;kn>j;j++)V[nord-kn+j]=Vt[j];
    for(j=0;nord>j;j++){Vq+=V[j]*V[j];
       Q[j]=Math.sqrt(Vq);
       uma[j]=[0];
       for(k=0;nord>k;k++)uma[j][k]=0;
       uso[0][j]=V[j];}
    Yn=1/Q[nord-1];
    uso[1][0]=Yn;
    uso[2][0]=kn;
    for(k=0;nord>k;k++)uma[0][k]=Yn*V[k];
    for(j=1;nord>j;j++){uma[j][j]=Q[j-1]/Q[j];
       uso[2][j]=uma[j][j];
       Yn=-V[j]/(Q[j-1]*Q[j]);
       uso[1][j]=Yn;
       for(k=0;j>k;k++)uma[j][k]=Yn*V[k];}
    for(k=0;nord>k;k++){uni[k]=[];
       for(j=0;kn>j;j++)uni[k][j]=uma[k][nord-kn+j];
       for(j=kn;nord>j;j++)uni[k][j]=uma[k][j-kn];}
    return [uni,uso]; }

In pratica basta calcolare solo i dati da usare per fare velocemente il prodotto della matrice ortogonale con un qualsiasi vettore. Qui calcola solo i dati necessari ( un triplice vettore ).

var fuso=function(Vt){
    var j,kn=0,Q=[0],Vm,Vq=0,Yn=-1;
    var uso=[[0],[0],[0]],nord=Vt.length;
    for(j=0;nord>j;j++){uso[0][j]=0;
       Vm=Math.abs(Vt[j]);
       if(Vm>Yn){kn=j;Yn=Vm}}
    for(j=kn;nord>j;j++)uso[0][j-kn]=Vt[j];
    for(j=0;kn>j;j++)uso[0][nord-kn+j]=Vt[j];
    for(j=0;nord>j;j++){Vq+=uso[0][j]*uso[0][j];
       Q[j]=Math.sqrt(Vq);}
    uso[1][0]=1/Q[nord-1];
    uso[2][0]=kn;
    for(j=1;nord>j;j++){
       uso[2][j]=Q[j-1]/Q[j];
       uso[1][j]=-uso[0][j]/(Q[j-1]*Q[j]);}
    return uso; }

Qui invece usa i dati precalcolati per fare velocemente il prodotto della matrice ortogonale per un qualsiasi vettore.

var fusa=function(uso,Vt){
    var j,kn,S=0,R=[0],V=[0],n=Vt.length;
    kn=uso[2][0];
    for(j=kn;n>j;j++)V[j-kn]=Vt[j];
    for(j=0;kn>j;j++)V[n-kn+j]=Vt[j];
    S=uso[0][0]*V[0]
    for(j=1;n>j;j++){
       R[j]=S*uso[1][j]+V[j]*uso[2][j];
       S+=uso[0][j]*V[j]}
    R[0]=S*uso[1][0];
    return R;}

Ecco l'applicazione pratica della matrice ortogonale : Fattorizza la generica matrice mat con le matrici unitarie producendo come risultato una triangolare superiore mn e l'insieme delle matrice unitarie da usare ossia uso.

var risol=function(mat){
    var j,k,v,uso=[],mn=[[]],n=mat.length;
    var fus=function(Vt){
       var j,kn=0,Q=[0],Vm,Vq=0,Yn=-1;
       var uso=[[0],[0],[0]],nord=Vt.length;
       for(j=0;nord>j;j++){uso[0][j]=0;
          Vm=Math.abs(Vt[j]);
          if(Vm>Yn){kn=j;Yn=Vm}}
       for(j=kn;nord>j;j++)uso[0][j-kn]=Vt[j];
       for(j=0;kn>j;j++)uso[0][nord-kn+j]=Vt[j];
       for(j=0;nord>j;j++){Vq+=uso[0][j]*uso[0][j];
          Q[j]=Math.sqrt(Vq);}
       uso[1][0]=1/Q[nord-1];
       uso[2][0]=kn;
       for(j=1;nord>j;j++){
          uso[2][j]=Q[j-1]/Q[j];
          uso[1][j]=-uso[0][j]/(Q[j-1]*Q[j]);}
       return uso;};
    var fua=function(uso,Vt){
       var j,kn,S=0,R=[0],V=[0],n=Vt.length;
       kn=uso[2][0];
       for(j=kn;n>j;j++)V[j-kn]=Vt[j];
       for(j=0;kn>j;j++)V[n-kn+j]=Vt[j];
       S=uso[0][0]*V[0]
       for(j=1;n>j;j++){
          R[j]=S*uso[1][j]+V[j]*uso[2][j];
          S+=uso[0][j]*V[j]}
       R[0]=S*uso[1][0];
       return R;};
    for(j=0;n>j;j++){
       mn[j]=[mat[j][0]];
       for(k=1;n>k;k++)mn[j][k]=mat[j][k];}
    for(j=0;n>j;j++){
       v=[];
       for(k=j;n>k;k++)v[k-j]=mn[k][j];
       uso[j]=fus(v);
       for(h=j;n>h;h++){
          v=[];
          for(k=j;n>k;k++)v[k-j]=mn[k][h];
          vu=fua(uso[j],v);
          for(k=j;n>k;k++)mn[k][h]=vu[k-j];} }
    return [mn,uso];}

Fattorizzata, con la precedente funzione risol, la matrice mat che ora chiamo mfa, la si può usare per trovare la soluzione del sistema lineare, usando il termine noto vt. 

var usolu=function(mfa,vt){
    var j,k,tn,x,v,vu,vs=[0],n=vt.length;
    var sup=function(mup,vt){
       var j,k,tn,x,vs=[0],n=mup.length;
       for(j=0;n>j;j++)vs[j]=vt[j];
       for(j=n-1;j>=0;j+=-1){
          tn=vs[j];
          x=tn/mup[j][j];
          for(k=0;j>k;k++)vs[k]+=-x*mup[k][j];
          vs[j]=x;}
       return vs;}; 
    var fua=function(uso,Vt){
       var j,kn,S=0,R=[0],V=[0],n=Vt.length;
       kn=uso[2][0];
       for(j=kn;n>j;j++)V[j-kn]=Vt[j];
       for(j=0;kn>j;j++)V[n-kn+j]=Vt[j];
       S=uso[0][0]*V[0]
       for(j=1;n>j;j++){
          R[j]=S*uso[1][j]+V[j]*uso[2][j];
          S+=uso[0][j]*V[j]}
       R[0]=S*uso[1][0];
       return R;};
    for(j=0;n>j;j++)vs[j]=vt[j];
    for(j=0;n>j;j++){
       v=[];
       for(k=j;n>k;k++)v[k-j]=vs[k];
       vu=fua(mfa[1][j],v);
       for(k=j;n>k;k++)vs[k]=vu[k-j];
       }
    v=sup(mfa[0],vs);
    return v; }

Trascrivo qui anche alcune funzioni di servizio usabili per fare varie cose, come stampare bene una matrice o fare la matrice trasposta, il prodotto matrice per vettore e il prodotto matrice per matrice etc... 

var stamat=function(mat){
    var j,k,sp,ss="[",xbr="\u003cbr/\u003e";
    for(j=0;mat.length>j;j++){
        ss+=xbr+"[";
        for(k=0;mat.length>k;k++){
           if(1e-8>Math.abs(mat[j][k])){sp=" 0.0";}
           else{sp=" "+mat[j][k];}
           if(sp.length>10)sp=sp.substr(0,10);
           if(mat.length-1>k)ss+=sp+","
           else{ss+=sp}}
           ss+="]";
        }
    return ss+"]";
    }
//
var trasp=function(mat){
    var j,k,mtr=[[0]];
    for(j=0;mat.length>j;j++){
       mtr[j]=[mat[0][j]];
       for(k=1;mat.length>k;k++)mtr[j][k]=mat[k][j];
       }
    return mtr;
    }
//
var matvet=function(ma,vet){
    var j,h,to,vr=[0];
    for(j=0;ma.length>j;j++){
       to=0;
       for(h=0;ma.length>h;h++)to+=ma[j][h]*vet[h];
       vr[j]=to;
       }
    return vr;
    }
//
var matmul=function(ma,mb){
    var j,k,h,to,mr=[[0]];
    for(j=0;ma.length>j;j++){
       mr[j]=[0];
       for(k=0;ma.length>k;k++){
           to=0;
           for(h=0;ma.length>h;h++)to+=ma[j][h]*mb[h][k];
           mr[j][k]=to;
           }
       }
    return mr;
    }
//
var diffev=function(va,vb){
    var j,vr=[0],n=va.length;
    for(j=0;n>j;j++)vr[j]=va[j]-vb[j];
    return vr;
    }
//
var diffu=function(ua,ub){
    var j,dd=0,n=ua.length;
    if(ub.length!=n)alert("Errore! dimensioni diverse "+n);
    for(j=0;n>j;j++){
       dd+=Math.abs(ua[0][j]-ub[0][j])+
           Math.abs(ua[1][j]-ub[1][j])+
           Math.abs(ua[2][j]-ub[2][j]);
       }
    return dd;
    }

...Attivare Javascript per vedere i calcoli!...

Anche se già l'ho messa in rete trascrivo qui questa agile, breve funzione che fa l'inversa in modo molto numericamente costoso ma questa funzione mi sembra utile per trovare la matrice inversa di una matrice di ordine piccolo...


Inversione fatta da una unica funzione
// Calcolo della matrice inversa della matrice mat.
// Come si vede e' un algoritmo molto compatto anche
// se il numero delle operazioni da fare cresce molto
// al crescere dell'ordine della matrice mat. 
//
var fainv=function(mat){
    var minore,xmino=function(mat,r,c){
       var j,jn=0,k,kn,id=[],ord=mat.length;
       for(j=0;ord>j;j++){
          if(j!=r){ kn=0;
             id[id.length]=[];
             for(k=0;ord>k;k++){
                if(k!=c){
                   id[jn][kn]=mat[j][k];
                   kn++;} }
             jn++;} };
       return id; }
    var xdet=function(mat){
       var j,ord=mat.length,pm=1,dd=0;
       if(ord==1){
          dd=mat[0][0];
          return dd;};
       for(j=0;ord>j;j++){
          mm=xmino(mat,0,j);
          dd+=pm*mat[0][j]*xdet(mm);
          pm=-pm;}
       return dd; }
    var j,k,m,d,id=[],ord=mat.length;
    d=1/xdet(mat);
    for(j=0;ord>j;j++){
       m=d*(1-2*(j%2));
       id[id.length]=[];
       for(k=0;ord>k;k++){
          minore=xmino(mat,k,j);
          id[j][k]=m*xdet(minore);
          m=-m;} }
    return id; }

☺ Costosa ma elegante questa function, NO ??? ☺