# MieLab: A Software Tool to Perform Calculations on the Scattering of Electromagnetic Waves by Multilayered Spheres

Authors: Ovidio Peña-Rodríguez; Pedro Pablo González Pérez; Umapada Pal
Journal: International Journal of Spectroscopy
 (2011)

---

## Abstract

In this paper, we present MieLab, a free computational package for simulating the scattering of electromagnetic radiation by multilayered spheres or an ensemble of particles with normal size distribution. It has been designed as a virtual laboratory, including a friendly graphical user interface (GUI), an optimization algorithm (to fit the simulations to experimental results) and scripting capabilities. The paper is structured in five different sections: the introduction is a perspective on the importance of the software for the study of scattering of light scattering. In the second section, various approaches used for modeling the scattering of electromagnetic radiation by small particles are discussed. The third and fourth sections are devoted to provide an overview of MieLab and to describe the main features of its architectural model and functional behavior, respectively. Finally, several examples are provided to illustrate the main characteristics of the software.


---

## Body

## 1. Introduction

The scattering and  absorption of electromagnetic (EM) radiation by small particles are important issues for a surprisingly large number of disciplines: solid-state physics, electrical engineering, meteorology, chemistry, biophysics, astronomy, and medicine are just a few of them [1, 2]. The optical properties of such small particles are closely related to their physical and geometrical parameters such as composition, size, shape, and embedding medium [3]. Hence, a quantitative knowledge of the EM scattering as a function of the particle parameters is required if we want to understand the optical phenomena occurring in them, as well as for developing remote sensing and laboratory techniques for particle characterization. However, obtaining such quantitative knowledge is not an easy task in many cases. Fortunately, in spite of the enormous amount of particle shapes and configurations which can be found in nature or in laboratory synthesis processes, it is often sufficient to study simpler systems not only for understanding the elementary phenomena occurring in them but also for practical applications.An excellent example of such kind of particle, which is both simple and of great interest for many applications, is the multilayered sphere (Figure1). In the field of atmospheric science, light absorption by aerosols has a heating effect in the atmosphere, which is of great interest to study several climatic effects. Their optical properties can be successfully modeled using a core-shell type particle [4, 5]. Moreover, in the field of nanotechnology, applications of multilayered spheres are multiple in catalysis [6], optoelectronics [7], biomedical diagnosis [8], and medicine [9–11], among others.Figure 1
(Color online) Schematic representation of the multilayered sphere geometries used for calculating the optical properties in MieLab.A problem that occurs while performing electromagnetic calculations, even when exact solutions exist, is that they are usually quite complicated as infinite series composed of combinations of badly behaved Bessel functions [1, 12] are involved. These difficulties make the integration of flexible and robust computational techniques (software applications) essential. Thus, it is not surprising that the number of existing calculation methods and computer codes related to this field is enormous [13] and continues to increase every day. Unfortunately, despite this variety, there is a marked trend (of course with some exceptions), where most of the freely available applications require not only the advanced knowledge of electromagnetism but also adequate programming skills to use them. Moreover, often the more user-friendly tools are the commercial ones, requiring licenses that are prohibitively expensive for many researchers, especially from the developing countries. Trying to help to alleviate this problem, we have set the ambitious goal of developing a tool capable of performing accurate calculations, which at the same time, would be freely available for everyone and would have a graphical user-friendly interface (GUI).In this paper, we present MieLab, a software tool released as free software, which is based on the implementation [14] of Yang’s algorithm for multilayered spheres [12] with several added features such as the possibility of calculating the optical response of an ensemble of independent (noninteracting) particles with a normal size distribution, a graphical user interface (GUI), an optimization algorithm (to fit the simulations to experimental results), and scripting capabilities.


## 2. Scattering of Electromagnetic  Radiation by  Multilayered Spheres

The methods used for computing the scattering properties of particles theoretically can be roughly divided in two groups: exact and approximated techniques. In general, the exact methods are more accurate than the approximated ones, but they are much more restricted with respect to the kind of particles they are suited for. Fortunately, in the case of spherically symmetric particles, the exact methods can be used to generate accurate results for essentially any size parameter and relative refractive index values. Therefore, in this paper, we will restrict ourselves to analyze only the exact methods in detail.All exact techniques for calculating electromagnetic scattering are based on solving the differential Maxwell equations or their integral counterparts in the time or frequency domain. To find an analytical solution, one has to solve the vector wave equations for the time-dependent electric fields outside and inside the particle, using the separation of variables technique. Then, the incident and internal fields are expanded in wave functions which are regular inside the scatterer, whereas the external scattered field is expanded in wave functions that behave as outgoing waves at infinity. Finally, using the requirement of continuity of the tangential component of the electric and magnetic fields at the particle boundary, the unknown coefficients in the internal and scattered field expansions can be determined from the known expansion coefficients of the incident field [12].Unfortunately, the separation of variables technique generates only manageable solutions in a few simple cases. The first of such exact solutions was developed by Mie [15] for spherical particles. This model was later extended first for the coated sphere [16] and then for the radially inhomogeneous (i.e., multilayered) sphere [17–20]. Similar solutions have been found for homogeneous infinite circular cylinders [21], infinite elliptical cylinders [22], and homogeneous and core-mantle spheroids [23].In the case of multilayered spheres, the expressions of the scattering coefficients, needed to compute the solution, are apparently simple and straightforward to implement. However, they involve several numerical difficulties, which make most of the existing algorithms [24, 25] inapplicable to several extreme cases [26, 27]. More recently, Yang [12] has developed an improved recursive algorithm that circumvents most of the numerical problems present in the previous algorithms. Based on Yang’s algorithm, we have developed a computer implementation which has been shown to be more accurate than the existing similar programs [14] and constitutes the basis of MieLab.


## 3. Overview of MieLab

### 3.1. Scattering Calculation

In this section, we will briefly discuss the most important aspects of the scattering algorithm implemented in MieLab. The interested reader can consult two additional papers: Yang [12] where the algorithm was originally presented and the paper containing the full description of our implementation [14]. In short, the analytical solution for the light scattering by a multilayered sphere can be found by expressing the EM field inside each layer of the sphere as a linear combination of the inward- and outward-travelling waves. Each layer is characterized by a size parameter xl=2πnmrl/λ=krl and a relative refractive index ml=nl/nm, l=1,2,…,L, where λ is the wavelength of the incident wave in vacuum, rl is the outer radius of the 1th layer, k is the propagation constant, and nm and nl are the refractive index of the medium outside the particle and its 1th component, respectively. For an L-layered sphere, there are 2L unknown field coefficients. Then, the scattering coefficients can be calculated by matching the tangential components of the electromagnetic fields at each interface and solving the obtained linear system of 2L independent equations. The final coefficients in the series are the scattering coefficients [12]
(1)an=anL+1=[Hna(mLxL)/mL+n/xL]ψn(xL)-ψn-1(xL)[Hna(mLxL)/mL+n/xL]ζn(xL)-ζn-1(xL),bn=bnL+1=[mLHnb(mLxL)+n/xL]ψn(xL)-ψn-1(xL)[mLHnb(mLxL)+n/xL]ζn(xL)-ζn-1(xL),
where ψn and ζn are the Riccati-Bessel functions, as defined by Kai and Massoli [27] and the determinants Hna and Hnb can be calculated by means of the following recursive expressions [12]
(2)Hna(m1x1)=Dn(1)(m1x1),Hna(mlxl)=G2Dn(1)(mlxl)-Qn(l)G1Dn(3)(mlxl)G2-Qn(l)G1,l=2,…,L,Hnb(m1x1)=Dn(1)(m1x1),Hnb(mlxl)=G̃2Dn(1)(mlxl)-Qn(l)G̃1Dn(3)(mlxl)G̃2-Qn(l)G̃1,l=2,…,L,Dn(1)(z)=ψn′(z)ψn(z),Dn(3)(z)=ζn′(z)ζn(z),Qn(l)=ψn(mlxl-1)/ζn(mlxl-1)ψn(mlxl)/ζn(mlxl),G1=mlHna(ml-1xl-1)-ml-1Dn(1)(mlxl-1),G2=mlHna(ml-1xl-1)-ml-1Dn(3)(mlxl-1),G̃1=ml-1Hnb(ml-1xl-1)-mlDn(1)(mlxl-1),G̃2=ml-1Hnb(ml-1xl-1)-mlDn(3)(mlxl-1).Once the scattering coefficients are obtained, we can calculate the extinction, scattering, radiation pressure, and radar backscattering efficiency factors (Qext, Qsca, Qpr, and Qbk) from the following relations [1]:(3)Qext=2xL2∑n=1∞(2n+1)Re{an+bn},Qsca=2xL2∑n=1∞(2n+1)(|an|2+|bn|2),Qpr=Qext-gQsca=Qext-4xL2{∑n=1∞[n(n+2)n+1Re{anan+1*+bnbn+1*}+2n+1n(n+1)Re{anbn*}]∑n=1∞[n(n+2)n+1Re{anan+1*+bnbn+1*}},Qbk=1xL2|∑n=1∞(2n+1)(-1)n(an-bn)|2.
Finally, from the previously calculated quantities, it is also possible to derive the absorption efficiency factor (Qabs), single scattering albedo (Λ) and asymmetry parameter (g)
(4)Qabs=Qext-Qsca,Λ=QscaQext,g=(Qext-Qpr)Qsca.


### 3.2. Optimization

The bounded, limited memory, Broyden-Fletcher-Goldfarb-Shanno (L-BFGS-B) algorithm [28, 29] was used in MieLab for the optimization process, together with the chi-square distribution, as the objective function
(5)χ2(r1,σ1,t2,σ2,…,tL,σL,nm,c)=1N∑i=1N(〈OD〉i(λi,r1,σ1,t2,σ2,…,tL,σL,nm,c)-ODexp⁡i)2ODexp⁡i,
where N is the number of points used for the fit, nm is the refractive index of the matrix, ODexpi are the values of the experimental optical density, while r1, ti, and σj are the radius of the first layer, the thickness of the ith layer (i=2,…,L), and the standard deviation of one of the two previous numbers for the jth layer (j=1,…,L), respectively. Then, starting with the given initial values of the variables to optimize (Xk=r1,σ1,t2,σ2,…,tL,σL,nm) and the experimental optical density, several optimization steps are performed iteratively. In each of the steps, the value of the theoretical optical density [30],  X2 and its gradient G are calculated, and subsequently, the L-BFGS-B algorithm is applied to obtain a new “optimized” value for all of the variables Xk. Lower and/or upper limits can be set for all the variables, in order to ensure that their optimized values have physical meaning (e.g., positive radii, etc.). The optimization process is finished when one of the following criteria is met:
(6)‖Gk‖≤εG,|(χ2)(i+1)-(χ2)(i)|≤εF⋅max⁡{(χ2)(i),(χ2)(i+1),1},|Xk(i+1)-Xk(i)|≤εX,
where εG, εF, and εX are positive numbers which define a precision of search, ∥·∥ means Euclidian norm, Gk, gradient projection onto a variable Xk, and the superindex represents the iteration number. Additionally, the optimization can be finished after a certain number of iterations. All these parameters controlling the fit can be tuned by the user, until (s)he finds the combination that fits best to her/his needs.


### 3.3. Optical Properties of Materials

At this point, MieLab allows two types of optical properties (more precisely, refractive index): a constant value and a table of values (“Element” in MieLab’s terminology). When a constant refractive index is selected, then the values ofn and k are entered directly into the text boxes corresponding to the layer being edited. On the other hand, if we choose to use the optical properties of an “Element”, then we have to select one of the available materials in the corresponding combo box. MieLab obtains the optical properties of the materials by reading, during the initialization, all the files with the extension “txt” found in the folder “RefIndex” (which in turn is located inside the folder where the program was installed). These are simply text files and the user can easily add a new or modify existing ones. The first lines of the file Ag.txt, which corresponds to “Silver”, are shown below, in order to exemplify the structure with which the data is stored:SilverPhys. Rev. B, 6 (1972) 437018.99 1.72331E-02 1.39e6 5.85e280.64 0.24 14.080.77 0.15 11.85.As can be seen, the first two lines of the file are the name of the “Element” (the same which is shown in the combo box) and a comment (in this case the reference to the source of the data), respectively. The third line can only take the values 0 or 1, where “1” means “apply size correction” [31] and “0” means the opposite. The next line depends on whether or not the size correction is to be applied; in the former case, it is expected to contain (in this order) the plasma frequency (ωp) in eV, the bulk damping constant (Γ) in eV, the Fermi velocity (vF) in m/s, and the atomic density in atoms/m3, while in the latter, only the last parameter is expected. The values of ωp, Γ and vF are used to apply the size correction [31], while the atomic density is necessary to calculate the optical density [30]. Finally, from the fifth line onwards, the actual values of the material refractive index are displayed. The first column presents the photon energy in eV, while the second and third ones present the real and imaginary parts of the complex refractive index, respectively.It must be said that although we believe that the options currently offered by MieLab for representing the optical properties of the materials are “good enough” for most applications, this is one of the areas we have selected to improve in future releases. In particular, we hope to add the possibility of using dispersion laws (Cauchy, Tauc-Lorentz, etc.) to represent the refractive index of materials.


### 3.4. Scripting Engine

RemObject Software’s Pascal Script (http://www.remobjects.com/ps.aspx), a free Pascal interpreter, which is compatible with both Delphi and Lazarus, has been used as the script engine in Mielab. Therefore, the scripts are written using a dialect of Object Pascal. When writing a script, in addition to the usual possibilities of Pascal, there are procedures to programmatically control almost all the MieLab’s functions. In addition, the object ActiveParticle (of type TLayeredSphere) is available to describe the geometry to be used in the calculations. A more detailed analysis of this functionality is beyond the scope of this paper. For more details, the reader can consult the sample script described below and MieLab’s Script Reference Guide, which is available on the website (http://scattering.sourceforge.net/ml_guide.html).


### 3.5. Software Architectural Model

As can be seen from the Figure2, MieLab architectural model is defined in terms of the following three main components:(i)
a shared data structure, namedworkspace, through which the components communicate with each other and the problem solution emerges,(ii)
a set of autonomous components, namedalgorithm objects, which contains specific knowledge and skills about the problem domain. Algorithm objects encapsulate the algorithms and models required for the elaboration of the data on the workspace,(iii)
a set of interactive components, namedinterface objects, for allowing the interaction of human users.Figure 2
(Color online) general overview of MieLab’s architectural model.Theworkspace contains all the information about the particle (e.g., composition, size, shape, refractive index, etc.), necessary to model the scattering of electromagnetic radiation and is organized into one or more hierarchies (depending on the application). The information at every level of a hierarchy represents partial solutions to the problem (e.g., scattering values, objective function values, results of the optimization algorithm, etc.). At every level, there is a “specific ontology” that describes the information. The knowledge and skills necessary for the scattering problem is distributed among the algorithm objects, which transform data from a hierarchy level to the same level or to the other. The algorithm objects execute these transformations using models, algorithmic procedures, or heuristic rules. Algorithm for calculating scattering, algorithm for calculating objective function, and optimization algorithm are some examples.The opportunistic reasoning for characterizing this organization is expressed in the following terms: at each step of the solution, thealgorithm object to be applied is determined dynamically, which produces an incremental generation of information in the workspace. The selection of a specific algorithm object is based on the status of the solution on the level and on the existence of other algorithm objects capable of improving the current state of the solution.As can be seen from the Figure2, the infrastructure also provides the tools for allowing the interaction with the system of human experts (interface objects) so that they can dynamically add knowledge on the scattering problem (e.g., modifying the particle parameters) or feedback on the results provided by the algorithms. In this way, the users can control the entire flow which gives rise to the final solution of the problem.


### 3.6. MieLab’s Graphical User Interface

Now, we describe MieLab’s graphical user interface, which is structured by four major components: a main menu, a toolbar, a drawing canvas, and a parameter panel. Figure3 shows MieLab’s graphical interface, integrating the previous components. The role of each of them is explained below.Figure 3
(Color online) Screenshot of MieLab’s graphical user interface.The main menu contains the “File” and “Calculation” submenus, among other options. The submenu “File” gives the user all needed functionalities for the creation, storage, and recovery of specific experiments. The functionality contained in the menu “File” is precisely what allows MieLab appearing to the users as a true virtual laboratory, where the experiments can be initiated, developed, stored, recovered, modified, restarted, and stopped. The menu “File” holds the following options: “New”, “Open”, “Save”, “Save As …”, “Open Experiment …”, “Export Theory”, and “Export Fit”. Through the “Calculation” submenu, we can run the theoretical calculations and fittings, as explained in Sections5.1 and 5.2, respectively.The toolbar component groups several toolbar buttons together and offers the user direct access to the same functionality offered by the main menu. The drawing canvas is the area in which the charts are created. Two types of chart are visualized in this area, the result of theoretical calculation and the optical responses for fitting.The “Parameters” panel component allows the user to set the parameter values characterizing the current problem, such as mean and standard deviation of the size distribution, concentration, and refractive index of the matrix. The provided parameter values can be entered manually, setting all those values in the “Parameters” panel (when option “New” in “File” submenu is selected) or via the reading of a given file (when option “Open” in “File” submenu is selected).


### 3.7. MieLab Workflow Behavior

Figure4 shows the flow chart of the major activities to be executed through MieLab while modeling the scattering of electromagnetic radiation. As depicted in this figure, the values of physical and geometrical parameters (e.g., composition, size, shape, refractive index, orientation, etc.) are set through the “Parameters” panel (1). The application runs the theoretical calculation of the optical properties of the particle or ensemble of particles, using the provided parameter values (2). The corresponding chart for the executed theoretical calculation is visualized on the drawing canvas (3). Now the fitting with an experimental file (data set) can be executed. The process implies the search of the physical parameters which give the best fit of the theoretical optical response to an experimental result. It comprises both the application of the algorithm for calculating the objective function and the optimization algorithm (4) and (5). The fitting results are improved (6). The fitted optical response is visualized on the drawing canvas (7).Figure 4
Schematic representation of MieLab’s workflow. The activity diagram describes the behavior of MieLab during a run.


## 3.1. Scattering Calculation

In this section, we will briefly discuss the most important aspects of the scattering algorithm implemented in MieLab. The interested reader can consult two additional papers: Yang [12] where the algorithm was originally presented and the paper containing the full description of our implementation [14]. In short, the analytical solution for the light scattering by a multilayered sphere can be found by expressing the EM field inside each layer of the sphere as a linear combination of the inward- and outward-travelling waves. Each layer is characterized by a size parameter xl=2πnmrl/λ=krl and a relative refractive index ml=nl/nm, l=1,2,…,L, where λ is the wavelength of the incident wave in vacuum, rl is the outer radius of the 1th layer, k is the propagation constant, and nm and nl are the refractive index of the medium outside the particle and its 1th component, respectively. For an L-layered sphere, there are 2L unknown field coefficients. Then, the scattering coefficients can be calculated by matching the tangential components of the electromagnetic fields at each interface and solving the obtained linear system of 2L independent equations. The final coefficients in the series are the scattering coefficients [12]
(1)an=anL+1=[Hna(mLxL)/mL+n/xL]ψn(xL)-ψn-1(xL)[Hna(mLxL)/mL+n/xL]ζn(xL)-ζn-1(xL),bn=bnL+1=[mLHnb(mLxL)+n/xL]ψn(xL)-ψn-1(xL)[mLHnb(mLxL)+n/xL]ζn(xL)-ζn-1(xL),
where ψn and ζn are the Riccati-Bessel functions, as defined by Kai and Massoli [27] and the determinants Hna and Hnb can be calculated by means of the following recursive expressions [12]
(2)Hna(m1x1)=Dn(1)(m1x1),Hna(mlxl)=G2Dn(1)(mlxl)-Qn(l)G1Dn(3)(mlxl)G2-Qn(l)G1,l=2,…,L,Hnb(m1x1)=Dn(1)(m1x1),Hnb(mlxl)=G̃2Dn(1)(mlxl)-Qn(l)G̃1Dn(3)(mlxl)G̃2-Qn(l)G̃1,l=2,…,L,Dn(1)(z)=ψn′(z)ψn(z),Dn(3)(z)=ζn′(z)ζn(z),Qn(l)=ψn(mlxl-1)/ζn(mlxl-1)ψn(mlxl)/ζn(mlxl),G1=mlHna(ml-1xl-1)-ml-1Dn(1)(mlxl-1),G2=mlHna(ml-1xl-1)-ml-1Dn(3)(mlxl-1),G̃1=ml-1Hnb(ml-1xl-1)-mlDn(1)(mlxl-1),G̃2=ml-1Hnb(ml-1xl-1)-mlDn(3)(mlxl-1).Once the scattering coefficients are obtained, we can calculate the extinction, scattering, radiation pressure, and radar backscattering efficiency factors (Qext, Qsca, Qpr, and Qbk) from the following relations [1]:(3)Qext=2xL2∑n=1∞(2n+1)Re{an+bn},Qsca=2xL2∑n=1∞(2n+1)(|an|2+|bn|2),Qpr=Qext-gQsca=Qext-4xL2{∑n=1∞[n(n+2)n+1Re{anan+1*+bnbn+1*}+2n+1n(n+1)Re{anbn*}]∑n=1∞[n(n+2)n+1Re{anan+1*+bnbn+1*}},Qbk=1xL2|∑n=1∞(2n+1)(-1)n(an-bn)|2.
Finally, from the previously calculated quantities, it is also possible to derive the absorption efficiency factor (Qabs), single scattering albedo (Λ) and asymmetry parameter (g)
(4)Qabs=Qext-Qsca,Λ=QscaQext,g=(Qext-Qpr)Qsca.


## 3.2. Optimization

The bounded, limited memory, Broyden-Fletcher-Goldfarb-Shanno (L-BFGS-B) algorithm [28, 29] was used in MieLab for the optimization process, together with the chi-square distribution, as the objective function
(5)χ2(r1,σ1,t2,σ2,…,tL,σL,nm,c)=1N∑i=1N(〈OD〉i(λi,r1,σ1,t2,σ2,…,tL,σL,nm,c)-ODexp⁡i)2ODexp⁡i,
where N is the number of points used for the fit, nm is the refractive index of the matrix, ODexpi are the values of the experimental optical density, while r1, ti, and σj are the radius of the first layer, the thickness of the ith layer (i=2,…,L), and the standard deviation of one of the two previous numbers for the jth layer (j=1,…,L), respectively. Then, starting with the given initial values of the variables to optimize (Xk=r1,σ1,t2,σ2,…,tL,σL,nm) and the experimental optical density, several optimization steps are performed iteratively. In each of the steps, the value of the theoretical optical density [30],  X2 and its gradient G are calculated, and subsequently, the L-BFGS-B algorithm is applied to obtain a new “optimized” value for all of the variables Xk. Lower and/or upper limits can be set for all the variables, in order to ensure that their optimized values have physical meaning (e.g., positive radii, etc.). The optimization process is finished when one of the following criteria is met:
(6)‖Gk‖≤εG,|(χ2)(i+1)-(χ2)(i)|≤εF⋅max⁡{(χ2)(i),(χ2)(i+1),1},|Xk(i+1)-Xk(i)|≤εX,
where εG, εF, and εX are positive numbers which define a precision of search, ∥·∥ means Euclidian norm, Gk, gradient projection onto a variable Xk, and the superindex represents the iteration number. Additionally, the optimization can be finished after a certain number of iterations. All these parameters controlling the fit can be tuned by the user, until (s)he finds the combination that fits best to her/his needs.


## 3.3. Optical Properties of Materials

At this point, MieLab allows two types of optical properties (more precisely, refractive index): a constant value and a table of values (“Element” in MieLab’s terminology). When a constant refractive index is selected, then the values ofn and k are entered directly into the text boxes corresponding to the layer being edited. On the other hand, if we choose to use the optical properties of an “Element”, then we have to select one of the available materials in the corresponding combo box. MieLab obtains the optical properties of the materials by reading, during the initialization, all the files with the extension “txt” found in the folder “RefIndex” (which in turn is located inside the folder where the program was installed). These are simply text files and the user can easily add a new or modify existing ones. The first lines of the file Ag.txt, which corresponds to “Silver”, are shown below, in order to exemplify the structure with which the data is stored:SilverPhys. Rev. B, 6 (1972) 437018.99 1.72331E-02 1.39e6 5.85e280.64 0.24 14.080.77 0.15 11.85.As can be seen, the first two lines of the file are the name of the “Element” (the same which is shown in the combo box) and a comment (in this case the reference to the source of the data), respectively. The third line can only take the values 0 or 1, where “1” means “apply size correction” [31] and “0” means the opposite. The next line depends on whether or not the size correction is to be applied; in the former case, it is expected to contain (in this order) the plasma frequency (ωp) in eV, the bulk damping constant (Γ) in eV, the Fermi velocity (vF) in m/s, and the atomic density in atoms/m3, while in the latter, only the last parameter is expected. The values of ωp, Γ and vF are used to apply the size correction [31], while the atomic density is necessary to calculate the optical density [30]. Finally, from the fifth line onwards, the actual values of the material refractive index are displayed. The first column presents the photon energy in eV, while the second and third ones present the real and imaginary parts of the complex refractive index, respectively.It must be said that although we believe that the options currently offered by MieLab for representing the optical properties of the materials are “good enough” for most applications, this is one of the areas we have selected to improve in future releases. In particular, we hope to add the possibility of using dispersion laws (Cauchy, Tauc-Lorentz, etc.) to represent the refractive index of materials.


## 3.4. Scripting Engine

RemObject Software’s Pascal Script (http://www.remobjects.com/ps.aspx), a free Pascal interpreter, which is compatible with both Delphi and Lazarus, has been used as the script engine in Mielab. Therefore, the scripts are written using a dialect of Object Pascal. When writing a script, in addition to the usual possibilities of Pascal, there are procedures to programmatically control almost all the MieLab’s functions. In addition, the object ActiveParticle (of type TLayeredSphere) is available to describe the geometry to be used in the calculations. A more detailed analysis of this functionality is beyond the scope of this paper. For more details, the reader can consult the sample script described below and MieLab’s Script Reference Guide, which is available on the website (http://scattering.sourceforge.net/ml_guide.html).


## 3.5. Software Architectural Model

As can be seen from the Figure2, MieLab architectural model is defined in terms of the following three main components:(i)
a shared data structure, namedworkspace, through which the components communicate with each other and the problem solution emerges,(ii)
a set of autonomous components, namedalgorithm objects, which contains specific knowledge and skills about the problem domain. Algorithm objects encapsulate the algorithms and models required for the elaboration of the data on the workspace,(iii)
a set of interactive components, namedinterface objects, for allowing the interaction of human users.Figure 2
(Color online) general overview of MieLab’s architectural model.Theworkspace contains all the information about the particle (e.g., composition, size, shape, refractive index, etc.), necessary to model the scattering of electromagnetic radiation and is organized into one or more hierarchies (depending on the application). The information at every level of a hierarchy represents partial solutions to the problem (e.g., scattering values, objective function values, results of the optimization algorithm, etc.). At every level, there is a “specific ontology” that describes the information. The knowledge and skills necessary for the scattering problem is distributed among the algorithm objects, which transform data from a hierarchy level to the same level or to the other. The algorithm objects execute these transformations using models, algorithmic procedures, or heuristic rules. Algorithm for calculating scattering, algorithm for calculating objective function, and optimization algorithm are some examples.The opportunistic reasoning for characterizing this organization is expressed in the following terms: at each step of the solution, thealgorithm object to be applied is determined dynamically, which produces an incremental generation of information in the workspace. The selection of a specific algorithm object is based on the status of the solution on the level and on the existence of other algorithm objects capable of improving the current state of the solution.As can be seen from the Figure2, the infrastructure also provides the tools for allowing the interaction with the system of human experts (interface objects) so that they can dynamically add knowledge on the scattering problem (e.g., modifying the particle parameters) or feedback on the results provided by the algorithms. In this way, the users can control the entire flow which gives rise to the final solution of the problem.


## 3.6. MieLab’s Graphical User Interface

Now, we describe MieLab’s graphical user interface, which is structured by four major components: a main menu, a toolbar, a drawing canvas, and a parameter panel. Figure3 shows MieLab’s graphical interface, integrating the previous components. The role of each of them is explained below.Figure 3
(Color online) Screenshot of MieLab’s graphical user interface.The main menu contains the “File” and “Calculation” submenus, among other options. The submenu “File” gives the user all needed functionalities for the creation, storage, and recovery of specific experiments. The functionality contained in the menu “File” is precisely what allows MieLab appearing to the users as a true virtual laboratory, where the experiments can be initiated, developed, stored, recovered, modified, restarted, and stopped. The menu “File” holds the following options: “New”, “Open”, “Save”, “Save As …”, “Open Experiment …”, “Export Theory”, and “Export Fit”. Through the “Calculation” submenu, we can run the theoretical calculations and fittings, as explained in Sections5.1 and 5.2, respectively.The toolbar component groups several toolbar buttons together and offers the user direct access to the same functionality offered by the main menu. The drawing canvas is the area in which the charts are created. Two types of chart are visualized in this area, the result of theoretical calculation and the optical responses for fitting.The “Parameters” panel component allows the user to set the parameter values characterizing the current problem, such as mean and standard deviation of the size distribution, concentration, and refractive index of the matrix. The provided parameter values can be entered manually, setting all those values in the “Parameters” panel (when option “New” in “File” submenu is selected) or via the reading of a given file (when option “Open” in “File” submenu is selected).


## 3.7. MieLab Workflow Behavior

Figure4 shows the flow chart of the major activities to be executed through MieLab while modeling the scattering of electromagnetic radiation. As depicted in this figure, the values of physical and geometrical parameters (e.g., composition, size, shape, refractive index, orientation, etc.) are set through the “Parameters” panel (1). The application runs the theoretical calculation of the optical properties of the particle or ensemble of particles, using the provided parameter values (2). The corresponding chart for the executed theoretical calculation is visualized on the drawing canvas (3). Now the fitting with an experimental file (data set) can be executed. The process implies the search of the physical parameters which give the best fit of the theoretical optical response to an experimental result. It comprises both the application of the algorithm for calculating the objective function and the optimization algorithm (4) and (5). The fitting results are improved (6). The fitted optical response is visualized on the drawing canvas (7).Figure 4
Schematic representation of MieLab’s workflow. The activity diagram describes the behavior of MieLab during a run.


## 4. MieLab as an Interactive Virtual Laboratory

The main features of the architectural model and functional behavior of MieLab are summarized below.(1) User-Friendly Interaction
MieLab provides a vast number of tools for allowing the interaction with the users so that they can dynamically modify the values of physical and geometrical parameters or give feedback on the results provided by the different algorithms. In this way, the users can control the entire flow which gives rise to the final solution of the scattering problem. Also, its modular construction enables to develop new parameters, tools or descriptors of the system under study.(2) Wide Range of Graphics for Interpreting Results
MieLab provides the user a wide variety of charts and tables to interpret the partial and final results obtained during execution. Absorption, scattering and extinction efficiencies, optical density, and so forth exist among the current charts.(3) Data Persistence and Recovery Experiments
The functionality and performance of MieLab make it appearing like a virtual laboratory to the users, where experiments can be started, developed, interrupted, stored, recovered, continued, modified, and restarted from any point of the run. Also, a more complex performance is intended in the future, where internal routines automatically optimize parameters or give feedback to the fitting function or optimization variables.(4) Extensibility and Reusability
The development of MieLab not only ensures the integrated use of the implemented models and algorithms for spherical particles using exact techniques but also it is flexible enough to incorporate new models and algorithms for computing the scattering properties of other type of particles using exact and approximated techniques. In this sense, MieLab architecture has been designed and implemented thinking about its extensibility and reusability.(5) Portability
MieLab has been implemented in Pascal (Delphi), so it is initially designed for MS Windows, but it is planned to port it to the Lazarus/Free Pascal development suite in order to make it available in a wider range of platforms, limited only for the processor speed and the available memory of the computer. Moreover, although this feature can hardly be called “portability”, MieLab can run now in other architectures such as Mac OSX or Linux using the WINE emulator (http://www.winehq.org/).(6) Task Automation through the Script Execution
MieLab includes a powerful scripting engine which exposes all the functionality of the internal algorithms in an easy to use programming language (a dialect of Object Pascal). This allows the end user to customize the application and to automate virtually all kind of complex/repetitive calculations without the need of dealing with the complication of the involved algorithms or recompiling the program.


## 5. Using MieLab

MieLab offers tools to perform several tasks related with the calculation of the optical properties of nanoparticles and particle ensembles, like the following:(i)
theoretical calculation of the optical properties of a nanoparticle or an ensemble of independent (noninteracting) nanoparticles,(ii)
finding of the physical parameters which give the best fit of the theoretical optical response to an experimental result,(iii)
automation of the previous tasks (or a combination of both of them) to perform unattended calculations.In this section, we try to provide examples of some possible (nontrivial) usages of MieLab. For each example, a detailed description of how the software was used is provided, so the reader could reproduce the examples seamlessly. We have restricted ourselves to relatively simple cases in order to illustrate the use of the program instead of devoting time to the examples themselves. However, it should be noted that we have used and tested extensively MieLab in our own research work. For instance, the reader can consult this work [30], where we have used an early incarnation of MieLab to calculate the size distributions of ensembles of metallic nanoparticles. As can be seen there, the values obtained from “real-world” spectra are equivalent to the size distributions calculated using other experimental techniques.


### 5.1. Theoretical Calculation

The execution of a theoretical calculation is relatively straightforward; it is, however, exemplified in order to use the result as an “experimental” file for a later test of the fitting capabilities. The calculations (Figure5, black line) were performed for the following parameters: a Gaussian distribution of radii with mean of 7.3 nm and standard deviation of 2.6 nm, a concentration of 4.8 atoms/cm2, and a refractive index of the matrix of 1.46. After setting all those values in the “Parameters” panel, at the right of the main window, the calculation can be performed (Calculation→Theory); the obtained optical response can then be saved using File→Export fit.Figure 5
(Color online) theoretical calculation and fittings performed in the described examples. The theoretical calculation was performed for a Gaussian distribution of radii with mean of 7.3 nm and standard deviation of 2.6 nm, a concentration of 4.8 atoms/cm2 and a refractive index of the matrix of 1.46. The fittings were done using as “experimental” file the previously generated optical response and the initial parameters shown in Table 1.


### 5.2. Fitting to an Experimental File

In order to illustrate the use of the fitting capability and its reliability, four fittings were performed using the previously generated optical response as “experimental” file. The initial parameters used for them are shown in Table1 and the corresponding optical response in Figure 5 (dotted lines). Once any of those initial configurations are ready, the parameters that will be fitted must be selected (mean radius, standard deviation, and concentration) using the checkbox at the left of the parameter value in the “Parameter” panel, and the experimental file must be loaded using File→Open experiment and selecting the file saved in the previous example.Table 1
Initial parameters used for the fitting.
Fit # 1Fit # 2Fit # 3Fit # 4Distribution of radiiNone (single value)GaussianNone (single value)GaussianMean radius3.0 nm3.0 nm20.0 nm20.0 nmStandard deviation—1.0 nm—10.0 nmConcentration2.0 atm/cm22.0 atm/cm220.0 atm/cm220.0 atm/cm2Refractive index of matrix1.461.461.461.46Once the above steps have been completed, we are ready to test the fitting, which is done usingCalculations→Fit; the parameters obtained from the fit are shown in Table 2 and the corresponding optical spectra in Figure 5 (continuous lines). The fitting is considerably faster if a single value of radius is used; in this case, however, the obtained mean radius is higher than the real one and, within the analyzed range, independent of the initial value. This happens because the intensity of the optical response of a particle is directly proportional to its volume, and therefore, the bigger ones contribute more to the optical spectrum. On the other hand, if Gaussian distributions are used, the obtained results are exact and equal to the real ones up to three significant digits, and independent of the initial values, which shows the reliability of the method, but the drawback in this case is that the calculations are considerably slower.Table 2
Final parameters, obtained from the fitting.
Fit # 1Fit # 2Fit # 3Fit # 4Distribution of radiiNone (single value)GaussianNone (single value)GaussianMean radius9.77 nm7.30 nm9.77 nm7.30 nmStandard deviation—2.60 nm—2.60 nmConcentration4.78 atm/cm24.80 atm/cm24.78 atm/cm24.80 atm/cm2Refractive index of matrix1.461.461.461.46


### 5.3. Scripting Capabilities

The final example is related to the scripting capabilities of MieLab. In order to access them, it is necessary to go to the “Scripts” window (Calculations→Scripts); once there, we can write and run (Run→Run) some code to automate some tasks. This capability of MieLab can be tested by running the following script, which will perform all the operations described in the previous examples in an unattended way. All the commands are describes in the code itself. For more details on the available commands the interested reader is encouraged to check MieLab’s Script Reference Guide, located at http://scattering.sourceforge.net/ml_guide.html.program Examples;varPath: String;begin//Assign the path  (the directory must exist)Path  := ‘C:∖Results∖';//Delete all layers in the active particleActiveParticle.DeleteAllLayers;//Add first layerActiveParticle.AddLayer;//Setup particle configurationActiveParticle.ThicknessType[0]:= ttDistribution;ActiveParticle.ThicknessValue[0]:= 7.3;ActiveParticle.ThicknessStdDev[0]:= 2.6;ActiveParticle.RefIndexType[0]:= ritElement;ActiveParticle.ElementName[0]:= 'Silver';ActiveParticle.RefIndexMatrix:= 1.46;ActiveParticle.Concentration:= 4.8;//Calculate ‘experimental' file//Update the changesNanoToText;//Do the theoretical calculationDoTheory;//Save the resultsSaveFit(Path + ‘Experimental.csv');//Save the projectSaveProject(Path + ‘Experimental.mlp');//Make first fitOpenExperimentFile(Path + ‘Experimental.csv');ActiveParticle.ThicknessType[0]:= ttSingleValue;ActiveParticle.ThicknessValue[0]:= 3.0;ActiveParticle.FitThicknessValue[0]:= True;ActiveParticle.Concentration:= 2.0;ActiveParticle.FitConcentration:= True;//Update the changesNanoToText;//Save the projectSaveProject(Path + ‘Fit1Initial.mlp');//Do the fitDoFit;//Save the resultsSaveFit(Path + ‘Fit1.csv');//Save the projectSaveProject(Path + ‘Fit1Final.mlp');//Make second fitOpenExperimentFile(Path + ‘Experimental.csv');ActiveParticle.ThicknessType[0]:= ttDistribution;ActiveParticle.ThicknessValue[0]:= 3.0;ActiveParticle.FitThicknessValue[0]:= True;ActiveParticle.ThicknessStdDev[0]:= 1.0;ActiveParticle.FitThicknessStdDev[0]:= True;ActiveParticle.Concentration:= 2.0;ActiveParticle.FitConcentration:= True;//Update the changesNanoToText;//Save the projectSaveProject(Path + ‘Fit2Initial.mlp');//Do the fitDoFit;//Save the resultsSaveFit(Path + ‘Fit2.csv');//Save the projectSaveProject(Path + ‘Fit2Final.mlp');//Make third fitOpenExperimentFile(Path + ‘Experimental.csv');ActiveParticle.ThicknessType[0]:= ttSingleValue;ActiveParticle.ThicknessValue[0]:= 20.0;ActiveParticle.FitThicknessValue[0]:= True;ActiveParticle.Concentration  := 20.0;ActiveParticle.FitConcentration  := True;//Update the changesNanoToText;//Save the projectSaveProject(Path + ‘Fit3Initial.mlp');//Do the fitDoFit;//Save the resultsSaveFit(Path + ‘Fit3.csv');//Save the projectSaveProject(Path + ‘Fit3Final.mlp');//Make fourth fitOpenExperimentFile(Path + ‘Experimental.csv');ActiveParticle.ThicknessType[0]:= ttDistribution;ActiveParticle.ThicknessValue[0]:= 20.0;ActiveParticle.FitThicknessValue[0]:= True;ActiveParticle.ThicknessStdDev[0]:= 10.0;ActiveParticle.FitThicknessStdDev[0]:= True;ActiveParticle.Concentration  := 20.0;ActiveParticle.FitConcentration  := True;//Update the changesNanoToText;//Save the projectSaveProject(Path + ‘Fit4Initial.mlp');//Do the fitDoFit;//Save the resultsSaveFit(Path + ‘Fit4.csv');//Save the projectSaveProject(Path + ‘Fit4Final.mlp');Writeln(‘Done!!!');end.


## 5.1. Theoretical Calculation

The execution of a theoretical calculation is relatively straightforward; it is, however, exemplified in order to use the result as an “experimental” file for a later test of the fitting capabilities. The calculations (Figure5, black line) were performed for the following parameters: a Gaussian distribution of radii with mean of 7.3 nm and standard deviation of 2.6 nm, a concentration of 4.8 atoms/cm2, and a refractive index of the matrix of 1.46. After setting all those values in the “Parameters” panel, at the right of the main window, the calculation can be performed (Calculation→Theory); the obtained optical response can then be saved using File→Export fit.Figure 5
(Color online) theoretical calculation and fittings performed in the described examples. The theoretical calculation was performed for a Gaussian distribution of radii with mean of 7.3 nm and standard deviation of 2.6 nm, a concentration of 4.8 atoms/cm2 and a refractive index of the matrix of 1.46. The fittings were done using as “experimental” file the previously generated optical response and the initial parameters shown in Table 1.


## 5.2. Fitting to an Experimental File

In order to illustrate the use of the fitting capability and its reliability, four fittings were performed using the previously generated optical response as “experimental” file. The initial parameters used for them are shown in Table1 and the corresponding optical response in Figure 5 (dotted lines). Once any of those initial configurations are ready, the parameters that will be fitted must be selected (mean radius, standard deviation, and concentration) using the checkbox at the left of the parameter value in the “Parameter” panel, and the experimental file must be loaded using File→Open experiment and selecting the file saved in the previous example.Table 1
Initial parameters used for the fitting.
Fit # 1Fit # 2Fit # 3Fit # 4Distribution of radiiNone (single value)GaussianNone (single value)GaussianMean radius3.0 nm3.0 nm20.0 nm20.0 nmStandard deviation—1.0 nm—10.0 nmConcentration2.0 atm/cm22.0 atm/cm220.0 atm/cm220.0 atm/cm2Refractive index of matrix1.461.461.461.46Once the above steps have been completed, we are ready to test the fitting, which is done usingCalculations→Fit; the parameters obtained from the fit are shown in Table 2 and the corresponding optical spectra in Figure 5 (continuous lines). The fitting is considerably faster if a single value of radius is used; in this case, however, the obtained mean radius is higher than the real one and, within the analyzed range, independent of the initial value. This happens because the intensity of the optical response of a particle is directly proportional to its volume, and therefore, the bigger ones contribute more to the optical spectrum. On the other hand, if Gaussian distributions are used, the obtained results are exact and equal to the real ones up to three significant digits, and independent of the initial values, which shows the reliability of the method, but the drawback in this case is that the calculations are considerably slower.Table 2
Final parameters, obtained from the fitting.
Fit # 1Fit # 2Fit # 3Fit # 4Distribution of radiiNone (single value)GaussianNone (single value)GaussianMean radius9.77 nm7.30 nm9.77 nm7.30 nmStandard deviation—2.60 nm—2.60 nmConcentration4.78 atm/cm24.80 atm/cm24.78 atm/cm24.80 atm/cm2Refractive index of matrix1.461.461.461.46


## 5.3. Scripting Capabilities

The final example is related to the scripting capabilities of MieLab. In order to access them, it is necessary to go to the “Scripts” window (Calculations→Scripts); once there, we can write and run (Run→Run) some code to automate some tasks. This capability of MieLab can be tested by running the following script, which will perform all the operations described in the previous examples in an unattended way. All the commands are describes in the code itself. For more details on the available commands the interested reader is encouraged to check MieLab’s Script Reference Guide, located at http://scattering.sourceforge.net/ml_guide.html.program Examples;varPath: String;begin//Assign the path  (the directory must exist)Path  := ‘C:∖Results∖';//Delete all layers in the active particleActiveParticle.DeleteAllLayers;//Add first layerActiveParticle.AddLayer;//Setup particle configurationActiveParticle.ThicknessType[0]:= ttDistribution;ActiveParticle.ThicknessValue[0]:= 7.3;ActiveParticle.ThicknessStdDev[0]:= 2.6;ActiveParticle.RefIndexType[0]:= ritElement;ActiveParticle.ElementName[0]:= 'Silver';ActiveParticle.RefIndexMatrix:= 1.46;ActiveParticle.Concentration:= 4.8;//Calculate ‘experimental' file//Update the changesNanoToText;//Do the theoretical calculationDoTheory;//Save the resultsSaveFit(Path + ‘Experimental.csv');//Save the projectSaveProject(Path + ‘Experimental.mlp');//Make first fitOpenExperimentFile(Path + ‘Experimental.csv');ActiveParticle.ThicknessType[0]:= ttSingleValue;ActiveParticle.ThicknessValue[0]:= 3.0;ActiveParticle.FitThicknessValue[0]:= True;ActiveParticle.Concentration:= 2.0;ActiveParticle.FitConcentration:= True;//Update the changesNanoToText;//Save the projectSaveProject(Path + ‘Fit1Initial.mlp');//Do the fitDoFit;//Save the resultsSaveFit(Path + ‘Fit1.csv');//Save the projectSaveProject(Path + ‘Fit1Final.mlp');//Make second fitOpenExperimentFile(Path + ‘Experimental.csv');ActiveParticle.ThicknessType[0]:= ttDistribution;ActiveParticle.ThicknessValue[0]:= 3.0;ActiveParticle.FitThicknessValue[0]:= True;ActiveParticle.ThicknessStdDev[0]:= 1.0;ActiveParticle.FitThicknessStdDev[0]:= True;ActiveParticle.Concentration:= 2.0;ActiveParticle.FitConcentration:= True;//Update the changesNanoToText;//Save the projectSaveProject(Path + ‘Fit2Initial.mlp');//Do the fitDoFit;//Save the resultsSaveFit(Path + ‘Fit2.csv');//Save the projectSaveProject(Path + ‘Fit2Final.mlp');//Make third fitOpenExperimentFile(Path + ‘Experimental.csv');ActiveParticle.ThicknessType[0]:= ttSingleValue;ActiveParticle.ThicknessValue[0]:= 20.0;ActiveParticle.FitThicknessValue[0]:= True;ActiveParticle.Concentration  := 20.0;ActiveParticle.FitConcentration  := True;//Update the changesNanoToText;//Save the projectSaveProject(Path + ‘Fit3Initial.mlp');//Do the fitDoFit;//Save the resultsSaveFit(Path + ‘Fit3.csv');//Save the projectSaveProject(Path + ‘Fit3Final.mlp');//Make fourth fitOpenExperimentFile(Path + ‘Experimental.csv');ActiveParticle.ThicknessType[0]:= ttDistribution;ActiveParticle.ThicknessValue[0]:= 20.0;ActiveParticle.FitThicknessValue[0]:= True;ActiveParticle.ThicknessStdDev[0]:= 10.0;ActiveParticle.FitThicknessStdDev[0]:= True;ActiveParticle.Concentration  := 20.0;ActiveParticle.FitConcentration  := True;//Update the changesNanoToText;//Save the projectSaveProject(Path + ‘Fit4Initial.mlp');//Do the fitDoFit;//Save the resultsSaveFit(Path + ‘Fit4.csv');//Save the projectSaveProject(Path + ‘Fit4Final.mlp');Writeln(‘Done!!!');end.


## 6. Conclusions

In this work, we have presented MieLab, a flexible virtual laboratory created by combining several algorithms, models, and techniques into a computational platform. Some of the most important capabilities of the application, which is freely available at the web addresshttp://scattering.sourceforge.net, have been described. This piece of software is expected to be useful for the students trying to learn the basics of light scattering, as well as the experts performing “serious” calculations. During the design and development of MieLab, a great deal of attention was paid to “hide” the inherent complexity of the rigorous and accurate algorithms used internally by means of a user-friendly graphical interface. The main features of MieLab are illustrated through several examples.However, the features described in the present paper are only a few of the possibilities that a tool like MieLab can offer, and the software in its present form is not a final one. In future, the software is planned to be improved in three important aspects. The first aspect which will be improved soon is the number of supported scattering algorithms. Currently, MieLab only allows the study of the scattering by a (multilayered) sphere. It is planned to add the T-matrix algorithm [2] in it to perform calculations on particles of arbitrary shapes having rotational symmetry. Possibly, the Discrete Dipole Approximation [32, 33] would be incorporated as well to improve the reliability of the obtained results on the particles of arbitrary shape. Inclusion of theoretical laws (Cauchy, Tauc-Lorentz, critical points, etc.) to describe the dielectric function of the materials is also under consideration for the future versions. Finally, we also want to increase the options of available optimization algorithms. Currently, MieLab uses the bounded limited memory Broyden-Fletcher-Goldfarb-Shanno (L-BFGS-B) [28, 29] multivariate optimization algorithm but a genetic algorithm is being developed to perform the optimization and hopefully will be included in MieLab soon.



---

*Source: 583743-2011-03-06.xml*