EUROPEAN JOURNAL OF PURE AND APPLIED MATHEMATICS
2025, Vol. 18, Issue 4, Article Number 6917
ISSN 1307-5543 – ejpam.com
Published by New York Business Global

A Linear Subdivision Scheme with Sixth-Order
Precision and C5 Smoothness

Fahad Sameer Alshammari1,∗, Pakeeza Ashraf2,∗, Saba Mehmood2, Ali Akgul3,4,5,6,7
1 Department of Mathematics, College of Science and Humanities in Alkharj,
Prince Sattam bin Abdulaziz University, Al-Kharj, Saudi Arabia
2 Department of Mathematics, Government Sadiq College Women University, Bahawalpur,
Pakistan
3 Department of Electronics and Communication Engineering, Saveetha School of
Engineering, SIMATS, Chennai, Tamilnadu, India
4 Department of Mathematics, Art and Science Faculty, Siirt University, Siirt, Turkey
5 Department of Computer Engineering, Biruni University, Topkapi, Istanbul, Turkey
6 Near East University, Mathematics Research Center, Department of Mathematics,
Near East Boulevard, Nicosia /Mersin 10, Turkey
7 Applied Science Research Center, Applied Science Private University, Amman, Jordan

Abstract. Subdivision schemes are a crucial component of geometric modeling, widely applied
in curve and surface design. Classical schemes such as the six-point interpolatory scheme and the
quintic B-spline are well known for their efficiency and smoothness. Yet, both have inherent draw-
backs, particularly with respect to approximation order and smoothness. This study introduces a
novel subdivision scheme that blends the strengths of the six-point interpolatory and the quintic
B-spline schemes. The proposed scheme achieves sixth-order approximation and C5 smoothness
while preserving the support size of the six-point scheme. A key feature of the scheme is the
introduction of a tension parameter, which provides flexibility to control the trade-off between
smoothness and approximation order. Moreover, the scheme preserves essential properties such as
monotonicity and convexity under mild conditions. Despite the other higher-order shape-preserving
schemes, which are non-linear and computationally complex, this scheme is linear and stationary.
Experimental results confirm that the proposed scheme consistently produces accurate and visually
elegant curves, outperforming existing schemes in both approximation order and smoothness.
2020 Mathematics Subject Classifications: 65D17, 65D10, 65D07, 65D05
Key Words and Phrases: Subdivision, continuity, monotonicity, convexity, precision

∗Corresponding author.
∗Corresponding author.
DOI: https://doi.org/10.29020/nybg.ejpam.v18i4.6917

Email addresses: f.alshammari@psau.edu.sa (F. S. Alshammari),
pakeeza@gscwu.edu.pk (P. Ashraf), msscholargscwu@gmail.com (S. Mehmood),
aliakgul00727@gmail.com (A. Akgul)

https://www.ejpam.com 1 Copyright: © 2025 The Author(s). (CC BY-NC 4.0)



F. S. Alshammari et al. / Eur. J. Pure Appl. Math, 18 (4) (2025), 6917 2 of 26

1. Introduction

Subdivision schemes (SSs) are iterative processes that refine an initial control polygon
or mesh to generate smooth curves or surfaces, using predefined refinement rules. These
schemes have become fundamental tools in computer graphics due to their extensive ap-
plications in areas such as medical imaging, engineering, visual computing, and image
analysis. The two major types of SSs are approximating and interpolating schemes. Inter-
polatory schemes ensure that the resulting limit curve passes through the initial control
points, while approximating schemes produce curves that may not necessarily interpolate
the original data points.
Mathematically, a stationary univariate subdivision scheme can be described using a re-
finement rule that updates the sequence of points at each level. Let gn = {gnm : m ∈ Z}
represent the set of points at level n. The new refined value at position i is given by:

(Spg
n)i =

∑
m∈Z

ci−2mgnm,

where gn = {gnm : m ∈ Z} represents the sequence of data at refinement level n, Sp denotes
the subdivision operator, which can be represented as a bi-infinite matric
Mc = (ci−2m : i,m ∈ Z) defined by the mask c = {cm : m ∈ Z}. The mask ′c′ is usually
compactly supported, which guarantees that the refinement process remains localized. A
SS is deemed uniformly convergent when the sequence of refined data gn converges uni-
formly to a continuous limiting function [1]. The convergence of the limit function is
influenced by the mask ′c′ and is typically evaluated using the Laurent polynomial [2],
denoted as C(v) =

∑
i∈Z civ

i, i ∈ C and joint spectral radius [3] which together help
determine the smoothness and stability of the scheme.
The design of novel SSs remains an important area of research in curve and surface mod-
eling. De Rham [4] and Chaikin [5] were renowned as the founding fathers of SSs. They
pioneered corner-cutting schemes, which yielded C1-limit curves. A significant advance-
ment came with the four-point binary interpolatory scheme developed by Dyn et al. [6],
which also achieved C1 continuity. Later, Deslauriers and Dubuc introduced a widely used
interpolatory scheme known as the DD scheme [7]. More recently, researchers explored
the idea of shape-controlling parameters. For instance, Mustafa et al. [8] introduced a
framework for constructing binary approximating schemes, allowing better control over
curve behavior. Zheng et al. [9] extended this idea by designing multi-parameter schemes
with high levels of continuity. Similarly, Tan et al. [10] proposed a four-point scheme that
incorporates two shape parameters, achieving C3 smoothness while preserving the shape
of the original data.
SSs are valued for key features such as smoothness and shape preservation. For deeper
discussion of these properties, readers are referred to works such as [11–19].
The quintic B-spline (QB-scheme) [20] and DD six-point interpolatory schemes (SD-
scheme) [21] are stationary SSs. Based on quintic polynomials, these schemes have distinct
advantages and limitations. The quintic B-spline maintains high-order continuity while
preserving the monotonicity and convexity of the initial sequence, having minimal support



F. S. Alshammari et al. / Eur. J. Pure Appl. Math, 18 (4) (2025), 6917 3 of 26

of the fundamental limit curve [3]. The scheme achieves an order of approximation limited
to ’two’ because it can only reproduce polynomials up to linear accuracy. However, the six-
point interpolatory scheme only provides C2-smoothness and, in general, lacks the shape-
preserving property, even though it provides sixth-order accuracy. Moreover, a significant
drawback of the SD-scheme is its interpolatory nature. Artifacts like self-intersections
and inflection points may arise if the edges of the control polygon differ substantially in
length. Furthermore, the refinement mask length, smoothness, and accuracy are essential
criteria for selecting a subdivision algorithm in practical applications. Multiple attempts
are made to increase the smoothness and accuracy by using non-uniform, non-linear, and
non-stationary SSs [22–24]. However, it increases the complexity and demands for a pre-
liminary smoothing process. To improve the smoothness of the interpolatory scheme, Choi
et. al.[25] introduced a family of quasi-interpolatory schemes by enhancing the support
length of the mask. Recently, Yang et al. [26] have presented a family of linear, stationary
SS with C2-smoothness and fourth-order accuracy. However, this scheme in [26] achieves
C2-smoothness; it still relies on a similar construction to the cubic B-spline, which may
limit its potential for further improving approximation accuracy.
Addressing the limitations mentioned above, this work introduces a novel class of station-
ary SSs, combining the benefits of QB and SD-schemes. The presented scheme features a
parameter µ0, bridging the gap between these two schemes and offering a versatile range
of subdivision rules to balance smoothness, precision, and shape-preservation. The main
advancements of our work are: First, by compromising on the interpolation property,
our scheme achieves enhanced smoothness of C5, having the equal support width [−5, 5]
as the SD-scheme. Second, despite only reproducing linear polynomials, our scheme can
attain sixth-order precision using the function that belongs to the Sobolev Space, such
as e ∈ W 6(R). Next, we prove that our proposed scheme maintains monotonicity and
convexity under certain mild conditions. A key advantage of our scheme is its linearity
and stationarity, which distinguishes it from other shape-preserving schemes that are often
non-linear, non-uniform, and computationally complex. To demonstrate the effectiveness
of the proposed SS, several numerical examples are presented.
The design of this paper is as follows. Section 2 discusses fundamental concepts and defi-
nitions. A novel class of binary 6-point SS (SP -scheme) is constructed in Section 3. The
convergence and accuracy of the SP -scheme are analyzed in Section 4. The monotonicity
and convexity preservation of the SP -scheme is discussed in Section 5. Lastly, in Section
6, we present various experimental results.

2. Fundamental concepts and definitions

This section provides fundamental concepts and essential background on stationary
SSs.

Definition 1. A stationary binary SS is a process that iteratively generates a sequence of
control point grids. The control points at level n + 1 are computed from an initial set as



F. S. Alshammari et al. / Eur. J. Pure Appl. Math, 18 (4) (2025), 6917 4 of 26

follows:
gn+1
j =

∑
i∈Z

cj−2ig
n
i , j ∈ Z,

where different masks are used for even and odd-indexed points. In contrast, non-stationary
SS apply recursive refinements using masks that may vary at different levels. The subdivi-
sion operator is assumed to be bounded and is characterized by the sequence c = {ci : i ∈ Z},
referred to as the subdivision mask, which is generally assumed to have finite support.

We define N0 as the set of non-negative integers. The notation
∏

i represents the space
of polynomials with a degree up to i. For a sequence g ∈ L∞(Z), the following difference
operators are defined:

(Dng
n)i = 2n(gni − gni−1), (D2

ng
n)i = 22n(gni−1 − 2gni + gni+1), i ∈ Z.

In this study, we consider a sequence gn = {gni : i ∈ Z} at level n ∈ N0 associated with
grid points 2−nZ.

Definition 2. A binary SS SP applied to a given initial sequence g0 ∈ C∞(Z) is called
uniformly convergent when it generates a continuous limit curve S∞

P g0 ∈ C(R) such that
for at least one nonzero initial data, the sequence of linear interpolants Gn, at grid points
2−nZ, satisfies:

limn→∞||Gn − S∞
P g0||L∞(N) = 0,

where N ⊂ R. If S∞
P g0 has continuous derivatives up to order e, the {SP } scheme is

termed Ce-convergent.
Specifically, as the provided data set δ = {δ0,i : i ∈ Z}, at zeroth level with kronecker delta
the fundamental limit function of SP -scheme is given as,

Ω = S∞
P δ.

Since the subdivision mask c has finite support, the basic limit function Ω also has finite
support. By the linearity of the refinement rule, it follows that

S∞
P g0 =

∑
i∈Z

gn0Ω(· − i).

A well-established condition for the convergence of an SS SP with mask c = {ci : i ∈ Z} is
that its even and odd values satisfy the partition of unity property:∑

i∈Z
cj−2i = 1, j = 0, 1.

Indeed, for an SS defined by mask c, its essential properties, including convergence, regu-
larity, and polynomial reproduction (or generation), are encapsulated in its symbol, defined
by the v-transform:

C(v) =
∑
i∈Z

civ
i, v ∈ C\{0}.



F. S. Alshammari et al. / Eur. J. Pure Appl. Math, 18 (4) (2025), 6917 5 of 26

Since c is finitely supported, the symbol C(v) forms a Laurent polynomial. The norm of
the associated SS SP is given by:

||SP ||∞ = max{
∑
i∈Z

|c2i|,
∑
i∈Z

|c1+2i|}.

Additionally, the Laurent polynomial corresponding to the m-iterated rule Sm
P is:

C [m](v) =
m−1∏
l=0

C(v2
l
) =

∑
i∈Z

c
[m]
i vi,

Where S0
P = SP and the norm of Sm

P is given by:

||Sm
P ||∞ = max{

∑
i∈Z

|c[m]
j+2mi| : j = 0, 1, ..., 2m − 1}.

A crucial property of an SS is its order of approximation, which is related to its ability to
reproduce polynomials.

Definition 3. Let Q0 = {q(m) : m ∈ Z} where q ∈
∏

e with e ∈ N0. We say that a
Stationary SS SP reproduces polynomials in

∏
e if SP is convergent and p = S∞

P Q0. Also,
the SS SP is said to be

∏
e-generating S∞

P Q0 ∈
∏

e [26].

Polynomial reproduction is a valuable property, as it guarantees that a convergent SS
SP with a compactly supported mask reproduces polynomials of degree e, ensuring order of
approximation e+1 [27]. However, generating polynomials of degree e doesn’t necessarily
imply the same level of approximation accuracy e + 1. For example, the QB-scheme
of degree e generates polynomials in

∏
e, but only reproduces linear polynomials (

∏
1),

yielding an order of approximation 2. According to Levin’s work [28], an SS that generates
polynomials of degree e can be modified to reproduce polynomials in

∏
e by applying

a specific operator to the initial sequence. Nevertheless, it is noted that the resulting
order of approximation is often not optimal. However, we introduce a new SS using a
tension parameter, offering improved smoothness C5 and sixth-order of approximation,
while reproducing linear polynomials and generating polynomials in

∏
5.

3. The proposed subdivision scheme

The following section introduces a new six-point binary subdivision schemes. We also
discuss the properties of polynomial generation.

3.1. Construction

We construct a novel class of stationary SSs that combining the advantages of the
QB-scheme and the SD-scheme via a tension parameter. Sequence of data is given by



F. S. Alshammari et al. / Eur. J. Pure Appl. Math, 18 (4) (2025), 6917 6 of 26

gk = {gni : i ∈ Z} at level n, sequence of refined data gn+1 = {gn+1
i : i ∈ Z} is computed

as a linear combination of previous value of gn in the following two rules:
gn+1
2i =

1

64
((µ0)g

n
i−2 + (16µ0)g

n
i−1 + (64− 34µ0)g

n
i + (16µ0)g

n
i+1 + (µ0)g

n
i+2),

gn+1
2i+1 =

1

256
((3− 3µ0)g

n
i−2 + (−25 + 49µ0)g

n
i−1 + (150− 46µ0)g

n
i

+ (150− 46µ0)g
n
i+1 + (−25 + 49µ0)g

n
i+2 + (3− 3µ0)g

n
i+3).

(1)

The construction of our subdivision mask relies crucially on the parameter µ0. The choice
of µ0 significantly affects the regularity, order of approximation, and shape-preservation of
the resulting scheme. Suppose that the sequence of initial data is sampled from a function
on the grid 2−n0Z, we investigate the optimal selection of µ0 as

µ0 = τ2−2n0 , τ ≥ 0. (2)

The key characteristics of our construction in (1) are as follows. First of all, the mask of
SD-scheme is

{
3

256 , 0,
−25
256 , 0,

75
128 , 1,

75
128 , 0,

−25
256 , 0,

3
256

}
, our SP -scheme refinement rule can

be interpreted as a linear modification of the SD-scheme presented below:

SP g
n = SDg

n + 2−2nµ0

64
D2(g

n). (3)

where linear operator is defined as D2
D2(g

n)2i = (D2gn)i−1 + 18(D2gn)i + (D2gn)i+1,

D2(g
n)2i+1 =

−3(D2gn)i−1 + 43(D2gn)i + 43(D2gn)i+1 − 3(D2gn)i+2

4
.

(4)

The new SP -scheme enhances the SD-scheme by integrating the discrete 2nd-derivative
information of gn. Additionally, the parameter µ0 is chosen based on the initial data
density, as specified in (2). This enables the new SP -scheme to attain sixth-order of
approximation, enhancing smoothness without expanding support length. Additionally,
because it exactly reproduces all polynomials up to degree five, i.e.,

∏
5, the scheme ensures

the presence of the smoothing factor (1 + v)6 in its symbol. For further use, the explicit
form of the scheme’s symbol is provided in equation (5):

C(v) =
(1 + v)6

27

{
1

2v
(3− 3µ0) +

1

v2
(−9 + 11µ0) +

1

2v3
(38− 38µ0)

+
1

v4
(−9 + 11µ0) +

1

2v5
(3− 3µ0)

}
. (5)

Remark 1. The suggested scheme is a generalization of the SD-scheme and the QB-
scheme. If µ0 = 0, the scheme is the SD-scheme. Furthermore, when µ0 = 1 it changes to
the quintic B-spline (represented by QB-scheme) with the mask { 1

64 ,
3
32 ,

1
4 ,

13
32 ,

15
32 ,

13
32 ,

1
4 ,

3
32 ,

1
64}.

Then, the quintic B-spline QB may also be used to illustrate the refinement rule in (1) as
follows:

SP g
n = QBg

n − 2−2n (1− µ0)

64
D2(g

n). (6)



F. S. Alshammari et al. / Eur. J. Pure Appl. Math, 18 (4) (2025), 6917 7 of 26

with the operator D2 in (4). This representation will help analyze our SS shape-preserving
characteristics, as we shall see later.

3.2. Property of polynomial generation

When the corresponding fundamental limiting function Ω of SP meets the Strang-
Fix condition of order of approximation e, it is equivalent to the SS symbol having the
smoothing parameter (1+v)e, e ∈ N. Up to degree e, it ensures that the suggested scheme
produces polynomials. At a later stage, we shall discover that the SP -scheme in (3) can
achieve sixth-order approximation even though it reproduces polynomials in

∏
1.

For convenience, the notation used is as follows: for q ∈
∏

5 and n ∈ N0,

qni = q(i2−n) and q(4),n+1 = D2(q
n).

where given operator D2 is defined in (4).

Theorem 1. Assume that SP is the subdivision scheme in (1) with parameter µ0. Each
of the polynomials q represents a part of

∏
5. Polynomial generation satisfies the following

property.
S∞
P q0 = q +

µ0

62
q(4)

Proof. By mathematical induction, to begin, let for any n ∈ N0,

Si+1
P q0 = qn+1 +

µ0

64
q(4), n+1

n∑
i=0

2−2i (7)

Through (3), this identity is true for n = 0 as SP -scheme reproducing polynomials in∏
5. Next, let that Eq.(7) is valid for n ∈ N0. Since q ∈

∏
5, clearly, (D2q

(4), n
i ) = 0.

Consequently, based on the expression in (3), by induction hypothesis, we derive the
equation

Sn+1
P q0 = SP (q

n +
µ0

64
q(4),n

i−1∑
i=0

2−2i)

= (SD + 2−2nµ0

64
D2)(q

n +
µ0

64
q(4),n

i−1∑
i=0

2−2i)

= qn+1 +
µ0

64
q(4), n+1

n∑
i=0

2−2i

So, mathematical induction has been proved. Then, we can easily find that,

q
(4), n+1
2i+m = q(4)(i2−n) +O(2−n) for m=0,1.

It precisely ends the proof.



F. S. Alshammari et al. / Eur. J. Pure Appl. Math, 18 (4) (2025), 6917 8 of 26

4. Evaluation of the proposed scheme

In this section, we discuss the convergence rate, smoothness, and order of approxima-
tion of the proposed scheme.

4.1. Evaluation of convergence and smoothness

The convergence and smoothness of a subdivision scheme, along with the length of
its refinement mask, are fundamental criteria for choosing a suitable subdivision scheme.
To examine such properties for the proposed scheme defined in equation (1), we begin
by introducing some essential concepts. Given an SP -scheme, let SP,1-scheme denotes
the scheme for the divided differences of the initial control points, which possesses the
following property.

Dgn+1 = SP,1Dgn,

Where gn = Sn
P g

0. The characteristic Laurent polynomial C1(v) for SP,1-scheme is given
as:

C1(v) =
2

1 + v
C(v),

Here, C(v) represents the symbol of the scheme SP . More, for the γ-th order difference
scheme SP,γ derived from SP , the associated symbol is given as

Cγ(v) =
∑
i∈Z

c
[γ]
i vi =

(
2

1 + v

)γ

C(v). (8)

Here, C0(v) = C(v) and c
[0]
i = ci. This result offers a standard approach for determining

whether a stationary SS converges.

Theorem 2. [26] Consider two stationary subdivision schemes, SP and SP,1, with respec-
tive symbols C(v) and C1(v). Then, the SP -scheme is uniformly convergence (i.e., C0)
of SP -scheme is the existence of a positive integer r satisfying

∥∥(1
2SP,1

)r∥∥
∞ < 1, where

1
2SP,1 is the scheme denoted by the symbol 1

2C1(v).

Moreover, a sufficient condition for a SP -scheme to be Cm, m ≥ 1, is given below.

Theorem 3. [26] Suppose that SP be a stationary SS with the symbol

C(v) =
1

2
(1 + v)c(v)

with a Laurent polynomial c(v). If the Sc-scheme corresponding with c(v) is Cm, then the
SP -scheme is Cm+1(R).

We will now demonstrate that the scheme in (1) achieves C5 smoothness for a specific
range of µ0 values. Note that when µ0 = 0, the mask in (1) reduces to the SD-scheme,
which is known to be C2. Therefore, we focus on the case where µ0 > 0.



F. S. Alshammari et al. / Eur. J. Pure Appl. Math, 18 (4) (2025), 6917 9 of 26

Theorem 4. Consider the SS SP using the refinement rule is given as (1), having the
tension parameter µ0 specified in (2). Then, SP -scheme generates limit functions of class
C5 for 9

11 < µ0 < 1.

Proof. We begin by investigating the uniform convergence of the SP -scheme. To
continue, we note that the Laurent polynomial C(v) from equation (5) corresponds to
the refinement mask c introduced in equation (1). Using the relation in equation (8), the
Laurent polynomial for the 1st-order difference SP,1-scheme is expressed as follows:

C1(v) =
(1 + v)5

64
{1
2
(3− 3µ0) +

1

v
(−9 + 11µ0) +

1

2v2
(38− 38µ0)

+
1

v3
(−9 + 11µ0) +

1

2v4
(3− 3µ0)}

The mask c[1] = (c
[1]
i : i ∈ Z) which is used in the SP,1-scheme is defined as follows:

c[1] =
1

128
{3− 3µ0,−3 + 7µ0,−22 + 42µ0, 22 + 22µ0, 128− 68µ0, 128− 68µ0,

22 + 22µ0,−22 + 42µ0,−3 + 7µ0, 3− 3µ0}.

Now, ‖1
2SP,1‖ < 1 for −39

49 < µ0 <
195
71 . This implies that, it has C0-continuity.

Next, we will check for C1 smoothness for SP , consider the 2nd-order difference SP,2-
scheme is,

C2(v) =
(1 + v)4

32
{v
2
(3− 3µ0) + (−9 + 11µ0) +

1

2v
(38− 38µ0)

+
1

v2
(−9 + 11µ0) +

1

2v3
(3− 3µ0)}

The mask c[2] = (c
[2]
i : i ∈ Z) which is used in the SP,2-scheme is defined as follows:

c[2] =
1

64
{3− 3µ0,−6 + 10µ0,−16 + 32µ0, 38− 10µ0, 90− 58µ0, 38− 10µ0,

16 + 32µ0,−6 + 10µ0, 3− 3µ0}.

Now, ‖1
2SP,2‖ < 1 for 0 < µ0 < 2. This implies that, it has C1-continuity.

Next, we will check for C2 smoothness for SP , consider the 3rd-order difference SP,3-scheme
is,

C3(v) =
(1 + z)3

16
{v

2

2
(3− 3µ0) + (v)(−9 + 11µ0) +

1

2
(38− 38µ0)

+
1

v
(−9 + 11µ0) +

1

2v2
(3− 3µ0)}

The mask c[3] = (c
[3]
i : i ∈ Z) which is used in the SP,3-scheme is defined as follows:

c[3] =
1

32
{3− 3µ0,−9 + 13µ0,−7 + 19µ0, 45− 29µ0, 45− 29µ0,−7 + 19µ0,

− 9 + 13µ0, 3− 3µ0}.



F. S. Alshammari et al. / Eur. J. Pure Appl. Math, 18 (4) (2025), 6917 10 of 26

Now, ‖1
2SP,3‖ < 1 for 0 < µ0 < 2. This implies that, it has C2-continuity.

Next, we will check for C3 smoothness for SP , consider the 4th-order difference SP,4-scheme
is,

C4(v) =
(1 + z)2

8
{v

3

2
(3− 3µ0) + (v2)(−9 + 11µ0) +

v

2
(38− 38µ0)

+(−9 + 11µ0) +
1

2v
(3− 3µ0)}

The mask c[4] = (c
[4]
i : i ∈ Z) which is used in the SP,4-scheme is defined as follows:

c[4] =
1

16
{3− 3µ0, 5 + 3µ0,−12 + 16µ0, 40− 32µ0,−12 + 16µ0, 5 + 3µ0,

3− 3µ0, }.

Now, ‖1
2SP,4‖ < 1 for 1

2 < µ0 <
3
2 . This implies that, it has C3-continuity.

Next, we will check for C4 smoothness for SP , consider the 5th-order difference SP,5-scheme
is,

C5(v) =
(1 + z)

4
{v

4

2
(3− 3µ0) + (v3)(−9 + 11µ0) +

v2

2
(38− 38µ0)

+(v)(−9 + 11µ0) +
1

2
(3− 3µ0)}

The mask c[5] = (c
[5]
i : i ∈ Z) which is used in the SP,5-scheme is defined as follows:

c[5] =
1

8
{3− 3µ0,−15 + 19µ0, 20− 16µ0, 20− 16µ0,−15 + 19µ0, 3− 3µ0}.

Now, ‖1
2SP,5‖ < 1 for 11

19 < µ0 <
27
19 . This implies that, it has C4-continuity.

Figure 1: Basic limit function for µ0 = 0.6, 0.7, 0.8, 0.9.



F. S. Alshammari et al. / Eur. J. Pure Appl. Math, 18 (4) (2025), 6917 11 of 26

Next, we will check for C5 smoothness for SP , consider the 6th-order difference SP,6-
scheme is,

C6(v) =
1

2
{v

5

2
(3− 3µ0) + (v4)(−9 + 11µ0) +

v3

2
(38− 38µ0)

+(v2)(−9 + 11µ0) +
v

2
(3− 3µ0)}

The mask c[6] = (c
[6]
i : i ∈ Z) which is used in the SP,6-scheme is defined as follows:

c[6] =
1

4
{3− 3µ0,−18 + 22µ0, 38− 38µ0,−18 + 22µ0, 3− 3µ0}.

Now, ‖1
2SP,6‖ < 1 for 9

11 < µ0 < 1. This implies that, it has C5-continuity.
Next, for C6 smoothness, the divided difference of 7th-order does not exist making SP -
scheme of C5-continuity.

Remark 2. The fundamental limit function of the SP -scheme is expressed as φ = S∞
P δ,

where δ = {δ0,i : i ∈ Z} is a sequence at level 0, and δi,0 denotes the Kronecker delta.
It is important to note that the support of this function is [-5,5], which is equivalent to
the SD-scheme. A visualization of the basic limit function of the SP -scheme for varying
values of the tension parameter µ0 is provided in Figure 1.

4.2. Order of approximation

We investigate the order of approximation of the SP -scheme with the mask in (1),
focusing on its connection to the parameter µ0. Specifically, our objective is to demonstrate
that the SP -scheme achieves sixth-order of approximation with an appropriate choice of
µ0. Subsequently, we will prove that if the initial data sequence is given by g0 = {g(i2−n0) :
i ∈ Z}, where g is a smooth function in W6

∞(R), then the scheme SP generates a limit
function S∞

P gn0 such that
‖g − S∞

P g0‖ ≤ kg2
−6n0 (9)

To simplify our presentation, we use the notation ĝn, for n ∈ Z+, to represent the sequence
obtained by sampling a smooth function g at the grid points 2−nZ. That is,

ĝn = {g(i2−n : i ∈ Z}, n ∈ Z+. (10)

Lemma 1. Let SP be the SPSS with the parameter µ0 = τ2−2n0 . Let ĝn, n ≥ n0 be a
sequence of data points obtained by sampling a function g ∈ W6

∞(R) with density 2−n.
Then

‖Spĝ
n − ĝn+1‖∞ ≤ kg2

−2(2n+n0)

for some constant kg > 0 depending on g .



F. S. Alshammari et al. / Eur. J. Pure Appl. Math, 18 (4) (2025), 6917 12 of 26

Proof. Let us first estimate the case of odd rule, i.e.,|SP ĝ
n
2i+1 − ĝn+1

2i+1|. To do this, the
proposed scheme is represented using SD-scheme:

(SP ĝ
n)2i+1 = SDĝ

n
2i+1 −

τ

64
2−2n0 .2−2n(D2g

n
2i+1) (11)

Clearly D2g
n
2i+1 = O(2−4n) for i ∈ Z because g is C5 and giv is bounded. It follows that

|(SP ĝ
n)2i+1 − ĝn+1

2i+1| ≤ |(SDĝ
n)2i+1 − ĝn+1

2i+1|+O(2−2(2n+n0)) = O(2−2(2n+n0)) (12)

Since, the SD-scheme provides the sixth-order approximation and g ∈ W6
∞(R), it readily

holds that
|(SDĝ

n)2i+1 − ĝn+1
2i+1| = O(2−6n).

Combining this with (4) verifies the lemmas claim,

‖Spĝ
n − ĝn+1‖∞ ≤ kg2

−2(2n+n0)

In a similar way, we can verify the even rule.

Theorem 5. Assume that the initial data sequence g0 := {g(2−n0i) : i ∈ Z} is obtained
by sampling a function g ∈ W 6

∞(R), i.e., g0 := ĝn0. Let {SP } denote the SPSS defined by
the refinement rule given in equation (2.1), and suppose the shape parameter µ0 satisfies
0 ≤ µ0 ≤ τ2−n0, for some constant τ > 0. Under these conditions, the following error
estimate holds:

‖S∞
P g0 − g‖∞ ≤ kg,τ2

−6n0 ,

where kg,τ > 0 is a constant that depends only on the function g and the parameter τ , but
is independent of the refinement level n0.

Proof. Let D = {i2−n : i ∈ Z, n ∈ Z+} be the set of dyadic points. Suppose the initial
data g0 is defined on the grid 2−n0Z, and let the iterated subdivision sequence be given
by gn = Sn−n0

P ĝn0 . Then the sequence gn satisfies the telescoping identity

gn = SP ĝ
n−1 +

n−2∑
l=n0

Sn−l−1(SP ĝ
l − ĝl+1).

Assume the proposed SP -scheme is stable i.e., there exists a constant k > 0 such that
||Sm

P ||∞ ≤ k for any m ∈ N. Then, for any n ≥ n0 and i ∈ Z, the following estimate holds:

|gni − g(i2−n)| ≤ ||SP ĝ
n−1 − ĝn||∞ +

n−2∑
l=n0

||Sn−l−1
P ||∞||SP ĝ

l − ĝl−1||∞

≤ k
n−1∑
l=n0

||SP ĝ
n − ĝl+1||∞.



F. S. Alshammari et al. / Eur. J. Pure Appl. Math, 18 (4) (2025), 6917 13 of 26

Our next step is to approximate the term ||Sĝl − ĝl+1||∞. By assumption,
0 ≤ µ0 ≤ τ2−2n0 for some τ > 0. Thus

||Sĝl − ĝl+1||∞ ≤ kg, τ2
−2(n0+l)

for some l ≥ g0 with a fixed value kg, τ > 0 depending on g and τ . It leads to the bound

n−1∑
l=n0

||SP ĝ
n − ĝl+1||∞ ≤ kg

n−1∑
l=n0

2−2(n0+l) ≤ kg2
−6n0 .

5. Shape-preserving properties

Smoothness, precision, and shape preservation, such as monotonicity and convexity, are
all important aspects of an SS in practical and commercial contexts. Thus, this section’s
goal is to examine our SS shape-preserving characteristics. To proceed, we express the
refinement rules given in equation (1) in the following rewritten form:

SP g
n = QBg

n − 2−2n−6(1− µ0)D2(g
n). (13)

Where D2 is given in (4) and QB is the quintic B-spline with the mask { 1
64 ,

3
32 ,

1
4 ,

13
32 ,

15
32 ,

13
32 ,

1
4 ,

3
32 ,

1
64}. If µ0 = 0, our scheme reduces to the SD-scheme. It is a known fact that the

SD-scheme has the inherent drawback of being interpolatory. If the initial control points
are irregularly distributed, this property can lead to undesirable artifacts in the resulting
curve, such as self-intersections. Thus, it is appropriate to examine the parameter µ0 away
from zero. However, if µ0 = 1, the scheme transforms into the QB-scheme. In essence,
the parameter µ0 acts as a balancing factor between the SD-scheme and QB-scheme. By
conducting numerical experiments on various µ0, we note that a suitable choice of µ0 is
in the range

9

11
< µ0 < 1. (14)

5.1. Preservation of monotonicity

We first demonstrate that SP -scheme in (13) preserves monotonicity under a mild
condition. Specifically, we show that (Dgn+1)i ≥ 0 whenever (Dgn)i ≥ 0 for all i ∈ Z and
n ∈ N0. It is important to note that the QB-scheme fulfills the following equations before
proceeding:

(DQBg
n)2i =

5

32
(Dgn)i−1 +

15

32
(Dgn)i +

11

32
(Dgn)i+1

1

32
(Dgn)i+2,

(DQBg
n)2i+1 =

1

32
(Dgn)i−1 +

11

32
(Dgn)i +

15

32
(Dgn)i+1

5

32
(Dgn)i+2.

(15)



F. S. Alshammari et al. / Eur. J. Pure Appl. Math, 18 (4) (2025), 6917 14 of 26

which means that QB-scheme is preserves the monotonicity . Then, given (13) and (15),
we can derive the following identities.

(Dgn)2i = (
5(Dgn)i−1 + 15(Dgn)i + 11(Dgn)i+1 + (Dgn)i+2

32
)

− 2−2n (1− µ0)

128
(−3(D3gn)i−1 + 36(D3gn)i + 7(D3gn)i+1),

(Dgn)2i+1 = (
(Dgn)i−1 + 11(Dgn)i + 15(Dgn)i+1 + 5(Dgn)i+2

32
)

− 2−2n (1− µ0)

128
(7(D3gn)i + 36(D3gn)i+1 − 3(D3gn)i+2).

(16)

Definition 4. (Condition-M ). Suppose given a number µ0 and gn0 be an initial sequence
defined on the grid 2−n0Z. Let the sequence gn0 satisfies ’condition-M’ when

(Dgn0)i ≥ 2−2n0
(1− µ0)

k
|(D3gn0)i|, (17)

Where k > 0.

Lemma 2. Let gn0 denote the initial data associated with the grid points 2−n0Z, and
suppose that gn+1 = SP g

n, n ≥ n0, with SP in (1). If gn0 fulfills the condition-M. Thus
for all n ≥ n0,

(Dgn)i ≥ 2−2n (1− µ0)

k
|(D3gn)i|, (18)

where k > 0.

Proof. Using mathematical induction, we establish this lemma. Given that the initial
sequence fulfills ’condition-M’, inequality (18) holds. We now assume that (18) is true for
n > n0 and strive to prove that

(Dgn+1)i ≥ 2−2n+1 (1− µ0)

k
|(D3gn+1)i|.

To prove the result, we consider the following two cases: (I) i = 2j and (II) i = 2j + 1.
We begin with Case (I), i = 2j, and we can readily obtain the expression using equation
(16)

(D3gn+1)2j = 22(n+1) 1

2
(−3(D3gn+1)2j−1 + 36(D3gn+1)2j + 7(D3gn+1)2j+1)

= (
27(D3gn)j−1 + 153(D3gn)j + 133(D3gn)j+1 + 7(D3gn)j+2

32
)− (1− µ0)

128
(417(D3gn)j−1 − 851(D3gn)j + 451(D3gn)j+1 − 17(D3gn)j+2)

Thus, it follows that

|(D3gn+1)2j | ≤
618

128
|(D3gn)j−1|+

2926

128
|(D3gn)j |+

162

128
|(D3gn)j+1|+

90

128
|(D3gn)j+2|.

This, together with (9), induces the relation



F. S. Alshammari et al. / Eur. J. Pure Appl. Math, 18 (4) (2025), 6917 15 of 26

(Dgn+1)2j − 2−2(n+1) (1−µ0)
k |(D3gn+1)2j | ≥ 5

32((Dgn)j−1 − 2−2(n) (1−µ0)
k |(D3gn)j−1|)

+15
32((Dgn)j − 2−2(n) (1−µ0)

k |(D3gn)j |)
+11

32((Dgn)j+1 − 2−2(n) (1−µ0)
k |(D3gn)j+1|)

+ 1
32((Dgn)j+2 − 2−2(n) (1−µ0)

k |(D3gn)j+2|).

Based on the inductive assumption for n, the nonnegativity of the final term follows. We
now turn our attention to the case i = 2j + 1, where the identities in (16) enable us to

|(D3gn+1)2j+1| = |( 7
16

(D3gn)j−1 +
133

16
(D3gn)j +

153

16
(D3gn)j+1 +

27

16
(D3gn)j+2)

+
(1− µ0)

128
(−17(D3gn)j−1 + 451(D3gn)j − 851(D3gn)j+1

+417(D3gn)j+2)|

≤ 90

128
|(D3gn)j−1|+

162

128
|(D3gn)j |+

2926

128
|(D3gn)j+1|+

618

128
|(D3gn)j+2|.

By the same argument as given above, a direct calculation using (16) shows that

(Dgn+1)2j+1 − 2−2(n+1) (1− µ0)

k
|(D3gn+1)2j+1| ≥ 0.

Where k > 0.

Theorem 6. Let the initial sequence gn0 be defined on the grid 2−n0Z, and consider the
subdivision process given by

gn+1 = SP g
n, for all n ≥ n0,

where SP is the subdivision scheme defined in equation (1). If the initial data gn0 satisfies
Condition-M, then the discrete first differences satisfy

(Dgn)i ≥ 0, for all n ≥ n0 and i ∈ Z.

Therefore, the scheme SP preserves monotonicity.

Proof. By lemma 2, we have already shown that (Dgn)i ≥ 2−2n (1−µ0)
k |(D3gn)i| ≥ 0 for

any n ≥ n0. Hence, the result follows directly.

5.2. Preservation of convexity

Next, we show that SP -scheme maintains convexity under a mild condition. For sim-
plicity, we introduce the following notation

∇gni = gni − gni−1 and ∆gni = gni−1 − 2gni + gni+1.



F. S. Alshammari et al. / Eur. J. Pure Appl. Math, 18 (4) (2025), 6917 16 of 26

Lemma 3. Let QB denote the quintic B-spline. Thus, for any i ∈ Z ,
(D2QBg

n)2i =
1

16
(D2gn)i−1 +

7

16
(D2gn)i +

7

16
(D2gn)i+1 +

1

16
(D2gn)i+2,

(D2QBg
n)2i+1 =

4

16
(D2gn)i +

8

16
(D2gn)i+1

4

16
(D2gn)i+2.

(19)

This shows that the QB-scheme preserves convexity.

Lemma 4. Assume SP is the SS using the refinement rule given in (1), suppose that µ0

is given in (2). When gn+1 = SP g
n for n ≥ n0. Next, we obtain

(D2gn+1)2i = (
1

16
(D2gn)i−1 +

7

16
(D2gn)i +

7

16
(D2gn)i+1 +

1

16
(D2gn)i+2)

−(1− µ0)

32
2−n(10∇(D3gn)i+1 + 10∇(D3gn)i+2),

(D2gn+1)2i+1 = (
4

16
(D2gn)i +

8

16
(D2gn)i+1

4

16
(D2gn)n+2)−

(1− µ0)

32
2−n

(−3∇(D3gn)i+1 + 26∇(D3gn)i+2 − 3∇(D3gn)i+3). (20)

Proof. We now concentrate on proving the first equation in (20). We note that QB

implies that the quintic B-spline, using (13) and Lemma 3 to obtain

(D2gn+1)2n = (
1

16
(D2gn)i−1 +

7

16
(D2gn)i +

7

16
(D2gn)i+1 +

1

16
(D2gn)i+2)

−(1− µ0)

32
(10(D2gn)i−1 − 10(D2gn)i − 10(D2gn)i+1 + 10(D2gn)i+2)

= (
1

16
(D2gn)i−1 +

7

16
(D2gn)i +

7

16
(D2gn)i+1 +

1

16
(D2gn)i+2)

−(1− µ0)

32
2−n(10∇(D3gn)i+1 + 10∇(D3gn)i+2).

By applying equation (13) and Lemma 3, we obtain the following expression through basic
calculation

(D2gn+1)2n+1 = (
4

16
(D2gn)i +

8

16
(D2gn)i+1 +

4

16
(D2gn)i+2)−

(1− µ0)

32
(−3D2gni−1 + 32D2gni − 58D2gni+1 + 32D2gni+2 − 3D2gni+3)

= (
4

16
(D2gn)i +

8

16
(D2gn)i+1

4

16
(D2gn)n+2)−

(1− µ0)

32
2−n

(−3∇(D3gn)i+1 + 26∇(D3gn)i+2 − 3∇(D3gn)i+3).

Lemma 5. Suppose Sp be the SS given by (1), acting on the initial data gn0 at grid points
2−n0Z. We then define the sequence gn+1 = Spg

n for n ≥ n0, we obtain

∇(D3gn+1)2i = (
−5

32
∇(D3gn)i−1 +

11

32
∇(D3gn)i +

5

16
∇(D3gn)i+1),



F. S. Alshammari et al. / Eur. J. Pure Appl. Math, 18 (4) (2025), 6917 17 of 26

∇(D3gn+1)2i+1 = (
3

64
∇(D3gn)i−1 −

24

64
∇(D3gn)i +

85

64
∇(D3gn)i+1

−32

64
∇(D3gni+2)). (21)

Proof. From Lemma 4, it follows by direct calculation that

2−(n+1)∇(D3gn+1)2i+1 = (−3(D2gn+1)2i−1 + 32(D2gn+1)2i − 58(D2gn+1)2i+1

+32(D2gn+1)2i+2 − 3(D2gn+1)2i+3)

= (
10

16
(D2gn)i−1 +

20

16
(D2gn)i −

40

16
(D2gn)i+1 −

20

16
(D2gn)i+2

+
30

16
(D2gn)i+3)−

(1− µ0)

32
2−n(130∇D3fn

i

−290∇D3fn
i+1 + 190∇D3gni+2 − 30∇D3gni+3)

∇(D3gn+1)2i+1 = (−30

32
∇(D3gn)i−1 −

40

32
∇(D3gn)i +

62

32
∇(D3gn)i+1

− 2

32
∇(D3gni+2)).

2−(n+1)∇(D3gn+1)2i = (10(D2gn+1)2i−1 − 10(D2gn+1)2i − 10(D2gn+1)2i+1

+10(D2gn+1)2i+2)

= (
30

16
(D2gn)i−1 −

20

16
(D2gn)i −

40

16
(D2gn)i+1 +

20

16
(D2gn)i+2

+
10

16
(D2gn)i+3)−

(1− µ0)

32
2−n(−30∇D3fn

i

+190∇D3fn
i+1 − 290∇D3gni+2 + 130∇D3gni+3)

∇(D3gn+1)2i = (− 2

32
∇(D3gn)i−1 +

62

32
∇(D3gn)i −

40

32
∇(D3gn)i+1

−30

32
∇(D3gn)i+1).

Definition 5. (Condition-C). Given a number µ0 as in (2), and an initial sequence gn0

defined on 2−n0Z, we say gn0 fulfills ’condition-C’ if

(D2gn0)j ≥
(1− µ0)

k
2−n0 |∇(D3gn)j+1|, ∀j ∈ Z. (22)

Where, k > 0.

Theorem 7. Suppose SP be the SS in (1) with µ0 in (14). Assume that the initial
sequence gn0 be connected to 2−n0Z, and suppose that it fulfills the condition-C. Then,
for any n > n0, (∆gn)i ≥ 0 for i ∈ Z, which implies that the SP -scheme preserves the
convexity.



F. S. Alshammari et al. / Eur. J. Pure Appl. Math, 18 (4) (2025), 6917 18 of 26

Table 1: Maximum errors in L∞-norm and corresponding order of approximations for g1 and g2 given in
Example 1.

Initial sampling interval g1 ∈ W 6
∞(R) g2 ∈ W 6

∞(R)
||.||∞ error Order ||.||∞ error Order

2−1 5.9981E−05 – 5.8756E−05 –
2−2 9.37E−07 6.0 9.32E−07 5.9
2−3 1.5E−08 6.0 1.5E−08 6.0
2−4 2.288E−10 6.0 2.3E−10 6.0
2−5 3.6E−12 6.0 3.7E−12 6.0
2−6 5.59E−14 6.0 5.6E−14 6.0

Note: The results confirm that both test functions g1 and g2 achieve sixth-order
convergence in the L∞-norm as the sampling interval decreases.

Proof. Applying the principle of mathematical induction, prove that for any n ≥ n0,

(D2gn)j ≥
(1− µ0)

k
2−n|∇(D3gn)j+1|, ∀ j ∈ Z. (23)

The initial sequence gn0 satisfies by hypothesis. Assuming the inequality is true for some
i. Now, we check if it is true for i + 1, starting with the case where j = 2i. By applying
Lemma 4 and 5, we deduce

(D2gn+1)2i − 2−n−1 (1−µ0)
k | ∇(D3gn+1)2i+1 |≥ 4

16((D
2gn)i−1 − 2−n (1−µ0)

k |∇(D3gn)i−1|)
+ 8

16((D
2gn)i − 2−n (1−µ0)

k |∇(D3gn)i|)
+ 4

16((D
2gn)i+1 − 2−n (1−µ0)

k |∇(D3gn)i+1|).

According to our inductive hypothesis, the final result is positive. We now turn j = 2i−1,
and use Lemma 5 to drive

(D2gn+1)j − 2−n−1 (1−µ0)
k | ∇(D3gn+1)j |≥ 1

16((D
2gn)i−1 − 2−n (1−µ0)

k |∇(D3gn)i−1|)
+ 7

16((D
2gn)i − 2−n (1−µ0)

k |∇(D3gn)i|)
+ 7

16((D
2gn)i+1 − 2−n (1−µ0)

k |∇(D3gn)i+1|)
+ 1

16((D
2gn)i+2 − 2−n (1−µ0)

k |∇(D3gn)i+2|).

Based on the assumption for the case n, we can conclude that the final term is positive,
thereby completing the proof.

6. Numerical Applications

This section presents numerical experiments to demonstrate the effectiveness of the
proposed subdivision schemes.



F. S. Alshammari et al. / Eur. J. Pure Appl. Math, 18 (4) (2025), 6917 19 of 26

Example 1. (Order of approximation)
To verify the order of approximation of SP -scheme given in equation (1), we consider the
following two smooth functions:

g1(x) = x6, g2(x) = x6 cos(x).

The initial data sequences are sampled from each function over the interval [0, 0.2] at
sampling densities 2−n0, where n0 = 0, . . . , 6. For every refinement level, the parameter
µ0 is selected as µ0 = 2−n0, corresponding to the choice τ = 1 in equation (2).

Table 1 presents the computed maximum errors, based on the L∞-norms, along with
the estimated rate of convergence. The functions g1 and g2 belong to the Sobolev space
and achieve the optimal order of approximation six. These numerical results validate the
theoretical order of approximations described in Theorem 5.

Example 2. (Tension parameter)
To show the effect of the parameter µ0 on the limit curves, we provide a numerical
example. The dotted lines represent the control polygon, and the control points are
marked with circles. Figure 2 presents the curves produced by the proposed scheme for
µ0 = 3/10, 5/10, 8/11, 10/11. As µ0 increases, the curves become smoother and deviate
more from the original control polygon. It shows how the tension parameter affects the
limit curve, resulting in different rates of convergence.

Figure 2: Subdivision limit curves produced with the parameter µ0 = 3/10, 5/10, 8/11 and 10/11, respectively.

Example 3. (Curve design)
This example illustrates the numerical results of our subdivision SP -scheme with its parent
schemes. Figure 3 and 4 show the limit curve when µ0 is set to zero, the curve corresponds
to the limit curve of the SD-scheme algorithm through this results in an undesirable artifact
and when µ0 increases toward 1 it becomes QB-scheme the curves deviate more significantly
from the control polygon.



F. S. Alshammari et al. / Eur. J. Pure Appl. Math, 18 (4) (2025), 6917 20 of 26

(a) SP -scheme (b) SD-scheme

(c) Zoomed area
Figure 3: 2D-model generated by SD-scheme and SP -scheme at µ0 = 10/11.



F. S. Alshammari et al. / Eur. J. Pure Appl. Math, 18 (4) (2025), 6917 21 of 26

(a) SP -scheme (b) QB-scheme

(c) Zoomed area
Figure 4: 2D-model generated by SP -scheme and QB-scheme at µ0 = 10/11.

Example 4. (Monotonicity preservation)
This example illustrates the effectiveness of the SS SP in Monotonicity preservation.
Suppose that initial data set gn0 be given on [-5,5] as follows:
gn0 = {0.007, 0.018, 0.047, 0.119, 0.269, 0.500, 0.731, 0.881, 0.952, 0.982, 0.993}.
A strictly monotone initial dataset gn0 is provided, and the scheme is applied with varying
tension parameter values µ0 = 0.5, 0.6, 0.7, 0.8, 0.9. The results demonstrate that the limit
curves generated by the scheme preserve the monotonicity of the initial data in Figure 5(a),
and their first derivatives in Figure 5(b) are also examined. The scheme’s performance is
showcased within the applicable domain, confirming its ability to maintain monotonicity.



F. S. Alshammari et al. / Eur. J. Pure Appl. Math, 18 (4) (2025), 6917 22 of 26

(a) (b)
Figure 5: Monotonicity preservation. Limit function for (a) µ0 = 0.5, 0.6, 0.7, 0.8, 0.9 and (b) their first
derivatives.

Example 5. (Convexity preservation)
This example demonstrates the effectiveness of the SS SP in Convexity preservation. Sup-
pose that initial data set gn0 be given on [-5,5] as follows:
gn0 = {0.00, 0.0816, 0.3265, 0.7347, 1.3061, 2.0408, 2.9388, 4.00}.
A strictly convex initial dataset gn0 is provided, and the scheme is applied with varying
tension parameter values mu0 = (0.5, 0.6, 0.7, 0.8, 0.9). The results illustrate that the limit
curves generated by the scheme preserve the Convexity of the initial data in Figure 6 (a),
and their second derivatives in Figure 6 (b) are also examined. The scheme’s performance
is showcased within the applicable domain, confirming its ability to maintain Convexity.

(a) (b)
Figure 6: Convexity preservation. Limit functions for (a) µ0 = 0.5, 0.6, 0.7, 0.8, 0.9 and (b) their second
derivatives.

Example 6. (Discontinuous function)
In this example, we analyze the performance of the proposed scheme on a discontinuous



F. S. Alshammari et al. / Eur. J. Pure Appl. Math, 18 (4) (2025), 6917 23 of 26

Table 2: Comparison with existing schemes.

Scheme Type Support Order Cn

Binary 3-point [29] Approximating 5 3 C3

Binary 4-point [30] Interpolating 6 4 C1

Binary 4-point [31] Approximating 7 4 C2

Binary 5-point [32] Approximating 9 2 C4

Binary 4-point [26] Approximating 6 4 C2

Binary 6-point [21] Interpolating 10 6 C2

SP -scheme Approximating 10 6 C5

function defined as:

f(x) =

{
1, if 0 ≤ x < 0.4

0, if 0.4 ≤ x ≤ 0.8

This step function is sampled uniformly at nine points from x=0 to x=0.8. The scheme is
applied using a fixed tension parameter µ0 = 9/11. In Figure 7, The results show that the
scheme handles the jump without introducing oscillations and preserves locality near the
discontinuity, demonstrating its robustness on non-smooth data.

]
Figure 7: Graph of the discontinuous function g(x).

7. Comparison and conclusion

In conclusion, we have proposed a new linear, stationary subdivision scheme that
unifies the strengths of the quintic B-spline and six-point interpolatory scheme. The pro-
posed scheme achieves C5-smoothness and sixth-order approximation with same support



F. S. Alshammari et al. / Eur. J. Pure Appl. Math, 18 (4) (2025), 6917 24 of 26

as of interpolatory scheme. The inclusion of a tension parameter allows for shape con-
trol, making the scheme flexible for preserving essential properties such as convexity and
monotonicity. A comparative analysis in Table 2 with existing schemes highlights that the
proposed scheme offers superior smoothness (C5 continuity) and higher approximation
order (sixth order), while maintaining same support of interpolatory scheme. Numerical
experiments further confirm its effectiveness in producing smooth, accurate, and visually
elegant curves, demonstrating its potential for practical applications in geometric model-
ing and computer-aided design. Despite these strengths, the scheme has some limitations.
At present, it is restricted to univariate curve subdivision, and the choice of the tension
parameter may vary depending on the data set or application. As future work, the scheme
can be extended to bivariate surfaces using the tensor product approach, which is expected
to further enhance smoothness and approximation order. Another promising direction is
the development of adaptive strategies for selecting the tension parameter automatically,
along with exploring applications in CAD, image processing, and scientific visualization.

Acknowledgements

The authors extend their gratitude to the Deanship of Scientific Research at Prince
Sattam Bin Abdulaziz University, Kingdom of Saudi Arabia.

References

[1] C. Conti, L. Romani, and J. Yoon. Approximation order and approximate sum rules
in subdivision. Journal of Approximation Theory, 207:380–401, 2016.

[2] N. Dyn. Subdivision schemes in cagd. Advances in Numerical Analysis, 2:36–104,
1992.

[3] N. Dyn and D. Levin. Subdivision schemes in geometric modelling. Acta Numerica,
11:73–144, 2002.

[4] G. de Rham. Sur une courbe plane. Journal de Mathématiques Pures et Appliquées,
35:25–42, 1956.

[5] G. M. Chaikin. An algorithm for high-speed curve generation. Computer Graphics
and Image Processing, 3(4):346–349, 1974.

[6] N. Dyn, D. Levin, and J. A. Gregory. A 4-point interpolatory subdivision scheme for
curve design. Computer Aided Geometric Design, 4(4):257–268, 1987.

[7] G. Deslauriers and S. Dubuc. Symmetric iterative interpolation processes. Construc-
tive Approximation, 5:49–68, 1989.

[8] G. Mustafa, P. Ashraf, and N. Saba. A new class of binary approximating subdivision
schemes. Jurnal Teknologi, 78(4-4), 2016.

[9] J. Tan, X. Zhuang, and L. Zhang. A new four-point shape-preserving c3 subdivision
scheme. Computer Aided Geometric Design, 31(1):57–62, 2014.

[10] H. C. Zheng, S. P. Huang, F. Guo, and G. H. Peng. Designing multi-parameter curve
subdivision schemes with high continuity. Applied Mathematics and Computation,
243:197–208, 2014.



F. S. Alshammari et al. / Eur. J. Pure Appl. Math, 18 (4) (2025), 6917 25 of 26

[11] I. Yadshalom. Monotonicity preserving subdivision schemes. Journal of Approxima-
tion Theory, 74(1):41–58, 1993.

[12] G. Mustafa and R. Bashir. Four-point n-ary interpolating subdivision schemes. In-
ternational Journal of Mathematics and Mathematical Sciences, page 893414, 2013.

[13] G. Mustafa, P. Ashraf, and J. Deng. Generalized and unified families of interpolating
subdivision schemes. Numerical Mathematics: Theory, Methods and Applications,
7(2):193–213, 2014.

[14] G. Akram, K. Bibi, K. Rehan, and S. S. Siddiqi. Shape preservation of 4-point inter-
polating non-stationary subdivision scheme. Journal of Computational and Applied
Mathematics, 319:480–492, 2017.

[15] H. Yang and J. Yoon. A shape preserving c2 non-linear, non-uniform subdivision
scheme with fourth-order accuracy. Applied and Computational Harmonic Analysis,
60:267–292, 2022.

[16] H. Yang, J. Kim, and J. Yoon. A shape preserving corner cutting algorithm with an
enhanced accuracy. Applied Mathematics Letters, 137:108487, 2023.

[17] C. Conti and K. Hormann. Polynomial reproduction for univariate subdivision
schemes of any arity. Journal of Approximation Theory, 163(4):413–437, 2011.

[18] N. Dyn, D. Levin, and J. Yoon. A new method for the analysis of univariate nonuni-
form subdivision schemes. Constructive Approximation, 40:173–188, 2014.

[19] Y. Kim, H. Yang, and J. Yoon. Shape-preserving subdivision scheme with the third-
order accuracy and c2 smoothness. Mathematics and Computers in Simulation,
235:160–174, 2025.

[20] J. Tan, Y. Yao, H. Cao, and L. Zhang. Convexity preservation of five-point bi-
nary subdivision scheme with a parameter. Applied Mathematics and Computation,
245:279–288, 2014.

[21] S. S. Siddiqi, W. us Salam, and K. Rehan. A new non-stationary binary 6-point
subdivision scheme. Applied Mathematics and Computation, 268:1227–1239, 2015.

[22] S. Zouaoui, S. Amat, S. Busquier, and J. Ruiz. On a nonlinear three-point subdivision
scheme reproducing piecewise constant functions. Mathematics, 10(15):2790, 2022.

[23] N. Dyn, K. Hormann, and C. Mancinelli. Non-uniform interpolatory subdivision
schemes with improved smoothness. Computer Aided Geometric Design, 94:102083,
2022.

[24] R. Hameed, A. Jamil, and J. Younis. Enhanced curve representations using piecewise
non-linear binary subdivision scheme. Arab Journal of Basic and Applied Sciences,
31(1):554–570, 2024.

[25] S. W. Choi, B. G. Lee, Y. J. Lee, and J. Yoon. Stationary subdivision schemes
reproducing polynomials. Computer Aided Geometric Design, 23(4):351–360, 2006.

[26] H. Yang, K. Kim, and J. Yoon. A family of c2 four-point stationary subdivision
schemes with fourth-order accuracy and shape-preserving properties. Journal of Com-
putational and Applied Mathematics, 446:115843, 2024.

[27] N. Dyn. Interpolatory subdivision schemes. In Tutorials on Multiresolution in Geo-
metric Modeling, pages 25–50. 2002.

[28] A. Levin. Polynomial generation and quasi-interpolation in stationary non-uniform



F. S. Alshammari et al. / Eur. J. Pure Appl. Math, 18 (4) (2025), 6917 26 of 26

subdivision. Computer Aided Geometric Design, 20(1):41–60, 2003.
[29] G. Mustafa and A. H. Randhawa. Complete analysis of 3-point binary subdivision

scheme. Journal of Pure and Applied Science, 24:33, 2014.
[30] C. Beccari, G. Casciola, and L. Romani. A non-stationary uniform tension controlled

interpolating 4-point scheme reproducing conics. Computer Aided Geometric Design,
24(1):1–9, 2007.

[31] N. Dyn, M. S. Floater, and K. Hormann. C2 four-point subdivision scheme with
fourth-order accuracy and its extensions. In Mathematical Methods for Curves and
Surfaces, pages 145–156, 2005.

[32] S. S. Siddiqi and N. Ahmad. A new five-point approximating subdivision scheme.
International Journal of Computational Mathematics, 85(1):65–72, 2008.