Implementing the Gram-Schmidt Algorithm

Let's implement the notions of orthogonality and the Gram-Schmidt algorithm.
This is a very nice example of the beauty and power of functional programming.

Orthogonality

Let's implement the notion of an orthogonal list of vectors.

In[1]:=

orthogonalQ :: Usage = "orthogonalQ[vs_?MatrixQ] returns True if the vectors in the list of vectors vs are mutually orthogonal, and returns False otherwise." ;

(* Inductive definition:
      A list containing one vector is orthogonal.
      A list containing more than one vector is orthogonal if
         the first vector in the list is orthogonal to the remaining vectors, and
         removing the first vector from the list results in a (shorter) list of vectors which are mutually orthogonal.
  *)

In[2]:=

orthogonalQ[vs_ ? MatrixQ] := If[Length[vs] 1, True, orthogonalQ[Rest[vs]] &&orthogonalQ[Rest[vs], First[vs]]] 

orthogonalQ[vs_ ? MatrixQ, v_ ? ListQ] := If[Length[vs] 1, orthogonal2Q[First[vs], v], orthogonalQ[Rest[vs], v] &&orthogonal2Q[First[vs], v]] 

orthogonal2Q[u_ ? ListQ, v_ ? ListQ] := If[u . v0, True, False]

Testing.

In[5]:=

v1 = {1, 0, 0} ;

v2 = {0, 1, 0} ;

v3 = {0, 0, 1} ;

v4 = {2, 3, 4} ; 

orthogonal2Q[v1, v2]

orthogonal2Q[v3, v4]

Out[9]=

True

Out[10]=

False

In[11]:=

orthogonalQ[{v1, v2, v3}]

orthogonalQ[{v1, v2, v4}]

Out[11]=

True

Out[12]=

False

Projection of a vector onto a subspace spanned by an orthogonal basis

Let's implement the notion of the projection of a vector onto a subspace which is spanned by an orthogonal basis.

In[13]:=

proj :: Usage = "proj[us,u] returns the projection of the vector u onto the subspace spanned by the orthogonal basis us." ;

In[14]:=

proj[us_ ? orthogonalQ, u_ ? ListQ] := Module[{u1 = u . us[[1]]/us[[1]] . us[[1]] us[[1]]}, If[Length[us] 1, u1, u1 + proj[Rest[us], u]]]

Testing.

In[15]:=

proj[{v1}, v4]

proj[{v2}, v4]

proj[{v3}, v4]

Out[15]=

{2, 0, 0}

Out[16]=

{0, 3, 0}

Out[17]=

{0, 0, 4}

In[18]:=

proj[{v1, v2}, v4]

proj[{v1, v3}, v4]

proj[{v2, v3}, v4]

Out[18]=

{2, 3, 0}

Out[19]=

{2, 0, 4}

Out[20]=

{0, 3, 4}

In[21]:=

proj[{v1, v2, v3}, v4]

Out[21]=

{2, 3, 4}

Gram-Schmidt

Let's implement the Gram-Schmidt algorithm.

In[22]:=

gramSchmidt :: Usage = "gramSchmidt[vs] returns an orthogonal list of vectors, us, with ... ctors in vs have the same span as the first k vectors in us, for 1<=k<=Length[vs]." ;

In[23]:=

gramSchmidt[vs_ ? MatrixQ] := If[Length[vs] 1, vs, Module[{us = gramSchmidt[Drop[vs, -1]], u = vs[[-1]]}, Append[us, u - proj[us, u]]]]

Testing.

In[24]:=

v1 = {1, 1/2, 0} ;

v2 = {1, 2, 0} ;

v3 = {1, 2, 3} ;

vs = {v1, v2, v3} ; 

us = {u1, u2, u3} = gramSchmidt[vs] ;

Display the results

The Arrow3D package can be downloaded from Wolfram Research.
showColorful3DVectors is based on a similar procedure written by Selwyn Hollis.

In[29]:=

<<Arrow3D`Arrow3D`

In[30]:=

showColorful3DVectors[vecs_, shaftColor_, headColor_, opts___] := Show[Graphics3D[Flatten[{T ... rheadColor],  {i, Length[vecs]}]}]], opts, ViewPoint {6, 2, 2} ]

In[31]:=

<<Graphics`Graphics`

In[34]:=

o = {0, 0, 0} ; 

DisplayTogether[vPlot = showColorful3DVectors[{{o, v1}, {o, v2}, {o, v3}}, M ... lue], uPlot = showColorful3DVectors[{{o, u1}, {o, u2}, {o, u3}}, Orange, Red]] ;

[Graphics:HTMLFiles/6_implementingGramSchmidt_45.gif]

Where is the third blue vector?

Created by Mathematica  (April 8, 2005) Valid XHTML 1.1!