key: cord-0761151-77p95mu3 authors: Daneh-Dezfuli, Alireza; Zarei, Mohammad Reza; Jalalvand, Mehdi; Bahoosh, Reza title: Simulation of time-fractional oxygen diffusion in cornea coated by contact-lens date: 2022-03-09 journal: Mech Time Depend Mater DOI: 10.1007/s11043-022-09545-0 sha: e09e9ed285cd63e68d677b3b2e08837bee055e0c doc_id: 761151 cord_uid: 77p95mu3 In this paper, the time-fractional oxygen diffusion has been simulated in a one-dimensional (1D) corneal-contact lens (CL) system. Different CLs have been employed as Balafilcon, thin- and thick-Polymacon. It is assumed that homogeneous and isotropic porous mediums of cornea and CL is saturated with compressible oxygen. The computations of the time-fractional derivations are done based on the Caputo method. The obtained results show that the fractional derivative order (FDO) severely affects pressure distribution in cornea and CL. Consequently, the magnitudes of post-lens-tear-film (PoLTF) pressure change due to diverse FDOs. Particularly, maximum changes have been observed in the results gained from the CLs with thicknesses more than 100 μm. The agreement of the results obtained from the time-fractional modeling with the experimental data compared to the standard diffusion modeling has been improved by more than 36%. Finally, it has been demonstrated that high-thickness CLs can cause exist anomalous diffusion process in cornea tissue. It is generally accepted that the standard topology of fluid particles diffused in a space has a linear relationship to time (Zayernouri 2015) . However, this point of mathematical view is sometimes unable to find suitable responses for experimental tasks. For the sake of this incapacitation, it is preferred on some issues to take into account the anomalous diffusion process. In this regard, the displacements of fluid particles can have a non-linear relationship to time. From the point of mathematical view, first-order derivations of time replace with fractional-time derivations. The application of the time-fractional calculus is nowadays increasing in computational analysis of diseases behavior such as COVID-19 transmission patterns (Awais et al. 2020) , epidemic childhood diseases (Baleanu et al. 2020) , damages in skin tissue caused by laser irradiation (Hobiny et al. 2020) , growth of tumor cells in the presence of chemotherapy (Kumar et al. 2021) and Alzheimer medical phenomena (Mohammad and Trounev 2021) . Fractional modeling of these medical problems allows finding appropriate curves fitted by the unpredicted distribution of experimental data. In this paper, the fractional modeling concept is extended to achieve a better understanding of oxygen diffusion in cornea tissue coated by CLs. Numerous studies have been accomplished related to the estimation of corneal oxygen diffusion. Fatt et al. (1974) studied pressure distribution for a corneal-CL system assuming steady-state conditions. The CLs included Softlens, Bionite, and Silicon rubber lenses. Harvitt and Bonanno (1999) considered the effect of acidosis on pressure distribution in cornea worn by the polymethyl methacrylate, Oxyflow f30, Fluorex 700, and Oxyflow 151 lenses. Brennan (2005) improved Harvitt and Bonanno's solution to eliminate negative pressure detections. In these three reviewed studies, the simplified governing equations have been solved employing the solution algorithms of ordinary differential equations. Alvord et al. (2007) simulated oxygen diffusion in a two-dimensional cornea geometry with average measures using the finite-element method. It is reported that the actual shape of any individual cornea can significantly have a deviation from the used dimensions. Although there is a layered pressure distribution same as 1D geometry in most cornea volumes, the peripheral part of the cornea contains spatial variation, especially in the closed-eye situation. Chhabra et al. (2009b) determined oxygen diffusion numerically for the cornea-CL system based on the second Fick's law. To this end, the reaction term is defined by the Monod kinetics model. Different values of oxygen solubility coefficient and oxygen permeability have been computed for a certain cornea when coated by the Balafilcon and thin-Polymacon CLs. Larrea and Büchler (2009) offered a new numerical simulation for transient diffusion of oxygen in light of the Michaelis-Menten-type relationship. Chhabra et al. (2009b) have developed a more complete model accounting for the role of reactive species such as glucose, bicarbonate, lactate, hydrogen, and carbon dioxide along with oxygen to define an index as oxygen deficiency factor in oxygenation of the cornea-CL system. Takatori and Radke (2012) prepared a quasi-two-dimensional model for respiration of the cornea to gauge corneal hypoxia. Del Castillo et al. (2015) repeated the methodology suggested by Chhabra et al., assuming constant cornea characteristics. It has been claimed that the maximum oxygen consumption rate (OCR) is unfixed when the oxygen pressure diminishes at the interface CL. Compañ et al. (2017) and Moreno et al. (2021) have been analyzed the application of the generalized Monod kinetics model to describe the human corneal OCR during various soft CLs wear. It has been observed that maximum OCR increases to oxygen pressure of 105 mmHg and then decreases. These variations can be related to limitations in all the OCR models. Therefore, a generalization of them should be performed to acquire a better description of the behavior of the cornea in humans. The common denominator of all the recent studies reviewed is that the standard form of Fick's law has been used for almost all simulations performed. Besides, the corneal-CL system has been taken into account as a 1D geometry. The cornea is thin compared to its area, so the effects of lateral diffusion are neglected (Harvitt and Bonanno 1999) . From the points of material and skeleton views, a certain cornea tissue must have constant characteristics such as oxygen solubility coefficient, oxygen diffusion coefficient, and oxygen permeability. The employment of different CLs can solely affect the value of the maximum OCR. This idea is utilized in the numerical simulations by Del Castillo et al. (2015) . However, the final outputs differed from the laboratory results by Bonanno et al. (2002) . To this end, anomalous diffusion simulations are presented again, considering the constant characteristics of corneal oxygen. In other words, the main focus of the current paper is on the development of the oxygen diffusion analysis in cornea based on the time-fractional derivative. Most of the fractional-time derivatives are modeled using the Riemann-Liouville (Fei et al. 2021; Wu et al. 2021; Zhou et al. 2021 ) and Caputo (Caputo 1967; Almeida 2017; Wu et al. 2020; Liu et al. 2021a,b) approaches. Unfortunately, the Riemann-Liouville approach leads to initial conditions containing the limit values of the Riemann-Liouville fractional derivatives at the lower terminal (Podlubny 1998). The main advantage of Caputo's approach is that the initial conditions for fractional differential equations with the Caputo derivatives take on the same form as for integer-order differential equations (Podlubny 1998). Therefore, the fractional-time derivative is modeled here employing the Caputo-based derivative. It is deserved to say that the current study explores the effect of various FDOs on the pressure distribution and PoLTF pressure of CL-equipped cornea as a novel experience. The obtained results can be used in the geometrical design of CLs. The presented simulations include a 2-layer system made of the cornea, c , and CL, l , mediums. The CLs are in the types of Balafilcon, thin-and thick-Polymacon materials. Mathematical modeling of the problem, as well as the results and discussions of time-fractional simulations, is described in the following sections with more details. Figure 1 depicts a corneal-CL system. The cornea geometry with a thickness of L c , which its interior surface contains endothelial cells. Epithelial cells of the cornea are coated by an isothermal post lens tear film. The CL geometry has a thickness of L l whose exterior surface is exposed to atmospheric air pressure at the open-eye situation and palpebral conjunctive at the closed-eye situation. Since the CL is porous, oxygen particles can penetrate and reach the cornea. Practical assumptions given below are made to simplify the present study: 1. The human cornea consists of three different porous mediums as Bowman's layer, corneal Storma, and Descemet's membrane (Mobaraki et al. 2019) , which transport metabolic species (oxygen, lactate ion, glucose, hydrogen ion, bicarbonate ion, and carbon dioxide) (Chhabra et al. 2009b) . Porous mediums of cornea and CLs are homogeneous, isotropic, and saturated with oxygen. 2. Characteristics of cornea and CLs are invariable. 3. Cornea and CL mediums have chemical equilibrium. 4. Cornea thickness is much less than its surface area, and thus, lateral diffusions have been neglected. 5. The spatial variation of oxygen pressure is 1D. 6. Thicknesses of endothelium, epithelium, and stroma are integrated into cornea thickness. 7. Thicknesses of pre-and post-lens tear films (∼3 µm) are ignored. 8. Oxygen is considered a compressible fluid. The standard governing equations are presented here to facilitate the creation of the timefractional equations represented in Sect. 2.2. In an open-eye situation, the governing equation describing anomalous oxygen diffusion in the cornea-CL system in Fig. 1 consists of the 1D and time-fractional Fick's law expressed as (Chhabra et al. 2009a ) where P represents the distribution of oxygen pressure at time t in x-direction. K and (Dk) are oxygen solubility coefficient and oxygen permeability for each medium, respectively. Q denotes the OCR, which can be computed based on the Monod kinetics model as shown below (Chhabra et al. 2009a) . where Q c,max is maximum OCR. K m denotes the Monod dissociation equilibrium constant. Initial conditions for the governing equation can be obtained from a closed-eye situation using the steady-state solution of Eq. (1) as (Chhabra et al. 2009a ) Both time-dependent and steady-state governing equations can be solved under boundary conditions of the Dirichlet type summarized in Table 1 on the basis of the pressures of the anterior chamber, P ac , air, P air , and palpebral conjunctiva, P pc . In addition, the employed characteristics of the cornea medium are represented in Table 2 Time-fractional modeling of the diffusion process in the cornea-CL system can be represented as shown below. where α is the order of fractional derivative, and ∂ α /∂t α is the fractional (non-integer) time derivative. The standard oxygen diffusion equation obtains considering α = 1. Therefore, Eq. (6) is the same as Eq. (1) when α = 1. Anomalous diffusion process can be observed with 0 < α < 1. Using the theory of operators, the time derivative term of Eq. (6) can be rewritten as below: where D denotes fractional derivative operator. For the numerical solution of the problem, the spaces of the cornea and CL have been decomposed with N c and N l points, respectively. The time scale is partitioned with n levels. Discretized form the time-fractional term presented in Eq. (7) can be computed on the basis of the Caputo methodology as (Caputo and Cametti 2008) : where i is spatial point indicator. denotes the Gamma function. The left-hand side of Eq. (8) can be computed numerically in the following manner (Toubaei et al. 2019) : where j is the time level indicator. t refers to time step, and Discretizing the diffusion term in Eq. (6) based on the central finite-difference approximation as well as the utilization of the result of Eq. (8), the implicit discretized time-fractional equation at nth time level can be presented as where x is the spatial step. The non-linear partial differential equations of Eqs. (3) and (11) can be solved by the Newton-Raphson method (Press et al. 1989 ). The simulation results of the time-fractional oxygen diffusion in a cornea coated by distinct CLs are first studied in this section and then discussed compared to standard oxygen transport. The utilized CLs are in the types of Balafilcon, thin-and thick-Polymacon; their characteristics are been given in Table 3 . It deserves to be mentioned that the lens characteristics are taken from the Chhabra et al. (2009b) data. The influence of various FDOs, α {0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1.}, on oxygen pressure distributions is shown in Fig. 2 after 60 , 40, and 70 s respectively inside the cornea separately coated by Balafilcon, thin-and thick-Polymacon CLs. As can be seen, the FDO value has a significant influence on the pressure distribution within the cornea-CL medium. Moreover, the usual linear distribution of pressure in the standard oxygen diffusion process, α=1 within the CL medium, has been converted to a non-linear behavior in the anomalous oxygen diffusion by adjusting different FDOs less than 1. The effect of the same FDOs on Fig. 2 The influence of the FDO on the oxygen pressure distribution in the cornea coated by different CLs Fig. 3 Effect of the FDO on the PoLTF pressure via various CLs wearing the post-lens tear film pressure is illustrated in Fig. 3 for the cornea while wearing different CLs. Again, it can be seen that different values of the FDOs strongly influence the process of the post-lens tear film pressure changes. To select the most appropriate amount of the FDO, the value of the L 2 relative error norm (L 2 REN) is calculated in terms of experimental, exp, and calculated, cal values based on the following form: where NOE refers to the number of those experimental data reported by Bonanno et al. (2002) . Figure 4 shows the value of the L 2 REN versus the FDOs for each of the CLs. It can be seen that the L 2 REN value for the Balafilcon and thick-Polymacon CLs has been obtained less than 10% at α = 0.9. The common denominator of both CLs belongs to their thicknesses which are more than 100 µm. These CLs can therefore reveal that the selection of high thicknesses for CLs can suggest the anomalous oxygen diffusion in the cornea medium. It can also be seen in the same figure that the L 2 REN of the thin-Polymacon CL in α = 1 is less than 10%. The thickness of this CL is 60 µm. Therefore, it can be claimed that a CL with thickness less than 100 µm can present standard oxygen diffusion in the corneal medium. Based on the appropriate values of the FODs calculated in Fig. 4 , the oxygen transfer simulation has been executed again, considering time-fractional modeling in the cornea covered by each of the CLs. To this end, steady-state simulations have been first accomplished using the iterative Newton-Raphson method to compute initial conditions. The convergence diagram of these simulations is shown in Fig. 5 . To achieve accurate initial conditions, the stop criterion is set to 10 −15 . This figure demonstrates that the initial pressure distributions have been obtained via fast iterative processes, finally up to 7 iterations per CL. The results of the PoLTF pressure distribution are plotted in Fig. 6 . The experimental Bonanno et al. (2002) data and the standard diffusion results of Del Castillo et al. (2015) are also given in this figure. It is obvious for the Balafilcon and thick-Polymacon CLs that the results of the time-fractional modeling rather than the standard diffusion results agree better with the experimental data. The improvement of the PoLTF pressure calculation obtained from the time-fractional modeling against the published numerical data of standard diffusion is presented in Table 4 based on the value of L 2 REN. It can be deduced that the timefractional method can improve the modeling of PoLTF pressure distribution by more than 36% for the CLs with thicknesses greater than 100 µm. In this paper, the process of anomalous oxygen diffusion has been simulated in a certain cornea coated by different CLs. The Balafilcon, thin-and thick-Polymacon CLs have been selected to investigate the results. It is assumed that compressible oxygen penetrates in 1D homogeneous and isotropic porous mediums of cornea and CL without any chemical reaction. The anomalous diffusion process has been modeled by the fractional-time method with various FDOs. To integrate non-integer time derivation, the Caputo approach based on an implicit discretization scheme is selected. The Newton-Raphson algorithm has been employed to solve non-linear discretized equations. The FDO has a significant effect on pressure distribution in a corneal-CL system. Linear pressure variation in the CL medium can be converted to a non-linear one. Moreover, minimum pressure magnitude inside the cornea changes in the fractional-time modeling. Similarly, PoLTF pressure manipulates based on the time-fractional modeling. Particularly, all significant changes are evident for the CLs with thicknesses greater than 100 µm. The agreement of the results obtained from time-fractional modeling with the experimental data demonstrates that high-thickness CLs lead to form anomalous diffusion process in cornea tissue. Funding This work was supported by Shahid Chamran University of Ahvaz under Grant No. SCU.EM1400. 26639. The authors confirm that the data supporting the findings of this study are available within the paper. A Caputo fractional derivative of a function with respect to another function Corneal oxygen distribution with contact lens wear Modeling and simulation of the novel coronavirus in Caputo derivative On modelling of epidemic childhood diseases with the Caputo-Fabrizio derivative by using the Laplace Adomian decomposition method Estimation of human corneal oxygen consumption by noninvasive measurement of tear oxygen tension while wearing hydrogel lenses Beyond flux: total corneal oxygen consumption as an index of corneal oxygenation during contact lens wear Linear models of dissipation whose Q is almost frequency independent-II Diffusion with memory in two cases of biological interest Diffusion and Monod kinetics to determine in vivo human corneal oxygen-consumption rate during soft contact-lens wear Modeling corneal metabolism and oxygen transport during contact lens wear Analysis of the application of the generalized Monod kinetics model to describe the human corneal oxygenconsumption rate during soft contact lens wear Diffusion and Monod kinetics model to determine in vivo human corneal oxygen-consumption rate during soft contact lens wear Oxygen tension distributions in the cornea: a re-examination A triaxial creep model for salt rocks based on variable-order fractional derivative Re-evaluation of the oxygen diffusion model for predicting minimum contact lens Dk/t values needed to avoid corneal anoxia The effect of fractional time derivative of bioheat model in skin tissue induced to laser irradiation Analysis of tumor cells in the absence and presence of chemotherapeutic treatment: the case of Caputo-Fabrizio time fractional derivative A transient diffusion model of the cornea for the assessment of oxygen diffusivity and consumption A Caputo fractional damage creep model and its experimental validation Nonlinear damage creep model based on fractional theory for rock materials Corneal repair and regeneration: current concepts and future directions Explicit tight frames for simulating a new system of fractional nonlinear partial differential equation model of Alzheimer disease A refined model on flow and oxygen consumption in the human cornea depending on the oxygen tension at the interface cornea/post lens tear film during contact lens wear Numerical Recipes A quasi-2-dimensional model for respiration of the cornea with soft contact lens wear Boundary functions determination in an inverse time fractional heat conduction problem New fractional variable-order creep model with short memory A new model of shear creep and its experimental verification Spectral and spectral element methods for fractional PDEs A viscoelastic-viscoplastic mechanical model of timedependent materials based on variable-order fractional derivative Publisher's Note Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations