Gaussian optics in Sympy: studying a lens

In the last post I presented the small module for gaussian optics that I’ve created for sympy. Here I’ll try to show a real problem solved with it and with sympy.

[update]: The api is still changing. Also here I’m using ABCD matrices for something that is more easily done with conjugation relations. It’s that way so I can show how the matrices are used in the module.

The problem is finding how a laser beam is modified by passing trough a lens. What is the new waist? Where is it obtained?

Making the boilerplate

Let us first introduce the variables for the input beam and the lens. Notice that as input for the creation of the BeamParameter the Rayleigh length z_r is used and not the waist w_0 = \sqrt{\frac{z_r \lambda}{\pi}}. That’s the way for now just because it gives simpler expressions:

In [1]: from sympy.physics.gaussopt import *

In [2]:

In [3]: z_r = Symbol("z_r", positive=True)

In [4]: l = Symbol("lambda", positive=True)

In [5]: z = Symbol("z", real=True)

In [6]: f = Symbol("f", real=True)

In [7]: z_r
Out[7]: zᵣ

In [8]: l
Out[8]: λ

Defining the parameter and the lens and then calculating the parameter for the output beam:

In [9]: q_in = BeamParameter(z, z_r, l)

In [10]: tl = ThinLens(f)

In [11]:

In [12]: q_out= tl*q_in

Let’s see what does this give for the waist of the output parameter:

In [13]: w_0_out = q_out.w_0

In [14]: z_out = q_out.z

In [15]: w_0_out
Out[15]:
                     ⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽
  ⎽⎽⎽   ⎽⎽⎽⎽        ╱         1
╲╱ λ ⋅╲╱ zᵣ ⋅      ╱  ──────────────────
                  ╱              2     2
                 ╱        2⋅z   z    zᵣ
                ╱     1 - ─── + ── + ───
               ╱           f     2     2
             ╲╱                 f     f
─────────────────────────────────────────
                ⎽⎽⎽
              ╲╱ π

In [16]: z_out
Out[16]:
                        2                        2
z                      z                       zᵣ
────────────────── - ────────────────────── - ──────────────────────
           2     2     ⎛           2     2⎞     ⎛           2     2⎞
    2⋅z   z    zᵣ      ⎜    2⋅z   z    zᵣ ⎟     ⎜    2⋅z   z    zᵣ ⎟
1 - ─── + ── + ───   f⋅⎜1 - ─── + ── + ───⎟   f⋅⎜1 - ─── + ── + ───⎟
     f     2     2     ⎜     f     2     2⎟     ⎜     f     2     2⎟
          f     f      ⎝          f     f ⎠     ⎝          f     f ⎠

In [17]: together(z_out) # to simplify it a bit
Out[17]:
  ⎛       2     2⎞
f⋅⎝f⋅z - z  - zᵣ ⎠
─────────────────────
 2            2     2
f  - 2⋅f⋅z + z  + zᵣ

Those are still mostly unreadable. To simplify the expressions further I’ll introduce the dimensionless parameters Z_i = \frac{z}{f} and Z_r = \frac{z_r}{f}. Also there is the obvious substitution for the input beam waist.

In [19]: Z_in = Symbol("Z_i", real=True)

In [20]: Z_r = Symbol("Z_r", real=True)

In [21]:

In [22]: Z_out = z_out/f

In [23]: Z_out = z_out.subs({z : Z_in*f, z_r : Z_r*f})

In [24]: Z_out = together(Z_out)

In [25]: Z_out
Out[25]:
  ⎛    2          2⎞
f⋅⎝- Zᵢ  + Zᵢ - Zᵣ ⎠
────────────────────
 2            2
Zᵢ  - 2⋅Zᵢ + Zᵣ  + 1

And same for the waist:

In [26]: w_0_in = Symbol("w_0^i", positive=True)

In [27]: w_0_in
Out[27]: wⁱ₀

In [28]: w_0_in_expression = q_in.w_0

In [29]: w_0_in_expression
Out[29]:
  ⎽⎽⎽   ⎽⎽⎽⎽
╲╱ λ ⋅╲╱ zᵣ
────────────
  ⎽⎽⎽
╲╱ π

In [30]:

In [31]: w_0_out = w_0_out.subs({w_0_in_expression : w_0_in, z : Z_in*f, z_r : Z_r*f})

In [32]: w_0_out
Out[32]:
          wⁱ₀
─────────────────────────
   ⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽
  ╱   2            2
╲╱  Zᵢ  - 2⋅Zᵢ + Zᵣ  + 1

Let’s study how the waist changes. First of all, there better be no solutions to the quadratic equations in the denominator. We want it always positive in order to have real waist. Let’s extract the expression.

In [41]: quad_equ = w_0_out.args[1].args[0]

In [42]: quad_equ
Out[42]:
 2            2
Zᵢ  - 2⋅Zᵢ + Zᵣ  + 1

Ok, obviously the discriminant is -Z_r^2 < 0 and we have no real solutions. But just as a crosscheck let’s pass it trough sympy’s solve.

In [44]: solve(quad_equ, Z_in)
 Out[44]: [-ⅈ⋅│Zᵣ│ + 1, ⅈ⋅│Zᵣ│ + 1]

Great. The waist is real and positive.

Some useful results

Until now we haven’t done anything really interesting from physical point of view. We have just derived some expressions. Here are some more interesting questions.

What is the maximal possible waist obtainable with a single lens?

To find it we minimize the denominator. I lost a few minutes searching for a maximize function in sympy but didn’t find one. But obviously the maximum is at Z_i = 1 and so w_{out}^{max} = \frac{w_{in}}{Z_r}=\frac{f \lambda}{\pi w_{in}}.

It’s quite interesting that the maximum is for Z_i = 1 (i.e. z_i = f). Let’s see what does that mean for a convergent lens (it’s mostly the same for a divergent one). The beam is converging to the focal point of the lens so after refraction the beam is parallel. And we know that in order to have parallel beam (i.e. minimize divergence) we must have big waist:

In [39]: q_in.divergence.subs({z_r : waist2rayleigh(w_0_in,l)})
Out[39]:
  λ
─────
π⋅wⁱ₀

What about the minimal waist?

The predicted waist goes to zero for very big Z_i. But the gauss beam is a solution to the paraxial equation which is invalid for great curvatures of the beam front. We have the following limit:

In [41]: q_out.waist_aproximation_limit
Out[41]:
2⋅λ
───
 π

It would be interesting to know for what distance is this obtainable. The sad thing is that currently sympy can solve only very explicit equations. To rewrite it I do the following:

In [48]: waist_limit = q_out.waist_aproximation_limit

In [49]: equation = (w_0_out.args[0] / waist_limit)**2 - (1/w_0_out.args[1])**2

In [50]: equation
Out[50]:
                          2    2
    2            2       π ⋅wⁱ₀
- Zᵢ  + 2⋅Zᵢ - Zᵣ  - 1 + ───────
                             2
                          4⋅λ

In [51]:

In [52]: solve(equation, Z_in)
Out[52]:
⎡       ⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽         ⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎽⎤
⎢      ╱       2  2    2    2         ╱       2  2    2    2 ⎥
⎢    ╲╱  - 4⋅Zᵣ ⋅λ  + π ⋅wⁱ₀        ╲╱  - 4⋅Zᵣ ⋅λ  + π ⋅wⁱ₀  ⎥
⎢1 - ─────────────────────────, 1 + ─────────────────────────⎥
⎣               2⋅λ                            2⋅λ           ⎦

There are many more interesting questions but that’s all for now.

Advertisements
Leave a comment

Leave a Reply

Fill in your details below or click an icon to log in:

WordPress.com Logo

You are commenting using your WordPress.com account. Log Out / Change )

Twitter picture

You are commenting using your Twitter account. Log Out / Change )

Facebook photo

You are commenting using your Facebook account. Log Out / Change )

Google+ photo

You are commenting using your Google+ account. Log Out / Change )

Connecting to %s