6
$\begingroup$

I have a requirement to integrate some large symbolic expressions containing products of trigonometric functions. A single integration is averaging 7 minutes, and I have 2628 expressions to integrate, so it will take 13 days to complete unless I can speed it up.

Here is an example:

expr = c14 (a0+a1 Cos[n x]+a3 Cos[2 n x]+a5 Cos[3 n x]+a7 Cos[4 n x]+a2 Sin[n x]+a4 Sin[2 n x]+a6 Sin[3 n x]+a8 Sin[4 n x]) (e1 Cos[n x]+e3 Cos[2 n x]+e5 Cos[3 n x]+e7 Cos[4 n x]+e2 Sin[n x]+e4 Sin[2 n x]+e6 Sin[3 n x]+e8 Sin[4 n x])+c44 (e1 Cos[n x]+e3 Cos[2 n x]+e5 Cos[3 n x]+e7 Cos[4 n x]+e2 Sin[n x]+e4 Sin[2 n x]+e6 Sin[3 n x]+e8 Sin[4 n x])^2+c66 (j0+j1 Cos[n x]+j3 Cos[2 n x]+j5 Cos[3 n x]+j7 Cos[4 n x]+j2 Sin[n x]+j4 Sin[2 n x]+j6 Sin[3 n x]+j8 Sin[4 n x])^2+c11 (1/2 (h0+h1 Cos[n x]+h3 Cos[2 n x]+h2 Sin[n x]+h4 Sin[2 n x])^2+3/2 (a0+a1 Cos[n x]+a3 Cos[2 n x]+a5 Cos[3 n x]+a7 Cos[4 n x]+a2 Sin[n x]+a4 Sin[2 n x]+a6 Sin[3 n x]+a8 Sin[4 n x])^2+1/2 (d0+d1 Cos[n x]+d3 Cos[2 n x]+d5 Cos[3 n x]+d7 Cos[4 n x]+d2 Sin[n x]+d4 Sin[2 n x]+d6 Sin[3 n x]+d8 Sin[4 n x])^2)+c14 ((b0+b1 Cos[n x]+b3 Cos[2 n x]+b5 Cos[3 n x]+b7 Cos[4 n x]+b2 Sin[n x]+b4 Sin[2 n x]+b6 Sin[3 n x]+b8 Sin[4 n x]) (d0+d1 Cos[n x]+d3 Cos[2 n x]+d5 Cos[3 n x]+d7 Cos[4 n x]+d2 Sin[n x]+d4 Sin[2 n x]+d6 Sin[3 n x]+d8 Sin[4 n x])+2 (a0+a1 Cos[n x]+a3 Cos[2 n x]+a5 Cos[3 n x]+a7 Cos[4 n x]+a2 Sin[n x]+a4 Sin[2 n x]+a6 Sin[3 n x]+a8 Sin[4 n x]) (e1 Cos[n x]+e3 Cos[2 n x]+e5 Cos[3 n x]+e7 Cos[4 n x]+e2 Sin[n x]+e4 Sin[2 n x]+e6 Sin[3 n x]+e8 Sin[4 n x])+(h0+h1 Cos[n x]+h3 Cos[2 n x]+h2 Sin[n x]+h4 Sin[2 n x]) (g0+g1 Cos[n x]+g3 Cos[2 n x]+g5 Cos[3 n x]+g7 Cos[4 n x]+g2 Sin[n x]+g4 Sin[2 n x]+g6 Sin[3 n x]+g8 Sin[4 n x]))+c12 (1/2 (b0+b1 Cos[n x]+b3 Cos[2 n x]+b5 Cos[3 n x]+b7 Cos[4 n x]+b2 Sin[n x]+b4 Sin[2 n x]+b6 Sin[3 n x]+b8 Sin[4 n x])^2+1/2 (e1 Cos[n x]+e3 Cos[2 n x]+e5 Cos[3 n x]+e7 Cos[4 n x]+e2 Sin[n x]+e4 Sin[2 n x]+e6 Sin[3 n x]+e8 Sin[4 n x])^2+1/2 (g0+g1 Cos[n x]+g3 Cos[2 n x]+g5 Cos[3 n x]+g7 Cos[4 n x]+g2 Sin[n x]+g4 Sin[2 n x]+g6 Sin[3 n x]+g8 Sin[4 n x])^2)+c13 (1/2 (c0+c1 Cos[n x]+c3 Cos[2 n x]+c2 Sin[n x]+c4 Sin[2 n x])^2+1/2 (f0+f1 Cos[n x]+f3 Cos[2 n x]+f5 Cos[3 n x]+f7 Cos[4 n x]+f2 Sin[n x]+f4 Sin[2 n x]+f6 Sin[3 n x]+f8 Sin[4 n x])^2+1/2 (j0+j1 Cos[n x]+j3 Cos[2 n x]+j5 Cos[3 n x]+j7 Cos[4 n x]+j2 Sin[n x]+j4 Sin[2 n x]+j6 Sin[3 n x]+j8 Sin[4 n x])^2);

Integrate[expr, {x, 0, 2 Pi}]

n is an integer, if that helps.

How can I compute these integrals more efficiently?

Take case with this example too:

a12 = Cos[n x] (c14 (a0+a1 Cos[n x]+a3 Cos[2 n x]+a5 Cos[3 n x]+a7 Cos[4 n x]+a2 Sin[n x]+a4 Sin[2 n x]+a6 Sin[3 n x]+a8 Sin[4 n x]) (e1 Cos[n x]+e3 Cos[2 n x]+e5 Cos[3 n x]+e7 Cos[4 n x]+e2 Sin[n x]+e4 Sin[2 n x]+e6 Sin[3 n x]+e8 Sin[4 n x])+c44 (e1 Cos[n x]+e3 Cos[2 n x]+e5 Cos[3 n x]+e7 Cos[4 n x]+e2 Sin[n x]+e4 Sin[2 n x]+e6 Sin[3 n x]+e8 Sin[4 n x])^2+c66 (j0+j1 Cos[n x]+j3 Cos[2 n x]+j5 Cos[3 n x]+j7 Cos[4 n x]+j2 Sin[n x]+j4 Sin[2 n x]+j6 Sin[3 n x]+j8 Sin[4 n x])^2+c11 (1/2 (h0+h1 Cos[n x]+h3 Cos[2 n x]+h2 Sin[n x]+h4 Sin[2 n x])^2+3/2 (a0+a1 Cos[n x]+a3 Cos[2 n x]+a5 Cos[3 n x]+a7 Cos[4 n x]+a2 Sin[n x]+a4 Sin[2 n x]+a6 Sin[3 n x]+a8 Sin[4 n x])^2+1/2 (d0+d1 Cos[n x]+d3 Cos[2 n x]+d5 Cos[3 n x]+d7 Cos[4 n x]+d2 Sin[n x]+d4 Sin[2 n x]+d6 Sin[3 n x]+d8 Sin[4 n x])^2)+c14 ((b0+b1 Cos[n x]+b3 Cos[2 n x]+b5 Cos[3 n x]+b7 Cos[4 n x]+b2 Sin[n x]+b4 Sin[2 n x]+b6 Sin[3 n x]+b8 Sin[4 n x]) (d0+d1 Cos[n x]+d3 Cos[2 n x]+d5 Cos[3 n x]+d7 Cos[4 n x]+d2 Sin[n x]+d4 Sin[2 n x]+d6 Sin[3 n x]+d8 Sin[4 n x])+2 (a0+a1 Cos[n x]+a3 Cos[2 n x]+a5 Cos[3 n x]+a7 Cos[4 n x]+a2 Sin[n x]+a4 Sin[2 n x]+a6 Sin[3 n x]+a8 Sin[4 n x]) (e1 Cos[n x]+e3 Cos[2 n x]+e5 Cos[3 n x]+e7 Cos[4 n x]+e2 Sin[n x]+e4 Sin[2 n x]+e6 Sin[3 n x]+e8 Sin[4 n x])+(h0+h1 Cos[n x]+h3 Cos[2 n x]+h2 Sin[n x]+h4 Sin[2 n x]) (g0+g1 Cos[n x]+g3 Cos[2 n x]+g5 Cos[3 n x]+g7 Cos[4 n x]+g2 Sin[n x]+g4 Sin[2 n x]+g6 Sin[3 n x]+g8 Sin[4 n x]))+c12 (1/2 (b0+b1 Cos[n x]+b3 Cos[2 n x]+b5 Cos[3 n x]+b7 Cos[4 n x]+b2 Sin[n x]+b4 Sin[2 n x]+b6 Sin[3 n x]+b8 Sin[4 n x])^2+1/2 (e1 Cos[n x]+e3 Cos[2 n x]+e5 Cos[3 n x]+e7 Cos[4 n x]+e2 Sin[n x]+e4 Sin[2 n x]+e6 Sin[3 n x]+e8 Sin[4 n x])^2+1/2 (g0+g1 Cos[n x]+g3 Cos[2 n x]+g5 Cos[3 n x]+g7 Cos[4 n x]+g2 Sin[n x]+g4 Sin[2 n x]+g6 Sin[3 n x]+g8 Sin[4 n x])^2)+c13 (1/2 (c0+c1 Cos[n x]+c3 Cos[2 n x]+c2 Sin[n x]+c4 Sin[2 n x])^2+1/2 (f0+f1 Cos[n x]+f3 Cos[2 n x]+f5 Cos[3 n x]+f7 Cos[4 n x]+f2 Sin[n x]+f4 Sin[2 n x]+f6 Sin[3 n x]+f8 Sin[4 n x])^2+1/2 (j0+j1 Cos[n x]+j3 Cos[2 n x]+j5 Cos[3 n x]+j7 Cos[4 n x]+j2 Sin[n x]+j4 Sin[2 n x]+j6 Sin[3 n x]+j8 Sin[4 n x])^2))
$\endgroup$
  • 1
    $\begingroup$ Upload your matrix (as a .nb file) to a repository and provide a link. A sub-matrix is good too. $\endgroup$ – QuantumDot Nov 15 '14 at 12:32
  • $\begingroup$ @QuantumDot, Hi, thanks for your help. But I didn't understand. Please explain it in detail. Thanks $\endgroup$ – Amandeep Nov 15 '14 at 13:06
  • $\begingroup$ You need to supply a couple of examples of the integrations, otherwise we don't know why it's so slow. $\endgroup$ – Simon Woods Nov 15 '14 at 13:08
  • $\begingroup$ @SimonWoods, N2matrix(1,1) = ((6*a02*c11 + 3*a12*c11 + 3*a22*c11 + 3*a32*c11 + 3*a42*c11 + 3*a52*c11 + 3*a62*c11 + 3*a72*c11 + 3*a82*c11 + 2*b02*c12 + b12*c12 + b22*c12 + b32*c12 + b42*c12 + b52*c12 + b62*c12 + b72*c12 + b82*c12 + 2*c02*c13 + c12*c13 + c13*c22 + c13*c32 + c13*c42 + 4*b0*c14*d0 + 2*c11*d02 + 2*b1*c14*d1 + c11*d12 + 2*b2*c14*d2 + c11*d22 + 2*b3*c14*d3 )*Pi)/8 This expression is too much long, I could not print the reuslts. This is the result of integration of one element. $\endgroup$ – Amandeep Nov 15 '14 at 13:12
  • 1
    $\begingroup$ "all elements contain quadratic expressions of trigonometric functions". Is there a pattern that can be generalized into a single or perhaps a few forms that can be integrated to reduce the number of integrations required. $\endgroup$ – Bob Hanlon Nov 15 '14 at 13:15
11
$\begingroup$

Your integrand(s), when expanded, consist of terms that are constants times zero, one or two factors of sine or cosine. The terms can be integrated symbolically and stored as rules. It will be convenient to find the average value and multiply by 2 Pi at the end. (This will automatically handle constant terms without a special rule.) Thanks to Simon Woods for a much simplified set of rules.

Clear[a, b];
iRules = {(Cos | Sin)[_]^2 -> 1/2, (Cos | Sin)[_] -> 0};

ia11 = 2 Pi Expand[a11] /. iRules; // AbsoluteTiming
(*  {0.010274, Null}  *)

It's a lot faster:

evalRules = {   (* assumes  n  is an  Integer *)
  Cos[(a_Integer: 1) n Pi] :> (-1)^a,
  Sin[(b_Integer: 1) n Pi] :> 0
  };
ja11 = Integrate[a11, {x, 0, 2 Pi}] /. evalRules; // AbsoluteTiming
(*  {177.270803, Null}  *)

And it gives the same result, provided the assumptions inherent in evalRules are correct:

ia11 == Expand[ja11] // Simplify
(*  True  *)

Update: Here's another way, about 20 times slower than iRules, but 1000 times faster than Integrate. It's a similar idea, but it should work on general products of sines and cosines. It assumes that sines and cosines are the only functions in the integrands.

ka11 = 2 Pi Expand[TrigToExp@a11] /. E^_ -> 0; // AbsoluteTiming
(*  {0.205734, Null}  *)

ka11 == Expand[ja11 /. evalRules] // Simplify
(*  True  *)

Or, sparked by a comment by the OP, more complicated combinations of sines and cosines can be fairly quickly simplified with TrigReduce.

2 Pi TrigReduce[a11] /. iRules; // AbsoluteTiming
(*  {0.086483, Null}  *)

Original rules

iRules =(*Dispatch@*)
 {Cos[(a_Integer: 1) n x] Sin[(a_Integer: 1) n x] ->
   Integrate[Cos[a n x] Sin[a n x], {x, 0, 2 Pi}]/(2 Pi),
  Cos[(a_Integer: 1) n x] Sin[(b_Integer: 1) n x] /; a == -b -> 
   Integrate[Cos[a n x] Sin[-a n x], {x, 0, 2 Pi}]/(2 Pi),
  Cos[(a_Integer: 1) n x] Sin[(b_Integer: 1) n x] -> 
   Integrate[Cos[a n x] Sin[b n x], {x, 0, 2 Pi}]/(2 Pi),
  Sin[(a_Integer: 1) n x] Sin[(b_Integer: 1) n x] -> 
   Integrate[Cos[a n x] Sin[b n x], {x, 0, 2 Pi}]/(2 Pi),
  Cos[(a_Integer: 1) n x] Cos[(b_Integer: 1) n x] -> 
   Integrate[Cos[a n x] Sin[b n x], {x, 0, 2 Pi}]/(2 Pi),
  Cos[(a_Integer: 1) n x]^2 -> 
   Integrate[Cos[a n x]^2, {x, 0, 2 Pi}]/(2 Pi),
  Sin[(b_Integer: 1) n x]^2 -> 
   Integrate[Sin[b n x]^2, {x, 0, 2 Pi}]/(2 Pi),
  Cos[(a_Integer: 1) n x] -> 
   Integrate[Cos[a n x], {x, 0, 2 Pi}]/(2 Pi),
  Sin[(b_Integer: 1) n x] -> 
   Integrate[Sin[b n x], {x, 0, 2 Pi}]/(2 Pi)};

Some remarks on the code: Note that the rule in iRules are Rule (->) not RuleDelayed (:>). This means that the integrals will be performed when iRules are defined instead of delayed until the application of iRules. Therefore the pattern symbols a and b must be cleared of any values. If that is not possible, then use other symbols or put in a Module (Module[{a, b}, iRules = {...}]).

Dispatch did not speed things up for me. That may change if the list of integration rules were to grow much larger.

$\endgroup$
  • 2
    $\begingroup$ +1 Same approach I took, though since n is an integer I didn't bother with Integrate at all. Just Expand[2 Pi expr]/.{(Cos|Sin)[_]^2 -> 1/2, (Cos|Sin)[_] -> 0} $\endgroup$ – Simon Woods Nov 15 '14 at 14:40
  • $\begingroup$ @SimonWoods Thanks. And duh. I started from the OP's comment that he had products of squares. The adjustment was through cut-and-paste, instead of rethinking the problem. Your approach is how I do integrals in spherical coordinates in my head when teaching vector calculus. (I teach it to the students, too. As Penn & Teller have shown, magic is even more fascinating you can see what's going on.) $\endgroup$ – Michael E2 Nov 15 '14 at 14:55
  • $\begingroup$ @SimonWoods That is a superior way, even if conceptually the same, and I felt the site deserved the best answer. So I included it and converted the answer to CW. (Or you could answer, and I'll delete.) $\endgroup$ – Michael E2 Nov 15 '14 at 15:19
  • $\begingroup$ No, please don't delete your original. I think it's good to keep it as it shows how you would approach similar problems in cases where the integrals aren't trivial. $\endgroup$ – Simon Woods Nov 15 '14 at 16:07
  • $\begingroup$ Thanks Simon woods and Michael E2. Thanks a lot $\endgroup$ – Amandeep Nov 15 '14 at 20:32

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.