In recent years, researchers have placed significant emphasis on high-speed short-range data transmission using plastic optical fiber (POF) [1,2]. POF has the advantages of a large core and an easy connection, making it a potentially best option for the home network. Various materials are used in the manufacturing of POFs, the most common being polymethyl methacrylate (PMMA) [3–16]. Because the POF material is flexible, it can be used to create POFs with different materials or specifications to suit different application requirements. PMMA has been the most widely utilized material for POF manufacture [17]. POF can typically be categorized as step-index (SI) [18] or graded-index (GI) [19], based on the distribution of the refractive index (RI), and as single-mode [20] and multimode [21], depending on the number of propagation modes. A type of POF known as GI POF has an RI distribution that steadily decreases from the core axis to the cladding. The GI distribution of RI can reduce intermodal dispersion, enhance the bandwidth, and extend the transmission range of the fiber. However, intricate doping techniques are needed to create GI POF.

The 1990s saw the successful proposal of photonic crystal fiber (PCF) [22]. The flexibility of the optical fiber is greatly increased by the microstructure of the PCFs [23–26]. The first PMMA polymer photonic crystal fiber (PPCF) was created by Argyros in 2001 [27]. As a result of its various applications, PPCF sparked research interest [28,29]. Figure 1 shows a PPCF with a core of air holes of different sizes (d), simulating a GI optical fiber. Greater control over air-hole sizes and pitch Λ is the advantage of the GI PPCF over traditional GI POF, as opposed to the latter’s requirement for complex doping procedures. Moreover, it has been found that GI PPCF outperforms conventional GI POF in terms of bandwidth and loss [30].

The performance of the GI PPCF is significantly influenced by mode coupling. Light scattering, which occurs when random anomalies in multimode optical fibers transfer power from one mode to another, is the main cause of mode coupling. As the fiber length increases, power distribution varies until an equilibrium mode distribution (EMD) is formed at “coupling length” L _c . The fiber length at which the highest-order guiding mode altered its distribution to m = 0 is indicated by the coupling length L _c at which EMD is attained. Light is evenly dispersed, and the coupling process is practically finished beyond L _c . Since the steady-state distribution (SSD) was developed, each distribution that is released has a distinct far-field pattern. In other words, length z _s indicates the fiber length at which the output angular power distribution becomes completely independent of the launch beam. Mode coupling reduces modal dispersion and increases the transmission bandwidth [30]. It is also noteworthy that mode coupling makes it impossible to precisely characterize bandwidth and attenuation unless the SSD is fully obtained. Therefore, knowing the fiber length at which an SSD is constructed is essential.

To date, it has not been possible to examine the transmission characteristics of multimode PCFs with commercial simulation software applications. To tackle this problem, we present in this paper an efficient application of the Langevin equation to the GI PPCF mode-coupling problem. Thus, the intrinsically stochastic problem of mode coupling could be mathematically described stochastically using the Langevin equation. To the best of our knowledge, this is the first time that the Langevin equation has been used for the investigation of mode coupling in multimode GI PPCFs. For multimode GI PPCF, we found lengths for obtaining the EMD and SSD using launch beam distributions with various launch beam radial offsets. It is assumed that the air holes in the core and cladding are spaced in a regular triangular pitch (see Figure 1).

A GI PPCF considered in this work is shown in Figure 1. Rings 1, 2,. . . , 6 represent the six air-hole rings on the GI PPCF. A triangular lattice with a pitch of Λ is used to hold the air holes, and a polymer is considered the fiber material. The diameter of the air holes in rings 5 and 6 is the same as that of ring 4 (d ₄ = d ₅ = d ₆).

The following Equation 1 is the RI profile of the GI optical fiber [31]: n r , λ = n c o λ   1 − 2 Δ λ r a g 1 / 2   0 ≤ r ≤ a n c o λ   1 − 2 Δ λ 1 / 2 = n c l λ   r > a . (1) Here, a is the core radius, g is the core index exponent, n _co (λ) is the core index measured at the fiber axis, n _cl (λ) is the cladding index, and ∆ = n c o − n c l / n c o is the relative index difference.

The time-independent power flow equation (TI PFE) for the GI optical fiber is [31] ∂ P m , λ , z ∂ z = D m ∂ P m , λ , z ∂ m + D ∂ P 2 m , λ , z ∂ m 2 , (2) where P m , λ , z is the mth principal mode power, z is the coordinate along the fiber axis, and D is a mode-coupling constant. The maximum principal mode number is given in Equation 3 [31], as follows: M λ = g Δ λ g + 2 a k n c o λ , (3) where k = 2π/λ.

The principal mode m excited at the input fiber end is given in Equation 4 [31], as follows: m M = Δ r a g + θ 2 2 Δ g + 2 / 2 g , (4) where ∆ r is the launch beam radial offset and θ is the launch beam angle.

We first approximate Equation 2 as follows: ∂ P m , λ , z ∂ z = − V ∂ P m , λ , z ∂ m + D ∂ 2 P m , λ , z ∂ m 2 . (5)

Equation 5 can be understood as a special Fokker–Planck equation [32]. One can compute the drift coefficient V using Equation 6: V = 1 K ∑ i = 1 K V i , (6) where V _i is the drift coefficient of the ith principal mode. Later in this article, an illustration of the drift coefficient determination technique is provided.

The discretized Langevin equation can be obtained by transforming the Fokker–Planck Equation 5 [32,33], where principal mode m _n+1 at fiber length z _n+1 is given as m n + 1 = m n + V k + D k ω n , (7) where k = z _f /N, z _f is the fiber length, N is the number of finite steps of length k, n = 0, 1, . ., N ^-1, and ω ₀, ω₁ , . . ., ω _N-1 are independent Gaussian random numbers, with properties <ω _n >=0 and <ω _n ω _n’ ≥2δ _nn’ . For m n = 0 , Equation 7 reduces to Equation 8 [33]: m n + 1 = D k ω n . (8)

Thus, one obtains m _N = m(z _f ). By calculating a large number of representations of ω _n and averaging in appropriate intervals Δm for 0 ≤ m ≤ M , one obtains < m(z _f )>. It should be noted that optical fiber perturbations are known to be random in nature. Examples of these perturbations include variations in diameter, stress-induced microscopic random bends, and defects in the fiber core. Thus, the stochastic process of energy redistribution in optical fiber produced by its perturbations is explicitly described and modeled using the Langevin equation because of its stochastic nature.

The following is the effective V parameter for GI PPCF: V = 2 π λ a e f f n 0 2 − n f s m 2 , (9) where a_eff = Λ / 3   [34] and n_fsm is the effective RI of different core and cladding layers, which is obtained from Equation 9, with the effective parameter V [35]: V λ Λ , d Λ = A 1 + A 2 1 + A 3 ⁡ exp A 4 λ / Λ , (10) where parameters A_i (i = 1 to 4) are given in Equation 11, as follows: A i = a i 0 + a i 1 d Λ b i 1 + a i 2 d Λ b i 2 + a i 3 d Λ b i 3 . (11)

The coefficients a _{i 0}–a _{i 3} and b _{i 1}–b _{i 3} (i = 1–4) are given in our previous work [35].

We applied our method to the GI PPCF with the following parameters: a core radius a = 4Λ = 16 μm, pitch Λ = 4 μm, fiber diameter b = 1 mm, n _co = 1.5220, and n _cl = 1.4920 [31,36]. We used M = 24 at λ = 633 nm, g = 2.0, and ∆ = n c o − n c l / n c o = 0.019711. The parameter D = 1482 1/m is a typical value for GI PPCF and conventional GI POFs [31,36], and V=( − 1.5 ±0.1) 1/m. We highlight that when modeling the GI PPCF, the typical values of D that characterize a standard GI POF can be utilized as the degree of mode coupling in both standard GI POFs and GI PPCFs correlates with the polymer core material. This assumption mirrors the approach taken in silica PCF modeling [37].

For Λ = 4 μm and air-hole diameters d ₁ = 0.6 μm, d ₂ = 0.7 μm, d ₃ = 1.3 μm, and d ₄ = 3.1 μm, the refractive indices n ₁ = 1.5201, n ₂ = 1.5145, n ₃ = 1.5050, and n ₄ = 1.4920, respectively, are calculated using Equations 9, 10. These parameters are chosen in order to enable the GI distribution (1) with g = 2, which then results in the best transmission properties (such as bandwidth) of the investigated PPCF. The diameter of the cladding air holes in rings 5 and 6 is d ₄ = d ₅ = d ₆ = 3.1 μm, which corresponds to the cladding refractive index n ₄ = n ₅ = n ₆ = n _cl = 1.4920. In Figure 2, the normalized output modal power distribution P(m,λ,z) obtained by solving the Langevin equation is compared to the numerical solutions of the TI PFE [36] at different fiber lengths. For these calculations, a Gaussian beam P(θ,z) launched with θ = 0^o and radial offsets ∆r = 0, 4, 8, and 12 µm is used. A good agreement between the solutions of the Langevin equation and the TI PFE is shown in Figure 2. The EMD (Figure 2D) is established at a coupling length of L _c = 18 m. The SSD is observed at z≡z _s = 60 m.

It should be mentioned that for standard GI POF, which we previously investigated in our study [31], a coupling length of L _c = 31 m is published. The coupling coefficient for this type of GI POF is D = 1,482 1/m. The coupling length in GI PPCF is shorter than that in traditional GI POF due to the smaller core radius and fewer propagating modes (the maximum principal mode number in conventional GI POF was M = 656). In other words, for a shorter distance, fewer propagating modes must couple together. Comparing silica PCFs to the GI PPCF, which was the focus of this investigation, showed that their mode coupling is noticeably weaker and, therefore, much longer lengths L _c ≃1.45 to 1.65 km at which an EMD is achieved and z _s ≃3.30 to 3.80 km at which an SSD are established [35].

In summary, we demonstrated that mode coupling in GI POFs may be effectively treated by solving the Langevin equation (stochastic differential equation), which explicitly acknowledges a stochastic nature of energy redistribution in optical fiber produced by its perturbations. It is applicable to all GI optical fibers. This is not an issue with the Langevin equation, in contrast to the Fokker–Planck equation and TI PFE, which call for extra care in the stability of their numerical solutions. In terms of effectiveness, speed of execution, and memory usage, the Langevin equation is preferable. The Langevin equation does not have this issue, in contrast to the Fokker–Planck equation and the TI PFE, where a very fine mesh in the finite difference approach is required in order to obtain a highly accurate numerical solution (high memory consumption). The Langevin equation integration algorithm and the explicit finite difference method algorithm for the numerical solution of the TI PFE are evaluated in terms of their time efficiency (speed of execution), space efficiency (memory consumption), and complexity (solution/algorithm structure). The Langevin equation and the TI PFE take 1.5 and 2 min, respectively, to execute on an Intel(R) Core(TM) i3 CPU 540@3.07 GHz computer for the longest examined fiber length of 60 m. When expressed in terms of a 2-dim array, the memory consumption for the TI PFE and the Langevin equation is 24 × 6.0 × 10⁶ and 24 × 6.0 × 10⁵, respectively. Compared to the solution of the Langevin equation, the explicit finite difference solution of the TI PFE is more complicated. It is worth noting that the experimental setup for future experiments with the GI PPCF investigated in this work would be similar to that used in our previous work with a standard GI POF [31].

The behavior of mode coupling plays a crucial role in determining the length-dependent bandwidth of GI PCFs. Notably, the bandwidth decreases inversely proportional to lengths shorter than the coupling length L _c . Beyond this coupling length L _c , it exhibits a z−^1/2 dependency. A shorter L _c results in a quicker transition to a phase of reduced bandwidth decrease [37,38]. Compared to conventional GI POFs, GI PPCFs require shorter lengths to establish EMD and SSD, leading to a faster improvement in bandwidth enhancement [39]. This characteristic suggests that GI PPCFs are more suitable for short-range telecommunications. The findings of this research have practical implications for various communication and sensing systems utilizing multimode GI PPCFs. Understanding the modal distribution of multimode GI PPCFs at specific lengths is crucial for their integration into optical fiber sensing systems. One should also mention that, on the other hand, single-mode optical fibers are successfully used as a part of various fiber optic sensor systems [40–42]. Although single-mode and multimode optical fibers are both used in fiber optic sensor systems, they have different characteristics and are usually suited for different sensing applications. In other words, the choice between these two types of optical fibers for a particular sensing system is governed by their core diameter, distance, bandwidth, and light source.

This paper presents the numerical solution of the Langevin equation, which was used to investigate the state of mode coupling along a GI PPCF. The Langevin equation recognizes and explicitly accounts for the stochastic nature of the intrinsic perturbation effects of the GI PPCF. The results show that the length required to develop SSD and the coupling length required to obtain the EMD are both low in this fiber due to a strong mode-coupling process that is typical of POFs. One explanation for such substantial mode coupling is the large intrinsic perturbation effects in the GI PPCF. The GI PPCF under investigation in this study achieves the EMD at a length L _c that is even shorter than that in the conventional GI POF (L _c = 18 m in GI PPCFs compared to L _c = 31 m in conventional GI POFs). This is a result of the GI PPCF having fewer propagating modes due to its smaller core radius. In particular, a shorter length is required to complete the mode-coupling process when there are fewer propagating modes. Thus, a shorter L _c would result in a quicker shift to the slower bandwidth regime decrease. The fiber characterization provided in this paper is important for its use in data transmission, sensing, and power supply systems.