Optimization

Advanced Optimization

Use global search and multi-configuration solvers for complex designs.

AdvancedOptimizationNumPy backend15 min read

Introduction

Building on tutorial 5a, this tutorial demonstrates further optimization functionalities including:

Various operand types (paraxial, real ray-based, aberrations)
Various variable types (radii, thickness, conic)
Global optimization (e.g. dual annealing, differential evolution, SHGO, basin-hopping)

Core concepts used

optimization.DifferentialEvolution(...)

A stochastic global optimizer that evolves a population of solutions, finding the global optimum across the entire search space.

real_y_intercept

Targets the exact coordinate where a ray strikes the image plane, used here to enforce a specific distortion profile.

conic

Allows the optimizer to change the surface conic constant, enabling aspheric shapes during the search.

workers=-1

Enables multi-core parallelization, drastically reducing the time required for a global optimization run.

Step-by-step build

Import numpy, analysis, optic, and optimization

python

import numpy as np

from optiland import analysis, optic, optimization

Define the doublet lens with aperture, fields, and wavelengths

Let's define a non-optimized doublet, which we will optimize. Our goal in optimization will be to:

Control distortion via real ray intersection operands
Minimize spherical aberration
Minimize the spot size for each field

Maintain a focal length of 50 mm

python

lens = optic.Optic()

 lens.surface_group.surfaces = []

 lens.add_surface(index=0, radius=np.inf, thickness=np.inf)
 lens.add_surface(index=1, radius=50, thickness=5, material="N-BK7")
 lens.add_surface(index=2, radius=-500, thickness=10)
 lens.add_surface(index=3, radius=500, thickness=5, material="SK16", is_stop=True)
 lens.add_surface(index=4, radius=-50, thickness=35)
 lens.add_surface(index=5)

 # add aperture
 lens.set_aperture(aperture_type="EPD", value=8)

 # add field
 lens.set_field_type(field_type="angle")
 lens.add_field(y=0)
 lens.add_field(y=5.6)
 lens.add_field(y=8)

 # add wavelength
 lens.add_wavelength(value=0.4861327)
 lens.add_wavelength(value=0.5875618, is_primary=True)
 lens.add_wavelength(value=0.6562725)

 lens.update_paraxial()

Draw the initial lens layout

python

lens.draw()

Add ray intercept, Seidel, RMS spot, and focal length operands

Let's define the optimization problem. First, let's add the operands:

python

problem = optimization.OptimizationProblem()

# add real ray operands - the intersection location is f*tan(theta), where f is the
# focal length and theta is the field half angle
input_data = {
 "optic": lens,
 "surface_number": 5,
 "Hx": 0,
 "Hy": 0.7,
 "Px": 0,
 "Py": 0,
 "wavelength": 0.5875618,
}
problem.add_operand(
 operand_type="real_y_intercept",
 target=4.903,
 weight=1,
 input_data=input_data,
)

input_data = {
 "optic": lens,
 "surface_number": 5,
 "Hx": 0,
 "Hy": 1,
 "Px": 0,
 "Py": 0,
 "wavelength": 0.5875618,
}
problem.add_operand(
 operand_type="real_y_intercept",
 target=7.027,
 weight=1,
 input_data=input_data,
)

# add operand to minimize spherical aberration - we will minimize the first
# seidel aberration (i.e. #1)
input_data = {"optic": lens, "seidel_number": 1}
problem.add_operand(operand_type="seidel", target=0, weight=1, input_data=input_data)

# add RMS spot size operand - let's minimize the spot size for each field at the
# primary wavelength we choose a 'uniform' distribution, so the number of rays actually
# means the rays on one axis. So, we acually trace ≈16^2 rays here.
for field in lens.fields.get_field_coords():
 input_data = {
     "optic": lens,
     "surface_number": 5,
     "Hx": field[0],
     "Hy": field[1],
     "num_rays": 16,
     "wavelength": 0.5875618,
     "distribution": "uniform",
 }
 problem.add_operand(
     operand_type="rms_spot_size",
     target=0,
     weight=10,
     input_data=input_data,
 )

# add focal length target
problem.add_operand(operand_type="f2", target=50, weight=1, input_data={"optic": lens})

Add thickness, radius, and conic variables and inspect the problem

Now, we can add the variables. Let's let one of the surfaces' conic constant vary. Lastly, we print the optimization problem information overview.

python

# thickness between the lenses and distance to image plane
problem.add_variable(lens, "thickness", surface_number=2, min_val=3, max_val=30)
problem.add_variable(lens, "thickness", surface_number=4, min_val=0, max_val=100)

# radii of all lenses
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.add_variable(lens, "radius", surface_number=3, min_val=-1000, max_val=1000)
problem.add_variable(lens, "radius", surface_number=4, min_val=-1000, max_val=1000)

# conic constant of first lens surface
problem.add_variable(lens, "conic", surface_number=1, min_val=-10, max_val=10)

# print optimization problem info
problem.info()

Run Differential Evolution global optimization

Let's optimize our system. There are many optimizers available, but we want to demonstrate a global optimization method: Differential Evolution. Optiland utilizes the scipy implementation of this optimization method.

This method can take some time. To speed up the calculation, we parallize the calculation and use all available CPU cores. We specify this by supplying -1 to the workers argument.

python

optimizer = optimization.DifferentialEvolution(problem)
# optimizer = optimization.SHGO(problem)  # note SHGO requires bounds
# optimizer = optimization.BasinHopping(problem)  # BasinHopping requires no bounds

res = optimizer.optimize(maxiter=256, disp=False, workers=-1)

Draw the optimized lens layout

Let's see how the lens has changed after optimization:

python

lens.draw()

Review optimization improvement with problem.info()

By viewing the problem info, we can see that there was a greater than 99% improvement in performance, based on the defined operands.

python

problem.info()

Plot the spot diagram after optimization

Let's perform a few analyses to check the final performance of our lens.

python

spot = analysis.SpotDiagram(lens)
spot.view()

Plot the grid distortion

python

grid = analysis.GridDistortion(lens)
grid.view()

Print the Seidel aberration coefficients

python

print("Seidel Aberrations:")
for k, seidel in enumerate(lens.aberrations.seidels()):
 print(f"\tS{k + 1}: {seidel:.3e}")

Show full code listing

python

import numpy as np

from optiland import analysis, optic, optimization

lens = optic.Optic()

lens.surface_group.surfaces = []

lens.add_surface(index=0, radius=np.inf, thickness=np.inf)
lens.add_surface(index=1, radius=50, thickness=5, material="N-BK7")
lens.add_surface(index=2, radius=-500, thickness=10)
lens.add_surface(index=3, radius=500, thickness=5, material="SK16", is_stop=True)
lens.add_surface(index=4, radius=-50, thickness=35)
lens.add_surface(index=5)

# add aperture
lens.set_aperture(aperture_type="EPD", value=8)

# add field
lens.set_field_type(field_type="angle")
lens.add_field(y=0)
lens.add_field(y=5.6)
lens.add_field(y=8)

# add wavelength
lens.add_wavelength(value=0.4861327)
lens.add_wavelength(value=0.5875618, is_primary=True)
lens.add_wavelength(value=0.6562725)

lens.update_paraxial()

lens.draw()

problem = optimization.OptimizationProblem()

# add real ray operands - the intersection location is f*tan(theta), where f is the
# focal length and theta is the field half angle
input_data = {
  "optic": lens,
  "surface_number": 5,
  "Hx": 0,
  "Hy": 0.7,
  "Px": 0,
  "Py": 0,
  "wavelength": 0.5875618,
}
problem.add_operand(
  operand_type="real_y_intercept",
  target=4.903,
  weight=1,
  input_data=input_data,
)

input_data = {
  "optic": lens,
  "surface_number": 5,
  "Hx": 0,
  "Hy": 1,
  "Px": 0,
  "Py": 0,
  "wavelength": 0.5875618,
}
problem.add_operand(
  operand_type="real_y_intercept",
  target=7.027,
  weight=1,
  input_data=input_data,
)

# add operand to minimize spherical aberration - we will minimize the first
# seidel aberration (i.e. #1)
input_data = {"optic": lens, "seidel_number": 1}
problem.add_operand(operand_type="seidel", target=0, weight=1, input_data=input_data)

# add RMS spot size operand - let's minimize the spot size for each field at the
# primary wavelength we choose a 'uniform' distribution, so the number of rays actually
# means the rays on one axis. So, we acually trace ≈16^2 rays here.
for field in lens.fields.get_field_coords():
  input_data = {
      "optic": lens,
      "surface_number": 5,
      "Hx": field[0],
      "Hy": field[1],
      "num_rays": 16,
      "wavelength": 0.5875618,
      "distribution": "uniform",
  }
  problem.add_operand(
      operand_type="rms_spot_size",
      target=0,
      weight=10,
      input_data=input_data,
  )

# add focal length target
problem.add_operand(operand_type="f2", target=50, weight=1, input_data={"optic": lens})

# thickness between the lenses and distance to image plane
problem.add_variable(lens, "thickness", surface_number=2, min_val=3, max_val=30)
problem.add_variable(lens, "thickness", surface_number=4, min_val=0, max_val=100)

# radii of all lenses
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.add_variable(lens, "radius", surface_number=3, min_val=-1000, max_val=1000)
problem.add_variable(lens, "radius", surface_number=4, min_val=-1000, max_val=1000)

# conic constant of first lens surface
problem.add_variable(lens, "conic", surface_number=1, min_val=-10, max_val=10)

# print optimization problem info
problem.info()

optimizer = optimization.DifferentialEvolution(problem)
# optimizer = optimization.SHGO(problem)  # note SHGO requires bounds
# optimizer = optimization.BasinHopping(problem)  # BasinHopping requires no bounds

res = optimizer.optimize(maxiter=256, disp=False, workers=-1)

lens.draw()

problem.info()

spot = analysis.SpotDiagram(lens)
spot.view()

grid = analysis.GridDistortion(lens)
grid.view()

print("Seidel Aberrations:")
for k, seidel in enumerate(lens.aberrations.seidels()):
  print(f"\tS{k + 1}: {seidel:.3e}")

Conclusions

We showed several different operands and variables that can be used in optimization.
We demonstrated the use of differential evolution for optimization
While the lens itself improved in performance, there is still room for improvement. Different operands, variables, or sequences of optimization steps would result in different, and possibly better, results.

Next tutorials

NextOptimization Case StudyFollow the end-to-end design of a real-world imaging system.RelatedSimple OptimizationThe foundational local-search optimization concepts used here.

Original notebook: Tutorial_5b_Advanced_Optimization.ipynb on GitHub · ReadTheDocs