## 2 Materials and Methods  
### Fiber-Gel Finite Element Model  
The fiber-gel model combines a viscoelastic ground substance (a gel) with directional fibers under tension in all layers of the vocal fold (superficial layer of the lamina propria (SLLP), ligament, and muscle) [42]. At the physiologic level of input, the control parameters are lung pressure and activation levels of five intrinsic laryngeal muscle pairs (if left-right symmetry is assumed), CT, TA, LCA, PCA, and IA. With these inputs, and stored anatomical dimensions of the larynx ([43], Chapter. 1) and the airways [44], pre-phonatory vocal fold and vocal tract shapes are produced. Figure 1 shows how all five intrinsic muscles are contributing to the vocal fold length, \(L\) and the position (\(\xi_{02}\), \(\psi_{02}\)) of the vocal process of the arytenoid cartilage in the horizontal plane. The plane is at the level of the vocal ligament. All muscles except the CT apply forces directly to the arytenoid cartilage. The CT has a forward and downward pull on the thyroid cartilage relative to the crioid cartilage, therewith elongating the vocal folds. The dotted-line vector for LCA indicates that the fibers run below those of the TA, while the dotted line for the PCA vector indicate an effective force direction due to curvilinear fiber attachment. The combined mechanics of the cartilages is solved with the second-order Runge-Kutta differential method [43]. Only right-side cartilages and forces were shown in Figure 1 for clarity. The mechanics of the cartilages on the left side are solved similarly. In the rostral-cranial (vertical) direction, the medial surface of the vocal folds is defined with three parameters, the vocal fold thickness \(T\), the lower adduction of the arytenoid cartilage \(\xi_{01}\), and a medial surface bulging factor \(\xi_{b}\). With these parameters, the glottal half-width at any point (\(y,z\)) along the surface of each vocal fold is defined as  
\[\xi_{0}(y,z)\ =\ (1-y/L)[\xi_{02}+(\xi_{01}-\xi_{02}-4\xi_{b}z/T)(1-z/T)] \tag{1}\]  
The parameters in these equations are determined by muscle activations with empirical rules  
\[L\ =\ L_{0}(1+\varepsilon) \tag{2}\]  
\[T\ =\ T_{0}/(1+\varepsilon) \tag{3}\]  
\[\xi_{01}\ =\ \xi_{02}+\max(0,0.1T(1-1.5a_{TA})) \tag{4}\]  
\[\xi_{b}\ =\ 0.05T(1-a_{TA}) \tag{5}\]  
Figure 1: Forces acting on arytenoid cartilage for cricoartytenoid movements ([PERSON], 2006 [43]). where \(\varepsilon\) is the vocal fold strain (determined by all muscles) and \(a_{TA}\) is the activation of the TA muscle. Here, \(L_{0}\) was set to be 1.6 cm, and \(T_{0}\) was set to be 0.7 cm for an average male vocal fold. Figure 2a shows the right vocal fold geometry as parameterized by the empirical rules given above. Also shown is the variable \(\xi_{m}\), the amplitude of vibration when self-sustained oscillation takes place. The left vocal fold is also governed by the same parameters. The geometric properties of the layers are modeled according to the histological data available from the literature [45; 46; 47]. The full details of vocal fold posturing with intrinsic laryngeal muscle activations is provided in ([43], Chapter 3). As tissue is deformed with muscle contraction (e.g., length, thickness, and medial surface of the vocal folds), viscoelastic parameters are modified with basic constitutive equations. The passive fibrous tissue properties in the y-direction in the three vocal fold layers and the other intrinsic laryngeal muscles are characterized with a combined linear and exponential stress function ([43], p. 77):  
\[\sigma_{y}\ =\ -\sigma_{0}(\varepsilon-\varepsilon_{1})+\sigma_{2}\Big{[}e^{B( \varepsilon-\varepsilon_{2})}-1-B(\varepsilon-\varepsilon_{2})\Big{]} \tag{6}\]  
In the above equation, \(\sigma_{0}\) and \(\sigma_{2}\) are scale factors, \(\varepsilon_{1}\) is the strain where the linear function begins and \(\varepsilon_{2}\) is the strain where the exponential term begins.