Discrete Normal Vector Field Approximation via Time Scale Calculus

The theory of time scales calculus have long been a subject to many researchers from different disciplines. Beside the unification and the extension aspects of the theory, it emerge as a powerful tool for mimetic discretization process. In this study, we present a framework to find normal vector fields of discrete point sets in R3 by using symmetric differential on time scales. A surface parameterized by the tensor product of two time scales can be analogously expressed as the vertex set of non-regular rectangular grids. If the time scales are dense, then the discrete grid structure vanishes. If the time scales are isolated, then the further geometric analysis can be executed by using symmetric dynamic differential. Moreover, we present an algorithmic procedure to determine the symmetric dynamic differential structure on the neighborhood of points in surfaces. Our results indicate that the method we present has good approximation to unit normal vector fields of parameterized surfaces rather than the Delaunay triangulation for some points.


Introduction
The mimetic discretization of differential operators is a process that maintains the fundamental properties of continuous differential operators. The main goal in this process is to ensure that the protected properties are maximized and, if not, to give up most of the properties. Geometric partial differential equations are very important to find discrete analogues of differential geometric operators such as mean curvature, Gaussian curvature and Laplace Beltrami, which are defined by surface normal. The de facto methods, such as finite differences and finite elements, are directly related to the discretization of the equation system. A disadvantage of this method is that the selected discretization process may have little connection to the underlying physical problem. However, mimetic methods start with discrete analogue of the continuous theory underlying the problem. Once discrete operators conform to extended physical laws, these mimetic operators can be applied to partial differential equations or integral equation systems. Consequently, it is possible to obtain the discretization of the boundary value problem, which is in accordance with the physical laws on the considered set of definitions. Mimetic methods are a fundamental tool for simulations that do not change their physical properties, whose solution is suddenly variable, irregular grid structures, or long-running simulations, and it is getting more and more important [1]. Mimetic methods are generally used in logically uniform grids [2][3][4], in regular or unstructured grids [5][6][7], in triangular grids [8][9][10], and in polygonal grids [11,12].
In [1], authors particularly point out that mimetic discretization is more effective in practice. The time scale calculus which is an efficient mathematical theory that unifies discrete and continuous calculus emerges as a powerful tool for mimetic discretization process. The first link between effective discretization and generalization of time scale theory and mimetic methods is presented in [13]. The theory of time scales calculus has also provided considerable development in recent years [14][15][16][17][18]. The first work on the geometric interpretation of the theory [19] provided the introduction of the concept of partial dynamic derivatives on time scales [20] and various geometric studies are introduced [21][22][23][24]28].
In this study, we present a framework to express the vertices of non-uniform rectangular grids as the points of surface parameterized on a tensor product of two time scales. The mimetic discretization aspect of the theory of time scales calculus in the geometric sense is first presented in [25] in the terms of symmetric differentiation. The rest of paper is as follows: in Section 2, we present the symmetric differentiation on time scales and the geometric aspect of the theory. We also explain how to compute normal vector fields of surfaces on time scales. In Section 3, we introduce a method to find normal fields of set of discrete points by using symmetric dynamic derivatives. To this end, we determine the set of discrete points in R 3 , namely point clouds, as the subset of a grid like structure emerge from a surface on time scales. In Section 4, we present the numerical results of our method and compare them with the results of well-known method Delaunay triangulation. Finally, in Section 5, we give the conclusions of our study.

Symmetric Calculus on Time Scales
An n-dimensional time scale is defined by the Cartesian product as where for i ∈ {1, . . . , n} the sets T i are time scales. A time scale T i is a non-empty closed subset of reals. We refer readers to [26,27] for details on the theory of single and multivariable dynamic calculus on time scales. Throughout this study we denote σ i and ρ i as the forward and backward jump operators of the time scale T i , respectively. Besides, we set that if T i has a left scattered maximum M and right scattered minimum m, then . . , A n independent of t ∈ T n such that for all t ∈ U δ (t 0 ) and i ∈ {1, . . . , n}, where δ is sufficiently small positive number, U δ (t 0 ) is the δ -neighborhood of t 0 , and α i = α i (t 0 ,t) defined on U δ (t 0 ) such that it is equal to zero for t = t 0 and lim t→t 0 α i = 0 for all i ∈ {1, . . . , n}.

Definition 2.
[25] Let f : T 1 × T 2 → R be a real valued function and (t 0 , s 0 ) ∈ (T 1 ) κ κ × (T 2 ) κ κ . For all ε 1 > 0, there is an open (relative to the topology of T 1 × T 2 ) neighborhood U 1 of (t 0 , s) such that for all (t, s) ∈ U 1 In ∆− and ∇− calculus on time scales, the differentiability of functions defined on time scales comes up with concepts called completely ∆− and completely ∇− differentiability [20,29]. Basically, completely differentiability hypotheses assume the equality of right and left hand side dynamic derivatives. If the point is left dense and right scatter the σ i -completely differentiability or if the point is right dense and left scatter the ρ i -completely differentiability concepts of the functions on T n emerge. However, these hypotheses have strong restriction to define geometric operators on time scales. Besides, geometric operators are not well-defined on isolated time scales. For instance, in [24], authors defined the curvature of curves on time scales by using ∆-derivatives. Analogously, the backward curvature can be defined by using ∇-derivatives. However, such definitions of the curvature are not in unified way, hence would not be useful tool to mimetic discretization process. The introduction of the symmetric dynamic calculus on time scales to overcome such geometric drawbacks are discussed in [25] in details.

Definition 4.
[25] Let S be a closed subset of R 3 . S is a surface if for each point P in S , there is a neighborhood A of P and a function ϕ : U → S where U is a closed set in R 2 and an open set in time scale topology satisfying the following conditions: i.e., ϕ is 3-regular.
The function ϕ : U → S is called a surface patch. S is called a 3-smooth surface if, for all points P in S there exists a surface patch such that P ∈ ϕ(U).
Proposition 1. Let U ⊂ T 2 and f be a 3-differentiable function. Then, the set determines a 3-smooth surface. Proof. Let {t, s} be the Euclidean coordinate system of T 2 . Since the coordinate functions and f are 3-differentiable, ϕ : U → S is also 3-differentiable. Jacobian matrix of ϕ respect to symmetric differentiation is rankJ(ϕ) = 2 for all (t, s) ∈ U, hence ϕ is 3-regular. Besides it is trivial that ϕ mapping is a homeomorphism.
In this study, we determine the metric tensor of a surface S on time scales by the partial 3-derivatives of ϕ.
Hence, mimetically, the surface normal vector can be compute as In the calculation of the normal vector fields of surfaces, two cases emerge on finiteness of time scales. If the time scales are isolated and infinite, then the normal vector fields can be computed directly by symmetric differentiation. However, if the time scales are isolated and finite, then we need to compute the vector fields on boundaries by using ∆ or ∇ differentiation. To see this, we give the following theorems: Theorem 2. Let U andŨ be nonempty closed subsets of R 2 and ϕ : U → S be a 3-regular surface patch. If φ :Ũ → U is diffeomorphism, then the functioñ By the chain rule, we obtain The coefficient at the right side of the equation (1) is equal to determinant of the jacobian matrix This completes the proof.
As a result of Theorem 2, we may remark that if T 1 × T 2 is an infinite isolated time scale, then it is possible to find a coordinate change which makes the point in ( T 1 ) κ κ × ( T 2 ) κ κ . Hence, normal vector field of a surface on time scales can be computed directly by symmetric differentiation. If T 1 × T 2 is a finite time scale, the natural boundaries of the surface arise at the end points of time scales. Let t m 0 and s m 0 be the minimum points and t M 0 and s M 0 be the maximum points of T 1 and T 2 , respectively. By Theorem 3, ∂ ϕ(t, s) Figure 1, we present a surface parameterized on a finite isolated time scale T 1 × T 2 , where T 1 = {t 1 ,t 2 ,t 3 ,t 4 } and T 2 = {s 1 , s 2 , s 3 , s 4 }. This figure also serves a good example of a non-regular rectangular grids which has vertices on a surface parameterized on an isolated time scale. In this figure, the arrows represent the normal vectors at respected points. T 1 and T 2 both have the usual time scale topology which is respect to partial order of the indices. Hence, the only non-boundary point is ϕ(t 2 , s 2 ) on the surface. Thus, the normal vectors can be computed as

Method
3D point clouds represent a discrete representation of the surfaces, namely discrete manifolds, present in the real world. If this sample is obtained from range sensors such as 3D scanners, noise is expected in the sample. Because of this type of noise, clear information about surface orientation and curvature can be lost. The normal vector estimate tries to reconstruct this information by creating a set of vectors perpendicular to the tangential plane of each surface [30][31][32]. In particular, the resulting normal vector is a reasonable procedure for integrating each data item into a feature space corresponding to a point of the point cloud. In this section, we present a method by using the concept of symmetric differential for predicting normal vectors on a discrete representation of surfaces existing in the real world.
Our approach is based on considering the points of point clouds as the points on a surface patch parameterized on two isolated time scales. However, the direct implementation of this approach has two major drawbacks. First drawback emerges when the points on discrete manifold are not aligned on grid-like structure. In this case it becomes impossible to determine the forward or backward jumps of the coordinate variable, therefore the dynamic differentiation even in symmetric sense is not well-defined. To cope with this drawback, we give the following proposition: Proposition 4. Let M ⊂ R 3 be the set of non-uniform points and S be a surface with an atlas on T 1 × T 2 . Then, M can be given a symmetric differentiable structure in such way that the inclusion i : M → S is an embedding. Proof. Let F : R 2 → M and S be a surface with an atlas on T 1 × T 2 with the coordinate chart (x, y). Let us first fix x ∈ M. Now choose the homeomorphism ξ x : F −1 x → S with x ∈ F −1 ∩ (T 1 × T 2 ). Since the symmetric covariant derivative 3 y ξ x is injective, we can choose a bijection γ : {1, 2, 3} → {1, 2, 3} such that the rows γ(1),γ(2),γ(3) of 3 y ξ x are linearly independent. Define π : S → S by π(x, y, z) → (x , y , z ). Then, the 3 y (π • ξ x ) is an isomorphism. Hence there are open sets in the time scale topology A ⊆ F −1 (M) and B ⊆ π(S ) and the map η x : B → A is the inverse of π • ξ x . Define Therefore, A is a symmetric differentiable atlas making M into a symmetric differentiable surface. The inclusion i : M → S is a homeomorphism, and for any patch we have Id • i • η −1 x , which has injective symmetric derivative.
The second drawback is subject to determining the parameterizations of the surface patches. When the unorganized points are the subject, there are several ways to interpolate them [35][36][37]. In this study, we use the quadratic surface fitting.
A quadratic surface passing through origin is It can be seen that the surface with the Equation 2 requires eight other given points (x i , y i , z i ), i = 1, . . . , 8, and is expressed by the system The accuracy of the fitting method is directly depending conditioning of the matrix in Equation 3. Let us assume that the patch in Equation 2 is expressed in vector form as From Proposition 1, we may conclude that the eight-point neighborhood of the R(0, 0) is a 3-smooth surface. Now, let us consider this neighborhood as in Figure 2. Hence, the partial symmetric dynamic derivatives are Therefore, accuracy of our method is order two with O(max{σ 2 1 (x), ρ 2 1 (x), σ 2 2 (y), ρ 2 2 (y)}). If T 1 × T 2 ≡ hZ × hZ or T 1 × T 2 ≡ R × hZ or T 1 × T 2 ≡ hZ × R, then the accuracy becomes of order two with O(h 2 ).
Given an non-uniform set of points M in R 3 , it is possible to determine the 3-smooth surface S containing M by parametrization of F : R 2 → M. We can say from Proposition 4 that the symmetrical differential structure of S is also in M by the embedding i : M → S . With the parametrization of M, the geometric structure of S , that is, the vertices of the non-regular rectangular grid structure, which contains the M, can be determined as S be the cardinally smallest 3-smooth surface due to reduce computational costs.
A quadratic surface patch to model an entire surface is not convenient but it is good for modeling neighborhood around a point. If an entire surface needs to be modelled, then we may implement our method to cubic splines [33,34]. However this is totally subject to another study. To determine neighborhood of the data point in M for fitting quadratic surface, we use the closeness relation of 3D data points. That is, we obtain an appropriate graph G = (V, E), where V is the set of points and E is the set of edges with In our study, the procedure in Table 1 is applied to achieve the bundle S = S . The computational complexity of the procedure is directly dependent on the size of M and the size of neighborhood of i-th point. Since quadratic fitting needs only eight points, this complexity can be reduced to O(|M|).
end for end for end for Output: Bundle S = S Table 1 The procedure to obtain the bundle of 3-smooth surfaces on M.

Results
In computer graphics and engineering analyses, the triangular meshing of surfaces plays a key role. A surface meshing can be achieved by mapping meshes in parameter space onto surfaces where the meshes can be triangular or rectangular grids [38,39]. A good looking mesh in the parameter space may have a problematic image on surface under these mappings since the transformation of geometry from the parameter space to the surface may be twisted along some directions. In many engineering applications including Finite Element analysis the triangular mesh is the most popular choice to get over this problem.
The de facto method to obtain triangular surface meshes is the Delaunay triangulation. Traditionally, the Delaunay triangulation is a cell complex that subdivides the convex hull of the discrete points in R 3 in which every circum-circle of a triangle is an empty circle [40]. The algorithmic complexity of the Delaunay triangulation of M ⊂ R 3 is O(|M| 2 ). However, if the points in M are well distributed on a smooth surface, then the Delaunay triangulation has the reduced complexity O(|M| log |M|) [41].
Let M ⊂ R 3 be the set of discrete points and the Delaunay triangulation of D(M) composed of the points P. For a local triangular mesh of a surface interpolating M with K triangles, the normal vector n at P ∈ M is equal to where n i and w i are the normal vector and weight of the i-th triangle P i PP i+1 , respectively. The weight w i of the i-th triangle can be computed as w i = | PP i × PP i+1 | 2 [31]. In this section, to make the comparative analysis, we give the computational results of the normal fields of our method and Delaunay triangulation. To see the comparison, we use the two parametric surfaces and measure the vector errors with infinity norm. The first surface has the parametrization as ϕ 1 (t, s) = (ts 2 ,t 2 s, 1 − √ t 2 + s 2 ) for t, s ∈ [−1, 1], and the second surface has the parametrization as ϕ 2 (t, s) = ((1, 16 s ) cos s(1 + cost), −1.16 s sin s(1 + cost), −2(1, 16 s )(1 + sint)) for t ∈ [0, 2π] and s ∈ [− 15,6]. The smooth surfaces and points of M sampled on them are presented in Figure 3. The Delaunay triangulations of M are also presented in Figure 4.  The graphs derived from the Delaunay triangulation are called the Delaunay graph of M. The procedure presented in Table 1 initially starts with a graph representation of M. Hence, to see the effectiveness of our method, we use the Delaunay graph of M as G = (V, E), where M is chosen as well distributed to reduce time complexity. Besides, the quadratic surface fit requires eight other given points. Therefore, if the 1-neighborhood N(i) does not composed of nine points, then we extend the neighborhood to k-neighborhood in which composed of optimally many points.
In order to measure the error in normal vectors, we need to measure the size or norm of the vectors || N − N D || ∞ and || N − N T S || ∞ at P ∈ M, where N is the unit normal vector of the parameterized surfaces, N D is the unit normal vector obtained by the Delaunay triangulation, and N T S is the unit normal vector obtained by our method. The measured errors are presented in Figures 5-6.
The error measurements show us that the present method based on symmetric dynamic derivatives yields better approximations for certain points of both ϕ 1 (t, s) and ϕ 2 (t, s). The numbers of sampled points on ϕ 1 (t, s) and ϕ 2 (t, s) are 2571 and 2295, respectively. The lesser error are obtained at 962 many points of ϕ 1 (t, s) and 113 many points of ϕ 2 (t, s) by using our method. Besides, the min{|| N − N T S || ∞ } = 9.24888 × 10  3.38129 × 10 −1 for ϕ 2 (t, s). The points where our method yields better approximation than the N D to unit normal vectors is presented in Figure 7. Fig. 7 The points where the unit normal is approximated better for φ 1 (t, s) and φ 2 (t, s).

Conclusions
The mimetic discretization of continuous operators yields us efficient way to model continuous theory underlying the physical problem. For the mimetic discretization process the theory of time scale calculus emerge as an efficient mathematical theory. If the nature of the problem involves the forward or backward discretization, then the delta or nabla dynamic differentiation may dominate the modelling and solutions. However, in geometric point of view, the symmetric dynamic differentiation mimics the discrete counterpart of the geometric modelling, since it ensures that the protected geometric properties such as curvatures are maximized.
In this study, we consider a surface on time scales as the vertex set of non-regular rectangular grids. Then, we determine the normal vector fields of surfaces parameterized by the tensor product of two time scales by using symmetric dynamic differentiation. If we have such a closed subset set of R 2 as the parameter domain of a surface, then the most basic mapping ϕ(t, s) = (t, s, f (t, s)) yields us a 3-smooth surface. Besides, if a surface patch has a finite geometry, then the symmetric differential acts as forward or backward dynamic derivatives on boundaries.
In real world applications, a set of discrete points in R 3 does not involve a regular geometric structure. In this paper, we also present an algorithmic procedure to approximate normal fields of such sets which are sampled on a smooth surface. Our procedure first start with the geometric closeness relation of discrete points. This relation is expressed as a finite graph. Then, by using the neighborhood of points in the graph, we fit a quadratic surface to approximate the parameter of underlying smooth surface. Afterwards, we assign a non-regular rectangular grid like structure to this parametrization and consider the embedding as a surface on time scales. We also show that such immersion of discrete points to surface on time scales is an embedding. Therefore, these discrete points also have the symmetric dynamic differential structure of a surface on time scales. Subsequently, the normal vector fields are determined by symmetric dynamic derivatives.