Explainable quantum neural networks for multi-material topology optimization


1 Introduction↩︎

Topology optimization (TO) distributes material within a prescribed design domain so that a structural objective, most commonly compliance, is optimized under equilibrium and resource constraints. It has matured through homogenization, density-based formulations, level-set methods, deformable mesh strategies, and broad computational implementations [1][8]. Multi-material topology optimization (MMTO) considers the assignment of material types in addition to the usual layout decision. This coupled layout–type decision creates material–material interfaces, void–material interfaces, and a larger combinatorial design space, motivating the development of specialized formulations for multi-material structural design [9][15].

Learning-assisted TO is attractive in this setting because repeated MMTO calculations are computationally demanding when boundary conditions, loading scenarios, mesh resolutions, or design requirements change. Convolutional neural networks, generative models, multi-stage predictors, reinforcement-learning strategies, and super-resolution based accelerators have demonstrated that near-optimized topologies can be inferred directly from physical fields, early optimizer states, or problem specifications [16][23]. While such learned predictors can substantially accelerate the design process, they also impose an interpretability requirement: producing a plausible material map is insufficient if the model cannot account for the signals that support its prediction. A mechanics surrogate should make clear whether its decision is driven by load-carrying regions, material interfaces, or features that separate material types, and these signals should remain accessible when the test case departs from the training distribution. This requirement aligns with a broader view of interpretable and physics-informed machine learning, where physical consistency is regarded as a modeling criterion [24][27]. Although generic post-hoc tools such as SHAP and gradient localization are useful in many domains [28], [29], their explanatory vocabulary in compact numerical surrogates is often inherited from the underlying architecture.

Parameterized quantum circuits (PQCs), often referred to as quantum neural networks (QNNs) in learning applications, provide a natural route to such explicitly accessible variables through measured expectation values of selected observables. A PQC encodes classical inputs into a quantum state, transforms the state through trainable unitary operations, and reads out the resulting state by bounded expectation values associated with selected qubits and measurement operators. In this work, each measured expectation value is termed a quantum observable channel; specifically, each channel corresponds to a qubit-wise Pauli-Z expectation, \(\langle Z_i\rangle\), obtained after feature encoding and trainable circuit transformation.

This observable-channel formulation is well suited to MMTO because the assignment of material types is not governed solely by scalar density values. It depends on competing material types, interface adjacent elements, and coupled mechanical descriptors. Compressing all measured qubit responses into a single scalar can obscure circuit-level information generated through superposition and entangling operations. Retaining qubit-wise observables instead preserves multiple bounded response channels, giving the readout layer independent degrees of freedom for multi-material classification while keeping the intermediate variables explicitly measurable.

This perspective extends earlier work on quantum-assisted topology optimization. PQCs have been shown to accelerate topology prediction by using optimizer-derived strain energy, sensitivity, and density descriptors from early TO iterations [30]. The binary-topology formulation used fixed early optimizer snapshots and converted a scalar quantum response into topology decisions through an empirical threshold. Here, we propose an explainable quantum neural network (XQNN) that retains the accelerator-oriented motivation of that approach but adapts it to MMTO. More specifically, fixed early snapshots are replaced by mechanically defined milestones, \(k_{\rm em}\) and \(k_{\rm em}+5\), to provide a richer post-emergence history, together with the prediction decision milestone \(k_{\rm pd}\); a scalar quantum response is replaced by qubit-wise \(Z\) observables followed by a softmax readout; and a Sobel boundary hint is introduced to provide interface sensitive information. Variational quantum latent encoding has also been explored for TO, where quantum circuits form compact latent representations of design information for optimization-oriented prediction [31]. These studies motivate the central question of this work: can a compact PQC generate measured observables that remain predictive after readout and can it also be interrogated in mechanics-facing terms?

2 Results↩︎

The learning task is element-wise material classification for MMTO layouts generated by the method of moving asymptotes (MMA) where each finite element is assigned to one of the available solid material types or to void. We denote the number of solid material types by \(M\) and study two-, three-, and four-material cases. Separately, the grayscale ratio \(G(k)\) is used only to quantify the discreteness of an evolving material field during MMA optimization and to select intermediate iterations, where \(k\) denotes the MMA iteration index. The problems and dataset are organized into non-overlapping roles (Tables S1 and S2). Training problems \(T_1\)\(T_3\) are used to fit the circuit and readout parameters, validation problems \(V_1\)\(V_3\) are used for model selection and early stopping, and out-of-distribution (OOD) test problems \(O_1\)\(O_4\) are used to evaluate the performance of the selected model under boundary/loading conditions unseen during training. The high-resolution transfer cases \(\mathit{HR}_1\)\(\mathit{HR}_3\), defined by progressively refined meshes, and the three-dimensional transfer cases \(\mathit{3D}_1\)\(\mathit{3D}_4\) are designed to demonstrate the generalization performance of XQNN. Across all material layout plots, colors represent the material type: blue, orange, yellow, and green for material types 1, 2, 3, and 4, respectively, with white for void.

2.1 Model architecture of the XQNN↩︎

The XQNN representation keeps the quantum output as a dictionary of measured observables rather than collapsing it into one scalar. A pure qubit can be written as \[|\psi(\vartheta,\varphi)\rangle= \cos\frac{\vartheta}{2}|0\rangle+ e^{i\varphi}\sin\frac{\vartheta}{2}|1\rangle, \label{eq:bloch95results}\tag{1}\] where \(\vartheta\) and \(\varphi\) are the polar and azimuthal angles on the Bloch sphere (Figure 1a). After feature encoding and ansatz transformation, each qubit is measured through a \(Z\) expectation, \[\langle Z_i\rangle=p_i(0)-p_i(1), \label{eq:z95expectation95results}\tag{2}\] which reports the population contrast of qubit \(i\). This bounded expectation value is the quantum observable channel used by the readout layer.

Hadamard and phase gates form the ZFeatureMap, and trainable \(R_Y\) rotations with CNOT entanglers form the RealAmplitudes ansatz (Figure 1b). The Hadamard and phase gates introduce descriptor-dependent quantum phases, whereas the \(R_Y\) rotations and CNOT entanglers mix the encoded state before the final \(Z\) measurements. The MMTO histories and preprocessing used to construct the descriptors are described first in Supplementary Notes 1 and 2, and the corresponding gate-level matrices and circuit ordering are given afterward in Supplementary Note 3.

Figure 1c defines the ten-qubit XQNN input. Each element is represented by ten descriptors: strain energy (\(E_e\)), sensitivity (\(S_e\)), and density (\(D_e\)) sampled at the emergence iteration \(k_{\rm em}\), five MMA steps later at \(k_{\rm em}+5\), and the prediction decision iteration \(k_{\rm pd}\), together with a Sobel boundary hint (\(B_e\)) computed from the prediction stage density field. The descriptor vector is then \[\begin{align} \mathbf{x}_e = \big[& E_e(k_{\rm em}),\, S_e(k_{\rm em}),\, D_e(k_{\rm em}),\\ & E_e(k_{\rm em}+5),\, S_e(k_{\rm em}+5),\, D_e(k_{\rm em}+5),\\ & E_e(k_{\rm pd}),\, S_e(k_{\rm pd}),\, D_e(k_{\rm pd}),\, B_e(k_{\rm pd}) \big]^T . \end{align} \label{eq:descriptor95results}\tag{3}\] Adapting the optimizer history strategy [30] to MMTO, the same relative MMA milestones are used in both training and prediction. Here, \(k_{\rm em}\) replaces a fixed first early snapshot as the first mechanically meaningful emergence state, \(k_{\rm em}+5\) replaces the next early snapshot as a richer post-emergence state, and \(k_{\rm pd}\) provides the early stopped prediction decision state. During inference, a new MMA run is stopped at \(k_{\rm pd}\), and descriptors from \(k_{\rm em}\), \(k_{\rm em}+5\), and \(k_{\rm pd}\) are passed to the trained XQNN.

The raw descriptor values are normalized with a training set min–max transform before quantum encoding. Ground-truth (GT) material type labels are refined with a \(3\times3\) mode filter, and the scaling transform is fitted independently for each material type setting because the physical ranges of strain energy, sensitivity, density, and Sobel boundary hint values differ across the two-, three-, and four-material simulations. For material setting \(M\) and descriptor/qubit index \(q\), \[\tilde{x}_{e,q}=\operatorname{clip}\!\left(\mathcal{S}_{M,q}(x_{e,q}),0,1\right), \qquad \alpha_{e,q}=\pi\tilde{x}_{e,q}. \label{eq:angle95clip95results}\tag{4}\] The pre-fitted transform is reused without refitting for OOD, high-resolution, and three-dimensional inference, and clipping keeps the encoded angle in the valid feature map range \([0,\pi]\).

The normalized descriptor is encoded by the ZFeatureMap with REPS=2, transformed by the RealAmplitudes ansatz with REPS=3, and measured as \[\mathbf{Q}_e=\big(\langle Z_0\rangle_e,\langle Z_1\rangle_e,\ldots,\langle Z_9\rangle_e\big)^T. \label{eq:observable95results}\tag{5}\] During training (Figure 1d), class balanced element samples from the selected MMA histories optimize the circuit parameters and readout weights. During prediction (Figure 1e), a new topology optimization run is stopped at \(k_{\rm pd}\), each element is encoded by the trained circuit, and the readout layer maps \(\mathbf{Q}_e\) to material type probabilities, \[\mathbf{z}_e=\mathbf{W}^T\mathbf{Q}_e+\mathbf{b}, \qquad p_{e,c}=\frac{\exp(z_{e,c})}{\sum_{r=0}^{C-1}\exp(z_{e,r})}, \qquad \hat{y}_e=\arg\max_c p_{e,c}. \label{eq:softmax95results}\tag{6}\] The readout has \(C=M+1\) classes for \(M\) solid material types plus void. The remaining implementation details, including the exact statevector evaluation, training schedule, and parameter counts, are given in Methods and Supplementary Note 3.

2.2 Model training and validation↩︎

We train the model using MMA histories of two-dimensional MMTO problems solved with a \(50\times25\) mesh. Three boundary/loading conditions are used to obtain the training dataset: (\(T_1\)) under a vertical point load at the right edge while fixed at the left edge, (\(T_2\)) under a vertical point load at the center of the top edge while fixed at both left and right edges, and (\(T_3\)) under a horizontal point load at the upper right corner while fixed at the bottom edge (Figure 2a). The optimizer fields sampled at \(k_{\rm em}\), \(k_{\rm em}+5\), and \(k_{\rm pd}\) are used for training after being converted into element-wise XQNN descriptors.

The validation set is constructed using three additional boundary/loading conditions: (\(V_1\)) under a vertical point load at the center of the top edge while simply supported at the bottom corners, (\(V_2\)) under a distributed vertical load on the central part of the top edge while fixed at the four corners, and (\(V_3\)) under two vertical point loads at the right corners while fixed at the left edge (Figure 2b). Validation is used only for model selection and early stopping. The validation predictions retain the dominant load paths and material type regions with high material type accuracy across \(V_1\)\(V_3\), showing that the selected intermediate fields carry useful information for boundary/loading conditions that are not used for fitting the hyperparameters. The two- and three-material predictions reproduce the main structural members well, with some local shifts near the void–solid or material–material boundaries. The four-material case is visually less clean, but the main load path and the largely clustered material regions remain recognizable. The accuracy ranges from 81.4% to 94.0%.

The model can be trained stably regardless of the number of material types used (Figure 2c). The training loss is evaluated on the class balanced random subset sampled at each epoch, while the validation accuracy is calculated by taking the mean of material type overlap for all elements in \(V_1\)\(V_3\). Nevertheless, convergence of the loss function becomes slower and the validation accuracy decreases as the number of material types increases, because it gets harder for the same compact observable dictionary to separate more material types and material–material interfaces.

Increasing the number of training designs further improves or stabilizes the validation behavior, but three training designs already provide sufficient validation accuracy (Figure S1 and Table S3). Training loss decreases as more optimization histories are supplied while the mean validation accuracy changes modestly only and can even appear highest when only \(k_{\rm pd}\) is used (Figure S2 and Table S4). More importantly, high-resolution transfer tests with altered stiffness settings in a two-material case show that using multiple milestones, \(k_{\rm em}\), \(k_{\rm em}+5\), and \(k_{\rm pd}\), produces more stable holes, material paths, and connector regions than using a single late snapshot. Varying the depth of the RealAmplitudes indicates that REPS=3 balances expressivity and stable prediction, while deeper circuits do not consistently improve the validation cases (Figure S3 and Table S5). The target grayscale threshold \(G_{\rm target}\) controls a trade-off between the computing time and solution accuracy (Figure S4 and Table S6). Lowering \(G_{\rm target}\) moves \(k_{\rm pd}\) closer to the final optimized layout and usually improves validation accuracy, but also reduces the benefit of early stopping. We select \(G_{\rm target}=0.20\), 0.15, and 0.10 for the two-, three-, and four-material cases, respectively, as a compromise that retains a clear early stopping benefit while providing sufficient boundary clarity for prediction. Feature space and prediction distribution diagnostics further show that confidence broadens and material type frequency imbalance becomes more pronounced in the four-material case, consistent with its more distributed errors (Figure S5).

2.3 Generalization performance↩︎

We evaluate the generalization performance of the trained model using two-dimensional OOD boundary/loading conditions, various mesh resolutions, and three-dimensional TO problems (Figure 3). In addition to pixel-wise accuracy of assigned material types, we calculate the normalized compliance, \(C_{\rm pred}/C_{\rm GT}\), indicating the relative stiffness of the inferred structure compared to that of the GT reference solution. The full GT–prediction layout pairs for these problems are provided in Figures S6–S8 for the two-, three-, and four-material settings, respectively.

OOD boundary/loading conditions We apply four additional boundary/loading conditions, denoted \(O_1\)\(O_4\), to the same \(50\times25\) mesh to assess OOD performance (Figure 3a). The scatter plot compares the accuracy of material types and normalized compliance for different numbers of materials across four cases. Overall, the two-material case exhibits relatively high accuracy and compliance values close to unity, indicating stable predictive performance. In contrast, the three- and four-material cases show larger problem-dependent variations, as reflected by the broader shaded regions. The case \(O_2\) in the two-material TO shows a noticeable deviation in normalized compliance while achieving high accuracy. This turns out to be due to a local shift at the material interface rather than a deficient structural layout (Figure S6). The other three cases with relatively high values of normalized compliance reveal a similar aspect in that the principal load-carrying skeleton is well predicted by XQNN but local errors in material type assignments particularly near the material interfaces deteriorate the structural stiffness. As the number of material types increases, the agreement between predicted and GT compliance becomes more sensitive to the problem because the number of material interfaces increases naturally.

Mesh resolution We test the generalization performance of the XQNN model by applying it, without additional training, to the \(V_1\) case solved using three higher-resolution meshes, \(\mathit{HR}_1\)(\(100\times50\)), \(\mathit{HR}_2\)(\(150\times75\)), and \(\mathit{HR}_3\)(\(200\times100\)) (Figure 3b). Note that this test evaluates the mesh-scale transfer of our element-wise observable model, rather than the fixed-domain mesh convergence. Across resolutions, the inferred structures retain arch-like load paths and central voids, indicating that XQNN is not dependent on the mesh resolution at which it is trained. For the three- and four-material cases, XQNN provides stable solutions with good material type accuracies (ranging from 73% to 91%) and normalized compliance values close to one. No clear mesh dependency is observed. However, the two-material cases become less stable, particularly when higher mesh resolutions are used. This is because wrong assignments of the material type to voids near material interfaces result in weakly connected or disconnected parts in the inferred structure, making the structure too flexible or statically unstable. This structural disconnection becomes more probable in a higher-resolution mesh as the number of elements increases while less probable with more material types. Therefore, this test clearly shows both the strength and the limitation of element-wise transfer. Global material composition and load paths can be preserved, whereas thin or disconnected members due to material type errors at the material interfaces remain the dominant source of mechanical sensitivity.

Three-dimensional MMTO We further apply the trained XQNN to three-dimensional TO problems without new training. Here, we use a \(20\times20\times10\) mesh with four different boundary/loading conditions, \(\mathit{3D}_1\)\(\mathit{3D}_4\), to assess its generalization performance in a three-dimensional setting (Figure 3c). The parameters learned from two-dimensional training data are applied voxel by voxel to the descriptors after a three-dimensional MMA run is stopped at \(k_{\rm pd}\). XQNN reproduces the main three-dimensional load transfer regions without additional three-dimensional labels or a new network architecture. The representative \(\mathit{3D}_2\) predictions show this trend for all material settings. The dominant three-dimensional and internal load paths are reproduced well for all cases, while the three- and four-material predictions lose finer material details around internal interfaces unlike the clean two-material prediction. Across the full \(\mathit{3D}_1\)\(\mathit{3D}_4\) sets (Figures S6–S8), the mean accuracy decreases from 89.3% to 79.6% and 73.6% as the number of material types increases. Interestingly, the corresponding mean normalized compliance values are 0.96, 0.98, and 0.83, respectively, indicating that the inferred three-dimensional structures are often stiffer than the GT references despite lower voxel-level accuracies. This is probably because the three-dimensional MMA references can be more prone to local minima than the two-dimensional training solutions, and the learned two-dimensional regularities may sometimes guide XQNN toward mechanically stiffer three-dimensional layouts. However, this effect should be interpreted together with the current limitation of XQNN in that it predicts material type labels element-wise and does not explicitly enforce the global volume constraints during inference.

2.4 Observable-level explanation and faithfulness test↩︎

The explainable QNN analysis is grounded in the final readout equation, which is linear in the measured observables. For element \(e\), the margin between the winning class \(c_1\) and the nearest competing class \(c_2\) is \[\Delta z_e=z_{e,c_1}-z_{e,c_2} =(b_{c_1}-b_{c_2})+\sum_{j=0}^{9}Q_{e,j}\left(W_{j,c_1}-W_{j,c_2}\right). \label{eq:margin95decomposition95results}\tag{7}\] Conditioned on the trained PQC output \(\mathbf{Q}_e\), the readout layer margin is therefore exactly decomposed across observable channels. We define the signed channel contribution and non-negative quantum observable attribution map (QOAM) magnitude as \[m_{e,j}=Q_{e,j}\left(W_{j,c_1}-W_{j,c_2}\right), \qquad a_{e,j}=|m_{e,j}|, \qquad A_e=\sum_{j=0}^{9}a_{e,j}. \label{eq:qoam95results}\tag{8}\] The displayed margin QOAM map is the normalized spatial field \(A_e\). High QOAM therefore has a precise meaning: the measured observables at that element make a large contribution to the local class margin.

The attribution maps concentrate on mechanically relevant and decision-sensitive regions (Figure 4a). High-QOAM regions appear along load-carrying members, material-transition zones, and interface adjacent regions where local material decisions strongly affect the class margin. The two-material case produces relatively localized attributions along the main load path and its boundaries, whereas the three- and four-material cases produce more distributed patterns. This broadening is mechanically plausible because additional material types introduce more material–material interfaces and more local label competition.

Specialist-observable analysis connects individual measured channels to mechanics-derived fields (Figure 4b). Four targets are used: strain energy, negative sensitivity, final material type label, and interface strength. Here, the interface strength denotes the local finite-difference magnitude of the discrete material type label, used as a proxy for material-boundary intensity rather than as a physical interface length. For each target, the observable channel with the largest absolute Pearson correlation in the validation set is selected and visualized on a representative OOD case. The selected channels and correlation coefficients are summarized in Table 1. In the two- and three-material cases, several observables behave as localized specialist detectors: strain energy, sensitivity, and material type labels can be associated with individual channels with high values of \(|R|\). This provides a direct mechanics-facing interpretation of part of the learned representation.

The same table also shows why a low single channel Pearson correlation should not be read as a failure of the quantum representation. As the number of material types increases, or when the target is the highly local and nonlinear field of interface strength, the largest individual \(|R|\) values become smaller. This trend suggests a transition from localized to distributed representation. The ZFeatureMap and CNOT-entangling RealAmplitudes ansatz project and mix the input descriptors across the multi-qubit state before measurement, so information about one mechanical quantity need not remain confined to one same index qubit. The final linear readout can then combine collectively entangled observable responses to recover nonlinear material type and interface relationships. The specialist maps are therefore not intended to be perfect reconstructions of the targets. They test which single observables remain visibly aligned with mechanics, while QOAM and perturbation tests examine how the full observable dictionary supports the decision.

Table 1: Validation-selected specialist observables and Pearson correlation coefficients. For each mechanics-derived target, the table reports the observable with the largest absolute Pearson correlation on the unified validation set.
Target 2mat 3mat 4mat
Strain energy \(\langle Z_0\rangle\) (\(R=-0.833\)) \(\langle Z_1\rangle\) (\(R=0.899\)) \(\langle Z_3\rangle\) (\(R=-0.554\))
Sensitivity \(\langle Z_1\rangle\) (\(R=-0.603\)) \(\langle Z_1\rangle\) (\(R=0.844\)) \(\langle Z_3\rangle\) (\(R=-0.518\))
Material type label \(\langle Z_8\rangle\) (\(R=0.938\)) \(\langle Z_9\rangle\) (\(R=-0.881\)) \(\langle Z_9\rangle\) (\(R=-0.847\))
Interface strength \(\langle Z_9\rangle\) (\(R=0.345\)) \(\langle Z_3\rangle\) (\(R=-0.186\)) \(\langle Z_9\rangle\) (\(R=-0.300\))

3pt

Faithfulness is tested by perturbing observables and inputs after training (Figure 5). Observable erasure replaces selected channels by neutral-reference values computed from a constant mid-range input vector. Masking high-QOAM channels produces a stronger erasure response than masking low-QOAM channels, with random masking between the two targeted cases (Figure 5a). This ordering is important because it tests the explanation, not only the prediction. Channels identified by the QOAM margin decomposition are the channels whose removal most strongly changes the classifier response. QOAM therefore highlights decision-relevant observables rather than visually plausible but unused features.

Post-hoc sparsity and input feature ablation provide complementary checks (Figure 5b). When only a small number of observables are retained, the channels ranked highly by readout-weight magnitude outperform random channel subsets, showing that the trained observable dictionary contains a compact decision core. The sparsity curves also explain why the response can plateau after several top-ranked channels have been retained: once the dominant margin-carrying observables are present, additional low-ranked channels mainly add redundant or case-specific information and may not change the predicted class for many elements. The four-material task requires a broader subset than the two- and three-material tasks, which is consistent with the larger number of material types and material–material interfaces that must be separated. Input-feature ablation further shows that strain energy, sensitivity, density, and the Sobel boundary hint play different roles. Energy and sensitivity support load-path recognition, density carries the evolving material layout, and the Sobel boundary hint contributes interface sensitive information. No single descriptor group explains all predictions, and the limited effect of removing some groups in selected cases is consistent with distributed encoding, where related mechanical information can be partially recovered from other observables. Together with the validation set specialist maps, these perturbation tests make the explanations causal at the readout level and mechanically interpretable at the input feature level.

3 Discussion↩︎

The results support an interpretation of XQNN: it is an observable-based surrogate for MMTO, not a black-box quantum classifier. The MMA solver supplies mechanics-derived element descriptors, the ZFeatureMap writes these descriptors into quantum phase-encoded qubit states, the RealAmplitudes ansatz transforms those states, and the readout layer uses measured \(Z\) observables to predict material type labels. Retaining qubit-wise observables is central to the method because the same quantities used for prediction can also be decomposed by QOAM, correlated with mechanics-derived targets, erased, compressed, and input-ablated.

Three conclusions follow from this design. First, a compact ten-qubit model can capture structural information that transfers across OOD two-dimensional boundary/loading conditions and mesh resolutions, and it can be applied to voxelized three-dimensional domains without additional three-dimensional training. Transfer is imperfect, especially in four-material cases, but the dominant load-carrying topology is usually preserved while errors concentrate near interfaces and thin material regions. The accuracy–compliance pairs are important here: pixel-wise material type accuracy can remain high even when a small connectivity error strongly affects compliance, so both metrics are needed to evaluate a mechanics surrogate. Second, retaining qubit-wise observables provides an auditable intermediate representation rather than a single compressed scalar response. The separate measured \(Z\) channels preserve independent readout degrees of freedom generated after superposition, phase encoding, and entangling transformations, allowing the classifier and the explanation analyses to operate on the same measurable variables. Third, the Sobel boundary hint contributes a complementary signal rather than a complete explanation. It gives the tenth qubit local edge-sensitive information, while the ablation results show that strain energy, sensitivity, density, and the Sobel boundary hint are all used, with their relative importance changing across the number of material types.

The specialist observable results clarify how the quantum circuit should be interpreted. A tempting but incorrect reading would be to associate each measured qubit directly with the same index input feature. The validation correlations show a different picture: after quantum phase encoding and CNOT-mediated mixing, observables become circuit-level responses to the descriptor vector. Strong single channel correlations in simpler cases provide localized interpretability, whereas lower single channel correlations for the four-material and interface-strength targets indicate that the relevant information is increasingly distributed across the observable dictionary. This distributed representation is useful because the readout can reconstruct complex material types and interface decisions from collective observable responses, but it also requires explanation methods that operate on the measured observables rather than on the raw inputs alone. QOAM, specialist alignment, and perturbation tests therefore play complementary roles: QOAM decomposes the decision margin, specialist alignment relates observables to mechanics-derived targets, and erasure/sparsity tests whether the identified channels are actually used by the classifier.

The main limitation of the current work is that all quantum evaluations use exact statevector simulation. Real devices would introduce shot noise, gate errors, and decoherence, and hence, it should be tested whether the observable-level explanations remain stable under those effects [32]. Nevertheless, the proposed framework is not limited to MMTO only. Any mechanics surrogate that converts local descriptors into point-wise or element-wise labels can use the same design principle: retain a compact dictionary of measured observables, train a readout layer, and interrogate the resulting decisions with margin decomposition, specialist alignment, and perturbation tests.

4 Methods↩︎

4.1 MMTO simulations and dataset organization↩︎

Ground-truth (GT) layouts are generated with an MMA-based multi-material topology optimization solver adapted from the MATLAB code of Zheng et al., which provides complete two-dimensional and three-dimensional implementations for compliance minimization of multi-material continuum structures using mapping-based interpolation, density filtering, projection, and sensitivity analysis [15]. In our dataset generation, the optimization minimizes compliance subject to finite element equilibrium, box constraints, and one volume constraint per material type. The material type settings, relative Young’s moduli, and volume constraints used for the main two-dimensional datasets are summarized in Table S2. The terminal GT/reference layouts are generated by running MMA until the grayscale ratio falls below 0.02 for the regular two-dimensional and three-dimensional simulations and below 0.01 for the high-resolution two-dimensional simulations. These terminal criteria define the converged reference layouts, whereas \(G_{\rm target}\) defines the earlier prediction decision state used for XQNN inference. Details of the material interpolation, finite element discretization, filtering, and three-dimensional problem settings are given in Supplementary Note 1. We use the disjoint training, validation, and final evaluation roles defined in Results and Figure 23. The same split is used throughout XQNN training and evaluation, and the boundary/loading definitions and mesh sizes are summarized in Table S1. This Methods section therefore specifies how the split is used in simulation, preprocessing, and training rather than repeating the full problem definitions.

4.2 Preprocessing and descriptor construction↩︎

The terminal MMTO material type label fields are refined with a \(3\times3\) mode filter, implemented as majority voting over the local label neighborhood, to suppress isolated high-frequency numerical artifacts and checkerboard-like patterns along material type boundaries. The input descriptor for each element consists of strain energy, sensitivity, and density sampled at \(k_{\rm em}\), \(k_{\rm em}+5\), and \(k_{\rm pd}\), together with a Sobel boundary hint evaluated from the prediction stage density field. The emergence iteration \(k_{\rm em}\) is selected using \(G_{\rm em}=0.5\), and \(k_{\rm pd}\) is selected using \(G_{\rm target}=0.20\), 0.15, and 0.10 for the two-, three-, and four-material cases, respectively. These \(G_{\rm target}\) values are grayscale-ratio thresholds and are not used as a synonym for the number of material types.

For each two-, three-, or four-material case, min–max scaling to \([0,1]\) is fitted independently and only on the corresponding training designs. The fitted transform is reused unchanged for validation, OOD, high-resolution, and three-dimensional inference. Scaled descriptors are clipped to \([0,1]\) and then mapped to quantum feature map angles in \([0,\pi]\). This training-only scaling prevents data leakage from evaluation distributions while preserving the valid input range of the ZFeatureMap. Detailed descriptor equations and Sobel definitions are provided in Supplementary Note 2.

4.3 Model training and analysis↩︎

The XQNN circuit uses ten qubits, a ZFeatureMap with REPS=2, and a RealAmplitudes ansatz with REPS=3. The ZFeatureMap uses Hadamard and phase gates, and the RealAmplitudes ansatz uses trainable \(R_Y\) rotations with CNOT entanglers, as shown in Figure 1b,c. Exact statevector simulation is used for all main results, and the measured qubit-wise \(Z\) observables are passed to a linear softmax readout. Training minimizes class balanced multi-class cross-entropy with Adam optimization. At each epoch, class balanced random sampling draws up to 20 elements per class from each of the three training designs, giving 180, 240, and 300 sampled elements per epoch for the two-, three-, and four-material tasks, respectively. The model with the highest mean validation accuracy on the full \(V_1\)\(V_3\) validation set is retained. Details of the training schedule, parameter shift gradients, sampling protocol, and configuration studies are provided in Supplementary Note 3.

All quantum evaluations use exact statevector simulation in Python 3.12.10 with Qiskit SDK v2.3.1 and Qiskit Aer v0.17.2 on an AMD Ryzen 7 9700X 8-Core Processor. The QOAM, specialist observable, erasure, sparsity, and input feature ablation protocols are summarized in Results and detailed in Supplementary Note 4.

5 Data availability↩︎

The data that support the plots within this paper and other findings of this study are available from the corresponding author upon reasonable request.

6 Code availability↩︎

The computer codes used to generate the results of this paper will be made available on request.

7 Acknowledgements↩︎

This work was supported by the National Research Foundation of Korea (NRF) grant funded by the Ministry of Science and ICT (RS-2024-00346176) and the Development of Core Technologies for Manufacturing Foundation Models Program (RS-2025-25458052) funded by the Ministry of Trade, Industry & Energy (MOTIE, Korea). D.J. also thanks Minwoo Kim at Seoul National University for fruitful discussions and helpful comments.

8 Author Contributions↩︎

D.J. conceived the project, designed the research, proposed the algorithm, and assembled input data; D.J. and D.-N.K. analyzed output data and wrote the initial manuscript; N.S. and K.T. provided the prior research code, guidance on the test direction, and result validation; D.-N.K. supervised the project; all authors contributed to the writing of the manuscript.

9 Competing Interests↩︎

The authors declare no competing interests.

Figure 1: XQNN workflow for multi-material topology prediction. a, Bloch-sphere representation of a single qubit. b, Gate primitives used by the implemented circuit: Hadamard and phase gates in the ZFeatureMap, and trainable R_Y rotations with CNOT entanglers in the RealAmplitudes ansatz. c, Ten-qubit circuit structure. The first nine qubits encode strain energy, sensitivity, and density descriptors from k_{\rm em}, k_{\rm em}+5 and k_{\rm pd}, while the tenth qubit encodes a Sobel boundary hint. The resulting state is measured through qubit-wise Z observables. d, Training path from selected MMA iterations to class balanced element-wise sampling and multi-class cross-entropy optimization. e, Prediction path in which a new MMA run is stopped at k_{\rm pd}, each element is encoded by the trained circuit, and the trained readout layer returns material type probabilities.
Figure 2: Training, validation, and convergence behavior. a, Training boundary/loading conditions T_1–T_3 and selected MMA fields at k_{\rm em}, k_{\rm em}+5, and k_{\rm pd} used as XQNN input descriptors. b, Validation conditions V_1–V_3 with GT and XQNN-inferred material type layouts; validation is used only for model selection and early stopping. c, Training loss and validation accuracy plotted with common axes; dashed lines denote the mean validation accuracy over V_1–V_3.
Figure 3: Generalization performance under OOD, high-resolution, and three-dimensional transfer. a, OOD tests O_1–O_4 on the same 50\times25 mesh and the corresponding scatter plot of material type accuracy versus normalized compliance C_{\rm pred}/C_{\rm GT} for all material settings. b, High-resolution transfer cases \mathit{HR}_1(100\times50), \mathit{HR}_2(150\times75), and \mathit{HR}_3(200\times100), with representative XQNN predictions for two-, three-, and four-material cases. c, Three-dimensional transfer conditions \mathit{3D}_1–\mathit{3D}_4 and representative \mathit{3D}_2 predictions on a 20\times20\times10 voxel mesh.
Figure 4: QOAM and mechanics-facing feature alignment. a, GT fields, XQNN predictions, and margin QOAM maps for two-, three-, and four-material cases. Margin QOAM decomposes the readout layer class margin into observable channel contributions and highlights regions where measured quantum observables affect the local material decision. b, Validation-selected specialist observables aligned with mechanics-derived target fields, including strain energy, sensitivity, material type label, and interface strength. Strong single channel correlations indicate localized specialist observables, whereas weaker correlations for more complex or interface dominated targets indicate distributed circuit-level representations.
Figure 5: Observable-level faithfulness, sparsity, and input feature ablation. a, Observable-level erasure response under high-QOAM, low-QOAM, and random channel masking. b, Post-hoc sparsity obtained by retaining only the highest-ranked observable channels; plateau regions indicate that the dominant margin-carrying observables have already been retained. The right panel shows input feature ablation obtained by removing strain energy, sensitivity, density, or the Sobel boundary hint.

References↩︎

[1]
M. P. Bendsøe and N. Kikuchi, “Generating optimal topologies in structural design using a homogenization method,” Computer Methods in Applied Mechanics and Engineering, vol. 71, no. 2, pp. 197–224, 1988, doi: 10.1016/0045-7825(88)90086-2.
[2]
M. P. Bendsøe and O. Sigmund, Topology optimization: Theory, methods, and applications. Berlin: Springer, 2003.
[3]
O. Sigmund, “A 99 line topology optimization code written in MATLAB,” Structural and Multidisciplinary Optimization, vol. 21, no. 2, pp. 120–127, 2001, doi: 10.1007/s001580050176.
[4]
O. Sigmund and K. Maute, “Topology optimization approaches: A comparative review,” Structural and Multidisciplinary Optimization, vol. 48, no. 6, pp. 1031–1055, 2013, doi: 10.1007/s00158-013-0978-6.
[5]
J. D. Deaton and R. V. Grandhi, “A survey of structural and multidisciplinary continuum topology optimization: Post 2000,” Structural and Multidisciplinary Optimization, vol. 49, no. 1, pp. 1–38, 2014, doi: 10.1007/s00158-013-0956-z.
[6]
N. P. van Dijk, K. Maute, M. Langelaar, and F. van Keulen, “Level-set methods for structural topology optimization: A review,” Structural and Multidisciplinary Optimization, vol. 48, no. 3, pp. 437–472, 2013, doi: 10.1007/s00158-013-0912-y.
[7]
G. Allaire, F. Jouve, and A.-M. Toader, “Structural optimization using sensitivity analysis and a level-set method,” Journal of Computational Physics, vol. 194, no. 1, pp. 363–393, 2004, doi: 10.1016/j.jcp.2003.09.032.
[8]
K. Jung and D.-N. Kim, “Density-based topology optimization using a deformable mesh,” Computers & Structures, vol. 316, p. 107879, 2025, doi: 10.1016/j.compstruc.2025.107879.
[9]
R. Tavakoli and S. M. Mohseni, “Alternating active-phase algorithm for multimaterial topology optimization problems: A 115-line MATLAB implementation,” Structural and Multidisciplinary Optimization, vol. 49, no. 4, pp. 621–642, 2014, doi: 10.1007/s00158-013-0999-1.
[10]
W. Zuo and K. Saitou, “Multi-material topology optimization using ordered SIMP interpolation,” Structural and Multidisciplinary Optimization, vol. 55, no. 2, pp. 477–491, 2017, doi: 10.1007/s00158-016-1513-3.
[11]
N. Vermaak, G. Michailidis, G. Parry, R. Estevez, G. Allaire, and Y. Bréchet, “Material interface effects on the topology optimization of multi-phase structures using a level set method,” Structural and Multidisciplinary Optimization, vol. 50, no. 4, pp. 623–644, 2014, doi: 10.1007/s00158-014-1074-2.
[12]
Y. Wang, Z. Luo, Z. Kang, and N. Zhang, “A multi-material level set-based topology and shape optimization method,” Computer Methods in Applied Mechanics and Engineering, vol. 283, pp. 1570–1586, 2015, doi: 10.1016/j.cma.2014.11.002.
[13]
A. M. Mirzendehdel and K. Suresh, “A pareto-optimal approach to multimaterial topology optimization,” Journal of Mechanical Design, vol. 137, no. 10, p. 101701, 2015, doi: 10.1115/1.4031088.
[14]
D. Li and I. Y. Kim, “Multi-material topology optimization for practical lightweight design,” Structural and Multidisciplinary Optimization, vol. 58, no. 3, pp. 1081–1094, 2018, doi: 10.1007/s00158-018-1953-z.
[15]
R. Zheng, B. Yi, X. Peng, and G.-H. Yoon, “An efficient code for the multi-material topology optimization of 2D/3D continuum structures written in Matlab,” Applied Sciences, vol. 14, no. 2, p. 657, 2024, doi: 10.3390/app14020657.
[16]
R. V. Woldseth, N. Aage, J. A. Bærentzen, and O. Sigmund, “On the use of artificial neural networks in topology optimisation,” Structural and Multidisciplinary Optimization, vol. 65, no. 10, p. 294, 2022, doi: 10.1007/s00158-022-03347-1.
[17]
Y. Yu, T. Hur, J. Jung, and I.-H. Jang, “Deep learning for determining a near-optimal topological design without any iteration,” Structural and Multidisciplinary Optimization, vol. 59, no. 3, pp. 787–799, 2019, doi: 10.1007/s00158-018-2101-5.
[18]
Z. Nie, T. Lin, H. Jiang, and L. B. Kara, TopologyGAN: Topology optimization using generative adversarial networks based on physical fields over the initial domain,” Journal of Mechanical Design, vol. 143, no. 3, p. 031715, 2021, doi: 10.1115/1.4049533.
[19]
D. Bielecki, D. Patel, R. Rai, and G. F. Dargush, “Multi-stage deep neural network accelerated topology optimization,” Structural and Multidisciplinary Optimization, vol. 64, no. 6, pp. 3473–3487, 2021, doi: 10.1007/s00158-021-03028-5.
[20]
L. Xue, J. Liu, G. Wen, and H. Wang, “Efficient, high-resolution topology optimization method based on convolutional neural networks,” Frontiers of Mechanical Engineering, vol. 16, no. 1, pp. 80–96, 2021, doi: 10.1007/s11465-020-0614-2.
[21]
N. K. Brown, A. P. Garland, G. M. Fadel, and G. Li, “Deep reinforcement learning for engineering design through topology optimization of elementally discretized design domains,” Materials & Design, vol. 218, p. 110672, 2022, doi: 10.1016/j.matdes.2022.110672.
[22]
Y. Gao, S. Zhou, and M. Li, “Structural topology optimization based on diffusion generative adversarial networks,” Engineering Applications of Artificial Intelligence, vol. 138, p. 109444, 2024, doi: 10.1016/j.engappai.2024.109444.
[23]
J. Lim, K. Jung, Y. Jung, and D.-N. Kim, “Accelerating topology optimization using deep learning-based image super-resolution,” Engineering Applications of Artificial Intelligence, vol. 133, p. 108370, 2024, doi: 10.1016/j.engappai.2024.108370.
[24]
C. Rudin, “Stop explaining black box machine learning models for high stakes decisions and use interpretable models instead,” Nature Machine Intelligence, vol. 1, no. 5, pp. 206–215, 2019, doi: 10.1038/s42256-019-0048-x.
[25]
R. Roscher, B. Bohn, M. F. Duarte, and J. Garcke, “Explainable machine learning for scientific insights and discoveries,” IEEE Access, vol. 8, pp. 42200–42216, 2020, doi: 10.1109/ACCESS.2020.2976199.
[26]
A. B. Arrieta et al., “Explainable artificial intelligence (XAI): Concepts, taxonomies, opportunities and challenges toward responsible AI,” Information Fusion, vol. 58, pp. 82–115, 2020, doi: 10.1016/j.inffus.2019.12.012.
[27]
G. E. Karniadakis, I. G. Kevrekidis, L. Lu, P. Perdikaris, S. Wang, and L. Yang, “Physics-informed machine learning,” Nature Reviews Physics, vol. 3, no. 6, pp. 422–440, 2021, doi: 10.1038/s42254-021-00314-5.
[28]
S. M. Lundberg and S.-I. Lee, “A unified approach to interpreting model predictions,” Advances in Neural Information Processing Systems, vol. 30, 2017.
[29]
R. R. Selvaraju, M. Cogswell, A. Das, R. Vedantam, D. Parikh, and D. Batra, “Grad-CAM: Visual explanations from deep networks via gradient-based localization,” Proceedings of the IEEE International Conference on Computer Vision, pp. 618–626, 2017, doi: 10.1109/ICCV.2017.74.
[30]
N. Sukulthanasorn and K. Terada, “Quantum neural network-assisted topology optimization: Concept and implementation with parameterized quantum circuits,” Computer Methods in Applied Mechanics and Engineering, vol. 448, p. 118411, 2026, doi: 10.1016/j.cma.2025.118411.
[31]
A. Tabarraei, “Variational quantum latent encoding for topology optimization,” Engineering with Computers, vol. 41, pp. 4549–4573, 2025, doi: 10.1007/s00366-025-02208-x.
[32]
M. Cerezo, G. Verdon, H.-Y. Huang, L. Cincio, and P. J. Coles, “Challenges and opportunities in quantum machine learning,” Nature Computational Science, vol. 2, no. 9, pp. 567–576, 2022, doi: 10.1038/s43588-022-00311-3.