Ismail Khalfaoui-Hassani
email: <ismail.khalfaoui-hassani@univ-tlse3.fr>.

Thomas Pellegrini

Timothée Masquelier

Dilated Convolution with Learnable Spacings: beyond bilinear interpolation

Keywords:

Dilated Convolution with Learnable Spacings (DCLS) is a recently proposed variation of the dilated convolution in which the spacings between the non-zero elements in the kernel, or equivalently their positions, are learnable. Non-integer positions are handled via interpolation. Thanks to this trick, positions have well-defined gradients. The original DCLS used bilinear interpolation, and thus only considered the four nearest pixels. Yet here we show that longer range interpolations, and in particular a Gaussian interpolation, allow improving performance on ImageNet1k classification on two state-of-the-art convolutional architectures (ConvNeXt and ConvFormer), without increasing the number of parameters. The method code is based on PyTorch and is available at github.com/K-H-Ismail/Dilated-Convolutionwith-

Introduction

Dilated Convolution with Learnable Spacings (DCLS) is an innovative convolutional method whose effectiveness in computer vision was recently demonstrated [START_REF] Khalfaoui-Hassani | Dilated convolution with learnable spacings[END_REF]. In DCLS, the positions of the non-zero elements within the convolutional kernels are learned in a gradient-based manner. The challenge of nondifferentiability caused by the integer nature of the positions is addressed through the application of bilinear interpolation. By doing so, DCLS enables the construction of a differentiable convolutional kernel.

DCLS is a differentiable method that only constructs the convolutional kernel. To implement the whole convolution, one can utilize either the native convolution provided by PyTorch or a more efficient implementation such as the "depthwise implicit gemm" convolution method proposed by [START_REF] Ding | Scaling up your kernels to 31x31: Revisiting large kernel design in CNNs[END_REF], which is suitable for large kernels.

The primary motivation behind the development of DCLS was to investigate the potential for enhancing the fixed grid structure imposed by standard dilated convolution in an input-independent way. By allowing an arbitrary number of kernel elements, DCLS introduces a free tunable hyperparameter called the "kernel count". Additionally, the "dilated kernel size" refers to the maximum extent to which the kernel elements are permitted to move within the dilated kernel (Fig. 1c). Both of these parameters can be adjusted to optimize the performance of DCLS. The positions of the kernel elements in DCLS are initially randomized and subsequently allowed to evolve within the limits of the dilated kernel size during the learning process. The main focus of this paper will be to question the choice of bilinear interpolation used by default in DCLS. We tested several interpolations and found in particular that a Gaussian interpolation with learnable standard deviations made the approach more effective.

To evaluate the effectiveness of DCLS with Gaussian interpolation, we integrate it as a drop-in replacement for the standard depthwise separable convolution in two state-ofthe-art convolutional models: the ConvNext-T model [START_REF] Liu | A convnet for the 2020s[END_REF] and the ConvFormer-S18 model [START_REF] Yu | Metaformer baselines for vision[END_REF]. In Section 5, we evaluate the training loss and the classification accuracy of these models on the ImageNet1k dataset [START_REF] Deng | Imagenet: A large-scale hierarchical image database[END_REF]. The remainder of this paper will present a detailed analysis of the methods, equations, algorithms and techniques regarding the application of the Gaussian interpolation in DCLS.

Related work

In the field of convolutional neural networks (CNNs), various approaches have been explored to improve the performance and efficiency of convolutional operations. Gaussian mixture convolutional networks have investigated the fit of input channels with Gaussian mixtures [START_REF] Celarek | Gaussian mixture convolution networks[END_REF], while [START_REF] Chen | Gaussian mask convolution for convolutional neural networks[END_REF] utilized Gaussian masks in their work. Additionally, continuous kernel convolution was studied in the context of image processing by [START_REF] Kim | Self-moving point representations for continuous convolution[END_REF]. Their approach is similar to the linear correlation introduced in [START_REF] Thomas | Flexible and deformable convolution for point clouds[END_REF]. The interpolation function used in the last two works corresponds to the DCLS-Triangle method described in 3.1. Romero et al. have also made notable contributions in learning continuous functions that map the positions to the weights (Romero et al., 2022a;b).

In the work by [START_REF] Jacobsen | Structured receptive fields in cnns[END_REF], the kernel is represented as a weighted sum of basis functions, including centered Gaussian filters and their derivatives. [START_REF] Pintea | Resolution learning in deep convolutional networks using scale-space theory[END_REF] extended this approach by incorporating the learning of Gaussian width, effectively optimizing the resolution. [START_REF] Shelhamer | Blurring the line between structure and learning to optimize and adapt receptive fields[END_REF] introduced a kernel factorization method where the kernel is expressed as a composition of a standard kernel and a structured Gaussian one. In these last three works the Gaussians are centered on the kernel.

Furthermore, the utilization of bilinear interpolation within deformable convolution modules has already shown its effectiveness. [START_REF] Dai | Deformable convolutional networks[END_REF], [START_REF] Qi | Deformable convolutional networks-coco detection and segmentation challenge 2017 entry[END_REF] and recently [START_REF] Wang | Internimage: Exploring large-scale vision foundation models with deformable convolutions[END_REF] leveraged bilinear interpolation to smoothen the non-differentiable regular-grid offsets in the deformable convolution method. Even more recently, in Kim et al. (2023), a Gaussian attention bias with learnable standard deviations has been successfully used in the positional embedding of the attention module of the ViT model [START_REF] Dosovitskiy | An image is worth 16x16 words: Transformers for image recognition at scale[END_REF] and leads to reasonable gains on ImageNet1k.

Methods

From bilinear to Gaussian interpolation

We denote by m ∈ N * the number of kernel elements inside the dilated constructed kernel and we refer to it as the "kernel count". Moreover, we denote respectively by s x , s y ∈ N * × N * , the sizes of the constructed kernel along the x-axis and the y-axis. The latter could be seen as the limits of the dilated kernel, and we refer to them as the "dilated kernel size".

The s x × s y matrix space over R is defined as the set of all s x × s y matrices over R, and is denoted M sx,sy (R). The real numbers w, p x , σ x , p y and σ y respectively stand for the weight, the mean position and standard deviation of that weight along the x-axis (width) and its mean position and standard deviation along the y-axis (height).

The mathematical construction of the 2D-DCLS kernel in [START_REF] Khalfaoui-Hassani | Dilated convolution with learnable spacings[END_REF] relies on bilinear interpolation and is described as follows :

f : R × R × R → M sx,sy (R) w, p x , p y → K (1)
where ∀i ∈ 1 .. s x , ∀j ∈ 1 .. s y :

K ij =            w (1 -r x ) (1 -r y ) if i = ⌊p x ⌋, j = ⌊p y ⌋ w r x (1 -r y ) if i = ⌊p x ⌋ + 1, j = ⌊p y ⌋ w (1 -r x ) r y if i = ⌊p x ⌋, j = ⌊p y ⌋ + 1 w r x r y if i = ⌊p x ⌋+1, j = ⌊p y ⌋+1 0 otherwise
(2) and where the fractional parts are:

r x = {p x } = p x -⌊p x ⌋ and r y = {p y } = p y -⌊p y ⌋ (3)
An equivalent way of describing the constructed kernel K in Equation 2 is:

K ij = w • g(p x -i) • g(p y -j) (4) with g : x → max(0, 1 -|x|) (5)
This expression corresponds to the bilinear interpolation as described in Dai et al. (2017, eq. 4).

In fact, this last g function is known as the triangle function (refer to Fig. 2 for a graphic representation), and is widely used in kernel density estimation. From now on, we will note it as

∀x ∈ R Λ(x) def = max(0, 1 -|x|) (6)
First, we consider a scaling by a parameter σ ∈ R + for the triangle function (the bilinear interpolation corresponds to σ = 1),

∀x ∈ R, ∀σ ∈ R + Λ σ (x) def = max(0, σ -|x|) (7)
We found that this scaling parameter σ could be learned by backpropagation and that doing so increases the performance of the DCLS method. As we have different σ parameters for the x and y-axes in 2D-DCLS, learning the standard deviations costs two additional learnable parameters and two additional FLOPs (multiplied by the number of the channels of the kernel and the kernel count). We refer to the DCLS method with triangle function interpolation as the DCLS-Triangle method.

Second, we tried a smoother function rather than the piecewise affine triangle function, namely the Gaussian function:

∀x ∈ R, ∀σ ∈ R * , G σ (x) def = exp - x 2 2σ 2 (8)
We refer to the DCLS method with Gaussian interpolation as the DCLS-Gauss method. In practice, instead of Equations 7 and 8, we respectively use:

∀x ∈ R, ∀σ ∈ R, Λ σ0+σ (x) = max(0, σ 0 + |σ| -|x|) (9) ∀x ∈ R, ∀σ ∈ R, G σ0+σ (x) = exp - 1 2 x 2 (σ 0 + |σ|) 2
(10) with σ 0 ∈ R * + a constant that determines the minimum standard deviation that the interpolation could reach. For the triangle interpolation, we take σ 0 = 1 in order to have at least 4 adjacent interpolation values (see Figure 1c). And for the Gaussian interpolation, we set σ 0 = 0.27.

Last, to make the sum of the interpolation over the dilated kernel size equal to 1, we divide the interpolations by the following normalization term :

A = ϵ + sx i=1 sy j=1 I σ0+σ x (p x -i) • I σ0+σ y (p y -j) (11)
with I an interpolation function (Λ or G in our case) and ϵ = 1e -7 for example, to avoid division by zero.

Other interpolations Based on our tests, other functions such as Lorentz, hyper-Gaussians and sinc functions have been tested with no great success. In addition, learning a correlation parameter ρ ∈ [-1, 1] or equivalently a rotation parameter θ ∈ [0, 2π] as in the bivariate normal distribution density did not improve performance (maybe because cardinal orientations predominate in natural images).

The 2D-DCLS-Gauss kernel construction algorithm

In the following, we describe with pseudocode the kernel construction used in 2D-DCLS-Gauss and 2D-DCLS-Triangle. I is the interpolation function (Λ or G in our case) and ϵ = 1e -7. In practice, w, p x , p y , σ x and σ y are 3-D tensors of size (channels out, channels in // groups, K count), but the algorithm presented here is easily extended to this case by applying it channel-wise.

Algorithm 1 2D-DCLS-interpolation kernel construction Require: w, p x , p y , σ x , σ y : vectors of dimension m Ensure: K : the constructed kernel, of size (s x × s y )

1: K ← 0 sx,sy {zero tensor of size s x , s y } 2: for k = 0 to m -1 do 3:

H ← 0 sx,sy 4:

p x k ← p x k + s x //2; p y k ← p y k + s y //2 5: σ x k ← |σ x k | + σ I 0 ; σ y k ← |σ y k | + σ I 0 6:
for i = 0 to s x -1 do 7:

for j = 0 to s y -1 do 8:

H[i, j] ← I σ x k (p x k -i) * I σ y k (p y k -j) 9:
end for 10:

end for 11:

H[:, :] ← H[:, :] /(ϵ + s x -1 i=0 s y -1 j=0 H[i, j]) 12: K ← K + H * w k
13: end for

Learning techniques

Having discussed the implementation of the interpolation in the DCLS method, we now shift our focus to the techniques employed to maximize its potential. We retained most of the techniques used in [START_REF] Khalfaoui-Hassani | Dilated convolution with learnable spacings[END_REF], and suggest new ones for learning standard deviations parameters. In Appendix C, we present the training techniques that have been selected based on consistent empirical evidence, yielding improved training loss and validation accuracy.

Results

We took two recent state-of-the-art convolutional architectures, ConvNeXt and ConvFormer, and drop-in replaced all the depthwise convolutions by DCLS ones, using the three different interpolations (bilinear, triangle or Gauss). Table 1 reports the results in terms of training loss and validation accuracy.

A first observation is that all the DCLS models perform much better than the baselines, whereas they have the same number of parameters. There are also subtle differences between interpolation functions. As Figure 3 shows, triangle and bilinear interpolations perform similarly, but the Gaussian interpolation performs significantly better.

Furthermore, the advantage of the Gaussian interpolation w.r.t. bilinear is not only due to the use of a larger kernel, as a 17x17 Gaussian kernel (5th line in Table 1) still outperforms the bilinear case (2nd line). Finally, the 6th line in Table 1 shows that there is still room for improvement by increasing the kernel count, although this slightly increases the number of trainable parameters w.r.t. the baseline.

Conclusion

In conclusion, this study introduces Gaussian and Λ interpolation methods as alternatives to bilinear interpolation in Dilated Convolution with Learnable Spacings (DCLS).

Evaluations on state-of-the-art convolutional architectures demonstrate that Gaussian interpolation improves performance of image classification task on ImageNet1k without increasing parameters. Future work could implement the Whittaker-Shannon interpolation instead of the Gaussian interpolation and search for a dedicated architecture, that will make the most of DCLS.

Figure 1 .

 1 Figure 1. (a) a standard 3 × 3 kernel. (b) a standard dilated 3 × 3 kernel. (c) a 2D-DCLS kernel using bilinear interpolation with 9 kernel elements and a kernel size of 9. (d) the same kernel as (c) with Gaussian interpolation. The numbers have been rounded in all figures and omitted in (d) for readability.

Figure 2 .

 2 Figure 2. 1D view of Gaussian and Λ functions with σ = 5.

  IDX is None or self.lim is None: 13 J = Parameter(torch.arange(0, self.dilated_kernel_size[0]), 14 requires_grad=False).to(device) 15 I = Parameter(torch.arange(0, self.dilated_kernel_size[1]), 16 requires_grad=False).to(device) 17 I = I.expand(self.dilated_kernel_size[0],-1) 18 J = J.expand(self.dilated_kernel_size[1],-1).t() 19 IDX = torch.cat((I.unsqueeze(0),J.unsqueeze(0)), 0) 20 IDX = IDX.expand(self.out_channels, self.in_channels//self.groups, .tensor(self.dilated_kernel_size).to(device) 24 self.lim = lim.expand(self.out_channels, 25 self.in_channels//self.groups, self.kernel_count, -1).permute(3

Table 1 .

 1 Classification accuracy on the validation set and training loss on ImageNet-1K. For the 17/34 bilinear, the 23/26 Triangle and Gaussian cases, the results have been averaged over 3 distinct seeds (the corresponding lines are highlighted in yellow).

	MODEL @ 224	KER. SIZE / COUNT	INTERPOLATION	# PARAM. TRAIN LOSS	TOP-5 ACC. TOP-1 ACC.
	CONVNEXT-T ▨	7 2 / 49		28.59M	2.828	96.05	82.08
	CONVNEXT-T ▩	17 2 / 34	BILINEAR	28.59M	2, 775	96.11	82.44
	CONVNEXT-T ⊙	23 2 / 26	TRIANGLE	28.59M	2.787	96.09	82.34
	CONVNEXT-T ⋆	23 2 / 26	GAUSSIAN	28.59M	2.762	96.18	82.44
	CONVNEXT-T	17 2 / 26	GAUSSIAN	28.59M	2.773	96.17	82.40
	CONVNEXT-T	23 2 / 34	GAUSSIAN	28.69M	2.758	96.22	82.60
	CONVFORMER-S18 ▨	7 2 / 49		26.77M	2.807	96.17	82.84
	CONVFORMER-S18 ▩	17 2 / 40	BILINEAR	26.76M	2.764	96.42	83.14
	CONVFORMER-S18 ⊙	23 2 / 26	TRIANGLE	26.76M	2.761	96.38	83.09
	CONVFORMER-S18 ⋆	23 2 / 26	GAUSSIAN	26.76M	2.747	96.31	82.99

Acknowledgments

This work was performed using HPC resources from GENCI-IDRIS (Grant 2021-[AD011013219]). Support from the ANR-3IA Artificial and Natural Intelligence Toulouse Institute is gratefully acknowledged. We would also like to thank the region of Toulouse Occitanie.

A. Code and reproducibility

The code of the method is based on PyTorch and available at https://github.com/K-H-Ismail/Dilated-Convolution-with-Learnable-Spacings-PyTorch. • Weight decay: No weight decay was used for positions. We apply the same for standard deviation parameters.

• Positions and standard deviations initialization: position parameters were initialized following a centered normal law of standard deviation 0.5. Standard deviation parameters were initialized to a constant 0.23 in DCLS-Gauss and to 0 in DCLS-Triangle in order to have a similar initialisation to DCLS with bilinear interpolation at the beginning.

• Positions clamping : Previously in DCLS, kernel elements that reach the dilated kernel size limit were clamped. It turns out that this operation is no longer necessary with the Gauss and Λ interpolations.

• Dilated kernel size tuning: When utilizing bilinear interpolation in ConvNeXt-dcls, a dilated kernel size of 17 was found to be optimal, as larger sizes did not yield improved accuracy. However, with Gaussian and Λ interpolations, there appears to be no strict limit to the dilated kernel size. Accuracy tends to increase logarithmically as the size grows, with improvements observed up to kernel sizes of 51. It is important to note that increasing the dilated kernel size does not impact the number of trainable parameters, but it does affect throughput. Therefore, a compromise between accuracy and throughput was achieved by setting the dilated kernel size to 23.

• Kernel count tuning: This hyper-parameter has been configured to the maximum integer value while still remaining below the baselines to which we compare ourselves in terms of trainable parameters. It is worth noting that each additional element in the 2D-DCLS-Gauss or 2D-DCLS-Triangle methods introduces five more learnable parameters: weight, vertical and horizontal position, and their respective standard deviations. To maintain simplicity, the same kernel count was applied across all model layers.

• Learning rate scaling: To maintain consistency between positions and standard deviations, we applied the same learning rate scaling ratio of 5 to both. In contrast, the learning rate for weights remained unchanged.

• Synchronizing positions: we shared the kernel positions and standard deviations across convolution layers with the same number of parameters, without sharing the weights. Parameters in these stages were centralized in common parameters that accumulate the gradients.

D. 1D and 3D convolution cases

For the 3D case, Equation 4 can be generalized as a product along spatial dimensions. We denote respectively by s x , s y , s z ∈ N * × N * × N * , the sizes of the constructed kernel along the x-axis, the y-axis and the z-axis. The constructed kernel tensor K 3D ∈ M sx,sy,sz (R) is therefore:

∀i ∈ 1 .. s x , ∀j ∈ 1 .. s y , ∀k ∈ 1 .. s z :

with I an interpolation function (Λ or G), σ 0 = 1 for the Λ interpolation and σ 0 = 0.27 for the Gaussian one. w, p x , σ x , p y , σ y , p z and σ z respectively representing the weight, the mean position and standard deviation of that weight along the x-axis (width), the mean position and standard deviation along the y-axis (height) and its mean position and standard deviation along the z-axis (depth).

The constructed kernel vector K 1D ∈ R sx for the 1D case is simply:

∀i ∈ 1 .. s x :

The Algorithm 1 as well as the Pytorch code B are readily adapted to these cases by following the above note.


domain: info math

text:
Carlotta Giannelli
email: carlotta.giannelli@unifi.it

Sofia Imperatore
email: sofia.imperatore@unifi.it

Angelos Mantzaflaris
email: angelos.mantzaflaris@inria.fr

Felix Scholz
email: felix.scholz@jku.at

BIDGCN: Boundary informed dynamic graph convolutional network for adaptive spline fitting of scattered data

Keywords: Geometric deep learning, Parameterization, Surface fitting, Adaptive spline approximation, THB-splines

Introduction

Geometric deep learning refers to the generalization of convolutional neural networks to non-Euclidean data such as discrete manifolds, graphs and general point clouds [START_REF] Bronstein | Geometric deep learning: Going beyond Euclidean data[END_REF]. While for data represented on a regular grid in a Euclidean space the convolution operator is uniquely defined for arbitrary dimensions, in the case of unstructured non-Euclidean data many different approaches for the definition of operators that mimic the behavior of standard convolution have been proposed [START_REF] Wu | A comprehensive survey on graph neural networks[END_REF].

Most current graph neural network architectures are classified as message passing neural networks [START_REF] Gilmer | Neural message passing for quantum chemistry[END_REF], meaning that hidden states at each vertex are computed by aggregating the output of message functions applied to the features of the vertex and its neighbors, the adjacent vertices directly connected by an edge.

While graph neural networks have been successfully applied to classification tasks such as shape recognition, segmentation, and registration [START_REF] Welling | Semi-supervised classification with graph convolutional networks[END_REF][START_REF] Qi | Pointnet: Deep learning on point sets for 3d classification and segmentation[END_REF][START_REF] Hanocka | Meshcnn: a network with an edge[END_REF][START_REF] Sharp | DiffusionNet: Discretization agnostic learning on surfaces[END_REF], regression tasks where the network 1 should predict a (potentially vector-valued) function on the input geometry are much less studied [START_REF] Berrone | Graphinformed neural networks for regressions on graph-structured data[END_REF][START_REF] Fan | Arbitrary surface data patching method based on geometric convolutional neural network[END_REF]. Moreover, graph neural networks that are used for processing discrete surface data have mostly been applied to closed surfaces, i.e. surfaces without boundary. This simplifies the design of suitable graph convolution operators significantly, since all vertices of the discrete manifold or point cloud can be handled in the same way and no explicit distinction between interior and boundary vertices needs to be made. However, in many applications in geometric modelling, geometry processing and numerical analysis, the input geometries are given as discrete surfaces with boundary. For a variety of problems boundary conditions are imposed on the corresponding vertices and often the solution to a problem is only uniquely determined up to these boundary conditions. For example this is the case for boundary value problems of elliptic partial differential equations on surfaces, as well as for the problem of scattered point cloud parameterization. In order to apply deep learning based methods to this kind of problem it is then necessary to devise a network architecture that takes into account the boundary conditions in addition to the standard vertex features of the discrete surface. Since boundary conditions can be regarded as additional features that are defined only on the boundary vertices but not on the interior vertices, this means that such a network architecture needs to be able to handle data with varying feature dimensions.

In this paper, we propose a new Boundary Informed Dynamic Graph Convolutional Network (BIDGCN) for processing scattered point clouds based on the dynamic edge convolution introduced in [START_REF] Wang | Dynamic graph CNN for learning on point clouds[END_REF]. The distinctive feature of this dynamic edge convolution is that the graph used for message passing is recomputed after each layer according to the new features, thereby enabling the network to predict the graph and to transport the features arbitrarily far in the point cloud. The key idea of the new BIDGCN is that we regard the boundary conditions as additional features at the boundary vertices of the point cloud. These features are propagated into the whole point cloud by a novel graph convolution operator that contains two separate trainable message functions: the first one for edges between interior and boundary vertices, and the second one for edges between interior vertices. In the subsequent hidden layers, the information stemming from the boundary conditions is further processed and used to predict the solution for the problem at hand. To summarize, the main advantages of the proposed network layer are: (i) the ability to incorporate boundary conditions in its prediction, (ii) the dynamic prediction of the graph used for message passing. The second property is directly inherited from the dynamic edge convolution approach [START_REF] Wang | Dynamic graph CNN for learning on point clouds[END_REF].

As a use case for our new network architecture, we consider its application to the problem of scattered point cloud parameterization. The parameterization problem is a key step in different geometric modeling tasks and its solution is usually the starting point of various data processing algorithms. Among others, it is essential for approximating a scattered three-dimensional point cloud with high-order geometries, such as B-splines, NURBS or locally refined (hierarchical) B-splines. Industrial applications include, for example, a variety of reverse engineering problems, as the adaptive re-construction of geometric models from a 3D scan [START_REF] Kiss | Adaptive CAD model (re-)construction with THB-splines[END_REF][START_REF] Bracco | Adaptive fitting with THB-splines: Error analysis and industrial applications[END_REF][START_REF] Bracco | THB-Spline Approximations for Turbine Blade Design with Local B-Spline Approximations[END_REF]. A class of methods for scattered point cloud parameterization was proposed in [START_REF] Floater | Meshless parameterization and surface reconstruction[END_REF] building upon the methods from [START_REF] Floater | Parametrization and smooth approximation of surface triangulations[END_REF] for the related problem of parameterizing triangulated surface data. In these methods, the one-dimensional boundary of the point cloud is parameterized separately, which is a much easier problem. The parameterization of the interior vertices is then found by assuming that the parameter for each interior vertex is a convex combination of the parameters of a certain vertex neighborhood. The specific choices of neighbors and coefficient values considered in the convex combinations determine the particular method from this class, see [START_REF] Floater | Meshless parameterization and surface reconstruction[END_REF]. Recently, machine learning tools have been applied to the problem of determining these coefficients [START_REF] Yavuz | A dynamic neural network model for accelerating preliminary parameterization of 3D triangular mesh surfaces[END_REF][START_REF] Giannelli | Learning meshless parameterization with graph convolutional neural networks[END_REF]. In particular, [START_REF] Yavuz | A dynamic neural network model for accelerating preliminary parameterization of 3D triangular mesh surfaces[END_REF] builds upon [START_REF] Floater | Parametrization and smooth approximation of surface triangulations[END_REF] and trains a network for a fixed triangular mesh, whereas [START_REF] Giannelli | Learning meshless parameterization with graph convolutional neural networks[END_REF] provides a more general meshless data-driven approach based on [START_REF] Floater | Meshless parameterization and surface reconstruction[END_REF]. While these methods lead to an improved performance compared to the standard choices of coefficients, they also inherit their limitations in terms of efficiency, robustness, and sensitivity to parameter-dependent graph connectivity.

In our approach, we assume that the boundary is already parameterized, as in the standard methods. We then use a network architecture based on the new BIDGCN model to predict the parameterization of the interior vertices. This leads to three key advantages compared to the classical methods for scattered point cloud parameterization: In order to achieve an effective scattered data approximation scheme, we suitably exploit the BIDGCN parameterization as starting point for adaptive fitting with truncated hierarchical Bsplines (THB-splines) [START_REF] Giannelli | THBsplines: The truncated basis for hierarchical splines[END_REF]. The performance of the approximation scheme are then enhanced by introducing a parameter correction step [START_REF] Hoschek | Intrinsic parametrization for approximation[END_REF] within each iteration of the adaptive strategy. A similar idea has been investigated in the case of structured data fitting with T-splines in [START_REF] Wang | Curvature-guided adaptive T-spline surface fitting[END_REF][START_REF] Shang | Effective re-parameterization and GA based knot structure optimization for high quality t-spline surface fitting[END_REF] and for template mapping with adaptive splines in [START_REF] Sajavičius | Template Mapping Using Adaptive Splines and Optimization of the Parameterization[END_REF].

The structure of the paper is as follows. Section 2 introduces the new boundary informed edge convolution. In Section 3, we summarize the point cloud parameterization problem, that is then addressed in Section 4 by developing a boundary informed dynamic graph convolutional network devoted to scattered data parameterization. Section 5 presents the adaptive fitting algorithm with moving parameters. A selection of numerical examples on synthetic and real data sets is presented in Section 6. In particular, the generalization capabilities of the proposed parameterization model are highlighted within the framework of an adaptive surface reconstruction scheme based on THB-splines. Finally, Section 7 concludes the paper.

Boundary informed edge convolution

In this section, after briefly describing the classic dynamic edge convolution operator, we present our new boundary informed dynamic convolutional network architecture.

Dynamic edge convolution

We assume to be given a point cloud P together with vertex features x i ∈ R n , where n ∈ N is the input feature dimension. The dynamic edge convolution operator defined in [START_REF] Wang | Dynamic graph CNN for learning on point clouds[END_REF] computes new features y i ∈ R m for all points in P, where m ∈ N is the output feature dimension. To this end, first the k-nearest neighbor graph G k is computed based on the input features. This is a directed graph where the existence of a directed edge (j, i) implies that j is among the k nearest neighbors of i with respect to the Euclidean distance

∥x i -x j ∥
of the input features. In the next step, edge features are computed for each directed edge in G k as follows. For all (j, i)

∈ G k e ji = h Θ (x i , x j -x i ),
is evaluated, where h Θ is a feed-forward neural network h Θ : R 2n → R m with trainable weights Θ. Finally, at each vertex of G k , the edge contributions are aggregated as

x ′ i = 2 (j,i)∈G k e ji ,
where 2 can be the sum, the mean, or the maximum value. In a dynamic edge convolution network, the k-nearest neighbor graph is recomputed after each layer with respect to the new features x ′ i , thereby allowing the network to pass information arbitrarily fast across the point cloud. The dynamic approach enables the network to automatically predict a suitable graph for message passing instead of using a fixed one. We refer to these graph neural networks characterized by dynamic edge convolution operators as DGCNN (Dynamic Graph Convolutional Neural Networks).

Finally, in the examples presented in this paper, we replace the k-nearest neighbor graph by the radius graph, meaning that instead of specifying the number of neighbors k, we specify a radius r > 0 and add directed edges (i, j) and (j, i) for all i, j ∈ V such that

∥x i -x j ∥ ≤ r.
By dynamically recomputing the radius graph instead of the k-nearest neighbor graph, we reduce the dependence of the network architecture on the local density of the input point clouds. In particular, the radius graph results in neighborhoods with a fixed maximum distance, while a k-nearest graph might associate points with very large distances if the point cloud is locally sparse.

BIDGCN: Boundary Informed Dynamic Graph Convolutional Network

Based on the dynamic edge convolution, we introduce a new boundary informed input layer that takes as input the point cloud with its vertex features as well as the boundary conditions. It then propagates the boundary conditions into the new features of the interior points. The output of this layer consists of all interior points together with new vertex features. We name the resulting neural network architecture Boundary Informed Dynamic Graph Convolutional Network (BIDGCN).

More precisely, we assume that the input point cloud is decomposed as P = P B ∪ P I into boundary and interior points. Moreover, we assume that each vertex i ∈ P I comes with features x i ∈ R n and every vertex j ∈ P B comes with features (x j , u j ) ∈ R n × R d , were we regard u j ∈ R d as boundary conditions.

We first compute two different radius graphs

G I→I and G B→I ,
where G I→I contains all directed edges (j, i) with i, j ∈ P I :

∥x i -x j ∥ ≤ r, while G B→I contains all directed edges (k, i) with i ∈ P I , k ∈ P B : ∥x i -x k ∥ ≤ r.
This means that even in G B→I , the edges do not depend on the boundary conditions u j . Note that vertices in P B only have outgoing edges and no incoming ones. We then compute edge contributions for all edges using two separate neural networks. For all (j, i) ∈ G I→I , we compute

e ji = h Θ (x i , x j -x i ),
where h is a feed-forward neural network h Θ : R 2n → R m with output feature size m and learnable weights Θ.

For edges (k, i) ∈ G B→I , the boundary conditions at the boundary vertices are concatenated with the vertex features and we compute

e ik = g Φ (x i , x k -x i , u k ),
where g is another feed-forward neural network g Φ : R 2n+d → R m with the same output feature size and independent learnable weights Φ.

Finally, the edge contributions are aggregated in the target vertices of the directed edges. By construction, all target vertices are contained in P I . For i ∈ P I , we have

x ′ i = 2 j:(j,i)∈GI→I e ji 2 2 k:(k,i)∈GB→I e ki ,
where the aggregation operator 2 can be the sum, the mean or the component-wise maximum value. The output of the layer consists of features of dimension m for the interior point cloud P I .

During training Θ and Φ are optimized simultaneously. The flow of the input layer is depicted in Fig. 1 and its application to an interior vertex is illustrated in Fig. 2.

Hidden and output layers

The input layer gives as output new features in R m for the interior point cloud P I .

i (x k0 , u k0 ) (x k1 , u k1 ) x k2 x k3 x k4 g Φ (x i , x k0 -x i , u k ) g Φ (x i , x k1 -x i , u k ) h Θ (x i , x k2 -x i ) h Θ (x i , x k3 -x i ) h Θ (x i , x k4 -x i )
Fig. 2: Edge contributions of the input layer for x i (blue) computed from the features at the neighboring boundary vertices (red) and interior vertices (green).

x ′ i ∈ R m . As a slight modification, while in [START_REF] Wang | Dynamic graph CNN for learning on point clouds[END_REF] the k nearest-neighbor graph was used, we propose to use the radius graph also in the hidden layers. This makes our network architecture more independent of the sampling density of the point cloud. Since the number of neighbors of the vertices is not constant, we use the mean value as the aggregation of the edge features, also motivated by the goal of achieving independence of the sampling density.

Trough the application of the different layers, the information that was transported by the input layer from the boundary vertices to their neighbors in G B→I is propagated further into the interior of the point cloud. Since the radius graph is recomputed after each layer, the information can travel arbitrarily far in each layer, as it is determined by the training process.

Scattered point cloud parameterization for surface fitting

In this section, we introduce the problem of scattered point cloud parameterization for approximation with smooth surfaces. Moreover, we shortly summarize the standard methods for this problem. In the following section, we will then report on how we use a network architecture based on BIDGCN to tackle this problem.

Problem formulation

For a given point cloud P such that each point i ∈ P is equipped with a feature x i ∈ R 3 , a parameterization is defined as a mapping

u : P → Ω,
where Ω ⊂ R 2 is called the parameter domain. We write u i = u(i) for the parameter of i ∈ P.

Finding a suitable parameterization of a point cloud is an important step in different geometry processing tasks. In particular, we consider the task of approximating a point cloud sampled from a surface with a smooth parametric surface

S : [0, 1] 2 → R 3 with S (u) := N j=0 c j B j (u) , for u ∈ Ω (1)
defined in the spline space spanned by the basis (B j ). In this paper, we consider either the tensorproduct Bernstein polynomials or (truncated hierarchical) B-splines as basis functions B j . The objective of surface fitting is to find the best approximating surface, i.e.

min c0,...,c N ,u1,...,u |P| |P| i=1 ∥S(u i ) -x i ∥ 2 .
Optimizing the parameterization u as well as the surface S at the same time is a complicated non-linear problem that can be simplified by separating the optimization with respect to the parameterization u i ∈ [0, 1] 2 and the optimization with respect to the control points c j ∈ R 3 . Once a parameterization u of P has been determined, the best-approximating surface S can be found by solving the linear least-squares problem min c0,...,c N ,

|P| i=1 N j=0 c j B j (u i ) -x i 2 .
(

) 2
using standard techniques from computational linear algebra.

Our goal is then to look for the parameterization u that minimizes the residual of (2).

Standard methods

In order to demonstrate the efficiency and accuracy of our proposed method, we will compare its performance with the three parameterization methods presented in [START_REF] Floater | Meshless parameterization and surface reconstruction[END_REF]: uniform parameterization, inverse distance parameterization and shape preserving (SP) parameterization. In these three methods the parameterization of the boundary of the point cloud is first determined separately. This can be done using standard methods such as uniform, chord length and centripetal parameterization [START_REF] Piegl | The NURBS Book[END_REF] and their improvements [START_REF] Fang | An improved parameterization method for B-spline curve and surface interpolation[END_REF][START_REF] Balta | Dynamic centripetal parameterization method for b-spline curve interpolation[END_REF], using iterative schemes [START_REF] Hoschek | Intrinsic parametrization for approximation[END_REF][START_REF] Saux | An improved Hoschek intrinsic parameterization[END_REF] as well as using feed-forward neural networks [START_REF] Laube | Deep Learning Parametrization for B-Spline Curve Approximation[END_REF][START_REF] Scholz | Parameterization for polynomial curve approximation via residual deep neural networks[END_REF].

Then, a graph of the input point cloud P = P I ∪ P B is constructed by choosing neighborhoods N i ⊂ P for all i ∈ P I . Once the parameterization weights λ ik > 0 with k∈Ni λ ik = 1 are chosen, the parameterization is determined by the system of equations

u i = k∈Ni λ ik u k , (3)
where the given parameters u j for j ∈ P B serve as boundary conditions.

For the uniform and the inverse distance parameterization, the neighborhoods N i are chosen as N i = {k ∈ P : ∥x k -x i ∥ ≤ r} for some radius r and the weights are set to

λ ik = 1 |N i | (uniform)
or to

λ ik = 1 ∥x k -xi∥ ℓ∈Ni 1 ∥x ℓ -xi∥ (inverse distance).
For the shape preserving weights, the neighborhoods N i are chosen by first projecting all points x k in a large ball neighborhood of x i onto their best approximating plane and then finding a Delaunay triangulation of the resulting planar point cloud. The neighbors of i in the Delaunay triangulation form the neighborhood N i and the weights λ ij are the shape preserving weights for triangulations presented in [START_REF] Floater | Parametrization and smooth approximation of surface triangulations[END_REF]. Their defining property is that they preserve the distribution of the angles of the triangulation around x i .

Point cloud parameterization using BIDGCN

In order to predict the optimal parameterization we propose a network architecture based on BIDGCN. As detailed in the previous section, the task for our trained neural network is to predict parameters u i ∈ [0, 1] 2 ⊂ R 2 for the given point cloud such that the resulting solution c 0 , . . . c N ∈ R 3 of (2) results in the smallest possible residual. The neural network should then approximate an operator that assigns to each sufficiently large set of input features a parameterization that is optimal with respect to [START_REF] Wu | A comprehensive survey on graph neural networks[END_REF]. Without fixing some of the parameter values a priori, this operator is not uniquely defined, making it difficult to directly train neural network for this task. Inspired by the standard methods summarized in Section 3.2, we fix the parameters u j for the boundary points using the standard one-dimensional parameterization algorithms. This determines a well-defined specific parameterization operator that takes as input the positions of all points, together with the boundary parameterization, and gives as output the unique optimal parameterization of the interior points. We train our neural network to approximate this operator.

In order to predict the optimal parameterization with respect to a specific choice of boundary parameters, the neural network needs to take into account the boundary parameters as boundary conditions. This motivates our choice to apply BIDGCN to the parameterization problem.

Network architecture

We design a neural network based on the new boundary informed input layer described in Section 2. More precisely, our network consists of the boundary informed input layer, four hidden dynamic edge convolution layers, and a multilayer perceptron as the output layer. In the input layer, we have two multi-layer perceptrons (MLP), one for the edges in G I→I and one for the edges G B→I . For G I→I we train an MLP with layer sizes {6, 64, 64}, while for G B→I we train a MLP with layer sizes {8, 64, 64}. Note that the input dimension corresponds to the edge features in the two different graphs. In all hidden dynamic edge convolution layers the MLP has size {128, 64}. The output feature dimension of the hidden layers is deliberately chosen to be moderate because a new radius graph is computed after each layer with respect to the output features, which can be costly for high dimensions. Finally, the output of the last hidden layer is concatenated with the output of all hidden layers and fed into a final MLP of size {320, 256, 256, 2}.

After each hidden layer, the ReLu activation function is applied. Finally, after the output layer, the sigmoid activation function is applied to enforce that the predicted parameters lie in [0, 1] 2 . Since all layers in our network are convolutional except for the last layer that is applied vertexwise, the network can be applied to any point cloud, independent of its size.

The overall number of trainable parameters in this network architecture is 192,130. The architecture of our point cloud parameterizing network is shown in Fig. 3.

Training data

In our training, we employ tensor-product Bézier surfaces, meaning that B j (u) = B ab (u 1 , u 2 ) are the tensor-product Bernstein polynomials of bidegree (p 1 , p 2 ), defined as

B ab (u 1 , u 2 ) = p 1 a p 2 b (1 -u 1 ) p1-a (u 1 ) a (1 -u 2 ) p2-b (u 2 ) b ,
where a(j), b(j) is an ordering of {0, . . . , p 1 } × {0, . . . , p 2 }. This choice of the discrete space has the advantage that for a sufficiently large point cloud the best approximation with respect to a fixed parameterization of the point cloud is well defined without further regularization.

In order to train the network for parameterizing scattered point cloud data, we generate a data set consisting of 100,000 point clouds. Each point cloud P = P I ∪ P B contains 1000 interior points and 34 boundary points and thus the ratio of interior points to boundary points is equivalent to a uniformly distributed point cloud. The point clouds are sampled from biquadratic tensorproduct Bézier surfaces whose control points c ij were randomly sampled from

C ij = i 2 - 1 4 , i 2 + 1 4 × j 2 - 1 4 , j 2 + 1 4 × [-1, 1]
for i, j = 0 . . . 2 according to the uniform distribution. This choice of sampling space ensures a large variation of complexity of the surfaces while avoiding self-intersections. In order to achieve rotation-invariance, we rotate each surface around a randomly sampled axis in R 3 by a random angle. Besides the vertex features x i ∈ R 3 for all interior and boundary points, we store the exact parameters u i ∈ R 2 for all boundary points. Note that the choice of point cloud size in the the training data set does not mean that the trained network is limited to point clouds of this size. The network described in the previous section can be applied to any point cloud and we will observe in the numerical experiments that it performs well for a large range of point cloud sizes.

Normalization

Before evaluating the network on a point cloud from the training, validation and test data sets or from a real-world sample, we normalize the vertex features x i ∈ R 3 by translating and scaling them so that they lie in [0, 1] 3 . Due to the affine invariance of Bézier and B-spline surfaces this does not affect the optimal choice of parameters u i ∈ R 2 .

Loss function

As we want to train the neural network to predict parameters that minimize (2), we follow an unsupervised strategy and use the predicted parameters for the interior points as well as the prescribed parameters for the boundary points to fit a biquadratic polynomial surface to the point cloud. More precisely, we assemble the collocation matrix

A ij = B j (u i ),
for i = 1, . . . , |P| and j = 0, . . . , 8. Here, B j are the biquadratic tensor-product Bernstein polynomials and u i ∈ [0, 1] 2 are the prescribed parameters if i ∈ P B and the parameters predicted by the neural network if i ∈ P I .

The right hand side b ∈ R |P|×3 is given by the features of the point cloud, i.e.

b i = x i .
We then have reformulated (2) into the linear least squares problem

min c∈R 9×3 ∥Ac -b∥ 2 (4)
that we solve using QR decomposition, keeping track of the gradients with respect to the learnable weights of the neural network. The loss for the predicted parameterization of the interior points is the residual of (4). We train the network using stochastic gradient descent. We start the training with learning rate 0.1 and decrease the learning rate whenever the loss plateaus.

Adaptive THB-spline fitting with moving parameters

The need of automatic representations in terms of free-form spline models which ensure sufficient flexibility, effective geometric modeling options, high approximation accuracy, while simultaneously preserving low computational costs, led to the developement of (truncated) hierarchical B-spline (THB-spline) adaptive approximation schemes, see e.g. [START_REF] Giannelli | THBsplines: An effective mathematical technology for adaptive refinement in geometric design and isogeometric analysis[END_REF] among others. The truncated basis for hierarchical B-splines was originally introduced in [START_REF] Giannelli | THBsplines: The truncated basis for hierarchical splines[END_REF] to recover the partition of unity property in the hierarchical spline setting and reduce the interaction of hierarchical B-splines inserted at different refinement levels.

We consider a hierarchical spline model S : [0, 1] 2 → R 3 of the form (1) that approximates the input point cloud P, where B j : [0, 1] 2 → R are bivariate THB-splines [START_REF] Giannelli | THBsplines: The truncated basis for hierarchical splines[END_REF]. State of the art of scattered point cloud fitting schemes address the data parameterization and the design of the hierarchical spline model separately, see [START_REF] Kiss | Adaptive CAD model (re-)construction with THB-splines[END_REF][START_REF] Bracco | Adaptive fitting with THB-splines: Error analysis and industrial applications[END_REF][START_REF] Bracco | THB-Spline Approximations for Turbine Blade Design with Local B-Spline Approximations[END_REF]. Given a point cloud P, a suitable parameterization {u 1 , . . . , u |P| } is initially computed and, subsequently, on such a fixed parameterization, the hierarchical spline space is iteratively refined within the adaptive approximation scheme. We here consider as a starting point the adaptive global THB-spline least squares approach, originally proposed in [START_REF] Kiss | Adaptive CAD model (re-)construction with THB-splines[END_REF]. Let the features x i ∈ R 3 be the Cartesian coordinates of the items i ∈ P, and u 1 , . . . , u |P| a suitable parameterization. By starting with the initial polynomial approximation, a THB-spline approximation is iteratively computed until the error satisfies a given tolerance or a maximum number of iterations is reached. At each iteration of the adaptive loop a set of elements is marked for refinement based on the error computation, and the hierarchical mesh is subsequently refined. More precisely, the iterative approximation scheme is characterized by four main steps. First, a solution of the current fitting problem is computed. In particular, for a fixed (TH)B-spline space, the control points in (1) are computed by solving the penalized least squares problem

min c0,...,c N 1 2 |P| i=1 ∥S (u i ) -x i ∥ 2 2 + λJ (c 0 , . . . , c N ) , (5)
where the regularization term J is the thin-plate functional, i.e.

J (c 0 , . . . , c N ) = [0,1] 2 ∂ 2 S ∂u 2 1 2 2 + ∂ 2 S ∂u 1 ∂u 2 2 2 + ∂ 2 S ∂u 2 2 2 2 du 1 du 2 ,
whose influence is tuned by the weight λ ≥ 0.

In the second step of each iteration, the THBspline approximant is evaluated on the parameter sites {u 1 , . . . , u |P| } related to the data points P and the pointwise distance ∥S(u i ) -x i ∥ 2 2 is measured for each i ∈ P. The pointwise error constitutes the error indicator on which the last two steps (marking and refinement) of the adaptive scheme are developed. The parameter values u i ∈ [0, 1] 2 associated to the points i ∈ P whose error exceeds a certain input threshold ϵ > 0 are marked to be refined, since they identify regions of the parametric domain that need to be refined. The hierarchical mesh elements that contain at least one marked parameter u i ∈ [0, 1] 2 are dyadically refined, together with a certain number of neighboring cells, see [START_REF] Kiss | Adaptive CAD model (re-)construction with THB-splines[END_REF] for detailed information. The iteration of these four steps is performed until the maximum point wise errors is within an input tolerance or a maximum number of iterations is reached.

Note that any adaptive fitting scheme of this kind is characterized by an iterative enrichment of the approximation space to progressively improve the accuracy of the result. Consequently, in our setting, even if the the parametric values associated to the input observations are initially (quasi-)optimal, there is no evidence that optimality is preserved when the hierarchical spline space is adaptively refined. We then introduce a parameter correction step [START_REF] Hoschek | Intrinsic parametrization for approximation[END_REF] at every iteration of the fitting scheme, after the (local) refinement step.

The parameter correction approach addresses an optimization problem to identify the optimal intrinsic parameterization of a given spline model. The points P are projected on the current spline model S(u) and for each i ∈ P, its footpoint S(û i ) is identified. Hence, if ∥S(û i ) -x i ∥ 2 < ∥S(u i ) -x i ∥ 2 , the new (corrected) parameter for the point i ∈ P corresponds to ūi = ûi , otherwise ūi = u i . Note that, for each i = 1, . . . , |P|, the projection is performed by considering the problem

ûi = min u∈[0,1] 2 ∥S(u) -x i ∥ 2 2 , (6)
which can be solved, for example, with a Newtonlike method. Once the updated parameter values {ū 1 , . . . , ū|P| } are computed, the surface is also updated, and a new projection can be performed. Very few iterations of parameter correction are usually enough to obtain higher accuracy. Since we here repetitively applied the parameter update within the adaptive loop, a single correction step at any iteration is considered to suitably update the fitting result with respect to the current hierarchical spline configuration. Note that, as shown in the last examples of the next section, the initial parameterization plays a key role to obtain high quality results.

Numerical results

We implemented our method using the PyTorch [START_REF] Paszke | PyTorch: An Imperative Style, High-Performance Deep Learning Library[END_REF] and PyG [START_REF] Fey | Fast graph representation learning with PyTorch Geometric[END_REF] libraries. We tested the deep learning model on unseen data sets, generated in the same way as the training data, as described in Section 4.2, as well as on a number of real-world data sets. Moreover, we also consider noisy data, where we added Gaussian noise with varying standard deviation. We start with an ablation study (Section 6.1) of our new boundary informed layer to demonstrate its fundamental role in tackling the parameterization problem. We then show the comparison of the proposed method with the three standard approaches described in Section 3.2 (see Sections 6.2, 6.3, 6.4, and 6.5). In order to further test the generalization capabilities of BIDGCN on complex data sets and adaptive spline constructions, we consider in Sections 6.6 and 6.7 the adaptive THB-spline approximation scheme with moving parameters introduced in the previous section.

The training, as well as all evaluations of the different methods, was performed on a standard workstation computer with an NVIDIA GeForce GTX 1060 GPU.

Ablation study

As previously described, the parameterization problem cannot be solved without taking into account the parameterization of the boundary curves. In order to demonstrate this necessity empirically and also to show that it is fulfilled by our new boundary informed layer, we first compare the performance of our network with a trained DGCNN, see Section 2, whose architecture is equal to the one of our point cloud parameterizing network shown in Fig. 3, but with a standard dynamic edge convolution layer as input layer instead of the boundary informed input layer. This means that for boundary vertices, this network only acts on the positional features in the same way as it acts on interior layers.

In Fig. 4, we show the results on two surfaces from our test data set. We observe that the parameterization predicted by DGCNN contains large voids and is far from the original. This is expected, since the optimal parameterization is not uniquely defined by the positional vertex features only. On the other hand, BIDGCN correctly takes the boundary conditions into account and predicts parameters that are very close to the original ones, leading to a much better approximation with a biquadratic surface. Quantitatively, DGCNN results in a mean squared error of 7.49 • 10 -4 for the surface on the left and 5.25 • 10 -4 for the surface on the right while BIDGCN results in a mean squared error one order of magnitude smaller, namely 5.26 • 10 -5 and 3.50 • 10 -5 , respectively.

When looking the average mean squared error over 100 point clouds that were sampled from quadratic surfaces in our test set, we observe that DGCNN resulted in a mean squared error of 2.93 • 10 -3 and BIDGCN resulted in a mean squared error of 6.38 • 10 -5 . This shows that the boundary informed layer results in a large improvement of the accuracy. This is in line with the theoretical consideration, that the optimal parameterization that is to be predicted is only uniquely defined when taking into account boundary conditions. As the network architectures of the two networks is similar, their evaluation times does not differ significantly.

Point clouds belonging to the test data set

In our first experiment, we study the performance of our neural network when applied to data from the same class as the training data, i.e., data sampled from biquadratic surfaces. As a benchmark, we use the three parameterization methods presented in [START_REF] Floater | Meshless parameterization and surface reconstruction[END_REF], parameterization with uniform, inverse distance and shape preserving weights. For these methods, we choose the radius for defining the local neighborhoods adaptively depending on the number of points per point cloud as

r = 3 |P| . ( 7
)
Choosing the optimal radius r for these methods is a complicated manual task that for general data cannot be solved efficiently. Here, the factor |P| is derived from the number of points on an axisaligned line in a uniformly distributed point cloud. The factor 3 was determined empirically to be close to optimal for the standard methods when parameterizing point clouds of 200 points. In order to compare the methods, we evaluate them on 100 point clouds and report the mean squared error as well as the total computation time needed to parameterize and approximate all 100 point clouds.

In particular, we study the dependence of the accuracy and the computation time on the number of points in each point cloud. Note that while our network was trained only on point clouds of size 1000, it is important to ensure that it gives good results for point clouds of any size. The left column of Fig. 5 shows the comparison of the methods on point clouds sampled from surfaces that were generated in the same way as the training data set for our neural network, described in Section 4.2.

All the plots are semi-log plots with a logarithmic scale used for the time and error axes. We observe that the mean squared errors produced by our method are much smaller than the ones resulting from the uniform and inverse distance parameterization. On the other hand, while the accuracy of our method on this synthetic data is similar to the one of the shape preserving parameterization, the computation time of the neural network based method is much lower than the time needed for all three other methods. In particular, SP is very costly, with a computation time which is about two orders of magnitude higher than the one of BIDGCN. Moreover we note that, even if the network was trained exclusively on data with 1000 points per point cloud, the performance of BIDGCN does not depend on the size of the scattered data set.

Evaluations on noisy data

Measured real-world data are always subjected to noise. For this reason, we study the behavior of the trained network and the benchmark methods when applied to noisy data. We evaluate all methods on 100 point clouds sampled from surfaces that were generated as described in Section 4.2.

We added Gaussian noise of standard deviation 0.005 and 0.01 to all surfaces. The middle column of Fig. 5 shows the behavior of all four methods when applied to data with Gaussian noise of standard deviation 0.005, while the right column of Fig. 5 shows their behavior on surfaces with noise of size 0.01. We observe that when applied to noisy data, BIDGCN performs much better than the three other methods. This means that the network is able to predict good parameterizations even for data that it was not trained on. As in the non-noisy data case, we observe that the evaluation of our method is much faster than the evaluations of the other methods. In particular, the speed-up with respect to SP is significant.

A further advantage of BIDGCN is that it is very robust with respect to the presence of noise, even if the noise becomes very large. For the standard methods suitably choosing the radius that determines the local neighborhood becomes near impossible when the data is non-regularly distributed. For a fixed choice of radius or our simple adaptive choice [START_REF] Sharp | DiffusionNet: Discretization agnostic learning on surfaces[END_REF], this means that these methods often fail due to neighborhoods that are too small. In particular, SP needs enough points in each neighborhood to generate a Delaunay triangulation. Fig. 6 shows how many surfaces out of 100 surfaces with Gaussian noise of standard deviation 0.05 could be successfully parameterized by the four methods. We observe that BIDGCN was always successful, while the number of surfaces that could be parameterized using SP strongly decreases with the number of points per point cloud. One possible way to tackle this problem could be by suitably choosing a different radius r i for each point i ∈ P; however, there is no obvious way how to realize such an adaptive local choice in a robust way. Moreover, a local choice, if one exists, would further increase the computational complexity and overall efficiency risks to deteriorate even further. Finally, we report the mean squared errors of the surfaces that were successfully parameterized in Fig. 7. Also in this case the neural network results in significantly smaller errors.

Generalization to higher degrees

In order to test the generalization properties of the neural network in terms of degree, we apply it to point clouds sampled from surfaces of higher bidegree, namely (p, p) with p = 3, 4, 5, with and without Gaussian noise. In doing so, we use tensor-product Bézier surfaces of the same degree for fitting the parameterized point clouds. Fig. 8 shows the mean squared errors as well as the timings for point clouds without noise. We observe that BIDGCN results in similar mean squared errors compared to the SP parameterization. However, the computation time for BIDGCN is more than one order of magnitude smaller with respect to SP.

Fig. 9 shows the results on point clouds with added Gaussian noise of standard deviation 0.01. We observe that in this case BIDGCN results in improved mean squared errors compared to the SP parameterization, while the computation time is still over one order of magnitude smaller.

Visual comparison of the methods

While the error in the previous examples was averaged over a large number of point clouds, we will now present particular examples that visually demonstrate the approximation power of the neural network. We applied the BIDGCN as well as SP to the first two point clouds of size 1000 in our test data set, with and without noise, sampled from biquadratic surfaces. For SP we compare two choices of the radius: r = 0.09 ≈ 3/ |P| and r = 0.2. Fig. 10 and 11 show the parameters, the evaluation of the resulting surface at the parameters as well as the fitted surface for both methods as well as the input data. While the point cloud in Fig. 10 does not have any noise, we added Gaussian noise of standard deviation 0.01 to the point cloud shown in Fig. 11. We observe that the neural network predicts parameters that are visually very close to the original parameters, both for the noisy and the non-noisy case. On the other hand, the behavior of SP largely depends on the choice of the radius r. For r = 0.09, SP parameterization performs well but appears to result in larger deviations from the original parameters compared to BIDGCN. For the choice r = 0.2, the parameterization contains large gaps in the interior of the parameter domain. This shows that for using SP effectively, one needs to carefully choose the radius r, while BIDGCN is able to predict the optimal graph in the first layers of the network architecture and therefore does not require any fine-tuning.

Reconstruction and comparison on real-world data

In this experiment we process real world point clouds of different size, representing a Nefertiti face model and we lead a comparison between the standard parameterization methods and the proposed BIDGCN, in terms of accuracy and computational time. As concerns the standard parameterization methods, a suitable choice of the radius r needs to be selected. To produce fair comparisons and avoid unreasonable parameter value distributions, as illustrated in Fig. 10 and 11, we execute an heuristic search for r = dataset by computing the polynomial approximation of bidegree (2, 2), and select the radius r giving the best mean squared error. By analyzing the value of the error with respect to r, we decide to fix the radius as defined in [START_REF] Sharp | DiffusionNet: Discretization agnostic learning on surfaces[END_REF] since it leads to the best approximation results for almost all the real data point clouds.

As first analysis, we collect 6 data sets of 532, 1043, 2062, 3253, 6024, 8818 points, acquired from the same Nefertiti face model. We then compute the parametric values for each Nefertiti point cloud with the standard and the proposed BIDGCN methods to subsequently construct a polynomial biquadratic least squares approximation. For each method we report the mean squared error associated with the reconstructed polynomial surface and the total computational time needed to parameterize and approximate each point cloud in Fig. 12. In line with the trends observed on the synthetic data investigated in the previous sections, the uniform and inverse distance parameterizations lead to mean squared errors that are, in general, much higher than SP and BIDGCN. On the other hand, the error values obtained with BIDGCN and SP are similar but our method scales favourably with the dimension of the point cloud always having the lowest computational time.

As a second analysis on the Nefertiti model, we consider the adaptive THB-spline fitting of the point cloud of size 8818, shown in Fig. 13 (left). In this case, we parameterize the input scattered data with the standard SP method as well as the proposed BIDGCN method and subsequently reconstruct the penalized least square adaptive spline model by performing the algorithm for bidegree [START_REF] Wu | A comprehensive survey on graph neural networks[END_REF][START_REF] Wu | A comprehensive survey on graph neural networks[END_REF]. Note that among all the listed values for the radius of the shape preserving method that we tested on this point cloud, the choice r = 3 √

|P|

gives the best results in terms of polynomial approximation. To properly compare the results, the refinement tolerances are chosen such that they result in (almost) the same number of degrees of freedom in the final hierarchical spline model. When the SP method is considered, we obtain a THB-spline model with 5492 degrees of freedom, which is 62% of the total number of scattered input items, characterized by a mean squared error of 9.32e-4. If we parameterize the Nefertiti data with the BIDGCN method, the THB-spline approximation has 5228 degrees of freedom, i.e. 59% of the number of points, and affects the final quality of the adaptive spline approximation.

Hierarchical spline reconstruction of different degrees

In this last example we investigate the generalization capabilities of our BIDGCN parameterization to spline configurations of different bidegrees for the reconstruction of different point clouds. We start by performing an adaptive THB-spline fitting of a human face model of 9283 points, shown on the right of Fig. 13. In particular, we reconstruct the THB-spline models of the input point cloud by performing adaptive penalized least squares approximation for bidegree [START_REF] Gilmer | Neural message passing for quantum chemistry[END_REF][START_REF] Gilmer | Neural message passing for quantum chemistry[END_REF] and [START_REF] Welling | Semi-supervised classification with graph convolutional networks[END_REF][START_REF] Welling | Semi-supervised classification with graph convolutional networks[END_REF]. In the first case, the final reconstructed geometry has 2687 degrees of freedom (29% of the number of points), that approximate the input point cloud registering a mean square error of 6.86e-5. As concerns the fitting for bidegree [START_REF] Welling | Semi-supervised classification with graph convolutional networks[END_REF][START_REF] Welling | Semi-supervised classification with graph convolutional networks[END_REF], the final reconstructed model is a THB-spline geometry with 2793 degrees of freedom (30% of the number of points) that approximates the input point cloud with a mean squared error of 6.08e-5. The reconstructed THB-spline models, their scaled error distributions on the point cloud, the computed parametric values, and the corresponding hierarchical meshes are illustrated in Fig. 14. We end by considering two additional point clouds of 9636 items collected from a ship hull and a wind turbine. In particular, we approximate the two scattered data sets by performing adaptive THB-spline fitting for bidegree [START_REF] Qi | Pointnet: Deep learning on point sets for 3d classification and segmentation[END_REF][START_REF] Wu | A comprehensive survey on graph neural networks[END_REF]. The point clouds, the reconstructed THB-spline model, the point wise error distribution on the scattered data and the hierarchical meshes are displayed in Fig. 16 for both the ship hull (top) and the wind turbine (bottom). In addition, the final THB-spline model for the ship hull geometry has 4092 degrees of freedom (i.e. 42% of the number of points), and a mean square error equal to 3.38e-6, whereas the THB-spline model for the wind turbine has 2449 degrees of freedom, equivalent in magnitude to 25% of the number of input data, and a mean squared error of 5.61e-6.

Conclusion

We introduced the novel graph convolutional neural network BIDGCN for learning on point clouds.

In particular, we developed a new input layer which takes into account boundary conditions and propagates this information into the interior of the point cloud. Further hidden layers can then incorporate the boundary information into the prediction.

We applied this network to the problem of point cloud parameterization for surface approximation. When compared to standard parameterization methods, the architecture based on BIDGCN parameterization scheme achieves higher accuracy with less computational time, especially in the presence of noise in the data. In addition, BIDGCN is parameter independent and robust, avoiding failure issues in the case of sparse neighborhoods. Finally, we demonstrate the suitability of the network parameterization results for adaptive surface reconstruction with a new THBspline fitting scheme with moving parameters. While we obtain valuable results on real-world data with BIDGCN trained on synthetic data, the performance of this new architecture can be highly optimized in real application settings by a suitable training on complex data sets in line with the specific test cases. This specialized training will be an important subject of our future work. Moreover, we plan to apply the new boundary informed graph convolutional layer to different computational tasks where boundary conditions play a role. A promising potential application of BIDGCN is the solution of partial differential equations on point cloud data. This could be done directly by adding Dirichlet or Neumann boundary conditions to the boundary vertices features and training a network to predict the solution of a specific PDE. Alternatively, a neural network based on BIDGCN can perform a supporting task for classical numerical methods for solving PDE, such as the mesh generation problem considered in [START_REF] Chan | Locally refined quad meshing for linear elasticity problems based on convolutional neural networks[END_REF], where adaptive mesh refinement for planar geometries was investigated.

Fig. 1 :

 1 Fig. 1: Boundary informed input layer.

  x

Fig. 3 :

 3 Fig. 3: Network architecture of the point cloud parameterizing network, n = 3 is the dimension of the point cloud and d = 2 is the parametric dimension.

Fig. 4 :Fig. 5 :

 45 Fig. 4: Comparison of our BIDGCN network with a pure dynamic edge convolution network trained for the point parameterization problem. Top row: Original parameters and parameters predicted by the networks. Middle row: Input data and the evaluation of the fitted surfaces at the predicted parameters. Bottom row: The original surface and the fitted biquadratic surfaces.

Fig. 6 :

 6 Fig. 6: Number of successfully parameterized surfaces out of 100 biquadratic surfaces with added Gaussian noise of size 0.05.

Fig. 7 :

 7 Fig.7: Mean squared error when applying the four methods to point clouds of different sizes sampled from 100 biquadratic surfaces with added Gaussian noise of standard deviation 0.05. All cases where a method did not succeed were removed from the mean.

3 √

 3 |P|, 0.05, 0.075, 0.1, . . . , 0.25, 0.275, 0.3 on each

Fig. 9 :

 9 Fig. 9: Mean squared error (top) and computation time (bottom) when applying the four methods to point clouds of different sizes sampled from 100 of bidegree (3, 3) (left), (4, 4) (middle) and (5, 5) (right), with added Gaussian noise of standard deviation 0.01.

Fig. 10 :

 10 Fig. 10: Visual comparison of the BIDGCN and SP for radius r = 0.2 and r = 0.09 on a point cloud of size 1000 sampled from a biquadratic surface without noise.

Fig. 11 :Fig. 12 :

 1112 Fig. 11: Visual comparison of the BIDGCN and SP for radius r = 0.2 and r = 0.09 on a point cloud of size 1000 sampled from a biquadratic surface with Gaussian noise of standard deviation 0.01.

Fig. 13 :

 13 Fig. 13: Nefertiti (left) and face (right) point clouds characterized by 8818 and 9283 points, respectively.

Fig. 14 :

 14 Fig. 14: Hierarchical spline model reconstruction of the Nefertiti point cloud shown on left of Fig. 13 for SP parameterization with 5492 degrees of freedom (top), as well al for the BIDGCN parameterization with 5288 degrees of freedom (bottom) with bidegree (2, 2). The reconstructed THB-spline models are shown together with the scaled error distributions on the point cloud, the computed parametric values, and the corresponding hierarchical meshes (from left to right).

Fig. 15 :

 15 Fig. 15: Hierarchical spline model reconstruction of the face point cloud shown on the right of Fig. 13 for bidegree (3, 3) (top) and (4, 4) (bottom) using the BIDGCN parameterization. The reconstructed THBspline models are shown together with the scaled error distributions on the point cloud, the computed parametric values, and the corresponding hierarchical meshes (from left to right).

Fig. 16 :

 16 Fig. 16: The point cloud, the hierarchical spline model reconstruction, the scaled error distribution and the hierarchical mesh for a ship hull (top) and a wind turbine (bottom) using the BIDGCN parameterization with bidegree (5, 2).

  Input data: P I (|P I | × n) and P B (|P B | × (n + d))

	Boundary informed input layer
	MLP sizes {6, 64, 64} and {8, 64, 64}
	ReLu
	|P I | × 64
	Dynamic edge convolution
	MLP size {128, 64}
	ReLu
	|P I | × 64
	Dynamic edge convolution
	MLP size {128, 64}
	ReLu
	|P I | × 64
	Dynamic edge convolution
	MLP size {128, 64}
	ReLu
	|P I | × 64
	Dynamic edge convolution
	MLP size {128, 64}
	ReLu
	|P I | × 64
	Concatenation
	MLP size {320, 256, 256, d}
	Sigmoid
	|P

I | × d

Data availability

The datasets generated during and/or analysed during the current study are available from the corresponding author on reasonable request.

Declarations

Competing interests

The authors have no relevant financial or nonfinancial interests to disclose.


Eddy Renaud
email: eddy.renaud@inrae.fr

David Brunet
email: david.brunet@swde.be

Eric Smit
email: eric.smit@swde.be

Yves Le Gat
email: yves.legat@inrae.fr

Applied research for multi-scale asset management of the Walloon Water Company's drinking water supply network

Keywords: Active leakage control, Drinking water networks, Sociology of work, Statistical modelling, User preferences, Water Infrastructure Asset Management

des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

This scientific program is based on a service-performance driven WIAM, mainly in terms of preservation of drinking water resources (control of water losses), continuity of service (control of network failures), quality of water, and service pressure. It adopts a multi-scale approach, which considers service performance from the local level (pipe segment) to the territorial level (District Metered Area (DMA), entire system) and from the short term (annual work programing) to the long term (investment strategy for the next decades). This requires being able to assess user vulnerabilities and knowing what elements of the network environment are affected by poor service performance. Furthermore, the methods and tools developed must take into account the operational organization and governance of the SWDE.

The GePaME project is structured around three interrelated subjects that are addressed through an interdisciplinary approach involving the methods listed in Table 1.

Results

The main scientific results of the GePaME project are the followings:

-Data processing and missing data completion methods and algorithms to feed the models -Statistical models making it possible to jointly predict networks breaks and detectable leaks, to simulate in the long term the probabilities of a pipe section to be renewed and to estimate the minimum night flow of a DMA according to its characteristics -A sociological analysis of the relationships between people / tools / decisions / actions at SWDE scale and from the particular angle of the implementation of the leak detection activity -An econometric analysis of users' willingness-to-pay for measures aimed at mitigating the effects of droughts on the water supply service

The operational exploitation of these scientific advances by SWDE concerns:

-Improving the interoperability of databases feeding the information system, rationalizing the recording of network failures and detected pipe leaks and the automation of data processing methods for asset management purposes -Improving the organization of leak detection (prioritization of DMAs, choice of types of action, information management) -The implementation of a new multi-criteria decision support tool for pipes renewal, making it possible to rank all the pipes in the network and taking into account the reliability of data and the vulnerability of users. -The experimentation of a simulation tool allowing to evaluate and compare the long-term performance of the system resulting from various investment scenarios

Conclusion

GePaME is a partnership, multi-scale and interdisciplinary applied research project.

Thanks to this partnership, INRAE had access to a large amount of structured data, benefited from the expertise of SWDE, and was thus able to renew its research questions in line with operational needs. At the same time, SWDE was able to engage in an improvement process by adjusting its practices to the project's innovations.

The multi-scale approach made it possible on the one hand to take into account the interactions between the organization, the territories and the infrastructures, and on the other hand, to reinforce the coherence between the short-term programming of pipes renewal and the long-term service performance strategy.

Finally, the joint mobilization of engineering sciences and social sciences, made it possible to take into account the human context (organization and internal practices of the SWDE, point of view of users) in the construction of the technical proposals.

The tools and methods developed within the framework of the project are a step towards a more integrated approach of WIAM practices and water loss control. After a necessary phase of appropriation and evaluation, they open up prospects to grasp a more global vision of the efficiency and performance of the service. This means, for example, taking into consideration disaster exposure, adaptation to climate change or the impacts of evolving technologies on skills needs and investment strategies.

Table 1 :

 1 Scientific fields and methods of the GePaME project An operational monitoring committee, bringing together INRAE engineers and researchers and SWDE experts, coordinates and adapts the Research and Development work through the framework of monthly working meetings.

	Topics	Scientific fields	Methods
	Transversal tasks	Data science	Statistical and spatial processing of SWDE data
		Sociology	Bibliography and interviews to analyse decision-making factors and governance
		Computer science	Implementation of methods in the SWDE' information system
	Leak process and	Water network engineering Multilinear modelling of DMAs' minimum night flows according to their characteristics
	effectiveness of	Sociology	Interviews and participative observations with leak detection teams
	leak detection	Statistics	Modelling of the leakage likelihood according to pipeline and DMA characteristics
	Multi-objective	Water network engineering Pipe failure prediction and multi-criteria analysis of their impacts on service performance
	and constraint	Environmental economics	Survey to assess user vulnerability to water droughts and their consequences in terms of
	selection of		water supply
	renewal works	Statistics	Coupled modelling of the risks of breakage and detectable leak according to pipes
			characteristics
	Long-term	Water network engineering Tuning of long-term performance indicators
	simulation of	Statistics	Development of a long-term risk-driven probabilistic algorithm integrating the different
	WIAM strategies		causes of renewal (constraint, opportunity, choice)

EFFICIENT 2023 -Pathways to sustainable development on Efficient Urban Water Management September, 13-15 2023. Bordeaux, France