MathGroup Archive 2011

[Date Index] [Thread Index] [Author Index]

Search the Archive

Re: Speeding Up Known Integrations

  • To: mathgroup at smc.vnet.net
  • Subject: [mg120983] Re: Speeding Up Known Integrations
  • From: "Sean Violante" <sean.violante07 at imperial.ac.uk>
  • Date: Sat, 20 Aug 2011 06:17:14 -0400 (EDT)
  • Delivered-to: l-mathgroup@mail-archive0.wolfram.com
  • References: <201108191034.GAA24006@smc.vnet.net> <4E4E906F.4030200@wolfram.com>

Hi Daniel
Thanks for your response.  Unfortunately you changed the integration
variables (w0t2 became w0t1) so that it became a 1-d integral - and 
thats why it sped up.

so my original integrand ( iterated ) takes 125 seconds
expanding and doing it as double integral takes 147 seconds
 
expanding and distributing the double integral (ie your method) takes 225s

What I am trying (see below) takes about 0.07 seconds, but I am not
confident in the results!


here is the integral again.  I am able to read this back into mathematica.

int = Integrate1[
  1/StdF1[t] E^(-t \[Kappa] +
     2 W0t2 \[Kappa]) (VC0 + E^(-c W0t2) VCC +
      E^(-W0t2 \[Kappa]) VCK)^2 Integrate1[
     E^(W0t1 \[Kappa]) (VC0 + E^(-c W0t1) VCC +
        E^(-W0t1 \[Kappa]) VCK) \[Rho] ((
        E^(-t \[Kappa] +
          W0t1 \[Kappa]) (VpC0 + E^(-c W0t1) VpCC) \[Kappa] Subscript[
         a, 1] Subscript[\[Xi], 2])/(
        c StdF1[t] - \[Kappa] StdF1[t]) + (
        E^(-c t +
          c W0t1) (VpC0 + E^(-c W0t1) VpCC) \[Kappa] Subscript[a, 1]
          Subscript[\[Xi], 2])/(-c StdF1[t] + \[Kappa] StdF1[t]) + (
        E^(-c t + c W0t1) (VpC0 + E^(-c W0t1) VpCC) Subscript[a, 2]
          Subscript[\[Xi], 2])/StdF1[t]), {W0t1, 0, W0t2}] Subscript[
    a, 1] Subscript[\[Xi], 1], {W0t2, 0, t}]

ReplaceExponential attempts to identify terms of the form t1^n e^{k t1} and
renames them e[n,k]
Then we replace each e[] term by its corresponding integral, Inte[,].

ReplaceExponential[
  WW0t1_] :=  {WW0t1 ^n_ E^z_ :>
   E^Coefficient[z, WW0t1, 0] e[n, Coefficient[z, WW0t1, 1]],
  WW0t1 E^z_ :>
   E^Coefficient[z, WW0t1, 0] e[1, Coefficient[z, WW0t1, 1]],
   E^z_ :> E^Coefficient[z, WW0t1, 0] e[0, Coefficient[z, WW0t1, 1]],
  WW0t1 ^n_. :> e[n, 0]}

Clear[Inte] (* NB definition has to be in this order, otherwise the \
for loop defn override for eg Inte[n,0,t] *)
Inte[n_, 0, t_] = t^(n + 1)/(n + 1)
Inte[0, 0, t_] = t
For[i = 0, i < 10, i++,
 Inte[i, k_, t_] = Integrate[t1^i E^(k t1), {t1, 0, t}]]
(* to avoid Gamma Function *)

Note that intB is the inner integrand and intA is the outer integrand
Timing[
 i1 = (Expand[intB] /. ReplaceExponential[W0t1]);
 i2 = i1 /. e[n_Integer, k_] :> Inte[n, k, W0t2];
 i3 = Expand[intA i2] /. ReplaceExponential[W0t2];
 i4 = i3 /. e[n_Integer, k_] :> Inte[n, k, t];]


On 19/08/2011 18:33, Daniel Lichtblau wrote:
On 08/19/2011 05:34 AM, Sean Violante wrote:

I have a large collection of iterated integrals to evaluate.
These are basically sums of terms of the form.

   t_1^n1 E^{ k1 t_1}]    t_2^n2 E^{ k2 t_2}  t_3^n3 E^{ k3 t_3}

Is there any way of speeding up this evaluation?  I was wondering
whether gathering up terms in the above form, and then applying the
integral to each would work out faster, or is this effectively what
Integrate is doing anyway?

If faster, how would you suggest gathering up terms?
I have tried a few simple things eg  multiplying by exp[n (k+c)W0t2]
( where n is the largest negative power of Exp[ k t] and Exp [c t] ie so 
expression  becomes polynomial)
and then using CoefficientList [ %,{W0t2, e^k W0t2, e^c W0t2} ]
However, I never quite manage to pull out all the exponential terms

Thanks for your help

Sean

below is an example ( where Integrate1 will be replaced by integrate).
It takes about 100 seconds on my machine.
   Integrate1[
     1/StdF1[t] \[ExponentialE]^(-t \[Kappa] +
        2 W0t2 \[Kappa]) (VC0 + \[ExponentialE]^(-c W0t2)
           VCC + \[ExponentialE]^(-W0t2 \[Kappa])
           VCK)^2 Integrate1[\[ExponentialE]^(
         W0t1 \[Kappa]) (VC0 + \[ExponentialE]^(-c W0t1) 
             VCC + \[ExponentialE]^(-W0t1 \[Kappa])
             VCK) \[Rho] ((\[ExponentialE]^(-t \[Kappa] +
             W0t1 \[Kappa]) (VpC0 + \[ExponentialE]^(-c W0t1)
                VpCC) \[Kappa] Subscript[a, 1] Subscript[\[Xi], 2])/(
           c StdF1[t] - \[Kappa] StdF1[t]) + (\[ExponentialE]^(-c t +
             c W0t1) (VpC0 + \[ExponentialE]^(-c W0t1)
                VpCC) \[Kappa] Subscript[a, 1] Subscript[\[Xi],
            2])/(-c StdF1[t] + \[Kappa] StdF1[
              t]) + (\[ExponentialE]^(-c t +
             c W0t1) (VpC0 + \[ExponentialE]^(-c W0t1) VpCC) Subscript[
            a, 2] Subscript[\[Xi], 2])/StdF1[t]), {W0t1, 0,
         W0t2}] Subscript[a, 1] Subscript[\[Xi], 1], {W0t2, 0, t}] \!
\*SubsuperscriptBox[\(\[Xi]\), \(1\), \(2\)] /.

This did not parse correctly when I tried to cut and paste it into
Mathematica. Below is what I think was intended, with Greek and subscripted
items renamed to similar in straight ascii. What seems to give better speed
is to expand the integrand, and map the integration over that sum.

In[86]:= integrand = Expand[E^(W0t1*kappa)*(VC0 + VCC/E^(c*W0t1) + VCK/E^(W0t1*kappa))*
         rho*((E^((-t)*kappa + W0t1*kappa)*(VpC0 + VpCC/E^(c*W0t1))*
         kappa*a1*
                 eps2)/(c*StdF1[t] - kappa*StdF1[t]) +
            (E^((-c)*t + c*W0t1)*(VpC0 + VpCC/E^(c*W0t1))*kappa*a1*
         eps2)/
              ((-c)*StdF1[t] + kappa*StdF1[t]) +
            (E^((-c)*t + c*W0t1)*(VpC0 + VpCC/E^(c*W0t1))*a2*eps2)/
              StdF1[t])*a1*eps1];

In[87]:= Timing[
 result = Map[
    Integrate[#, {W0t1, 0, W0t2}, {W0t2, 0, t},
      Assumptions -> {0 < W0t2 < t}] &, integrand];]

Out[87]= {5.88, Null}

If you have conditions on the various variables (positivity, etc.) then
adding those to the Assumptions might improve speed further. Or maybe not,
but worth experimenting with if possible.

Daniel Lichtblau
Wolfram Research






  • Prev by Date: Re: Simple Doppler effect
  • Next by Date: time polarized em wave interferometer
  • Previous by thread: Re: Speeding Up Known Integrations
  • Next by thread: CDF files now on Gyre&Gimble (math blog)