SARS-CoV-2 rate of spread in and across tissue, groundwater and soil: A meshless algorithm for the fractional diffusion equation.

The epidemiological aspects of the viral dynamic of the SARS-CoV-2 have become increasingly crucial due to major questions and uncertainties around the unaddressed issues of how corpse burial or the disposal of contaminated waste impacts nearby soil and groundwater. Here, a theoretical framework base on a meshless algorithm using the moving least squares (MLS) shape functions is adopted for solving the time-fractional model of the viral diffusion in and across three different environments including water, tissue, and soil. Our computations predict that by considering the α (order of fractional derivative) best fit to experimental data, the virus has a traveling distance of 1 m m in water after 22, regardless of the source of contamination (e.g., from tissue or soil). The outcomes and extrapolations of our study are fundamental for providing valuable benchmarks for future experimentation on this topic and ultimately for the accurate description of viral spread across different environments. In addition to COVID-19 relief efforts, our methodology can be adapted for a wide range of applications such as studying virus ecology and genomic reservoirs in freshwater and marine environments.

Introduction

SARS-CoV-2 virus, responsible for covid-19 pandemic, caught the whole world rapidly and by surprise. Since its initial reported positive cases in 2019, the global effort has been put in three major directions: implementation of effective mitigation and prevention protocols, finding effective vaccinations and treatments and finally predicting the rate and pattern of the virus [1], [2], [3], [4], [5], [6], [7]. The latter has become increasingly more crucial as no one is yet immune to virus. This includes how fast the virus spreads in various environments [8], [9], [10], [11]. The main environmental factors includes but not limited to corpse burial or disposing of contaminated waste in soil and ground water. The latter case is of particular importance as even in a developed country such as USA, nearly 65% of the outbreak of aquatic diseases stems from the entry and transmission of the virus into the groundwater aquifers [12].

The World Health Organization (WHO) guidelines for burying deaths from Ebola include 12 steps to maximize safety and prevent transmission of the virus [13]. Although researchers have designed to improve the condition of current corpse bags to prevent the spread of the Coronavirus to the soil [14], due to the likelihood of oversight, it is still imperative to have information about the pattern, rate, and rate of virus transmission from a contaminated corpse to the water or soil that transmits the virus to locations far beyond the burial site. The same possibility applies to the accessories and disposals from those infected with covid with or without symptoms. As a result, he recent focus has been on developing computational and experimental framework to determine the pattern of virus spread in different media.

While most studies are currently focused on how Covid- is transmitted through the air and airborne particles [15], [16], [17], [18], the epidemiological aspect of the bulk virus spread through other environments has been left relatively under studied. Here, we describe a novel framework based on moving least squares (MLS) method for solving the time fractional equation of viral diffusion. In this methodology, we provide a precise and flexible mathematical model based on fractional space to simultaneously consider different virus diffusion coefficients and spread rate in 4 different environments, including tissuewater, tissuesoil, tissuesolidwater and finally tissuewatersoil. This allows us to investigate how the virus is spread and the contamination patterns across different environments.

The origin of fractional derivatives and integrals (i.e. the fractional calculus) as an extension of the ordinary derivatives and integrals (i.e. the classical calculus) can be attributed to Leibniz in 1965 [19]. After Leibniz, from 1730 to 1976, fractional calculus came to the attention of many mathematicians, including Euler, Lagrange, Laplace, Fourier, Liouville, Riemann, Grunwald, Letnikov and Caputo. The Grunwald–Letnikov, Riemann–Liouville and Caputo definitions of the fractional derivative are the well-known fractional derivatives [19]. Fractional derivatives are able to model dynamical systems more accurately due to their non-local property and a greater degree of freedom versus the classical derivatives. Recently, applications of fractional derivatives have been reported in various problems, including anomalous diffusion [20], diffusion–reaction problems [21], electrochemistry [22], chaotic problems [23], viscoelasticity [24], image processing [25], fluid mechanics [26], teletraffic [27], oscillators problems [28], Gaussian noise [29], vibrators [30], biology [31] and medicinal modeling [32]. Generally, analytical and numerical methods are two classes of methods for solving fractional differential equations. The analytical methods include Laplace transform, Fourier transform, Adomain decomposition method, Mellin transform, Green’s function method, etc. For more details, see [33] and references therein. However, in most cases it is often impossible to analytical solve fractional differential equations. In recent years, scholars have developed various numerical methods for solving fractional differential equations including finite difference/spectral method [34], Petrov–Galerkin method [35] Sinc-Chebyshev collocation method [36], etc.

Fractional diffusion equations are one of the most important categories of fractional partial differential equations which have been widely studied in recent years. These equations have been utilized in modeling viscoelastic materials [21], turbulent flow [34], biology [37], chaotic dynamics of classical conservative systems [38] and other problems [39]. We recall that analytical solution of such problems is extremely difficult and in most cases is impossible. This caused that in recent years, several numerical methods have been developed for solving fractional diffusion equations. In [40], a multigrid method developed for spatial fractional diffusion equations with variable coefficients. Iterated fractional Tikhonov regularization method used in [41] for the spherically symmetric backward time fractional diffusion equation. In [33], neural networks based on Legendre polynomials have been developed to solve space and time fractional diffusion equations. Non-standard finite difference and Chebyshev collocation methods have been used in [42] for solving time fractional diffusion equation. In [43], Bayrak et al. have proposed a Chebyshev collocation scheme to solve time fractional diffusion equation. The authors of [44] have applied a numerical method for solving time fractional nonlinear reaction–diffusion equations. An analytic algorithm has been utilized in [45] for finding approximation solution of nonlinear time fractional reaction–diffusion equation. In [46], Cheng et al. used a novel linearized compact alternating direction implicit scheme for Riesz space fractional nonlinear reaction–diffusion equations.

Meshless (or meshfree) approaches are the most important numerical methods for finding the solution of high dimensional (in most cases with complex geometries) differential equations [47]. During last years, meshless approaches provided using the shape functions of moving least squares (MLS) have been extensively applied for divers problems, such as 2D elliptic interface problems [48], fractional telegraph problem [49], fractional version of advection–diffusion problem [50], fractional form of reaction–diffusion equation [51] and integral equations systems [52].

In this study, a mesh free algorithm regarding the shape functions of MLS is developed for finding the solution of the time fractional diffusion equation of the coronavirus with different diffusion coefficients in tissue, soil and water environments. The presented algorithm contains these steps: Applying the finite differences technique (accomplished with weight) for approximating the fractional derivative and consequently making a recurrence algorithm. Next, approximating the problem solution, as well as its partial derivatives using the MLS functions, and inserting them into the main equation. Eventually, extracting an algebraic system of equations which its solution should be found at each time level.

The rest of this work is as follows: We review the MLS approximation in Section 2. Sections 3, 4 explain the proposed meshless method and materials, respectively. The obtained results and related discussion are provided in Section 5. Eventually, in Section 6, the conclusion of this work is provided.

The MLS approximation

Using the shape functions of MLS [51], we able to represent any real function as follows: where and are respectively the vector of basis functions and coefficients. The following cases of are often utilized for two dimensional problems [53]: The vector in relation (2.1) is as follows: The vector is evaluated by minimizing the relation where and is the number of the nodes in the neighborhood of , which satisfy the condition . Despite different weight functions, in this study, we apply the Gaussian weight functions [54] where and is the radius of the influence domain of the node . By considering relation (2.3) and putting , we get where , (the weight moment matrix) is given by and the matrix is as follow: Regarding relation (2.5), we conclude that has the following structure: Relations (2.1), (2.8) result in where the functions () are called the shape functions, and can be computed as So, the shape functions vector can be defined as follows: The second order partial derivatives of the vector can be expressed as follows [51], [53]: where

The MLS meshless method

Different versions of the diffusion equation have been adopted for modeling the transmission of viruses, cells, biomolecules [1], [4], modeling of tumor growth [16], formation and diffusion of amyloid-beta plaques in Alzheimer’s disease [17], and transmission of viruses in unsaturated sand columns [18]. In this study, in order to have a more compatibility of real viral behavior of virus, the time fractional diffusion equation (expressed in relation (3.1)) is considered and solved by using the MLS meshless method. Since in the case of disposal of the infected corpse or personal protective equipment (PPE) in soil or water, it is assumed that there are no nutritional resource of the host (or living cell) available for reproducing the virus, the time fractional diffusion equation is considered without the source term as follows: where accompanied by the initial and boundary conditions in which is the unknown solution and is a given constant. Moreover, , , and are given functions. Herein, we have used the Caputo fractional derivative, which is given in [19] as follows:

Approximation of the fractional derivative

Regarding relation (3.5) and the finite difference scheme, we can discretize the fractional derivative expressed in relation (3.1) as follows (more details are given in [55]): where , for , and

Establishing a recursive algorithm

Regarding the finite difference technique (accomplished with weight), we can write the following relation for Eq. (3.1): So, a recursive relation is obtained by putting relation (3.6) (without error term) into relation (3.7) as follows:

The MLS approximation

Utilizing the shape functions of MLS, the function expressed in relation (3.8) can be represented as follows: where are expansion coefficients which should be computed. To find the values of the coefficients , we should compute the function in the points as follows: Relation (3.10) can be rewritten by so that and are -column vectors, and is an matrix. This matrix is expressed in the form of with Regarding relations (2.12), (3.9), we approximate and as follows: So, regarding the above relation, we get By defining the column vectors , and the -order square matrices we rewritten relation (3.15) as follow: Using the internal points, the matrices and are defined as follows: We derive the following relation by putting relations (3.11), (3.17) into relation (3.8) and employing the boundary condition (3.4): with Here, , , and are -column vectors as follows: Relation (3.19) ables us to compute for . On the other hand, using the initial condition provided in (3.3), we obtain So, we can find as follows:

In this study, we have utilized the following parameters for doing the numerical simulations: Moreover, all simulations have been performed using MATLAB R2019b software.

Materials

The governing equation is solved for exploring the diffusion of the virus of four cases with different combinations of tissue, soil, and water domains which are introduced in Fig. 1 schematically. In all cases under the study, it is assumed that only the central region (tissue) is infected by the virus with a normalized concentration of 1 and other regions are initially virus free. In the current analysis, the Dirichlet boundary conditions are used for the problem. However, due to the physics of the problem, the outer domain of all models are large enough to eliminate the effect of the boundary condition on the diffusion of the virus through the inner region of the mediums. As the diffusion coefficient of the virus are different in the tissue, soil, and water, we had to obtained appropriate values for analysis the virus propagation in these different mediums.

Fig. 1

Four understudy cases for exploring the virus diffusion pattern in different mediums.

Due to the lack of the diffusion coefficient of SARS-CoV-2 virus in water, soil, and tissue, one possibility is using the Stokes–Einstein equation to estimate the required coefficients based on the available experimental data of other viruses. This approach gives a good approximation for spherical shape viruses because the Stokes–Einstein equation says the diffusion coefficient of special particles, , is reversely proportional to the virus diameter, : where is the Boltzmann constant, is the temperature in Kelvin, and is the viscosity of the medium. In the case of soil, the tortuosity coefficient should be considered to modify the diffusion coefficient of the virus (please see Anders et al. [18] for more details).

So, using (4.1) and considering the available experimental data for Influenza, Brome Mosaic virus, and Escherichia virus MS2 and their size and molecular weight analogy, the diffusion coefficients of SARS-CoV-2 virus in water, soil, and tissue, can be estimated accurately. It should be noted that although the high rate of Covid-19 infection and the specific conditions that have arisen may be related to the possible specific structure of the virus’ spike proteins, high mutation rates of the virus or its membrane fusion mechanism in other cells. However, knowledge of the spread of the virus from one environment to another, based on available virus diffusion data, provides useful information for authorities to develop more accurate health protocols.

Among member viruses of the Coronaviridae family, the Coronaviruses are often spherical with a diameter of 80–220 nm [56]. Although some researchers have considered the same diffusion coefficient for all classes of viruses, considering the specific values for viruses based on their shape and size leads to increases in the accuracy of the analysis results. The specific diffusion coefficient of Coronavirus in water, soil, and tissue which are estimated with help of Eq. (4.1) and the available experimental data for Influenza [1], Brome Mosaic Virus [57], Escherichia virus MS2, and Bacterial Virus PRD1 [18]. These values are calculated for environment temperature (298 K) and presented in Table 1.

Table 1

The estimated diffusion coefficients of the SARS-CoV-2 virus in different mediums based on the Stokes–Einstein equation and available experimental data for other viruses.

Medium		Diffusion coefficient (cm^2/s)		Diffusion coefficient (cm^2/s)
		SARS-CoV-2		for other viruses
		(radius: 40–100nm)
Tissue		2.27×10^−8		3.18×10^−8 for Influenza virus (radius: 40−60nm) [1].
Water		3.10×10^−8		15.5×10^−8 for Brome Mosaic Virus (radius: 14nm) [1].
Soil		3.80×10^−8		7.78×10^−8 for PRD1 Virus (radius: 33nm) [18].
				19.5×10^−8 for MS2 Virus (radius: 13.5nm) [18].

Results and discussion

The governing Eqs. (3.1)–(3.4) are solved for exploring four cases of virus spread of the virus in and across the tissue, soil, and water (Fig. 1). Among these four cases, Case I and Case III are more important in terms of available experimental data and comprehensiveness of the studied media, respectively. Results of cases I and III are presented in the main text and the results of the second and fourth cases are given in the supporting information (SI) file.

Fig. 2, Fig. 3 show the normalized virus concentration as a function of distance from the initial infected point and time from tissuewater (case I) and from tissuesoilwater (case III), respectively. Similar graphs for tissuesoil (case II) and tissuewatersoil (case IV) have been illustrated in Figs. and in the supporting information (SI) file, respectively. As indicated in the former section, all the tissues are considered as the center of infection with a normalized concentration of 1. As can be seen from Fig. 2, Fig. 3, as the virus spreads to a larger area over time, the concentration of the virus decreases in infected regions. The changes in the virus concentration at different radial distances can also be monitored over time. In this case, the virus concentration in the remote areas is first zero and then gradually increases over time until it reaches its peak value and after that, the concentration of the virus in the area decreases (see Fig. 4 for case I and case III and also see Fig. for case II and case IV).

Fig. 2

The virus spread as a function of time and distance for tissuewater (case I). The virus spread has been estimated by solving the time-fractional diffusion equation. The computations show over time, the concentration of the virus decreases, but wider areas become infected.

Fig. 3

The virus spread as a function of time and distance for tissuesoilwater (case III). The virus spread has been estimated by solving the time-fractional diffusion equation. The computations show over time, the concentration of the virus decreases, but wider areas become infected.

Fig. 4

The propagating pattern of the normalized virus concentration in different radial distance from the focal infection medium (tissue) during the time for A: case I (tissue water) and B: case III (tissuesoilwater).

The percentage of infected areas over time is also shown in Fig. 5 for all the cases. According to this graph, the rate of environmental pollution increases nonlinearly over time. Therefore, it is not possible to consider a constant velocity for the spread of the virus and predict a linear relationship across different areas.

Fig. 5

The percentage of fully infected area during the time for all cases.

Given the values of diffusion coefficients for the tissue, soil, and water medium (Table 1), it is clear that the rate of virus progression in the soil environment is higher than the other one. For example, in Case 4, as the virus reaches the soil area, the spread rate increases relative to the same area in other cases, and this increase can be easily seen from the slope at the tail-end of the curve in Fig. 5.

It should be noted that the virus diffusion coefficients for the aquatic environment are considered without considering the dynamics or mass transform of the water molecules. Therefore, the coefficient obtained for porous soil is reported to be higher than static water mass (see Table 1). Hence the data of Fig. 5 somehow confirms the result of solving problem for the values of the given coefficients.

To investigate fractional version of diffusion equation and to be comparable with the available experimental data, Eq. (3.1) has been solved with various fractional derivative order of . The results of cases I and III are shown in Fig. 6.

Fig. 6

The percentage of fully infected area in A: case I (tissue water) and B: case III (tissuesoilwater) during the time for different fractional order of time derivative.

As can be seen from Fig. 6, as the derivative order of the diffusion equation decreases, the equation describes the virus behavior with a higher propagation rate. So that the time required to contaminate 90% of the area for is about times the for . It should be noted that the change in slope observed in the graphs at different time intervals is related to the change of media from tissue to water, tissue to soil, or from soil to water, which occurs as a result of differences in virus diffusion coefficients in these environments. Which of the derivatives gives us more accurate results requires further investigation and comparison with the experimental data, which will be discussed in the discussion.

The result of the virus concentration at different radial distances over time for Case I and III are shown in Fig. 7, Fig. 8, Fig. 9. The virus diffusion rate increases due to the change in the order of the fractional derivative (Fig. 7, Fig. 8, Fig. 9). In these figures, which are plotted for different radial distances from the central infected site, it can be seen that in the early times of propagation ( distance from the infection site), the diffusion rate of all derivatives is close to each other. However, at longer distances, the difference in equation behavior becomes more apparent for the order of the fractional derivatives. For example, the required time to reach the peak of contamination to a radius of for is approximately 2.2, 3.0, and 6.5 times to that of , and , respectively. Therefore, the rate of virus propagation strongly depends on the order of the derivative of the governing equation.

Fig. 7

The virus concentration in A: case I (tissue water) and B: case III (tissuesoilwater) at the radial distance of (at the tissuewater and tissuesoil interface respectively).

Fig. 8

The virus concentration in A: case I (tissue water) and B: case III (tissuesoilwater) at the radial distance of .

Fig. 9

The virus concentration in A: case I (tissuewater) and B: case III (tissuesoilwater) at the radial distance of .

It is important to note that although the rate of virus spread depends on the derivative order, the results indicate that the maximum concentration in this region is the same for all the fractional derivative orders. In other words, the maximum concentration of the virus at different distances will not be a function of the fractional derivative order.

The fact that the description of the actual diffusion behavior by Eq. (3.1) is determined by choosing which can be determined by comparing the data of this analysis with the experimental results of virus diffusion.

In Fig. 10, we compared our computational predictions with available values of virus spread for different alpha values. We compared them when the virus reaches its maximum concentration at distances of , , and for case I (tissuewater) [58]. For ease of display, the distance traveled by the virus (horizontal axis) is illustrated logarithmically. As can be seen from the comparison of our computational results with the experimental results (Fig. 10), it is obvious that the virus diffusion behavior in the transition from tissue to water at the intervals of diffusion is completely dependent on the order of the fractional derivative. In such a way that in initial propagation times (at distances ) choosing of lower orders and in times of more than one hour (distance ) assigning higher order of the fractional derivative (close to ) will produce more compatible results with the experimental data.

Fig. 10

The comparison of diffusion distance of the SARS-COV-2 virus with the experiment data [58] for different time fractional derivative order. Since the time required for traveling distances of 1 and are in the range of s, for better display and comparison of the results, the initial part of the distance traveled is given as an inset in the main diagram.

It is worth noting that in the experimental study, the dynamics and mass transfer of water molecules have certainly been effective in increasing or decreasing the spread of the virus, which, as we have emphasized, is not included in the governing Eq. (3.1). However, these results show that fractional equations in the early stages of virus propagation have far higher predictive accuracy than the integer-order derivative () equation.

The data obtained in Fig. 10 show that for distances of 10 to , the order of the fractional derivative that is most consistent with the experimental data is . Fig. 11 shows the results of the relative concentration of the virus in radius of for all cases with . One extrapolation of our data can be used to estimate the diffusion time of the virus from a contaminated body until it reaches one meter from the center of the contamination. According to our model, considering the optimal alpha of 0.9 for the diffusion equation, the estimated diffusion time of the virus from a contaminated tissue buried in the soil at a distance of will be about 925 days. While this time will be more than 1300 days if one uses the usual form of diffusion equation (). Therefore, there should be less concern about corps burial at such distances as the chances that the virus does not biologically survive due to various reasons such as temperature, humidity, and UV light.

Fig. 11

The virus concentration in radial distance of for the optimum fractional order of time derivative () for all cases.

Conclusion

The time-fractional equation of diffusion of the virus was solved to investigate the viral dynamics of the SARS-CoV-2 virus. We adopted a moving least squares meshless method to investigate virus diffusion in and across several mediums including tissue, water, and soil. The time fractional diffusion equation was solved for different values of fractal derivative orders and based on the Stokes–Einstein equation and available experimental data. Our results indicated that in spite of the virus concentration in the infected area which is independent of the time-fractional orders (), the velocity and propagation pattern of the virus is completely different for various fractal derivative orders in all combinations of tissue/water and soil mediums. A comparison of the obtained mathematical result with the available experimental data reveals that adapting a fractal derivative orders form of diffusion equation leads to having a more precise and realistic description of the viral behavior. As a result, by considering the fractal derivative order of 0.9 () the virus has a traveling distance of in water medium after 22 h, which is consistent with the experimental values. By extrapolation of the time-diffusion curve of the virus, one can easily estimate the contamination time of the environment of the focal infection. As an example, the virus needs approximately 925 days for traveling 1 m far away from a buried corpse in the soil. Based on these results, by considering the isolation protocols and removing the surrounding potential nutrients for surviving the virus, the territory of virus diffusion can be limited to a few centimeters. Regardless of the exact estimate of the diffusion coefficient, the focus of the present study is to provide a methodological framework to solve the fractional problem, and comparing them with the results obtained using the conventional solution of the diffusion equation. We showed the effect of the fractional order derivative on the actual propagation behavior is considerable. Therefore, we believe that considering a time-varying alpha (variable-order fractional derivatives) seems to be the best solution for accurately predicting the behavior of known viruses over large time and dimensional scales. These observations reshape our understanding of the dynamics of COVID-19 spread in and across various environments and ultimately, our work should be used to take precautions to avoid environmental cross-contamination. In addition to pandemic relief efforts, our methodology can be adopted for a wide range of applications, such as studying virus ecology and genomic reservoirs in freshwater and marine environments.

CRediT authorship contribution statement

O. Bavi: Conceptualization, Methodology, Validation, Writing – original draft, Visualization. M. Hosseininia: Methodology, Formal analysis, Software, Validation. M.H. Heydari: Supervision, Project administration, Writing, Review and editing. N. Bavi: Conceptualization, Validation, Review and editing.

Declaration of Competing Interest

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.