I am looking for a more effcient method to perform the following analysis. I did my best, but running the following code on my machine takes 9.766s. I need to perform this analysis several times.


Jfunction = 
Compile[{{NOunknownPar, _Integer, 0}, {dt, _Real, 
0}, {Kgeneral, _Real, 2}, {Jacob1, _Real, 2}}, 
Module[{Ident, sj1, sj2, sj3, sf},
Ident = IdentityMatrix[2*NOunknownPar]; 
sj1 = Ident - (0.5*dt)*Kgeneral; sj2 = Ident + (0.5*dt)*Kgeneral; 
sj3 = Dot @@ {sj2, Jacob1}; sf = LinearSolve[sj1, sj3]; sf], 
CompilationTarget -> "C"];

NoVar = 79;
NoRep = 6001;
deltaT = 0.0001;
Jmat = IdentityMatrix[2*NoVar];
KK = ConstantArray[0, NoRep - 1];
Do[KK[[i]] = RandomReal[1, {2*NoVar, 2*NoVar}], {i, NoRep - 1}];

Do[Jmat = Jfunction[NoVar, deltaT, KK[[jJa]], Jmat], {jJa, 
NoRep - 1}] // Timing
  • 4
    $\begingroup$ There's not much use putting LinearSolve[] inside a compiled function. Also, sj3 = sj2.Jacob1 is much cleaner-looking. But, more importantly, how about stepping back from the code, and try to elaborate what you're trying to do here? $\endgroup$ – J. M. is away Jun 8 '15 at 9:53
  • 1
    $\begingroup$ @Guesswhoitis Thanks for your comments. Probably it is not easy to explain exactly why I need this code, but generally I am trying to perform a non-linear dynamic analysis and in order to compute the Jacobian of the system with respect to the initial conditions, I need to solve a similar matrix equation at every time-step. $\endgroup$ – mak maak Jun 8 '15 at 10:05
  • $\begingroup$ That last loop can be more cleanly done with FoldList[]. Also, KK = RandomReal[1, {NoRep - 1, 2*NoVar, 2*NoVar}] is a better initialization. $\endgroup$ – J. M. is away Jun 8 '15 at 10:19

Your Answer

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

Browse other questions tagged or ask your own question.