- GiNaC-list - ginac.de

Indefinite integration support
by Jorge Cardona 16 Jul '09

16 Jul '09

Hi, I recently wrote a code that will let define a entity that represents a indefinite integral, i called iintegral, basically i took the integral class, and remove the integration range, then i add to the function.pl extra code to define a new option on function, called iintegral_func, with that is possible to define the indefinite integral of a function when we define the function itself. I want to know if you like that idea and if you have some extra comment on this part. Look a definition of a function on inifcns_trans : static ex log_iinteg(const ex & x, unsigned iinteg_param) { GINAC_ASSERT(iinteg_param==0); // int {log(x)} -> x * (log (x) - 1) return x * (log (x) - _ex1); } .... REGISTER_FUNCTION(log, eval_func(log_eval). evalf_func(log_evalf). derivative_func(log_deriv). iintegral_func(log_iinteg). series_func(log_series). real_part_func(log_real_part). imag_part_func(log_imag_part). latex_name("\\ln")); ... Right now i call the iintegral constructor with 2 arguments, the integrand expression and the integration variable (the integration variable remains after the indefinite integral, is not as dummy as in the definite integral), and this iintegral instance remains until the eval_integ is called, which in a future (near future i hope) will have a better algorithm to integrate the expression, but right now just have the code that integrate a polynomial expression and a code that check if the expression is a functions, and if that function has an integral defined it also check if the argument is equal to the integration variable, and then call that function, if any of this don't match it return the iintegral object. I want to know what do you think about that way, or if you have some extra ideas, basically i want to implement the integration algorithm, but don't want to have some code that wont be added to GiNaC, so, please let me know your comments. Code on function::iintegral : ex function::iintegral(const symbol & s) const { GINAC_ASSERT(serial<registered_functions().size()); const function_options &opt = registered_functions()[serial]; // No iintegral defined? Then return abstract iintegral object if (opt.iintegral_f == NULL) return GiNaC::iintegral(s,this->eval()); size_t iinteg_param = -1; size_t num = seq.size(); for (size_t i=0; i<num; i++) { ex a = s; ex b = seq[i]; if(a.is_equal(b)) { iinteg_param = i; break; } } current_serial = serial; if (opt.iintegral_use_exvector_args) return ((iintegral_funcp_exvector)(opt.iintegral_f))(seq, iinteg_param); if(iinteg_param < 0) throw(std::logic_error("Composed functions on integral.")); switch (opt.nparams) { // the following lines have been generated for max. 14 parameters case 1: return ((iintegral_funcp_1)(opt.iintegral_f))(seq[1-1],iinteg_param); case 2: return ((iintegral_funcp_2)(opt.iintegral_f))(seq[1-1], seq[2-1],iinteg_param); case 3: return ((iintegral_funcp_3)(opt.iintegral_f))(seq[1-1], seq[2-1], seq[3-1],iinteg_param); case 4: return ((iintegral_funcp_4)(opt.iintegral_f))(seq[1-1], seq[2-1], seq[3-1], seq[4-1],iinteg_param); ...(Go with more arguments) Code on iintegral::eval_integ(): ex iintegral::eval_integ() const { if (!(flags & status_flags::expanded)) return this->expand().eval_integ(); if (f==x) return x*x/2; if (is_a<power>(f) && f.op(0)==x) { if (f.op(1)==-1) return log(x); if (!f.op(1).has(x)) { ex primit = power(x,f.op(1)+1)/(f.op(1)+1); return primit; } } if (is_a<function>(f)) { function fun = ex_to<function>(f); // TODO: Change the signature. return fun.iintegral(ex_to<symbol>(x)); } return *this; } There is a patch at: http://wiki.scilab.org/Gsoc2009Symbolic?action=AttachFile&do=get&target=Ext… This is a little code test: #include "ginac/ginac.h" using namespace std; using namespace GiNaC; int main() { symbol x("x"); ex f,g,h; f = 2*x; cout << f << endl; g = iintegral(x,f); cout << g << endl; h = g.eval_integ(); cout << h << endl; f = 2*sin(x) + exp(x) + log(x); cout << f << endl; g = iintegral(x,f); cout << g << endl; h = g.eval_integ(); cout << h << endl; } That returns: jcardona@terminus:~/stuff/personal/symbolic/spikes/Integrals$ ./test 2*x iintegral(x,2*x) x^2 exp(x)+log(x)+2*sin(x) iintegral(x,exp(x)+log(x)+2*sin(x)) exp(x)-2*cos(x)+(-1+log(x))*x Bye. -- Jorge Eduardo Cardona jorgeecardona(a)gmail.com jorgeecardona.blogspot.com ------------------------------------------------ Linux registered user #391186 Registered machine #291871 ------------------------------------------------

2 1

Release GiNaC 1.5.2
by Jens Vollinga 11 Jul '09

11 Jul '09

Hi everybody, GiNaC 1.5.2 is out and available. The changes are: * Fixed documentation of return_type_info(). * Fixed cstdio include issue. * GiNaC now requires at least version 1.2.2 of cln. As always, this release can be downloaded from http://www.ginac.de/Download.html Best wishes, Jens _______________________________________________ GiNaC-list mailing list GiNaC-list(a)ginac.de https://www.cebix.net/mailman/listinfo/ginac-list

1 0

legendre and chebyshev polynomials
by Jorge Cardona 10 Jul '09

10 Jul '09

Hi, i'm working in three classes to generate the legendre polynomial, and the chebyshev of first and second kind, in chebyshev i define the derivative in terms of chebyshev polynomials, but it get messy when is expanded because the polynomials at numerator and denominator are not factorized, and then the expression is complex that it should. What do you think is better, to let the derivatives of the polynomials in terms on the same polynomials, and expaned as like the users, or always expand at the derivative and return the polynomial expanded at the derivative (as in the legendre polynomial)? This is the code: / Chebyshev polynomial of first kind class chebyshev_1 : public basic { GINAC_DECLARE_REGISTERED_CLASS(chebyshev_1, basic) public: // Constructor of chebyshev of first kind chebyshev_1(numeric n, symbol x); ex derivative(const symbol & s) const; ex expand(unsigned options) const; protected: void do_print(const print_context & c, unsigned level = 0) const; private: numeric n; symbol x; }; // Chebyshev polynomial of second kind class chebyshev_2 : public basic { GINAC_DECLARE_REGISTERED_CLASS(chebyshev_2, basic) public: // Constructor of chebyshev of first kind chebyshev_2(numeric n, symbol x); ex derivative(const symbol & s) const; ex expand(unsigned options) const; protected: void do_print(const print_context & c, unsigned level = 0) const; private: numeric n; symbol x; }; // Legendre polynomial class legendre : public basic { GINAC_DECLARE_REGISTERED_CLASS(legendre, basic) public: // Constructor of chebyshev of first kind legendre(numeric n, symbol x); ex derivative(const symbol & s) const; ex expand(unsigned options) const; protected: void do_print(const print_context & c, unsigned level = 0) const; private: numeric n; symbol x; }; /*************************************************************** * Chebyshev polynomial of first kind ****************************************************************/ GINAC_IMPLEMENT_REGISTERED_CLASS_OPT(chebyshev_1, basic, print_func<print_context>(&chebyshev_1::do_print)) // Default constructor chebyshev_1::chebyshev_1() {} // Print the function void chebyshev_1::do_print(const print_context & c, unsigned level) const { c.s << "Tn(" << n << "," << x << ")"; } int chebyshev_1::compare_same_type(const basic & other) const { const chebyshev_1 &o = static_cast<const chebyshev_1 &>(other); if (x.is_equal(o.x) and (n == o.n)) return 0; return 1; } ex chebyshev_1::expand(unsigned options) const { if (n == 0) return 1; if (n == 1) return x; exvector e; numeric m = n / 2; if (n.is_odd()) m = (n-1)/2; e.reserve(m.to_long()); ex tmp = pow(2 * x,n) / 2; e.push_back(tmp); for(numeric i = 1; i <= m; i++) { tmp *= -(n - 2*i + 2) * (n - 2*i + 1); tmp /= (4 * i * pow(x,2) * (n-i)); e.push_back(tmp); } return add(e); } ex chebyshev_1::derivative(const symbol & s) const { if(s.is_equal(x) and (n!=0)) return n * chebyshev_2(n-1,x); return 0; } chebyshev_1::chebyshev_1(numeric n, symbol x) : n(n), x(x){} /*************************************************************** * Chebyshev polynomial of second kind ****************************************************************/ GINAC_IMPLEMENT_REGISTERED_CLASS_OPT(chebyshev_2, basic, print_func<print_context>(&chebyshev_2::do_print)) // Default constructor chebyshev_2::chebyshev_2() {} // Print the function void chebyshev_2::do_print(const print_context & c, unsigned level) const { c.s << "Un(" << n << "," << x << ")"; } int chebyshev_2::compare_same_type(const basic & other) const { const chebyshev_2 &o = static_cast<const chebyshev_2 &>(other); if (x.is_equal(o.x) and (n == o.n)) return 0; return 1; } ex chebyshev_2::expand(unsigned options) const { if (n == 0) return 1; if (n == 1) return 2*x; exvector e; numeric m = n / 2; if (n.is_odd()) m = (n-1)/2; e.reserve(m.to_long()); ex tmp = pow(2 * x,n); e.push_back(tmp); for(numeric i = 1; i <= m; i++) { tmp *= -(n - 2*i + 2) * (n - 2*i + 1); tmp /= (4 * i * pow(x,2) * (n - i + 1)); e.push_back(tmp); } return add(e); } ex chebyshev_2::derivative(const symbol & s) const { if(s.is_equal(x) and (n!=0)) return ((n+1) * chebyshev_1(n+1,x) - x * chebyshev_2(n,x)) / (pow(x,2) - 1); return 0; } chebyshev_2::chebyshev_2(numeric n, symbol x) : n(n), x(x){} /*************************************************************** * Legendre polynomial ****************************************************************/ GINAC_IMPLEMENT_REGISTERED_CLASS_OPT(legendre, basic, print_func<print_context>(&legendre::do_print)) // Default constructor legendre::legendre() {} // Print the function void legendre::do_print(const print_context & c, unsigned level) const { c.s << "Pn(" << n << "," << x << ")"; } int legendre::compare_same_type(const basic & other) const { const legendre &o = static_cast<const legendre &>(other); if (x.is_equal(o.x) and (n == o.n)) return 0; return 1; } ex legendre::expand(unsigned options) const { if (n == 0) return 1; if (n == 1) return x; exvector e; numeric m = n / 2; if (n.is_odd()) m = (n-1)/2; e.reserve(m.to_long()); ex tmp = factorial(2*n) * pow(x,n) /(pow(2,n) * pow(factorial(n),2)); e.push_back(tmp); for(numeric i = 1; i <= m; i++) { tmp *= -(n - 2 * i + 1) * (n - 2 * i + 2); tmp /= (4 * n -4 * i + 2) * i * pow(x,2); e.push_back(tmp); } return add(e); } ex legendre::derivative(const symbol & s) const { return (this->expand(0)).diff(s); //if(s.is_equal(x) and (n!=0)) // return ((n+1) * chebyshev_1(n+1,x) - x * legendre(n,x)) / (pow(x,2) - 1); //return 0; } legendre::legendre(numeric n, symbol x) : n(n), x(x){} -- Jorge Eduardo Cardona jorgeecardona(a)gmail.com jorgeecardona.blogspot.com ------------------------------------------------ Linux registered user #391186 Registered machine #291871 ------------------------------------------------

2 2

GiNaC and Scilab
by Jorge Cardona 10 Jul '09

10 Jul '09

Hi, My name is Jorge Eduardo Cardona, i'm working rigth now on GSoC with Scilab, my project consist in add a symbolic module to it. In order to do that i had choose GiNaC library for it, and your ToDo list has two topics in which i want to work: Limits and symbolic integration. Well, about integrals i check the source and found a little piece of code that expands integrals (integrals.cpp) and resolve integrals on pow types, and with it solve polynomial integrals. I was thinking in start there to add more methods, but i want to know if there is some work already done. One way could be to add an argument to REGISTER_FUNCTION and pass an integral_func(...) as the extra argument, and we have then integrals for all the functions registered, and we can then start to use more heuristic methods on integrals.cpp, well thats just an initial idea, i want to know all the previous ideas on it to start to work. With limits, i don't even look at the code for it, but this night i will check it. Well, please let me know what do you think about it, you can see a wiki with information in: http://wiki.scilab.org/Gsoc2009Symbolic, and the code on: http://gitweb.scilab.org/?p=symbolic;a=summary. I'm reading some papers on ACM Digital Library about symbolic integration, specially Symbolic integration: the stormy decade, and the one that is linked at the TODO list(that link is actually broken: http://www-sop.inria.fr/cafe/Manuel.Bronstein/publications/issac98.pdf) Well, i really want to add these features to GiNaC and then use it on scilab, please let me know everything about it, what do you expect of it? what basic guidelines on code writing, basic ideas, anything. Bye. -- Jorge Eduardo Cardona jorgeecardona(a)gmail.com jorgeecardona.blogspot.com ------------------------------------------------ Linux registered user #391186 Registered machine #291871 ------------------------------------------------

2 3

Non commutative symbols ?
by ahabrard 30 Jun '09

30 Jun '09

Hello, We are new to the Ginac community and we want to manipulate non commutative symbols. We didn't find how to declare non commutative symbols in the tutorial and we found this on the web: const unsigned tinfo=2; symbol x("x", return_types::noncommutative,&tinfo) But we don't know if this declaration is correct ? We will appreciate any documentation on using non commutative symbols. We tried the following code: const unsigned tinfo=2; symbol x("x", return_types::noncommutative,&tinfo); symbol y("y",return_types::noncommutative,&tinfo); cout << x*y+y*x+x*y << endl; cout << x/x << endl; cout << x*pow(x,-1) << endl; The first expression works fine, however the last two expressions are not simplified to the neutral element. Is there anyone who knows the reason ? Another problem with matrix of non commutative symbols, consider the following code: matrix m(2,2); m = x, y, y, x; cout << m.inverse() << endl; Ginac stops with the following message: "Terminate called after throwing an instance of 'std::logic_error' what(): add::eval(): sum of non-commutative objects has non-zero numeric term Abort trap" Anyone knows the reason ? Thanks in advance for your help. Best regards, Amaury -- View this message in context: http://www.nabble.com/Non-commutative-symbols---tp24257628p24257628.html Sent from the Ginac - General mailing list archive at Nabble.com.

3 2

parse own function
by Jan.M 14 Jun '09

14 Jun '09

hello, i'm new to ginac and got some problems with parsing own functions: #include <iostream> #include <ginac/ginac.h> using namespace std; using namespace GiNaC; DECLARE_FUNCTION_3P(myfct) REGISTER_FUNCTION(myfct,dummy()) int main(int argc, char **argv){ symbol x("x"),y("y"); lst all; all = x,y; string str = "myfct(x,y,2)"; ex e(str, all); cout << e << endl; return 0; } ---> terminate called after throwing an instance of 'GiNaC::parse_error' what(): GiNaC: parse error at line 0, column 0: no function "myfct" with 3 arguments [GiNaC::ex GiNaC::parser::parse_identifier_expr()(parser/parser.cpp:66)] have i forgot something? thanks in advance. jan -- View this message in context: http://www.nabble.com/parse-own-function-tp24001248p24001248.html Sent from the Ginac - General mailing list archive at Nabble.com.

2 1

Re: [GiNaC-list] GiNaC::marix::solve(..)
by Lisa Maletzki 04 Jun '09

04 Jun '09

> Message: 1 > Date: Tue, 02 Jun 2009 18:37:08 +0200 > From: Lisa Maletzki <l.maletzki(a)tu-bs.de> > Subject: [GiNaC-list] GiNaC::marix::solve(..) problem > To: ginac-list(a)ginac.de > Message-ID: <20090602183708.y2zyinropw4o840k(a)webmail.tu-bs.de> > Content-Type: text/plain; charset=ISO-8859-1; DelSp="Yes"; > format="flowed" > > Hey, > > I have a smallish problem with the method mentioned above > (matrix::solve(..)). I have two nxn numeric matrices A and b (from the > equation Ax=b) and of course want the matrix x in the end but I keep > on getting the error that the first argument must be a matrix of > symbols. Now my question how can I do that? I tried it with ex > symbolic_matrix(..) but I couldn't get it to work since it is of the > type ex and not matrix. > > Kind regards, > > Lisa > > > > ------------------------------ > > Message: 2 > Date: Tue, 02 Jun 2009 19:35:08 +0200 > From: Jens Vollinga <jensv(a)nikhef.nl> > Subject: Re: [GiNaC-list] GiNaC::marix::solve(..) problem > To: GiNaC discussion list <ginac-list(a)ginac.de> > Message-ID: <4A2562CC.2020405(a)nikhef.nl> > Content-Type: text/plain; charset=ISO-8859-1; format=flowed > > Hi, > > Lisa Maletzki schrieb: >> I have a smallish problem with the method mentioned above >> (matrix::solve(..)). I have two nxn numeric matrices A and b (from the >> equation Ax=b) and of course want the matrix x in the end but I keep on >> getting the error that the first argument must be a matrix of symbols. >> Now my question how can I do that? I tried it with ex >> symbolic_matrix(..) but I couldn't get it to work since it is of the >> type ex and not matrix. > > take the following example as a guide: > > matrix A(2,2); > A = 1, 2, 3, 4; > symbol x1("x1"), x2("x2"); > matrix x(2,1); > x = x1, x2; > matrix b(2,1); > b = 5, 6; > cout << A.solve(x, b) << endl; > > As the manual warns, you will get an exception thrown if the system is > overdetermined. A more robust version is shown in the following example > using no matrices: > > symbol x1("x1"), x2("x2"); > lst equations; > equations.append(1*x1 + 2*x2 == 5); > equations.append(3*x1 + 4*x2 == 6); > cout << lsolve(equations, lst(x1,x2)) << endl; > > Regards, > Jens > Thank you Jens. But the problem I have is, that I do not know how big my matrix will be. In the smallest case it is a 18x18 matrix, so getting all those symbols would be rather complicated, or am I just thinking that? Kind regards, Lisa

3 6

GiNaC::marix::solve(..) problem
by Lisa Maletzki 02 Jun '09

02 Jun '09

Hey, I have a smallish problem with the method mentioned above (matrix::solve(..)). I have two nxn numeric matrices A and b (from the equation Ax=b) and of course want the matrix x in the end but I keep on getting the error that the first argument must be a matrix of symbols. Now my question how can I do that? I tried it with ex symbolic_matrix(..) but I couldn't get it to work since it is of the type ex and not matrix. Kind regards, Lisa

2 1

Use of STL iterators in GiNaC
by Michael Goffioul 24 Apr '09

24 Apr '09

Hi, I'm trying to compile ginac-1.5.1 with VC++ under Windows. I know "I should get myself a better compiler", but that being said, I've encountered runtime problems due to usage of STL iterators that might be wrong, essentially on std::vector. Problems I've noticed are: dereferencing end(), incrementing end(), decrementing begin()... I'm currently tracking down these issues by running the test suite in the source package, but this might not catch all of them. If you are interested, I can send a patch once I succeed to run the test suite completely. Bye, Michael.

8 25

Re: [GiNaC-list] [CLN-list] msvc support patch
by Richard B. Kreckel 20 Apr '09

20 Apr '09

Hi! Bruno Haible wrote: > Michael Goffioul wrote: >> - struct/class, const/non-const: at some places, CLN foward declares >> objects as being "class", while there are actually "struct"; VC++ generates >> different mangled names for struct and class, so you end up with unresolved >> references. > > What does the C++ standard say about this? Is CLN right, or is MSVC++ right? I don't think it is worth checking what the standard says about this. Even if it permits such sloppiness, they should be avoided. I'll prepare a patch for CLN. -richy. -- Richard B. Kreckel <http://www.ginac.de/~kreckel/>

1 0