1
$\begingroup$

This question is followed up from this Question. I would like to thank Dr. Hintze and I_Mariusz for the comments and help. I am pretty new to mathematica ( I just learned it 4 days ago) so I would like to ask if there is a way to accelerate the following codes : 

Clear["Global`*"]
eq = (2*S0/(y*sigma^2))^
nu*(Gamma[nu + (2*mu)/sigma^2]/Gamma[2*nu + (2*mu)/sigma^2])*
Hypergeometric1F1[nu, 2*nu + (2*mu)/sigma^2, -2*S0/(y*sigma^2)];
int = Integrate[eq, {y, K, Infinity}, GenerateConditions -> False]
int2 = int /. nu -> (-vu/2 + Sqrt[vu^2 + 8*alpha/sigma^2]/2) /. 
 vu -> (2*mu/sigma^2 - 1)

mu = 15/100;
sigma = 5/100;
S0 = 100;
K = 95;
T = 1;
F[alpha_] = int2/alpha;

A = 18.4;
n = 40;
m = 50;
S = ConstantArray[0, m + 1];
B = 1;
For[k = 0, k <= m, k++,
   S[[k + 1]] = Exp[A/2]/(2*B)*Re[F[A/(2*B)]];
   For[j = 1, j <= n + k, j++,
    S[[k + 1]] = 
S[[k + 1]] + Exp[A/2]/B*(-1)^j*Re[F[(A + 2*j*Pi*I)/(2*B)]];
]    ]
f = 0
   For[k = 0, k <= m, k++,
   f = f + Binomial[m, k]*S[[k + 1]]*2^(-m);
   ]
f

Currently, it takes me hours to run this program. Thank you so much for your time. I truly appreciate it.

Improve this question
$\endgroup$
8
  • 1
    $\begingroup$ I believe that if you change F[alpha_] = int2/alpha to F[alpha_] := int2/alpha this will run quite quickly. Using = in a function assignment doesn't actually assign it as a pattern (see the documentation on Set (=) versus SetDelayed (:=)). $\endgroup$
    – nben
    Commented Dec 15, 2015 at 16:30
  • $\begingroup$ Thank you ! please let me try . $\endgroup$
    – KNN
    Commented Dec 15, 2015 at 16:31
  • $\begingroup$ If you make that a delayed definition you need to also make int2 an explicit function of alpha. The other thing I'd suggest is to convert your For loop to a Table[ Sum[] ] structure. $\endgroup$
    – george2079
    Commented Dec 15, 2015 at 16:49
  • $\begingroup$ george2079 : Thanks ! , let me understand how Table is used. I've have never used it before. $\endgroup$
    – KNN
    Commented Dec 15, 2015 at 16:51
  • 1
    $\begingroup$ The other thing you should do here is note you are frequently recomputing F for the same argument. Define F so it remembers: F[alpha_] := F[alpha] = int2[alpha]/alpha; ( or restructure your loop so you aren't repeating calculations , that inner For loop just adds one new term for each k.) $\endgroup$
    – george2079
    Commented Dec 15, 2015 at 16:57

1 Answer 1

Reset to default
3
$\begingroup$

This takes about a minute.

(* with int2 redefined as : int2[alpha_] = .. *)
F[alpha_] := F[alpha] = int2[alpha]/alpha;
terms = Monitor[ 
  Exp[A/2]/B Table[ (-1)^j*Re[F[(A + 2*j*Pi*I)/(2*B)]] ,
       {j, 0, n + m + 1}], {j}]
S = Accumulate[terms][[n ;;]]
f = Sum[Binomial[m, k]*S[[k + 1]]*2^(-m), {k, 0, m}]

I'd encourage you to study this carefully to be sure it is faithful to your calculation. (I'm sure its close but I could have shifted an index or some such doing it quickly)

Improve this answer
$\endgroup$
0

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

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