Hi Doug,

thanks for your reply. A cartoon of what I would like to do is as follows: I want to compute symbolically a double integral whose integrand is

N_1(u)*e_1(u) + N_2(u)*e_2(u) +...+N_6(u)*e_6(u)   (EQ 1)

u(x,y) is a function defined on a square. N_i(u) and e_i(u) are polynomial functions of u (in general N_i(u)*e_i(u) is a polynomial of order 4 in u). What I need to do is:

1. introduce a Fourier-modes decomposition for u(x,y), for instance if  u(x,y) = \sum \sum u_{ij}* cos(m_i *x) * cos(m_j*y)
2. plug the decomposition in the integrand (EQ 1)
3. compute symbolically the integral of (EQ 1)
4. the integral found will be a function F(u_00,u_10,u_01,u_11,...,u_ij)
5. write F on a matlab file
6. compute the jacobian of F and write it to a Matlab file

In the post linked above, I upload a .mw file in which I do that. u(x,y) is the deformation of a plate, and cos(m_i*x) the corresponding modes. By increasing the number of Fourier modes, I should converge to a better result. The computation gets stuck at step 3.

Hope this is clearer now.

Best.

What's upsetting is that it's hard to believe a sum of integrals of sin and cos powers can last so long, and I don't know how to speed this up