5

I am not very used to do numerical simulations on Mathematica. Do you have any ideas how to improve i.e. speed up my code?

f4[a_?NumericQ, b_?NumericQ, c_?NumericQ, delta_?NumericQ, 
K_?NumericQ, d_?NumericQ] := Module[
{kk1, kk2, out, sigma, beta, rho},
kk1 = NDSolve[{xt'[t] == 10 (yt[t] - xt[t]),
  yt'[t] == xt[t] (28 - zt[t]) - yt[t],
  zt'[t] == xt[t] yt[t] - 8/3 zt[t],
  xt[0] == yt[0] == zt[0] == 1},
 {xt, yt, zt}, {t, 15}];
sigma = {13.25, 7, 6.5};
rho = {19, 18, 38};
beta = {3.5, 3.7, 1.7};
kk2 = 
NDSolve[Join[
  Table[x[i]'[t] == 
    sigma[[i]] (y[i][t] - x[i][t]) + 
     a Sum[x[j][t] - x[i][t], {j, 1, 3}], {i, 1, 3}],
  Table[
   y[i]'[t] == 
    x[i][t] (rho[[i]] - z[i][t]) - y[i][t] + 
     b Sum[y[j][t] - y[i][t], {j, 1, 3}], {i, 1, 3}],
  Table[
   z[i]'[t] == 
    x[i][t] y[i][t] - beta[[i]] z[i][t] + 
     c Sum[z[j][t] - z[i][t], {j, 1, 3}], {i, 1, 3}],
  Table[x[i][0] == 1, {i, 1, 3}],
  Table[y[i][0] == 1, {i, 1, 3}],
  Table[z[i][0] == 1, {i, 1, 3}]], 
 Join[Table[x[i], {i, 1, 3}], Table[y[i], {i, 1, 3}], 
  Table[z[i], {i, 1, 3}]], {t, 15}];
outt[tt_] := {xt[tt + 5], yt[tt + 5], zt[tt + 5]} /. kk1;
ti = Table[i + 5, {i, 0, (K - 1)*d, d}];
out[tt_] := {1/3 (Sum[x[i][tt + 5], {i, 3}]), 
  1/3 (Sum[y[i][tt + 5], {i, 3}]), 
  1/3 (Sum[z[i][tt + 5], {i, 3}])} /. kk2;
dist[tt_] := ((out[tt][[1, 1]] - 
      outt[tt][[1, 1]])^2 + (out[tt][[1, 2]] - 
      outt[tt][[1, 2]])^2 + (out[tt][[1, 3]] - 
      outt[tt][[1, 3]])^2)*0.4^tt;
FC = 1/(K delta ) Sum[
  NIntegrate[dist[tt], {tt, ti[[i]], ti[[i]] + delta}], {i, 1, K}];
FC];
f4[1, 1, 1, 1, 10, 0.2] // AbsoluteTiming

Thanks a lot!!

HarryH
  • 53
  • 1
  • 3
  • 3
    If you could, I don't know, give a short description of what this code of yours is actually trying to do, we might be able to suggest a better approach... – J. M.'s missing motivation Jul 27 '12 at 11:49
  • Basically 3 oscillator systems are coupled to each other via a, b and c. In the end I want to minimize the function FC with regard to these parameters, but this takes quite too long with the current code... :-( Just the general question whether my code is very cumbersome?! – HarryH Jul 27 '12 at 12:10

1 Answers1

23

If you insert a bunch of commands like "Print@First@AbsoluteTiming[...]" in the middle of your function you will see that literally all time is spent on the last line with the Sum[NIntegrate[...]].

To solve your problem insert the following commands into NIntegrate:

Method -> {Automatic, "SymbolicProcessing" -> 0}

Caution: I have personally encountered situations where this option impressively reduces the computational burden with no change to the result whatsoever; conversely, I have also seen situations where it gives the wrong answer. I am not completely aware of how this function works so my advice is: check your answer without it every now and then.

Gabriel Landi
  • 2,590
  • 1
  • 19
  • 22
  • 1
    Great!!! Thanks a lot, this is exactly was I was hoping for!! – HarryH Jul 27 '12 at 12:41
  • 4
    @Gabriel Landi, if you do find such cases where Method -> {Automatic, "SymbolicProcessing" -> 0} give incorrect results please send them to support@wolfram.com for investigation. Thanks. –  Jul 27 '12 at 13:47
  • @ruebenko Ok, sure thing! :) – Gabriel Landi Jul 27 '12 at 13:47