July 2003 - GiNaC-list

Re: C++ code generation: any way of deflating large expressions?
by hugo_lec＠club-internet.fr 21 Jul '03

21 Jul '03

Hello, I also worked on such a program : I have to convert huge expressions in optimized C++ code (to automate design of finite element code). I tried to make something which does not calculate two times the same thing. There's an algorithm to use the less temporary variables as possible (in order to help the compiler to optimize). It also looks for factorizations with * and + (A=x*y; B=sin(x*y*z) -> B=sin(A*z)... ) and for the cost of operations (A=x/y; B=z/y; gives double tmp0=1/y; A=x*tmp0; B=z*tmp0)... It can generate code for huge expressions in a few seconds. It is a bit crappy because it does not use the internal representation given by GiNaC. But if someone wants to clean it... For now it is usable enough for me. A simple example : Optimized_ex_output mm("float"); symbol x ("x"), y ("y"), z ("z"), u("u"), v("v"); mm.add ("A{1]", x/u*v, true); mm.add ("B", (x+y)/u*v+cos((x+y+z)*y), false); mm.write(std::cout,2); // 2 spaces before each line mm.display_tree(); // graph output produces : float reg0 = v; float reg1 = z; float reg2 = u; float reg3 = x; float reg4 = y; float reg5 = reg3; reg3 += reg4; reg2 = 1 / reg2; float reg6 = reg2; reg2 *= reg3; reg2 *= reg0; reg0 *= reg6; reg3 += reg1; reg4 *= reg3; reg0 *= reg5; A{1] += reg0; reg4 = cos( reg4 ); reg4 += reg2; float B = reg4; (and a nice dot graph) In the near future, I've planed to allow generation of vectorial code to produce 4x the same calculi on different data. I've also planed to generate code for GPUs but not for the near future. Hugo LECLERC.

3 2

Changing symbol string names on-the-fly
by Dr. Vassilis S. Vassiliadis 11 Jul '03

11 Jul '03

Is it possible to change on-the-fly the string names used in a symbol? For example, I am considering programs with several thousands of variables and wish to output the code for the expressions in csrc-format. The names used internally utilize the object-partitioning of my code e.g. : symbol x("system[1].composition.x[2]") The output in raw C++ form is desired to be a simple counting scheme e.g. for the above symbol: "var[101]" will be substituted. I have tried to use an auxiliary vector list with new symbols and then substitute the new symbols for the old, using .subs(). However, when the number of substitutions is beyond say a couple of hundred even a 2.5 GHz machine cannot run the problem efficiently. I proceeded to test further abstract examples of substitutions and found that the .subs() method can indeed be extremely slow (inefficient?) in general. Without resorting to brute force outputing and replacing of strings, I wonder if there is an elegant and economical way to do the above in GiNaC. Vassilis. ----------------------------------------------------------------- Dr. Vassilis S. Vassiliadis, Lecturer, UNIVERSITY OF CAMBRIDGE, DEPARTMENT OF CHEMICAL ENGINEERING Pembroke Street, Cambridge CB2 3RA, UK. Tel: (reception): +44 1223 334777 Fax: (departmental): +44 1223 334796 WWW: www.cheng.cam.ac.uk/~vsv20 e-mail: vsv20(a)cheng.cam.ac.uk -----------------------------------------------------------------

3 5

how to subs t.i with i*t.i
by Yong Xiao 11 Jul '03

11 Jul '03

Hi, it seems that the following lines ex e1 = indexed(t,i); e1.subs(indexed(t,idx(i_,3))==idx(i_,3).get_value()*indexed(t,idx(i_,3))); ,where i and i_ are idx and wildcard objects respectively, doesn't substitute t.i with i*t.i. t.i still remains it's original form. Is there any way I can do such subs: t.i-->i*t.i? Thank you. Yong

2 2

The speed of subs()
by Christian Bauer 10 Jul '03

10 Jul '03

Hi! Prompted by the recent discussion I've done some timings and experiments with subs(). Two benchmarks were used: a) "Binary tree" constructs a binary expression tree (alternating sums and products) of depth N with symbols as the leaves (e.g. N=3 leads to (a1+a2)*(a3+a4)+(a5+a6)*(a7+a8)), then substitutes all 2^N symbols by a set of new symbols (using one subs() call). b) "Linear polynomial" constructs an expression of the form "1+a1+a1*a2+ a1*a2*a3+...+a1*...*aN" of length N and also substitutes all N symbols by new ones. The time measured was only the time for the subs() step, excluding the expression and list set-up time. All tests were done on an Athlon XP 2000+ machine with 768MB of RAM. Results: GiNaC 1.0 --------- With default options: Binary tree, 2^N symbols.... N: 8 9 10 11 time/s: 0.17 1.33 13.269 129.4 Linear polynomial, N symbols...... N: 100 200 300 400 500 600 time/s: 0.13 1.39 6.32 18.539 49.429 109.64 With no_pattern = true: Binary tree, 2^N symbols.... N: 8 9 10 11 time/s: 0.17 1.31 12.29 129.629 Linear polynomial, N symbols...... N: 100 200 300 400 500 600 time/s: 0.08 1.03 4.879 15.429 41.77 99.43 GiNaC 1.0 wastes its time in the inefficient lst::op(). no_pattern isn't of much help here. GiNaC 1.1 --------- With default options: Binary tree, 2^N symbols..... N: 8 9 10 11 12 time/s: 0.01 0.05 0.17 0.739 5.06 Linear polynomial, N symbols...... N: 100 200 300 400 500 600 time/s: 0.05 0.4 1.34 3.149 6.15 10.92 With subs_options::subs_no_pattern: Binary tree, 2^N symbols..... N: 8 9 10 11 12 time/s: 0 0 0.029 0.17 2.25 Linear polynomial, N symbols...... N: 100 200 300 400 500 600 time/s: 0 0.04 0.1 0.25 0.479 0.81 Much better. subs_no_pattern significantly speeds up the second case. GiNaC 1.2 (i.e. what's currently in CVS) --------- Performs about 5% better than GiNaC 1.1 but shows the same scaling behavior, so I've omitted the results. GiNaC 1.2 with basic::subs() using a map<ex,ex,ex_is_less> instead of two lists --------- Default options: Binary tree, 2^N symbols..... N: 8 9 10 11 12 time/s: 0.01 0.029 0.13 0.55 3.899 Linear polynomial, N symbols...... N: 100 200 300 400 500 600 time/s: 0.059 0.39 1.29 3.12 5.919 10.46 With subs_options::subs_no_pattern: Binary tree, 2^N symbols..... N: 8 9 10 11 12 time/s: 0 0.01 0.029 0.13 1.53 Linear polynomial, N symbols...... N: 100 200 300 400 500 600 time/s: 0 0.02 0.04 0.059 0.1 0.149 A significant speedup, especially with subs_no_pattern. In the general case, subs() still has to iterate through the list of substitutions to apply match() to every single one. Only with subs_no_pattern can it use map::find(). But wait, there's more: I've noticed that expairseq::subschildren() scans the substitution list and uses a different algorithm depending on whether the list contains products or not. This step should be done only once, and the result passed on as a subs_options flag. Preliminary tests indicate that doing so speeds up the "binary tree" test by two orders of magnitude, and the "linear polynomial" test by 10%). Should these changes go into GiNaC 1.2? basic::normal() should then probably also be altered to take a map<> instead of two lists. What about match(), to_rational() and to_polynomial()? Also, this will break the API for class implementors (again), and Stefan Weinzierl will beat me. :-) (ok, it's broken in 1.2 anyway since basic::copy() and basic::destroy() are gone) Bye, Christian -- / Physics is an algorithm \/ http://www.uni-mainz.de/~bauec002/

1 0

C++ code generation: any way of deflating large expressions?
by Dr. Vassilis S. Vassiliadis 10 Jul '03

10 Jul '03

I wonder if someone has run into the following and has found a easy way out: In generating C-source code from symbolic expressions it is always possible to run into trouble with large size expression output. For example, if one stores in an ex variable something like: ex SUM; `` SUM = \sum_{i=1}^{10000} x[i] ; '' then printing this out as csrc will generate a huge line. The best way is to see if there is an easy algorithmic construct (perhaps within GiNaC, or hints provided for its efficient implementation) so that at each level of the tree all operation-containing nodes are assigned a dummy variable (e.g. those of arithmetic operations or implementing a function). Then upon exhaustion, to reverse the procedure from the leaves of the tree to the root and thus generate a C++ code fragment with all necessary steps & intermediate variables to finally get the value of the root. ----------------------------------------------------------------- Dr. Vassilis S. Vassiliadis, Lecturer, UNIVERSITY OF CAMBRIDGE, DEPARTMENT OF CHEMICAL ENGINEERING Pembroke Street, Cambridge CB2 3RA, UK. Tel: (reception): +44 1223 334777 Fax: (departmental): +44 1223 334796 e-mail: vsv20(a)cheng.cam.ac.uk -----------------------------------------------------------------

2 1

symbolic expressions: machine precision calculation
by Dr. Vassilis S. Vassiliadis 10 Jul '03

10 Jul '03

I wonder if there is a possibility to use within GiNaC the C++ native precision arithmetic rather than that of CLN. In particular, following symbolic generation and manipulation of a mathematical model, it would be desirable in many cases to switch to machine arithmetic for numerical computations based on the tree of each constructed expression of interest. Such a possibility would obviate the need to spool external C++ source code, to recompile externally from the application built on GiNaC in the first place, and would lead to usability of GiNaC in limitless domains where eventual efficient numerical computation is paramount. Mathematica, although in a rudimentary and rather cumbersome way, has this option for example (the so-called compiled version of a function). Vassilis. ----------------------------------------------------------------- Dr. Vassilis S. Vassiliadis, Lecturer, UNIVERSITY OF CAMBRIDGE, DEPARTMENT OF CHEMICAL ENGINEERING Pembroke Street, Cambridge CB2 3RA, UK. Tel: (reception): +44 1223 334777 Fax: (departmental): +44 1223 334796 e-mail: vsv20(a)cheng.cam.ac.uk -----------------------------------------------------------------

1 0

nestedsums for GiNaC 1.1
by Christian Bauer 03 Jul '03

03 Jul '03

Hi! Due to the popularity of Stefan Weinzierl's nestedsums library (see http://www.fis.unipr.it/~stefanw/nestedsums/) I have taken the liberty to create a version of the library that works with GiNaC 1.1. Sources and pre-packaged RPM files can be found in ftp://ftpthep.physik.uni-mainz.de/pub/nestedsums/ Changes to the original nestedsums 1.1.1 include: - adapted to GiNaC 1.1 API - testsuite can be run with "make check" - includes a nestedsums.m4 that provides an AM_PATH_NESTEDSUMS macro for autoconf scripts - includes a spec file for RPM generation Happy nesting, Christian -- / Physics is an algorithm \/ http://www.uni-mainz.de/~bauec002/

1 0

Multivariate polynomial factorization
by poulard 02 Jul '03

02 Jul '03

Dear all, I have a great challenge to manage and a deep thinking before beginning the developpement. First before entering into the real subject I wish tell you (if it has not been before) that GiNaC+CLN compile well on W2K under cygwin with gcc (gcc version 3.2 20020927 (prerelease)) All chech and tests are ok! Great job. NOW THE REAL SUBJECT OF THE MAIL -------------------------------- For one of our project WE DEFINITLY need a package for multivariate polynomial factorization under Z. For this I will have rapidely the man.power to begin this great challenge (next week). My great wish would be to work within open-source framework. As this subject is one of the ToDo list of GiNac, and a very difficult task. Before compiling GiNaC, I beleived that I could use it as a native library un Win32. But it seems that not, unless modification. Perhaps someone has did it ? I have a constraint that I can't modify : I MUST BE ABLE to use and link with a classic native win32 application. So, unless someone tell me the modification to do to be able to compile cln + ginac under VC6++ I cannot use ginac as a basis for further developpment. ----------------------------------------------------------- The solution for me today is the use of NTL (which compile well with VC6) as univariate basis. ----------------------------------------------------------- BUT the developpement of this package under NTL could serve as a basis for GiNac also ! So I need you help for this developpement at a information level. --------------------------------------------------------- What is today the best method for multivariate polynomial factorization under Z ? --------------------------------------------------------- I have the old algorithm from Berlekamp & Hensel for univariate (p-adique method) and only an extension for multivariate, but the algorithm complexity seems to be very high, and I guess that new algorithms have been found far now ? Someone could help me to indicate what is today the best algorithm for this task ? REMARK = The rapidity and accuracy of the result of such task made by Maple5 indicates that there should be a quit rapid algorithm ! For finish this - call to contribution - , I must be honest and indicate that the polynom I have to factorize are little special : - they are very very very long in term of monomial : the number of monomial could be 100, 1000 !!! - they could be a lot of variable (when I say a lot, it could be 10 and more) !!!! - THEY ARE homogeneous : the total degree of monomial are always the same - The coefficient of each monomial is : 1, 2 or not very high - IT SEEMS, after lot experiment under Maple5 that they are squarre free (BUT we have not demonstrate this assertion) UNDER this condition I can tell you that the kernel of Maple5 (in Matlab), for a HUGE HUGE polynomial like I describe, result is almost immediate !!!!! Very very speed. I hope that it will interest you to collaborate and help me to develop this piece of software under open-source. Herve PS = If you the modification to do to GiNaC to be able to be a native win32 library is not to complex I promiss to collaborate to do this task. Thank you, if you have read all :-) _/_/_/_/_/_/_/_/_/_/_/_/_/_/_/_/_/_/_/_/_/_/_/_/_/_/_/_/_/_/_/_/ _/ Dr Herve Poulard | ACTIA - ACTIELEC Hold. _/ _/ Research Manager | BP 4215 _/ _/ Tel : (33) 5.61.17.68.79 | 25 Chemin de Pouvourville _/ _/ Fax : (33) 5.61.55.42.31 | 31432 Toulouse Cedex 04 _/ _/ | FRANCE _/ _/ E-Mail | _/ _/ - poulard(a)actia.fr | http://www.actia.com/ _/ _/ - poulard(a)laas.fr | http://www.laas.fr/~poulard/ _/ _/ - poulard.herve(a)wanadoo.fr | http://www.actielec.com/ _/ _/ (plz cc for security) | _/ _/_/_/_/_/_/_/_/_/_/_/_/_/_/_/_/_/_/_/_/_/_/_/_/_/_/_/_/_/_/_/_/

2 1