[comp.graphics] Spline coefficients

bls@texhrc.UUCP (Brian L. Sumner) (11/16/89)

I am tring to determine a set of spline coefficients, and to do
so I need to solve the linear system
    Ax = b
where A can be a very large (positive definite?) banded matrix.
This matrix looks like

    [ A B . . . ]
    [ B C D . . ]
    [ . D E F . ]
    [ . . F G H ] 
    [ . . . H I ]

where each of A,B, ... I, is itself a (positive definite?) banded matrix
with the same number of bands, e.g.,

    [ a b . . . ]
    [ b c d . . ]
    [ . d e f . ]
    [ . . f g h ] 
    [ . . . h i ]

where each of a,b, ... , i, is again a positive definite?) banded matrix
with the sme number of bands, ...

The number of levels this continues depends on the dimensionality
of the problem, but in practice will be 1,2,3, or 4.

I do not believe that I should just throw such a matrix into a banded
solver.  The regular structure of this matrix begs for a special (faster)
technique.  If you are aware of a method to attack this problem with
I would appreciate it.

Brian Sumner
...!nuchat!convex!texhrc!bls

woody@rpp386.cactus.org (Woodrow Baker) (11/17/89)

In article <352@texhrc.UUCP>, bls@texhrc.UUCP (Brian L. Sumner) writes:
> 
> I am tring to determine a set of spline coefficients, and to do
> so I need to solve the linear system
>     Ax = b
> where A can be a very large (positive definite?) banded matrix.
> This matrix looks like
> 
>     [ A B . . . ]
>     [ B C D . . ]
>     [ . D E F . ]
>     [ . . F G H ] 
>     [ . . . H I ]
> 
> where each of A,B, ... I, is itself a (positive definite?) banded matrix
> with the same number of bands, e.g.,
> 
>     [ a b . . . ]
>     [ b c d . . ]
>     [ . d e f . ]
>     [ . . f g h ] 
>     [ . . . h i ]
> 
> where each of a,b, ... , i, is again a positive definite?) banded matrix
> with the sme number of bands, ...
> 
> The number of levels this continues depends on the dimensionality
> of the problem, but in practice will be 1,2,3, or 4.
> 
> I do not believe that I should just throw such a matrix into a banded
> solver.  The regular structure of this matrix begs for a special (faster)
> technique.  If you are aware of a method to attack this problem with
> I would appreciate it.
> 
> Brian Sumner
> ...!nuchat!convex!texhrc!bls

I am currently struggling with the same problem.  Depending on what you are
trying to do, there are several ways to approach it.  Dr. Dobbs had an article
published a while back, written by a man named Ashdown that dealt explicitly
with this problem.  A deeper article on pur curvetracing using splines,
"Curve fitting with piecewise splines" or somthing like that by Micheal Plass
and M. Stone from Xerox Palo Alto Reasearch center was published in Siggraph
in '83 or '84.  I have a copy at home, as well as some nice reference work
on splines.  There are several nice books on the subject.  My project involvs
being able to take a scanned image, and curve trace it.  I have a rudimentry
spline drawing program that uses a mouse and does splines sort of.  Computing
the coefficients of the spline given the x and y coordinates of the control
points, and the end points, is fairly straight forward, and then evaluating
the spline is very straight forward.  My problem is that I'm lousy at math, and
could not even solve the redbook equations and re-arrange them such that
I could solve for the coefficients.  That has been done now, and I'll be
glad to post them.  Don Lancaster has published them in last months (or
maybe this months) computer shopper.  If you are trying to fit a curve to data
points, then the job becomes much harder.   Basically, you are on the right
track.  You can use a Gauss-Jordan elimination to solve the matrix.  Be sure
to use a piviot row modification however.  Alan Ashdons' article handles
the solving of the matrix nicely.  Plass's paper details an algorithm
that iterativly moves the data points closer to the line.  If you know the
total line length, then you can get a close approximation for t by taking
the distance between each line using a simple difrence squared and doing
a square root, to get delta x and delta y.  Then you can compute the length
of the line between the data points, and travel the curve, adding them up.
This will give you an APPROXIMATION of    t   and you can then divide any
individual segment by this value to get an approximate t for a given x and
y coordinate.  Then you do a least squares fit of the line, and iterativly
drive the x and y coordinates to the closes t value to them along the line.

This is sort of a off the top of the head crude overview of the Plass stuff.
I'd be glad to send a copy to anyone who needs it, just include a stamp....

Cheers

Woody

rhbartels@watcgl.waterloo.edu (Richard Bartels) (11/20/89)

In article <352@texhrc.UUCP> bls@texhrc.UUCP (Brian L. Sumner) writes:
>
>I am tring to determine a set of spline coefficients, and to do
>so I need to solve the linear system
>    Ax = b
>where A can be a very large (positive definite?) banded matrix.
>
>[   ... description deleted ... ]
>

The description sounds suspiciously like the sort of matrix that arises
from interpolating or fitting over a rectilinear mesh (a tensor-product grid).
If this is the case, then deBoor's book, "A Practical Guide to Splines"
and the following two articles might be instructive.  Briefly, a tensor-product
problem in n dimensions can be reduced to n successive one-dimensional
problems.  With B-splines, a each one-dimensional problem interpolation
problem, for example, typically involves a k-banded matrix, where k is the
degree of the spline; e.g.,  cubics involve tridiagonal matrices.
Such systems can be solved in linear time.  deBoor's book offers code
(also available on NETLIB), and the articles show how the tensor-product
insight extends to least squares fitting (as well as being more readible
than the explanation given by deBoor).

-Richard
==========================================================================
%A P. Dierckx
%L Dierckx77-
%T An algorithm for least-squares fitting of cubic spline surfaces
to functions on a rectilinear mesh over a rectangle
%J Journal of Computational and Applied Mathematics
%V 3
%N 2
%D 1977
%P 113-129
%l journal-article

%A P. Dierckx
%L Dierckx81-
%T An Algorithm for Surface-Fitting with Spline Functions
%J IMA Journal of Numerical Analysis
%V 1
%D 1981
%P 267-283
%l journal-article

woody@rpp386.cactus.org (Woodrow Baker) (11/21/89)

In article <12349@watcgl.waterloo.edu>, rhbartels@watcgl.waterloo.edu (Richard Bartels) writes:
> In article <352@texhrc.UUCP> bls@texhrc.UUCP (Brian L. Sumner) writes:
> >
> >I am tring to determine a set of spline coefficients, and to do
> >so I need to solve the linear system
> >    Ax = b
> >where A can be a very large (positive definite?) banded matrix.
> >
> >[   ... description deleted ... ]
> >
> 
> The description sounds suspiciously like the sort of matrix that arises
> from interpolating or fitting over a rectilinear mesh (a tensor-product grid).
> If this is the case, then deBoor's book, "A Practical Guide to Splines"
> and the following two articles might be instructive.  Briefly, a tensor-product
> problem in n dimensions can be reduced to n successive one-dimensional
> problems.  With B-splines, a each one-dimensional problem interpolation
> problem, for example, typically involves a k-banded matrix, where k is the
> degree of the spline; e.g.,  cubics involve tridiagonal matrices.
> Such systems can be solved in linear time.  deBoor's book offers code
> (also available on NETLIB), and the articles show how the tensor-product
> insight extends to least squares fitting (as well as being more readible
> than the explanation given by deBoor).
> 
> -Richard
> ==========================================================================
> %A P. Dierckx
> %L Dierckx77-
> %T An algorithm for least-squares fitting of cubic spline surfaces
> to functions on a rectilinear mesh over a rectangle
> %J Journal of Computational and Applied Mathematics
> %V 3
> %N 2
> %D 1977
> %P 113-129
> %l journal-article
> 
> %A P. Dierckx
> %L Dierckx81-
> %T An Algorithm for Surface-Fitting with Spline Functions
> %J IMA Journal of Numerical Analysis
> %V 1
> %D 1981
> %P 267-283
> %l journal-article

How do you access NETLIB.  I'd like specific details, as I am sort of a
novice at using the net at this time.  Can post or email, or call.  Thanks

Woody Baker (512)-837-8317

peter@mit-amt.MEDIA.MIT.EDU (Peter Schroeder) (11/22/89)

to acces netlib use the address
	netlib@anl-mcs.arpa

This is an automatic system. Start with a one line message
	send index
and you'll be able to figure it out from there.

It has lots of very good codes and I highly recommend it.

Peter

rhbartels@watcgl.waterloo.edu (Richard Bartels) (11/22/89)

In article <1064@mit-amt.MEDIA.MIT.EDU> peter@media-lab.media.mit.edu
(Peter Schroeder) writes:
>to acces netlib use the address
>	netlib@anl-mcs.arpa
>

In fact, there are two addresses to NETLIB, and the one Peter gives
may be in a state of flux, since its manager, Jack Dongarra, is moving
from "anl-mcs" (Argonne National Laboratories - Math and Computer Science)
to the University of Tenn.  The other address is at AT&T/Bell Labs
and is under the management of Eric Grosse.  Together Jack and Eric provide
a valuable service to the world and are due a round of apple sauce.
Below is more information on NETLIB, with the AT&T address included.
It is the header of the file you will receive if you mail the request
"send index".

-Richard

===== How to use netlib =====

This file is the reply you'll get to:
	mail netlib@research.att.com
	send index
Here are examples of the various kinds of requests.
*  get the full index for a library
	send index from eispack
*  get a particular routine and all it depends on
	send dgeco from linpack
*  get just the one routine, not subsidiaries
	send only dgeco from linpack
*  get dependency tree, but excluding a subtree
	send dgeco but not dgefa from linpack
*  just tell how large a reply would be, don't actually send the file
	send list of dgeco from linpack
*  search for somebody in the SIAM membership list:
	who is gene golub
*  keyword search for netlib software  (somewhat out of date)
	find cubic spline
*  keyword search in the approximation catalog
	find schumaker from approximation
*  notes on the Fortran-to-C converter  (netlib@research.att.com only)
	send index for f2c

The Internet address "netlib@research.att.com" refers to a gateway
machine at AT&T Bell Labs in Murray Hill, New Jersey.  This address
should be understood on all the major networks.  For systems having
only uucp connections, use the address uunet!research!netlib.  In this
case, someone will be paying for long distance 1200 baud phone calls,
so keep your requests to a reasonable size!  An excellent guide to the
mysteries of networks and address syntax is: Donnalyn Frey and Rick
Adams (1989) "!%@:: A Directory of Electronic Mail Addressing and
Networks", O'Reilly & Associates, Inc, 632 Petaluma Ave, Sebastopol CA
95472.  Background about netlib is in Jack J. Dongarra and Eric Grosse,
Distribution of Mathematical Software Via Electronic Mail, Comm. ACM
(1987) 30, 403--407.

Bugs reports, comments, and annual lists of recipients will be
forwarded to the code authors when possible.  Many of these codes are
designed for use by professional numerical analysts who are capable of
checking for themselves whether an algorithm is suitable for their
needs.  One routine can be superb and the next awful.  So be careful!

[ further information deleted (long) ]