Tutorial 3a: Simple Optimization
This tutorial shows how to optimize optical systems in Optiland. In this example, we will optimize a singlet.
[1]:
import numpy as np
from optiland import analysis, optic, optimization, wavefront
We first define a simple singlet. We will make it double convex and will define a few fields and wavelengths.
[2]:
lens = optic.Optic()
# add surfaces
lens.surfaces.add(index=0, thickness=np.inf)
lens.surfaces.add(index=1, thickness=5, radius=100, is_stop=True, material="N-SF11")
lens.surfaces.add(index=2, thickness=59, radius=-1000)
lens.surfaces.add(index=3)
# add aperture
lens.set_aperture(aperture_type="EPD", value=25)
# add field
lens.fields.set_type(field_type="angle")
lens.fields.add(y=0.0)
lens.fields.add(y=0.7)
lens.fields.add(y=1.0)
# add wavelength
lens.wavelengths.add(value=0.4861)
lens.wavelengths.add(value=0.5876, is_primary=True)
lens.wavelengths.add(value=0.6563)
lens.update_paraxial()
[3]:
lens.draw()
[3]:
(<Figure size 1000x400 with 1 Axes>, <Axes: xlabel='Z [mm]', ylabel='Y [mm]'>)
Now, we need to define the optimization problem, which contains all information about the variables and operands.
In particular, we will optimize for minimal optical path difference for each field. We will also add an operand to specify the focal length at 100 mm. As for variables, we will let the lens radii vary and the distance from the lens to the image plane. Once defined, we can view all information related to the optimization problem using the ‘info’ method.
There are many operand types available, which can be seen in optiland.operand. Variable options can be seen in optiland.variable,
[4]:
problem = optimization.OptimizationProblem()
for wave in lens.wavelengths.get_wavelengths():
for Hx, Hy in lens.fields.get_field_coords():
input_data = {
"optic": lens,
"Hx": Hx,
"Hy": Hy,
"num_rays": 3,
"wavelength": wave,
"distribution": "gaussian_quad",
}
problem.add_operand(
operand_type="OPD_difference",
target=0,
weight=1,
input_data=input_data,
)
problem.add_operand(
operand_type="f2",
target=100,
weight=10,
input_data={"optic": lens},
)
problem.add_variable(lens, "thickness", surface_number=2, min_val=0, max_val=1000)
problem.add_variable(lens, "radius", surface_number=1, min_val=-1000, max_val=1000)
problem.add_variable(lens, "radius", surface_number=2, min_val=-1000, max_val=1000)
problem.info()
╒════╤════════════════════════╤═══════════════════╕
│ │ Merit Function Value │ Improvement (%) │
╞════╪════════════════════════╪═══════════════════╡
│ 0 │ 209255 │ 0 │
╘════╧════════════════════════╧═══════════════════╛
╒════╤════════════════╤══════════╤══════════════╤══════════════╤══════════╤═════════╤═════════╤════════════════╕
│ │ Operand Type │ Target │ Min. Bound │ Max. Bound │ Weight │ Value │ Delta │ Contrib. [%] │
╞════╪════════════════╪══════════╪══════════════╪══════════════╪══════════╪═════════╪═════════╪════════════════╡
│ 0 │ OPD difference │ 0 │ │ │ 1 │ 253.838 │ 253.838 │ 30.79 │
│ 1 │ OPD difference │ 0 │ │ │ 1 │ 84.642 │ 84.642 │ 3.42 │
│ 2 │ OPD difference │ 0 │ │ │ 1 │ 84.678 │ 84.678 │ 3.43 │
│ 3 │ OPD difference │ 0 │ │ │ 1 │ 216.741 │ 216.741 │ 22.45 │
│ 4 │ OPD difference │ 0 │ │ │ 1 │ 72.277 │ 72.277 │ 2.5 │
│ 5 │ OPD difference │ 0 │ │ │ 1 │ 72.306 │ 72.306 │ 2.5 │
│ 6 │ OPD difference │ 0 │ │ │ 1 │ 196.482 │ 196.482 │ 18.45 │
│ 7 │ OPD difference │ 0 │ │ │ 1 │ 65.523 │ 65.523 │ 2.05 │
│ 8 │ OPD difference │ 0 │ │ │ 1 │ 65.548 │ 65.548 │ 2.05 │
│ 9 │ f2 │ 100 │ │ │ 10 │ 116.082 │ 16.082 │ 12.36 │
╘════╧════════════════╧══════════╧══════════════╧══════════════╧══════════╧═════════╧═════════╧════════════════╛
╒════╤═════════════════╤═══════════╤═════════╤══════════════╤══════════════╕
│ │ Variable Type │ Surface │ Value │ Min. Bound │ Max. Bound │
╞════╪═════════════════╪═══════════╪═════════╪══════════════╪══════════════╡
│ 0 │ thickness │ 2 │ 59 │ 0 │ 1000 │
│ 1 │ radius │ 1 │ 100 │ -1000 │ 1000 │
│ 2 │ radius │ 2 │ -1000 │ -1000 │ 1000 │
╘════╧═════════════════╧═══════════╧═════════╧══════════════╧══════════════╛
We now need to define an optimizer, which will take the optimization problem as an input. We choose the standard “generic” optimizer, but there are many other options, including a least squares optimizer.
Note that the generic optimizer simply utilizes scipy.optimize.minimize.
[5]:
optimizer = optimization.OptimizerGeneric(problem)
# We can also use a least squares optimizer as follows:
# optimizer = optimization.LeastSquares(problem)
We can now simply call the optimizer.
[6]:
res = optimizer.optimize()
If we now view the optimization problem info, we can see that we’ve significantly improved the system (>99.9%)
[7]:
problem.info()
╒════╤════════════════════════╤═══════════════════╕
│ │ Merit Function Value │ Improvement (%) │
╞════╪════════════════════════╪═══════════════════╡
│ 0 │ 102.608 │ 99.951 │
╘════╧════════════════════════╧═══════════════════╛
╒════╤════════════════╤══════════╤══════════════╤══════════════╤══════════╤═════════╤═════════╤════════════════╕
│ │ Operand Type │ Target │ Min. Bound │ Max. Bound │ Weight │ Value │ Delta │ Contrib. [%] │
╞════╪════════════════╪══════════╪══════════════╪══════════════╪══════════╪═════════╪═════════╪════════════════╡
│ 0 │ OPD difference │ 0 │ │ │ 1 │ 6.608 │ 6.608 │ 42.55 │
│ 1 │ OPD difference │ 0 │ │ │ 1 │ 2.224 │ 2.224 │ 4.82 │
│ 2 │ OPD difference │ 0 │ │ │ 1 │ 2.247 │ 2.247 │ 4.92 │
│ 3 │ OPD difference │ 0 │ │ │ 1 │ 3.046 │ 3.046 │ 9.04 │
│ 4 │ OPD difference │ 0 │ │ │ 1 │ 0.997 │ 0.997 │ 0.97 │
│ 5 │ OPD difference │ 0 │ │ │ 1 │ 0.979 │ 0.979 │ 0.93 │
│ 6 │ OPD difference │ 0 │ │ │ 1 │ 5.568 │ 5.568 │ 30.21 │
│ 7 │ OPD difference │ 0 │ │ │ 1 │ 1.84 │ 1.84 │ 3.3 │
│ 8 │ OPD difference │ 0 │ │ │ 1 │ 1.823 │ 1.823 │ 3.24 │
│ 9 │ f2 │ 100 │ │ │ 10 │ 100.011 │ 0.011 │ 0.01 │
╘════╧════════════════╧══════════╧══════════════╧══════════════╧══════════╧═════════╧═════════╧════════════════╛
╒════╤═════════════════╤═══════════╤═══════════╤══════════════╤══════════════╕
│ │ Variable Type │ Surface │ Value │ Min. Bound │ Max. Bound │
╞════╪═════════════════╪═══════════╪═══════════╪══════════════╪══════════════╡
│ 0 │ thickness │ 2 │ 96.0069 │ 0 │ 1000 │
│ 1 │ radius │ 1 │ 85.5752 │ -1000 │ 1000 │
│ 2 │ radius │ 2 │ -922.283 │ -1000 │ 1000 │
╘════╧═════════════════╧═══════════╧═══════════╧══════════════╧══════════════╛
And when we view the lens, we see that it properly focuses the light on our image plane.
[8]:
lens.draw()
[8]:
(<Figure size 1000x400 with 1 Axes>, <Axes: xlabel='Z [mm]', ylabel='Y [mm]'>)
Let’s look at a few other common analyes:
[9]:
fan = analysis.RayFan(lens)
fan.view()
[9]:
(<Figure size 1000x999 with 6 Axes>,
[<Axes: title={'center': 'Hx: 0.000, Hy: 0.000'}, xlabel='$P_y$', ylabel='$\\epsilon_y$ (mm)'>,
<Axes: title={'center': 'Hx: 0.000, Hy: 0.000'}, xlabel='$P_x$', ylabel='$\\epsilon_x$ (mm)'>,
<Axes: title={'center': 'Hx: 0.000, Hy: 0.700'}, xlabel='$P_y$', ylabel='$\\epsilon_y$ (mm)'>,
<Axes: title={'center': 'Hx: 0.000, Hy: 0.700'}, xlabel='$P_x$', ylabel='$\\epsilon_x$ (mm)'>,
<Axes: title={'center': 'Hx: 0.000, Hy: 1.000'}, xlabel='$P_y$', ylabel='$\\epsilon_y$ (mm)'>,
<Axes: title={'center': 'Hx: 0.000, Hy: 1.000'}, xlabel='$P_x$', ylabel='$\\epsilon_x$ (mm)'>])
[10]:
spot = analysis.SpotDiagram(lens)
spot.view()
[10]:
(<Figure size 1200x400 with 3 Axes>,
[<Axes: title={'center': 'Hx: 0.000, Hy: 0.000'}, xlabel='X (mm)', ylabel='Y (mm)'>,
<Axes: title={'center': 'Hx: 0.000, Hy: 0.700'}, xlabel='X (mm)', ylabel='Y (mm)'>,
<Axes: title={'center': 'Hx: 0.000, Hy: 1.000'}, xlabel='X (mm)', ylabel='Y (mm)'>])
[11]:
opd = wavefront.OPD(lens, field=(0, 1), wavelength=0.55)
opd.view(projection="2d", num_points=512)
[11]:
(<Figure size 700x550 with 2 Axes>,
<Axes: title={'center': 'OPD Map: RMS=2.100 waves'}, xlabel='Pupil X', ylabel='Pupil Y'>)
