Resonant transmission line method for unconventional fibers

We provide a very general review of the resonant transmission line method for optical fiber problems. The method has been found to work seamlessly for a variety of difficult problems including elliptical and eccentric core fibers as well as “holey”, photonic crystal fibers. A latest version has been shown to offer great versatility with respect to cases of unconventional, inhomogeneous index profiles.


Introduction
The method of resonant transmission lines (RTL) has evolved from some original observations around the generic theme of transverse resonance in direct association with the standard theory of the solutions of the general telegrapher's equation in transmission lines.The particular advantage of the method is that it offers the possibility of reducing even full three-dimensional problems into one its dimensional equivalent, at least wherever a full variable separation is possible due to the particular symmetries of a given coordinate system.A full historical review of the use of transmission lines for modeling electromagnetic problems would require too much space but we can provide some important milestones in chronological order in brief.
The transmission line theme can be traced as far back as in the old treatise by Schelkunoff [1].Further physical justification is offered by a seminal paper by Marcuvitz and Schwinger [2] where we find a first version of the transverse resonance condition as a means to satisfy certain boundary conditions via appropriate impedance matching.Extensive use of the telegrapher's equation is also found in Ramachandran [3], Gallawa [4], Clarricoats and Oliner [5,6], Yoneyama [7], Borneman and Arndt [8], Tao [9], Shigesawa and Tsuji [10] as well as Dahl et al. [11].More recently, Moshonas et al. [12] also used transmission line methods for the study of photovoltaic cells.An important comparison between transverse stationary modes and their similarity with quantum wells is mentioned by Bialynicki-Birula [13] and has led to the application of the RTL method by some of the authors in Schrödinger type and general Sturm-Liouville problems [14][15][16].An early application of similar ideas in optical fibers can be traced in a paper by Yeh and Lindgren [17].Tamir also discussed guided waves based on impedance matching in Reference [18], while later Carlin and Zmuda [19] applied transfer matrix methods for inhomogeneous fibers.More recently, Mencarelli and Rossi [20] used transverse resonance for studying multilayered photonic crystals.
The modern RTL method was originally introduced by Papageorgiou and Boucouvalas [21][22][23][24][25][26] in simple cases of standard step index and radial index profiles.Lately, it was also applied successfully in the case of superlattices with periodic potentials [27] which are similar to a certain class of photonic crystal fibers (PCFs) or Bragg fibers in general.Several other cases of unusual refractive index profiles have been successfully treated with the same method [28][29][30][31].In the present work, we expose more recent developments in an effort to expand the validity of the RTL method in all unconventional fiber models.
We shall heretofore refer to all the cylindrical optical fibers as the class of conventional optical fibers (COFs) for as long as they can be separated in a set of n very thin successive cylindrical layers of average radius r with uniform refractive index values given as n i .Indices of successive layers will be allowed in general to be different for each step varying as n i = n(r i ).The resulting, total index profile n(r) variation from the center of the fiber up to the limit outer air medium, completely defines the propagation properties inside the fiber.
We shall also separate another class of unconventional optical fibers (UOFs) as the ones in which at least in some of its successive thin cylindrical layers present an additional variation of the total index profile along its radial coordinate ϕ in the form η(r,ϕ).Such cases include elliptic core, non-symmetric or eccentric core fibers, and in general, all cases of not strictly circular core fibers.In these cases, any discretization into thin cylindrical layers that "cuts" through the core and cladding, results in a variation of the refractive index along ϕ.In Figure 1, we show such a case of a thin, circular cylindrical layer cutting an elliptic core fiber.As is evident in the schematic, the variation is analogous to the arc length of the circular sector cut by the ellipse for each discretization step.
Appl.Sci.2019, 9, x 2 of 24 transfer matrix methods for inhomogeneous fibers.More recently, Mencarelli and Rossi [20] used transverse resonance for studying multilayered photonic crystals.The modern RTL method was originally introduced by Papageorgiou and Boucouvalas [21][22][23][24][25][26] in simple cases of standard step index and radial index profiles.Lately, it was also applied successfully in the case of superlattices with periodic potentials [27] which are similar to a certain class of photonic crystal fibers (PCFs) or Bragg fibers in general.Several other cases of unusual refractive index profiles have been successfully treated with the same method [28][29][30][31].In the present work, we expose more recent developments in an effort to expand the validity of the RTL method in all unconventional fiber models.
We shall heretofore refer to all the cylindrical optical fibers as the class of conventional optical fibers (COFs) for as long as they can be separated in a set of n very thin successive cylindrical layers of average radius r with uniform refractive index values given as i n .Indices of successive layers will be allowed in general to be different for each step varying as . The resulting, total index profile n(r) variation from the center of the fiber up to the limit outer air medium, completely defines the propagation properties inside the fiber.
We shall also separate another class of unconventional optical fibers (UOFs) as the ones in which at least in some of its successive thin cylindrical layers present an additional variation of the total index profile along its radial coordinate φ in the form η(r,φ).Such cases include elliptic core, nonsymmetric or eccentric core fibers, and in general, all cases of not strictly circular core fibers.In these cases, any discretization into thin cylindrical layers that "cuts" through the core and cladding, results in a variation of the refractive index along φ.In Figure 1, we show such a case of a thin, circular cylindrical layer cutting an elliptic core fiber.As is evident in the schematic, the variation is analogous to the arc length of the circular sector cut by the ellipse for each discretization step.Additionally, photonic fibers made of silica with a set of small air holes around their centers will have many cylindrical layers of varying refractive indexes along θ.In Figure 2, we show an example of a PCF with a hexagonal lattice of 6 rows of air holes.
Appl.Sci.2019, 9, x 3 of 24 Additionally, photonic fibers made of silica with a set of small air holes around their centers will have many cylindrical layers of varying refractive indexes along θ.In Figure 2, we show an example of a PCF with a hexagonal lattice of 6 rows of air holes.Again, radial discretization with thin, circular cylindrical layers results in alternating cuts through air holes and silica of a PCF as depicted in the schematic.The main aim in the present work is to develop the original RTL method so as to be able to include all such UOF cases via an appropriate transformation to mathematically equivalent COF cases.
In the following sections, we first introduce the exact mathematical formulation of the discretization process for Maxwell's equations in the context of the RTL method in the following section.We also show how to perform an algebraic decoupling of the resulting transmission line equations and the transfer function method for computing the actual fields from the eigenvalues obtained.In Section 3, the same tools are applied in four more difficult cases of UOF with either, asymmetry or eccentricity of the central core.We also perform the same analysis for elliptic and rectangular cores as well as for the PCF case.

Mathematical Equivalence of Homogeneous Circular Cylindrical Layers to Electric Transmission Lines
The basis for the application of the previously introduced RTL method is the radial discretization of all cylindrical fibres via a separation into a succession of thin cylindrical layers, each one with its own constant refractive index n.These layers can be made to extend outside of the cladding in order to take into consideration the effect of the surrounding air (n = 1).Each thin cylindrical layer could have thickness δr proportional to each average radius r which means that given discrete steps as  =  2 −  1 with  =  1 +  2 2 one has Again, radial discretization with thin, circular cylindrical layers results in alternating cuts through air holes and silica of a PCF as depicted in the schematic.The main aim in the present work is to develop the original RTL method so as to be able to include all such UOF cases via an appropriate transformation to mathematically equivalent COF cases.
In the following sections, we first introduce the exact mathematical formulation of the discretization process for Maxwell's equations in the context of the RTL method in the following section.We also show how to perform an algebraic decoupling of the resulting transmission line equations and the transfer function method for computing the actual fields from the eigenvalues obtained.In Section 3, the same tools are applied in four more difficult cases of UOF with either, asymmetry or eccentricity of the central core.We also perform the same analysis for elliptic and rectangular cores as well as for the PCF case.

Mathematical Equivalence of Homogeneous Circular Cylindrical Layers to Electric Transmission Lines
The basis for the application of the previously introduced RTL method is the radial discretization of all cylindrical fibres via a separation into a succession of thin cylindrical layers, each one with its own constant refractive index n.These layers can be made to extend outside of the cladding in order to take into consideration the effect of the surrounding air (n = 1).Each thin cylindrical layer could have thickness δr proportional to each average radius r which means that given discrete steps as δr = r 2 − r 1 with r = r 1 +r 2 2 one has Appl.Sci.2019, 9, 270 4 of 21 For any such circular cylindrical layer Maxwell equations (for a constant wavelength, i.e., constant frequency "ω") can be written in their standard form as Taking into consideration the relations ωµ 0 = k 0 z 0 and ωε 0 = k 0 z 0 where k 0 = ω c , z 0 = 120π and replacing z 0 In circular cylindrical geometry of coordinates (r, ϕ, z) the following set of three partial differential equations can be derived by the first vector Maxwell equation as Applying a Fourier Transform along "z" and "ϕ" with wave numbers "β" and "l", where l is integer (because along "ϕ" we have Fourier series of period 2π), the set (4) becomes: In (5) we use new variables E r , E ϕ , E z , H r , H ϕ , H z to denote the Fourier Transforms of the respective electromagnetic field components.Furthermore, replacing β and r by their reduced variables according to the following relations: Following a similar approach, the second Maxwell vector Equation (3) can be written in the form Furthermore, following a cumbersome analysis, it is possible to prove that the system of Equations ( 6) and ( 7) can be transformed in a set of four differential Equation (8), relating the equivalent "voltage' and "current" functions V M , I M , V E , I E defined as follows: where we use the notation F = (βr) In (8) we introduced the total propagation factor γ 2 = l 2 r 2 + β 2 − n 2 and the auxiliary function At this point it is noticed that V M , I M , V E , I E are continuous functions at the boundaries because the tangential components of electric and magnetic fields E z on the cylindrical surface are continuous functions passing the boundaries' of the cylindrical layer.Using the previous relations, the Fourier Transforms of the electromagnetic field components along (r, l, β) can be expressed as functions of their equivalent "voltages" and "currents" functions with the auxiliary relations It becomes evident by inspection that the final Equation ( 8) represents two coupled electric transmission lines.

Decoupling the Transmission Line Equations
The prescribed set of Equation ( 8) constitutes a homogeneous set of ordinary differential equations of r and considering that the all vectors [V M , I M , V E , I E ] can be turned into exponential functions of r given by V M = V M e ξr , I M = I M e ξr , V E = V E e ξr , I E = I E e ξr , where V M , I M , V E , I E are constants, i.e., not functions of r.Thus, the system (8) can be transformed in an algebraic set of the following four equations This then leads to the eigenvalue equations From the standard form of the eigenvalue problem we obtain through the determinant differential equations as follows Hence the system has two eigenvalues and two mutually excluded or "normal" eigenvectors.The eigenvectors will be found by replacing ξ 2 by its value.Thus, for Since the auxiliary M function has the sign of l, the set (V s , I s ), for l = −l becomes equal to the set (V d , I d ).Thus, we can consider as a unique solution for the set (V s , I s ) and the integer 'l 'varies −∞ to + ∞, and of course Furthermore, V s , I s should be continuous functions at their boundaries although n(r) varies from layer to layer.This is achieved via the adjustment Another option for achieving continuity is to consider the functions V ss = V M n + V E and I ss = nI M + I E .In this case, V ss = 2V E and I ss = 2I E are also continuous.Thus Thus, the set of two coupled transmission lines ( 9) is equivalent to two independent transmission lines (12) and (13).
The two waves represented by the equations of transmission lines ( 12) and ( 13), are geometrically normal because the first is related to the magnetic field and the second to the electric field that are geometrically normal for transmitted EM waves.This property is an inherent property of EM modes in optical fibers related to birefringence phenomena.However, the β respective values, for any mode, are always found to be very close and can be considered as practically equal.

Equivalent Circuits for Cylindrical Layers, Boundary Conditions, and Birefringence
Taking into consideration the transmission line theory, it can be proved that each layer of infinitesimal thickness δr is equivalent to a T-circuit as the one shown in Figure 3    for any mode, are always found to be very close and can be considered as practically equal.

Equivalent Circuits for Cylindrical Layers, Boundary Conditions, and Birefringence
Taking into consideration the transmission line theory, it can be proved that each layer of infinitesimal thickness δr is equivalent to a T-circuit as the one shown in Figure For  ≪ 1 the impedances can be approximated by the equivalent relations If  2 > 0, both   ,   are "capacitive" reactances, for  2 < 0 however   becomes "inductive" reactance.For (  ,   ) the approximate respective impedances of the T-circuit are given as For ξδr 1 the impedances can be approximated by the equivalent relations If ξ 2 > 0, both Z B , Z p are "capacitive" reactances, for ξ 2 < 0 however Z B becomes "inductive" reactance.For (V E , I E ) the approximate respective impedances of the T-circuit are given as As previously stated, the functions (V M , I M ) of each layer are continuous at the cylindrical boundaries of the layer, thus if we divide the fiber (including a sufficient number of air layers) in successive thin layers and replace them by their equivalent T-circuits, an overall lossless transmission line is formed with only reactive elements.For given "l", the "β" values that lead to the resonance of the overall transmission line are the eigenvalues of the whole optical fiber.
When a transmission line is in resonance, at any arbitrary point r 0 of the line, the sum of reactive impedances arising from the successive T-circuits on the left and right sides of r 0 should be equal to zero, thus the equation giving the eigenvalues of the transmission line is the following: .
Z L.r 0 , .Z R.r 0 are the overall reactive impedances of successive T-circuits on the left and right of r 0 , using Equations ( 14) or (15).The value of r 0 is usually given by the core radius.For the same "l" the Equations ( 14) and (15) give usually slightly different values of 'β'.This phenomenon is called "Birefringence".For circular step index fibers, the birefringence is negligible; however, for elliptic fibers and fibers of any other non-circular cores, the birefringence phenomenon could be not negligible.
In order to calculate the overall reactive impedances on the left and right of r 0 we should find the impedances for r → 0 and for r → ∞ .As we proceed to 0 or to ∞, the remaining piece of transmission line becomes "homogeneous", i.e., its overall reactive impedance is equal to its characteristic impedance given by Z = ξ jF or ξ jn 2 F .Then we must have For l = 0 Z r→0 = ∞ (open circuit at the center of the equivalent transmission line) It is useful to notice that there is an equivalence between our formulation and the classic formulation modes of optical fibers.In particular, for l = 0, the modes (V M ,I M ) are the TM modes, while the modes (V E ,I E ) are the TE modes.For l > 0, the modes (V M ,I M ) are the HE modes, while the modes (V E ,I E ) are their HE birefringence modes.For l <0 the modes (V M ,I M ) are the EH modes, while the modes (V E ,I E ) are their EH birefringence modes.For any given l, using the resonance technique the β values of the two birefringence modes can be calculated.The Equation ( 16) is given as a MATLAB function in Appendix A. Let us consider for example a step-index fiber of n 1 = 1.54, n 2 = 1.47 the V M , V E , fundamental modes for V = 3.3, can be calculated and their β/k 0 values are respectively 1.518934962534846 and 1.518340184686295, hence their birefringence is equal to 0.0004947 or 0.0391%.The β/k 0 value for the equivalent mode V eq was also calculated and was equal to 1.518638548412019 (that is approximately equal to the mean value of the previous β/k 0 values), while the β/k 0 value calculated classically by Bessel functions is equal to 1.518642063686336.These β values are very close differing only 0.0002315%.
In the following Figure 4, the normalized birefringence of the step-index fibers for n 1 = 1.54, n 2 = 1.47, and of n 1 = 1.475, n 2 = 1.47 as functions of V are shown.We notice that for any V, the normalized birefringence is almost proportional to Δn 2 = (n1 − n2) 2 , thus the birefringence of step-index fibers of very small Δn is negligible.For instance, for a value of V = 2.4, and Δn = 1.54 − 1.47 = 0.07, the birefringence is found to be 0.168 × 0.0049 = 0.0008232 or ~0.055% on the average β, while for Δn = 1.475 − 1.47 = 0.005, the birefringence becomes 0.168 × 0.000025 = 0.000042 or ~0.0028% on the average β.What is remarkable is that our method is sensitive and calculates it.

Calculating "Voltages" VM, VE and "Currents" IM, IE and Resulting Fields
For any given , using the resonance technique the β values of the two birefringence modes can be calculated.These β values are practically the same, thus we can consider them as equal or we can consider as the proper value of β the mean value of the two modes.Taking VM = 1 at the center point of the fiber (r = 0), the respective value of IM at the same point can be calculated by the respective terminal impedance.Using the matrix relations between input-output for the equivalent successive T-circuits, the values of   and   at the rest thin cylindrical layers can be calculated.In fact, from the general theory of the telegrapher's equation we know that the inputs and outputs are associated via a transfer matrix as follows  We notice that for any V, the normalized birefringence is almost proportional to ∆n 2 = (n 1 − n 2 ) 2 , thus the birefringence of step-index fibers of very small ∆n is negligible.For instance, for a value of V = 2.4, and ∆n = 1.54 − 1.47 = 0.07, the birefringence is found to be 0.168 × 0.0049 = 0.0008232 or ~0.055% on the average β, while for ∆n = 1.475 − 1.47 = 0.005, the birefringence becomes 0.168 × 0.000025 = 0.000042 or ~0.0028% on the average β.What is remarkable is that our method is sensitive and calculates it.

Calculating "Voltages" V M , V E and "Currents" I M , I E and Resulting Fields
For any given l, using the resonance technique the β values of the two birefringence modes can be calculated.These β values are practically the same, thus we can consider them as equal or we can consider as the proper value of β the mean value of the two modes.Taking V M = 1 at the center point of the fiber (r = 0), the respective value of I M at the same point can be calculated by the respective terminal impedance.Using the matrix relations between input-output for the equivalent successive T-circuits, the values of V M and I M at the rest thin cylindrical layers can be calculated.In fact, from the general theory of the telegrapher's equation we know that the inputs and outputs are associated via a transfer matrix as follows In Equation ( 17), the characteristic impedance should be taken as Z(r) = ξ(r)/jF(r) to fit with the previous analysis.Using the relations nV E = V M and nI M = I E the respective values of their birefringence partners can also be calculated for every thin cylindrical layer r i .Finally, we obtain the actual fields via the relations Fr I M + βV E A very useful field component for optical fibers is the value of the overall electric field at any thin cylindrical layer of average radius r that can be calculated by the formula: After some algebra, this leads to the formula In the next section, we extend our analysis in certain UOF cases.

Unconventional Fibers
The refractive index n(r,ϕ) of the fiber with a UOF profile in general can be described as a function of both r and ϕ.Each cylindrical layer of an average radius r is considered to have a local value η (ϕ) for r = r 1 +r 2 2 .Again, we make use of the generic form of Maxwell equations as in Equation (2) of the previous section.Fourier Transforming the first vector equation gives the same set of equations as in (6), since there is no difference with the UOF case, however the second Maxwell vector equation due to the presence of the general n function should now be written as Appl.Sci.2019, 9, 270 10 of 21 The symbol ⊗ means convolution arising by the product of two functions of the variable ϕ.In the following paragraphs it will be shown how to escape this mathematical difficulty for the usual unconventional optical fibers.

UOF with Non-Circular, Non-Symmetric, or Eccentric Cores
For unconventional fibers of non-circular cores there are a set of circular layers where the refractive index varies between the inner and outer core and cladding values, respectively.In any such case, the function n(ϕ) 2 is a sum of a steady component n 2 and a periodic function of ϕ of period 2π thus can be written as a Fourier series n(ϕ) 2 = n 2 + ∑ +∞ −∞ N k exp(jkϕ).Taking into consideration that the convolution of the product of an exponential function exp(jkϕ) with any function A(ϕ) of a Fourier Transform A(l) is equal to A(l + k), i.e., the convolution generates "harmonics".The function n(ϕ) 2 is in a set of cylindrical thin layers, a sum of step functions alternating between the values n 2 1 and n 2 2 , where n 1 and n 2 are refractive indexes of core and cladding.Considering that in optical fibers the refractive indices of core and cladding are very close, one effectively has that (n 1 − n 2 )/n 1 1.As a result, any harmonic factors N k of the function n(ϕ) 2 are negligible in comparison to its steady component n 2 and can be omitted.
As an example, the harmonics become maximal for equal alternation steps.In this case, the first harmonic, that has the maximum value of all harmonics, is equal to π , while the steady component n 2 equals (n 2 1 + n 2 2 )/2 and we may make an approximation as 1. Thus, for optical fibers we can always assume that n(l) 2 n 2 .Then the system (18) will become equivalent to the following We can then follow the analysis that we did with the conventional fibers, where n 2 is the average value of the η 2 (ϕ) of each layer along ϕ in the [0, 2π] interval.

Application to Elliptic Core Fibers
The method was applied in the calculation of fundamental modes of a fiber of elliptic core of a and b major and minor semi-axis, respectively, with refractive index n 1 = 1.54, and a cladding value of n 2 = 1.47 (Figure 5) for various wavelengths (defined by various V factor values V = 2πb n 2 1 − n 2 2 ) and four ellipticity ratios a/b = 1.1, 1.3, 1.5, and 2.0.Results are presented in tabulated format (Tables 1-4) compared with previous results calculated with Mathieu functions together with differences and relative differences showing a deviation which goes as only 0.01/0.123%on the average.Note the results in Tables 1-4 are not normalised the same way as in Figure 4, they are simply the odd modes, b 11 = β o /k o .We can see that for small ellipticities especially, the results compare very well for all V values quite well with the Mathieu Functions results.As the ellipticity becomes very large as in Table 4, the results begin to differ.However, the accuracy of the Mathieu Functions used in this case is not known precisely, so the trend is correct and the actual difference could be debated.
The steady component of the refractive index for the calculations for each radius r is defined as n1 for r < b, n2 for r > a, and as ( 1  1 +  2  2 )/ when b < r < b, where φ1, φ2 are the arcs of the circle of radius r, inside and outside the ellipse in the upper semi ellipse.The steady component of the refractive index for the calculations for each radius r is defined as n 1 for r < b, n 2 for r > a, and as (n 1 ϕ 1 + n 2 ϕ 2 )/π when b < r < b, where ϕ 1 , ϕ 2 are the arcs of the circle of radius r, inside and outside the ellipse in the upper semi ellipse.
In the following Figure 6, the β diagram of the fundamental even mode of an elliptic fiber with semi axis ratio aa/bb = 2, n 1 = 1.54 and n 2 = 1.47 and variable factor defined by: ) is shown for completeness in agreement with [32].

Application to a Rectangular Core Fiber
The method was applied also in the calculation of fundamental modes of a fiber with an rectangular core, Figure 7

Application to a Rectangular Core Fiber
The method was applied also in the calculation of fundamental modes of a fiber with an rectangular core, Figure 7, with aa and bb semi-sides, with refractive index n 1 = 1.54, and a cladding of refractive index n 2 = 1.47, for various wavelengths, defined by various V factor values and four ratios a/b = 1.1, 1.3, 1.5, and 2.0.

Application to a Rectangular Core Fiber
The method was applied also in the calculation of fundamental modes of a fiber with an rectangular core, Figure 7, with aa and bb semi-sides, with refractive index n1 = 1.54, and a cladding of refractive index n2 = 1.47, for various wavelengths, defined by various V factor values  =  •  0 • √ 1 2 −  2 2 and four ratios a/b = 1.1, 1.3, 1.5, and 2.0.Birefringence results between the elliptical and rectangular waveguides are compared in tabular form in Tables 5-8 with equal semi axis.The constant component of the refractive index for calculating with each radius r is defined as where ϕ 1 , ϕ 2 being the arcs of the circle of radius r, inside and outside the rectangle in the upper half plane.For the results in the Tables 5-8 we have defined birefringence as simply the difference in the orthogonal fundamental modes normalized to k 0 , for the rectangular and elliptical waveguides.We can see the elliptical waveguide has greater birefringence in all cases of ellipticity and for all wavelengths (represented by the V value).
To make a direct comparison, the average refractive indexes as functions of r, for an elliptic core fiber and for a rectangular core fiber of the same aa and bb and aa/bb = 2, are shown in Figure 8.In Appendix B, the MATLAB codes for the average refractive indexes of an elliptic and a rectangular core fiber are given.In Figure 8 we observe that the average refractive index of the elliptical waveguide is dropping off faster than the rectangular and this facilitates the explanation of the stronger birefringence for this case.

The PCF Case
Ιn the case of a holey fiber, we separate the whole fiber circular cross-section into a set of thin cylindrical layers, variable η along φ extending beyond the cladding to take into account the surrounding air with η = 1.Each layer's thickness is δr = r1 − r2.We can then approximate n(r,φ) as n(φ) for the average r <r> = r + δr/2.The refractive index can be written as a Fourier series, i.e., as Taking into account the properties of the Fourier Transform we see that so that the expressions in the second terms of Equation ( 18) spread around a spectrum of harmonics.This is also to be understood as a result of successive scatterings from the bored air holes.We can now use the natural geometry of the usual hexagonal lattice to see that for each set of holes we can have either 6k harmonics.For the fundamental harmonic of l′ = 1, the derived harmonics passing through a layer of 6k holes should be 6k + 1.Thus, for the fundamental wave crossing the successive layers it "sees" a different set of periodic rectangle functions that will be shown rigorously to contribute a different number of harmonics (7, 13, 19, …).For a common harmonic to pass through, one must then take an integer product which leads to higher and higher harmonics, thus cutting out the entire spectrum apart from the last highest frequency.We conclude that for holey optical fibers, the approximation for any of its cylindrical thin layers, 2 ( , ) ( ) r l r η η η ≈ = suffices for further analysis of the resulting equations.Thus, the original system (18) becomes For the usual hexagonal pattern of holes, we may utilize elementary analytical geometry to derive the two separate regions where the refractive index alternates between the air refractive index n = 1 value and the higher value of the crystal material.We assume that along each separate layer a large circle corresponding to each cylindrical shell of radius r from the center of the fiber to the center of a smaller hole of radius r << r0 is cut while moving clockwise along the large circle.

The PCF Case
In the case of a holey fiber, we separate the whole fiber circular cross-section into a set of thin cylindrical layers, variable η along ϕ extending beyond the cladding to take into account the surrounding air with η = 1.Each layer's thickness is δr = r 1 − r 2 .We can then approximate n(r,ϕ) as n(ϕ) for the average r <r> = r + δr/2.The refractive index can be written as a Fourier series, i.e., as Taking into account the properties of the Fourier Transform we see that FT(exp(jlφ) • f (ϕ)) = f (l + l ) so that the expressions in the second terms of Equation ( 18) spread around a spectrum of harmonics.This is also to be understood as a result of successive scatterings from the bored air holes.We can now use the natural geometry of the usual hexagonal lattice to see that for each set of holes we can have either 6k harmonics.For the fundamental harmonic of l = 1, the derived harmonics passing through a layer of 6k holes should be 6k + 1.Thus, for the fundamental wave crossing the successive layers it "sees" a different set of periodic rectangle functions that will be shown rigorously to contribute a different number of harmonics (7, 13, 19, . . . ).
For a common harmonic to pass through, one must then take an integer product which leads to higher and higher harmonics, thus cutting out the entire spectrum apart from the last highest frequency.We conclude that for holey optical fibers, the approximation for any of its cylindrical thin layers, η(r, l) 2 ≈ η 2 (r) = η 2 suffices for further analysis of the resulting equations.Thus, the original system (18) becomes For the usual hexagonal pattern of holes, we may utilize elementary analytical geometry to derive the two separate regions where the refractive index alternates between the air refractive index n = 1 value and the higher value of the crystal material.We assume that along each separate layer a large circle corresponding to each cylindrical shell of radius r from the center of the fiber to the center of a smaller hole of radius r << r 0 is cut while moving clockwise along the large circle.
Prescribing a set of circles of successive radii r for each of which we can find the air holes (in 1/6 angle of the PCF) which are cut by the particular radius each time.Each arc is computed inside its respective air hole and the total sum of them divided by π/3 expresses the average squared refractive index.As a matter of fact, the square of the refractive index in this sum is equal to one, while the refractive index in the rest arc is the square of the silica refractive index.Hence the average refractive index can be easily calculated along r.In Figures 9 and 10, we show the average refractive index and the electric field of a hexagonal PCF, of n = 1.46 with a twelve layers lattice, as functions of the reduced distance from the center of the fiber for the fundamental mode.Figure 8 was generated by a MATLAB code for air-hole diameter equal to 0.8 of the air-hole distance and the air hole diameter was 3.14 times the transmitted wave length (V = 1.607847).We also notice the parametrization used as = (Λ − d/2) × 2 × π n 2 1 − n 2 2 ), n 1 for the silica refractive index, n 2 minimum refractive index = 1.123,Λ for the reduced air hole distance, d for the reduced air hole diameter, and Λ − d/2 = reduced inner core of PCF.
the electric field of a hexagonal PCF, of n = 1.46 with a twelve layers lattice, as functions of the reduced distance from the center of the fiber for the fundamental mode.Figure 8 was generated by a MATLAB code for air-hole diameter equal to 0.8 of the air-hole distance and the air hole diameter was 3.14 times the transmitted wave length (V = 1.607847).We also notice the parametrization used as = ( − /2) × 2 × � 1 2 −  2 2 ), n1 for the silica refractive index, n2 ≃ minimum refractive index = 1.123,Λ for the reduced air hole distance, d for the reduced air hole diameter, and  − /2 = reduced inner core of PCF.
For completeness we also show Figure 11, which is the β/V diagram of the even fundamental for a PCF.This shows how easily we can produce some very useful results for a series of unconventional fibers using the transmission line method we have developed.
In Appendix B, the MATLAB code for the average refractive index of a hexagonal photonic crystal fiber is also given.
In Appendix B, the MATLAB code for the average refractive index of a hexagonal photonic crystal fiber is also given.For completeness we also show Figure 11, which is the β/V diagram of the even fundamental for a PCF.This shows easily we can produce some very useful results for a series of unconventional fibers using the transmission line method we have developed.In Appendix B, the MATLAB code for the average refractive index of a hexagonal photonic crystal fiber is also given.
In Appendix B, the MATLAB code for the average refractive index of a hexagonal photonic crystal fiber is also given.

Conclusions
The presented resonance technique can be used for the study of unconventional fibers, i.e., fibers with cores of any shape, as long as the difference between core and cladding refractive indices is sufficiently small which holds true for almost all the monomode and holey fibers.The unconventional case is proven reducible to the same technique of conventional fibers, where for each mode order l we can approximate by a set of two, independent and non-homogeneous, resonant transmission lines (RTLs), each one representing one mode of the birefringence.
The simulation of unconventional fibers with RTLs gives a new, simple, and effective method for computing the eigenvalues of the RTLs representing the various modes of the holey fibers.Furthermore, for each eigenvalue, the average values of E.M. fields for every thin cylindrical layer of radius r of the unconventional fiber is directly computable from the relevant eigenfunctions of the RTLs.

Figure 1 .
Figure 1.Example of the alternating character of the local index value from a thin shell radial discretization.

Figure 1 .
Figure 1.Example of the alternating character of the local index value from a thin shell radial discretization.

Figure 2 .
Figure 2. Schematic depiction of the alternating step index resulting from radial discretization in a standard photonic crystal fiber (PCF).

Figure 2 .
Figure 2. Schematic depiction of the alternating step index resulting from radial discretization in a standard photonic crystal fiber (PCF).

H
to have the same units (V/m), Maxwell equations become then

Figure 3 .
Figure 3.The equivalent quadrupole for each cylindrical sector.

Figure 3 .
Figure 3.The equivalent quadrupole for each cylindrical sector.

Figure 4 .
Figure 4. Normalized birefringence of two step-index fibers with different refractive indexes as functions of their parameters V.

Figure 4 .
Figure 4. Normalized birefringence of two step-index fibers with different refractive indexes as functions of their parameters V.

Figure 5 .
Figure 5. Elliptic fiber with three indicative elliptic thin layers.Inside the ellipse r < b (n = n 1 ), outside the ellipse r > a (n = n 2 ) and partly outside b < r < a (n 1 > n > n 2 ).

Table 1 .
Comparison of ellipse for the fundamental modes.

Figure 6 .
Figure 6.β-V diagram of the even fundamental mode an elliptic fiber of semi axis ratio aa/bb = 2 and core refractive index 1.54 and cladding index 1.47.

Figure 6 .
Figure 6.β-V diagram of the even fundamental mode an elliptic fiber of semi axis ratio aa/bb = 2 and core refractive index 1.54 and cladding index 1.47.

Figure 6 .
Figure 6.β-V diagram of the even fundamental mode an elliptic fiber of semi axis ratio aa/bb = 2 and core refractive index 1.54 and cladding index 1.47.

Figure 8 .
Figure 8.Average refractive indexes of circular thin layers of elliptic and rectangular core fibers.

Figure 8 .
Figure 8.Average refractive indexes of circular thin layers of elliptic and rectangular core fibers.

Figure 10 .
Figure 10.Average electric field of the even fundamental mode along the radial distance of the PCF.

Figure 10 .
Average electric field of the even fundamental mode along the radial distance of the PCF.

Figure 11 .
Figure 11.β/V diagram of the even fundamental for a PCF.

Figure 11 .
Figure 11.β/V diagram of the even fundamental for a PCF.

Table 1 .
Comparison of ellipse for the fundamental modes.

Table 2 .
Comparison of ellipse for the fundamental modes.

Table 3 .
Comparison of ellipse for the fundamental modes.

Table 4 .
Comparison of ellipse for the fundamental modes.

Table 5 .
Comparison of birefringence for elliptical and orthogonal core fibers.

Table 6 .
Comparison of birefringence for elliptical and orthogonal core fibers.

Table 7 .
Comparison of birefringence for elliptical and orthogonal core fibers.

Table 8 .
Comparison of birefringence for elliptical and orthogonal core fibers.