Low dimensional matrices, vectors, quaternions and a cubic equation solver

Andrei Alexandrescu SeeWebsiteForEmail at erdani.org
Fri Apr 16 19:41:54 PDT 2010


On 04/16/2010 04:25 PM, Gareth Charnock wrote:
> Okay, here goes. I've collected together the basic functionality that
> would probably make a good starting point. As I've mentioned my code is
> very messy and has bits missing (e.g. I never had a use for the cross
> product but it's pretty important in general). I guess a good way to
> begin would be to write the pubic interfaces then start on the
> implementation.
>
> Cubic Solvers:
> General complex cubic solver with two algorithms (one requiring a
> complex cosine and and arcosine one using only +-*/ and roots). A
> special case cubic solver for the reduced cubic x^^3 + px - q = 0.

This sounds a candidate for std.numeric. How popular/needed are cubic 
solvers?

> Quaternions:
> opAdd, opSub, opMult(quaternion), opMult(vector), opDiv, Normalise,
> Normalized, conjugate, conjugated, toEulerAngles*, fromEulerAngles,
> this(real,i,j,k), this(angle,axis), getAngle(), getAxis()

Sounds good, but you'd need to convert the code to the new overloaded 
operators approach.

> *Currently I've only got a quaternion->euler angles routine that works
> in the ZYZ convention but I have read a paper that generalises my method
> to all valid axis conventions. Will probably impliment as something like:
> toEulerAngles(string convention="XYZ")()
>
> Vectors:
> opAdd, opSub, opMult(scalar), opMult(vector)*, cross**, Normalise,
> Normalized, Length

What is the representation of vectors? I'm afraid the design above would 
be too limited for what we need.

> * dot product. Would this be better named as dot()?

We already have dot product and normalization routines that work with 
general ranges.

http://www.digitalmars.com/d/2.0/phobos/std_numeric.html#dotProduct
http://www.digitalmars.com/d/2.0/phobos/std_numeric.html#normalize

Generally I'd strongly suggest making operations free generic functions 
instead of members.

> ** 3D vectors only. Perhaps defining a cross product on the
>
> Matrices:
> opAdd, opSub, opMult(scalar), opMult(vector), opMult(matrix)**, Invert,
> Inverted, Orthogonalize, Orthogonalized, Reorthogonalize***,
> Reorthogonalized***, Det, Transpose, Transposed, Dagger*, Daggered*,
> Eigenvalues****, Eigenvectors****

What is the representation of matrices?

> *The hermitian conjugate/conjugate transpose. Reduces to the transpose
> for a real matrix

Transposition should also be handled in a static manner, e.g. define a 
transposed view of a matrix that doesn't actually move elements.

> ** Matrix-matrix multiplication doesn't commute. Could this be a problem
> when using operator notation?

Should be fine.

> *** Othogonalize assuming the matrix is nearly orthogonal already
> (possibly using some quick, approximate method such as a Taylor series)
> **** I have a eigenvalue/vector solver for 3x3 matrices which seems
> reasonably stable but needs more testing.
>
> Free functions:
> MatrixToQuaternion
> QuaternionToMatrix
> + code to allow easy printing to stdout/streams and such.

Sounds encouraging.

I think a good next step is to go through a community scrutiny process 
by dropping the code somewhere on the Web so people can review it.


Andrei



More information about the Digitalmars-d mailing list