The nematic liquid crystal (NLC), which is one of the strong nonlocal nonlinear media (SNNM), has been widely used due to its cheap price. We choose it as the transmission media, and the even and odd Laguerre Gaussian beams as the incident beams. The even and odd Laguerre Gaussian beam is gaven by Eq. (5). Then, the split-step Fourier method is used to simulate the evolution of the odd Laguerre Gaussian beam in NLC. For convenience, we set the fundamental mode beam width a0 = 1 and the characteristic length σ = 20 in all numerical simulation. In addition, the windows of the system on x and y are − 30 to 30, which is discretized to 513 points.

**3.1 The unbounded state in mode ***p = 0, m ≥ 1.*

We choose the radial mode number p = 0, and then change the azimuthal mode number m. It is found that the beam profile of the odd Laguerre Gaussian beam in mode(0,1) maintian well, and can keep quasi periodic oscillation transmission for a long distances, as shown in Fig. 1(a). Meanwhile, the beam width of mode (0,1) do not change, as shown in Fig. 2(a). Therefore, we call it quasi-periodic oscillating solitons. Then, we observe the evolution of the beams in mode (0,2), (0,5). It is found that the beams oscillate quasi-periodically at the initial stage, then break up and evolve into a set of single-hump profiles to emit energy due to the scattering, as shown in Fig. 1(b) and Fig. 1.(c). The radiation waves arrive to infinity results in the increasing of the beam width as seen in Fig. 2(a). It is because the scattering of mode (0,5) is weaker than mode (0,2), the energy of mode (0,5) lost slower than mode (0,2) during the propagation. The beam width of mode (0,5) in Fig. 2(a) increases at a lower rate than that of mode (0,2). The beam width of mode (0, 6) also increases slower than that of mode (0,2). Therfore, we conclude that the larger the azimuthal mode number m is, the weaker of the scattering is, the energy of beam is not easy to be lost during the propagation.

Comparing Fig. 1 and Fig. 2(a), we can see that the beam of mode (0,1) is stable and can propagates as soliton in NLC. The beams of mode (0,2), (0,5), (0,6) are unstable in the unbounded state in fact. The proper input power, which is regard as the critical power Ppm, only can make the beams (*ψ*mp) propagate stably in initial stage. The beams break up rapidly after a short distance. Therefore, we use cross-correlation function to observe the change of the waveform. The cross-correlation function is as follows:

$$c(z)=\frac{1}{P}\int_{{ - \infty }}^{\infty } {\int_{{ - \infty }}^{\infty } {\psi (x,y,0){\psi ^*}(x,y,z)dxdy} }$$

7

where the superscript * denotes the conjugate complex, the value of |C| starts at 1.

Figure 2 (b) shows that the waveform does not deform along z for *m* = 1. Whereas for the case of *m* = 2, 5, 6, the waveform can not maintain and deform sinificantly. The position of the deformation can be seen clearly.

**3.2 The chaoticon state with ***p ≥ 1, m ≥ 0.*

Similarly, we increase the value of p and m to study the evolution characteristics of the beams in mode(1,m). We fix p = 1, then increase the value of m from 1 to 7 to study their beam width, waveform and beam profile. It can be seen from Fig. 3(a) and Fig. 3(b) that the beam width of mode (1,0) become wider during the propagation process, and its waveform has obvious deformation in a short distance after the initial stage. The broken beam profile of it in Fig. 4(a) lost energy quickly due to the scattering in unbounded state. The beam of mode (1,0) is unstable, and the evolution of it is very different with the beam of mode (0,1). Although the waveform deformation of modes (1,4) and (1,5) is more obvious, the beam width of them in Fig. 3 widen at a slower rate. Figure 4(b) and Fig. 4(c) show that their beam profile break up after a short quasi periodic oscillation, the energy of the broken beams are localized well and can be confined to a certain range during the propagation. These conclusions are consistent with those in Fig. 3. By comparing of Fig. 3 and Fig. 4, we conclude that a lager m is conductive to slowing down the beam broadening.

We continue to study the beams of modes (2,m) and (3,m), and find that the waveforms deformation will be more servere; but if m > 4, with the increase of m, the beam width becomes wider more slowly. Especially, it does not get wider at all in some higher modes, such as mode (1,6), (2,7), (3,8). The beam width of modes (1,6), (2,7), (3,8) can maintain well over longe distances, as shown in Fig. 5(a). The beam profile of them in Fig. 6 show that after a short stable quasi periodic oscillation, the energy of the broken beam are localized well in a certian range, and propagate over long distances with very little scattering. We analyze these beams as chaoticons with soliton-like properties, which usually appear in chaotic system.

The chaotic properties of a system are usually determined by the Lyapunov exponent. Therefore, we have studied the Lyapunov exponent of modes (1,m), (2,m), (3,m). The maximal Lyapunov exponent can be numerically obtained as[11, 12]

$$L=\mathop {\lim }\limits_{{\varepsilon \to 0}} \mathop {\lim }\limits_{{z \to \infty }} \frac{1}{z}\ln \frac{{d({\psi _1},{\psi _2};z)}}{{d({\psi _1},{\psi _2};0)}}$$

8

$$d({\psi _1},{\psi _2};z)={\left[ {\int_{{ - \infty }}^{\infty } {\int_{{ - \infty }}^{\infty } {{{\left| {{\psi _1}(x,y,z) - {\psi _2}(x,y,z)} \right|}^2}dxdy} } } \right]^{1/2}}$$

9

Where d(*ψ*1, *ψ*2,; z) is the distance between two functions *ψ*1(x, y; z) and *ψ*2(x, y; z) in the Hilbert space (the L2 norm in the Hilbert space), the initial values *ψ*1(x ,y; 0)=ψ(x, y; 0) and *ψ*2(x, y; 0)=ψ(x,y;0) + r(x, y), and *r(x, y)* is a random perturbation function (as small as machine precision allows, e.g., in the order of 10− 8).

It is found that the maximal Lyapunov exponents are all positive for the condition of p ≥ 1, m ≥ 0, and the maximal Lyapunov exponents of mode (1,6), (2,7), (3,8) are shown in Fig. 7. The positive maximal Lyapunov exponents prove the chaotic properties of the beams, and further confirm correctness of our analysis. Therefore, we confirm that the beams of mode (1,6), (2,7), (3,8) propagate as chaoticons in NLC. Furthermore, Fig. 7 shows the maximal Lyapunov exponents in higher modes is lager than in lower modes. Interestingly, some modes may have the same maximal Lyapunov exponents, just like mode (1,6), (2,7).

As is well known, for a partial differential system with respect to time and space, there are two kinds of chaotic states: the temporal chaos with spatial coherence and the spatiotemporal chaos with spatial decoherence[11–13]. Hence, for making clear the chaoticon properties, we have calculated the spatial cross correlation of the beams in mode (0,1), (1,6), (2,7), (3,8). The spatial cross correlation function of two long enough wave-amplitude series at locations (x1, y1) and (x2, y2) is as follows:

$$\mu [({x_1},{y_1}),({x_2},{y_2})]=\mathop {\lim }\limits_{{z \to \infty }} \frac{{\int_{0}^{{{z_0}}} {\psi ({x_1},{y_1},z){\psi ^*}({x_2},{y_2},z)dz} }}{{\sqrt {\int_{0}^{{{z_0}}} {|\psi ({x_1},{y_1},z){|^2}dz\int_{0}^{{{z_0}}} {|\psi ({x_2},{y_2},z){|^2}dz} } } }}$$

10

Here, the point (x1, y1) represents a arbitrary point on the beams, and is marked by a red asterisk. The spatial cross correlation of this point with other points (x2, y2) in space are shown in Fig. 8.

The first row show that |µ| is equals to 1 for the beam of mode (0,1). This means the beam of mode (0,1), which propagates as soliton, is spatial coherence. The second to fourth row show that |µ| decrease for the beams of mode (1,6), (2,6), (3,8), and drop faster along radial direction than along angular direction. These means the beams of mode (1,6), (2,6), (3,8), which propagate as chaoticons, are spatial decoherence.