Re: [GiNaC-list] Term ordering and compiling C++ code

14 May 2010

      I see we have similar problems. I'm using GiNaC to do Multibody
simulations, with some applications
in realtime control.

May you are interested on my work related to GiNaC

http://www.imac.unavarra.es/3d_mec/pages/lib3d_mec-ginac.php

Reference: J.Ros et al. LIB3D MEC–GINAC, A LIBRARY FOR SYMBOLIC
MULTIBODY
DYNAMICS. MULTIBODY DYNAMICS 2007, ECCOMAS Thematic Conference. Milano,
Italy. 2007.

The implementation of what I call atomization of expresions, is not
optimal at all,
and uses a completely different approach to the one commented here.

(Beware that the version in the we is not the last one)

Javier

On Thu, 2010-05-13 at 17:25 -0700, Doug wrote:
...
Javier,
The output code might be a tad simpler, but the code that generates
the code might be less simple, because now you need a special case for
the leaves/symbols.  Compilers are really good at detecting a = my_a;
and e1 = a;, so that wouldn't worry me at all.
In my case, the output code is already so complex that it doesn't much
matter.
I agree that a feature like this would make Ginac more attractive for
other people in my company who use similar symbolic code for computing
jacobians and hessians.  For the record, the code Ginac prints gets
evaluated at several hertz, like 10 or 20 Hz, so speed is critical.
Support NPR 20 seconds at a time. www.twentysecondsatatime.org
--- On Thu, 5/13/10, jros <jros@unavarra.es> wrote:
From: jros <jros@unavarra.es>
        Subject: Re: [GiNaC-list] Term ordering and compiling C++ code
        To: "GiNaC discussion list" <ginac-list@ginac.de>
        Date: Thursday, May 13, 2010, 7:13 PM
Yes, thats the idea. I suppose that can be implemented.
        An obvious simplification, is  to avoid declarations like
double e1 = a;
for symbols, as declaration
double a = my_a;
makes the symbol already available.
I would consider such a functionality really important for
        GiNaC. Lots of people is using GiNaC, to
        output code (thats my case, matlab/octave/ and c), and this
        would produce highly optimal
        code in terms of space and speed as lots of computations can
        be avoided.
It would be nice if some power developer could say something
        about the proper way to go, and
        impressions about possibilities to make it a part of GiNaC.
Thank you,
Javier
On Thu, 2010-05-13 at 11:47 -0700, Doug wrote:
> Javier, All,
        > 
        > I think what we need can be done with two print functions,
        > as described below.  It seems straightforward enough that I
        > can implement this myself, but I wonder if it already exists
        > somewhere?  Is this something the developers would consider
        > adding to Ginac?
        > 
        > Thanks,
        > Doug
        > 
        > ---------------------
        > 
        > Perhaps this is a clean solution.  One could write a print
        > function that outputs C++ code that explicitly names and
        > declares sub-expressions.  You also have a companion print
        > function that just emits the name of the top-level
        > expression.  Call them print_declarations() and
        > print_variable_name().  
        > 
        > For example, suppose you have an expression like this: 
        >      ex myex = a + b*pow(c+d+e,2.0) + f/(g*h + i*j);
        > For simplicity's sake, let each symbol just be a variable
        > (e.g. symbol a("my_a")), but of course they can be anything.
        > print_declarations() can emit these declarations below
        > during a post-order traversal of the expression tree.  
        > 
        > // Symbols
        > double a = my_a;
        > double b = my_b; 
        > [...]
        > double j = my_j;
        > // Expressions
        > double e1 = a;
        > double e2 = b;
        > double e3 = c;
        > double e4 = d;
        > double e5 = e;
        > double e6 = e3 + e4 + e5;
        > double e7 = pow(e6, 2.0);
        > double e8 = e2 * e7;
        > double e9 = f;
        > double e10 = g; 
        > double e11 = h;
        > double e12 = e10*e11;
        > double e13 = i;
        > double e14 = j;
        > double e15 = e13*e14;
        > double e16 = e12 + e15;
        > double e17 = e9 / e16;
        > double e18 = e1 + e8 + e17;
        > 
        > print_variable_name() can then just emit e18 for the example
        > expression, to be used however you would normally use C++
        > code that you print.
        > 
        > Then, assuming ginac is smart about not having separate
        > expression tree nodes for duplicated sub-expressions, you
        > will declare each sub-expression just once.  The order ginac
        > uses for expressions does not matter at all.
        > 
        > Support NPR 20 seconds at a time.
        > www.twentysecondsatatime.org
        > 
        > --- On Thu, 5/13/10, jros <jros@unavarra.es> wrote:
        > 
        >         
        >         From: jros <jros@unavarra.es>
        >         Subject: Re: [GiNaC-list] Term ordering and
        >         compiling C++ code
        >         To: "GiNaC discussion list" <ginac-list@ginac.de>
        >         Date: Thursday, May 13, 2010, 10:31 AM
        >         
        >         Hello,
        >         
        >         Although I don't have a solution for your problem,
        >         as I'm myself addressing similar problems
        >         matching common subexpressions to variables, in down
        >         top manner, I think that such a functionality
        >         is implicitly implemented in GiNaC.
        >         
        >         If I understand GiNaC internal structure correctly,
        >         subexpressions common to two expresions,
        >         are frequently  shared internally, to save memory.
        >         
        >         So it must be possible to write a print algorithm
        >         that goes trough an/some expression/s tree, and that
        >         replaces
        >         every shared subexpression (let say sum product)
        >         with a variable, that again is assigned a expression
        >         that would be printed
        >         in the same way using recursion. 
        >         
        >         Probably allowing/disallowing some kind automatic
        >         simplifications (so that subexpression sharing
        >         expected value increases) can probably help to
        >         obtain improved results.
        >         
        >         I wonder what do the developers think about this.
        >         
        >         Cheers,
        >         
        >         Javier
        >         
        >         
        >         On Thu, 2010-05-13 at 07:11 -0700, Doug wrote:
        >         
        >         > Hi,
        >         > 
        >         > I'm running into a term ordering problem like the
        >         > one mentioned in this post.
        >         > http://permalink.gmane.org/gmane.comp.mathematics.ginac.general/1107
        >         > 
        >         > Relevant details about my application are below,
        >         > for the curious.  The short of it, though, is that
        >         > the C++ code I am generating with ginac is so
        >         > large (50MB) that the compiler can't handle it,
        >         > and moreover, the values being represented by
        >         > ginac symbols are array accesses that the compiler
        >         > can't optimize away.  So I need to create
        >         > temporary values and do sub-expression elimination
        >         > on the C++ code just so the compiler can have a
        >         > chance.  I have successfully done this using
        >         > search-and-replace on the C++ code for other
        >         > applications, but I guess I was just lucky that
        >         > ginac printed terms in a consistent order.
        >         > 
        >         > Is this a problem other ginac users have
        >         > confronted?  From the answer to the above post, it
        >         > sounds like I need to write my own print function
        >         > that can ensure expressions are printed in a
        >         > predictable order.  It looks like I can also
        >         > inspect the expression tree and possibly do a
        >         > proper job of common sub-expression elimination.
        >         > But if there is an easier approach, I'd be very
        >         > grateful for a pointer or two.
        >         > 
        >         > Thanks,
        >         > Doug Baker
        >         > 
        >         > -------------------
        >         > Background Details
        >         > -------------------
        >         > 
        >         > My application is related to visual odometry.  I
        >         > am computing the gradient and hessian for a
        >         > function of 10 variables that involves converting
        >         > quaternions to rotation matrices, applying those
        >         > matrices to 3D points and projecting those points
        >         > onto a 2D plane.
        >         > 
        >         > I am generating code with ginac that is initially
        >         > 50MB of C++ code. (And thank you for ginac being
        >         > able to handle that monstrosity with such grace.)
        >         > The input symbols are all array elements (e.g.
        >         > symbol x("state[t+X]")) or function calls (e.g.
        >         > symbol rTi11("rTi.getValue<1,1>()")).
        >         > 
        >         > If I take the printed output from ginac as is, the
        >         > C++ compiler can't properly optimize it because it
        >         > doesn't know that state[t+X] is a constant or that
        >         > getValue is constant.  Moreover, the expressions
        >         > have huge repeated expressions, and these bog the
        >         > compiler down, too.
        >         > 
        >         > I've been able to use a simple search-and-replace
        >         > approach in other problems, and on this one that
        >         > same approach has gotten me from 50MB to 6MB or
        >         > so, but the term ordering issue is biting me now.
        >         > For example, I have one sum (a+b+c+d) that occurs
        >         > 232608 times.  But there are 4! = 24 permutations
        >         > of that sum that I'd have to replace.  It's worse
        >         > than that, though, because that's not the only
        >         > expression like it.  Some are bigger, like this:
        >         > (-2.0*sK/sKSIJ+4.0*sK*sS2/pKSIJ2
        >         > +-4.0*sI*sJ*sS/pKSIJ2). This has something like 72
        >         > permutations, not even counting the permutations
        >         > of pKSIJ and pKSIJ2, which are themselves
        >         > expressions like (a+b+c+d).  And not only are
        >         > there variations on this due to term ordering,
        >         > similar expressions occur with different variables
        >         > (e.g. with sJ instead of sK, etc.)  So when all
        >         > these expressions can be output in various orders,
        >         > it's no longer feasible to generate all the
        >         > permutations for search-and-replace.
        >         > 
        >         > 
        >         > Support NPR 20 seconds at a time.
        >         > www.twentysecondsatatime.org
        >         > 
        >         > 
        >         > _______________________________________________
        >         > GiNaC-list mailing list
        >         > GiNaC-list@ginac.de
        >         > https://www.cebix.net/mailman/listinfo/ginac-list
        >         
        >         
        >         -----Inline Attachment Follows-----
        >         
        >         _______________________________________________
        >         GiNaC-list mailing list
        >         GiNaC-list@ginac.de
        >         https://www.cebix.net/mailman/listinfo/ginac-list
        >         
        >
-----Inline Attachment Follows-----
_______________________________________________
        GiNaC-list mailing list
        GiNaC-list@ginac.de
        https://www.cebix.net/mailman/listinfo/ginac-list
_______________________________________________
GiNaC-list mailing list
GiNaC-list@ginac.de
https://www.cebix.net/mailman/listinfo/ginac-list