In recent years, a lot of researches have been devoted to the nonconventional applications of optical computing  - . In fact, different binary patterns have been used in many optical devices to implement intriguing phase profiles and for generating the complex amplitudes of the well-known catastrophe functions. The application of spatial light modulators has permitted the experimental realization of some artificial catastrophes in optics, e.g., the fundamental catastrophe beams such as the Airy (Fold) and Pearcey caustic modes   . It is well established in the catastrophe theory      that the catastrophe functions can be described by seven elementary structures which are classified according to their codimensions K. The caustic beams emerging as diffraction catastrophes are particular solutions of the paraxial Helmholtz equation, and they can be expressed in the typical integral representation     .
where s represents the state variable, Cj are the n-2 control parameters and is the catastrophe phase polynomial of n-th order that depends on the state variable s and the control parameters Cj with the codimension K = n-2.
The three members of the hierarchy catastrophe structures that we will examine in this paper are given in Table 1, where we have listed the considered diffraction catastrophes with the corresponding phase polyomials and codimensions K. For the sake of comparison we have reported, the Berry’s  and Stewart’s  unfolding phase polyomials.
By using the change of the variable , one obtains the relationship
between the control parameters C and in the Refs.  and . For K = 1, the phenomenon is called Fold diffraction caustic, we have
and the structure function is proportional to the well-known Airy function
Table 1. The elementary diffraction catastrophes of odd codimensions : with Berry and Stewart unfoldings.
where the Airy function is defined as ,
For K = 3, the phenomenon is called Swallowtail diffraction catastrophe     , so we have
and the caustic function reads
The case K = 5 corresponds to the so-called Wigwam diffraction catastrophe  and we obtain the following identities
It is worth noting that, apart from the case K = 1, the caustic functions above haven’t been expressed, to the best of our knowledge, in terms of well-known mathematical functions. In the remainder of the paper, we will show the relationship between these catastrophe functions and the Airy-Hardy integrals, and propose an optical method by using tailored binary screens, for the creation of the elementary optical catastrophes of odd codimensions (K = 1, 3 and 5).
2. Caustic Functions and Their Connection with the Airy-Hardy Integrals
One can note that for particular values of the control parameters, the phase functions in Equations (6) and (7) are reduced to odd polynomials, and thus the associated catastrophe functions are proportional to the Airy-Hardy cosine integrals . In fact, if one takes in the case K = 3, or for K = 5, the associated phase polynomials and will be odd functions, then it will be simpler to express the caustic functions (K = 1, 3 or 5) as
These last integrals belong to the Airy-Hardy integrals class and then they can be expressed in the canonical form  
where denotes the Airy-Hardy function and
, for (11b)
where is the hypergeometric function defined by 
where is the Pochhammer symbol.
The first odd-order functions are given by
It is well-known that the Airy-Hardy functions are solutions of the following differential equation  
Hobbs et al.  have expressed these solutions in terms of Bessel functions of fractional orders as
where and are the υ-th order Bessel and the modified Bessel functions of the first kind, respectively.
For , the approximate value of reads
, for n = 1, 2 and 3 (15)
By putting and using the similarity between the Equations (11b)
and (16a), we obtain
The integral expression of Equation (9) can be rewritten
If we take and , the last diffraction integral reduces
The Wigwam catastrophe integral can also be expressed as
By putting , this leads to
We inferred from the above procedure that the caustics of odd codimensions K = 1, 3 and 5 are connected to the Hardy-Airy functions of orders n = 3, 5 and 7, respectively. In the forthcoming section, we will be interested in displaying in optical way these non-conventional functions.
3. Generation of the Elementary Optical Catastrophes with K = 1, 3 and 5
Based on the fact that the complex caustic structures are expressed in terms of polynomial phase functions, we adopt an optical technic similar to that employed by Lohmann et al.   for implementing Airy, Bessel and Laguerre functions. The proposed setup is a Fourier transformer whose schematic diagram is given in Figure 1: A plane wave front illuminates a binary screen which is transparent only along the curve on the XY plane. Thus, the transmittance function in the binary screen can be expressed as
where is the curve profile and is the Dirac function.
The amplitude distribution of the Fraunhofer diffraction pattern in the back focal plane of the lens (L) is obtained by taking the Fourier transform of Equation (19). This leads to
Figure 1. Diagram of the optical setup: at the input plane XY a binary screen with the transmittance illuminated by a plane wave. Detector for output pattern at the back focal plane of the lens L.
Substituting from Equation (19) into Equation (20a) and using the integral property of the Dirac function yields
It is worth noting that Equation (20b) is a typical integral representation for many special functions employed is lasers physics, e.g., the Airy, Bessel, Laguerre functions  . In particular, if F(x) is odd function, the Equation (20b) can be rewritten as
In the following, we propose the transmission profiles that may create oscillating integral patterns encountered in many fields, e.g, in physics, chemistry and biology.
3.1. The Fold Catastrophe Pattern
First, let us consider a mask transmission profile of the form , where is an arbitrary constant (Figure 2(a)).
The Fourier transform of this input amplitude can be written as
with . The irradiance distribution of the Fraunhofer pattern in the
output plane is depicted in Figure 2(b).
By placing a narrow slit in the output plane along the υ-axis, at just behind the slit, the output diffracted amplitude reads
It follows from this last equation that by displacing the slit to arbitrary values
Figure 2. (a) Input binary pattern with , (b) associated output Fraunhofer pattern, and (c) square modulus of Airy-Hardy Ch3(x) with = 2.
of , one can produce the Hardy-Airy function along the υ-axis. Figure 2(c) shows the irradiance of which is consistent with the well-known Airy function, with an arbitrary .
3.2. The Swallowtail Catastrophe Pattern
Now, by considering the curve function of the form (Figure 3(a))
where p and q are arbitrary constants, the Fourier transform of the associated input amplitude reads
Making the change of the variable , with , this leads to
Putting and , and after some algebraic operations, Equation (25) is reduced to
Reciprocally, one can writes
The irradiance distribution of the Fraunhofer pattern of the input amplitude is shown in Figure 3(b). By placing a narrow slit along the υ-axis at , just behind the slit, the resulted amplitude reads
Hence, one can visualize the profile of the Hary-Airy function along the υ-axis. Figure 3(c) shows the irradiance profile proportional to with = 2.
3.3. The Wigwam Catastrophe Pattern
By considering the mask transmission function Fw
where p, q and r are arbitrary constants, its Fourier transform is given by
Now, making the change of variable , with , and by taking
Figure 3. (a) Input pattern with , (b) associated output Fraunhofer pattern and (c) square modulus of Airy-Hardy Ch5(x) with = 2.
, Equation (30) can be expressed as
Following the same procedure as described in Sections (3.1) and (3.2) and using the binary mask of Equation (29), we have displayed the irradiance distribution of the Fraunhofer pattern associated with in Figure 4(c). It is
Figure 4. (a) Input pattern with , (b) associted output Fraunhofer pattern and (c) square modulus of Airy-Hardy Ch7(x) with = 2.
worth noting that the patterns of Figure 2(c), Figure 3(c) and Figure 4(c) can be controlled by adjusting the parameters σ, p, q, r; e.g., a negative value of σ will give the symmerical profile of Figure 2(c).
The parameter σ can be regarded as a scaling length for controlling the Fraunhofer pattern structure of Figure 2(b), and in the same way, different mapping caustics for Swallowtail and Wigwam optical catastrophes (Figure 3(c) and Figure 4(c)) can be performed, globally, by varying the scaling parameters p, q and r of the binary amplitude transmittance in the input plane This procedure may permit application opportunities in micromachining on tailored curves and light guiding paths.
In summary, based on the fact that the integral diffraction catastrophes can be reduced in the case of odd polynomials phase functions into Airy-Hardy cosine integrals, we obtained the closed-form expressions for the related caustic functions. We showed that by the use of an optical Fourier Transformer device with tailored binary screens one can obtain the displaying of caustic beams such as the Fold, the Swallowtail and the Wigwam catastrophes. Furthermore, it is shown that the mapping caustics for Swallowtail and Wigwam caustics can be controlled by the scaling parameters p, q and r of the binary mask transmittance. The result of this work can be useful in shaping lasers with catastrophe amplitude function and may have application opportunities in micromachining and light guiding paths. This study can be extended to higher oscillatory integral functions.
 Lohmann, A.W., O-Castaneda, J. and Ibarra, J.G. (1994) Airy Function and Laguerre Polynomials: Optical Display and Processing. Optics Communications, 109, 361-367. https://doi.org/10.1016/0030-4018(94)90482-0
 Siviliglou, G.A., Broky, J., Dogariu, A. and Chritodoulides, D.N. (2007) Observation of Accelerating Airy Beams. Physical Review Letters, 99, Article ID: 213901.
 Ez-Zariy, L., Nebdi, H., Boustimi, M. and Belafhal, A. (2014) Transformation of a Two-Dimensional Finite Energy Airy Beam by an ABCD Optical System with a Rectangular Annular Aperture. Physical and Chemical News, 73, 39-49.
 Ez-Zariy, L., Hennani, S., Nebdi, H. and Belafhal, A. (2014) Propagation Characteristics of Airy-Gaussian Beams Passing through a Misaligned Optical System with Finite Aperture. Optics and Photonics Journal, 4, 325-336.
 Ebrahim, A.A.A., Ez-Zariy, L. and Belafhal, A. (2015) Propagation of Finite Airy-Hermite Hollow Gaussian Beams through a Paraxial ABCD Optical System. International Journal of Advancement in Earth and Environmental Sciences, 3, 11-20.
 Ez-Zariy, L., Hricha, Z. and Belafhal, A. (2015) Novel Finite Airy Beams Generated from Gaussian Array Beams Illuminating an Optical Airy Transform System. Progress in Electromagnetics, 49, 41-50. https://doi.org/10.2528/PIERM16041909
 Ez-Zariy, L., Boufalah, F., Dalil-Essakali, L. and Belafhal, A. (2018) Conversion of the Hyperbolic-Cosine Gaussian Beam to a Novel Finite Airy-Related Beam Using an Optical Airy Transform System. Optik, 171, 501-506.
 Ouahid, L., Dalil-Essakali, L. and Belafhal, A. (2018) Effect of Light Absorption and Temperature on Self-Focusing of Finite Airy-Gaussian Beams in a Plasma with Relativistic and Ponderomotive Regime. Optical and Quantum Electronics, 50, 216.
 Ouahid, L., Dalil-Essakali, L. and Belafhal, A. (2018) Relativistic Self-Focusing of Finite Airy-Gaussian Beams in Collisionless Plasma Using the Wentzel-Kramers- Brillouin Approximation. Optik, 154, 58-66.
 Yaalou, M., Hricha, Z., El Halba, E.M. and Belafhal, A. (2019) Production of Airy-Related Beams by an Airy Transform Optical System. Optical and Quantum Electronics, 51, 308. https://doi.org/10.1007/s11082-019-2023-5
 Yaalou, M., Hricha, Z. and Belafhal, A. (2019) Propagation Properties of Finite Cosh-Airy Beams through an Airy Transform Optical System. Optical and Quantum Electronics, 51, 356. https://doi.org/10.1007/s11082-019-2075-6
 Ring, J.D., Lindberg, J., Mourka, A., Mazilu, M., Dholakia, K. and Denis, M.R. (2012) Autofocusing and Self-Healing of Pearcey Beams. Optics Express, 20, 18955-18966. https://doi.org/10.1364/OE.20.018955
 Boufalah, F., Dalil-Essakali, L., Nebdi, H. and Belafhal, A. (2016) Effect of Turbulent Atmosphere on the On-Axis Average Intensity of Pearcey-Gaussian Beam. Chinese Physics B, 25, Article ID: 064208. https://doi.org/10.1088/1674-1056/25/6/064208
 Hobbs, C.A., Connor, J.N.L and Kirk, N.P. (2007) Theory and Numerical Evaluation of Oddoids and Evenoids: Oscillatory Cuspoid Integrals with Odd and Even Polynomial Phase Functions. Journal of Computational and Applied Mathematics, 207, 192-213. https://doi.org/10.1016/j.cam.2006.10.079
 Zannotti, A., Diebel, F., Boguslawski, M. and Denz, C. (2017) Optical Catastrophes of the Swallowtail and Butterfly Beams. New Journal of Physics, 19, Article ID: 053004. https://doi.org/10.1088/1367-2630/aa6ecd
 Cheng, K., Lu, G. and Zhong, X. (2017) The Poynting Vector and Angular Momentum Density of the Swallowtail-Gauss Beams. Optics Communications, 396, 58-65. https://doi.org/10.1016/j.optcom.2017.03.038
 Cheng, K., Lu, G., Zhou, Y., Yao, N. and Zhong, X. (2018) The Poynting Vector and Angular Momentum Density of the Autofocusing Butterfly-Gauss Beams. Optics and Laser Technology, 105, 23-34. https://doi.org/10.1016/j.optlastec.2018.02.029