## 3 Discrete Formulation of Topology Optimization Problem Statements and Sensitivities
### Stress Functional
\begin{split}\frac{\partial\sigma^{PN}(\rho)}{\partial\rho}= \frac{\partial\sigma^{PN}(\rho)}{\partial\sigma_{a}^{pM}}\frac{\partial \sigma_{a}^{pM}(\rho)}{\partial\sigma_{a}}\frac{\partial\sigma_{a}(\rho)}{ \partial\rho}=\\ =\frac{\partial\sigma^{PN}(\rho)}{\partial\sigma_{a}^{pM}}\frac{ \partial\sigma_{a}^{pM}(\rho)}{\partial\sigma_{a}}\frac{\left(1\over 2\right.} \rho^{-\frac{1}{2}}\mathbb{C}_{0}B_{a}u(\rho)\\ -\rho^{\frac{1}{2}}\mathbb{C}_{0}B_{a}K^{-1}(\rho)\frac{\partial K (\rho)}{\partial\rho}u(\rho)\end{split}  \tag{53}\]  
Renaming these two terms as \(\partial\sigma_{spf}^{PN}(\rho)\) (Equation (54)), referring to the component of the derivative that is dependent on the stress penalization factor, and \(\partial\sigma_{u}^{PN}(\rho)\) (Equation (55)), referring to the component of the derivative that is dependent on the nodal displacement:  
\[\partial\sigma_{spf}^{PN}(\rho)=\frac{\partial\sigma^{PN}(\rho)}{\partial \sigma_{a}^{pM}}\frac{\partial\sigma_{a}^{pM}(\rho)}{\partial\sigma_{a}} \left(\frac{1}{2}\rho^{-\frac{1}{2}}\mathbb{C}_{0}B_{a}u(\rho)\right) \tag{54}\]  
\[\partial\sigma_{u}^{PN}(\rho)=\frac{\partial\sigma^{PN}(\rho)}{\partial\sigma _{a}^{pM}}\frac{\partial\sigma_{a}^{pM}(\rho)}{\partial\sigma_{a}}\times \left(-\rho^{\frac{1}{2}}\mathbb{C}_{0}B_{a}K^{-1}(\rho)\frac{\partial K(\rho) }{\partial\rho}u(\rho)\right). \tag{55}\]  
Notice that \(\partial\sigma_{spf}^{PN}(\rho)\) can be easily determined, since the only information required is the stress vector \((\sigma_{a}(\rho)=\mathbb{C}_{0}B_{a}u(\rho))\) and the design densities. \(\frac{\partial\sigma^{PN}(\rho)}{\partial\sigma_{a}^{pM}}\) can be determined with a simple summation and \(\frac{\partial\sigma_{a}^{pM}(\rho)}{\partial\sigma_{a}}\) through the product of two vectors and a matrix. The term \(\partial\sigma_{u}^{PN}(\rho)\), on the other hand, requires the explicit definition of the matrices \(B_{a}\), \(K^{-1}(\rho)\), and \(\frac{\partial K(\rho)}{\partial\rho}\), which are dependent on the element formulation used in the numerical model. Furthermore, for topology optimization problems where the number of design variables is larger than the number of constraints, the most-efficient way to determine \(\partial\sigma_{u}^{PN}(\rho)\) is using an adjoint model, defining the adjoint variable as:  
\[ \begin{split} K(\rho)\lambda&=\frac{\partial \sigma^{PN}(\rho)}{\partial\sigma_{a}^{MM}}\frac{\partial\sigma_{a}^{\
u M}( \rho)}{\partial\sigma_{a}}\mathbb{C}_{0}B_{a}\Leftrightarrow\\ \Leftrightarrow\lambda&=\frac{\partial\sigma^{PN}( \rho)}{\partial\sigma_{a}^{\
u M}}\frac{\partial\sigma_{a}^{\
u M}(\rho)}{ \partial\sigma_{a}}\mathbb{C}_{0}B_{a}K^{-1}(\rho).\end{split}  \tag{56}\]  
Therefore, the adjoint variable can be extracted from a finite element analysis, where the load applied on each node is equal to \(\frac{\partial\sigma^{PN}(\rho)}{\partial\sigma_{a}^{\
u M}}\frac{\partial \sigma_{a}^{\
u M}(\rho)}{\partial\sigma_{a}}\mathbb{C}_{0}B_{a}\), and \(\lambda\) is equal to the node displacement. Introducing the adjoint variable in Equation (54) then leads to:  
\[\partial\sigma_{u}^{PN}(\rho)=\frac{\partial\sigma^{PN}(\rho)}{\partial\sigma _{a}^{\
u M}}\frac{\partial\sigma_{a}^{\
u M}(\rho)}{\partial\sigma_{a}}\bigg{(} -\rho^{\frac{1}{2}}\lambda\frac{\partial K(\rho)}{\partial\rho}u(\rho)\bigg{)}. \tag{57}\]  
Finally, since \(K\) can be defined as \(K=AB_{a}^{T}\mathbb{C}(\rho)B_{a}\), since \(B_{a}\lambda\) and \(B_{a}u\) are equal to the deformation vectors of the adjoint and state models (\(\xi_{adj}\) and \(\xi_{s}\), respectively), Equation (57) can be simplified to:  
\[\partial\sigma_{u}^{PN}(\rho)=\frac{\partial\sigma^{PN}}{\partial\sigma_{a}^{ \
u M}}\frac{\partial\sigma_{a}^{\
u M}}{\partial\sigma_{a}}\Big{(}-\rho^{ \frac{1}{2}}\xi_{adj}^{T}(\rho)\mathbb{C}^{\prime}(\rho)\xi_{s}(\rho)\Big{)}. \tag{58}\]  
With this information, the first three steps of the generic optimization process can be rewritten as follows:  
* Find \(u\) such that: \(K(\rho)u=F(\rho)\).