Numerical analysis of eikonal equation

by   D. S. Kulyabov, et al.
RUDN University

The Maxwell equations have a fairly simple form. However, finding solutions of Maxwell's equations is an extremely difficult task. Therefore, various simplifying approaches are often used in optics. One such simplifying approach is to use the approximation of geometric optics. The approximation of geometric optics is constructed with the assumption that the wavelengths are small (short-wavelength approximation). The basis of geometric optics is the eikonal equation. The eikonal equation can be obtained from the wave equation (Helmholtz equation). Thus, the eikonal equation relates the wave and geometric optics. In fact, the eikonal equation is a quasi-classical approximation (the Wentzel-Kramers-Brillouin method) of wave optics. This paper shows the application of geometric methods of electrodynamics to the calculation of optical devices, such as Maxwell and Luneburg lenses. The eikonal equation, which was transformed to the ODE system by the method of characteristics, is considered. The resulting system is written for the case of Maxwell and Luneburg lenses.



page 1

page 2

page 3

page 4


A geometric optics ansatz-based plane wave method for two dimensional Helmholtz equations with variable wave numbers

In this paper we develop a plane wave type method for discretization of ...

Relative Schwarz Method

This paper proposes Relative Schwarz Method, which is a new domain decom...

Solitary water wave interactions for the Forced Korteweg-de Vries equation

The aim of this work is to study solitary water wave interactions betwee...

Numerical study of the stabilization of 1D locally coupled wave equations

In this paper, we study the numerical stabilization of a 1D system of tw...

A Flawed Dataset for Symbolic Equation Verification

Arabshahi, Singh, and Anandkumar (2018) propose a method for creating a ...
This week in AI

Get the week's most popular data science and artificial intelligence research sent straight to your inbox every Saturday.

I Introduction

In this article, we consider the approach to transform the eikonal equations to the ODE system. The first part of the article describes in detail all the mathematical calculations. In the second part we briefly describe Maxwell and Luneberg lenses, and explain the approach to their numerical modeling, which allows to obtain ray trajectories and wave fronts from sources of different shapes.

Ii Application of the characteristics method to the eikonal equation solution

ii.1 The eikonal equation

The eikonal equation can be obtained from Maxwell’s equations, written for the regions free of currents and charges, and under the condition of a time-changing harmonic electromagnetic field in a nonconducting isotropic medium born-wolf:principles_optics::en; stratton:1948::en; ll:2::en; ll:8::en

. In general, the eikonal equation is written as a partial differential equation of the first order:


is radius-vector,

is the boundary condition, is the refractive index of the medium. The function is the real scalar function with a physical meaning of time. It is also often called the eikonal function bruns:1895; klein:1901:eikonal.

For visualization of lens modeling results, we will consider their projection on the plane. In this case, the eikonal equation is reduced to the following two-dimensional form:


Using the method of characteristics, the eikonal equation can be transformed into an ODE system that can be solved by standard numerical methods.

ii.2 Characteristics of the eikonal equation

Let us briefly describe the method of characteristics jeong:2007:eikonal_parallel; kimmel:1998:geodesic; zhao:2004:fsm; beliakov:1996:numerical_luneburg; gremaud:2006:fsm_eikonal; bak:2010:fast-sweeping-method and the application of this method to the eikonal equation.

The partial differential equation of the following form is considered:


where , , and are sufficiently smooth functions. This equation is equivalent to the statement that a vector field with components , , is tangent to the surface , which has a normal vector with components . Accordingly, for this equation one can write the system of ODE, called equations of characteristics. It has the following form:

This ODE system reduces the solution of the partial differential equation of the first order to the solution of the ODE system of the first order.

To get the equations of characteristics for the eikonal equation one has to perform two steps. At the first stage, the equation should be converted to the form (2), and after that the ODE system may be written down. For the two-dimensional case the conversion of the eikonal equation to (2) is performed by replacing

In this case, the equation itself is converted to form:

A number of changes should be made.

After that the following system of equations is obtained:



Using this equality our expressions may be converted in

As a result:

Thus, the goal is achieved — the equation (1) is transformed in two equations of the form (2).

The characteristics for each equation may be written down:

So the ODE system of four equations with four functions: , , , , is derived:

The initial conditions:

Constants and are bonded by following relation . These constants may be presented as and . The initial conditions give a mathematical description of the source of the rays. For example, to model a point source, we need to fix the initial coordinates , and change the angle , which will set the angle of the beam exit from the source-point. To simulate the radiating surface, on the contrary, it is necessary to fix the angle and change the coordinates and .

Let us to find the relation between the parameter and the function . Since:



Due to the fact that , we obtain:

The solution of the equation is the function which implies that the parameter has a physical meaning of the signal propagation time from the point to the point

In polar coordinates, the eikonal equation has the following form::

and the corresponding system of ODEs will have the form:

The initial conditions:

Iii Numerical simulation of Luneburg and Maxwell lenses

Let’s consider the examples of lenses kulyabov:2017:sfm:geometrization_maxwell; kulyabov:2018:sfm:lens-calculations.

iii.1 Luneburg lens

Luneburg lens luneburg:1964; morgan:1958:luneberg_lens; lock:2008:luneburg_ray; lock:2008:luneburg_wave is a spherical lens of radius with the center at point (consider the projection on the plane Oxy) with a refractive index of the following form

where is the distance from the center of the lens to an arbitrary point in the plane. The formula implies that the coefficient continuously varies from to starting from the center of the lens and ending with its boundary. The refractive index of the medium outside the lens is constant and is equal to . Usually is equal to .

To solve the eikonal equation by the method of characteristics it is necessary to find partial derivatives of the function . For the case of Luneburg lens the partial derivatives are:

Outside the lens region derivatives are equal to .

iii.2 Maxwell fish eye lens

Maxwell fish eye lens maxwell:1854:fish-eye is also a spherical lens of radius with the center at point (consider the projection on the plane Oxy) with a refractive index of the following form:

To solve the eikonal equation by the method of characteristics it is necessary to find partial derivatives of the function . For the case of Maxwell lens partial derivatives have the form:

iii.3 Description of the numerical modeling

Julia programming language joshi:book:learning-julia is used to simulate the trajectories of rays through the Maxwell and Luneburg lenses. We use classical Runge–Kutta methods with constant step to solve the ODE system.

We carry on numerical modeling for lenses with a radius , the refractive index of the external medium , the center of the lens was placed in the point , the boundary region was set as the rectangle , , and . The point source was placed on the lens boundary at . values of the parameter have been taken from the interval , which allowed to simulate rays trajectories from a point source within an angle slightly smaller than . The parameter was changed within the interval.

Each parameter value sets new initial conditions for the ODE system. The process of numerical simulation consists in multiple solution of this system for different initial conditions. The numerical solution of the ODE system for a particular initial condition gives us a set of points , approximating the trajectories of a particular beam. After performing calculations for all the selected initial conditions, we obtain a set of rays. To visualize the rays, it is enough to depict each of the obtained numerical solutions. The result of the simulation can be seen in the Fig. 2 and Fig. 3 (the trajectories of the rays) and Fig. 2 and Fig. 4 (the wavefronts).

Figure 1: The trajectories of the rays in case of Maxwell’s lens for a point source and
Figure 2: The wavefronts for in case of Maxwell’s lens for a point source and
Figure 3: The trajectories of the rays in case of Luneburg lens for a point source and .
Figure 4: The wavefronts in case of Luneburg lens for a point source and .

To visualize the wave fronts with the resulting numerical data it is necessary to carry out additional recalculations. From each numerical solution, we must select points that correspond to a specific point in time .

The use of a numerical method with a fixed step gives an advantage, since each numerical solution will be obtained for the same uniform grid .

Iv Conclusion

The paper presents the description of the numerical solution of the eikonal equation for the case of Luneburg and Maxwell lenses. The results are visualized as trajectories of rays passing through lenses and as fronts of electromagnetic waves.

The publication has been prepared with the support of the ‘‘RUDN University Program 5-100’’ and funded by Russian Foundation for Basic Research (RFBR) according to the research project No 19-01-00645.