DOI: 10.24193/subbmath.2020.2.08

## Constrained visualisation using

## Shepard-Bernoulli interpolation operator

### Teodora C˘ atina¸s

Abstract. We consider Shepard-Bernoulli operator in order to solve a problem of interpolation of scattered data that is constrained to preserve positivity, using the technique described by K.W. Brodlie, M.R. Asim and K. Unsworth (2005).

We also give some numerical examples.

Mathematics Subject Classification (2010):41A29, 41A05, 41A25, 41A35.

Keywords:Scattered data, Shepard operator, constrained interpolation.

### 1. Introduction

The interpolation operators and the radial basis functions are the usual tools used for approximating scattered data. Sometimes we may have data that have to pre- serve some constraints, subject to certain physical laws (e.g., the densities, percentage mass concentrations in a chemical reaction, volume and mass are always positive, see [1], [2]); such cases require to impose some special conditions to the interpolants.

The Shepard method is a well suited method for multivariate interpolation of very large scattered data sets, but it does not guarantee to preserve positivity.

In [3] and [4] there have been introduced some combined Shepard operators of Bernoulli type which diminish the drawbacks of the Shepard operator. In [4] the com- bined operators are obtained using the classical and the modified Shepard operator, introduced, respectively, in [12] and [8]. They preserve the advantages and improve the reproduction qualities, have better accuracy and better computational efficiency.

We recall some results from [6]. Bernoulli polynomials are defined by:

B0(x) = 1,

B_{n}^{0}(x) =nB_{n−1}(x), n≥1,
Z 1

0

B_{n}(x)dx= 0, n≥1.

(1.1)

Forf ∈C^{m}[a, b], the univariate Bernoulli interpolant is given by
(B_{m}f)(x) :=B_{m}[f;a, b] =f(a) +

m

X

i=1

S_{i}

x−a h

h^{i−1}

i! ∆_{h}f^{(i−1)}(a), (1.2)
whereh=b−aand

S_{i}

x−a h

=B_{i}

x−a h

−B_{i}, i≥1, (1.3)

∆hf^{(i−1)}(a) =f^{(i−1)}(b)−f^{(i−1)}(a), 1≤i≤m.

LetX = [a, b]×[c, d]. Denoteh:=b−a, k:=d−c and consider the operators:

∆(h,0)f(x, y) :=f(x+h, y)−f(x, y), (1.4)

∆_{(0,k)}f(x, y) :=f(x, y+k)−f(x, y),

∆_{(h,k)}f(x, y) := ∆_{(h,0)}∆_{(0,k)}f(x, y) = ∆_{(0,k)}∆_{(h,0)}f(x, y).

Forf ∈C^{m,n}(X),the Bernoulli interpolant on the rectangle is [6]:

(B_{m,n}f)(x, y) :=f(a, c) +

m

P

i=1

∆_{(h,0)}f^{(i−1,0)}(a, c)h^{i−1}
i! S_{i}

x−a h

(1.5) +

n

P

j=1

∆_{(0,k)}f^{(0,j−1)}(a, c)k^{j−1}
j! Sj

y−c k

+

m

P

i=1 n

P

j=1

∆_{(h,k)}f^{(i−1,j−1)}(a, c)h^{i−1}k^{j−1}
i!j! S_{i}

x−a h

S_{j}

y−c k

, where Sk, k >1 are given in (1.3). The polynomial from (1.5) satisfies the following interpolation conditions [6]:

(B_{m,n}f)(a, c) =f(a, c), (1.6)

(∆_{(h,0)}B_{m,n}f)^{(i,0)}(a, c) = ∆_{(h,0)}f^{(i,0)}(a, c), 0≤i≤m−1,
(∆_{(0,k)}B_{m,n}f)^{(0,j)}(a, c) = ∆_{(0,k)}f^{(0,j}^{)}(a, c), 0≤j≤n−1,

(∆_{(h,k)}B_{m,n}f)^{(i,j)}(a, c) = ∆_{(h,k)}f^{(i,j)}(a, c), 0≤i≤m−1,0≤j ≤n−1.

The Shepard method, introduced in [12], is a well suited method for multivariate interpolation of very large scattered data sets. The bivariate Shepard operator is given by

(Sf) (x, y) =

N

X

i=0

A_{i,µ}(x, y)f(x_{i}, y_{i}), (1.7)
where

A_{i,µ}(x, y) =

N

Q

j=0 j6=i

r_{j}^{µ}(x, y)

N

P

k=0 N

Q

j=0 j6=k

r_{j}^{µ}(x, y)

, (1.8)

withµ >0 andri(x, y) are the distances between (x, y) and the given points (xi, yi), i= 0,1, ..., N.

In [4] we have introduced the bivariate Shepard-Bernoulli operator that preserve the advantages and improve the reproduction qualities, have better accuracy and computational efficiency:

(SBf)(x, y) =

N

X

i=0

Ai,µ(x, y)(B^{i}_{m,n}f)(x, y), µ >0, (1.9)
whereB_{m,n}^{i} f denotes the Bernoulli interpolantB_{m,n}[f; (x_{i}, y_{i}),(h_{i}, k_{i})] in the rectan-
gle with opposite vertices (xi, yi),(xi+1, yi+1), given by (1.5), havinghi =xi+1−xi,
k_{i}=y_{i+1}−y_{i}, i= 0, ..., N.

Remark 1.1. The operatorSB has the following interpolation properties:

(S_{B}f)(x_{p}, y_{p}) =f(x_{p}, y_{p}), 0≤p≤N;µ > m+n−2
and its degree of exactness is (m, n).

There are flat spots at each data point and the accuracy tends to decrease in the areas where the interpolation nodes are sparse. This can be improved using the local version of Shepard interpolation, introduced by Franke and Nielson in [8] and improved in [7], [10], [11]:

(Sf) (x, y) =

N

P

i=0

Wi(x, y)f(xi, yi)

N

P

i=0

Wi(x, y)

, (1.10)

with

Wi(x, y) =

(Rw−ri)+

R_{w}r_{i}
^{2}

,

whereR_{w}is a radius of influence about the node (x_{i}, y_{i}) and it is varying withi.This
is taken as the distance from node i to the jth closest node to (xi, yi) for j > Nw

(N_{w} is a fixed value) and j as small as possible within the constraint that the jth
closest node is significantly more distant than the (j−1)st closest node (see, e.g. [10]).

As it is mentioned in [9], this modified Shepard method is one of the most powerful software tools for the multivariate interpolation of large scattered data sets.

With these assumptions, for f ∈ C^{(m,n)}(X) and distinct points (x_{i}, y_{i}) ∈ X,
i= 0, ..., N,the Shepard-Bernoulli operator, given in (1.9), becomes (see [4]):

(S_{B}^{w}f) (x, y) :=

N

P

i=0

Wi(x, y) (B_{m,n}^{i} f)(x, y)

N

P

i=0

W_{i}(x, y)

. (1.11)

### 2. Constraints of the Shepard-Bernoulli operator

There are two most important classes of interpolation methods of very large scattered data sets: radial basis functions and Shepard type methods. Both are widely used in practice.

There could be cases when the data are inherently positive. We will make the modified Shepard-Bernoulli operator to preserve positivity by forcing the quadratic basis functions to be positive, using the method introduced in [2].

The modified Shepard-Bernoulli operator preserves the advantages of the classi- cal Shepard operator and improves the reproduction qualities, have better accuracy and computational efficiency. There are cases when we have additional information to take into account in reconstruction by interpolation as the case when the information are the subject to certain physical laws that constrain their behavior. In [2] there are mentioned the case when the information refer to some densities and the case when data values are specified as fractions of a whole. In the first case the reconstruction must be positive and in the second must be within [0,1] to be realistic.

The classical Shepard operatorS,given in (1.7) satisfies the following property:

min{f(x_{i}, y_{i})} ≤(Sf)(x, y)≤max{f(x_{i}, y_{i})}, i= 0, ..., N. (2.1)
A consequence of this property is that a positive interpolant is guaranteed if the data
values are positive.

The modified Shepard operator, given in (1.10), has superior qualities but it does not satisfy the property (2.1).

For a functionf ∈C^{(m,n)}(X), X= [a, b]×[c, d] and a set ofN+1 distinct points
(x_{i}, y_{i})∈X, i= 0, ..., N,we consider Shepard-Bernoulli operator given by (1.9). We
will impose constraints to positivity, using the method from [2].

We will use as a basis function a linear transformation of the old one, namely the function

(C_{B}^{i}f)(x, y) =α(B^{i}_{m,n}f)(x, y) +β, i= 0, ..., N (2.2)
whereαandβ are chosen as

α= f(xi, yi)

f(xi, yi)− min

(x,y)∈[xi,x_{i+1}]×[yi,y_{i+1}]{B_{m,n}^{i} f(x, y)} ∈[0,1]

β= (1−α)f(x_{i}, y_{i}), fori= 1, ...N.

Remark 2.1. B_{m,n}^{i} f, i = 0, ..., N could have negative values but the constrained
functionC_{B}^{i} f, i= 0, ..., N have just positive values.

Theorem 2.2. If(xA, yA)and(xB, yB)are two points such that
(B_{m,n}^{i} f)(xA, yA)≤(B_{m,n}^{i} f)(xB, yB)
then

(C_{B}^{i}f)(x_{A}, y_{A})≤(C_{B}^{i} f)(x_{B}, y_{B}).

Proof. The proof follows directly taking into account the form (2.2).

Remark 2.3. The method can be extended to handle other types of constraints, for example, in the interval [0,1] or, furthermore, in any arbitrary interval [a, b], a > b, a, b∈R.

T he constrained Shepard-Bernoulli operator of first kind is given by
(S_{B}^{c}f)(x, y) =

N

X

i=0

A_{i,µ}(x, y)(C_{B}^{i}f)(x, y), µ >0, (2.3)
withAi,µ(x, y) given in (1.8).

Theorem 2.4. For f ∈ C^{(m,n)}(X) the operator SB has the following interpolation
properties:

(S_{B}^{c}f)(xp, yp) =f(xp, yp),
for0≤p≤N andµ > m+n−2.

Proof. We have

(S_{B}^{c}f)(xp, yp) =

N

X

i=0

Ai,µ(xp, yp)(C_{B}^{i} f)(xp, yp)

=α

N

X

i=0

Ai,µ(xp, yp)(B_{m,n}^{i} f)(xp, yp) +β

N

X

i=0

Ai,µ(xp, yp)

Taking into account that

N

P

i=0

Ai,µ(xp, yp) = 1, we get

(S^{c}_{B}f)(xp, yp) =α

N

X

i=0

Ai,µ(xp, yp)(B_{m,n}^{i} f)(xp, yp) +β

=α(SBf)(x, y)(xp, yp) +β

and by the interpolation properties ofSB (given in [4]) the conclusion follows.

Theorem 2.5. The degree of exactness of the operatorS_{B}^{c} is(m, n).

Proof. The proof follows considering the form of

(C_{B}^{i}f)(x, y) =α(B_{m,n}^{i} f)(x, y) +β

and the property that degree of exactness of the operator SB is (m, n), (as it was

proved in [4]).

We consider also the modified Shepard-Bernoulli operator given by (1.11). We will keep the benefits of the modified Shepard-Bernoulli interpolation and impose constraints, using the previous method (see [2]).

The constrained Shepard-Bernoulli operator of second kind is given by

(S_{B}^{c,w}f) (x, y) :=

N

P

i=0

W_{i}(x, y) (C_{B}^{i}f)(x, y)

N

P

i=0

Wi(x, y)

, (2.4)

Remark 2.6. If f(xi, yi) = 0 for any i then α= 0 and β = f(xi, yi) therefore the interpolants becomes the classical Shepard interpolants.

### 3. Numerical examples

We consider the following test functions ([7], [10], [11]):

f1(x, y) = exp

−81

16((x−0.5)^{2}+ (y−0.5)^{2})

/3, (Gentle) f2(x, y) =p

64−81((x−0.5)^{2}+ (y−0.5)^{2})/9−0.5 (Sphere)

Table 1 contains the maximum errors for approximating by the Shepard,
Shepard-Bernoulli, the modified Shepard-Bernoulli and the coresponding constrained
methods, (2.3) and (2.4), considering 52 random generated nodes in the unit square,
m=n= 2 andµ= 2. By numerical experiments we have obtained that for these data
the optimal value forN_{w}isN_{w}= 8. We compare the obtained numerical results with
some combined Shepard operators known in the literature, namely with the combined
Shepard operators of Lagrange, Hermite and Taylor type, denoted respectively byS_{L},
SH andST.

In Figures 1 and 2 we plot the graphs offi, SBfi, S_{B}^{w}fi, S_{B}^{c}fi, S_{B}^{c,w}fi,fori= 1,2.

Table 1.Maximum approximation errors.

f_{1} f_{2}
Sf 0.1870 0.2374
SBf 0.0905 0.0274
S_{B}^{w}f 0.0628 0.0187
S_{B}^{c}f 0.2975 0.3563
S_{B}^{c,w}f 0.3085 0.3767
S_{L}f 0.5353 0.5798
S_{H}f 0.4646 0.0883
STf 0.2110 0.3635

0 0.2

0.4 0.6

0.8 1

0 0.2 0.4 0.6 0.8 1 0 0.05 0.1 0.15 0.2 0.25 0.3 0.35

f1

0 0.2 0.4 0.6 0.8 1

0 0.2 0.4 0.6 0.8 1 0.08

0.1 0.12 0.14 0.16 0.18 0.2 0.22

Sf1

0 0.2

0.4 0.6

0.8 1

0 0.2 0.4 0.6 0.8 1 0 0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4

SBf1

0 0.2

0.4 0.6

0.8 1

0 0.2 0.4 0.6 0.8 1

−0.05 0 0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4

S_{B}^{w}f_{1}

0 1 0.1

1 0.2

0.8 0.3

0.5 0.6

0.4

0.4 0.2

0 0

S^{c}_{B}f1

0 1 0.1

1 0.2

0.8 0.3

0.5 0.6

0.4

0.4 0.2

0 0

S_{B}^{c,w}f1

Figure 1. Graphs for the functionf1.

0 0.2

0.4 0.6

0.8 1

0 0.2 0.4 0.6 0.8 1 0 0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4

f2

0 0.2 0.4 0.6 0.8 1

0 0.2 0.4 0.6 0.8 1 0.2 0.22 0.24 0.26 0.28 0.3 0.32 0.34

Sf2

0 0.2

0.4 0.6

0.8 1

0 0.2 0.4 0.6 0.8 1 0 0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4

SBf2

0 0.2

0.4 0.6

0.8 1

0 0.2 0.4 0.6 0.8 1 0 0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4

S_{B}^{w}f2

0 1 0.1

1 0.2

0.8 0.3

0.5 0.6

0.4

0.4 0.2

0 0

S_{B}^{c}f_{2}

0 1 0.1

1 0.2

0.8 0.3

0.5 0.6

0.4

0.4 0.2

0 0

S_{B}^{c,w}f_{2}
Figure 2. Graphs for the functionf2.

### References

[1] Asim, M.R., Mustafa, G., Brodlie, K.W.,Constrained Visualization of 2D Positive Data using Modified Quadratic Shepard Method, WSCG, 2004.

[2] Brodlie, K.W., Asim, M.R., Unsworth, K.,Constrained Visualization Using the Shepard Interpolation Family,Computer Graphics Forum,24(2005), 809-820.

[3] Caira, R., Dell’Accio, F.,Shepard-Bernoulli operators,Math. Comp.,76(2007), 299-321.

[4] C˘atina¸s, T., The bivariate Shepard operator of Bernoulli type,Calcolo, 44(2007), 189- 202.

[5] Coman, Gh.,The remainder of certain Shepard type interpolation formulas, Stud. Univ.

Babe¸s-Bolyai Math.,32(1987), no. 4, 24-32.

[6] Costabile, F.A., Dell’Accio, F., Expansion Over a Rectangle of Real Functions in Bernoulli Polynomials and Applications, BIT,41(2001), no. 3, 451-464.

[7] Franke, R.,Scattered data interpolation: tests of some methods,Math. Comp.,38(1982), 181-200.

[8] Franke, R., Nielson, G., Smooth interpolation of large sets of scattered data, Int. J.

Numer. Meths. Engrg.,15(1980), 1691-1704.

[9] Lazzaro, D., Montefusco, L.B.,Radial basis functions for multivariate interpolation of large scattered data sets,J. Comput. Appl. Math.,140(2002), 521-536.

[10] Renka, R.J.,Multivariate interpolation of large sets of scattered data,ACM Trans. Math.

Software,14(1988), 139-148.

[11] Renka, R.J., Cline, A.K.,A triangle-basedC^{1} interpolation method,Rocky Mountain J.

Math.,14(1984), 223-237.

[12] Shepard, D.,A two dimensional interpolation function for irregularly spaced data,Proc.

23rd Nat. Conf. ACM, 1968, 517-523.

Teodora C˘atina¸s

“Babe¸s-Bolyai” University

Faculty of Mathematics and Computer Sciences 1, Kog˘alniceanu St.

400084 Cluj-Napoca, Romania e-mail:[email protected]