Another service from Omega

Orthogonal Projections, exp, leastsquares


*****

Orthogonal Projections:

Compute the orthogonal projection of the vector (2,1,-1,3) onto the plane generated by (-1,2,1,1) and (2,-1,-1,1).

Answer:

Orthonormalize the basis for the plane and then project.

> u1 := vector([-1,2,1,1]); u2 := vector([2,-1,-1,1]); u3 := vector([2,1,-1,3);

                             u1 := [ -1, 2, 1, 1 ]

                             u2 := [ 2, -1, -1, 1 ]

                             u3 := [ 2, 1, -1, 3 ]
> v1 := evalm(u1/sqrt(7));
                              1/2       1/2       1/2       1/2
               v1 := [ - 1/7 7   , 2/7 7   , 1/7 7   , 1/7 7    ]
> v2 := evalm( u2 - innerprod(u2,v1)*v1 );
                        v2 := [ 10/7, 1/7, -3/7, 11/7 ]

Now normalize to make it of length 1.

> v2 := evalm(v2/sqrt(sum(v2[i]^2,i=1..4)));

                  10    1/2           1/2            1/2          1/2
         v2 := [ --- 231   , 1/231 231   , - 1/77 231   , 1/21 231    ]
                 231

finally project u3 onto v1 and v2 to get the answer:

> sol := innerprod(u3,v1)*v1 + innerprod(u3,v2)*v2;

                                  1/2       19        1/2
                      sol := 2/7 7    v1 + ---- v2 231
                                            77
> sol := evalm(sol);
                                  24
                        sol := [ ----, 9/11, -5/11, 3 ]
                                  11

Now watch this:

> A := augment(u1,u2);

                                     [ -1   2 ]
                                     [        ]
                                     [  2  -1 ]
                                A := [        ]
                                     [  1  -1 ]
                                     [        ]
                                     [  1   1 ]
> SOL := evalm( A &* leastsqrs(A,u3) );
                                  24
                        SOL := [ ----, 9/11, -5/11, 3 ]
                                  11

Hey did you see that? It can be done in only one line!

The exponential of a matrix.

Problem:

Compute the exponential of the matrix B where:

> B := matrix(3,3,[1,0,0,0,1,2,0,2,1]);

                                     [ 1  0  0 ]
                                     [         ]
                                B := [ 0  1  2 ]
                                     [         ]
                                     [ 0  2  1 ]

Solution:

You first need to diagonalize B. We get a diagonal matrix G and e non singular matrix P such that, B = P G P^(-1)

> eigensys := eigenvects(B);

   eigensys :=

       [-1, 1, {[ 0, -1, 1 ]}], [3, 1, {[ 0, 1, 1 ]}], [1, 1, {[ 1, 0, 0 ]}]
> G := diag(-1,3,1);
                                    [ -1  0  0 ]
                                    [          ]
                               G := [  0  3  0 ]
                                    [          ]
                                    [  0  0  1 ]
> P := augment([0,-1,1],[0,1,1],[1,0,0]);
                                    [  0  0  1 ]
                                    [          ]
                               P := [ -1  1  0 ]
                                    [          ]
                                    [  1  1  0 ]
> evalm(P &* G &* inverse(P));
                                  [ 1  0  0 ]
                                  [         ]
                                  [ 0  1  2 ]
                                  [         ]
                                  [ 0  2  1 ]

We got 'em.

> ExpB := evalm(P &* diag(exp(-1),exp(3),exp(1)) &* inverse(P));

           [ exp(1)               0                           0             ]
           [                                                                ]
   ExpB := [    0     1/2 exp(-1) + 1/2 exp(3)   - 1/2 exp(-1) + 1/2 exp(3) ]
           [                                                                ]
           [    0    - 1/2 exp(-1) + 1/2 exp(3)   1/2 exp(-1) + 1/2 exp(3)  ]
> evalf(");
                   [ 2.718281828       0            0      ]
                   [                                       ]
                   [      0       10.22670818  9.858828739 ]
                   [                                       ]
                   [      0       9.858828739  10.22670818 ]

Link to the commands in this file
Carlos Rodriguez <carlos@math.albany.edu>
Last modified: Thu Jun 19 11:46:28 EDT 1997