1. Introduction
This article will use bone fractal operators as the backdrop, the error function as the link, and combine the symmetry evolution laws of functional fractal trees in the physical fractal space to reveal the fractional-order correlation between different special functions.
In recent years, our research on the self-similar structures of biomaterials and fractional-order mechanics has involved a series of symmetry and symmetry-breaking issues. For instance, when establishing the fractional-order constitutive model of tendon/ligament fibers, we found that the symmetry breaking of fractal cells introduces a time modulation factor in their relaxation response, further extending the relaxation time and successfully revealing the inevitable correlation between the fractional-order viscoelastic emergent behavior of biological fibers and their fractional-dimension self-similar structures [
1]. When establishing the fractional-order impulse characteristics of neurons, we discovered that the symmetry breaking of step cells introduces an exponential modulation term in the dynamic operator, which also provides a functional explanation for why nerve fibers have spiny dendrites [
2]. When studying the infinite-level physical fractal circuit model and fractional-order time-varying dynamic response of arterial blood flow, we found that the symmetry breaking of fractal cells is a necessary condition for the formation of the reflux time-varying response phenomena in real arterial blood flow [
3].
In studying the aforementioned issues, we found that as the symmetry of the functional fractal tree evolves in the physical fractal space, the properties of the fractal operators become more enriched, as shown in
Figure 1. Specifically,
Figure 1a represents the abstract process of fractal cells in nerve fibers, forming a resistance-capacitance network with an infinitely expanding tree-like topology that is structurally perfectly symmetrical [
2].
Figure 1b represents the abstract process of fractal cells in neuronal dendrites. Its assembly of components is consistent with that of nerve fibers, and it maintains a macroscopic binary bifurcation structure [
2]. However, at this point, due to the introduction of an additional physical component, the functional fractal tree undergoes its first symmetry breaking. This phenomenon is caused by the unequal number of physical components at the two ends of the bifurcation. Physiologically, this multi-scale fractal tree model of dendrites can also be used to explain the electrophysiological changes associated with aging and neurodegenerative diseases [
4,
5]. If the topological structure of the functional fractal tree and the number of the physical components at the two ends of the bifurcation are changed simultaneously, a second symmetry breaking will occur, forming a new fractal step topology.
Figure 1c represents the abstract process of fractal cells in compact bone fibers [
6]; this soft/hard alternating multi-level self-assembly pattern is not unique to bone and is widely found in biological materials such as shell [
7] and nacre [
8]. In our study of the operator algebra method for viscoelastic theory, we found that fractal tree and fractal step operators share commonalities: both are quadratic radical operators, or in other words, both are non-rational operators. There are also differences: the fractal tree operator is extremely simple, while the fractal step operator is quite complex. The simplicity of the former comes from the symmetry of the fractal tree topology, while the complexity of the latter arises from the symmetry breaking of the fractal step topology. It is important to note that the methods for breaking symmetry are not unique.
Figure 1b achieves symmetry breaking by changing the number of physical components at the two ends of the functional fractal tree. If the number of physical components is not changed, but only their functional properties are altered, symmetry breaking can also be achieved. Therefore, if the topological structure of the functional fractal tree and the functional properties of the physical components at the two ends of the bifurcation are changed simultaneously, the fractal step topology shown in
Figure 1c can also be formed. We found that as the number and degree of symmetry breaking increase, the properties of the fractal operators become more enriched. For a perfectly symmetrical functional fractal tree (
Figure 1a), it is difficult to characterize the various mechanical properties of complex organisms, such as the high modulus and toughness of bone [
9,
10,
11,
12,
13], with just two types of physical components and the classical spring-dashpot model. The introduction of symmetry breaking increases the designability of fractal operators, enriches their basic functions, and broadens the applicability of the physical fractal space.
In studying the biomechanical and biophysical issues, we abstracted the “physical fractal space” and established the mechanics in the “physical fractal space [
6]”, which is different from the traditional “geometric fractal model [
14,
15,
16]”. Surprisingly, the mechanics in the physical fractal space exhibit common features: non-rational fractional-order fractal operators, mechanics governed by non-elementary functions, non-integer-order mechanics, and non-localized mechanics. These fractional-order fractal operators all possess rich invariance properties. Based on the invariant properties of fractal operators, we confirmed that the special functions from different disciplines and research backgrounds indeed have fractional-order correlations. With the help of fractal operators, the deep dependency relationships between special functions can be revealed. The fractal operators in the physical fractal space provide a new perspective for understanding non-elementary functions and fractional-order calculus [
6].
In studying the convolutional kernel function of bone fractal operators, we confirmed that the error function is its core component [
6]. We found that the fundamental laws in the physical fractal space can all be characterized by fractal operators and the convolution kernel functions of each type of the operators are related to the special functions. Despite the completely different objects of study, an astonishing consistency is displayed. How can we understand this consistency? Is there a correlation between these special functions? This study attempts to provide answers, including the following contents: (a) an analysis of the intrinsic expressions and characteristics in the convolutional kernel function of bone fractal operators; (b) a demonstration of the correlation between the convolutional kernel function of bone fractal operators and special functions under various input signals (external load, external stimulus); and (c) a revelation of the fractional-order correlation between different special functions based on bone fractal operators. The definition and operations of traditional fractional-order calculus often involve certain special functions, such as the Gamma function. In other words, the emergence of special functions has provided the theoretical foundation for the establishment and development of fractional-order calculus [
17,
18,
19]. Now, inspired by the symmetry evolution laws of functional fractal trees in the physical fractal space, with the aid of bone fractal operators, we established a series of correlations between special functions, which represents a novel attempt to understand the commonalities across different disciplinary fields.
2. The Intrinsic Expression in Convolutional Kernel Function of Bone Fractal Operators
The previous article abstracted the physical fractal space from the compact bone fibers [
6]. This section mainly discusses the intrinsic expression in the convolutional kernel function of bone fractal operators. For convenience, we first provide the algebraic expressions of bone fractal operators. When stress
is applied to bone fibers and strain
is generated, then we have
where
is the fractal operator on the physical fractal space of the bone. The differential operator
is defined as follows:
Unlike the conventional approach, our research object is selected as a microfiber. Force and deformation are transmitted along the fiber chain, from which we can abstract a multi-level self-similar component tree. If the number of levels in the multi-level structure tends toward infinity, the self-similar component tree becomes a fractal component tree, as shown in
Figure 1c. The “=” sign means that the intrinsic invariance in the physical fractal space is the equivalence of the fractal cell and fractal component [
6].
Therefore, the bone fractal operator
needs to satisfy the following algebraic equation:
In mechanics, the right-hand side of Equation (3) can be regarded as the stiffness of the fractal component, and the left-hand side can be regarded as the stiffness of the fractal cell. It is further deduced that
Equation (4) is the algebraic equation for the fractal operator
. This is a quadratic equation with a radical solution:
So far, the fractal operator of compact bone is determined, as shown in Equation (5).
If the physical component operators are taken as
where both
and
are elastic elements and
is a viscous element, then the main body of the fractal operator in Equation (5) becomes
where
is the physical constant, equivalent to the characteristic frequency [
6].
Interestingly,
operators appeared in both hemodynamics and bone mechanics [
6]. In hemodynamic research, Peng et al. [
3] considered the operator
as a whole. Unlike Peng et al., in bone mechanics, Jian et al. separated the quadratic radical term
in the operator
and studied it separately [
6]:
The quadratic radical operator
is decomposed into the product of two 1/2-order differential operators
and
. In Courant’s work [
20], the operator
acting on the Heaviside unit step function
generates a convolutional kernel function
:
With the help of the displacement theorem [
21], the expression for the kernel function
can be written as follows:
In Equation (7), is the characteristic function, is the characteristic integral of the convolutional kernel function of bone fractal operators, where is the exponential term brought by the displacement theorem. comes from the convolutional kernel function of the operator .
If the input signal is
and the fractal operator
acts on
, then the output signal can be expressed as
Obviously, the operation of bone fractal operator
inevitably involves the integrals such as
and
:
we call it the weighted integration of the characteristic function
. Obviously, the weighted integration is the core component of the 1/2-order fractal operator
. Moreover, the weighted integration is closely related to the error function. The previous article [
6] has exported
namely,
In Equation (9), integration
can be considered as the 0 moment of the characteristic function
. Therefore, the kernel function can be written as
and
Combining the differential properties of error function,
the final expression of the kernel function can be exported
and
Equation (12) shows that the core component of characteristic integration is the error function . It can also be said that the core component of the 1/2-order fractal operator is the error function . Therefore, the following analysis primarily focuses on the invariance properties of weighted integration and the error function.
3. The Moment of the Characteristic Function and Its Correlation with Error Function
This section considers the typical cases of the general input signal . Common input signals are typically elementary functions, such as the unit slope function , the unit acceleration function , the sine function , the cosine function , the exponential function , and the logarithmic function , among others.
This section attempts to answer the following basic question: when the input signal takes different types of elementary functions, is still related to the error function? If the answer is yes, can this correlation reveal the correlation between the error function and different special functions? It should be emphasized that there are theoretically an infinite variety of elementary functions and their combinations. This study, however, will only discuss the following typical elementary functions.
When
is taken as a polynomial function, special cases can be considered first. The unit function
, corresponding to the 0 moment of the characteristic function, which has been discussed in the previous section. Here, when the input signal is the unit slope function
,
In Equation (13), integration can be regarded as the first moment of the characteristic function
When the input signal takes the unit acceleration function
,
In Equation (14), integration can be regarded as the second moment of the characteristic function
From Equations (13) and (14), the integration of the weighted characteristic function is related to the error function. To export more general expressions, let
take the derivative of the parameter variable
at both sides of Equation (15):
is the
n-th moment of the characteristic function. It can be inferred that
A very simple recursive relationship and general formula can be written using Equation (18):
According to Equations (9) and (15), the 0 moment
can also be written as
It can be inferred that the n-th moment is related to the error function. Therefore, we give the following proposition:
Proposition 1. When the input signal
takes any n-degree polynomial function,
is related to the error function.
5. The Integral Property of Error Function
This section will confirm the following proposition.
Proposition 3. When the input signal
takes different elementary functions, the integration operation results with the error function still possess correlation with the error function.
In the case of
taking the Heaviside unit function
, using the partial integration method and combining it with the differential properties of the error function (Equation (11)), we can derive
In Equation (25), is the 0 moment of error function. The 0 moment of the error function is still related to the error function.
The input signal takes the unit slope function—that is,
. By using the partial integration method and combining with the differential properties of the error function (Equation (11)), we can derive
In Equation (26), is the first moment of error function. The first moment of the error function is still related to the error function.
The input signal takes the unit acceleration function—that is,
. By using the partial integration method and combining the differential properties of the error function (Equation (11)), we can derive
In Equation (27), is the second moment of the error function. The second moment of the error function is still related to the error function.
The 0, 1st, and 2nd moments of the error function are all related to the error function. We speculate that this phenomenon should be universal. There should be the following proposition: the n-th moment of the error function is still related to the error function.
Now we confirm the above proposition holds. We directly analyze the general situation. By using the partial integration method, we can obtain the
n-th moment of the error function:
The last term in Equation (28) involves the (
n + 1)-order moment of the negative exponential function. The definition of the lower incomplete Gamma function is as follows [
22]:
the (
n + 1)-order moment of the negative exponential function can be written as
The general formula for the
n-th moment of error function can be derived from Equation (29):
At this point, we confirm the following general proposition: the integration of the polynomial weighted error function is not only related to the error function but also to the Gamma function.
Let us take a different perspective and look at the general formula for the
n-th moment. We focus on the recursive properties of the
n-th moment. Equation (30) gives the moment of (
n − 1)-order as follows:
namely,
Based on Equation (32), recalculating the
n-th moment, we can derive the recursive formula:
The recursive Formula (33) shows that the
n-th moment of the error function can be expressed as an integral of the (
n − 1)-order moment. Similarly, the moment of (
n − 1)-order can be expressed as an integral of the (
n − 2)-order moment. By repeating the above operations, we can conclude that
By combining Equations (30) and (34), we can conclude that
The last equation of Equation (35) is summarized as follows:
namely,
Equation (36) indicates that the high-order integration of the error function is not only related to the error function but also to the high-order integration of the Gamma function. We can write the error function and Gamma function in Equation (36) on both sides of the equal sign:
In Equation (37), the left end only has an error function, and the right end only has a Gamma function. Two unrelated functions can be given identity relationships through complex integration operations.
We can also understand Equation (37) from the perspective of integral equations. Consider the following integral equation:
It is obvious that this integral equation has the following special solutions:
That is to say, the error function is the solution of the nonlinear integral Equation (38).
It is important to note that the first term
of the last equal sign in Equation (35) is closely related to the weighted feature integration formula
. According to the differential properties of the error function (Equation (11)), the weighted feature integration can be expressed as
The convolution in the last term of Equation (40) can be calculated as follows (the proof can be found in
Appendix A):
Therefore, the weighted feature integration can be transformed into
Equation (42) indicates that when taking , is still related to the n-th integration of the error function.
For ease of comparison, Equation (42) can be rewritten as
Comparing Equations (36) and (43),
Simplifying the above equation yields
Equation (44) indicates that when
taking
,
is not only related to error function but also to the
n-th integration of Gamma function. Furthermore, according to Equation (44), we have
By combining Equations (45) and (46) and eliminating the error function, we can obtain the recursive relationship:
When the input signal
takes the sine
and cosine
functions, respectively, we can derive that
When the input signal
takes the exponential function
, we can derive that
We can introduce the 0 moment:
Taking the derivative of the parameter variable
at both ends yields, we can obtain the first moment:
Taking the derivative of the parameter variable
at both ends yields again, we can obtain the second moment:
The
n-th moment is defined as
Therefore, there are the recursive formulas and general formulas:
By exporting the 0-order moment , any n-th moment can be easily exported.
The above equations can also deduce the recursive formula for the
n-th moment of the error function. Let
:
Therefore, there exists an extremely simple recursive relationship:
For the first item, the following equation holds:
where
is the upper incomplete Gamma function
. Equation (59) can also be written as
and
satisfy the following relationship [
23,
24]:
Equation (60) is a special case of Equation (30). Taking , Equation (30) can degenerate into Equation (60).
When takes any n-degree polynomial function, the weighted integration of the error function is not only related to the error function but also to the Gamma function.
8. Fractional-Order Correlation between Error Function and Struve Function
When studying hemodynamics, Zhou et al. derived the kernel functions of fractal operators
and
as follows [
26]:
Zhou et al. also confirmed that the kernel function at the right end of Equations (65) and (66) can be integrated into an explicit form, namely [
26]:
wherein
and
are 0-order and 1st-order Struve functions, respectively.
and
are 0th-order and 1st-order modified Struve functions, respectively [
27]. The above equations show that there is a profound intrinsic correlation between the Bessel function and the Struve function, as well as between the modified Bessel function and the modified Struve function, and Equations (67) and (68) have similar forms, reflecting good symmetry.
In
Section 7, we established the fractional-order correlation between the error function and Bessel function using the operator
. If we make the following transformation to Equation (66), we have
It is important to note the relationship
and substitute Equation (12) into Equation (69):
According to Equations (66) and (70), we can conclude that
Equation (71) indicates that there is a fractional-order correlation between the error function and modified Struve function. Therefore, we give the following proposition.
Proposition 5. The correlation between the modified Struve function and the error function can be established through a 1/2-order fractional differential operator .
The above equation can be regarded as a 1/2-order fractional differential equation:
Equation (72) must possess the following special solution:
By combining Equations (62) and (70), we can conclude that
By combining Equations (66) and (74), we can conclude that
Equation (75) shows that there is a fractional-order correlation between the Gamma function and modified Struve function. Therefore, we give the following proposition.
Proposition 6. The correlation between the modified Struve function and the Gamma function can be established through a 1/2-order fractional differential operator .
Equation (75) can also be regarded as a 1/2-order fractional differential equation:
and the above equation must possess the following special solution:
It is important to note that the Struve function has the following recursive properties:
The modified Struve function has similar recursive properties, which will not be further elaborated here. Equations (77) and (78) show that the correlation between the Struve function and Gamma function is diverse. As mentioned above, in Equations (71) and (75) is a 1/2-order operator, which characterizes the fractional-order correlations. In Equation (77), is the first derivative, which characterizes the integer-order correlations. There is no derivative in Equation (78), which characterizes the algebraic correlations.
10. Fractional-Order Correlation between Error Function and Generalized Hypergeometric Function
In
Section 5, we provide the expression for the
n-th moment of the characteristic function when
takes any n-degree polynomial function (Equation (40)). In fact, when
takes
, the
n-th moment of the characteristic function can also be expressed in the following form (the proof process can be found in
Appendix D):
In Equation (80),
is the generalized hypergeometric function and
is the binomial coefficient. In mathematics, hypergeometric functions are defined using hypergeometric series, and many special functions are their special cases or limits [
27]. The
n-th moment at the left end of Equation (80) is related to the fractional-order operator. Therefore, the generalized hypergeometric function at the right end characterizes the fractional-order effect.
We can rewrite Equation (80) as
we can conclude that
Therefore, we give the following proposition.
Proposition 7. There is an integer-order integral correlation between the hypergeometric series, error function, and Gamma function.
This proposition generally holds, but there are also special cases, such as those that follow.
Taking
in Equation (81), we can obtain the following special case:
The second equation of Equation (82) is transformed into
Equation (83) indicates that the error function can be decomposed into the product of an exponential function and a hypergeometric function. This indicates that there is an algebraic identity between the hypergeometric function and error function. It is important to note that hypergeometric functions constitute a vast set. Equation (83) shows that an element in this set, i.e., , can be represented by the error function.
Taking
in Equation (81), we can obtain the following special case:
The second equation of Equation (84) indicates that there is also a simple algebraic transformation relationship between the error function, Gamma function, and hypergeometric function.
By combining Equations (83) and (84), we can conclude that
Equation (85) is organized as follows:
Equation (86) indicates that the Gamma function can be expressed as the product of a power function, an exponential function, and a hypergeometric function.
Similarly, Tremblay and Choi et al. [
28,
29] obtained a series of transformation formulas for hypergeometric functions through a fractional-order operator, providing a new perspective for establishing the correlations between special functions.
Thus, this article uses bone fractal operators as the backdrop, the error function as the link, and combines the symmetry evolution laws of functional fractal trees in the physical fractal space to reveal the fractional-order correlation between different special functions. Moreover, in the study of the calculus properties of special functions, a notable phenomenon of symmetry was also discovered. The network of interrelationships between special functions, with the error function as its core, is shown in
Figure 2.