## Abstract

Chromatic dispersion characteristics of nonlinear photonic crystal fibers are, for the first time to our knowledge, theoretically investigated. A self-consistent numerical approach based on the full-vector finite-element method in terms of all the components of electric fields is described for the steady-state analysis of axially-nonsymmetrical nonlinear optical fibers. Electric fields obtained with this approach can be directly utilized for evaluating nonlinear refractive index distributions. To eliminate nonphysical, spurious solutions and to accurately model curved boundaries of circular air holes, curvilinear hybrid edge/nodal elements are introduced. It is found from the numerical results that under high optical intensity, chromatic dispersion characteristics become different from those of the linear state due to optical Kerr-effect nonlinearity, especially in short wavelength region.

© 2003 Optical Society of America

## 1. Introduction

It is well-known that there are many interesting phenomena in optical fibers arising from nonlinear optical effects, such as optical soliton transmission, optical pulse compression, and generation of ultrashort pulse trains by a modulational instability. Chromatic dispersion is one of the most important parameters ruling these phenomena and under high optical intensity, it becomes different from that of the linear state due to optical Kerr-effect nonlinearity [1–3].

Photonic crystal fibers (PCFs) [4], also called holey fibers (HFs), consisting of a central defect region surrounded by multiple air holes running parallel to the fiber length have been one of the most interesting recent developments in fiber optics. Such fibers possess numerous unusual properties such as a wide single-mode wavelength range [5], a bend-loss edge at short wavelengths [5], a very large [6] or small [7] effective core area, anomalous group-velocity dispersion at visible and near-infrared wavelengths [8], and strong wavelength-dependent beam divergence [9]. Furthermore, as the effective core area can be controlled freely, extensive applications are expected, such as soliton transimission in the visible wavelength region [10], generation of super-continuum with short fiber lengths [11], and so on. But, also in HF, chromatic dispersion may change under high optical intensity, as in conventional optical fiber.

In this paper, the influence of the optical Kerr-effect on chromatic dispersion characteristics of PCFs are, for the first time to our knowledge, theoretically investigated. A self-consistent numerical approach based on the full-vector finite-element method (VFEM) for the steady-state analysis of axially-nonsymmetrical nonlinear optical fibers is newly formulated. It is found from the numerical results that under high optical intensity, chromatic dispersion characteristics become different from those of the linear state due to optical Kerreffect nonlinearity, especially in short wavelength region.

## 2. Full-vector finite element method

Various methods have been proposed for the steady-state modal analysis of nonlinear optical fibers. In early days, approximation methods based on the Gaussian approximation and the equivalent-step-index approximation have been proposed [12–16]. Although these methods are very simple and result in quick solutions, errors are inevitable when the trial fields deviate from the true modal fields. Hence, powerful numerical methods have been required for modeling nonlinear optical fibers and under the weakly guiding approximation, a finite element method (FEM) with one-dimensional, line elements has been adopted to circularly symmetric nonlinear optical fibers [1,3]. Unfortunately, this type of FEM cannot be applied to PCFs because of their axially-nonsymmetrical geometry. Furthermore, a full-vector model [17–21] is crucial for accurately predicting the modal properties, owing to the large refractive-index difference between silica materials and air holes. Although various approaches based on the VFEM have been effectively utilized for investigating the unique properties of PCFs [17–21], their applications are limited only to linear one.

Not only in PCFs but in optical channel waveguides, the guided waves are truly hybrid and for the steady-state analysis of the nonlinear optical channel waveguides, various types of VFEM in terms of all the magnetic or electric field components [22–28] have been developed using the conventional, in usual, triangular elements so far. In these approaches, to eliminate spurious solutions, the penalty function method with an artificial coefficient (penalty coefficient) has been introduced and thus, the accuracy of solutions depends on it. Recently, Obayya *et al.* [29] have proposed an imaginary distance finite-element beam propagation method in terms of transverse magnetic fields. In this approach, line integrals are involved in the weak form corresponding to the partial difference equation for the slowly-varying transverse magnetic fields. The evaluation of the line integrals seems to be very cumbersome, because the line integral occurs on every interface between elements with different media. Furthermore, when using magnetic fields as the working variables [22–27,29], the electric fields necessary for evaluating nonlinear refractive-index distributions cannot be directly calculated. In [24,28], although the working variables are the electric fields, complicated boundary conditions between different media must be imposed because the conventional nodal elements are used for discretizing the waveguide cross section.

Here, we newly formulate a VFEM in terms of all the components of electric fields for the steady-state analysis of nonlinear optical waveguides with arbitrarily shaped cross section. In the present formulation, by using hybrid edge/nodal elements [30], spurious solutions are eliminated, artificial coefficients [22–28] and line integral terms [29] are not included, and complicated boundary conditions [24,28] described above are not necessary. Electric fields obtained with this approach can be directly utilized for evaluating nonlinear refractive index distributions.

## 3. Basic equation

We consider a nonlinear optical waveguide with arbitrary cross section in the *xy* (transverse) plane and assume that the structure is uniform along the propagation direction (*z* axis).

From Maxwell’s equations the following vectorial wave equation is derived

where ** E** is the electric field vector,

*k*

_{0}is the free-space wavenumber, and

*n*is the refractive index. Because of the unifomity of the waveguide, we can write the electric fields

**as**

*E*where β is the propagation constant. Also, the refractive index *n* depends on the electric field intensity and is given by

## 4. Self-consistent approach

We divide the waveguide cross section into high-order hybrid edge/nodal elements [30] as shown in Fig. 1 and expand the electric fields as

$$\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}=\left({\mathbf{i}}_{x}{\left\{U\right\}}^{T}+{\mathbf{i}}_{y}{\left\{V\right\}}^{T}\right){\left\{{e}_{t}\right\}}_{e}\mathrm{exp}\left(-j\mathit{\beta z}\right)\mathbf{+}{\mathbf{i}}_{z}j\mathit{\beta}{\left\{N\right\}}^{T}{\left\{{e}_{z}\right\}}_{e}\mathrm{exp}\left(-j\mathit{\beta z}\right)$$

within each element, where the superscript^{T} denotes a transpose and

{*e _{t}*}

_{e}and {

*e*}

_{z}_{e}edge and nodal electric fields for each element, respectively;

{*U*} and {*V*} shape function vectors for edge elements;

{*N*} shape function vector for nodal elements.

Applying the finite element technique to Eq. (1), we obtain

with

$$\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}-\frac{\partial \left\{V\right\}}{\partial x}\frac{\partial {\left\{V\right\}}^{T}}{\partial x}+{k}_{0}^{2}{n}^{2}\left\{U\right\}{\left\{U\right\}}^{T}+{k}_{0}^{2}{n}^{2}\left\{V\right\}{\left\{V\right\}}^{T}]$$

$$\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}\phantom{\rule{.2em}{0ex}}=\sum _{e}{\iint}_{e}\left[\left\{U\right\}\frac{\partial {\left\{N\right\}}^{T}}{\partial x}+\left\{V\right\}\frac{\partial {\left\{N\right\}}^{T}}{\partial y}\right]\mathit{dx}\mathit{dy}$$

where {0} is a null vector, [0] is a null matrix, and Σ_{e} extends over all different elements. The optical power *P* is evaluated as

$$\phantom{\rule{.2em}{0ex}}=\frac{\beta}{2{k}_{0}{Z}_{0}}\left({\left\{{e}_{t}\right\}}^{T}\left[{M}_{\mathit{tt}}\right]\left\{{e}_{t}\right\}+{\left\{{e}_{t}\right\}}^{T}\left[{M}_{\mathit{tz}}\right]\left\{{e}_{z}\right\}\right)$$

where *Z*
_{0} is the free-space impedance and * denotes complex conjugate. To obtain the intensity-dependent refractive index in Eq. (3), it is necessary to compute the actual electric field {*e*}. The relation between the actual field {*e*} and the eigenvector {*e*’} of Eq. (5) can be written as

with

where {*e _{t}*’} and {

*e*’} are the edge and nodal electric fields obtained by solving Eq. (5), respectively, and the dash is used for discriminating the solutions of Eq. (5) from the actual ones related to the optical power

_{z}*P*in Eq. (9).

Equation (5) is a nonlinear eigenvalue problem whose eigenvalue and eigenvector correspond to β^{2} and {*e*’}, respectively. Hence, one can solve it self-consistently using the following iterative scheme:

- Specify the refractive index
*n*, the wavelength λ=2π/*k*_{0}, and optical power*P*as input data. - To assign intial values to β and {
*e*’}, solve the corresponding linear problem. - To obtain a new set of β and {
*e*’}, solve Eq. (5). - Iterate the above procedures 3) and 4) until the solution converges within the desired criterion.

## 5. Chromatic dispersion characteristics of nonlinear holey fibers

#### 5.1 Nonlinear optical fibers

To evaluate chromatic dispersion, the second order derivative of the propagation constant with respect to frequency is needed and here, numerical differentiation is introduced. Chromatic dispersion of a circular-core optical fiber can be exactly calculated but under high nonlinearity, there is no analytical solution, even though the circular core is assumed. To show the validity and reliability of our approach, we consider a step-index fiber with nonlinear circular core, the group velocity dispersion of which is given by

where *V* and *b* are, respectively, the well-known normalized frequency and normalized propagation constant. For comparison, the fiber parameters are chosen to have the same values as in [1–3], namely, the core radius *a*=2.5 µm, the cladding index *n _{cl}*=1.47, and
$a\sqrt{{n}_{\mathit{co}}^{2}-{n}_{\mathit{cl}}^{2}}=0.22\phantom{\rule{.2em}{0ex}}\mathrm{\mu m}$
with

*n*being the linear part of core index. Also,

_{co}*n*=6.4×10

_{2}P^{-14}m

^{2}is assumed, where

*n*

_{2}[m

^{2}/W] is the nonlinear coefficient,

*P*[W] is the optical power in Eq. (9), and the refractive index of nonlinear core

*n*is given by

The actual optical power *P* considered here is 200 kW when the fiber is made of silica glass with *n*
_{2}=3.2×10^{-20} m^{2}/W and is less than 1 W when the fiber is made of organic material with *n*
_{2}≈10^{-14} m^{2}/W [1,2].

Figure 2 shows the group velocity dispersion *g*(*V*). The calculated results for the linear case agree well with the exact solutions. Also, under high nonlinearity, our results are in good agreement with those of FEM for circulary symmetric fibers [3], indicating that our approach is appropriate for investigating chromatic dispersion characteristics of various nonlinear optical fibers.

#### 5.2 Nonlinear holey fibers

We consider a PCF, also called HF, as shown in Fig. 3, where *d* is the hole diameter, Λ is the hole pitch (center-to-center distance between the holes) of triangular lattice structure, and the background nonlinear refractive index *n* is given by

with *n _{L}* being the linear part of the background index. Because of the symmetry nature of the system, only one-quarter of the fiber cross section is divided into curvilinear edge/nodal hybrid elements shown in Fig. 1. Here, the background material is assumed to be silica and the material dispersion is taken into account by using the three-term Sellmeier equation.

Figures 4(a), (b), (c), and (d) show the chromatic dispersion of nonlinear holey fibers with *d*/Λ=0.9 for different hole pitches of Λ=1.0 µm, 1.5 µm, 2.0 µm, and 2.5 µm, respectively, where the values of *n*
_{2}
*P* is taken as a parameter. The influence of optical Kerr-effect on chromatic dispersion is larger in the shorter wavelength region and the zero-dispersion wavelength shifts to the longer wavelength region. For the shorter wavelength region, the light is strongly confined into the core region more and more with increasing optical power and so, the effect of air holes on the waveguide dispersion becomes smaller. In the longer wavelength region, on the other hand, the light confinement is not so strong, in other words, the effective core area becomes large, and so, the chromatic dispersion tends to be insensitive to the change of nonlinearity, compared to that in the shorter wavelength region. To confirm these facts, the effective core area *A _{eff}* defined as

is shown in Fig. 5, where Λ=1.5 µm and *d*/Λ=0.9. It can be seen clearly that in the longer wavelength region, *A _{eff}* becomes large, resulting in reducing the effective nonlinearity. As the optical power increases, because of stronger field confinement,

*A*becomes smaller and the effective nonlinearity is enhanced.

_{eff}Figures 6(a), (b), (c), and (d) show the chromatic dispersion of nonlinear holey fibers with Λ=1.5 µm for different ratios of hole diameter to hole pitch of *d*/Λ=0.5, 0.6, 0.7, and 0.8, respectively, where the values of *n*
_{2}
*P* is taken as a parameter. As the value of *d*/Λ increases, because of stronger field confinement, the effective nonlinearity is enhanced and the change of zero-dispersion wavelength becomes large. Figures 7(a) and (b) show the zero-dispersion wavelength as a function of hole pitch Λ(*d*/Λ=0.9) and of ratio of hole diameter to hole pitch *d*/Λ(Λ=1.5 µm), respectively, where the value of *n*
_{2}
*P* is taken as a parameter. We can see that the zero-dispersion wavelength shifts greatly to the longer wavelength region with increasing optical power, especially in smaller values of Λ and in larger values of *d*/Λ.

## 6. Conclusion

We have investigated, for the first time, chromatic dispersion characteristics of nonlinear PCFs theoretically. To analyze nonlinear PCFs accurately, a self-consistent VFEM in terms of all the components of electric fields is newly formulated for nonlinear optical waveguides with arbitrary cross section. To eliminate nonphysical, spurious solutions and to accurately model curved boundaries of circular air holes, curvilinear hybrid edge/nodal elements are introduced. It is found from the numerical results that under high optical intensity, chromatic dispersion characteristics become different from those of linear state due to optical Kerr-effect nonlinearity, especially in short wavelength region. As the hole pitch is smaller and the ratio of hole diameter to hole pitch is larger, because of stronger field confinement, the effective core area becomes very small, resulting in enhanced effective nonlinearity.

## References and links

**1. **K. Okamoto and E.A.J. Marcatili, “Chromatic dispersion characteristics of fibers with optical Kerr-effect,” J. Lightwave Technol. **7**, 1988–1994 (1989). [CrossRef]

**2. **R.A. Sammut and C. Pask, “Group velocity and dispersion in nonlinear-optical fibers,” Opt. Lett. **16**, 70–72 (1991). [CrossRef]

**3. **H.Y. Lin and H.-C. Chang, “An efficient method for determining the chromatic dispersion characteristics of nonlinear single-mode optical fibers,” J. Lightwave Technol. **10**, 1188–1195 (1992). [CrossRef]

**4. **J.C. Knight, T.A. Birks, P.St.J. Russell, and D.M Atkin, “All-silica single-mode optical fiber with photonic crystal cladding,” Opt. Lett. **21**, 1547–1549 (1996). [CrossRef]

**5. **T.A. Birks, J.C. Knight, and P.St.J. Russell, “Endlessly single-mode photonic crystal fiber,” Opt. Lett. **22**, 961–963 (1997). [CrossRef]

**6. **J.C. Knight, T.A. Birks, R.F. Cregan, P.St.J. Russell, and J.-P. de Sandro, “Large mode area photonic crystal fiber,” Electron. Lett. **34**, 1347–1348 (1998). [CrossRef]

**7. **N.G.R. Broderick, T.M. Monro, P.J. Bennett, and D.J. Richardson, “Nonlinearity in holey optical fibers: measurement and future opportunities,” Opt. Lett. **24**, 1395–1397 (1999). [CrossRef]

**8. **M.J. Gander, R. McBride, J.D.C. Jones, D. Mogilevtsev, T.A. Birks, J.C. Knight, and P.St.J. Russell, “Experimental measurement of group velocity in photonic crystal fiber,” Electron. Lett. **35**, 63–64 (1999). [CrossRef]

**9. **M.J. Gander, R. McBride, J.D.C. Jones, T.A. Birks, J.C. Knight, P.St.J. Russell, P.M. Blanchard, J.G. Burnett, and A.H. Greenaway, “Measurement of the wavelength dependence of beam divergence for photonic crystal fiber,” Opt. Lett. **24**, 1017–1019 (1999). [CrossRef]

**10. **W.J. Wadsworth, J.C. Knight, A. Ortigosa-Blanch, J. Arriaga, E. Silvestre, and P.St.J. Russell, “Soliton effects in photonic crystal fibres at 850 nm,” Electron. Lett. **36**, 53–55 (2000). [CrossRef]

**11. **W.J. Wadsworth, A. Ortigosa-Blanch, J.C. Knight, T.A. Birks, T.-P. Martin Man, and P.St.J. Russell, “Supercontinuum generation in photonic crystal fibers and optical fiber tapers: a novel light source,” J. Opt. Soc. Am. B **19**, 2148–2155 (2002). [CrossRef]

**12. **A.W. Snyder, Y. Chen, L. Poladian, and D.J. Mitchell, “Fundamental mode of highly nonlinear fibers,” Electron. Lett. **26**, 643–644 (1990). [CrossRef]

**13. **R.A. Sammut and C. Pask, “Variational approach to nonlinear waveguides-gaussian approximations,” Electron. Lett. **26**, 1131–1132 (1990). [CrossRef]

**14. **R.A. Sammut and C. Pask, “Gaussian and equivalent-step-index approximations for nonlinear waveguides,” J. Opt. Soc. Am. B **8**, 395–402 (1991). [CrossRef]

**15. **Y. Chen, “Nonlinear fibers with arbitrary nonlinearity,” J. Opt. Soc. Am. B **8**, 2338–2341 (1991). [CrossRef]

**16. **M.J. Holmes, D.M. Spirit, and F.P. Payne, “New gaussian-based approximation for modeling non-linear effects in optical fibers,” J. Lightwave Technol. **12**, 193–201 (1994). [CrossRef]

**17. **F. Brechet, J. Marcou, D. Pagnoux, and P. Roy, “Complete analysis of the characteristics of propagation into photonic crystal fibers,” Opt. Fiber Technol. **6**, 181- (2000). [CrossRef]

**18. **M. Koshiba and K. Saitoh, “Numerical verification of degeneracy in hexagonal photonic crystal fibers,” IEEE Photon. Technol. Lett. **13**, 1313–1315 (2001). [CrossRef]

**19. **M. Koshiba, “Full-vector analysis of photonic crystal fibers using the finite element method,” IEICE Trans. Electron. **E85-C**, 881–888 (2002).

**20. **K. Saitoh and M. Koshiba, “Full-vectorial imaginary-distance beam propagation method based on a finite element scheme: application to photonic crystal fibers,” IEEE J. Quantum Electron. **38**, 927–933 (2002). [CrossRef]

**21. **A. Cucinotta, S. Selleri, L. Vincetti, and M. Zoboli, “Holey fiber analysis through the finite element method,” IEEE Photon. Technol. Lett. **14**, 1530–1532 (2002). [CrossRef]

**22. **K. Hayata and M. Koshiba, “Full vectorial analysis of nonlinear-optical waveguides,” J. Opt. Soc. Am. B **5**, 2494–2501 (1988). [CrossRef]

**23. **R.D. Ettinger, F.A. Fernandez, B.M.A. Rahman, and J.B. Davies, “Vector finite element solution of saturable nonlinear strip-loaded optical waveguides,” IEEE Photon. Technol. Lett. **3**, 147–149 (1991). [CrossRef]

**24. **X.H. Wang and G.K. Cambrell, “Full vectorial simulation of bistability phenomena in nonlinear-optical channel waveguides,” J. Opt. Soc. Am. B **10**, 1090–1095 (1993). [CrossRef]

**25. **X.H. Wang and G.K. Cambrell, “Simulation of strong nonlinear effects in optical waveguides,” J. Opt. Soc. Am. B **10**, 2048–2055 (1993). [CrossRef]

**26. **S. Selleri and M. Zoboli, “An improved finite element method formulation for the analysis of nonlinear anisotropic dielectric waveguides,” IEEE Trans. Micorwave Theory Tech. **43**, 887–892 (1995). [CrossRef]

**27. **M. Zoboli, F.Di Pasquale, and S. Selleri, “Full-vectorial and scalar solutions of nonlinear optical fibers,” Opt. Comuun. **97**, 11–15 (1993). [CrossRef]

**28. **X.H. Wang and G.K. Cambrell, “Vectorial simulation and power-parameter characterization of nonlinear planar optical waveguides,” J. Opt. Soc. Am. B **12**, 265–274 (1995). [CrossRef]

**29. **S.S.A. Obayya, B.M.A. Rahman, K.T.V. Grattan, and H.A. El-Mikati, “Full-vectorial finite-element solution of nonlinear bistable optical waveguides,” IEEE J. Quantum Electron. **38**, 1120–1125 (2002). [CrossRef]

**30. **M. Koshiba and Y. Tsuji, “Curvilinear hybrid edge/nodal elements with triangular shape for guided-wave problems,” J. Lightwave Technol. **18**, 737–743 (2000). [CrossRef]