Cross-sectional Structural Analysis for 3D Printing Optimization
Nobuyuki Umetani∗,
Ryan Schmidt†
Autodesk Research

Abstract
We propose a novel cross-sectional structural analysis techniquethat efficiently detects critical stress inside a 3D object. We slice
the object into cross-sections and compute stress based on bending
momentum equilibrium. Unlike traditional approaches based on
finite element methods, our method does not require a volumetricmesh or solution of linear systems, enabling interactive analysis
speed. Based on the stress analysis, the orientation of an object is
optimized to increase mechanical strength when manufactured with
3D printing.

|T|=10kt = 0.13sec
max
weaknessminweakness(a)(b)
pulling right fore leg

1Introduction

Democratized digital manufacturing devices such as desktop 3D
printers enable non-professionals to casually create physical objects. Various software packages such as MeshMixer [Schmidt
and Singh 2010] and Autodesk 123D [2010] provide end-users
with intuitive 3D modeling interfaces for 3D printing. However,novice users often create structurally-unsound shapes which caneasily break. Analyzing structural weakness is a significant challenge. Even simple rules-of-thumb may involve considering a large
number of possible force location/magnitude combinations, and itis hard to estimate the effect of such forces on complex 3D shapes.
Umetani et al. [2012] showed that interactive analysis significantly
helps novice users for plank-based furniture design. Our goal is
to apply a similar concept for free-form 3D objects. Traditionally,
structural validity of 3D objects was evaluated using the Finite Element Method (FEM). However, FEM involves time-consuming 3D
mesh generation and the solution of large linear systems. As a result, FEM is generally not integrated into the real-time visualizations of interactive tools, making it difficult for users to consider
structural robustness during incremental trial-and-error design.
To achieve fast shape validation, we propose cross-sectional structural analysis. We compute a number of cross sections of an object and compute forces on the cross sections. We extend the wellknown Euler-Bernoulli assumption [Timoshenko 1970] to greatly
reduce the complexity of the problem. The advantage of this assumption is that it focuses purely on the geometric relationship between a cross section and external load, ignoring deformation outside the cross section. Based on this simplification, we can find theweak parts of a 3D object at interactive rates.

pulling right fore leg

broke at
0.4kgf

CR Categories: I.3.5 [Computer Graphics]: Computational Geometry and Object Modeling—Physically based modeling;
Keywords: structural analysis, optimization, 3D printing

optimalprinting
direction(c)broke at3.7kgf(d)Figure 1: (a) our system visualizes structural weakness at interactive rates during mesh editing. (b) A computed optimal printing
direction. (c) A naı̈ve printed model easily fractures under external
force, while (d) the optimized print bears much larger forces.
Our manufacturing target is 3D printing with the Fused Deposition
Method (FDM). In this process, objects are constructed in layers,
with each layer composed of a thin filament of melted plastic [Gibson et al. 2009]. This layer-by-layer construction introduces significant structural anisotropy – in particular, vertical bonding between
layers is much weaker than the in-layer bonds. In this paper, we optimize the orientation of models for 3D printing to maximize theirmechanical strength, based on our novel cross-sectional structural
analysis.

2

Related WorkFEM has been applied for structural weakness detection to find theworst loading scenario [Zhou et al. 2013], to strengthen [Stava et al.
2012] and to partition [Luo et al. 2012] an input 3D object. We analyze structural weakness to optimize printing direction using novel
cross-sectional structural analysis.
Make It Stand [Prévost et al. 2013] deforms objects so that theystand upright without toppling by optimizing momentum equilibrium under gravity. We similarly consider momentum equilibrium,
but on cross sections under many different external forces.
Structural soundness has also been studied in graphics for masonry
design (e.g. [Panozzo et al. 2013]). Whiting et al. [2012] propose
a structural soundness evaluation function which considers bending
force between bricks. In the same spirit, we consider bending force
on cross sections to analyze structural soundness.

∗ e-mail:n.umetani@gmail.com† e-mail:ryan.schmidt@autodesk.com

3

MethodObjects manufactured via layer-by-layer 3D
printing techniques such as FDM exhibit significant anisotropy;
they can bear considerable tensile force in the direction of filament
strands, while being rather fragile when forces perpendicular to thefilament strands are applied. Because FDM accumulates layers ina certain printing direction Θ ∈ R3 , the printed object will fracFracture criteria.

real-time session
stress analysis
on cross
sections
slicingclustering

layers of cross sections

input shapesampling
multiple
directions
detected weakcross sections

clustered cross sections

optimized printing direction(shown with red arrow)

Figure 2: Overview of our algorithm.
ture much more easily in this direction. Within a layer, filament isusually oriented in many directions, so fracturing perpendicular toΘ is significantly less likely. We denote σ̂L as the maximum tensile stress under which a layer of filament separates. We conducteda simple experiment to break a cantilever beam and measured thevalue σ̂L as 60MPa.
An overview of our algorithm is shown in Fig. 2. In detecting possible breaking scenarios, we first assume that an object is uniformly
filled with an isotropic material which can bear maximum tensilestress σ̂L . Next, we find a set of cross sections for which external
forces can produce a stress that exceeds σ̂L inside that cross section.
Then, we determine a printing direction Θ such that the normal of
the weakest cross sections are as perpendicular to Θ as possible.

3.1

Cross-sectional Structural Analysis

Before describing our method, we
briefly overview the Euler-Bernoulli (EB) assumption, on which
our cross-sectional analysis is based. Assume a beam which has
axis direction z ∈ R3 , a cross section Ω of the beam that is perpendicular to z, and a line inside the beam – the neutral axis –which goes through the centroid of Ω (the neutral axis has directionz). For a beams bending deformation, the assumption is threefold(see Fig. 3). First, there is no in-plane deformation and the crosssection Ω remains planar. Second, the deformed Ω remains (approximately) perpendicular to the deformed neutral axis. Finally,
the neutral axis does not stretch. Based on these assumptions, theinfinitesimal strain tensor  ∈ R3×3 arising from bending is written
as  = z ⊗ z, where  ∈ R is a linear function over Ω and is zeroat its centroid. This assumption is widely-used in the field of engineering. In the computer graphics literature, PriMo [Botsch et al.
2006] used a similar assumption of the rigidity of cross sections.
Euler-Bernoulli assumption.

ΩΩ
neutral axisz(a) rest shapeto occur, and hence we can still use this assumption for general 3D
shapes for the purpose of critical stress detection.
Virtual cross-section. To implement
the EB assumption on 3D objects, we firstx
Π’need to construct a cross section and neug’tral axis (Ω, z) such that z is perpendicz
ular to Ω and passes through the centroidg ynof Ω. Although the medial axis [Miklos
ΠΩ
et al. 2010] can be used to find candidates,∆nthe medial axis is costly to generate. Instead we introduce a fast virtual cross section technique to approximately generate (Ω, z). First, take an arbitrary cross section of the 3D object Γ which has a normal directionn ∈ R3 . We consider another section Γ0 at a small distance ∆nthat has the same normal direction n. For sections Γ and Γ0 , wecompute the centroids g and g0 . Next, we compute a neutral axis
direction Γ as z = (g0 − g)/|g0 − g| with infinitesimally-small
distance ∆n → 0. Finally, we get a virtual cross-section Ω byprojecting Γ perpendicular to z so that the centroid of Ω remains g.

We consider the local coordinate x and y on the virtual cross section Ω taking the origin at centroid g. We define axis y such that itis on the original section Γ. If the angle between n and x – denotedθ – is large, the virtual section is no longer a good approximation ofan actual section in the direction of z. Hence, we discard sections
where θ > 45◦ . We consider a fixed set of sample directions for n,
specifically the 13 unique directions of the form (A, B, C)T , where
each of A, B, C takes the values (−1, 0, +1).
Given a virtual cross section and its
neutral axis direction (Γ, n) → (Ω, z), we now consider bending
moment equilibrium to compute the stress on Ω. We assume theexternal force f is applied at a point p on the object’s surface. The
resulting external bending moment τ ex ∈ R3 which rotates thevirtual cross section can be written as
External bending moment.

τ ex = (I − z ⊗ z) {(p − g) × f } ,(b) deformed shape

Figure 3: Deformation of a beam under the Euler-Bernoulli assumption; The cross section Ω is rigid and remains perpendicular
to the unstretched neutral axis.
While the EB assumption is a widely-used engineering model, its
application is traditionally limited to beam structures. Our challenge is to extend its application to free-form 3D objects. In thiscase, the assumption is reasonable if the representative length of across section Ω (e.g. circumradius) is much longer than the distance from Ω to the location the force is applied (i.e. moment arm).
Hence, the assumption becomes more accurate as the cross-section
thins and the moment-arm lengthens. Our key observation is thatthese are precisely the conditions in which fracture is most likely(1)where I is an identity matrix. Note that τ ex is in the plane of Ω. To
achieve momentum equilibrium, the virtual cross section must generate a bending moment which balances out this external moment.
Using the EByηassumption, we can denote z-directionalstrain as a linear function of xy-coordinateηmax(x, y) = ¯X x + ¯Y y. From Hooke’s law, zξ
directional stress is given as σ = E whereψ
E is the Young’s modulus. We define orx
0
thogonal axes ξ and η such that strain isΩzero (¯X x + ¯Y y = 0) on the axis ξ. Letψ be the angle between ξ and x. With thisaxis transformation, the strain can be denoted in ξη-coordinate as(ξ, η) = ¯ξ η, where ¯ξ = ¯X cos ψ + ¯Y sin ψ. The InternalInternal bending moment.

bending moment vector then becomes:Zτ in =
η(E)dΩξ = ¯ξ EIξ ξ,(2)clusteringΩ

where Iξ is the second moment of area: Iξ =RΩ

η 2 dΩ.

Considering the moment equilibriumτ ex + τ in = 0, we obtain ξ = τ ex /|τ ex | and ¯η = |τ ex |/(EIξ ).
The maximum stress on the virtual cross section can be written
as σmax = |τ ex |/Zξ , where Zξ = Iξ /|η|max and |η|max is the
maximum of the absolute value of coordinate η inside this virtualcross section. In the engineering literature, Zξ is called the section modulus, relating the bending moment and maximum stress.
When σmax exceeds σ̂L , this section may fracture. Note that the
maximum stress does not depend on the Young’s modulus or othermaterial parameters; based on the cross-section’s shape, force position, and its direction, we can estimate stress inside the object in apurely geometric manner.
Moment equilibrium.

Figure 5: Left: two cross sections are adjacent if triangles areconnecting them without crossing other sections. Right: sectionsin the same plane are grouped based on adjacency information. The
neutral axis (black line) is computed for each section group.

3.2

Optimized 3D Printing Direction

The result of the analysis technique described above is a large number of combinations of cross section and external force. Next wemust detect which section/force pairs are most likely to break.
Structural weakness of a cross section isdefined considering two factors; the magnitude of force needed tobreak a section, and the area over which the force is distributed.
Let fmin (Γ, p, ψ) be the minimal force at position p on the surface needed to produce the critical stress on section Ω with bending
direction ψ. From Eq.(1) and Eq.(2), we have
Weakness evaluation.

So far, we have described a way
to compute maximum stress on a cross section, given an external
force. However, in a complex shape, there may be multiple disjoint cross-sections in a single slice, and the distribution of forces
depends on the boundary conditions (i.e. where the object is fixed,
and where the external force is applied). Hence, for a group of planar cross sections, we divide the object into two connected surfaces
S + and S − , with S + having smaller area than S − . Then, we apply
the force somewhere on S + and fix the object at a point somewhere
on S − . The motivation for this heuristic is that we usually holdthe larger part of an object when applying a force. Fig. 4 illustrates this process. First, we compute all cross-sections of an object {γ1 , γ2 , . . .} in a given plane. Then, we greedily cluster these
sections into groups Γ := {Γ1 , Γ2 , . . .} such that each group’s sections Γi ∈ Γ divide the object into a S + and S − . Each group isformed by picking an initial section γi , and then for remaining sections γj , checking if a path from γi to γj exists on each side of thesection plane. If so, γj is clustered with γi . We iterate this processuntil all the sections have been separated into groups, and then foreach group Γi , we compute the virtual cross section Ωi .
Clustering cross-sections.

+

γ1=Γa

( γ2 γ3 ) =Γb

γ4=Γc

Sa

Sa
(a)

Sb

+

Sb

S+
c

Sc

(b)

(c)

Figure 4: Three alternatives for separation of an example object
into surfaces S + (red) and S − (blue). We apply force at somewhere on S + . The cross-sections are grouped as {Γa , Γb , Γc } ={γ1 , (γ2 , γ3 ), γ4 }. As shown in (b), γ2 and γ3 are in the same group
because it is connected from both side.

To efficiently determine if sections areconnected, we use a connectivity graph. We first compute sectionsof the object at interval h along n (we use h = 32). Next we construct a section adjacency graph, where sections are connected if
there is a path of triangles between them that does not cross othersections (Fig. 5-left). We can then quickly determine if sections areaccessible from both sides using this adjacency information (Fig. 5right). Finally, for each section Γi we denote the set of other sections that lie on Si+ as Γ+i , which will be used for approximate
integration in the next section.
Clustering Computation.

fmin (Γ, p, ψ) = |(I − ξ ⊗ ξ)(p − g)| σ̂L Zξ η

(3)

Next we define a function B(Γ, p, ψ) which indicates the relative
ease with which a section/force configuration can be broken. B
is defined as 1/fmin if |fmin | < fˆmax , and 0 otherwise, wherefˆmax is a maximum external force (we choose 0.2 kgf). We thendefine the relative weakness of a section Γi by integrating B overthe potential force point locations on S +Z
Z
W (Γi ) =B(Γi , p(s), ψ)dψds,
(4)
s∈Si+

2π

where p(s) is the position of force point at a surface point s. We
approximate this costly surface integral by computing only on thesections. Since the input shape is a triangle mesh, its cross sections
are polylines. Thus, we compute the integration by aggregating B
multiplied by lm , which is half the sum of the two neighboring edge
lengths for the vertex sm on a polyline. Furthermore, ψ is sampled
at intervals J (we used J = 16):JX
2πh X X
2πj
W (Γi ) 'lm
B Γi , p(sm ),,JJ+ s ∈Γ
j=1
Γk ∈Γim(5)k

We further accelerate computation by culling sections Γk ∈ Γ+iwhich cannot break the cross section Γi . This is accomplished bycomputing the bounding volume of Γk (we used the circumcircle
as a bounding volume) and checking that the minimum breaking
force |fmin | is larger than the maximum-force threshold fˆmax for
all points on the section.
Using the machinery describedthus far, we can proceed to optimize the object orientation for
printing. For a set of given normal directions n, we slice theobject at intervals h and compute our weakness metric W (Γi )for all sections, and then store the identified “weak” sections as
Γw := {Γj ∈ Γw |W (Γj ) 6= 0}. Once the weak cross-sections
are detected, we need only to select an optimal printing direction.
As we previously mentioned, objects printed with FDM are mostly
Print Orientation Optimization.

likely to break between layers. Since the bending moment excites
force in the direction of neutral axis z, the object should be oriented
such that z is as vertical as possible. Hence, we use our weakness metric to weight the neutral axis directions zj corresponding
to identified weak sections Γj , and then obtain the optimal printing
direction Θopt asX
Θ · W (Γj )zj .
(6)
Θopt = arg minΘ∈R3 ,|Θ|=1

Γj ∈Γw

Eq. (6) can be solved
P by eigen-analysis of a 3x3 symmetric positive
definite matrix Γj ∈Γw W (Γj )zj ⊗ zj . The smallest eigenvector gives the strongest printing direction Θopt . Note that in somecases no printing direction can strengthen all weak cross sections.
However we have observed that such cases are rare.

4

Results

In Fig. 6, we compare our method with a standard second-orderlinear FEM stress analysis. The force distribution on two different
virtual cross sections is displayed. As expected, when the initial
section is relatively orthogonal to the local surface, so is the virtual section, and the force magnitude/directions agree well with theground-truth FEM. As the virtual section diverges from the original
section, the directions necessarily diverge as well, but the maximum
magnitudes remain comparable. Unlike FEM, our method is not affected by the quality of the surface mesh. Degenerate and “skinny”
sliver triangles are common in CAD meshes, but they do not affectour computation (see Fig. 7-rightmost).

|T|=34kt = 0.31s

40mm|T|=105kt = 0.90s
25mm

max
weakness|T|=10kt = 0.13s

23mm

min
weakness

Figure 7: More examples. |T | is the triangle count, and t is thetotal time taken to compute our cross-sectional analysis (using anIntel XeonTM 2.5 GHz CPU machine). Red arrows show the optimal
printing directions.

shape editing. Based on this analysis we optimize the orientation of
the model for 3D printing to increase structural soundness.

References
AUTODESK, 2010. Autodesk 123d. http://www.123dapp.com/.
B OTSCH , M., PAULY, M., G ROSS , M., AND KOBBELT, L. 2006.
PriMo: coupled prisms for intuitive surface modeling. In Proc.
of the 4th Eurographics Symp. on Geom. Proc., SGP ’06, 11–20.

L UO , L., BARAN , I., RUSINKIEWICZ , S., AND M ATUSIK , W.
2012. Chopper: partitioning models into 3D-printable parts.
ACM Trans. Graph. 31, 6 (Nov.), 129:1–129:9.

fixed b.c.
(b) forces on the cross section

Figure 6: Comparison with FEM. Cyan vectors show FEM forces,
purples are ours, blue loops are original input cross section, andgreen loops are virtual cross sections.
We also implemented our technique in a 3D mesh modeling
tool [Schmidt and Singh 2010]. The accompanying video showsinteractive modeling sessions with our structural analysis runningon a background thread, and updating real-time visualizations whenthe model is edited. Fig. 7 shows weak regions, optimal printing directions, and computation times for various objects. Our algorithm
scales well up to detailed meshes. In the video, we see the visualization update for a mesh with 650k triangles in a few seconds.
We 3D printed several models shown in the figures, and applied
forces to evaluate their strength (see Fig. 1 and accompanying
video). We used the MakerBot Replicator 2 and MakerWare v2.2
for printing, with the objects printed as solid (100% infill) to maximize strength. The material we used was polylactic acid (PLA) filament. Our algorithm assumes the object breaks with bending force,hence it cannot handle complex fracture such as buckling. However, we observed that the locations where objects fracture agreed
well with the weakest cross-sections our algorithms identified, andthat our optimal orientation makes objects significantly stronger.

5

65mm

G IBSON , I., ROSEN , D. W., AND S TUCKER , B.
2009.
Additive Manufacturing Technologies: Rapid Prototyping to
Direct Digital Manufacturing, 1st ed. Springer Publishing, Inc.

external
force(a) boundary conditions|T|=17kt = 0.18s

Conclusion

We present a framework for structural weakness detection of 3D objects which is fast enough to provide real-time feedback during 3D

M IKLOS , B., G IESEN , J., AND PAULY, M. 2010. Discrete scale
axis representations for 3D geometry. ACM Trans. Graph. 29(July), 101:1–101:10.
PANOZZO , D., B LOCK , P., AND S ORKINE -H ORNUNG , O. 2013.
Designing unreinforced masonry models. ACM Trans. Graph.
32, 4 (July), 91:1–91:12.
P R ÉVOST, R., W HITING , E., L EFEBVRE , S., AND S ORKINE H ORNUNG , O. 2013. Make it stand: balancing shapes for 3d
fabrication. ACM Trans. Graph. 32, 4 (July), 81:1–81:10.
S CHMIDT, R., AND S INGH , K. 2010. MeshMixer: an interface for
rapid mesh composition. ACM SIGGRAPH 2010 Talks.
S TAVA , O., VANEK , J., B ENES , B., C ARR , N., AND M ĚCH , R.
2012. Stress relief: improving structural strength of 3D printableobjects. ACM Trans. Graph. 31, 4 (July), 48:1–48:11.
T IMOSHENKO , S. 1970. Theory of Elasticity, 3 ed. Mcgraw-Hill
College, 6.
U METANI , N., I GARASHI , T., AND M ITRA , N. J. 2012. Guided
exploration of physically valid shapes for furniture design. ACM
Trans. Graph. 31, 4 (July), 86:1–86:11.
W HITING , E., S HIN , H., WANG , R., O CHSENDORF, J., AND D U RAND , F. 2012. Structural optimization of 3D masonry buildings. ACM Trans. Graph. 31, 6 (Nov.), 159:1–159:11.
Z HOU , Q., PANETTA , J., AND Z ORIN , D. 2013. Worst-case structural analysis. ACM Trans. Graph. 32, 4 (July), 137:1–137:12.