Multi-Objective Design Optimization of a Shape Memory Alloy Flexural Actuator ^†

Abstract

This paper presents a computational model and design optimization strategy for shape memory alloy (SMA) flexural actuators. These actuators consist of curved SMA wires embedded within elastic structures; one potential application is positioning microcatheters inside blood vessels during clinical treatments. Each SMA wire is shape-set to an initial curvature and inserted along the neutral axis of a straight elastic member (cast polydimethylsiloxane, PDMS). The elastic structure preloads the SMA, reducing the equilibrium curvature of the composite actuator. Temperature-induced phase transformations in the SMA are achieved via Joule heating, enabling strain recovery and increased bending (increased curvature) in the actuator. Actuator behavior is modeled using the homogenized energy framework, and the effects of two critical design parameters (initial SMA curvature and flexural rigidity of the elastic sleeve) on activation curvature are investigated. Finally, a multi-objective genetic algorithm is utilized to optimize actuator performance and generate a Pareto frontier, which is subsequently experimentally validated.

1. Introduction

Liver cancer is the second leading cause of cancer-related death worldwide [1]. Surgical resection is the preferred, and most effective, treatment approach for primary and metastatic liver malignancies. Despite the large number of patients affected by these diseases, current treatment options are less than ideal, as 80–90% of patients are not good candidates for surgery [2], 50% of patients with solid tumors have full resistance to chemotherapy [2], and liver tissue is more sensitive to external radiation than the tumor itself [3].

Selective internal radiation therapy (SIRT) is an emerging treatment for liver cancer during which the blood vessels supplying tumors are embolized with radioactive microspheres. These microspheres contain the radioactive isotope yttrium-90 (⁹⁰Y) which delivers high-energy, low-penetrating radiation that destroys tumor tissue while limiting adverse effects to surrounding healthy liver tissue [4]. Currently, SIRT uses a single-lumen microcatheter (1.0 mm diameter) and a manually operated syringe to deliver the ⁹⁰Y microspheres into the hepatic artery. Recent studies have documented significant increases in patient survival (29.4 months with SIRT and systemic chemotherapy vs. 12.8 months with chemotherapy alone [5]). Despite the documented advantages of SIRT, the inability to directly target tumor sites has limited its widespread adoption. The tortuous and patient-variable arterial anatomy limits the direct targeting of tumors, causing damage to healthy liver tissue and other organs via ischemia or radiation [6].

Computational fluid dynamic (CFD) simulations have shown that precise tumor targeting can be achieved by controlling the catheter tip location within the cross-section of the common hepatic artery (CHA) (mean diameter 6.3 ± 0.9 mm) [7], as shown in Figure 1. This finding has been validated experimentally by Richards et al. [6], who concluded that the development of an actuated microcatheter is necessary for direct tumor targeting to be realized. Therefore, to position a microcatheter anywhere within the cross section of the CHA, a miniature actuator needs to be developed that can produce a stroke of at least 7 mm with minimal disruption to the blood flow around the catheter.

Because shape memory alloys (SMAs) exhibit high energy densities, they are attractive options for microscale actuation [8]. One actuator option consists of shape-set SMA tendons embedded within an elastic structure to generate large bending moments. The elastic structure of this “flexural actuator” produces a stress-induced martensitic phase transformation in the SMA, which can be recovered by Joule heating. While flexural actuators can also be designed by offsetting the SMA from the neutral axis, this offset increases the diameter of the actuator, making it less suitable for size-critical applications. One proposed active catheter concept, shown in Figure 2, uses eight shape-set flexural actuators arranged radially. With this configuration, the microcatheter’s release point within the cross section of the CHA could be precisely regulated while keeping the catheter parallel with the surrounding blood flow.

In this paper, we build on our previous research [9] to present a design optimization strategy for shape-set flexural actuators based on the homogenized energy framework [10,11]. Previous research on flexural actuators has focused on modeling and understanding the influences of flexible beam composition. Lagoudas, et al. successfully modeled a flexible rod with a single embedded off-neutral axis SMA fiber [12]. Wang, et al. presented numerical simulations and an experimental validation of a similar actuator composed of a flexible inner cylindrical beam, a SMA tendon, and an outer elastic cylindrical beam [12]. These designs, however, are difficult to miniaturize because the offset SMA increases the diameter of the actuator and introduces complexity into the manufacturing process.

Baz, et al. modeled shape-set SMA strips embedded on the neutral plane in a fiberglass beam [13]. The model showed good correlation with experimental data, proving that shape-setting SMA is a viable way to control the shape of a flexural actuator. However, the model did not take into account stress non-linearities across the cross-section of the SMA. Ryu, et al. [14] investigated the relationship between composite beam stiffness and bending performance.

The remainder of this paper is organized as follows. The flexural actuator system model is presented first. Next, the experimental setup is presented. The SMA model parameters are then estimated using a single-objective genetic algorithm. A multi-objective genetic algorithm is utilized to optimize actuator design parameters. Finally, modeled and experimental data are compared.

2. System Model

Each flexural actuator design consists of two components: a curved SMA tendon embedded within a straight elastic sleeve. The SMA tendon (diameter a, straight length L, and connection length x as in Figure 3a) is shape-set (or thermally trained) to a particular curvature κ₀, where

$${\kappa }_0=\frac{1}{R}$$

(1)

and R is the radius of the curvature, shown in Figure 3a. The SMA tendon is embedded within a rectangular elastic sleeve (height h, width w, and length L, as in Figure 3b) which produces a strain-induced phase transition to detwinned martensite, resulting in a reduced deactivated equilibrium curvature κ_d (Figure 3c). As the tendon is heated to temperatures above the austenitic finish temperature A_f, the material undergoes a phase transformation to austenite and recovers some of this strain, establishing a larger activated equilibrium curvature κ_a, such that κ_d < κ_a < κ₀.

To accurately model the deactivated and activated displacements, the equilibrium curvature must be found for each SMA state. Curvature is related to tip displacement by

$$y=\frac{1}{\kappa }-\frac{1}{\kappa }\cos \left(\kappa L\right)$$

(2)

To enable precise catheter navigation and positioning, the deactivated displacement y_d must be minimized while the activated displacement y_a must be maximized (Figure 4).

2.1. Homogenized Energy Model of SMA Flexural Actuators

The homogenized energy model (HEM) is used to describe the nonlinear hysteretic relationships between stress, strain, and temperature [10]. The HEM is a macroscopic model that incorporates material inhomogeneities and interaction effects. In [10], the macroscopic strain is modeled as a function of applied stress. Here, we invert the relationship to obtain the macroscopic stress

$$\sigma (\epsilon ,T)={\displaystyle {\int }_{-\infty }^{\infty }{\displaystyle {\int }_0^{\infty }{\nu }_R}}({\sigma }_R){\nu }_I({\sigma }_I)\overline{\sigma }(\epsilon ,T)d{\sigma }_{R}d{\sigma }_I,$$

(3)

where the local stress

$$\overline{\sigma }(\epsilon ,T)=\frac{\epsilon -{\epsilon }_T\left(x_{M+}-x_{M-}\right)}{\frac{x_A}{E_A}+\frac{x_{M+}+x_{M-}}{E_M}}$$

(4)

depends on the SMA phase fractions: austenite x_A, martensite plus x_M+, and martensite minus x_M-. The relative stress σ_R and interaction stress σ_I are assumed to be manifestations of underlying densities that affect the equilibrium phase fractions. In Equation (4), E_A is the austenitic elastic modulus, E_M is the martensitic elastic modulus, and ε_T is the maximum recoverable strain.

2.2. Equilibrium Equations

At equilibrium, the composite actuator’s internal forces and moments must sum to zero. Neglecting gravity and external forces, the actuator’s equilibrium position is solely determined by the moment balance between the SMA tendon and the elastic sleeve (Figure 5a). The elastic moment is given by

$$M_e\left(\kappa \right)=EI\kappa ,$$

(5)

where E is the sleeve’s elastic (Young’s) modulus, I is its area moment of inertia, and κ is its curvature. The SMA moment is given by

$$M_{sma}\left(\kappa ,T\right)={\displaystyle \int {\displaystyle \int y}}\sigma (\epsilon ,T)dA,$$

(6)

where dA is the differential area of the tendon cross-section, as shown in Figure 5b. Note that the stress depends non-linearly on strain ε and temperature T.

We assume that strain varies linearly over the cross-sectional area according to

$$\epsilon =\left({\kappa }_0-\kappa \right)y.$$

(7)

Therefore, at the initial (shape-set) curvature κ₀, the SMA has zero strain over its cross-sectional area.The moment integral (6) is discretized using rectangular elements [15], as shown in Figure 6. Assuming the stress distribution is symmetric about the neutral axis, this discretization yields

$$M_{sma}=2{\displaystyle \sum _{i=1}^{N_i}2}\sqrt{a^2-y_i^2}\sigma \left(\epsilon ,T\right)y_i\Delta y,$$

(8)

where a is the SMA tendon radius, y_i is the midpoint of the $i^{th}$ rectangle, and $2\sqrt{a^2-y_i^2}$ is the width of the rectangle.

The equilibrium deactivated curvature minimizes

$${\kappa }_d=\underset{\kappa }{\arg \min }{\left(M_e(\kappa )-M_{sma}(\kappa ,T_{\infty })\right)}^2,$$

(9)

where T_∞ is the ambient temperature. Equation (9) can be solved using a variety of nonlinear optimization algorithms; here, we employ the golden section search method [16]. After finding the deactivated curvature (and associated equilibrium phase fractions), the equilibrium activated curvature minimizes

$${\kappa }_a=\underset{\kappa }{\arg \min }{\left(M_e(\kappa )-M_{sma}(\kappa ,T_{ss})\right)}^2,$$

(10)

where T_ss is the SMA’s activated steady-state temperature.

To determine the ideal number of discretization layers, simulations were conducted for a varying initial curvature κ₀ and sleeve flexural rigidity EI. The nonlinear relationship between the initial and final curvature in Figure 7 is due to the nonlinear stress–strain characteristics of SMA, which are included in our HEM model. At very low initial strains, the SMA is in its austenite phase, and the stress–strain relationship is relatively linear. This characteristic is evident in Figure 7a for initial curvatures up to approximately 100 m⁻¹. At larger strains, the austenite to martensite transformation causes a plateauing of the stress–strain relationship. This characteristic is evident in Figure 7a for initial curvatures between approximately 100 and 700 m⁻¹. When the phase transformation is complete (to fully detwinned martensite), a linear stress–strain relationship resumes. This characteristic is evident in Figure 7a for initial curvatures above approximately 700 m⁻¹. As shown in Figure 7, the accuracy of computed actuator curvature improves with the number of discretization layers (N). However, because computational time is directly proportional to N, there exists a tradeoff between model accuracy and computational time. Based on the simulation results, N = 8 was found to adequately balance this tradeoff, as it produced curvature estimates within 2% of the full-order model (N = 128) with only 6.25% of the computational burden.

3. Experimental Setup

To constrain the SMA tendons during the shape setting process, a custom fixture (Figure 8) was designed and fabricated from a 304 stainless steel plate (0.305 m × 0.102 m × 0.013 m). This fixture features eight precisely machined curvature profiles ranging from 50 to 333 (m⁻¹), with bolts and washers used to secure the SMA tendon during shape setting (Figure 9). This range of initial curvatures was determined from simulations to provide the largest difference in the activated and deactivated curvature (Figure 10b). To standardize the electrical resistance of each 0.31 mm diameter SMA tendon (Dynalloy, Inc., Irvine, CA, USA), each curvature profile was machined to provide the same arc length and straight length. Constant tendon tension was maintained by suspending a 400 g mass vertically from each SMA tendon during the fixture tightening process.

To shape-set the SMA tendons, a computer-controlled furnace (Thermolyne 1500) was preheated to 500 °C, and the loaded fixture was placed inside for 20 min [17]. Immediately following heating, the fixture was quenched in room-temperature water. After cooling, each SMA tendon was removed from the fixture, resulting in stress-induced deformations. To experimentally quantify each tendon’s initial curvature κ₀, a programmable power supply (Agilent E3615A) was used to activate the specimen (maintaining 1.4 amps of DC current for 5.0 s, resulting in thermally-induced austenitic phase transformation). The curvature of the actuated specimen was measured by comparison to a printed scale.

To fabricate elastic sleeves with specific flexural rigidities (EI) in the bending plane, molds were designed with rectangular cross-sections (Figure 11). By specifying the cross-sectional height to be three times the width, the sleeve is 27 times more compliant in the bending plane than out-of-plane; the effects of self-weight are thus minimized during experimental testing. The range of flexural rigidities (6 × 10⁻⁶ to 3 × 10⁻⁴ Nm²) was determined from simulations to provide the largest change in activated and deactivated curvature (Figure 10). The molds were fabricated using rapid prototyping in ABSplus material (Dimension Elite, Stratasys, Eden Prairie, MN, USA). Nylon monofilament (0.31 mm diameter) was tensioned along the neutral axis of each rectangular mold to provide space for SMA tendon insertion. Each mold was cast with polydimethylsiloxane (PDMS, Dow Corning SYLGARD 184) and allowed to cure for two days.

Each actuator was assembled by inserting a shape-set SMA tendon into a molded elastic sleeve. The experimental set up consisted of a cantilevered actuator connected to a programmable power supply (E3615A, Agilent Technologies) by 36 AWG magnet wire. Actuation was achieved by passing 1.4 A of DC current through the SMA tendon for 30 s. The actuator was then deactivated for one minute, cooling via heat transfer to the ambient air. High-resolution photographs (Canon DS126191 EOS Rebel 10.1MP) were taken immediately before and after each actuation cycle. Digital images were acquired and stored for subsequent evaluation of κ_a and κ_d. The process of activating, deactivating and image acquisition was automated using a custom LabView program (National Instruments, Austin, TX, USA).

To compute κ_a and κ_d, an image processing code was developed using MATLAB’s Image Processing toolbox (Mathworks, Inc., Natick, MA, USA). The code applied golden-section optimization to find the bending radius R (Figure 12) that minimizes the difference between predicted (L_pred) and actual bending length (L_act) according to

$$L_{pred}=\theta R,$$

(11)

$$L_{act}=R{\cos }^{-1}\left(\frac{2R^2-C^2}{2R^2}\right),$$

(12)

where C is the corresponding chord length and θ is the included angle between the linear sections of the sleeve (Figure 12). To evaluate the curvature along the actuator’s neutral axis, half of the known width of the elastic sleeve was added to the computed inner radius (R). This computation was repeated for 110 consecutive actuations, averaging the last 10 cycles to eliminate variance and drift associated with low cycle fatigue [18].

Design optimization was conducted using multi-objective genetic algorithm (MOGA) techniques. Genetic algorithms (GAs) are optimization methods based on the principles of natural selection. These algorithms minimize (or maximize) one or more objective functions and are reliable and efficient in optimizing non-linear, multi-modal and complex systems. An excellent overview of GAs can be found in [19,20].

4. Results

4.1. HEM Parameter Optimization

A single-objective genetic algorithm was utilized to minimize

$$J={\displaystyle \sum _i{\left({\kappa }_{a_i}-{\widehat{\kappa }}_{a_i}\right)}^2},$$

(13)

the sum of squared errors (SSE) between predicted $\left({\widehat{\kappa }}_{a_i}\right)$ and experimentally measured $\left({\kappa }_{a_i}\right)$ actuator curvatures. The 11 HEM stress density coefficients (α_k,m and β_k) and seven SMA constitutive model parameters (E_A, E_M, σ_L, Δσ_T, ε_T, h and T_max) were chosen as design variables. The design variables and their associated bounds are listed in

Table 1. Single-objective design variables and bounds.

Variable	Description	Lower Bound	Upper Bound	Units
${\alpha }_{k,m}$	Relative stress density coefficients	0	3.0	-
${\beta }_k$	Interaction stress density coefficients	0	3.0	-
$E_A$	Austenitic elastic modulus	10	100	GPa
$E_M$	Martensitic elastic modulus	10	100	GPa
${\sigma }_L$	Martensitic transition stress at 348 K	100	400	MPa
$\Delta {\sigma }_T$	Hysteresis loop’s temperature dependence	1.0	15.0	MPa/K
${\epsilon }_T$	Maximum recoverable strain	1.0	7.0	%
h	Convection Coefficient	0.01	1.6	W/(m2K)
Tmax	Maximum SMA Temperature	353	403	K

.

Experimental data was obtained for actuators with a fixed flexural rigidity (EI of 3 × 10⁻⁵ Nm²) and six different initial curvatures (κ₀ = 100, 125, 167, 200, 250, and 333 m⁻¹). Additional data was obtained for a fixed initial curvature (κ₀ of 200 m⁻¹) and four different values of flexural rigidity (EI = 6 × 10⁻⁶, 1.5 × 10⁻⁵, 3 × 10⁻⁵, 6 × 10⁻⁵ Nm²). The GA function in the MATLAB optimization toolbox was employed for single-objective design optimization. The design variables were real value encoded and 36 individuals were systematically populated within initial design variable bounds. The bounds were chosen from knowledge of the design space to speed convergence. Heuristic crossover was utilized with a crossover fraction of 1.2. The convergence criteria consisted of a minimum objective function gradient (1 × 10⁻⁶) and a generation limit (200). The initial and optimized model parameters are shown in

Table 2. Comparison of initial and optimized homogenized energy model (HEM) parameters.

Variable	Initial Value	Optimized Value	Units
$E_A$	30.7	59.02	GPa
$E_M$	26.0	16.83	GPa
${\sigma }_L$	295	204.17	MPa
$\Delta {\sigma }_T$	9.2	9.96	MPa
${\epsilon }_T$	4.4	7	%
h	0.8	0.4279	W/(m2K)
Tmax	380	381.4	K

and in Figure 13. The initial parameters were taken from [10].

GA convergence was monitored for five simulations with different initialization parameters; all five simulations converged based on the minimum function gradient criteria and the simulation average SSE decreased from 7263 to 514.2 with a standard deviation of 7.1. Though the designs did not converge to the same values in the design space, all designs converged to approximately the same SSE. The optimized model parameter set that reduced SSE by 93% (

Table 3. Sum of squared errors (SSE) for initial (non-optimized) and optimized models.

	SSE Deactivated	SSE Activated	SSE Sum
Initial Model	4792	2483	7263
Optimized Model	478.1	36.7	514.2
Percent Improvement	90%	99%	93%

) was selected for subsequent multi-objective optimization. Comparisons of the experimental data and initial (unoptimized) parameters are presented in Figure 14; the dependence of the actuator curvature on sleeve flexural rigidity is shown in Figure 14a, and its dependence on initial curvature is shown in Figure 14b. Differences between the two are primarily attributed to material parameter uncertainties (modulus of elasticity, etc.) and unmodeled dynamics (sliding friction between the SMA and sleeve, heat transfer modes, etc.). Figure 15 shows similar comparisons for the optimized model. The optimized HEM accurately characterizes the nonlinear dependence of the equilibrium curvature on sleeve flexural rigidity (Figure 15a) and the dependence of the actuator curvature on initial curvature (Figure 15b).

Comparing the optimized HEM parameters to the imposed parameter bounds reveals that the maximum recoverable strain lies on the upper bound. This indicates that the equilibrium model may need to take into account additional variables to more accurately predict equilibrium curvature: internal friction, non-linear elastic modulus of PDMS at high strains or the effect contact has on temperature distribution across the tendon. Despite these modeling limitations, the optimized model can still be used for design optimization.

4.2. Multi-Objective Optimization

To optimize the performance of the actuator, the deactivated displacement needs to be minimized while the activated displacement needs to be maximized. To accomplish this, a MOGA was utilized with κ₀ and EI as the design variables. The optimized HEM and equilibrium curvature model are used to evaluate κ_a and κ_d, which are related to y_a and y_d via Equation (2). The actuator length (L = 0.011 m) was chosen to achieve 0.007 m displacement; the actuator is $\theta ={90}^{\circ }$(Figure 12). To improve numerical conditioning, log₁₀ (EI) was used to normalize the individuals for reproduction. The gamultiobj function in MATLAB was employed for the MOGA. Most of the GA settings were identical to the single-objective optimization case, with the exceptions of population (50) and minimum objective function gradient (1 × 10⁻⁴). Upper and lower parameter bounds were utilized for the population initialization (UB_i and LB_i, respectively), and values for subsequent populations (UB and LB) are shown in

Table 4. Upper and lower parameter bounds for the initial population (UB_i and LB_i, respectively) and subsequent populations (UB and LB).

	${\mathit{\kappa }}_0$	EI
UBi	600	1 × 10−5
LBi	0	1 × 10−2
UB	600	1
LB	0	0

. The MOGA converged based on the minimum objective function gradient.

The Pareto frontier of Figure 16 shows the non-dominated designs from the optimization, and clearly reveals the inherent design tradeoff between y_d and y_a: no design combination of EI and κ₀ perfectly minimizes y_d (y_d = 0 m) while simultaneously maximizing y_a (y_a ≥ 0.007 m). The results of Figure 17a indicate that the Pareto-optimal displacements y_d and y_a have strong non-linear dependencies on EI (as was the case with κ_a and κ_d in Figure 15a). Furthermore, the correlation between y_d, y_a and κ₀ in Figure 17b resembles the quasi-linear dependence of κ_a and κ_d on κ₀ in Figure 15b. It is also clear from Figure 17 that the change in activated displacement occurs predominantly over initial curvatures ranging from 50 to 375 (m⁻¹) and flexural rigidities, EI, ranging from 1 × 10⁻⁵ to 5 × 10⁻⁴ Nm².

The interaction of Pareto-optimal EI and κ₀ in Figure 18 illustrates how these design variables relate to actuation range: generally, this range is proportional to the initial curvature and inversely proportional to flexural rigidity. An additional performance metric, the actuation ratio (y_a/y_d), helps illustrate the design tradeoff between y_a and y_d, as shown in Figure 19. Arguably, the “best” design balances both performance metrics, such as the design (y_a = 0.004 m, y_d = 0.0019 m) indicated in Figure 19.

4.3. Experimental Validation of Pareto Frontier

To experimentally validate the Pareto frontier, six Pareto-optimal designs (

Table 5. Selected Pareto-optimal designs for experimental validation.

Design	Initial Curvature (m−1)	Sleeve EI (Nm2)	Deactivated Displacement (m)	Activated Displacement (m)
1	100	3.05 × 10−5	0.0010	0.0025
2	125	2.23 × 10−5	0.0016	0.0036
3	167	1.76 × 10−5	0.0025	0.0049
4	200	1.58 × 10−5	0.0031	0.0057
5	250	1.47 × 10−5	0.0037	0.0065
6	333	1.16 × 10−5	0.0049	0.0075

) were selected, fabricated, experimentally evaluated, and compared to the simulated Pareto frontier. The results of this comparison are displayed in Figure 20.

The error between the modeled and experimental actuator displacements was quantified via Euclidean norms and normalized with respect to maximum actuator displacement (0.008 m). The average error for all six selected Pareto-optimal designs was 6.4%. The average activated displacement error (2.5%) was significantly lower than average deactivated displacement error (5.6%). Figure 15b provides a likely explanation for this effect: the optimized HEM underestimates ${\kappa }_d$ at low ${\kappa }_0$ and overestimates ${\kappa }_d$ at high ${\kappa }_0$. This also explains why certain experimental evaluations of representative designs out-perform the Pareto frontier: because ${\kappa }_d\propto y_d$, the experimental evaluations of y_d at low κ₀ (100, 125 m⁻¹) were higher than predicted and high κ₀ (200, 250, 333 m⁻¹) were less than predicted.

5. Conclusions

This paper demonstrates a method for optimizing the design of shape-set SMA flexural actuators. Optimized HEM and bending equilibrium models reveal the non-linear dependencies of the actuator bending curvature on flexural rigidity and initial curvature. Experimental results validate the accuracy of these models as predictors of activated and deactivated curvature. A Pareto frontier based on these models reveals inherent design tradeoffs between activated and deactivated displacement. While the goal here is optimizing the performance of an SMA flexural actuator, the methods and results can easily be extended to other smart materials and actuator topologies. Future work will focus on weighting functions and “strength of preference” curves to guide the selection of the “best” design. The model will also be refined to include internal friction, non-linear sleeve elastic modulus and non-uniform temperature distribution.