[update] The API shown here will most definitely change.
During my internship at LKB, Paris I’m expected to prepare a number of optical setups. I thought that it would be interesting to do the work in SymPy. That’s why I started working on the gaussopt module in sympy.physics (for the foreseeable future it won’t be upstream). I plan to make it small and basic (containing only the most basic functionality).
For the moment the only things on the blueprint are ray transfer matrices, geometric and gaussian beams.
To define a ray transfer matrix for free space or for a thin lens you just use the build in classes (a general matrix is also available):
In [1]: from sympy.physics.gaussopt import *
In [2]: distance = Symbol('distance', positive=True)
In [3]: focus = Symbol('focus', real=True)
In [4]: fs = FreeSpace(distance)
In [5]: tl = ThinLens(focus)
In [6]: fs
Out[6]:
⎡1 distance⎤
⎢ ⎥
⎣0 1 ⎦
In [7]: tl
Out[7]:
⎡ 1 0⎤
⎢ ⎥
⎢ -1 ⎥
⎢───── 1⎥
⎣focus ⎦
The module supports geometric beams but gaussian beams are more interesting, so I’ll show only them. The formalism that supports them is based on the complex beam parameter. Here is how you define a beam at its waist:
In [8]: waist = Symbol('waist', positive=True)
In [9]: wavelen = Symbol('wavelen', positive=True)
In [10]: z_rayleigh = waist**2 * pi / wavelen
In [11]: p = BeamParameter(0, z_rayleigh, wavelen)
In [12]: p.q
Out[12]:
2
ⅈ⋅π⋅waist
──────────
wavelen
In [13]: p.w
Out[13]: waist
In [14]: p.w_0
Out[14]: waist
In [15]: p.R
Out[15]: 0
In [16]: p.z
Out[16]: 0
In [17]: p.z_r
Out[17]:
2
π⋅waist
────────
wavelen
For the different properties shown here the names should be suggestive enough. The pretty_printer supports Greek letters for the wavelen but this is not so important here.
Now I would like to calculate how a telescopic system will act on the beam. For example
In [18]: tl1 = ThinLens(focus/10)
In [19]: p_out = tl1*( fs*( tl*p))
The last line takes about 60 seconds and produces very UNsimplified expression. I suppose that I should do some implicit simplification. Also notice the parentheses. I should check if they are really required.
simplify(p_out) takes too much time so I’ll cheat a bit (from symbolics to numerics):
In [20]: p_out1 = p_out.subs({focus:1, distance:1*11/10,
wavelen:532e-9, waist:5e-3})
In [21]: N(p_out1.q)
Out[21]: -0.110000000001061 + 1.47631233721549⋅ⅈ
In [22]: N(p_out1.w_0)
Out[22]: 0.000500000000000380
In [23]: N(p_out1.z)
Out[23]: -0.110000000001061
So with balanced telescopic system of power 10 the waist diminished 10 times and will be obtained at 11 cm after the second lens. Nice.
It would be nice to post some more realistic problems later on.