I am trying to do transfer matrix calculations (optics), where I have a layered sample and want to calculate how light propagates through it, depending on some parameters:

enter image description here

The charm of the method is that you can do a matrix multiplication for every subsequent layer, no matter how many layers you have; you just need a transition matrix T between every layer and a propagation matrix P for passing through the layer, which can be easily defined. An exemplary transfer matrix of a single layer with air on both sides would look like this:
M[some properties] := T[air properties,layer properties].P[layer properties].T[layer properties,air properties]

What I want to do now is a script, where you can give a number of layers such as n=3 and the transfer matrix will be automatically generated (of course you have to give a list of parameters as well). I thought of doing this recursively with a for-loop, such as:

For[i=1, i<n+1, i++, M[..] := M[..].P[properties i].T[properties i, properties i+1]];

but this does not seem to work, probably because I am declaring a function M and then overwrite it (M is initialized as a unity matrix before this function).

The thing is, for a higher number of layers the transition matrix becomes a bit unhandy, and even more so if you want to edit it manually. I was hoping to generate this function dynamically from the input of n and a list of material properties and then store it in a file, to be able to transfer it to a different script later.

Is there a way to do this? Where did I go wrong? I know that := is used to delay the calculation of the value (effectively making the thing left of it a function) and that does not seem suited for recursive use... but I also tried different variations and nothing worked. Also, the input parameters of the function in the brackets M[...] will probably be a problem, because it will be more parameters for each layer. Can I bridge this? I know that I can input anything into a function but I need to name it, because I need the values at specific points of the function.

Somehow I can't believe that Mathematica does not support the dynamic generation of functions depending on input, it seems like such an obvious thing to me. I tried to search for solutions but the name of the Dynamic[] function interferes with my search and I couldn't think of a better word to describe the type of function that I want.

EDIT: so here is some real code:

P[K_, e_, Q_, d_] := ...
T[K_, ea_, eb_, Q_] := ...

M[K_,Q_,e1_,d1_,e2_] := T[K,e0,e1[K],Q].P[K,e1[K],Q,d1].T[K,e1[K],e2[K],Q];

M[K_,Q_,e1_,d1_,e2_,d2_,e3_] := T[K,e0,e1[K],Q].P[K,e1[K],Q,d1].T[K,e1[K],e2[K],Q].P[K,e2[K],Q,d2].T[K,e2[K],e3[K],Q];

M[K_,Q_,e1_,d1_,e2_,d2_,e3_,d3_,e4_] := T[K,e0,e1[K],Q].P[K,e1[K],Q,d1].T[K,e1[K],e2[K],Q].P[K,e2[K],Q,d2].T[K,e2[K],e3[K],Q].P[K,e3[K],Q,d3].T[K,e3[K],e4[K],Q];

these three examples are:

  1. air - layer 1 - substrate
  2. air - layer 1 - layer 2 - substrate
  3. air - layer 1 - layer 2 - layer 3 - substrate

I tried it this way now:

P[K_, e_, Q_, d_] := ...;
T[K_, ea_, eb_, Q_] := ...;

M[K_, Q_, e_, d_, m_] :=  T[K, 1, e[[1, K]], Q].For[i = 1, i < m + 1, i++, P[K, e[i, K], Q, d[i]].T[K, e[i, K], e[i + 1, K], Q].];

e = {e1, e2, e3, e4};
d = {d1, d2, d3};
n = 3;

M[Kvalue, Qvalue, e, d, n]

So if only this Matrix multiplication would work with a loop, I would be fine :/ but unfortunately, it doesn't seem to work this way.

Any help would be greatly appreciated!

  • $\begingroup$ Why couldn't you use Piecewise[]? $\endgroup$ May 27, 2020 at 11:48
  • $\begingroup$ I don't really see how I could use Piecewise[] to get the desired result. Would you test for n and then return a different function for each n? If so, I would have to put in all options manually as far as I can see. Also, can I put a function definition with arguments in there? I didn't find an example of that in the documentation... $\endgroup$ May 27, 2020 at 13:03
  • $\begingroup$ Could you give us a complete example for, say, two layers? at the moment, the question is a bit fuzzy and, although I understand your physical problem in general, I do not know exactly what you need to do in Mathematica so I can't help much. $\endgroup$
    – MarcoB
    May 27, 2020 at 14:52
  • $\begingroup$ You are trying to modify a function so it generates different matrices based on the layer correct? Couldn't you just write a second function that calls M with the different parameters for each layer? Also, are you aware of Function? $\endgroup$
    – flinty
    May 27, 2020 at 15:02
  • $\begingroup$ @MarcoB I edited the post to include real code (even thought I changed some things, it should represent it adequately). $\endgroup$ May 27, 2020 at 17:48

1 Answer 1


You don't need to modify functions, you just need to define the e[i, K] and d[i, K] and call this:

M[K_, Q_, n_] := 
 Fold[Dot, T[K, e0, e[1, K], Q], 
  P[K, e[#, K], Q, d[#, K]].T[K, e[#, K], e[# + 1, K], Q] & /@ Range[n]]

M[K, Q, 1] 
(* returns T[K, e0, e[1, K], Q].P[K, e[1, K], Q, d[1, K]].T[K, e[1, K], e[2, K], Q] *)

M[K, Q, 2] 
(* returns T[K, e0, e[1, K], Q].P[K, e[1, K], Q, d[1, K]].T[K, e[1, K], e[2, K], 
  Q].P[K, e[2, K], Q, d[2, K]].T[K, e[2, K], e[3, K], Q] *)

M[K, Q, 3] 
(* returns T[K, e0, e[1, K], Q].P[K, e[1, K], Q, d[1, K]].T[K, e[1, K], e[2, K], 
  Q].P[K, e[2, K], Q, d[2, K]].T[K, e[2, K], e[3, K], Q].P[K, e[3, K],
   Q, d[3, K]].T[K, e[3, K], e[4, K], Q] *)
  • $\begingroup$ I didn't know you could use Dot like this... that looks very elegant and exactly what I wanted to do! I will try it out tomorrow and let you know if this works like I want it to. Thanks! $\endgroup$ May 27, 2020 at 20:12
  • $\begingroup$ Did this work for you and is my answer acceptable? $\endgroup$
    – flinty
    May 31, 2020 at 21:16
  • $\begingroup$ Sorry, as I simplified the problem a bit for this question I needed some time to implement it properly. It does exactly what I want it to, thanks! $\endgroup$ Jun 3, 2020 at 21:32

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.