# Neural Slicer for Multi-Axis 3D Printing

Tao Liu  0000-0003-1016-4191 The University of ManchesterManchesterUnited
Kingdom ,  Tianyu Zhang  0000-0003-0372-0049 The University of
ManchesterManchesterUnited Kingdom ,  Yongxue Chen  0000-0001-6236-4158 The
University of ManchesterManchesterUnited Kingdom ,  Yuming Huang
0000-0001-5900-2164 The University of ManchesterManchesterUnited Kingdom and
Charlie C. L. Wang  0000-0003-4406-8480 The University of
ManchesterManchesterUnited Kingdom

###### Abstract.

We introduce a novel neural network-based computational pipeline as a
representation-agnostic slicer for multi-axis 3D printing. This advanced
slicer can work on models with diverse representations and intricate topology.
The approach involves employing neural networks to establish a deformation
mapping, defining a scalar field in the space surrounding an input model.
Isosurfaces are subsequently extracted from this field to generate curved
layers for 3D printing. Creating a differentiable pipeline enables us to
optimize the mapping through loss functions directly defined on the field
gradients as the local printing directions. New loss functions have been
introduced to meet the manufacturing objectives of support-free and strength
reinforcement. Our new computation pipeline relies less on the initial values
of the field and can generate slicing results with significantly improved
performance.

curved slicing, neural network, field optimization, multi-axis motion, 3D
printing

††ccs: Computing methodologies Shape modeling††ccs: Computing methodologies
Machine learning approaches

(a)(b)(c)(d)

Figure 1. Our Neural Slicer for multi-axis 3D printing can compute curved
layers for a Bunny Head model (a), incorporating diverse representations
within the same model. These representations encompass a tetrahedral mesh
(depicted in blue) as solid, an open surface (depicted in green) as the kernel
of a shell, the skeletons (depicted in thin black lines) of struts as
cylindrical solids, and also the skeletons (depicted in bold black lines) of
tubular solids. (b) The curved layers are generated through the optimization
of a neural network-based mapping according to the objectives of support-free
(SF) and strength reinforcement (SR), which are formulated as the loss
functions in terms of local printing directions (LPDs). (c) The red bands in
the histograms depict the required angle ranges between LPDs and other
reference vectors. (d) Curved layers satisfying these requirements can be
effectively generated to supervise the physical fabrication taken on a robotic
system with eight degree-of-freedoms (DOFs). The mapping is numerically
computed with the assistance of a caging mesh that is independent of the
discrete representations of the input model. In this example, the genus
numbers of the input model and the caging mesh are $g=22$ and $g=0$
respectively.

##  1\. Introduction

The additional degree-of-freedoms (DOFs) provided by multi-axis 3D printing
offer many advantages over planar layer based printing, including the relieved
need for support structures (Dai et al., 2018), the improved surface
smoothness (Etienne et al., 2019), and the enhanced mechanical strength (Fang
et al., 2020). For an input model $\mathcal{M}$, these methods always define a
scalar field $G(\mathbf{x})$ inside the input model $\mathcal{M}$ and extract
the isosurfaces of $G(\mathbf{x})$ as curved layers for multi-axis 3D
printing. Different strategies have been applied to compute an optimized field
$G(\mathbf{x})$.

(a)(b)

Figure 2. An illustration of the distortion caused by indirect optimization in
$S^{3}$-Slicer: (a) the layers generated in the deformed space, and (b) curved
layers in the model space obtained by the mapping. In the deformed space, the
support-free requirement has been fulfilled – i.e., the angle between the
printing direction and surface normal is less than $135^{\circ}$. However, the
angle expands when mapped back to the model space, resulting in a larger
overhang area (see also the 3D printing result shown in Fig.7).

The most recent work of $S^{3}$-slicer (Zhang et al., 2022) can achieve
multiple objectives on the same model by deforming the tetrahedral mesh of
$\mathcal{M}$ via nonlinear optimization driven by rotations. Height values on
the deformed mesh are mapped back to the input model as $G(\mathbf{x})$ and
therefore obtain the curved layers. Although it can generate curved layers
satisfying multiple manufacturing objectives, $S^{3}$-Slicer also has the
following problems:

  * •

The computation requires a tetrahedral mesh with high quality – i.e., the mesh
can become dense and challenging to obtain when applied to models with complex
geometry and topology;

  * •

Objective functions of optimization are indirectly defined as rotations of
elements rather than directly on the curved layers – as a result, there is a
risk of violating manufacturing requirements due to the distortion between the
deformed space and the model space (see Fig.2 for an example);

  * •

The success of nonlinear optimization is significantly dependent on the
initial pose of an input model, which may give sub-optimal results if the
initial orientation is not well-selected (see Fig.9 for an example).

We propose a computational pipeline based on neural networks (NN) to address
the challenges associated with curved slicers in multi-axis 3D printing. For
instance, as shown in Fig.1, our slicer generates curved layers for a Bunny
Head model with diverse representations and complex topology. The curved
layers are optimized in accordance with manufacturing objectives of support-
free (SF) and strength reinforcement (SR), formulated as the requirements for
local printing directions (LPDs). The resultant curved layers for the Bunny
Head model have been verified through physical fabrication conducted on a
robotic system with 8-DOFs.

###  1.1. Our Method

For all discrete representations, an input model $\mathcal{M}$ can always be
evaluated as an implicit function $H(\mathbf{x})$, where the value defines
whether a query point $\mathbf{x}$ is inside (i.e., $H(\mathbf{x})<0$) or
outside the solid of $\mathcal{M}$ (i.e., $H(\mathbf{x})>0$). The smaller the
absolute value $|H(\mathbf{x})|$ returned, the closer the query point
$\mathbf{x}$ is to the boundary surface of $\mathcal{M}$, approximated as the
zero level-set of $H(\mathbf{x})$. The function $H(\mathbf{x})$ only needs to
exhibit this monotonic property near the model’s boundary, rather than
strictly adhering to being a distance field. We compute a continuous mapping
$\lambda:\mathbf{x}\mapsto\mathbf{y},\forall\mathbf{x},\mathbf{y}\in\mathbb{R}^{3}$
that defines the scalar field $G(\mathbf{x})$ for layer generation by the
$z$-component of $\mathbf{y}$. The optimization of this mapping occurs within
the model $\mathcal{M}$ (i.e.,
$\forall\mathbf{x}\in\mathbb{R}^{3},H(\mathbf{x})\leq 0$), and it is tailored
to meet various manufacturing demands specific to multi-axis 3D printing.

We parameterize the mapping with two continuous functions,
$\mathbf{q}(\mathbf{x})$ and $\mathbf{s}(\mathbf{x})$, represented as
$\lambda(\mathbf{q}(\mathbf{x}),\mathbf{s}(\mathbf{x}))$, where
$\mathbf{q}(\mathbf{x})$ and $\mathbf{s}(\mathbf{x})$ specify the local
deformation’s quaternion and scaling ratios for all
$\mathbf{x}\in\mathbb{R}^{3}$. The functions $\mathbf{q}(\mathbf{x})$ and
$\mathbf{s}(\mathbf{x})$ are modeled by neural networks with coefficients
$\theta$, which are the weights associated with activation functions. The
mapping $\lambda(\cdot)$ is computed numerically as a deformation applied to a
volumetric mesh $\mathcal{C}$ caging the input model $\mathcal{M}$ with our
effort to make the deformation differentiable. Then, the scalar field is
defined as

(1)
$G(\mathbf{x}):=\mathsf{proj}_{z}\mathbf{y}=\mathsf{proj}_{z}\lambda_{\theta}(\mathbf{x}).$

Different fields are defined by different network coefficients $\theta$.

With this new formulation, we develop a neural network-based computational
pipeline to generate the optimized curved layers (i.e., the scalar field
$G(\mathbf{x})$), utilizing $\theta$ as variables. This new pipeline
effectively addresses all the previously mentioned challenges:

  * •

The highly nonlinear nature of NN can help model more complicated deformation
than those defined as piecewise functions in tetrahedral elements (i.e., the
optimization is less blocked by the topology of meshes employed);

  * •

By establishing a differentiable neural pipeline for evaluating
$G(\mathbf{x})$, we can directly define loss functions based on
$G(\mathbf{x})$ and its gradients $\nabla G(\mathbf{x})$ (e.g., the LPD
requirements for SF and SR) at any point in the computation domain;

  * •

With our novel neural pipeline, LPDs can be easily adjusted during the
optimization of the scalar field $G(\mathbf{x})$, facilitated by the
robustness of modern stochastic NN solvers to different initial guesses.

We employ NN to represent continuous functions for $\mathbf{q}(\mathbf{x})$
and $\mathbf{s}(\mathbf{x})$ so that the mapping
$\lambda(\mathbf{q}(\mathbf{x}),\mathbf{s}(\mathbf{x}))$ can be effectively
and efficiently optimized by borrowing the computational power of the modern
machine learning pipeline.

The volumetric mesh employed for caging an input model $\mathcal{M}$ is
independent of $\mathcal{M}$’s discrete representation. It is an intermediate
representation used in numerical computation to evaluate the mapping
$\lambda$. This volumetric mesh is also employed for extracting the
isosurfaces of $G(\mathbf{x})$. The resultant curved layers are obtained by
trimming these isosurfaces using the implicit solid defined by
$H(\mathbf{x})$.

###  1.2. Contributions

We present a neural network-based computational pipeline as a representation-
agnostic slicer for multi-axis 3D printing, which makes the following
technical contributions:

  1. (1)

Formulation of the curved slicing problem as an optimization task for two
continuous functions, $\mathbf{q}(\mathbf{x})$ and $\mathbf{s}(\mathbf{x})$,
defining a mapping $\lambda(\mathbf{q}(\mathbf{x}),\mathbf{s}(\mathbf{x}))$,
and consequently, the scalar field $G(\mathbf{x})$ for slicing.

  2. (2)

Development of a differentiable pipeline of neural networks for optimization,
with loss functions directly based on $\nabla G(\mathbf{x})$ (representing
real LPDs), which reduces dependency on initial guesses.

  3. (3)

Derivation of loss functions within the neural pipeline to address
manufacturing objectives for multi-axis 3D printing, including support-free
and strength reinforcement.

To the best of our knowledge, this is the first curved slicer that can handle
models with diverse representations and complex topology. The manufacturing
demands on LPDs are directly optimized, leading to improved slicing solutions,
as validated through physical fabrication experiments.

##  2\. Related work

###  2.1. Multi-axis 3D Printing

Conventional planar layer based 3D printing has relatively simpler software
implementation, which however introduces many limitations such the weak
mechanical strength (Ahn et al., 2002) and the poor surface quality
(Chakraborty et al., 2008). Supporting structures are needed for large
overhangs, leading to challenges such as hard to remove, surface damage and
material waste (Zhang et al., 2015). With the help of additional DOFs in
motion, multi-axis 3D printing offers opportunities to improve these aspects
substantially. Recent developments enable the function of support-free
material accumulation (Mitropoulou et al., 2020; Dai et al., 2018; Wu et al.,
2016; Huang et al., 2016), enhance the mechanical strength of 3D printed parts
(Tam and Mueller, 2017a; Fang et al., 2020), and reduce the staircase effect
(Etienne et al., 2019). Many of these approaches were implemented on robotic
hardware (e.g., (Li et al., 2022; Wu et al., 2020; Bhatt et al., 2022)) or
those modified from a multi-axis CNC machine (e.g., (Fang et al., 2020; Zhong
et al., 2023; Zhang et al., 2021)), which have caught a lot of attention in
the community of computational fabrication (e.g., (Duenser et al., 2020;
Bartoň et al., 2021)).

Scalar fields have been widely employed to define curved layers inside a solid
by the isosurfaces. For instance, Dai et al. (2018) sought the help of voxel
representation to progressively compute a material accumulation field while
keeping a convex-front as a conservative method to prevent potential
collisions. To improve the surface quality of 3D printed models, Etienne et
al. (2019) developed a height-field like deformation to generate slightly
curved layers for 3-axis machines. Fang et al. (2020) presented a pipeline
exploiting the anisotropy of mechanical properties in fused filament
fabrication with the help of vector fields. They first optimized a vector
field to align with the minimal stresses, and later computed a scalar field
with its gradients approximating the vector field. The gradients of scalar
fields are not directly optimized. Similar indirect problem exists in the
recent work of $S^{3}$-Slicer (Zhang et al., 2022). Moreover, many of these
methods require a volumetric mesh with high quality, which is challenging for
models with complex geometry.

###  2.2. Field-based Optimization

The challenge of achieving different objectives in multi-axis 3D printing can
be formulated as a problem of optimizing LPDs at different points, essentially
forming vector fields. A possible approach to model LPDs is to extend the
techniques for processing vector fields defined on mesh surfaces (e.g., (De
Goes et al., 2016)) to volumetric meshes. Another option is frame-based field
optimization (Huang et al., 2012; Li et al., 2012; Ray et al., 2016), which
has been used to tackle other problems of computational fabrication (e.g.,
(Maestre et al., 2023; Mitra et al., 2023)). However, these methods can be
very time-consuming, as evidenced by the work presented in (Arora et al.,
2019) for volumetric Michell trusses.

Recently, neural network-based implicit neural representation (INR) has become
very popular because of its capability to effectively solve problems in
reconstruction and optimization (De Luigi et al., 2023). INRs are usually
formulated as neural networks in a differentiable pipeline of computation and
determined by self-learning based optimization using a NN solver (e.g., Adams
(Kingma and Ba, 2015)). Example works include DeepSDF (Park et al., 2019) that
approximates continuous signed distance functions for surface reconstruction
and SIRENs (Sitzmann et al., 2020) that introduces periodic activation
functions for modeling complex natural signals and solving partial
differential equations.

With the help of INR, recent effort has been made to overcome the limitations
caused by the discrete nature of meshes in deformation. Methods such as
(Groueix et al., 2018, 2019; Huang et al., 2020; Jiang et al., 2020; Yang et
al., 2021) utilized the implicit vector fields to represent the deformation,
which can be applied to models in different representations. The optimized
deformation of a template model can be computed by optimizing INR according to
different objectives and constraints. This strategy has been adopted to solve
problems of non-rigid registration (Deng et al., 2021; Sundararaman et al.,
2022), surface reconstruction (Williams et al., 2019) and parameterization
(Morreale et al., 2021; Aigerman et al., 2022). However, the translation based
nature of implicit vector fields does not provide rotational-invariant
deformation. In this paper, we employ NNs as continuous functions for
quaternions and scaling ratios to determine optimized curved layers for multi-
axis 3D printing.

###  2.3. Mesh-based Deformation and Parameterization

Deforming a given model is usually formulated as an optimization problem
(e.g., (Sorkine and Botsch, 2009)). Earlier approaches mainly defined the
geometric structures according to the local connectivity, where a surface mesh
or a volumetric mesh with high quality is needed. To overcome this challenge,
cage-based deformation (e.g., (Yifan et al., 2020; Hu et al., 2015)) has been
employed to handle wild geometry such as scattered points or polygon soups.

We proposed to construct a ‘cage’-like tetrahedral mesh, which encloses the
input model but has completely different connectivity, as the intermediate
discretization for computing the mapping. The scale-controlled deformation
proposed in (Zhang et al., 2022) is employed in our work to compute the
mapping. Specifically, the as-rigid-as-possible (ARAP) energy (Sorkine and
Alexa, 2007; Sorkine et al., 2004) is defined as the function to be minimized
while allowing to scale with specified ratios along three orthogonal
directions. The ARAP energy can be efficiently minimized by a local/global
scheme (Bouaziz et al., 2012; Liu et al., 2008). Although the effort of
differentiation has been made in (Bächer et al., 2014; Romain et al., 2013),
ARAP deformation has not been employed in the NN-based computational pipeline
– especially with respect to the fields of quaternions and scaling ratios.

(a)(b)(c)(d)(e)(f)

Figure 3. An overview of our neural slicer to generate curved layers for
multi-axis 3D printing. (a) The input Yoga model $\mathcal{M}$ with its
implicit solid $H(\mathbf{x})$ and the distribution of principal stresses
obtained from voxel-based FEA. (b) A volumetric mesh $\mathcal{C}$ caging the
input model $\mathcal{M}$ is constructed, serving as the intermediate
representation in numerical computation. (c) Two continuous functions
$\mathbf{q}(\mathbf{x})$ and $\mathbf{s}(\mathbf{x})$ specify the quaternion
and the scaling ratios of local deformation for all
$\mathbf{x}\in\mathbb{R}^{3}$, which are represented as neural networks (NN)
to be optimized. (d) The function values of $\mathbf{q}(\mathbf{x})$ and
$\mathbf{s}(\mathbf{x})$ are sampled to drive a differential deformation to
obtain a deformed caging mesh $\mathcal{C}^{d}$ therefore also the mapping
$\lambda(\mathbf{q}(\mathbf{x}),\mathbf{s}(\mathbf{x}))$. (e) The scalar field
$G(\mathbf{x})$ is obtained from the mapping, with its gradient $\nabla
G(\mathbf{x})$ serving as local printing directions (LPDs) that play a pivotal
role in formulating the loss functions utilized for optimization. (f) After
computing an optimized $G(\mathbf{x})$, its isosurfaces are extracted on the
caging mesh $\mathcal{C}$ and trimmed by the implicit solid $H(\mathbf{x})$ to
obtain the curved layers.

Our work to compute an optimized mapping shares certain similarities with the
volumetric parameterization problems (Patané et al., 2013). Specifically, the
task of computing a mapping for the SF and SR demands can be viewed as a
specialized volumetric parameterization, incorporating orientation
constraints. Different from the large literature of parameterization that
mainly focus on resolving foldovers to ensure local injectivity (ref. (Du et
al., 2020; Fu et al., 2015; Garanzha et al., 2021; Kovalsky et al., 2015; Liao
et al., 2021; Schüller et al., 2013; Su et al., 2019; Rabinovich et al., 2017)
), local injectivity is not particularly required in our problem.

##  3\. Neural Optimization for Slicing

We introduce the computational pipeline of our Neural Slicer in this section.
First of all, a neural network-based formulation is introduced for the mapping
$\lambda(\cdot)$ to be optimized. After that, the steps of our slicing
algorithm are presented. An overview illustration is given in Fig.3.

###  3.1. Differentiable Mapping as Neural Network

The basic idea of our approach is to construct a differentiable mapping
$\lambda(\cdot)\in\mathbb{R}^{3}$ by neural networks. The scalar field for
layer generation can be defined by the $z$-component of $\lambda(\cdot)$ as
$G(\mathbf{x})$, where the function $\lambda(\cdot)$ is parameterized with the
network coefficients $\theta$. Therefore, the scalar field $G(\mathbf{x})$ for
curved layer generation can be optimized by self-learning according to the
manufacturing requirements directly defined on $G(\mathbf{x})$ and $\nabla
G(\mathbf{x})$ as loss functions.

Simply defining the mapping $\lambda(\cdot)$ by a neural network without any
geometric meaning leads to a learning process that converges very slowly. We
propose to formulate $\lambda(\cdot)$ as
$\lambda(\mathbf{q}(\mathbf{x}),\mathbf{s}(\mathbf{x}))$ by two continuous
functions $\mathbf{q}(\mathbf{x})\in\mathbb{R}^{4}$ and
$\mathbf{s}(\mathbf{x})\in\mathbb{R}^{3}$ that define the local rotation (as
quaternion) and the local scaling ratios at any point
$\mathbf{x}\in\mathbb{R}^{3}$. $\mathbf{q}(\mathbf{x})$ and
$s(\mathbf{x})\in\mathbb{R}^{3}$ are represented by two neural networks with
coefficients $\theta_{q}$ and $\theta_{s}$ respectively. $\lambda(\cdot)$ are
then computed from $\mathbf{q}(\mathbf{x})$ and $\mathbf{s}(\mathbf{x})$ in a
discrete manner. The space $\Omega$ around the input model $\mathcal{M}$ is
discretized into a volumetric mesh $\mathcal{C}$ with tetrahedral elements,
where $\mathcal{C}$ encloses the solid of $\mathcal{M}$ as a cage. The scale-
controlled ARAP deformation (Zhang et al., 2022) is employed111Note that the
deformation is computed in a space $\Omega$ instead of the input model to
determine the mapping. to compute a deformed caging mesh $\mathcal{C}^{d}$ so
that the quaternions and the scaling ratios can be approximated at every
element of $\mathcal{C}$. Here we use the centers of elements as the sample
points. With the help of barycentric coordinates, every point
$\mathbf{x}\in\Omega$ is mapped to a new position $\mathbf{y}\in\Omega^{d}$ in
the new space $\Omega^{d}$ defined around $\mathcal{C}^{d}$. A discrete form
of the mapping $\lambda(\cdot)$ is obtained, the $z$-component of which is the
scalar field $G(\mathbf{x})$. By deriving the formulas to compute $\partial
G/\partial\theta_{q}$ and $\partial G/\partial\theta_{s}$ w.r.t. the positions
of vertices on $\mathcal{C}^{d}$ (details can be found in Sec. 5.1), a mapping
$\lambda(\cdot)$ can be optimized via the differentiable pipeline of neural
networks.

###  3.2. Slicing Algorithm

The algorithm consists of three major stages: pre-processing (Steps 1-5),
mapping optimization (Steps 6-10), and post-processing (Step 11). The steps
are presented as follows (see the illustration in Fig.3):

  1. (1)

Construct the solid $\mathcal{H}$ as an implicit function $H(\mathbf{x})$ for
the input model $\mathcal{M}$;

  2. (2)

Apply voxel-based FEA to compute stress field inside the solid
$H(\mathbf{x})\leq 0$;

  3. (3)

Generate a tetrahedral mesh $\mathcal{C}$ as the intermediate discretization
of the space $\Omega$ enclosing $\mathcal{H}$ (Sec. 5.3);

  4. (4)

Build a set of sample points as $\mathcal{B}$ on the surface boundary of
$\mathcal{M}$, which is employed to evaluate the loss functions for SF
requirements;

  5. (5)

Build a set of sample points as $\mathcal{T}$ at the center of voxels for
evaluating the strength reinforcement loss in the relevant region – i.e., the
region of the top 10% maximal stresses;

  6. (6)

Initialize the neural networks for $\mathbf{q}(\mathbf{x})$ and
$\mathbf{s}(\mathbf{x})$;

  7. (7)

Compute a deformed mesh of $\mathcal{C}$ as $\mathcal{C}^{d}$ by the scale-
controlled ARAP deformation;

  8. (8)

Evaluate the differentiation of the deformation to prepare for backpropagation
(Sec. 5.1);

  9. (9)

Update the networks of $\mathbf{q}(\mathbf{x})$ and $\mathbf{s}(\mathbf{x})$
by the NN solver according to the loss functions (Sec. 4);

  10. (10)

Go back to step (7) until the learning converges;

  11. (11)

Extract the isosurfaces of $G(\mathbf{x})$ from the volumetric mesh
$\mathcal{C}$ and trim them by the implicit solid $H(\mathbf{x})\leq 0$ (Sec.
5.4).

As the mapping is computed on an intermediate mesh $\mathcal{C}$, our slicing
algorithm does not require a perfect volumetric mesh as the discretization of
the input model $\mathcal{M}$.

##  4\. Loss Functions

In order to optimize a mapping function that leads to curved layers satisfying
different manufacturing requirements, loss functions are defined in this
section for strength reinforcement, support-free, and collision avoidance. All
loss functions are defined according to the LPDs – i.e.,
$\mathbf{d}_{p}={\nabla G(\mathbf{x})}/{\|\nabla G(\mathbf{x})\|}$ w.r.t. the
scalar field $G(\mathbf{x})$ to be optimized. Harmonic losses of
$\mathbf{q}(\mathbf{x})$ and $\mathbf{s}(\mathbf{x})$ are introduced for the
purpose of regularization.

###  4.1. Strength Reinforcement

Anisotropic mechanical behavior is a prominent characteristic observed in
models fabricated through filament-based 3D printing, primarily due to the
limited adhesion between deposited filaments (Riddick et al., 2016). Under
specific load conditions, the mechanical strength of printed models can be
significantly enhanced when the filaments are aligned along the distribution
of principal stresses (Tam and Mueller, 2017b). Specifically, the 3D printed
model can have a reinforced strength when the LPDs are nearly perpendicular to
directions of the maximal principal stresses $\mathbf{\tau}_{\max}$ as
$|\mathbf{d}_{p}\cdot\mathbf{\tau}_{\max}|\leq\sin\beta$ (ref. (Fang et al.,
2020)). According to this criterion, we define the following loss function
evaluated at the sample points $\mathbf{p}\in\mathcal{T}$ (i.e., the centers
of voxels in the relevant region) for the purpose of strength reinforcement.

(2)
$\mathcal{L}_{SR}:=\sum_{\mathbf{p}\in\mathcal{T}}|V_{e}|~{}\sigma\left(k_{SR}(|\mathbf{d}_{p}\cdot\mathbf{\tau}_{\max}(\mathbf{p})|-\sin\beta)\right)$

where $\sigma(\cdot)$ is the sigmoid function and $|V_{e}|$ is the volume of
each voxel element for FEA. $\beta$ is a parameter that introduces a level of
tolerance and is related to material properties and printing temperature.
$k_{SR}$ is the logistic growth rate to tune the ‘steepness’ of the loss
function. By experimental tests, $\beta=10^{\circ}$ and $k_{SR}=15$ are
employed for all examples in this paper.

###  4.2. Support-Free

Two different criteria are considered to define the loss functions to optimize
the curved layers to enable support-free (SF) 3D printing based on 1) surface
normal and 2) point extremity respectively. All are evaluated on the sample
points $\mathbf{p}\in\mathcal{B}$. Note that $\mathcal{B}$ can exclude points
in the regions that have been naturally supported by the printing platform
(e.g., those on the bottom plate of the Bunny Head model).

Given an implicit function $H(\mathbf{x})$, the surface normal at a sample
point $\mathbf{p}\in\mathcal{B}$ can generally be evaluated by the implicit
function as $\mathbf{n}(\mathbf{p})=-\nabla H(\mathbf{p})/\|\nabla
H(\mathbf{p})\|$. Specifically, $\mathbf{n}(\mathbf{p})$ can be precisely
obtained from the input model or the skeleton of the convolution solid
(including the shell, the cylindrical solid, and the tubular solid) as shown
above. In the context of multi-axis 3D printing, the SF requirement for LPDs
and surface normals is defined according to the local self-supporting angle
$\alpha$ as $-\mathbf{n}(\mathbf{p})\cdot\mathbf{d}_{p}\leq\sin\alpha$ (ref.
(Zhang et al., 2022)). Therefore, we define the following loss function.

(3)
$\mathcal{L}_{SF}:=\sum_{\mathbf{p}\in\mathcal{B}}|A_{\mathbf{p}}|~{}\sigma\left(k_{SF}(-\mathbf{n}(\mathbf{p})\cdot\mathbf{d}_{p}-\sin\alpha)\right)$

where $\sigma(\cdot)$ is the sigmoid function and $k_{SF}=30$ was used as the
logistic growth rate for this loss function. $|A_{\mathbf{p}}|$ is a weight
reflecting the surface area covered by a sample point, which can be estimated
by the average squared distance between $\mathbf{p}$ and its $k$-nearest
neighbors. We use $k=10$ for all our examples.

Similar to the scenario of point overhang discussed in (Vanek et al., 2014),
we need to prevent making a surface point as local minimal with reference to
the LPD $\mathbf{d}_{p}$. Specifically, giving $\mathcal{N}_{\mathbf{p}}$ as
the set of $\mathbf{p}$’s $k$-nearest neighbors, $\mathbf{p}$ is a local
minimum (i.e., point overhang) when

$(\mathbf{p}_{j}-\mathbf{p})\cdot\mathbf{d}_{p}>0$
($\forall\mathbf{p}_{j}\in\mathcal{N}_{\mathbf{p}}$).

The following loss is defined to avoid generating LPDs that lead to local
‘minimum’ on boundary surfaces (i.e., point overhang).

(4)
$\mathcal{L}_{PO}:=\sum_{\mathbf{p}\in\mathcal{B}}|A_{\mathbf{p}}|\max(0,\min_{\mathbf{p}_{j}\in\mathcal{N}_{\mathbf{p}}}((\mathbf{p}_{j}-\mathbf{p})\cdot\mathbf{d}_{p}))$

where $\min(\cdot)$ function is implemented by the min pooling layer and
$\max(0,\cdot)$ can be easily realized by the ReLU function. Again,
$|A_{\mathbf{p}}|$ is the area weight of a sample point.

###  4.3. Collision Avoidance

As more DOFs of motion are enabled for the process of multi-axis 3D printing,
the collision can happen locally between the working surface and the printer
head if the curved layers are not designed properly. As shown on the right,
the geometry of a printer head forms a cone shape with the apex angle
$\varphi$. When using this cone to touch a point on the working surface, local
collision will occur in any concave region with a dihedral angle less than
$\varphi$.

When the scalar field $G(\mathbf{x})$ is defined in a piecewise linear manner
on the elements of the caging mesh $\mathcal{C}$, the normal of isosurfaces in
an element $e$ is a constant vector as $\mathbf{n}_{e}=\nabla
G(\mathbf{x})/\|\nabla G(\mathbf{x})\|$. The concavity of isosurfaces in two
neighboring elements can be evaluated by
$(\mathbf{n}_{L}\times\mathbf{n}_{R})\cdot\mathbf{h}$, where $\mathbf{n}_{L}$
and $\mathbf{n}_{R}$ are the normal vectors in the left and the right elements
$e_{L}$ and $e_{R}$ and h is the unit vector obtained from the isoline of
$G(\mathbf{x})$ defined on the triangle shared by two elements. Borrowing the
study of gouging (i.e., local collision) in multi-axis CNC machining (Bartoň
et al., 2021), the local concavity needs to satisfy the following condition to
be collision-free

$\left\\{\begin{array}[]{rl}-\sin{\varphi}\leq(\mathbf{n}_{L}\times\mathbf{n}_{R})\cdot\mathbf{h}<0&(\varphi\leq\frac{\pi}{2}),\\\
(\mathbf{n}_{L}\times\mathbf{n}_{R})\cdot\mathbf{h}\leq-\sin{\varphi}&(\varphi>\frac{\pi}{2}).\end{array}\right.$

This can be further formulated into a loss as

(5) $\displaystyle\mathcal{L}_{CA}:=\sum_{(e_{L},e_{R})\in\mathcal{N}_{e}}$
$\displaystyle\max(0,(\mathbf{n}_{L}\times\mathbf{n}_{R})\cdot\mathbf{h})$
$\displaystyle+\max(0,-(\mathbf{n}_{L}\times\mathbf{n}_{R})\cdot\mathbf{h}-\sin(\varphi))$

for a relative sharp printer head with $\varphi\leq\frac{\pi}{2}$. Or we
define the loss function as follows

(6)
$\mathcal{L}_{CA}:=\sum_{(e_{L},e_{R})\in\mathcal{N}_{e}}\max(0,(\mathbf{n}_{L}\times\mathbf{n}_{R})\cdot\mathbf{h}+\sin(\varphi))$

for a relatively flat printer head with $\varphi>\frac{\pi}{2}$. Again,
$\max(0,\cdot)$ is realized by the ReLU function in our implementation.

Different from previous loss functions, $\mathcal{L}_{CA}=0$ needs to be
processed as a hard constraint to ensure collision-free. Details can be found
in Sec. 5.2. This collision avoidance loss only considers the local collision.
Global collision is checked and resolved during the motion planning process of
the physical realization (ref. (Zhang et al., 2021; Ezair et al., 2018)),
which is not the focus of this paper.

###  4.4. Harmonics

Harmonic losses are employed to control the smoothness of
$\mathbf{q}(\mathbf{x})$ and $\mathbf{s}(\mathbf{x})$ to be learned. Let
$\mathcal{N}_{e}$ define the set of neighboring element pairs in the caging
mesh $\mathcal{C}$, the differences of scale ratios and quaternions are
measured on pairs of neighboring elements.

First, a loss is introduced to minimize the difference of scaling ratios
$\mathbf{s}(\cdot)$ so that avoids the radical change of layer thickness in
neighboring regions.

(7)
$\mathcal{L}_{HS}:=\sum_{(e_{i},e_{j})\in\mathcal{N}e}\frac{1}{2}(|V_{e_{i}}|+|V_{e_{j}}|)\|\mathbf{s}(\mathbf{c}_{e_{i}})-\mathbf{s}(\mathbf{c}_{e_{j}})\|^{2}$

where $\mathbf{c}_{e_{i}}$ and $\mathbf{c}_{e_{j}}$ are the centers of the
neighboring elements $e_{i},e_{j}$, and $|V_{e_{i}}|$ and $|V_{e_{j}}|$ are
the volumes of the elements. The scaling ratios are evaluated on the NN-based
continuous function $\mathbf{s}(\cdot)$.

Furthermore, it is important to generate curved layers with smoothly varied
normals to ensure smooth and continuous motion of the printer head during
fabrication. We control this property of resultant curved layers by
introducing a loss that evaluates the difference of quaternions in neighboring
regions.

(8)
$\mathcal{L}_{HQ}:=\sum_{(e_{i},e_{j})\in\mathcal{N}e}\frac{1}{2}(|V_{e_{i}}|+|V_{e_{j}}|)\left(1-\frac{\mathbf{q}(\mathbf{c}_{e_{i}})}{\|\mathbf{q}(\mathbf{c}_{e_{i}})\|}\cdot\frac{\mathbf{q}(\mathbf{c}_{e_{j}})}{\|\mathbf{q}(\mathbf{c}_{e_{j}})\|}\right)^{2}$

where $\mathbf{q}(\cdot)$ is the NN-based continuous function of quaternions.

Among all metrics for 3D rotations (Huynh, 2009), we select one with a simple
form that is normalized to bound the loss value.

###  4.5. Total Loss

In summary, the curved layer generation problem can be formulated as a
learning-based optimization process to minimize the total loss:

(9)
$\mathcal{L}=w_{1}\mathcal{L}_{SF}+w_{2}\mathcal{L}_{SR}+w_{3}\mathcal{L}_{OP}+\mathcal{L}_{HS}+\mathcal{L}_{HQ},$

where $w$’s are the balancing weights for different loss terms. Considering
the collision avoidance loss as a hard constraint, that is

(10) $\arg\min_{\theta_{q},\theta_{s}}\mathcal{L}\quad
s.t.\;\mathcal{L}_{CA}=0.$

##  5\. Implement Details

###  5.1. Differentiable Deformation as Mapping

In our current implementation, the mapping $\lambda(\cdot)$ is computed on the
volumetric caging mesh $\mathcal{C}$ by using the scale-controlled ARAP
deformation (Zhang et al., 2022). That is to determine the new position
$\mathbf{v}^{d}$ for every vertex $\mathbf{v}\in\mathcal{C}$ so that the
quaternions and the scale ratios given by the functions $\textbf{q}(\cdot)$
and $\textbf{s}(\cdot)$ are satisfied in the elements of $\mathcal{C}$ via an
optimization formulation. We provide the details of this deformation below and
then analyze its differentiability.

For every element $e\in\mathcal{C}$, we determine its rotation matrix
$\mathbf{R}_{e}$ by the quaternion $\mathbf{q}(\mathbf{c}_{e})$ at the
element’s center $\mathbf{c}_{e}$. Its scale matrix
$\mathbf{S}_{e}=\mathrm{diag}(\mathbf{s}^{x}(\mathbf{c}_{e}),\mathbf{s}^{y}(\mathbf{c}_{e}),\mathbf{s}^{z}(\mathbf{c}_{e}))$
is also computed by $\textbf{s}(\cdot)$ at the center $\mathbf{c}_{e}$. Then,
the locally scaled and rotated element $e$ can be computed by
$\mathbf{R}_{e}\mathbf{S}_{e}(\mathbf{NV}_{e})^{\mathrm{T}}$ with
$\mathbf{V}_{e}$ being a position matrix formed by the coordinates of $e$’s
four vertices. $\mathbf{N}$ is used to transfer an element’s center to the
origin as $\mathbf{N}=\mathbf{I}_{4\times 4}-\frac{1}{4}\mathbf{1}_{4\times
4}$. We compute the deformed caging mesh $\mathcal{C}^{d}$ as

(11)
$\arg\min_{\mathcal{C}^{d}}\underbrace{\sum_{e\in\mathcal{C}}\|(\mathbf{NV}_{e}^{d})^{\mathrm{T}}-\mathbf{R}_{e}\mathbf{S}_{e}(\mathbf{NV}_{e})^{\mathrm{T}}\|_{F}^{2}}_{\mathrm{Position-
Compatibility}}+\underbrace{\gamma\sum_{\mathbf{v}\in\mathcal{C}}\|\mathbf{v}^{d}-\mathbf{v}\|^{2}}_{\mathrm{Regularization}}$

where $\|\cdot\|_{F}$ is the Frobenius norm, and the positions of vertices in
$\mathcal{C}^{d}$ to be determined are kept in the position matrix
$\mathbf{V}_{e}^{d}$. The regularization term is added to improve the
conditional number of the linear system for solving this problem.

Denoting $\mathbf{\xi}^{0}=\\{\mathbf{v}\\}$ and
$\mathbf{\xi}=\\{\mathbf{v}^{d}\\}$ as the vectors containing the positions of
all vectors before and after deformation, the above equation is actually in
the regularized least-squares form as

(12)
$\arg\min_{\mathbf{\xi}}\|\mathbf{A}\mathbf{\xi}-\mathbf{b}\|^{2}+\gamma\|\mathbf{\xi}-\mathbf{\xi}^{0}\|^{2}.$

with $\mathbf{A}$ and $\mathbf{b}$ derived from the the position compatibility
term of Eq.(11). Then, the positions of the deformed mesh are obtained as

(13)
$\mathbf{\xi}=(\mathbf{A}^{T}\mathbf{A}+\gamma\mathbf{I})^{-1}(\mathbf{A}^{T}\mathbf{b}+\gamma\mathbf{\xi}^{0})$

with $\mathbf{I}$ as an identity matrix.

Given a query point $\mathbf{x}\in\mathbb{R}^{3}$ contained by an element of
$\mathcal{C}$ with four vertices $\mathbf{v}_{i}$ ($i=1,\ldots,4$), we can
determine $\mathbf{x}$’s barycentric coordinate as $a_{i}(\mathbf{x})$ with
$\sum_{i}a_{i}(\mathbf{x})\equiv 1$ and $a_{i}(\mathbf{x})\in[0,1]$. As a
result, the mapping $\lambda(\cdot)$ is obtained as
$\mathbf{y}=\sum_{i}a_{i}(\mathbf{x})\mathbf{v}^{d}_{i}$. The scalar field
$G(\mathbf{x})$ as the $z$-components of $\mathbf{y}$ is

(14)
$G(\mathbf{x})=\mathrm{diag}(0,0,1)\sum_{i}a_{i}(\mathbf{x})\mathbf{v}^{d}_{i}.$

This scalar field $G(\mathbf{x})$ determined by the NN-based mapping is
differentiable with reference to the network coefficients $\theta_{q}$ and
$\theta_{s}$. Specifically, we have

(15) $\frac{\partial G}{\partial\theta_{q}}=\sum_{i}\frac{\partial
G}{\partial\mathbf{v}^{d}_{i}}\left(\frac{\partial\mathbf{v}^{d}_{i}}{\partial\mathbf{q}}\frac{d\mathbf{q}}{d\theta_{q}}\right)\;\mathrm{and}\;\frac{\partial
G}{\partial\theta_{s}}=\sum_{i}\frac{\partial
G}{\partial\mathbf{v}^{d}_{i}}\left(\frac{\partial\mathbf{v}^{d}_{i}}{\partial\mathbf{s}}\frac{d\mathbf{s}}{d\theta_{s}}\right),$

among which ${d\mathbf{q}}/{d\theta_{q}}$ and ${d\mathbf{s}}/{d\theta_{s}}$
can be directly obtained from the neural networks representing
$\mathbf{q}(\mathbf{x})$ and $\mathbf{s}(\mathbf{x})$. By using the calculus
of matrix, we can also have

(16)
$\frac{\partial\mathbf{\xi}}{\partial\mathbf{s}}=(\mathbf{A}\mathbf{\xi}^{0})^{T}\otimes((\mathbf{A}^{T}\mathbf{A}+\mathbf{\gamma\mathbf{I}})^{-1}\mathbf{A}^{T}\mathbf{R})$
(17)
$\frac{\partial\mathbf{\xi}}{\partial\mathbf{q}}=\frac{\partial\mathbf{\xi}}{\partial\mathbf{R}}\frac{\partial\mathbf{R}}{\partial\mathbf{q}}=(\mathbf{s}\mathbf{A}\mathbf{\xi}^{0})^{T}\otimes((\mathbf{A}^{T}\mathbf{A}+\gamma\mathbf{I})^{-1}\mathbf{A}^{T})\frac{\partial\mathbf{R}}{\partial\mathbf{q}}$

where $\otimes$ denotes the Kronecker product and $\mathbf{R}$ is the rotation
matrix form of $\mathbf{q}$.
${\partial\mathbf{v}^{d}_{i}}/{\partial\mathbf{s}}$ and
${\partial\mathbf{v}^{d}_{i}}/{\partial\mathbf{q}}$ can then be obtained from
${\partial\mathbf{\xi}}/{\partial\mathbf{s}}$ and
${\partial\mathbf{\xi}}/{\partial\mathbf{q}}$ from the corresponding location
of $\mathbf{v}^{d}$ in $\xi$.

By the chain rule, the total derivatives from the loss function $\mathcal{L}$
to the network coefficients $\theta_{q}$ and $\theta_{s}$ are given as follows

(18) $\frac{d\mathcal{L}}{d\theta_{q}}=\frac{\partial\mathcal{L}}{\partial
G}\left(\sum_{i}\frac{\partial
G}{\partial\mathbf{v}^{d}_{i}}\left(\frac{\partial\mathbf{v}^{d}_{i}}{\partial\mathbf{q}}\frac{d\mathbf{q}}{d\theta_{q}}\right)\right)+\frac{\partial\mathcal{L}}{\partial\mathbf{q}}\frac{d\mathbf{q}}{d\theta_{q}},$
(19) $\frac{d\mathcal{L}}{d\theta_{s}}=\frac{\partial\mathcal{L}}{\partial
G}\left(\sum_{i}\frac{\partial
G}{\partial\mathbf{v}^{d}_{i}}\left(\frac{\partial\mathbf{v}^{d}_{i}}{\partial\mathbf{s}}\frac{d\mathbf{s}}{d\theta_{s}}\right)\right)+\frac{\partial\mathcal{L}}{\partial\mathbf{s}}\frac{d\mathbf{s}}{d\theta_{s}}.$

The backpropagation of our NN-based computational pipeline can be realized. In
practical implementation, the total derivatives can be computed by automatic
differentiation (Paszke et al., 2019).

###  5.2. Network Architecture and Constrained Optimization

An architecture of Sinusoidal Representation Network (SIREN) is adopted for
$\mathbf{q}(\mathbf{x})$ and $\mathbf{s}(\mathbf{x})$ that contains $10$
hidden layers, each of which is formed by $512$ neurons with the periodic
activation functions. The hyperparameters of this network are chosen by the
strategy introduced in (Sitzmann et al., 2020). The input of the networks is
$\mathbf{x}\in\mathbb{R}^{3}$, and the output is a four-dimensional vector for
$\mathbf{q}(\mathbf{x})\in\mathbb{R}^{4}$ and a three-dimensional vector as
scaling ratios for $\mathbf{s}(\mathbf{x})\in\mathbb{R}^{3}$.

When computing the network coefficients according to the total loss, the Adam
optimizer (Kingma and Ba, 2015) is employed due to its efficacy in handling
large-scale data and complex computation. The initial learning rate is set to
1.0e-3, with a minimum learning rate threshold established at 1.0e-6. To
adaptively adjust the learning rate during the optimization process, we
utilize the ‘ReduceLROnPlateau‘ learning rate scheduler that can dynamically
alter the learning rate based on the optimization progress.

To implement the optimization with hard constraints in our NN-based
computational pipeline, we adopt the DC3 framework as introduced in (Donti et
al., 2021). Specifically, we implement a gradient-based correction step before
the loss calculation and back-propagation. This involves taking gradient steps
in the parameter space $(\theta_{q},\theta_{s})$ towards the feasible region
with $\mathcal{L}_{CA}=0$, thereby aligning solutions closer to the feasible
region.

Table 1. Computational statistics for Neural Slicer

Model | Fig. | Objectives |  Model Genus |  Cage Genus |  # of Cage Elements |  Pre-processing Time† (sec.) |  Optimization Time (sec.) |  Post-processing Time (sec.) |  Slicer Comp. Time (sec.) |  Toolpath Gen. Time (sec.)
---|---|---|---|---|---|---|---|---|---|---
Bunny Head | 1 | SF+SR | 22 | 0 | $23,139$ |  $93.5$ ($60.1$) | $553.2$ | $124.1$ | $677.3$ | $52.7$
Yoga | 3 | SF+SR | 2 | 1 | $6,279$ |  $55.1$ ($34.1$) | $312.3$ | $42.5$ | $354.8$ | $20.6$
Shelf | 6 | SF+SR | 30 | 26 | $42,150$ |  $72.4$ ($52.3$) | $647.5$ | $104.1$ | $65.1$ | $65.1$
Ring | 7 | SF | 1 | 1 | $5,264$ | $33.4$ | $73.4$ | $15.4$ | $88.8$ | $2.3$
Tubes | 8 | SF | 5 | 4 | $16,551$ | $21.4$ | $153.0$ | $45.7$ | $3.5$ | $49.2$
Spiral Fish | 12 | SF | 3 | 3 | $43,147$ | $25.3$ | $223.5$ | $62.0$ | $285.5$ | $29.4$
Bridge | 9 | SR | 21 | 21 | $47,520$ |  $52.1$ ($30.6$) | $174.3$ | $54.8$ | $229.1$ | $35.8$

† Including the computing time for voxel-based FEA as given in the brackets.

###  5.3. Pre-processing: Sampling and Cage Generation

As discussed above, a tetrahedral mesh $\mathcal{C}$ caging the input model
$\mathcal{M}$ is constructed in our approach as the intermediate
representation to determine the mapping $\mathbf{\lambda}(\cdot)$ according to
$\mathbf{q}(\mathbf{x})$ and $\mathbf{s}(\mathbf{x})$. This representation for
numerical computation is selected because of its good balance between the
computational efficiency and the fidelity of deformation, making it highly
suitable for our NN-based slicing framework. Meanwhile, the loss function of
support-free is evaluated on a set of surface sampling points $\mathcal{B}$.
The methods for generating $\mathcal{B}$ and $\mathcal{C}$ are introduced
below.

(a)(b)(c)

Figure 4. The illustration of cage generation: (a) the implicit surface
$H(\mathbf{x})$ (gray) of the input model $\mathcal{M}$ in an abstract
representation (e.g., as a convolution surface of the skeletons illustrated in
black curves), (b) the polygonal mesh of $H(\mathbf{x})=0$ (yellow) and the
caging surface mesh (black), (c) the tetrahedral mesh (purple) caging the
input model.

(a)(b)(c)$G(\mathbf{x})$$H(\mathbf{x})\leq
0$$\mathcal{G}_{i}$$\mathcal{P}_{i}$

Figure 5. The illustration of slicing process to generate the curved layers:
(a) a scalar field $G(\mathbf{x})$ that is optimized on the volumetric caging
mesh $\mathcal{C}$, (b) the isosurfaces of $G(\mathbf{x})$ as polygonal
surface meshes $\\{\mathcal{G}_{i}\\}$ are extracted from $\mathcal{C}$’s
tetrahedra – the colors on an isosurface visualize the function values of
$H(\mathbf{x})$, and (c) the curved layers $\\{\mathcal{P}_{i}\\}$ are
obtained by trimming $\\{\mathcal{G}_{i}\\}$ with the implicit solid
$H(\mathbf{x})\leq 0$.

After converting all different representations of the input model
$\mathcal{M}$ into an implicit function $H(\mathbf{x})$ (see Fig.4(a)), a
polygonal mesh $\tilde{\mathcal{H}}$ of its zero level-set surface
$H(\mathbf{x})=0$ is extracted by the Marching Cubes algorithm (Lorensen and
Cline, 1987). Surface points are sampled on $\tilde{\mathcal{H}}$ and
projected onto $H(\mathbf{x})=0$ by using the gradient $\nabla H(\mathbf{x})$.
Sample points in $\mathcal{B}$ are obtained.

A surface caging mesh can then be generated by using the Nested Cage approach
(Sacht et al., 2015), which encloses all regions with $H(\mathbf{x})\leq 0$
(see Fig.4(b)). Lastly, the volumetric caging mesh $\mathcal{C}$ is
constructed from the surface caging mesh by Tetgen (Si, 2015) (see Fig.4(c)).
The mapping and the curved layers are computed on this mesh $\mathcal{C}$.

###  5.4. Post-processing: Slicing on Cage

The resultant curved layers are computed from the scalar field $G(\mathbf{x})$
that is defined on the volume caging mesh $\mathcal{C}$ (see Fig.5(a)). The
isosurfaces of $G(\mathbf{x})$ are firstly extracted from the tetrahedra of
$\mathcal{C}$ as piece-wise linear polygonal meshes $\\{\mathcal{G}_{i}\\}$
(ref. (Treece et al., 1999)). The function values of $H(\mathbf{x})$ can then
be evaluated at every point on a surface mesh $\mathcal{G}_{i}$ (see the
values visualized as colors in Fig.5(b)). Lastly, the polygons of
$\\{\mathcal{G}_{i}\\}$ are trimmed by the implicit solid $H(\mathbf{x})\leq
0$ to form the curved layers $\\{\mathcal{P}_{i}\\}$. Toolpaths are generated
on each $\mathcal{P}_{i}$ and converted into trajectories of robot motion by
the method of $S^{3}$-Slicer (Zhang et al., 2022).

##  6\. Result and Discussion

Our computational pipeline is implemented in Python together with C++. PyTorch
(Paszke et al., 2019) is employed for constructing neural networks, automatic
differentiation, and solving linear equation systems. PyVista (Sullivan and
Kaszynski, 2019) is employed for mesh processing. Time-consuming steps such as
the point-to-surface distance query by (Gottschalk et al., 1996), the surface
cage generation by (Sacht et al., 2015), and the toolpath generation by (Zhang
et al., 2022) are based on C++ implementation. The source code of our Neural
Slicer will be released upon the acceptance of this paper.

(a)(b)(c)(d)DeformedCage $\mathcal{C}^{d}$

Figure 6. The result generated by our Neural Slicer vs. the result by a planar
slicer on the Shelf model: (a) the stress field under the given forces (shown
as the arrows) and the model’s cage used in the computation, (b) the curved
layers (bottom) generated from the mapping determined by a deformed cage
(top), (c) the histograms for evaluating the quality of results in terms of
(top) the angles between LPDs and surface normals for the SF requirement and
(bottom) the angles between LPDs and the maximal stresses for the SR
requirement, (d) the results of FEA simulation by using anisotropic material
orientations defined according to the LPDs for planar layers (top) and our
curved layers (bottom).

###  6.1. Computational Experiments

All the computational experiments are conducted on a desktop PC with an
Intel(R) Core(TM) i5-12600K CPU (10 cores @3.6GHz), NVIDIA RTX 4080 GPU, and
32GB RAM, running Ubuntu 20.04 LTS (Focal Fossa).

####  6.1.1. Examples and computational statistics.

We have tested our Neural Slicer on a variety of models with complicated
geometry and topology. The first model is a Bunny Head model in hybrid
representations (see Fig.1(a)), including $15,635$ tetrahedra for solid
(blue), $3,495$ triangles for shell (green), $41$ skeletons (thin black lines)
for the struts as cylindrical solids and $36$ skeletons (bold black lines) for
the tubular solids. The genus number of this input model is $g=22$. The
implicit function $H(\textbf{x})$ of this model is computed with the help of
the signed distance function and the convolution surfaces. Both the support-
free (SF) and the strength reinforcement (SR) requirements are applied to this
model. Curved layers can be successfully computed by our Neural Slicer (see
Fig.1(b)). The second example is the Yoga model with both the SF and SR
objectives required as shown in Fig.3, where the input is a mesh surface with
$13,938$ triangles. The third example is a Shelf model generated by topology
optimization. It is a model with a high genus number as $g=30$ and again both
the SF and the SR requirements are applied. The resultant curved layers are
shown in Fig.6. We can find that the curved layers follow the directions of
maximal stresses very well.

We also tested our Neural Slicer on a few other models by applying only the SF
requirement, including the Ring model (Fig.7), the Tubes model (Fig.8) and the
Spiral Fish model (Fig.12). The Ring model is selected to compare with the
curved layers generated by $S^{3}$-Slicer. As illustrated in Fig.2,
$S^{3}$-Slicer defined the SF objective in the deformed space. Layers that are
self-supported in the deformed space can become overhangs in the model space
caused by the mapping distortion. Differently, our Neural Slicer directly
evaluates the SF loss in the model space so that the overhang can be better
prevented (see the comparison shown in Fig.7). The Tubes model is an implicit
solid generated by convolution surfaces from the representation of $186$ line
segments as skeletons. Our slicer can generate curved layers directly from the
implicit solid. The Spiral Fish model is employed to demonstrate that our
approach is robust to different initial guesses. Details will discussed later
in Sec. 6.1.5.

Lastly, we have tested this Neural Slicer on the Bridge example and compared
our approach with the $S^{3}$-Slicer (see Fig.9). Only the SR requirement is
applied in this example. Computational statistics of all examples are given in
Table 1, where the computation of our Neural Slicer can be completed within 15
minutes. It can be observed that the genus number of a caging mesh can be
different from an input model. In other words, we can use a caging mesh with
simple topology as the computational domain to generate curved layers for a
model with complicated topology (see the Bunny Head model and its cage in
Fig.1 for an example). The toolpaths on the curved layers are generated by the
method presented in (Zhang et al., 2022), the computing time of which has also
been reported in Table 1.

(a)(b)

Figure 7. Comparisons of the results generated by (a) $S^{3}$-Slicer and (b)
our Neural Slicer for the Ring model.

####  6.1.2. Statistics in LPDs

As both the SF and the SR requirements are defined according to the LPDs, we
evaluate the quality of results generated by our Neural Slicer as the
histograms in terms of i) the angles between LPDs and the surface normals at
$H(\textbf{x})=0$ for the SF requirement and ii) the angles between LPDs and
the maximal stresses for the SR requirement. The histograms are generated by
LPDs (i.e., $\nabla G(\mathbf{x})/\|\nabla G(\mathbf{x})\|$) evaluated at all
the sample points in the set $\mathcal{B}$ for SF and the set $\mathcal{T}$
for SR. We have given these statistical results as histograms for the Bunny
Head model in Fig.1 and the Shelf model in Fig.6. It can be observed from the
histograms that both the SF and the SR requirements have been well achieved on
the curved layers generated by our Neural Slicer. Despite this, supporting
structures are still needed to help mount the Shelf model tightly on the
platform during the printing process.

The statistics for the curved layers generated by the $S^{3}$-slicer (Zhang et
al., 2022) are also given on the Tubes model (Fig.8) and the Bridge model
(Fig.9) for the purpose of comparison. It is found that our Neural Slicer can
further reduce the remaining overhang region by $95\%$ on the Tubes model. The
results are also visualized by the histogram shown in Fig.8(c). Note that the
result of the Tubes model is computed by $S^{3}$-slicer on a very dense
tetrahedral mesh with 430k elements while our Neural Slicer only employs a
caging mesh with 16.5k elements.

(a)(b)(c)

Figure 8. The results of the Tubes model that is generated by convolution
surface from a skeleton representation (see Fig.4(a)): the curved layers
generated by (a) our method and (b) $S^{3}$-Slicer. (c) The histograms of the
angles between LPDs and the surface normals are used to visualize the quality
of the SF objective achieved.

(a)(b)(c)(d)(e)(f)(g)(h)Max. Strain: 7.9e-2Max. Strain: 4.6e-2

Figure 9. Comparing our results with the result of $S^{3}$-Slicer on the
Bridge model for the SR requirement: (a) the input model and the stress field
under the given loading, (b) the histogram to visualize the results in terms
of angles between LPDs and maximal principal stresses, (c & e) the deformed
shape obtained by $S^{3}$-Slicer and its resultant curved layers, (d & f) the
deformed caging mesh generated by our Neural Slicer and the corresponding
curved layers. Our approach can find a very ‘smart’ solution that has nearly
planar layers that align with the distribution of maximal stresses better. To
further verify the mechanical strength, the anisotropic FEA is conducted to
generate strain distribution for the results of (g) $S^{3}$-Slicer and (h) our
Neural Slicer.

When applying our Neural Slicer on the Bridge model for the SR requirement
(see Fig.9), we can generate a result of curved layers that outperform those
obtained from the $S^{3}$-slicer (Zhang et al., 2022). Instead of the highly
curved layers, our slicer automatically generates a very ‘smart’ solution with
nearly planar layers (see Fig.9(d & f)), which surprisingly aligns better with
the maximal stresses – see the histogram given in Fig.9(b).

####  6.1.3. Verification by FEA

To verify the mechanical strength of models to be fabricated by different
layers, we conducted the FEA simulation with anisotropic material properties
by assigning different Young’s modulus along different directions at the
element level. Specifically, $Y_{1}=3.5~{}\mathrm{GPa}$ is used as the
strongest modulus and is assigned to the toolpath’s tangential direction. The
weakest modulus is assigned both to the surface normal direction of each layer
and the third orthogonal direction as $Y_{2}=Y_{3}=1.2~{}\mathrm{GPa}$.

FEA results of the Shelf model have been given in Fig.6(d), where the model
with curved layers generated by our Neural Slicer has reduced the maximal
strain by $43.3\%$ compared to the model using planar layers. For the FEA
simulation conducted on the Bridge model, the Neural Slicer can reduce the
maximal strain by $40.5\%$ than the $S^{3}$-Slicer (see Fig.9(g & h)). In
Fig.10, we compare the FEA results on the Bunny Head model for the results
with only the SF requirement vs. with both the SF and the SR requirements.
Adding the SR requirement can reduce the maximal strain by $36.8\%$.

(a)(b)(c)(d)Max. Strain: 3.8e-2Max. Strain: 2.4e-2

Figure 10. The results of FEA simulation with anisotropic material properties
on the Bunny Head model using curved layers generated by our Neural Slicer
with (a & b) only the SF requirement and (c & d) the SF + SR requirements.

(a)(b)(c)(d)

Figure 11. Ablation study of the point overhang loss $\mathcal{L}_{PO}$ on the
Bunny Head model: (a & c) illustrate the deformed cages $\mathcal{C}^{d}$
without / with $\mathcal{L}_{PO}$, (b & d) show the printing order of slices
by colors (from pink to gray). Printing sequence that starts from the ear
region is not manufacturable as the materials will be deposed onto the ‘air’.

####  6.1.4. Ablation Study for Point Overhang Loss

To demonstrate the effectiveness of the point overhang loss $\mathcal{L}_{PO}$
(Eq.(4)), we conducted an ablation study on the Bunny Head model imposing both
the SF and the SR requirements. The only difference is whether the point
overhang loss is omitted. When skipping the point overhang loss, the
optimization will lead to a deformed caging mesh $\mathcal{C}^{d}$ as shown in
Fig.11(a). Compared to the result with the full set of losses (as Fig.11(c)),
this omission leads to an overhang at the tip of the right ear. As a result,
the printing sequence determined according to $\mathcal{C}^{d}$ in Fig.11(b)
will require additional support at the tip of the right ear. Differently, the
printing sequence given in Fig.11(d) can be realized in a completely support-
free way.

####  6.1.5. Robustness to Initial Guess

One major issue of previous approaches based on nonlinear optimization is how
to obtain good initial guesses for LPDs. As can be found in the example of the
Tubes model shown in Fig.8, the $S^{3}$-Slicer computes results by using the
distance field generated by the heat method (Crane et al., 2017) as an initial
guess. Apparently, the resultant curved layers are less smoother than those
generated by our Neural Slicer (see the zoom views of Fig.8). Another example
is the Bridge model as shown in Fig.9. When posing the input model into an
orientation as Fig.9(d) and using the height field as the initial guess for
LPDs, the $S^{3}$-Slicer can obtain a result similar to our Neural Slicer.
However, a ‘good’ initial guess is never easy to be obtained in general.

(a)(b)(c)(d)(e)

Figure 12. Study conducted on a Spiral-Fish model by using different initial
guesses for our Neural Slicer: (a) planar layers as a height field, (b) curved
layers from a field generated by the heat method (Crane et al., 2017) and (c)
curved layers as the result of $S^{3}$-Slicer (Zhang et al., 2022). The
optimization as a network learning process (see the learning curves in (d))
can always converge to the result with curved layers as shown in (e).

Differently, our Neural Slicer is robust to the initial guess of the mapping
$\lambda(\cdot)$. A study has been conducted on the Spiral Fish model as shown
in Fig.12 to demonstrate this advantage of our approach. The computing
processes starting from different initial guesses are tested, which include:

  1. (1)

Height Field – We first test the initial guess as planar layers, which is in
fact an identity mapping for $\lambda(\cdot)$ that can be easily realized by
setting $\mathbf{q(\mathbf{x})}$ and $\mathbf{s(\mathbf{x})}$ with neither
rotation nor scaling (see Fig.12(a)).

  2. (2)

Heat Transfer – The other initial field is obtained by using the heat method
(Crane et al., 2017) to generate an approximated geodesic distance field to
the bottom of the model (see Fig.12(b)). The network coefficients $\theta_{q}$
and $\theta_{s}$ for $\mathbf{q(\mathbf{x})}$ and $\mathbf{s(\mathbf{x})}$ are
pre-trained to make $G(\textbf{x})$ fit this field by minimizing the mean
squared error.

  3. (3)

$S^{3}$-Slicer – The scalar field obtained from the $S^{3}$-Slicer (see
Fig.12(c)) is employed to pre-train the initial network coefficients
$\theta_{q}$ and $\theta_{s}$ for $\mathbf{q(\mathbf{x})}$ and
$\mathbf{s(\mathbf{x})}$.

When applying our Neural Slicer to compute curved layers by using the above
different initial guesses, the learning curves of the NN optimizer are as
shown in Fig.12(d). All converge to the same result of curved layers as shown
in Fig.12(e). The area of overhang on our result can be further reduced by
$94.2\%$ w.r.t. the result from $S^{3}$-Slicer. This study demonstrates the
robustness of our framework with different initial guesses.

###  6.2. Physical Experiments

The curved layers generated by our framework are tested to fabricate models on
a multi-axis 3D printing hardware system as shown in Fig.1(d), which is
composed of an ABB IRB 2600 robotic arm with 6-DOFs and an ABB A250 positioner
with 2-DOFs. The repeatability of this robotic system is $\pm 0.05$mm and the
TCP/IP protocol is used for the communication between the extrusion system,
the robotic arms, and the PC. With the kinematic redundancy provided by this
system, materials can be deposed onto the in-process model along arbitrary
directions in the model’s space while still controlling the physical
orientation and the speed / acceleration of the printer head. The extruder
system is controlled by a Duet3D board and the nozzle dimension is $1.0$mm.
Polylactic Acid (PLA) and Polyvinyl Alcohol (PVA) filaments with 1.75mm
diameter are used in our physical fabrication, where all models are printed by
PLA filaments. PVA is only used to print the support structures for the Bridge
model. With the help of speed control, layers with thickness in the range of
$[0.4\mathrm{mm},1.0\mathrm{mm}]$ can be reliably produced by this hardware
system.

Figure 13. The results of all models that are fabricated using the support-
free curved layers generated by our Neural Slicer. Table 2. Statistics of
physical fabrication.

|  Curved Layers† |  Planar Layers†
---|---|---
Model | #Layer | Weight | Time | #Layer |  Weight‡ | Time
Bunny Hd. | $250$ |  $512.9$g | 21.3h | $243$ |  $659.4$g |  $24.6$h
Yoga | 150 |  $140.8$g | 18.7h | / | / | /
Shelf | 320 |  $521.0$g | 39.3h | / | / | /
Ring | $99$ |  $52.3$g | 2.4h | $87$ |  $100.1$g |  $4.5$h
Tubes | $150$ |  $324.7$g | 22.3h | $132$ |  $514.8$g |  $20.9$h
Spiral Fish | $500$ |  $326.6$g | 24.8h | $340$ |  $445.5$g |  $19.5$h
Bridge (Neu) | 350 |  $769.1$g | 26.2h | / | / | /
Bridge ($S^{3}$) | 600 |  $921.1$g | 44.8h | / | / | /

† The thicknesses for curved layers are in the range of [$0.4$mm, $1.0$mm],
and the thickness for planar layers is $0.8$mm.

‡ The weight of a model fabricated by planar layers includes the weight of its
supporting structures.

(a)(b)

Figure 14. Fabrication results of the Bridge model using (a) the layers
generated by our Neural Slicer and (b) the layers generated by the
$S^{3}$-Slicer, where the building orientations have been specified by the
arrows.  Figure 15. Results of the Bridge model under 3-Point bending test,
where the force is applied at the top-middle of the Bridge and the bottom of
the model is fixed. Force-displacement curves are generated to study the
mechanical strength of models fabricated from different curved layers. The
breaking force has been doubled ($\uparrow 101.9\%$) compared to the result of
$S^{3}$-Slicer (Zhang et al., 2022).  Figure 16. The results of compression
tests taken on two specimens of the Bunny Head model, where the force is
applied on the right ear of the model. The specimens are fabricated from
layers with only the SF requirement and the SF + SR requirements – both are
generated by our Neural Slicer. After imposing the SR requirement, the model’s
breaking force can be increased by $30.6\%$.

The models fabricated by our system using support-free curved layers are shown
in Fig.13. The statistics of physical fabrication are given in Table 2. For
the purpose of comparison, a few models are also fabricated by planar layers,
where additional supporting structures need to be printed therefore the models
are heavier. To have a larger search space, curved layers of the Bridge model
are computed by only imposing the SR requirement. Two Bridge models are
fabricated using layers generated by the $S^{3}$-Slicer and our Neural Slicer
(see Fig.14).

Mechanical tests have been conducted to test the strength of 3D printed models
on an INSTRON tensile machine 5960 with a capacity of 10kN. First of all, the
3-point bending tests are conducted on two Bridge models shown in Fig.14 after
removing the PVA supports. Force-displacement curves are measured during the
bending tests (see Fig.15). An improvement of $101.9\%$ in the breaking force
can be observed on the specimen fabricated from the layers generated by our
Neural Slicer, which is consistent with the simulation results as given in
Fig.9. We also tested the mechanical strength of the 3D printed Bunny Head
model on the same tensile machine. The bottom of the model is fixed by a
specially designed fixture and the compressing force is applied on the right
ear of the model (see Fig.16). Two specimens are tested in our experiment. One
is fabricated from the layers generated by imposing only the SF requirement
while the other is from the layered with both the SF and the SR requirements.
From the force-displacement curves, we can observe a $30.6\%$ increase of the
breaking force on the specimen with SR requirement. The broken region is
changed from the root of the ear to the hole of the ear (see the zoom-views in
Fig.16). The results of mechanical tests are also consistent with the
simulation results given in Fig.10, where the maximal strain on the model with
the SF + SR requirements is reduced by $36.5\%$ compared to the one with only
the SF requirement.

###  6.3. Discussion

One major limitation of our current approach is that the stress field employed
to define the strength reinforcement loss $\mathcal{L}_{SR}$ is obtained from
FEA using isotropic material properties (i.e., $\mathbf{\tau}_{\max}$ is given
and unchanged during the optimization). This simplification ignores the change
of $\mathbf{\tau}_{\max}$ caused by the anisotropic material properties
introduced by curved layers. We study this influence by studying the change of
angles between LPDs and the maximal stress directions when using anisotropic
FEA instead of the input distribution obtained from isotropic FEA. The
anisotropic FEA in this study employs the same parameters as those used in
Sec. 6.1.3. The statistics of angle changes on two models have been given in
Fig.17. It can be found that this change indeed has a certain level of
influence but is not very significant. Moreover, considering the change of
$\mathbf{\tau}_{\max}$ in curved slicing needs to incorporate the anisotropic
FEA in the loop of optimization, which will further increase the computational
complexity. How to make the procedure of principal stress analysis
differentiable is also a problem to be further explored. We plan to consider
this in our future research.

Figure 17. Statistics of the angles between LPDs and the directions of maximal
stresses obtained by isotropic FEA vs. anisotropic FEA.

The second limitation is that the computation of mapping is based on an
intermediate representation – the caging mesh $\mathcal{C}$. We conducted the
computational tests on the Spiral Fish model by using $\mathcal{C}$ in three
different resolutions. As shown in Fig.18, the computation of our Neural
Slicer successes on all cages although the converging speed and the level of
final loss are different. Caging meshes in different resolutions may have
different genus numbers.

(a)(b)(c)

Figure 18. The learning curves for the Spiral Fish model by using the caging
mesh in different resolutions – i.e., with (a) $43,147$ elements ($g=3$), (b)
$35,453$ ($g=2$) elements and (c) $16,153$ elements ($g=0$), where $g$
indicates the genus number of a caging mesh. To conduct a fair comparison,
only $\mathcal{L}_{SF}+\mathcal{L}_{PO}$ are visualized in the learning
curves.

The third limitation comes from how global collision is handled. When the
collision between a printer-head and the model / platform occurs, we employ
the same strategy as $S^{3}$-slicer (Zhang et al., 2022) to increase the
weight of the harmonic term to give more flat layers. In extreme cases, it
will give planar layers that can guarantee no collision while other
manufacturing objectives are compromised. The collision with robotic arms can
be prevented in the downstream steps of motion planning (Zhang et al., 2021;
Dai et al., 2020).

(a)(b)

Figure 19. Study of computational efficiency by using a scalar field, a vector
field, and the quaternion-scaling field to parameterize the mapping, where the
tests are conducted on (a) the Bridge model and (b) the Spiral Fish model.
Clearly, the computation on quaternion-scaling field converges faster.

In our Neural Slicer, the mapping $\lambda(\cdot)$ is ‘parameterized’ as a
quaternion field $\mathbf{q}(\mathbf{x})$ and a scaling field
$\mathbf{s}(\mathbf{x})$ (simply denoted as Q-S fields). In literature,
shifting-vector based fields and scalar fields are widely used to represent a
deformation field for non-rigid registration (ref. (Deng et al., 2021;
Sundararaman et al., 2022; Park et al., 2019). We argue that the deformation
parameterized on the quaternion and the scaling fields is more efficient for
computing a mapping for curved slicing. However, we are not able to prove this
by theory. Comparisons are conducted on example models to study the
computational efficiency by using different fields (see Fig.19).

To demonstrate the effectiveness of our NN-based method to optimize the
deformation as a mapping for slicing, we have conducted tests of directly
using stochastic gradient descent (SGD) and Adaptive Moment Estimation (Adam)
(Kingma and Ba, 2015) to optimize the quaternions and scales at the center of
tetrahedral elements and compared with our NN-based approach using different
SIREN layers (5 & 10 layers) on both the Bridge model and the Spiral Fish
model (see Fig.20 for the comparison). Significant differences in convergence
speed can be observed – especially when SIREN with 10 layers is employed. This
is because such a neural network structure can capture and realize the
required deformation more easily.

Explicitly imposing the requirement of injectivity in ARAP is not necessary in
our framework. This is due to the fact that we directly evaluate manufacturing
objectives as loss functions defined in the model space but not in the
deformed space. The resultant $G(\mathbf{x})$ is a function mapping from
$\mathbb{R}^{3}$ to $\mathbb{R}$ which allows the local rotation around the
$z$-axis in the deformed space. This gives more flexibility to achieve optimal
slicing results. The quality of resultant layers as isosurfaces are softly
controlled by the harmonic terms. In examples with challenging geometry (e.g.,
the Spiral Fish model), we find that the allowed self-intersection is helpful
to obtain results satisfying the manufacturing objectives better.

(a)(b)

Figure 20. Study of computational efficiency in directly using stochastic
gradient descent (SGD) and adaptive moment estimation (Adam) (Kingma and Ba,
2015) to optimize the quaternions and scales at the center of tetrahedral
elements. The comparison with our approach using different SIREN layers (5 &
10 layers) is also given on (a) the Bridge and (b) the Spiral Fish models.

##  7\. Conclusion

This paper introduces a representation-agnostic slicer for multi-axis 3D
printing, utilizing a computational pipeline based on neural networks. Our
slicer is designed to work on models with diverse representations and complex
topology. We leverage a deformation field to compute the mapping, determining
a scalar field in the space surrounding an input model. Isosurfaces of this
scalar field are then extracted to generate curved layers for 3D printing. In
contrast to previous approaches, where manufacturing objectives are indirectly
optimized, our approach enables direct optimization of the scalar field by
using loss functions directly based on local printing directions. This
optimization is achieved through the implementation of a differentiable
computation pipeline. Our approach demonstrates the successful handling of
models with intricate topology, facilitated by the utilization of a volumetric
mesh as a numerical computation cage. Notably, this allows for the input model
and the cage to possess (i.e., different genus numbers). Leveraging the
capabilities of a robust neural network solver, our approach allows for the
flexible adjustment of local printing directions within a model, eliminating
the need of good initial guesses for the optimization process. Our Neural
Slicer presented in this paper can generate results with significantly
enhanced performance, which has been verified in a variety of examples and in
physical experiments.

###### Acknowledgements.

The project is supported by the chair professorship fund at the University of
Manchester and UK Engineering and Physical Sciences Research Council (EPSRC)
Fellowship Grant (Ref.#: EP/X032213/1).

## References

  * (1)
  * Ahn et al. (2002) Sung-Hoon Ahn, Michael Montero, Dan Odell, Shad Roundy, and Paul K. Wright. 2002.  Anisotropic material properties of fused deposition modeling ABS.  _Rapid Prototyping Journal_ 8, 4 (2002), 248 – 257.
  * Aigerman et al. (2022) Noam Aigerman, Kunal Gupta, Vladimir G. Kim, Siddhartha Chaudhuri, Jun Saito, and Thibault Groueix. 2022.  Neural Jacobian Fields: Learning Intrinsic Mappings of Arbitrary Meshes.  _ACM Transactions on Graphics_ 41, 4, Article 109 (2022).
  * Arora et al. (2019) Rahul Arora, Alec Jacobson, Timothy R. Langlois, Yijiang Huang, Caitlin Mueller, Wojciech Matusik, Ariel Shamir, Karan Singh, and David I.W. Levin. 2019.  Volumetric michell trusses for parametric design & fabrication.  _Proceedings: SCF 2019 - ACM Symposium on Computational Fabrication_ (2019).
  * Bächer et al. (2014) Moritz Bächer, Emily Whiting, Bernd Bickel, and Olga Sorkine-Hornung. 2014.  Spin-it: optimizing moment of inertia for spinnable objects.  _ACM Transactions on Graphics_ 33, 4 (2014), 96:1–96:10.
  * Bartoň et al. (2021) Michael Bartoň, Michal Bizzarri, Florian Rist, Oleksii Sliusarenko, and Helmut Pottmann. 2021.  Geometry and tool motion planning for curvature adapted CNC machining.  _ACM Transactions on Graphics_ 40, 4 (2021).
  * Bhatt et al. (2022) Prahar M. Bhatt, Ashish Kulkarni, Rishi K. Malhan, Brual C. Shah, Yeo Jung Yoon, and Satyandra K. Gupta. 2022.  Automated Planning for Robotic Multi-Resolution Additive Manufacturing.  _Journal of Computing and Information Science in Engineering_ 22, 2 (2022).
  * Bouaziz et al. (2012) Sofien Bouaziz, Mario Deuss, Yuliy Schwartzburg, Thibaut Weise, and Mark Pauly. 2012.  Shape-up: Shaping discrete geometry with projections.  _Computer Graphics Forum_ 31, 5 (2012), 1657 – 1667.
  * Chakraborty et al. (2008) Debapriya Chakraborty, B. Aneesh Reddy, and A. Roy Choudhury. 2008.  Extruder path generation for Curved Layer Fused Deposition Modeling.  _CAD Computer Aided Design_ 40, 2 (2008), 235 – 243.
  * Crane et al. (2017) Keenan Crane, Clarisse Weischedel, and Max Wardetzky. 2017.  The heat method for distance computation.  _Commun. ACM_ 60, 11 (2017), 90 – 99.
  * Dai et al. (2020) Chengkai Dai, Sylvain Lefebvre, Kai-Ming Yu, Jo M. P. Geraedts, and Charlie C. L. Wang. 2020.  Planning Jerk-Optimized Trajectory with Discrete Time Constraints for Redundant Robots.  _IEEE Transactions on Automation Science and Engineering_ 17, 4 (2020), 1711 – 1724.
  * Dai et al. (2018) Chengkai Dai, Charlie C.L. Wang, Chenming Wu, Sylvain Lefebvre, Guoxin Fang, and Yong-Jin Liu. 2018.  Support-free volume printing by multi-axis motion.  _ACM Transactions on Graphics_ 37, 4 (2018).
  * De Goes et al. (2016) Fernando De Goes, Mathieu Desbrun, and Yiying Tong. 2016.  Vector field processing on triangle meshes.  _ACM SIGGRAPH 2016 Courses, SIGGRAPH 2016_ (2016).
  * De Luigi et al. (2023) Luca De Luigi, Adriano Cardace, Riccardo Spezialetti, Pierluigi Zama Ramirez, Samuele Salti, and Luigi Di Stefano. 2023.  Deep learning on implicit neural representations of shapes. In _The Eleventh International Conference on Learning Representations_.  https://openreview.net/forum?id=OoOIW-3uadi
  * Deng et al. (2021) Yu Deng, Jiaolong Yang, and Xin Tong. 2021.  Deformed Implicit Field: Modeling 3D Shapes with Learned Dense Correspondence.  _Proceedings of the IEEE Computer Society Conference on Computer Vision and Pattern Recognition_ (2021), 10281 – 10291.
  * Donti et al. (2021) Priya L. Donti, David Rolnick, and J. Zico Kolter. 2021.  DC3: A Learning Method for Optimization with Hard Constraints.  _ICLR 2021 - 9th International Conference on Learning Representations_.
  * Du et al. (2020) Xingyi Du, Noam Aigerman, Qingnan Zhou, Shahar Z. Kovalsky, Yajie Yan, Danny M. Kaufman, and Tao Ju. 2020.  Lifting simplices to find injectivity.  _ACM Transactions on Graphics_ 39, 4 (2020).
  * Duenser et al. (2020) Simon Duenser, Roi Poranne, Bernhard Thomaszewski, and Stelian Coros. 2020.  RoboCut: Hot-wire Cutting with Robot-controlled Flexible Rods.  _ACM Transactions on Graphics_ 39, 4 (2020).
  * Etienne et al. (2019) Jimmy Etienne, Nicolas Ray, Daniele Panozzo, Samuel Hornus, Charlie C.L. Wang, Jonàs Martínez, Sara McMains, Marc Alexa, Brian Wyvill, and Sylvain Lefebvre. 2019.  Curvislicer: Slightly curved slicing for 3-axis printers.  _ACM Transactions on Graphics_ 38, 4 (2019).
  * Ezair et al. (2018) Ben Ezair, Saul Fuhrmann, and Gershon Elber. 2018.  Volumetric covering print-paths for additive manufacturing of 3D models.  _CAD Computer Aided Design_ 100 (2018), 1 – 13.
  * Fang et al. (2020) Guoxin Fang, Tianyu Zhang, Sikai Zhong, Xiangjia Chen, Zichun Zhong, and Charlie C. L. Wang. 2020.  Reinforced FDM: Multi-axis filament alignment with controlled anisotropic strength.  _ACM Transactions on Graphics_ 39, 6 (2020).
  * Fu et al. (2015) Xiao-Ming Fu, Yang Liu, and Baining Guo. 2015.  Computing locally injective mappings by advanced MIPS.  _ACM Transactions on Graphics_ 34, 4 (2015).
  * Garanzha et al. (2021) Vladimir Garanzha, Igor Kaporin, Liudmila Kudryavtseva, François Protais, Nicolas Ray, and Dmitry Sokolov. 2021.  Foldover-free maps in 50 lines of code.  _ACM Transactions on Graphics_ 40, 4 (2021).
  * Gottschalk et al. (1996) Stefan Gottschalk, Ming C Lin, and Dinesh Manocha. 1996.  OBBTree: A hierarchical structure for rapid interference detection.  _Proceedings of the ACM SIGGRAPH Conference on Computer Graphics_ , 171 – 180.
  * Groueix et al. (2018) Thibault Groueix, Matthew Fisher, Vladimir G. Kim, Bryan C. Russell, and Mathieu Aubry. 2018.  3D-CODED: 3D Correspondences by Deep Deformation.  _Lecture Notes in Computer Science (including subseries Lecture Notes in Artificial Intelligence and Lecture Notes in Bioinformatics)_ 11206 LNCS (2018), 235 – 251.
  * Groueix et al. (2019) Thibault Groueix, Matthew Fisher, Vladimir G. Kim, Bryan C. Russell, and Mathieu Aubry. 2019.  Unsupervised cycle-consistent deformation for shape matching.  _Computer Graphics Forum_ 38, 5 (2019), 123 – 133.
  * Hu et al. (2015) Kailun Hu, Shuo Jin, and Charlie C.L. Wang. 2015.  Support slimming for single material based additive manufacturing.  _CAD Computer Aided Design_ 65 (2015), 1 – 10.
  * Huang et al. (2020) Jingwei Huang, Chiyu Max Jiang, Baiqiang Leng, Bin Wang, and Leonidas Guibas. 2020.  Meshode: A robust and scalable framework for mesh deformation.  _arXiv preprint arXiv:2005.11617_ (2020).
  * Huang et al. (2012) Jin Huang, Yiying Tong, Hongyu Wei, and Hujun Bao. 2012.  Boundary aligned smooth 3D cross-frame field.  _Jisuanji Fuzhu Sheji Yu Tuxingxue Xuebao/Journal of Computer-Aided Design and Computer Graphics_ 24, 2 (2012), 137 – 139.
  * Huang et al. (2016) Yijiang Huang, Juyong Zhang, Xin Hu, Guoxian Song, Zhongyuan Liu, Lei Yu, and Ligang Liu. 2016.  FrameFab: Robotic Fabrication of Frame Shapes.  _ACM Transactions on Graphics_ 35, 6 (2016).
  * Huynh (2009) Du Q. Huynh. 2009.  Metrics for 3D rotations: Comparison and analysis.  _Journal of Mathematical Imaging and Vision_ 35, 2 (2009), 155 – 164.
  * Jiang et al. (2020) Chiyu Max Jiang, Jingwei Huang, Andrea Tagliasacchi, and Leonidas Guibas. 2020.  ShapeFlow: Learnable Deformations among 3D Shapes. In _Proceedings of the 34th International Conference on Neural Information Processing Systems_.
  * Kingma and Ba (2015) Diederik P. Kingma and Jimmy Lei Ba. 2015.  Adam: A method for stochastic optimization.  _3rd International Conference on Learning Representations, ICLR 2015 - Conference Track Proceedings_ (2015).
  * Kovalsky et al. (2015) Shahar Z. Kovalsky, Noam Aigerman, Ronen Basri, and Yaron Lipman. 2015.  Large-scale bounded distortion mappings.  _ACM Transactions on Graphics_ 34, 6 (2015).
  * Li et al. (2022) Yamin Li, Dong He, Shangqin Yuan, Kai Tang, and Jihong Zhu. 2022.  Vector field-based curved layer slicing and path planning for multi-axis printing.  _Robotics and Computer-Integrated Manufacturing_ 77 (2022).
  * Li et al. (2012) Yufei Li, Yang Liu, Weiwei Xu, Wenping Wang, and Baining Guo. 2012.  All-hex meshing using singularity-restricted field.  _ACM Transactions on Graphics_ 31, 6 (2012).
  * Liao et al. (2021) Wentao Liao, Renjie Chen, Yuchen Hua, Ligang Liu, and Ofir Weber. 2021.  Real-time locally injective volumetric deformation.  _ACM Transactions on Graphics_ 40, 4 (2021).
  * Liu et al. (2008) Ligang Liu, Lei Zhang, Yin Xu, Craig Gotsman, and Steven J. Gortler. 2008.  A local/global approach to mesh parameterization.  _Computer Graphics Forum_ 27, 5 (2008), 1495 – 1504.
  * Lorensen and Cline (1987) William E. Lorensen and Harvey E. Cline. 1987.  Marching cubes: A high resolution 3D surface construction algorithm.  _Proceedings of the 14th Annual Conference on Computer Graphics and Interactive Techniques, SIGGRAPH 1987_ (1987), 163 – 169.
  * Maestre et al. (2023) Juan Montes Maestre, Yinwei Du, Ronan Hinchet, Stelian Coros, and Bernhard Thomaszewski. 2023.  Differentiable Stripe Patterns for Inverse Design of Structured Surfaces.  _ACM Transactions on Graphics_ 42, 4 (2023).
  * Mitra et al. (2023) Rahul Mitra, Liane Makatura, Emily Whiting, and Edward Chien. 2023.  Helix-Free Stripes for Knit Graph Design.  _Proceedings - SIGGRAPH 2023 Conference Papers_ (2023).
  * Mitropoulou et al. (2020) Ioanna Mitropoulou, Mathias Bernhard, and Benjamin Dillenburger. 2020.  Print Paths Key-framing: Design for non-planar layered robotic FDM printing.  _Proceedings - SCF 2020: ACM Symposium on Computational Fabrication_.
  * Morreale et al. (2021) Luca Morreale, Noam Aigerman, Vladimir Kim, and Niloy J. Mitra. 2021.  Neural Surface Maps.  _Proceedings of the IEEE Computer Society Conference on Computer Vision and Pattern Recognition_ (2021), 4637 – 4646.
  * Park et al. (2019) Jeong Joon Park, Peter Florence, Julian Straub, Richard Newcombe, and Steven Lovegrove. 2019.  Deepsdf: Learning continuous signed distance functions for shape representation.  _Proceedings of the IEEE Computer Society Conference on Computer Vision and Pattern Recognition_ 2019-June (2019), 165 – 174.
  * Paszke et al. (2019) Adam Paszke, Sam Gross, Francisco Massa, Adam Lerer, James Bradbury, Gregory Chanan, Trevor Killeen, Zeming Lin, Natalia Gimelshein, Luca Antiga, Alban Desmaison, Andreas Köpf, Edward Yang, Zach DeVito, Martin Raison, Alykhan Tejani, Sasank Chilamkurthy, Benoit Steiner, Lu Fang, Junjie Bai, and Soumith Chintala. 2019.  PyTorch: An imperative style, high-performance deep learning library.  _Advances in Neural Information Processing Systems_ 32 (2019).
  * Patané et al. (2013) Giuseppe Patané, David Xianfeng Gu, and Xin Shane Li. 2013.  SIGGRAPH Asia 2013 - Course ”surface- and volume-based techniques for shape modeling and analysis”.  _SIGGRAPH Asia 2013 Courses, SA 2013_ (2013).
  * Rabinovich et al. (2017) Michael Rabinovich, Roi Poranne, Daniele Panozzo, and Olga Sorkine-Hornung. 2017.  Scalable locally injective mappings.  _ACM Transactions on Graphics_ 36, 2 (2017).
  * Ray et al. (2016) Nicolas Ray, Dmitry Sokolov, and Bruno Lévy. 2016.  Practical 3D frame field generation.  _ACM Transactions on Graphics_ 35, 6 (2016).
  * Riddick et al. (2016) Jaret C. Riddick, Mulugeta A. Haile, Ray Von Wahlde, Daniel P. Cole, Oluwakayode Bamiduro, and Terrence E. Johnson. 2016.  Fractographic analysis of tensile failure of acrylonitrile-butadiene-styrene fabricated by fused deposition modeling.  _Additive Manufacturing_ 11 (2016), 49 – 59.
  * Romain et al. (2013) Prévost Romain, Whiting Emily, Lefebvre Sylvain, and Sorkine-Hornung Olga. 2013.  Make It Stand: Balancing Shapes for 3D Fabrication.  _ACM Transactions on Graphics_ 32, 4 (2013), 81:1–81:10.
  * Sacht et al. (2015) Leonardo Sacht, Etienne Vouga, and Alec Jacobson. 2015.  Nested cages.  _ACM Transactions on Graphics_ 34, 6 (2015).
  * Schüller et al. (2013) Christian Schüller, Ladislav Kavan, Daniele Panozzo, and Olga Sorkine-Hornung. 2013.  Locally injective mappings.  _Computer Graphics Forum_ 32, 5 (2013), 125 – 135.
  * Si (2015) Hang Si. 2015.  TetGen, a delaunay-based quality tetrahedral mesh generator.  _ACM Trans. Math. Software_ 41, 2 (2015).
  * Sitzmann et al. (2020) Vincent Sitzmann, Julien N.P. Martel, Alexander W. Bergman, David B. Lindell, and Gordon Wetzstein. 2020.  Implicit neural representations with periodic activation functions.  _Advances in Neural Information Processing Systems_ 2020-December (2020).
  * Sorkine and Alexa (2007) Olga Sorkine and Marc Alexa. 2007.  As-rigid-as-possible surface modeling. In _Proceedings of the Fifth Eurographics Symposium on Geometry Processing_ , Vol. 4. Citeseer, 109–116.
  * Sorkine and Botsch (2009) Olga Sorkine and Mario Botsch. 2009.  Interactive Shape Modeling and Deformation.. In _Eurographics (Tutorials)_. 11–37.
  * Sorkine et al. (2004) O. Sorkine, D. Cohen-Or, Y. Lipman, M. Alexa, C. Rössl, and H.-P. Seidel. 2004.  Laplacian surface editing.  _ACM International Conference Proceeding Series_ 71 (2004), 175 – 184.
  * Su et al. (2019) Jian-Ping Su, Xiao-Ming Fu, and Ligang Liu. 2019.  Practical Foldover-Free Volumetric Mapping Construction.  _Computer Graphics Forum_ 38, 7 (2019), 287 – 297.
  * Sullivan and Kaszynski (2019) C. Bane Sullivan and Alexander Kaszynski. 2019.  PyVista: 3D plotting and mesh analysis through a streamlined interface for the Visualization Toolkit (VTK).  _Journal of Open Source Software_ 4, 37 (2019), 1450.
  * Sundararaman et al. (2022) Ramana Sundararaman, Riccardo Marin, Emanuele Rodolà, and Maks Ovsjanikov. 2022.  Reduced Representation of Deformation Fields for Effective Non-rigid Shape Matching.  _Advances in Neural Information Processing Systems_ 35 (2022).
  * Tam and Mueller (2017a) Kam-Ming Mark Tam and Caitlin T. Mueller. 2017a.  Additive Manufacturing Along Principal Stress Lines.  _3D Printing and Additive Manufacturing_ 4, 2 (2017), 63 – 81.
  * Tam and Mueller (2017b) Kam-Ming Mark Tam and Caitlin T. Mueller. 2017b.  Additive Manufacturing Along Principal Stress Lines.  _3D Printing and Additive Manufacturing_ 4, 2 (2017), 63 – 81.
  * Treece et al. (1999) G.M. Treece, R.W. Prager, and A.H. Gee. 1999.  Regularized marching tetrahedra: Improved iso-surface extraction.  _Computers and Graphics (Pergamon)_ 23, 4 (1999), 583 – 598.
  * Vanek et al. (2014) J. Vanek, J.A.G. Galicia, and B. Benes. 2014.  Clever support: Efficient support structure generation for digital fabrication.  _Computer Graphics Forum_ 33, 5 (2014), 117 – 125.
  * Williams et al. (2019) Francis Williams, Teseo Schneider, Claudio Silva, Denis Zorin, Joan Bruna, and Daniele Panozzo. 2019.  Deep geometric prior for surface reconstruction.  _Proceedings of the IEEE Computer Society Conference on Computer Vision and Pattern Recognition_ 2019-June (2019), 10122 – 10131.
  * Wu et al. (2020) Chenming Wu, Chengkai Dai, Guoxin Fang, Yong-Jin Liu, and Charlie C. L. Wang. 2020.  General Support-Effective Decomposition for Multi-Directional 3-D Printing.  _IEEE Transactions on Automation Science and Engineering_ 17, 2 (2020), 599 – 610.
  * Wu et al. (2016) Rundong Wu, Huaishu Peng, François Guimbretière, and Steve Marschner. 2016.  Printing arbitrary meshes with a 5DOF wireframe printer.  _ACM Transactions on Graphics_ 35, 4 (2016).
  * Yang et al. (2021) Guandao Yang, Serge Belongie, Bharath Hariharan, and Vladlen Koltun. 2021.  Geometry Processing with Neural Fields.  _Advances in Neural Information Processing Systems_ 27 (2021), 22483 – 22497.
  * Yifan et al. (2020) Wang Yifan, Noam Aigerman, Vladimir G. Kim, Siddhartha Chaudhuri, and Olga Sorkine-Hornung. 2020.  Neural cages for detail-preserving 3D deformations.  _Proceedings of the IEEE Computer Society Conference on Computer Vision and Pattern Recognition_ (2020), 72 – 80.
  * Zhang et al. (2021) Tianyu Zhang, Xiangjia Chen, Guoxin Fang, Yingjun Tian, and Charlie C. L. Wang. 2021.  Singularity-Aware Motion Planning for Multi-Axis Additive Manufacturing.  _IEEE Robotics and Automation Letters_ 6, 4 (2021), 6172 – 6179.
  * Zhang et al. (2022) Tianyu Zhang, Guoxin Fang, Yuming Huang, Neelotpal Dutta, Sylvain Lefebvre, Zekai Murat Kilic, and Charlie C. L. Wang. 2022.  S3-Slicer: A General Slicing Framework for Multi-Axis 3D Printing.  _ACM Transactions on Graphics_ 41, 6 (2022).
  * Zhang et al. (2015) Xiaoting Zhang, Xinyi Le, Athina Panotopoulou, Emily Whiting, and Charlie C.L. Wang. 2015.  Perceptual models of preference in 3D printing direction.  _ACM Transactions on Graphics_ 34, 6 (2015).
  * Zhong et al. (2023) Fanchao Zhong, Haisen Zhao, Haochen Li, Xin Yan, Jikai Liu, Baoquan Chen, and Lin Lu. 2023.  VASCO: Volume and Surface Co-Decomposition for Hybrid Manufacturing.  _ACM Transactions on Graphics_ 42, 6 (2023).