> 
> read("/afs/eos.ncsu.edu/users/k/kaltofen/www/courses/LinAlgebra/Maple/
> LsqPkg/InitPkg.mpl");

  libname := /afs/eos.ncsu.edu/users/k/kaltofen/www/courses/LinAlg\

        ebra/Maple, "/afs/bp/dist/maple60/lib"

> with(LsqPkg);

                            [Gram_Schmidt]

> ?Gram_Schmidt
> with(LinearAlgebra):
> A :=<<1,1,1,1>|<1,2,3,4>|<1,4,9,16>>;

                              [1    1     1]
                              [            ]
                              [1    2     4]
                         A := [            ]
                              [1    3     9]
                              [            ]
                              [1    4    16]

> Q := Gram_Schmidt(A, 'R');

                             [1    -3/2     1]
                             [               ]
                             [1    -1/2    -1]
                        Q := [               ]
                             [1    1/2     -1]
                             [               ]
                             [1    3/2      1]

> R;

                          [1    5/2    15/2]
                          [                ]
                          [0     1      5  ]
                          [                ]
                          [0     0      1  ]

> Q . R;

                            [1    1     1]
                            [            ]
                            [1    2     4]
                            [            ]
                            [1    3     9]
                            [            ]
                            [1    4    16]

> Diag := Transpose(Q) . Q;

                                [4    0    0]
                                [           ]
                        Diag := [0    5    0]
                                [           ]
                                [0    0    4]

> b := Vector([0,4,9,16]);

                                   [ 0]
                                   [  ]
                                   [ 4]
                              b := [  ]
                                   [ 9]
                                   [  ]
                                   [16]

> yhat := Diag^(-1) . Transpose(Q) . b;

                                    [29/4]
                                    [    ]
                                    [ 53 ]
                            yhat := [ -- ]
                                    [ 10 ]
                                    [    ]
                                    [3/4 ]

> bhat := Q . yhat;

                                    [1/20]
                                    [    ]
                                    [ 77 ]
                                    [ -- ]
                                    [ 20 ]
                                    [    ]
                            bhat := [183 ]
                                    [--- ]
                                    [20  ]
                                    [    ]
                                    [319 ]
                                    [--- ]
                                    [20  ]

> LinearSolve(R,yhat);

                                [-9/4]
                                [    ]
                                [ 31 ]
                                [ -- ]
                                [ 20 ]
                                [    ]
                                [3/4 ]

>