3
$\begingroup$

I wish to implement the following recursion efficiently.

n = 999;
mat1 = RandomReal[{0, 1}, {3, n}, WorkingPrecision -> 20];

mat2 = RandomReal[{0, 1}, {2, n}, WorkingPrecision -> 20];

ct = 1;
While [ct <= n - 1,

 mat1[[1, ct + 1]] = (const1*mat1[[1, ct]] + const2*mat1[[2, ct]])*mat2[[1, ct + 1]];

 mat1[[2, ct + 1]] = (const3*mat1[[1, ct]] + const4*mat1[[2, ct]])*mat2[[2, ct + 1]];

ct++;

];

I was thinking of using FoldList, but I have no idea how...

Edit: I've changed the number of constants.

Any help would be appreciated.

$\endgroup$
2
$\begingroup$

Here is how you can do it using NestList and some code vectorization

mat1[[1 ;; 2]] = Module[{ct = 1}, 
 NestList[({{const1, const2}, {const3, const4}}.#) mat2[[1 ;; 2, ++ct]] &, 
   mat1[[1 ;; 2, 1]], n - 1] //Transpose];

Assuming that the variables const1 to const4 have numerical values close enough to 1, e.g.

const1 = 1.0;
const2 = 0.9;
const3 = 0.95;
const4 = 1.2;

and therefore mat1 will only contain machine size numbers, one can use a compiled function:

fC0 = Compile[{{m1, _Real, 1}, {m2, _Real, 1}},
 ({{const1, const2}, {const3, const4}}.m1) m2, 
 CompilationOptions -> {"InlineExternalDefinitions" -> True}, CompilationTarget -> "C" ];

and then

mat1[[1 ;; 2]] = Module[{ct = 1}, 
 NestList[fC0[#, mat2[[1 ;; 2, ++ct]]] &, mat1[[1 ;; 2, 1]], n - 1] // Transpose];

One can even include the relevant parts of mat2 in a compiled function

fC = Compile[{{m, _Real, 1}, {counter, _Integer, 0}},
 ({{const1, const2}, {const3, const4}}.m) mat2[[1 ;; 2, counter]], 
 CompilationOptions -> {"InlineExternalDefinitions" -> True}, CompilationTarget -> "C"];

and then use

mat1[[1 ;; 2]] = Module[{ct = 1}, 
 NestList[fC[#, ++ct] &, mat1[[1 ;; 2, 1]], n - 1] // Transpose];
$\endgroup$
  • $\begingroup$ I still haven't had the time to study this. This weekend the latest I'll confirm that it works and accept it. Until then, +1. Many thanks for your effort. ;) $\endgroup$ – An old man in the sea. Oct 16 '14 at 20:39
  • $\begingroup$ @Anoldmaninthesea. If you have trouble understanding this code, I can break it down into smaller pieces. But I don't want to spoil your learning experience. $\endgroup$ – Karsten 7. Oct 19 '14 at 4:03
  • $\begingroup$ I've detected a mistake in my question. I had to increase the number of variables 'const'. So, I changed a bit your answer to : amat1[[1 ;; 2]] = Module[{ct = 1}, NestList[({{const1, const2}, {const3, const4}}. #) mat2[[1 ;; 2, ++ct]] &, mat1[[1 ;; 2, 1]], n - 1] // Transpose]; The thing is the gain in efficiency is only of 0.002 or 0.001 seconds with n=999. Is there a way to speed things up? $\endgroup$ – An old man in the sea. Oct 19 '14 at 11:11
  • $\begingroup$ And thanks for the learning experience. $\endgroup$ – An old man in the sea. Oct 19 '14 at 11:15
  • $\begingroup$ @Anoldmaninthesea. Using a compiled function reduces the computation time by a factor two on my PC. $\endgroup$ – Karsten 7. Oct 19 '14 at 17:41

Your Answer

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

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