Take the 2-minute tour ×
Mathematica Stack Exchange is a question and answer site for users of Mathematica. It's 100% free, no registration required.

I have found an alternative way to express an old solution:

PwD = (1/(4*Sqrt[Pi])) Integrate[(Exp[-(0.25/t)]/t^1.5)*
    Sum[Exp[-((n^2*Pi^2*t)/40000)] 
        Integrate[
            Exp[-((z^2*Tan[45*Degree]^2 + (200*n - z)^2)/(4*t))] + 
            Exp[-((z^2*Tan[45*Degree]^2 + (100 + 200*n - z)^2)/200)], 
            {z, -50, Plus[50]}
        ], {n, -4, 4}
    ], {t, 0, y}  
];

Solving a partial solution, I have:

PwD = (Exp[-(0.25/t)]/t^1.5)*Sum[Exp[-((n^2*Pi^2*t)/40000)]*
    Integrate[
        Exp[-((z^2*Tan[45*Degree]^2 + (200*n - z)^2)/(4*t))] + 
        Exp[-((z^2*Tan[45*Degree]^2 + (100 + 200*n - z)^2)/200)], 
        {z, -50, Plus[50]}
    ], {n, -4, 4}
];

Now, for plotting the indefinite integral in time (t) I had used the following codes:

function[y_] := (1/(4 Sqrt[π]))*NIntegrate[PwD, {t, 0.01, y}];
tabulate = Table[{y, function[y]}, {y, 0, 100000, 1}];
finalfunction = Interpolation[tabulate];
LogLogPlot[{finalfunction[y], y*finalfunction'[y]}, {y, 1, 1000000},  
  PlotRange -> {0.01, 10}, Frame -> True]

I used this same code on other similar solutions and had positive results. But for this specific solution I couldn't get it to work. Or it's seems that would take a very long time for it.

Could someone provide an alternative way to help Mathematica along?

share|improve this question
3  
As we mentioned last time: Please only post input text (right click, copy as -> input text), and make sure that was you post can be copied back into Mathematica. This is important otherwise people will need to clean up the question manually. –  Szabolcs Mar 13 '12 at 13:56
    
Have you separated out each line to determine which is running slowly? –  rcollyer Mar 13 '12 at 14:34
add comment

1 Answer

up vote 7 down vote accepted

If you break it up as below it will be considerably faster. The idea is to use symbolic integration on the inner integral (Why? Because we can, and it makes this faster).

In[263]:= 
i1[n_, t_] = 
 Integrate[
  Exp[-((z^2 + (200*n - z)^2)/(4*t))] + 
   Exp[-((z^2 + (100 + 200*n - z)^2)/200)], {z, -50, 50}, 
  Assumptions -> {t > 0, -4 <= n <= 4}]

Out[263]= (5*Sqrt[Pi]*(-Erf[10*n] + Erf[10*(1 + n)]))/
  E^(25*(1 + 2*n)^2) + 
   (Sqrt[Pi/2]*Sqrt[t]*(Erf[(25*Sqrt[2]*(1 - 2*n))/Sqrt[t]] + 
           Erf[(25*Sqrt[2]*(1 + 2*n))/Sqrt[t]]))/E^((5000*n^2)/t)

In[264]:= i2[t_] = Sum[i1[n, t], {n, -4, 4}];

In[265]:= 
pwd[y_] := (1/(4*Sqrt[Pi]))*
  NIntegrate[Exp[-1/(4*t)]/t^(3/2)*i2[t], {t, 0, y}]

The tabulation for the interpolation should complete in something approximating reasonable time. Well, within several hours, in any case, if counting by hundreds is any indication.

In[272]:= Timing[tabulate = Table[{y, pwd[y]}, {y, 0, 100000, 100}];]

Out[272]= {150.71, Null}
share|improve this answer
    
Thanks, it really helped! However, would be better (faster) with the introducing of a "counting variable"? With the purpose of taking just a few points on each Log cicle, for eventual interpolation... i'm trying to reach this, but my programming skills isn't so good... Could give me an ideia? Thanks again –  Bruno Rangel Mar 22 '12 at 13:25
add comment

Your Answer

 
discard

By posting your answer, you agree to the privacy policy and terms of service.

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