I'm just starting to use Maxima in order to simplify some big equations
and generate optimised code. I'm going to use C++ ultimately but I think
the fortran output converted automatically to C will be suitable.
Anyway, when I OPTIMIZE some code I noticed that it isn't fully
optimized, and I'm having problems generating the code after I optimize.
I'd be grateful if you could answer either of my two questions:
(i) How do you generate fortran from a BLOCK expression?
(ii) Can I encourage maxima to identify more common sub-expressions, or
would Mathematica or something else do a better job?
As a small example, here's one of my matrices unoptimized and converted
to fortran:
Y(1,1) = p(30,1)**2+p(30,0)**2+(-1.0)/2.0
Y(1,2) = p(30,0)*p(30,3)+p(30,1)*p(30,2)
Y(1,3) = p(30,1)*p(30,3)-p(30,0)*p(30,2)
Y(1,4) = 0
Y(1,5) = 0
Y(1,6) = 0
Y(2,1) = p(30,1)*p(30,2)-p(30,0)*p(30,3)
Y(2,2) = p(30,2)**2+p(30,0)**2+(-1.0)/2.0
Y(2,3) = p(30,2)*p(30,3)+p(30,0)*p(30,1)
Y(2,4) = 0
Y(2,5) = 0
Y(2,6) = 0
Y(3,1) = p(30,1)*p(30,3)+p(30,0)*p(30,2)
Y(3,2) = p(30,2)*p(30,3)-p(30,0)*p(30,1)
Y(3,3) = p(30,3)**2+p(30,0)**2+(-1.0)/2.0
Y(3,4) = 0
Y(3,5) = 0
Y(3,6) = 0
Y(4,1) = 0
Y(4,2) = -0.06*(p(30,1)*p(30,3)-p(30,0)*p(30,2))
Y(4,3) = 0.06*(p(30,0)*p(30,3)+p(30,1)*p(30,2))
Y(4,4) = p(30,1)**2+p(30,0)**2+(-1.0)/2.0
Y(4,5) = p(30,0)*p(30,3)+p(30,1)*p(30,2)
Y(4,6) = p(30,1)*p(30,3)-p(30,0)*p(30,2)
Y(5,1) = 0
Y(5,2) = -0.06*(p(30,2)*p(30,3)+p(30,0)*p(30,1))
Y(5,3) = 0.06*(p(30,2)**2+p(30,0)**2+(-1.0)/2.0)
Y(5,4) = p(30,1)*p(30,2)-p(30,0)*p(30,3)
Y(5,5) = p(30,2)**2+p(30,0)**2+(-1.0)/2.0
Y(5,6) = p(30,2)*p(30,3)+p(30,0)*p(30,1)
Y(6,1) = 0
Y(6,2) = -0.06*(p(30,3)**2+p(30,0)**2+(-1.0)/2.0)
Y(6,3) = 0.06*(p(30,2)*p(30,3)-p(30,0)*p(30,1))
Y(6,4) = p(30,1)*p(30,3)+p(30,0)*p(30,2)
Y(6,5) = p(30,2)*p(30,3)-p(30,0)*p(30,1)
Y(6,6) = p(30,3)**2+p(30,0)**2+(-1.0)/2.0
When I optimized this matrix I get some reductions:
BLOCK([%1, %2, %3, %4, %5, %6, %7, %8, %9, %10, %11, %12, %13],
2 2 1
%1 : p , %2 : p + %1 - -, %3 : p p ,
30, 0 30, 1 2 30, 1 30, 2
%4 : p p + %3, %5 : p p , %6 : %5 - p p ,
30, 0 30, 3 30, 1 30, 3 30, 0 30, 2
2 1
%7 : %3 - p p , %8 : p + %1 - -, %9 : p p ,
30, 0 30, 3 30, 2 2 30, 2 30, 3
%10 : %9 + p p , %11 : %5 + p p , %12 : %9 - p p ,
30, 0 30, 1 30, 0 30, 2 30, 0 30, 1
[ %2 %4 %6 0 0 0 ]
[ ]
[ %7 %8 %10 0 0 0 ]
[ ]
2 1 [ %11 %12 %13 0 0 0 ]
%13 : p + %1 - -, [ ])
30, 3 2 [ 0 - 0.06 %6 0.06 %4 %2 %4 %6 ]
[ ]
[ 0 - 0.06 %10 0.06 %8 %7 %8 %10 ]
[ ]
[ 0 - 0.06 %13 0.06 %12 %11 %12 %13 ]
Sorry that's a bit of a mess... anyway, there are some common sub
expressions identified.
However, when I try to convert this to fortran using fortran( M ); where
I previously assigned M:optimize(ev(sybolic expression for the matrix));
I get:
(C205) fortran(M);
BLOCK([%1,%2,%3,%4,%5,%6,%7,%8,%9,%10,%11,%12,%13],%1:p(30,0)**2,%
1 2:p(30,1)**2+%1+(-1.0)/2.0,%3:p(30,1)*p(30,2),%4:p(30,0)*p(30,3
2 )+%3,%5:p(30,1)*p(30,3),%6:%5-p(30,0)*p(30,2),%7:%3-p(30,0)*p(3
3 0,3),%8:p(30,2)**2+%1+(-1.0)/2.0,%9:p(30,2)*p(30,3),%10:%9+p(30
4 ,0)*p(30,1),%11:%5+p(30,0)*p(30,2),%12:%9-p(30,0)*p(30,1),%13:p
5 (30,3)**2+%1+(-1.0)/2.0,MATRIX([%2,%4,%6,0,0,0],[%7,%8,%10,0,0,
6 0],[%11,%12,%13,0,0,0],[0,-0.06*%6,0.06*%4,%2,%4,%6],[0,-0.06*%
7 10,0.06*%8,%7,%8,%10],[0,-0.06*%13,0.06*%12,%11,%12,%13]))
which isn't fortran!
The other problem is that not all CSEs are identified, e.g., %10:%9 +
p[30,0]*p[30,1] and %12:%9-p[30,0]*p[30,1] both have p[30,0]*p[30,1] as
a common term.
Thanks for any help!
Best,
Josh.