# Numerical Study of an Indicator Function for Front-Tracking Methods

Authors: Darae Jeong; Seokjun Ham; Junxiang Yang; Youngjin Hwang; Soobin Kwak; Haobo Hua; Xuan Xin; Junseok Kim
Journal: Mathematical Problems in Engineering
 (2022)

---

## Abstract

In this paper, we present a detailed derivation and numerical investigation of an indicator function for front-tracking methods. We use the discrete Dirac delta function to construct an indicator function from a set of Lagrangian points and solve the resulting discrete Poisson equation with the zero Dirichlet boundary condition using an iterative method. We present several computational tests to investigate the effect of parameters such as distance between points, uniformity of the distance, and types of the Dirac delta functions on the indicator function.


---

## Body

## 1. Introduction

Dirac delta function plays an important role in various scientific and industrial problems. One of them is on the numerical simulations of fluid flows including multiphase flows and fluid solid interaction (FSI) problems. The immersed boundary methods (IBM) are on the front burner for the numerical simulations of these problems and one of the key points lay on the methods to accurately separate the moving multiphase or solid boundaries. The phase field method (PHM) [1–4], level set method (LSM) [5–7], and volume of fluid methods (VOF) [8–10], express the different phases using particular scalar functions. These functions are designed to be advected by fluid flows and one phase changes to another one continuously to avoid numerical stability problem. It should be noted that the interphase is diffused and implicitly captured in the three methods. However, the immersed boundary method is different for its explicit boundary expression. Unlike PHM, LSM, and VOF, distinguish of different phase can be denoted directly by scalar phase functions; in IBM, a so-called indicator function needs to be constructed for these phase distinctions indirectly. In this process, if a Dirac delta function is used for the construction of the indicator function, the diffused interaction forms as PHM, LSM, and VOF. From this point of view, these four methods are all belonged to diffused models in simulations of multiphase flows. In IBM, the fluid interface is explicitly represented by several Lagrangian points and the fluid flow is solved on a stationary Eulerian mesh. Furthermore, when solving Richards’ equation [11, 12] for layered soils, the layers of soils can be expressed using an indicator function.Multiphase fluid flows are important in various scientific and industrial problems. There are many mathematical models and numerical methods for the multiphase fluid flows such as the immersed boundary method (IBM) [13–16], level set method [5–7], phase-field method [1–4], volume of the fluid method [9, 10], and lattice Boltzmann method [17–20], to name a few.In the front-tracking method [21] such as IBM, the fluid-fluid interface is represented by a Lagrangian grid and the flow is solved on a stationary Eulerian mesh. The Lagrangian points are represented by a set of marker points and move with the fluid flow defined at the Eulerian grid. To represent viscosity or density difference between different fluids, and indicator function is calculated based on the marker points [21–23]. Therefore, the indicator function has an important role in not only theoretical issues coupling with the hydrodynamics equations [24], but also in practical approaches in recent studies. The density field was determined by the indicator function derived from an irrotational discrete delta vector in the droplet simulations with a large density ratio in [25]. To incorporate the viscoelasticity effect of Oldroyd-B fluid, the smoothed Heaviside indicator function was used for simulating the dynamics of Newtonian vesicle in [26]. The role of the indicator function was represented by the solid fraction in the diffuse-interface IB framework for conjugate-heat-transfer problems proposed by [27]. Moreover, the indicator function was used to identify the spatial distribution of physical quantities to predict the structure and refractive index profile of fiber-optic components in [28]. In [21], the authors presented the indicator function Ix,y based on a continuous version as follows:(1)Ix,y=∫Aδx−x′δy−y′da′,where A is an area. Then, symbolically, we have the following equation:(2)∇Ix,y=∫A∇x−x′δy−y′da′,=−∫A∇′δx−x′δy−y′da,=−∮∂Aδx−x′δy−y′n′da′,=−∫∂Aδx−x′δy−y′n′da′,where n′ is the outward unit normal vector at the domain boundary ∂A.However, it is not straightforward to connect the relationship between the numerical scheme and the continuous version of the indicator function. Therefore, the objective of this paper is to present a detailed derivation and numerical investigation of the indicator function for front-tracking methods.The paper is organized in the following manner. In Section2, we present the detailed derivation of the indicator function. In Section 3, various numerical experiments are performed to investigate the effect of parameters such as distance between points, uniformity of the distance, and types of the Dirac delta functions on the indicator function. In Section 5, conclusions are drawn.


## 2. Derivation of the Indicator Function

Letδx be a smoothed Dirac delta function and satisfy the following equation:(3)∫−∞∞δx−x0dx=1for anyx0∈ℝ.For example, a 4-point delta function [29–32] is given as follows:(4)δx=1hϕxh,where h is the space step size in the discrete equation and given by the following equation:(5)ϕr=3−2r+1+4|r|−4r28,ifr≤1,5−2r−−7+12|r|−4r28,if1≤r≤2,0,otherwiseNow, we consider a series of the shifted delta functionδx−X0 for some −h<X0<h and it is shown in Figure 1 for some X0=0.7 and h=1, here we take h=1 for better visualization. We describe the h effect in Section 3.1.Figure 1
Schematic of delta functionδx and shifted delta function δx−X0.Letxi=ih for i=±1,±2,…. From the definition (4), we have the following equation:(6)∑i=−∞∞hδxi−X0=hδx−2−X0+hδx−1−X0+hδx0−X0+hδx1−X0+hδx2−X0=1.For anyX0∈ℝ, ∑i=−∞∞hδxi−X0=1 holds.Let us consider a pointx=X0 at the interface in the one-dimensional space. We want to construct an indicator function which increases monotonically from zero to one across the point X0 as x changes its position from left to right of the point, i.e.,(7)Ix=∫−∞xδt−X0dt.Figure2 shows the delta function δx (solid line) and the indicator function Ix (dashed line) with h=1.Figure 2
Delta functionδx (solid line) and indicator function Ix (dashed line) with h=1.Taking a second derivative of equation (7) with respect to x, then we have the following equation:(8)Ixxx=δxx−X0.It implies that the indicator function can be found by solving the Poisson equation (8). Let us consider a 4-point delta function with X0=0 and a discretization for equation (8) as follows:(9)Ii+1−2Ii+Ii−1=δi−1−δi+12,where we have used the unit space step size h=1 and Ii=Ixi=Ii, δi=δxi=δi. Let I−4=0 and I4=I3, then we can rewrite equation (9) for i=−3,−2,…,3 as a matrix vector format:(10)−21000001−21000001−21000001−21000001−21000001−21000001−1I−3I−2I−1I0I1I2I3=δ−2−δ−42δ−1−δ−32δ0−δ−22δ1−δ−12δ2−δ02δ3−δ12δ4−δ22.By multiplying both sides by the inverse of the coefficient matrix, we have the following equation:(11)I−3I−2I−1I0I1I2I3=−1−1−1−1−1−1−1−1−2−2−2−2−2−2−1−2−3−3−3−3−3−1−2−3−4−4−4−4−1−2−3−4−5−5−5−1−2−3−4−5−6−6−1−2−3−4−5−6−7δ−2−δ−42δ−1−δ−32δ0−δ−22δ1−δ−12δ2−δ02δ3−δ12δ4−δ22.From equation (11), we obtain the following results:(12)I−3=12δ−4+δ−3−δ3−δ4=0,Ii=12δi+∑j=−2i−1δjfori=−2,−1,0,1,2,I3=δ−2+δ−1+δ0+δ1+δ2=1,where we have used δ−4=δ−3=δ3=δ4=0. Since δx≥0, the indicator function Ix is a monotonically increasing function, i.e., 0=I−4=I−3≤I−2≤…≤I3=I4=1. Figure 3(a) shows discrete delta function δ−4=δ−3=δ−2=δ2=δ3=δ4=0,δ−1=δ1=0.25, δ0=0.5 and discrete indicator function I−4=I−3=I−2=0, I−1=0.125, I0=0.5, I1=0.975, I2=I3=I4=1 obtained from equation (12). Figure 3(b) shows discrete delta and indicator functions with various point positions, −1<X0<1.Figure 3
(a) Delta functionδx (solid line) and indicator function Ix (dashed line) with X0=0. (b) Discrete delta and indicator functions with various point positions, −1<X0<1.
(a)(b)Next, let us consider a hat shaped indicator function as shown in Figure4 and the sum of the corresponding Dirac delta functions is δx+4−δx−4. If we solve the following equation for the indication function with zero Dirichlet boundary condition, then we have the result in Figure 4.(13)Ixxx=δxx+4−δxx−4.Figure 4
Sum of the Dirac delta functionsδx+4−δx−4 and the associated indicator function Ix.Now, we consider this in the two-dimensional space. Let a computational domainΩ=a,b×c,d be partitioned in Cartesian geometry. Let Nx and Ny be the number of cells in the x- and y-directions, respectively. We assume a uniform mesh with mesh spacing h=b−a/Nx=d−c/Ny. Let xij=xi,yj=a+ih,c+jh for i=0,…,Nx and j=0,…,Ny. Let Xl=Xl,Yl for l=1,…,M be a set of M Lagrangian points, see Figure 5.Figure 5
Eulerian pointsxij, Lagrangian points Xl , and normal vector n for Lagrangian points in the computational domain Ω=a,b×c,d.BecauseIx is an indicator function, the following equation holds at the transition layer of the function:(14)∇Ix∇Ix=−nx,where nx is the outward normal vector at x. Next, let us rewrite equation (14) as follows:(15)∇Ix=−nx∇Ix.We replacenx∇Ix in equation(14). by ∑l=1MnXlδx−XlΔsl, then we have the following equation:(16)∇Ix=−∑l=1MnXlδx−XlΔsl,where δx−Xl=δx−Xlδy−Yl and Δsl=0.5Xl+1−Xl+Xl−1−Xl. After taking the divergence operator to equation (16), we have the following Poisson’s equation:(17)ΔIx=∇⋅∇Ix=−∇⋅∑l=1MnXlδx−XlΔsl.We solve equation (17). with appropriate Dirichlet boundary condition. In this study, we use the Gauss–Seidel method.


## 3. Numerical Experiments

Let us consider two-dimensional indicators on the computational domain,Ω=−4,4×−4,4 with h=1. In Figures 6 and 7, we consider uniform distance change of y-direction, Δs for fixed X0=0. In the immersed boundary method, there are many types of delta functions. Most commonly used delta functions are 2-point [31–33], 3-point [31, 32, 34], 4-point [29–32], 4-point cosine [30, 32], and 6-point [29, 31, 35, 36] functions. These functions are defined as follows:(i)
2-point delta function [31–33]:(18)ϕr=1−r,ifr≤1,0,otherwise.(ii)
3-point delta function [31, 32, 34]:(19)ϕr=131+−3r2+1,ifr≤0.5,165−3r−−31−r2+1,if0.5<r≤1.5,0,otherwise.(iii)
4-point cosine delta function [30, 32]:(20)ϕr=141+cosπr2ifr≤2,0otherwise.(iv)
6-point delta function [29, 31, 35, 36]:(21)ϕr=61112−1142r−1156r2+112r3+3336243+1584r−748r2−1560r3+500r4+336r5−112r6ifr≤1,2116+712r−78r2+16r3−32ϕr−11/2if1<r≤2,98−2312r+34r2−112r3+12ϕrif2<r≤3,0otherwise.Figure 6
From top to the bottom row, 2-point, 3-point, 4-point, 4-point cosine, and 6-point delta functions: (a) 2D delta functionδx−X0=δx−X0δy−Y0 with X0=0 and Y0=0, (b) 2D delta functions δx−Xl=δx−Xlδy−Yl with Xl=0 and various point positions with one distance, (c) sum of 2D delta functions ∑lδx−Xl, and (d) 2D indicator function Ix.
(a)(b)(c)(d)Figure 7
From top to bottom row, 2D delta functionsδx−Xl=δx−Xlδy−Yl with Xl=0 and Yl=±Δs,±2Δs,…, sum of 2D delta functions ∑lδx−Xl, and 2D indicator function Ix: (a) Δs=h, (b) Δs=2h, and (c) Δs=3h.
(a)(b)(c)Figure6 shows 2-point, 3-point, 4-point, 4-point cosine, and 6-point delta functions from top to bottom row with Δs=1. Figures 6(a)–6(d) are 2D delta function δx−X0=δx−X0δy−Y0 with X0=0 and Y0=0; 2D delta functions δx−Xl=δx−Xlδy−Yl with Xl=0 and various point positions with one distance; sum of 2D delta functions ∑lδx−Xl; and 2D indicator function Ix, respectively. We obtain very similar results of indicator functions Ix for all cases except the 2-point delta function. Figure 7 shows the result of delta functions, the sum of delta functions, and indicator functions according to different Δs=h,2h,3h, respectively. As shown in Figure 7, in cases of Δs=h and Δs=2h, the indicator functions are similar. On the other hand, when Δs=3h, the indicator function seems unstable compared to Figures 7(a) and 7(b). For this reason, it is recommended to use Δs which is at least smaller than 2h for a stable indicator function construction.Let us consider an indicator on the computational domain,Ω=−15,15×−15,15 with h=1 using the 4-point delta function. Figure 8 shows immersed boundary positions as a circle with different uniform distance Δs, 2D delta functions δx−Xl=δx−Xlδy−Yl, l=1,2,…,M, the sum of 2D delta functions ∑lδx−Xl, and 2D indicator function Ix from top to bottom row. We can observe that if Δs is too large, then the result of the indicator function is not good.Figure 8
(a) Immersed boundary positions with different distanceΔs, (b) 2D delta functions of at the boundary positions, (c) sum of 2D delta functions, and (d) 2D indicator function Ix with different uniform distance Δs.
(a)(b)(c)(d)The first and second columns in Figure9 show the construction of the indicator functions with randomly distributed points ranging 0.5h<Δsl<1.5h and 2.5h<Δsl<4h, respectively. We can observe that the indicator function is well reconstructed from randomly distributed points ranging 0.5h<Δsl<1.5h.Figure 9
(a) Randomly distributed boundary positions, (b) 2D delta functions of boundary positions, (c) sum of 2D delta functions, and (d) 2D indicator functionIx with different nonuniform distance Δsl.
(a)(b)(c)(d)


### 3.1. Effects ofh

In this section, we describe the effects of spatial grid sizeh. Figure 10 shows 4-point delta functions and corresponding indicator functions for different spatial step h=1,0.5, and 0.25. We observe that as h is smaller, the delta function becomes sharper and the interface transition layer of the indicator function becomes shorter.Figure 10
Delta and indicator functions with different spatial stepsh=1,0.5, and 0.25.In computational domainΩ=−15,15×−15,15, we consider uniform immersed boundary points for a circle radius of 10 with different h=1,0.5, and 0.25 as shown in Figure 11(a). In Figure 11(b), we solve the indicator function, accordingly the points in Figure 11(a). It has a more stiff interface as h is smaller. Figure 11(c) shows Ik−Ik−12 for the Gauss–Seidel iteration with a tolerance tol=1.0e−6, where k is the number of iterations and I0 is a zero vector. We can observe that the number of iterations is getting larger as h approaches zero.Figure 11
(a)-(b) From left to right column, different spatial stepsh=0.25,0.5, and 1: (a) Lagrangian point Xl and (b) indicator functions. (c) l2-error of Ik and Ik−1 for the number of iterations k using the gauss–seidel method with different spatial steps h=1, 0.5, and 0.25.
(a)(b)(c)


## 3.1. Effects ofh

In this section, we describe the effects of spatial grid sizeh. Figure 10 shows 4-point delta functions and corresponding indicator functions for different spatial step h=1,0.5, and 0.25. We observe that as h is smaller, the delta function becomes sharper and the interface transition layer of the indicator function becomes shorter.Figure 10
Delta and indicator functions with different spatial stepsh=1,0.5, and 0.25.In computational domainΩ=−15,15×−15,15, we consider uniform immersed boundary points for a circle radius of 10 with different h=1,0.5, and 0.25 as shown in Figure 11(a). In Figure 11(b), we solve the indicator function, accordingly the points in Figure 11(a). It has a more stiff interface as h is smaller. Figure 11(c) shows Ik−Ik−12 for the Gauss–Seidel iteration with a tolerance tol=1.0e−6, where k is the number of iterations and I0 is a zero vector. We can observe that the number of iterations is getting larger as h approaches zero.Figure 11
(a)-(b) From left to right column, different spatial stepsh=0.25,0.5, and 1: (a) Lagrangian point Xl and (b) indicator functions. (c) l2-error of Ik and Ik−1 for the number of iterations k using the gauss–seidel method with different spatial steps h=1, 0.5, and 0.25.
(a)(b)(c)


## 4. Fluid Application

In this section, we give two examples for the applications of indicator function in multiphase flows and fluid-structure interaction (FSI) problem. First, we consider the two-phase flows problem as the dynamics of a droplet suspended in the simple shear flow. The droplet is between two parallel plates as shown in Figure12. Let ux,t=ux,t,vx,t be fluid velocity at x=x,y and time t. The Lagrangian variable, denoted as Xs,t, where s is arc length parameter, expresses the immersed boundary. The governing equations on the domain Ω are as follows:(22)∂ux,t∂t+ux,t⋅∇ux,t=−∇px,t+1ReΔux,t+1Wefx,t,(23)∇⋅ux,t=0,fx,t=∫TFXs,tδ2x−Xtds,Fs,t=∂2Xs,t∂s2,∂Xt∂t=Uxt,UXs,t=∫Ωux,tδ2x−Xs,tdx,where Γ is the interface between the two fluids. We define the Reynolds number in equation (22) as Re=ρ2R2γ˙, We=Re⋅Ca, and Ca=ηRγ˙/σ, where ρ is fluid density, η fluid viscosity, R droplet radius, γ˙ the shear rate, and σ the surface tension coefficient. We define a gradient field,(24)∇Ix,t=∫ΓnXs,tδ2x−Xs,tds.Figure 12
Schematic illustration of a liquid droplet.To compute the indicator function, the following equation is solved:(25)ΔIx,t=∇⋅∫ΓnXs,tδ2x−Xs,tds.Then, the variable densityρ and viscosity μ are defined as follows:(26)ρIx,t=ρ2+ρ1−ρ2Ix,t,μIx,t=μ2+μ1−μ2Ix,t,where ρi and μi for i=1,2 are density and viscosity of fluid i, respectively.For simplicity, the viscosity and density of the droplet is set equal to that of the matrix phase. Our focus is on the generations of indicator functions from the unit normal of the interface in the drop deformation process. We refer to [37] for the details of numerical algorithms. We consider the simulations of drop deformation on Ω=0,2×0,2 with a space mesh size 128×128. A droplet with a radius 0.5 is put into the center of the domain with initially three sets of different Lagrangian points for comparisons. Delta function (5), Reynolds number Re=1 , and capillary number Ca=0.3 are used in simulations. We compared three different Δs settings to investigate the effect of Δs on the formation of indicator function at the end time t=2.5. The initial Δsini values are set to 1.0h, 5.0h and 10.0h with the number of boundary points 201, 38, and 20, separately, see the first column of Figure 13. While the drop deforms, the interface is remeshed by adding and deleting points with the condition Δsmax/Δsmin≤2.0 for different Δsini leading to the boundary points being more evenly distributed. As the time evolves, the indicator function could be well resolved for Δsini=1.0h as shown in the second row of Figure 13. However, increasing Δs will lead the indicator function a slightly serrated border as shown in the second and sixth rows of Figure 13 (case Δsini=5.0h and 10.0h), which is unacceptable to separate two immiscible fluids. Moreover, we could see the drop interface shrink in the normal direction quickly by setting the Δs larger while drop deforms in the shear flow, as shown in the odd row of Figure 13. Table 1 lists the area losses as the drop evolves at the time t=2.5. It shows that the larger Δs, the more the drop area loss.Figure 13
From top to bottom, the first, third, and fifth rows are the boundary interfaces with different initialΔsini=1.0h, 5.0h and 10.0h (attached with normal vectors). The even rows are corresponding to theirs’ indicator functions. Column (a)–(c), drop evolutions at different times.
(a)(b)(c)Table 1
Drop area loss for the differentΔsini at time t=2.5.
Δsini1.0 h5.0 h10.0 hInitial area0.78520.78220.7725Final area0.78520.68560.2070Area loss (%)0.012.3573.20Next, we take the direct-forcing immersed boundary method as an example application in the FSI problem of the indicator function. The nondimensional governing equations are similar to the multiphase flows and read as follows:(27)∂ux,t∂t+ux,t⋅∇ux,t=−∇px,t+1ReΔux,t+fx,t,∇⋅ux,t=0,fx,t=∫TFXsδ2x−Xtds,∂Xs∂t=UXs,UXs=∫Ωux,tδ2x−Xtdx,where u is the velocity, p is the pressure, u is the body force, and Re is the Reynolds number. From equation (27), we can get the following equation:(28)f=∂u∂t+u⋅∇u+∇p−1ReΔu=∂u∂t−rhs,∂u∂t−rhs=un+1−unΔt−rhs,where n and n+1 are time steps, and rhs=−u⋅∇u+∇p−1/ReΔu. The forcing exerted on the Lagrangian points at the immersed boundary can be written as follows:(29)FXs=Un+1Xt−U∗XsΔt,where Un+1Xs is the desired velocity and U∗Xs is a temporary velocity of the boundary. The algorithm in detail can be found in [38]. The difference of immersed boundary to handle the multiphase flows and FSI lies on the how to evaluate the boundary force.For a moving rigid body immersed on the fluid, the imposed force term, the authors in [39] pointed out it consists of two parts that can be written as the following equations (30) and (31):(30)fx=∮∂ΩFxXsds−∬Ω∂u∂tdV≈∑j=1mFxXjΔVj−∑i=1Nx∑j=1Nyuijn+1−uijnΔtIijΔVij,(31)fy=∮∂ΩFyXsds−∬Ω∂v∂tdV≈∑j=1mFyXjΔVj−∑i=1Nx∑j=1Nyvijn+1−vijnΔtIijΔVij,where ∂Ω is the solid boundary, Ω the solid phase, Δt the time step, Iij the indicator function, uij, vij the velocity in Eulerian gird, Nx, Ny the mesh grid numbers, and Vij the mesh volume in ij grid point.One force in the first terms of right hand sides of equations (30) and (31) is that of a submersed body acting on the external fluid, which can be determined simply by integrating the immersed boundary force densities by equation (29) on the solid boundary that is already evaluated by the direct forcing. Another one that contributes to the unsteady flow inside the solid phase, is the so-called internal or virtual flow force [39] in the second terms of the right hand side of equations (30) and (31). Usually, if the object is fixed in the fluid, we impose the f=fx,fy term to satisfy the no-slip boundary condition on the solid surface leading to the second term equals zero. However, when the solid boundary is moving, we should not neglect the effect of the virtual flow, even though its contribution is small [39]. For a complexly enclosed geometry, the indicator function can be used to calculate to its inside virtual flow by letting its value 1, and outside 0 of the moving solid object.For this moving solid problem in fluid, we consider a 2D NACA0012 airfoil as the representative of a fish body periodically undulating in the incoming flow. The study of this kind of fish kinematics can be found in [40, 41]. Our focus is on the fish shape expressed by the indicator function. As in the first example of drop deformation in the shear flow, we investigate the effect of ds on the indicator function. The computational domain is 0,4×0,2, with mesh size 256 for y direction. The airfoil is initially represented by 202, 102, and 42 Lagrangian points leading to the point distance Δs is about 0.918h, 2.117h, and 6.134h, where h is the Eulerian grid size. The odd rows of Figure 14 shows the fish undulations with outer normal vectors attached on the boundary points for different Δs settings. Because the boundary deforms not too much, we do not need to remesh it. As shown in Figure 14, by increasing Δs, the indicator functions get zigzag distribution on their boundaries between 0 and 1, which could lead to the calculations of fx and fy in equations (30) and (31) a little oscillations. However, the Δs≈1.0h is good enough to resolve the smooth distribution of indicator function as in the first drop deformation example. It should be noted that the diffused fish interface is not well resolved at the high curvature section of the fishtail for larger Δs as shown in the third case of Figure 14. Letting Δs≈1.0h resolves it quite well.Figure 14
From top to bottom, the first, third, and fifth rows are fish boundaries with different initialΔs=0.918h, 2.117h and 6.134h (attached with normal vectors), corresponding to 202, 102, and 42 boundary points. The even rows are corresponding to theirs’ indicator functions. Column (a)–(d), fish undulates at different times in a period. (a) t = 0. (b) t = 0.15. (c) t = 0.3. (d) t = 0.45.
(a)(b)(c)(d)


## 5. Conclusion

In this article, we presented a detailed derivation and numerical investigation of an indicator function together with its fluid applications in the immersed boundary method. We used the discrete Dirac delta function to construct an indicator function from a set of Lagrangian points. We solved the resulting discrete Poisson equation with the zero Dirichlet boundary condition using an iterative method. We presented several computational tests to study the effect of parameters such as distance between points, uniformity of the distance, and types of the Dirac delta functions on the indicator function. Two applications indicate that the construction of the indicator function are effective, multipurpose, and easy-to-use in separating the changing phases. We note that the indicator function can be extended directly and readily to the three-dimensional space. Based on the above-given numerical investigations, as future works, we believe that investigation of the 3D indicator function has more practical values.



---

*Source: 7381115-2022-07-31.xml*