[GiNaC-list] collecting on similar series expansion

Charlls Quarra charlls_quarra at yahoo.com.ar
Mon Jan 29 23:27:29 CET 2007


--- Sheplyakov Alexei <varg at theor.jinr.ru> escribió:


> Could you please post the complete example which
> illustrates this?
> 

> I tnink this should do the job:
> 
> ex rord = jacobian.expand().collect(lst(F(q),
> F(q).diff(q),
>           F(q).diff(q, 2), F(q).diff(q, 3),
> F(q).diff(q, 4)),
> 					true /* treat the expression as Z[x,y,...] 
> 					        (as opposed to Z[x][y]...) */ );
> 

now seeing that example (and trying it out) i realize
what was my error; i was calling jacobian.expand() and
assumed that was enough to change the jacobian
expression, but when i called collect later on the
(still unexpanded) jacobian ex it wasn't doing the job
at all.

thanks

this is what i was doing

 #include <iostream>
     #include <ginac/ginac.h>
     using namespace std;
     using namespace GiNaC;
    
DECLARE_FUNCTION_1P(F)

static ex F_eval(const ex & x)
{
	return F(x).hold();
}

static ex F_evalf(const ex & x)
{
	return F(x).hold();
}



   REGISTER_FUNCTION(F, eval_func(F_eval).
                            evalf_func(F_evalf).
                            latex_name("\\F"));
 
     int main()
     {
	cout << latex;
         symbol x("x"), y("y");
 //        ex F;
         ex f = pow( x , 5 )*sin( x )+ F(x);
 
         cout << " f is " << f << endl;
         cout << " D(f) is " << f.diff(x) << endl;
         cout << " D(f,2) is " << f.diff(x,2) << endl;


         cout << " integrator " << endl;
         symbol q("q") , p("p");
         symbol m("m"), dt("dt");
         symbol alpha("alpha") , beta("beta") ,
gamma("gamma") , delta("delta"), epsilon("epsilon") ,
nu("nu");
	 symbol dq_iota("dq_iota");
	symbol dq_lambda("dq_lambda");
         ex dq = p/m - F(q)*dt/m -
F(q).diff(q)*dq_iota*pow(dt,2)/(2*m) -
gamma*(F(q).diff(q,2)*pow(dq_iota,2) + F(q).diff(q)*
dq_iota.diff(q) )*pow(dt,3)/(factorial(3) * m);
	 ex dp = -F(q) -F(q).diff(q)*dq_lambda*dt -
delta*(F(q).diff(q,2)*pow(dq_lambda,2) +
F(q).diff(q)*dq_lambda.diff(q) )*pow( dt , 2
)/factorial(3);
	cout << "dq is " << dq << endl; 
	cout << "dp is " << dp << endl;
        exmap s_map;
	s_map[ dq_iota ] = alpha*(p/m - beta*F(q)*dt/m);
	s_map[ dq_lambda ] = epsilon*(p/m - nu*F(q)*dt/m);
	ex dq_new = dq;
	ex dp_new = dp;
	dq_new = dq_new.subs( s_map );
	dp_new = dp_new.subs( s_map );
	cout << "now dq is " << dq_new << endl;
	cout << "now dp is " << dp_new << endl;
	ex jacobian = diff( q + dq_new*dt , q )*diff( p +
dp_new*dt , p ) - diff( q + dq_new*dt , p )*diff( p +
dp_new*dt , q );
	cout << " the jacobian of the integrator is " <<
expand(jacobian) << endl;
        exmap diff_map , inv_diff_map;
	symbol foo0 , foo1 , foo2 , foo3 , foo4;
	lst params;
	params = foo0 , foo1 , foo2 , foo3 , foo4;
	diff_map[ F(q) ] = foo0; inv_diff_map[ foo0 ] = F(q);
	diff_map[ F(q).diff(q) ] = foo1; inv_diff_map[ foo1 ]
= F(q).diff(q);
	diff_map[ F(q).diff(q,2) ] = foo2; inv_diff_map[ foo2
] = F(q).diff(q,2);
	diff_map[ F(q).diff(q,3) ] = foo3; inv_diff_map[ foo3
] = F(q).diff(q,3);
	diff_map[ F(q).diff(q,4) ] = foo4; inv_diff_map[ foo4
] = F(q).diff(q,4);
        ex rord = collect( jacobian.subs( diff_map ) ,
params );
	rord = rord.subs( inv_diff_map );
	cout << "reordered we obtain " << rord << endl;

	ex eqns = jacobian == 1;
	ex vars = alpha;
	cout << lsolve( eqns , vars ) << endl;	
};



	

	
		
__________________________________________________ 
Preguntá. Respondé. Descubrí. 
Todo lo que querías saber, y lo que ni imaginabas, 
está en Yahoo! Respuestas (Beta). 
¡Probalo ya! 
http://www.yahoo.com.ar/respuestas 



More information about the GiNaC-list mailing list