Scalable Stress-Constrained Topology Optimization: A Super-resolution and Simultaneous Approach Based on Machine Learning

Xiangyu Sun; Shirin Hosseinmardi; Ramin Bostanabad

cs.LGOct 20, 2026

Scalable Stress-Constrained Topology Optimization: A Super-resolution and Simultaneous Approach Based on Machine Learning

Authors: Xiangyu Sun, Shirin Hosseinmardi, Ramin Bostanabad

Organizations: Department of Mechanical and Aerospace Engineering, University of California, Irvine, CA, USA · Department of Mechanical Engineering, Northwestern University, Evanston, IL, USA · Department of Civil and Environmental Engineering, University of California, Irvine, CA, USA

Abstract

Large-scale topology optimization (TO) remains computationally demanding, especially for stress-constrained problems. To address this challenge, we propose the simultaneous neural topology optimization (SNTO) framework where we represent state, adjoint, and design variables with custom neural networks (NNs) whose parameters are jointly optimized within a single optimization loop that leverages continuous adjoint analysis and automatic differentiation. This optimization minimizes a multi-component loss function that combines the design objective, an aggregation of the local stress constraints, and the potential energy functionals. Our NNs have a unique architecture that accelerates convergence while enabling the networks to capture sharp gradients and boundaries which are critical to stress-constrained TO. To achieve scalability to millions of degrees of freedom (DoF), we develop a hardware-aware training scheme in which the design domain is partitioned into subdomains that are distributed across multiple GPUs. To control peak GPU memory usage, we further slice each subdomain into smaller regions whose loss gradients are evaluated sequentially and accumulated to update the NN parameters. Unlike domain decomposition approaches, our partitioning naturally ensures field continuity across the subdomains and does not require the introduction of additional terms into the loss function. For validation, we first evaluate the proposed framework on various 2D stress-constrained benchmarks and compare our results to a conventional finite element-based approach, showing that SNTO produces super-resolution designs with sharper boundaries, smaller gray areas, and lower volume fractions. We then consider large-scale 3D problems with up to eight million collocation points and distribute the training over up to eight GPUs.

Figures & tables

Panel of Figure 1 — Figure 1: SNTO framework for stress-constrained TO: (a) Simultaneous optimization pipeline. Three independent PGCANs generate the raw outputs $\bm{m}_{\bm{u}}(\bm{x};\bm{\theta}_{u})$ , $\bm{m}_{\bm{u}}^{(a)}(\bm{x};\bm{\theta}_{u}^{(a)})$ , and $m_{\rho}(\bm{x};\bm{\theta}_{\rho})$ , which are analytically transformed via $A_{\bm{u}}$ , $A_{\bm{u}}^{(a)}$ , and $A_{\rho}$ into displacement, adjoint, and density fields that satisfy their respective constraints by construction. The potential energies $\mathcal{L}_{M}(\cdot)$ and $\mathcal{L}_{M}^{(a)}(\cdot)$ , the stress constraint penalty $\mathcal{P}^{(a)}(\cdot)$ , and the objective of volume fraction $C_{V}(\cdot)$ are combined into a single loss that is minimized simultaneously with respect to $\bm{\theta}_{u}$ , $\bm{\theta}_{u}^{(a)}$ , and $\bm{\theta}_{\rho}$ . The individual loss terms are listed on the right. (b) PGCAN architecture. The trainable feature tensor $\bm{F}_{0}$ is convolved to obtain $\bm{F}_{c}$ , which is interpolated at each query point $\bm{x}$ ; the interpolated features are then passed to a shallow gated decoder.

Panel of Figure 2 — Figure 2: Memory-efficient training scheme for large-scale stress-constrained TO: (a) Single-GPU setup, in which the centroid CPs $\bm{x}_{C}$ and the associated nodal CPs $\bm{x}_{N}$ of the design domain are partitioned into $S$ sub-slices that are passed to the GPU sequentially; (b) multi-GPU setup, in which the design domain is sliced into $R$ subdomains $\Omega_{1},\dots,\Omega_{R}$ , one per GPU, and each subdomain is further divided into $S$ sub-slices that are swept sequentially on its GPU, so that the autograd graph never holds more than one sub-slice. The parameter gradients are accumulated over the sub-slices and synchronized across the GPUs before every update of the shared NNs.

Panel of Figure 3 — Figure 3: Design domains, boundary conditions, geometries and applied external loads for 2D benchmarks: (a) the L beam example with a re-entrant corner, (b) the crack example with a constrained edge segment, (c) the hook example with distributed loading along the inner contact surface, and (d) the portal frame example with the half design domain.

Panel of Figure 4 — Figure 4: Large-scale 3D stress-constrained TO benchmarks used to evaluate the scalability of the proposed framework: (a) the L beam example with a re-entrant corner; (b) the portal frame example with a horizontally constrained side face, a supported lower surface, and a distributed top load surface.

Example	$H$ ( $\mathrm{mm}$ )	$L$ ( $\mathrm{mm}$ )	$W$ ( $\mathrm{mm}$ )	$R_{1}$ ( $\mathrm{mm}$ )	$R_{2}$ ( $\mathrm{mm}$ )	$L_{s}$ ( $\mathrm{mm}$ )	$f$ ( $\mathrm{N}$ )
L beam 2D	1	1	0.1	–	–	0.05	0.2
Crack 2D	2	1	0.1	–	–	0.1	0.5
Hook 2D	1.5	1	0.1	0.1	0.19	0.597	0.5
Portal frame 2D	0.6	0.6	0.1	–	–	0.05	0.15
L beam 3D	1	1	0.1	–	–	0.05	0.2
Portal frame 3D	0.6	0.6	0.1	–	–	0.05	0.15

Table 1: Geometry and loading of the benchmark examples:

H

and

L

are the height and width of the design domain, and

W

is the out-of-plane dimension, which is the plane stress thickness in 2D and the extrusion depth in 3D. For the hook,

R_{1}

and

R_{2}

are the radii of the fixed hole at the top and the loaded hole at the bottom, respectively; the upper half of the

R_{1}

rim is clamped and the lower half of the

R_{2}

rim carries the load.

L_{s}

is the length of the loaded boundary segment, which for the hook is the arc length

\pi R_{2}

of the loaded half-rim, and

f

is the total external force distributed uniformly over

L_{s}

in 2D and over

L_{s}\times W

in 3D.

Example (2D)	PolyStress		SNTO 2D 125K		SNTO 2D 250K		SNTO 2D 500K		SNTO 2D 1M		$\operatorname{Res}$
Example (2D)	Nodal CPs	Centroid CPs	Nodal CPs	Centroid CPs	Nodal CPs	Centroid CPs	Nodal CPs	Centroid CPs	Nodal CPs	Centroid CPs	$(N_{x}^{e},N_{y}^{e})$
L beam	251251	251000	125153	124270	250891	249641	502704	500934	1006008	1003503	(24, 24)
Crack	251695	250632	125900	125148	251525	250462	501501	500000	1002510	1000387	(24, 48)
Hook	495406	250000	125219	124123	251262	249704	501870	499666	1001114	997998	(24, 36)
Portal frame	495557	250000	125000	124261	250790	249743	500823	499343	1000785	998693	(24, 24)
Example (3D)	SNTO 3D 500K		SNTO 3D 1M		SNTO 3D 2M		SNTO 3D 4M		SNTO 3D 8M		$\operatorname{Res}$

Table 2: Mesh and feature-grid resolutions: Numbers of nodal and centroid CPs of the meshes used for the 2D examples (top) and the 3D examples (bottom); the nominal sizes are given in the column headers. The number of centroid CPs equals the number of local stress constraints (

N_{c}=n_{e}

). PolyStress uses a single polygonal mesh per 2D example, whereas SNTO uses structured quadrilateral meshes in 2D and structured hexahedral meshes in 3D. All SNTO 2D cases are trained on a single GPU; the 125K, 250K, and 500K meshes are processed as one slice, while the 1M meshes are divided into two sub-slices. The last column lists the resolution

\operatorname{Res}

of the PGCAN feature grid.

Figure 5: Comparison of optimized 2D topologies obtained by PolyStress and SNTO: The left column shows the problem geometry for each benchmark, followed by the PolyStress solution (the second column) and three SNTO solutions from the five runs with 250K nodal CPs corresponding to the maximum, median, and minimum final volume fraction. The density field is plotted from void (blue) to solid material (red), with the same color scale used across all cases. The comparison shows that SNTO recovers physically interpretable geometries containing corners and cracks while producing designs with reduced gray regions compared to PolyStress. With SNTO, the density field can be visualized at any resolution, which is set here to the centroids corresponding to the 250K nodal CPs used in training.

Example	Metric	PolyStress 2D	SNTO 2D 125K	SNTO 2D 250K	SNTO 2D 500K	SNTO 2D 1M
L beam 2D	Volume fraction	0.365	0.341	0.349	0.349	0.343
	Gray fraction	0.087	0.027	0.035	0.026	0.023
	Runtime (sec)	2435	1432	2722	5221	10705
Crack 2D	Volume fraction	0.258	0.244	0.247	0.241	0.248
	Gray fraction	0.054	0.021	0.017	0.012	0.015
	Runtime (sec)	2559	1276	2341	4423	8758

Table 3: Median volume fraction, gray area fraction, and wall-clock runtime of the optimized 2D topologies: For each example, the volume fraction of the final design, the gray area fraction (fraction of centroid CPs with

0.1<\rho<0.9

), and the wall-clock runtime in seconds are reported for PolyStress and for SNTO with 125K, 250K, 500K, and 1M nodal CPs. Each value is the median over 3 independent runs for PolyStress and 5 for SNTO.

Panel of Figure 6 — Figure 6: Evolution of the volume fraction and the selected topology during SNTO training: Volume fraction histories and density snapshots are shown for the (a) L beam, (b) crack, (c) hook, and (d) portal frame benchmarks with 250K nodal CPs. Red markers indicate the epochs at which the intermediate designs are sampled; the first marker lies in the warm-up stage, the second at the peak of the volume fraction, and the last at the end of training. The snapshots show the transition from a uniform initial design to a diffuse layout to a crisp binary structure as the volume fraction decreases and settles.

Panel of Figure 7 — Figure 7: Training histories of SNTO for the L beam example with 250K nodal CPs: (a) active element fraction and wall-clock time per epoch; (b) volume fraction and gray area fraction; (c) stress-violation measures $SV_{sd}$ and $SV_{al}$ ; (d) adjoint scaling factor $\alpha_{adj}$ of Equation 43 ; (e) twice the strain energy of the displacement field compared with the external work; (f) twice the strain energy of the adjoint field compared with the body work. At equilibrium the two quantities in (e) and (f) coincide (Clapeyron’s theorem), so their agreement is a direct check that both networks satisfy Equations 5 and 41 throughout the optimization.

Figure 8: Mesh-resolution study for optimized 2D topologies: Final density fields of the run with the median volume fraction are shown for the L beam, crack, hook, and portal frame benchmarks with 125K, 250K, 500K, and 1M nodal CPs. The consistent topologies across refinements indicate that SNTO is essentially mesh independent, preserving the main structural features as the number of local stress constraints increases by nearly an order of magnitude. The density scale is shared across all panels, with blue denoting void and red denoting solid material. The SNTO density fields are plotted at the centroids of the nodal CPs used in training, although they can be evaluated at any resolution.

Panel of Figure 9 — Figure 9: Sensitivity of evolution histories with respect to number of nodal CPs using the Crack benchmark: Evolution histories are compared for the runs with 125K, 250K, 500K, and 1M nodal CPs. The panels report (a) active element fraction, (b) volume fraction, (c) gray area fraction, (d) stress violation over the realized solid material, (e) adjoint scaling factor, and (f) body work contribution in the adjoint problem. The overlapping trends show that the training dynamics and constraint response remain consistent under mesh refinement.

Panel of Figure 10 — Figure 10: Stress-constraint verification of optimized 2D designs: Distributions of normalized von Mises stress, $\sigma_{\mathrm{VM}}/\sigma_{\mathrm{y0}}$ , with the corresponding normalized principal-stress states are shown for SNTO with 250K nodal CPs and PolyStress on the left and right, respectively. The stress is computed with FEniCS on the thresholded and remeshed designs. In the scatter plots, each point is the principal-stress state of one element, the black ellipse is the von Mises yield envelope, and the colors indicate the local stress ratio.

Panel of Figure 11 — Figure 11: Optimized topologies (SNTO with 2M nodal CPs): (a) the L beam and (b) the portal frame, each shown from an isometric view and a front view. The front views closely resemble the 2D designs in Fig. 5 , while the isometric views reveal 3D features such as the smoothly rounded load path at the re-entrant corner of the L beam and the fillet through the thickness where the diagonal strut of the portal frame meets the upper chord.

Example	Metric	SNTO 3D 500K	SNTO 3D 1M	SNTO 3D 2M	SNTO 3D 4M	SNTO 3D 8M
L beam 3D	GPUs $\times$ sub-slices	$1\times 2$	$2\times 2$	$2\times 4$	$4\times 4$	$4\times 10$
	Volume fraction	0.342	0.347	0.350	0.345	0.351
	Gray fraction	0.015	0.013	0.009	0.015	0.015
	Runtime (sec)	7661	8912	17784	22037	50532
Portal frame 3D	GPUs $\times$ sub-slices	$1\times 2$	$2\times 2$	$2\times 4$	$4\times 4$	$4\times 10$
	Volume fraction	0.234	0.224	0.228	0.248	0.258

Table 4: Volume fraction, gray area fraction, and wall-clock runtime of the optimized 3D benchmarks: For each example, the values are reported for the five mesh resolutions listed in Table 2 . Each value is the mean over two independent runs, except for the 8M cases, which were run once. All cases are trained for 30,000 epochs on NVIDIA A100 GPUs. The number of GPUs and the number of sub-slices per GPU are listed in the first row of each example.

Figure 12: Mesh-resolution study for optimized 3D topologies: Final designs are shown for the L beam (top) and portal frame (bottom) benchmarks with progressively more nodal CPs. The main structural features are preserved as the number of local stress constraints grows by more than an order of magnitude, which shows that SNTO is essentially mesh independent even without any filtering.

Panel of Figure 13 — Figure 13: Sensitivity of evolution histories with respect to number of nodal CPs using the 3D L beam benchmark: Same quantities and panel layout as Fig. 9 , for the runs with 500K, 1M, 2M, 4M, and 8M nodal CPs.

Panel of Figure 14 — Figure 14: Stress-constraint verification for the optimized 3D designs. Normalized von Mises stress distributions are shown for the L beam in the left column and the portal frame in the right column using 500K, 1M, 2M, 4M, and 8M nodal CPs. Each panel presents the normalized von Mises stress $\sigma_{\mathrm{VM}}/\sigma_{\mathrm{y0}}$ over the optimized geometry, a transparent rendering in which regions with $\sigma_{\mathrm{VM}}>\sigma_{\mathrm{y0}}$ are highlighted in dark green, and the corresponding normalized principal-stress states $(\sigma_{1},\sigma_{2},\sigma_{3})/\sigma_{\mathrm{y0}}$ . The black curve denotes the von Mises limit, and the scatter colors indicate the local stress ratio.

Panel of Figure 15 — Figure 15: Computational cost and memory scaling of the 3D examples: (a) Sub-slicing sweep on the L beam with 1M nodal CPs using two GPUs. The training time (left axis, dark red circles) and the peak GPU memory usage (right axis, dark blue diamonds) are plotted against the number of sub-slices per GPU (2, 4, 8, and 16). Finer sub-slicing lowers the peak memory at the cost of a longer training time. (b) Multi-GPU sweep on the portal frame with 2M nodal CPs. The same two quantities are plotted against the number of GPUs (1, 2, 4, and 8), with the number of sub-slices per GPU reduced accordingly (8, 4, 2, and 1) so that the total number of sub-slices and the peak memory usage stay roughly constant. The training time decreases nearly linearly with the number of GPUs. The optimized topologies are shown next to the corresponding data points.

Appendix figures & tables4 assets

Supplementary material from the paper’s appendix.

Appendix

Panel of Figure 16 — Figure 16: Sensitivity of evolution histories with respect to number of nodal CPs using the L beam benchmark: Same quantities and panel layout as Fig. 9 , for the runs with 125K, 250K, 500K, and 1M nodal CPs.

Panel of Figure 17 — Figure 17: Sensitivity of evolution histories with respect to number of nodal CPs using the hook benchmark: Same quantities and panel layout as Fig. 9 , for the runs with 125K, 250K, 500K, and 1M nodal CPs.

Panel of Figure 18 — Figure 18: Sensitivity of evolution histories with respect to number of nodal CPs using the portal frame benchmark: Same quantities and panel layout as Fig. 9 , for the runs with 125K, 250K, 500K, and 1M nodal CPs.

Panel of Figure 19 — Figure 19: Sensitivity of evolution histories with respect to number of nodal CPs using the 3D portal frame benchmark: Same quantities and panel layout as Fig. 9 , for the runs with 500K, 1M, 2M, 4M, and 8M nodal CPs.

Topology Optimization

Topology optimization finds material layouts that meet performance and manufacturing constraints across engineering systems.

See all 34 papers

Sep 23, 2026

KATOsuper: Surrogate-accelerated neural topology optimization with sensitivity-consistent Fourier neural operators

Topology optimization (TO) remains computationally intensive due to repeated finite element analysis (FEA) evaluations required at each iteration. While neural network-based surrogates offer potential acceleration, existing approaches often suffer from gradient inconsistency between predicted objectives and sensitivities, leading to optimization instability. This work presents KATOsuper, an objective-agnostic framework that couples neural-reparameterized topology optimization with a Sensitivity-Consistent Fourier Neural Operator (SC-FNO). The framework employs the forward_split architecture, which derives deployed sensitivities via automatic differentiation through the predicted objective field and thereby preserves consistency between the predicted objective and the gradient used for optimization. The case studies include three 2D benchmark problems and three 3D structures considering compliance or stress minimization. A physics-informed multi-channel input encoding with Fourier position embedding enables resolution-invariant learning, supporting zero-shot extrapolation beyond the training resolution, with useful performance at moderate scaling factors and topology-preserving exploration at up to 64x without retraining. The framework extends to 3D through KATO3D, featuring novel KANConv3D blocks with learnable B-spline activations. KATOsuper demonstrates 15--110x deployment-time speedup over MATLAB baselines while maintaining competitive optimality, with the clearest gains observed in complex 3D and stress-optimization cases. The insight that sensitivity direction matters more than magnitude enables robust optimization even with approximate physics evaluation, extensible to other differentiable physics-driven design objectives.

Shengyu Yan, Jasmin Jelovica

Department of Mechanical Engineering, The University of British Columbia, 6250 Applied Science Ln, Vancouver, V6T 1Z4, Canada · Department of Civil Engineering, The University of British Columbia, 6250 Applied Science Ln, Vancouver, V6T 1Z4, Canada

Neural Surrogate Modeling Topology Optimization Surrogate Modeling Fourier Neural Operator Surrogate-Assisted Optimization

Ask PDF

Sep 18, 2026

Bilevel Optimization of Topology and Hyperparameters (BOTH)

Topology optimization (TO) represents a significant step towards automating the design process: given a working simulation, TO can produce a viable prototype at the press of a button by differentiating the simulation and iteratively improving the design. In practice, however, TO is riddled with

magic numbers''---hyperparameters whose tuning significantly affects the outcome. Finding the right values typically requires not only deep problem-specific knowledge but also extensive trial-and-error. While practitioners can use surrogate-assisted hyperparameter optimization as an alternative, this approach requires strictly limiting the number of hyperparameters through careful problem formulation. Here, we propose differentiating TO itself using automatic differentiation. This yields

hypergradients'' that allow us to tune these hyperparameters in tandem with the primary optimization. We show that evaluating just one or two steps of TO is sufficiently informative and that the method scales favorably to thousands of hyperparameters at an expense comparable to only a few standard TO runs. We demonstrate this approach on stress-constrained and compliance problems, with the latter utilizing a neural parameterization of the density field.

Suryanarayanan Manoj Sanu, Miguel Anibal Bessa, Alejandro Marcos Aragón

Faculty of Mechanical Engineering Delft university of technology Delft The Netherlands · School of engineering Brown University Providence United States of America

Hypergradient Optimization Topology Optimization Hyperparameter Optimization Bilevel Optimization Differentiable Optimization Constrained Optimization

Ask PDF

Jun 18, 2026

eCNNTO: A Highly Generalizable ConvNet for Accelerating Topology Optimization

This work proposes an element-based Convolutional Neural Network (CNN) to accelerate density-based Topology Optimization (TO), termed eCNNTO. TO generally undergoes a large number of iterations, where finite element analysis is performed in every iteration, leading to the efficiency bottleneck especially when dense meshes are used to achieve high-resolution designs. To address this limitation, eCNNTO is proposed to build upon Kallioras et al. (2020), where a Deep Belief Network (DBN) was trained for every element to predict its near-optimal density from its early history, thereby skipping the great majority of iterations and significantly accelerating the TO procedure. However, the method lacks spatial correlations among neighboring elements and may lead to disconnected features in the final structure. The proposed method employs CNN with residual connections to address this issue. On top of it, a novel training strategy is introduced to further enhance the optimization efficiency, where the training dataset consists of the final stage density histories rather than early ones. This change can also help reduce the required training data size. eCNNTO requires only a small dataset to train and yet it can be generalized to problems with largely different boundary conditions, loading cases, design domain geometries, mesh resolutions, as well as non-design domains. In the end, the generalization capabilities and efficiency of eCNNTO are demonstrated through a variety of examples in two and three dimensions, achieving up to 90% and 97% reduction of iterations, respectively.

Shengbiao Lu, Xiaodong Wei

Global college, Shanghai Jiao Tong University, 800 Dongchuan Road, Shanghai, 200240, Shanghai, China

Topology Optimization Deep Learning Optimization Convolutional Neural Networks

Ask PDF

Abstract

Topology Optimization

Explore similar work

KATOsuper: Surrogate-accelerated neural topology optimization with sensitivity-consistent Fourier neural operators

Bilevel Optimization of Topology and Hyperparameters (BOTH)

eCNNTO: A Highly Generalizable ConvNet for Accelerating Topology Optimization