Hi I am very interested in GiNaC, but I also need symbolic integration, and GiNaC only compute basic integrals. I wrote a simple recursive routine that applies a basic integration rules, matching the input function. It is a good way to implement integration algoritm? Thanks int(c,x) = c*x int(x,x) = x^2/2 int(f+g,x) = int(f,x) + int(g,x) int(c*f,x) = c*int(f,x) int( f^c*1,x) = f^(c+1)/(c+1) or log(f) int( f^c*f',x) = f^(c+1)/(c+1) or log(f) some examples: cout << "int(1+x,x) = " << myinteg( 1 + x, x ) << endl; // 1/2*x^2+x cout << "int(3*x^2+x+6,x) = " << myinteg( 3*pow(x,2)+x+6, x ) << endl; // 1/2*x^2+x^3+6*x cout << "int((3*x+5)^2*3,x) = " << myinteg( pow(3*x+5,2)*3, x ) << endl; // 1/3*(5+3*x)^3 cout << "int((3*x+5)^(-1)*3,x) = " << myinteg( pow(3*x+5,-1)*3, x ) << endl; // log(5+3*x) cout << "int((3*x+5)^b*3,x) = " << myinteg( pow(3*x+5,b)*3, x ) << endl; // (1+b)^(-1)*(5+3*x)^(1+b) the code (with debug messages): ex myinteg( const ex &f, const symbol &x ){ exmap repls; cout << endl << "f(x) = " << f << endl; //int(c,x) c*x if( !has( f, x ) ){ cout << "match int(c,x)" << endl; return f*x; } //int(x,x) x^2/2 if( f.is_equal(x) ){ cout << "match int(x,x)" << endl; return pow(x,2)/2; } //int(f+g,x) int(f,x)+int(g,x) if( f.match( wild(0)+wild(1), repls ) ){ cout << "match int(f+g,x)" << endl; cout << "wild(0)+wild(1) = " << repls << endl; ex f, g; for( exmap::const_iterator i=repls.begin(); i!=repls.end(); ++i){ if(i->first.is_equal(wild(0)) ){ f = i->second; }else if(i->first.is_equal(wild(1))){ g = i->second; } } cout << "f = " << f << endl; cout << "g = " << g << endl; return myinteg( f, x ) + myinteg( g, x ); } //int(c*f,x) c*int(f,x) if( f.match( wild(0)*wild(1), repls ) ){ cout << "match int(c*f,x)" << endl; cout << "wild(0)*wild(1) = " << repls << endl; ex c, f; for( exmap::const_iterator i=repls.begin(); i!=repls.end(); ++i){ if(i->first.is_equal(wild(0)) ){ c = i->second; }else if(i->first.is_equal(wild(1))){ f = i->second; } } if( !has( c, x ) && has( f, x )){ return c*myinteg( f, x ); }if( has( c, x ) && !has( f, x )){ return f*myinteg( c, x ); } } //int( pow( f, c ), x ) pow( f, c+1 ) * pow( c+1, -1 ) if c != -1, or ln( f ) if c == -1 if( f.match( pow(wild(0),wild(1)), repls ) ){ cout << "match int( pow( f, c ), x )" << endl; cout << "pow(wild(0),wild(1)) = " << repls << endl; ex f, c; for( exmap::const_iterator i=repls.begin(); i!=repls.end(); ++i){ if(i->first.is_equal(wild(0)) ){ f = i->second; }else if(i->first.is_equal(wild(1))){ c = i->second; } } cout << "f = " << f << endl; cout << "c = " << c << endl; if( has( f, x ) && !has( c, x ) ){ ex dft = f.diff(x,1); ex df = 1; cout << "df/dft = " << df/dft << endl; if( !has(df/dft,x) ){ cout << "!has(df/dft,x)" << endl; if( c.is_equal(-1) ){ return log( f ) * (df/dft); }else{ cout << "df = " << df << endl; cout << "dft = " << dft << endl; cout << "(df/dft) = " << (df/dft) << endl; cout << "(1/dft) = " << (1/dft) << endl; cout << "pow(dft,-1) = " << pow(dft,-1) << endl; cout << "f = " << f << endl; cout << "c = " << c << endl; cout << "(pow( f, c+1 )/( c+1 )) = " << (pow( f, c+1 )/( c+1 )) << endl; cout << "pow( f, c+1 )*pow( c+1, -1) = " << pow( f, c+1 )*pow( c+1, -1) << endl; cout << "(pow( f, c+1 )/( c+1 ))*(df/dft) = " << (pow( f, c+1 )/( c+1 ))*(df/dft) << endl; cout << "pow( f, c+1 )*pow( c+1, -1) = " << pow( f, c+1 )*pow( c+1, -1)*df*pow( dft, -1) << endl; return (pow( f, c+1 )/( c+1 )) * (df/dft); } } } } //int( pow( f, c ) * f', x ) pow( f, c+1 ) * pow( c+1, -1 ) if c != -1, or ln( f ) if c == -1 if( f.match( wild(0)*pow(wild(1),wild(2)), repls ) ){ cout << "match int( pow( f, c ) * f', x )" << endl; cout << "wild(0)*pow(wild(1),wild(2)) = " << repls << endl; ex df, f, c; for( exmap::const_iterator i=repls.begin(); i!=repls.end(); ++i){ if(i->first.is_equal(wild(0)) ){ df = i->second; }else if(i->first.is_equal(wild(1))){ f = i->second; }else if(i->first.is_equal(wild(2))){ c = i->second; } } cout << "f' = " << df << endl; cout << "f = " << f << endl; cout << "c = " << c << endl; if( has( f, x ) && !has( c, x ) ){ ex dft = f.diff(x,1); cout << "df/dft = " << df/dft << endl; if( !has(df/dft,x) ){ cout << "!has(df/dft,x)" << endl; if( c.is_equal(-1) ){ return log( f ) * (df/dft); }else{ return pow( f, c+1 )*pow( c+1, -1) * (df/dft); } } } } return 0; }