Inverse of 3x3 Matrix with Complex Numbers

[deetee]

Hello!

My next calculator should be able to work with 3x3 matrices with complex numbers (determinant, invert, add, substract, multiply, ...). That for I want to impement fast explicit formulae (found on Wikipedia) rather than some algorithms.

The calculator itself works with two double stacks for real and imaginary parts of numbers (stack[4] and stacki[4]). Basic stack operations are done with _add(), _sub() and _mult().
For matrix operations I use matrix A (ma[3][3] and mai[3][3]) and B (mb[3][3] and mbi[3][3]) and a temporary matrix C (mc[3][3] and mci[3][3]).

So far I implemented the following code which works fine:

Code:

static void _mdet(void) { // MATRIX DET A
  mdet = mdeti = 0.0;
  stack[2] = ma[0][0]; stacki[2] = mai[0][0];
  stack[1] = ma[1][1]; stacki[1] = mai[1][1];
  stack[0] = ma[2][2]; stacki[0] = mai[2][2];
  _mult(); _mult(); mdet += stack[0]; mdeti += stacki[0];
  stack[2] = ma[0][1]; stacki[2] = mai[0][1];
  stack[1] = ma[1][2]; stacki[1] = mai[1][2];
  stack[0] = ma[2][0]; stacki[0] = mai[2][0];
  _mult(); _mult(); mdet += stack[0]; mdeti += stacki[0];
  stack[2] = ma[0][2]; stacki[2] = mai[0][2];
  stack[1] = ma[1][0]; stacki[1] = mai[1][0];
  stack[0] = ma[2][1]; stacki[0] = mai[2][1];
  _mult(); _mult(); mdet += stack[0]; mdeti += stacki[0];
  stack[2] = ma[2][0]; stacki[2] = mai[2][0];
  stack[1] = ma[1][1]; stacki[1] = mai[1][1];
  stack[0] = ma[0][2]; stacki[0] = mai[0][2];
  _mult(); _mult(); mdet -= stack[0]; mdeti -= stacki[0];
  stack[2] = ma[2][1]; stacki[2] = mai[2][1];
  stack[1] = ma[1][2]; stacki[1] = mai[1][2];
  stack[0] = ma[0][0]; stacki[0] = mai[0][0];
  _mult(); _mult(); mdet -= stack[0]; mdeti -= stacki[0];
  stack[2] = ma[2][2]; stacki[2] = mai[2][2];
  stack[1] = ma[1][0]; stacki[1] = mai[1][0];
  stack[0] = ma[0][1]; stacki[0] = mai[0][1];
  _mult(); _mult(); mdet -= stack[0]; mdeti -= stacki[0];
  stack[0] = mdet; stacki[0] = mdeti;
}
static void _minv(void) { // MATRIX INVERSE A
  _mdet();
  stack[3] = ma[1][1]; stacki[3] = mai[1][1]; // C00
  stack[2] = ma[2][2]; stacki[2] = mai[2][2];
  stack[1] = ma[2][1]; stacki[2] = mai[2][1];
  stack[0] = ma[1][2]; stacki[3] = mai[1][2];
  _mult(); _rot(); _mult(); _rotup(); _sub();
  mc[0][0] = stack[0]; mci[0][0] = stacki[0];
  stack[3] = ma[0][2]; stacki[3] = mai[0][2]; // C01
  stack[2] = ma[2][1]; stacki[2] = mai[2][1];
  stack[1] = ma[0][1]; stacki[2] = mai[0][1];
  stack[0] = ma[2][2]; stacki[3] = mai[2][2];
  _mult(); _rot(); _mult(); _rotup(); _sub();
  mc[0][1] = stack[0]; mci[0][1] = stacki[0];
  stack[3] = ma[0][1]; stacki[3] = mai[0][1]; // C02
  stack[2] = ma[1][2]; stacki[2] = mai[1][2];
  stack[1] = ma[0][2]; stacki[2] = mai[0][1];
  stack[0] = ma[1][1]; stacki[3] = mai[1][2];
  _mult(); _rot(); _mult(); _rotup(); _sub();
  mc[0][2] = stack[0]; mci[0][2] = stacki[0];
  stack[3] = ma[1][2]; stacki[3] = mai[1][2]; // C10
  stack[2] = ma[2][0]; stacki[2] = mai[2][0];
  stack[1] = ma[1][0]; stacki[2] = mai[1][2];
  stack[0] = ma[2][2]; stacki[3] = mai[2][0];
  _mult(); _rot(); _mult(); _rotup(); _sub();
  mc[1][0] = stack[0]; mci[1][0] = stacki[0];
  stack[3] = ma[0][0]; stacki[3] = mai[0][0]; // C11
  stack[2] = ma[2][2]; stacki[2] = mai[2][2];
  stack[1] = ma[0][2]; stacki[2] = mai[0][0];
  stack[0] = ma[2][0]; stacki[3] = mai[2][2];
  _mult(); _rot(); _mult(); _rotup(); _sub();
  mc[1][1] = stack[0]; mci[1][1] = stacki[0];
  stack[3] = ma[0][2]; stacki[3] = mai[0][2]; // C12
  stack[2] = ma[1][0]; stacki[2] = mai[1][0];
  stack[1] = ma[0][0]; stacki[2] = mai[0][0];
  stack[0] = ma[1][2]; stacki[3] = mai[1][2];
  _mult(); _rot(); _mult(); _rotup(); _sub();
  mc[1][2] = stack[0]; mci[1][2] = stacki[0];
  stack[3] = ma[1][0]; stacki[3] = mai[1][0]; // C20
  stack[2] = ma[2][1]; stacki[2] = mai[2][1];
  stack[1] = ma[1][1]; stacki[2] = mai[1][1];
  stack[0] = ma[2][0]; stacki[3] = mai[2][0];
  _mult(); _rot(); _mult(); _rotup(); _sub();
  mc[2][0] = stack[0]; mci[2][0] = stacki[0];
  stack[3] = ma[0][1]; stacki[3] = mai[0][1]; // C21
  stack[2] = ma[2][0]; stacki[2] = mai[2][0];
  stack[1] = ma[0][0]; stacki[2] = mai[0][0];
  stack[0] = ma[2][1]; stacki[3] = mai[2][1];
  _mult(); _rot(); _mult(); _rotup(); _sub();
  mc[2][1] = stack[0]; mci[2][1] = stacki[0];
  stack[3] = ma[0][0]; stacki[3] = mai[0][0]; // C22
  stack[2] = ma[1][1]; stacki[2] = mai[1][1];
  stack[1] = ma[0][1]; stacki[2] = mai[0][1];
  stack[0] = ma[1][0]; stacki[3] = mai[1][0];
  _mult(); _rot(); _mult(); _rotup(); _sub();
  mc[2][2] = stack[0]; mci[2][2] = stacki[0];
  for (byte j = 0; j < 3; j++) for (byte i = 0; i < 3; i++) { // Aij = 1/det * Cij
      stack[0] = mdet; stacki[0] = mdeti;
      stack[1] = mc[j][i]; stacki[1] = mci[j][i];
      _div();
      ma[j][i] = stack[0]; mai[j][i] = stacki[0];
    }
}

Unfortunately this code "eats up" (too) many valuable bytes. I'm sure it can be optimized and shrinked - at least with some loops. But I must resign to manage such loop indices.

Does anyone of you (very smart guys) see how to automate this code in an efficient way?
Any idea and hint is highly welcome.
Regards deetee