T-splines computational membrane–cable structural mechanics with continuity and smoothness: I. Method and implementation

We present a T-splines computational method and its implementation where structures with different parametric dimensions are connected with continuity and smoothness. We derive the basis functions in the context of connecting structures with 2D and 1D parametric dimensions. Derivation of the basis functions with a desired smoothness involves proper selection of a scale factor for the knot vector of the 1D structure and results in new control-point locations. While the method description focuses on C0\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$C^0$$\end{document} and C1\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$C^1$$\end{document} continuity, paths to higher-order continuity are marked where needed. In presenting the method and its implementation, we refer to the 2D structure as “membrane” and the 1D structure as “cable.” It goes without saying that the method and its implementation are applicable also to other 2D–1D cases, such as shell–cable and shell–beam structures. We present test computations not only for membrane–cable structures but also for shell–cable structures. The computations demonstrate how the method performs.


Introduction
The isogeometric analysis (IGA), with the superior accuracy it offers, brought fluid and solid mechanics computations to a new level [1][2][3][4]. Being able to use the IGA basis functions also in time in the context of space-time (ST) computational analysis expended the scope of the IGA and led to the introduction of the ST-IGA [5][6][7]. The terminology "ST-IGA" implies, depending on the context, discretization with IGA basis functions in space or time or both. In the 2D test cases reported in [5], the computation of flow past an airfoil was with IGA basis functions in space, and the advection com-putations with IGA basis functions in both space and time. The advection computations, accompanied by a stability and accuracy analysis for the pure equation, showed the advantages of using higher-order basis functions, not only in space, but also in time. Related to that, keeping in mind that the increased accuracy the ST-IGA with IGA basis functions in space brings is attained with fewer control points, the effective element sizes will be larger. With that, larger time steps can be taken while still keeping the Courant number at or below the levels we target for good accuracy.
Moving to solid and structural mechanics computations with IGA basis functions in space, it was pointed out as early as in 2007 (see [61]) that the image-based geometries used in patient-specific arterial FSI computations are not for the zero-stress state (ZSS) of the artery and that a ZSS estimation method is needed. The ZSS estimation methods introduced in and after 2016 [51,[62][63][64][65] stand on the IGA basis functions in space, and so does the related hyperelastic shell analysis [66,67]. The IGA basis functions in space have been a part of quite a few advanced computational methods targeting design and structural analysis, those reported in [68][69][70][71][72][73][74][75][76][77] are examples of that, and turbine blades and heart valves are among the examples.
In IGA discretization, specifying Dirichlet boundary conditions could be challenging. That is because the basis functions are generally not interpolatory. In the case of differential equations allowing Dirichlet conditions also on the derivatives, specifying conditions on the derivatives could be easier in IGA discretization than it is in finite element discretization. For example, specifying the slope can be accomplished by constraining the motion of the nearest interior point to the horizontal line passing through the boundary point.
Our objective in the work presented here is to address a related challenge. That is, computational structural analysis where structures with different parametric dimensions are connected with continuity and smoothness. In structural analysis of a membrane and cable, for example, connecting the two with C 0 or C 1 continuity involves challenges similar to those associated with specifying Dirichlet boundary conditions in IGA discretization. Even if the two structures have the same parametric dimensions, connecting them with C 1 continuity could still be challenging.
There has been some earlier related work. They include the bending-strip method for shell-shell [78] and shell-beam [79] structures, where the beam is actually a bending-stabilized cable, penalty formulation [80], and techniques based on Nitsche's method [81]. They also include using extra mesh refinement along the membrane edge [45] to attain C 0 continuity in both the edge direction and the other direction. Using T-splines [82] to attain C 1 continuity when the two structures have the same parametric dimensions is also among the earlier related work.
In our work here, we derive the basis functions that give us the desired smoothness between structures with 2D and 1D parametric dimensions. The derivation involves proper selection of a scale factor for the knot vector of the 1D structure and results in new control-point locations. While our method description will focus on C 0 and C 1 continuity, paths to higher-order continuity will be marked where needed. In presenting the method and its implementation, we will refer to the 2D structure as "membrane" and the 1D structure as "cable." The method and its implementation will, of course, be applicable also to other 2D-1D cases, such as shell-cable and shell-beam structures. When the membrane and cable are connected with smoothness, the strain and rotational freedoms are transferred between the two structures. For easy and efficient implementation of the method, we introduce the Bézier extraction row operators to be used in obtaining the basis functions.
As the limitations of the method we will be presenting, we can mention two points. The smoothness can only be achieved along the cable and its parametric-line continuation in the membrane. The starting membrane and cable meshes need to meet some additional requirements if they are to be connected with a continuity higher than C 1 .
We present test computations not only for membranecable structures but also for shell-cable structures. We use four structure models: membrane-cable, membranebending-stabilized-cable, shell-cable, and shell-bendingstabilized-cable. We use meshes with C 0 continuity and C 2 continuity. In total, we compute eight test cases.
In Sect. 2, we introduce the method for connecting a membrane and cable by T-splines. The test computations are presented in Sect. 3, and the concluding remarks are given in Sect. 4. In Appendix A, we provide the details on the mesh examples, and in Appendix B, we give the derivation of the smoothness constraints used in the method.

Connecting a membrane and cable by T-splines
In IGA, connecting a 1D structure, such as a cable, to a 2D structure, such as a membrane, is not that straightforward.
For example, in [45], in connecting the cables to the ram-air parachute, extra mesh refinement was used along the membrane edge to attain C 0 continuity in both the edge direction and the other direction. The extra refinement, because of the The membrane and cable elements after connecting them knot insertion, reduces not only the computational efficiency but also the continuity. Especially with isogeometric shells and bending-stabilized cables [79], the continuity is essential.
Here we describe methods to connect a cable to a membrane at any location along its edge, and we also describe methods to attain smoothness along the cable and its parametric-line continuation in the membrane.

Connecting the cable to the membrane
We assume that the membrane is of rectangular shape. For illustration purposes, it is made of 2 × 3 quadratic B-spline elements, and the cable consists of one quadratic B-spline element. Figure 1 shows the membrane and cable (for more details on the mesh, see Appendix A.1), and Fig. 2 shows the mesh after connecting them (for more details on the mesh, see Appendix A.2). We will explain how we obtain this mesh. We represent the local basis functions in the parametric space −1 ≤ ξ α ≤ 1, where α = 1, . . . , n pd , and n pd is the number of parametric dimensions. For a global basis function index a or b in element e, we denote the local basis functions as M e a (ξ 1 , ξ 2 ) for the membrane elements and as L e b (ξ 1 ) for the cable elements. They are expressed as and with the index mapping e, a → (a 1 (e, a), a 2 (e, a)) and e, b → b 1 (e, b), and we are dropping "(e, a)" and "(e, b)" not to crowd the notation. Here, N e,α k represent 1D functions identified by element number e, direction α, and local index k. Although it is not included, the polynomial order is assumed to be p = p e,α , and we may omit the superscripts for notational convenience. The index a α denotes the local index in α direction, and there is no unique mapping from a α to a. In fact, quite often, a α corresponds to multiple a. The symbols M e and L e represent, for the membrane and cable, the sets of functions with nonzero value in element e.

Remark 1
All the T-splines used in this article can be expressed in this product form. Therefore, Eq. (1) is applicable even after connecting the cable to the membrane. Some of the notation may not be general enough as in NURBS, but they should be straightforward. We are giving up some generality so that we do not further complicate the notation. and Remark 2 Depending on the parametrization, i) one of ξ A,α c is either −1 or 1, and ii) ξ B,1 is either −1 or 1.
We note that N A,1 a 1 (ξ A,1 c ) is either 1 or 0 at the membrane edge. We define the set of functions with value 1 at the membrane With that, the connection point as represented by the membrane becomes Similarly, at the cable end, is either 1 or 0, and, in fact, there is only one function that has value 1, and that is b 1 = c. With that, the connection point as represented by the cable becomes To connect the membrane and cable, we need c L , and this can be enforced at the basis-function level. For that, we first exclude the basis function and control point associated with c from the cable element (and actually from the mesh). Then, we add control points x a , where a ∈ M A c , to the cable element. The corresponding basis functions, indexed by b = a, are and the overbar is for distinguishing between the basis functions before and after connecting the cable to the membrane. With that, the cable end is always on the membrane edge. The number of functions in the cable element becomes where |X | represents the cardinality of a set X .

Remark 3
The control points 7, 11, and 15 are distinct because of the membrane.

Remark 4
If the membrane and cable had rotational freedoms, they would not be transferred between the two.

Representation of the connected cable with Bézier extraction operators
We now describe how connecting the cable is implemented by using Bézier extraction operators (see the notation in [83,84] where k = 0, . . . , p, and p k = k! k!( p−k)! are the binomial coefficients. The representation of the 1D basis functions we have in Eqs. (1) and (2) will then be in the form for l = 0, . . . , p. The set of coefficient matrices are the Bézier extraction operators. However, in a T-spline element, the number of unique N A,α functions may not be equal to p A,α + 1, and the order of the rows has no significance. Therefore, we represent the Bézier extraction operator as a set of row operators, denoted by C C C A,α l ∈ R 1×( p A,α +1) , with l being the index that identifies the unique function number in an element A, in direction α.
With that, the membrane local basis functions can now be expressed as and the cable local basis functions can be expressed as From that and Eq. (2), Eq. (7) can be written as The scalar term can be expressed in a matrix form. For that, we first define a column matrix with the evaluated Bézier functions as its components: Then the scalar is written as We note that this scalar depends on a 2 . With that, the Bézier extraction operator C C C B,1 b associated with Eq. (12) becomes Now, for the examples in Figs. 1 and 2, we express the Bézier extraction row operators for the cable element. The original cable element has C C C 6,1 C C C 6,1 After the connection, we have C C C 6,1 C C C 6,1 C C C 6,1

Connecting the membrane and cable with smoothness
The smoothness can play a role if the bending stiffness of the cable will influence the membrane structure. For that purpose, we describe the case with higher-order continuity in the cable direction at the connection point. In this method, the smoothness is along the cable and its parametric-line continuation in the membrane. Again, we use Fig. 1 as an example, and Fig. 6 shows the mesh after connecting the cable and membrane with C 1 continuity along the cable line (for more details on the mesh, see Appendix A.3). The smoothness is along the cable and the parametric line ξ 2 = ξ A,2 c in the membrane element. We use the knot removal technique to have the desired continuity in the function space. We will next explain the process.
In the B-spline mesh (Fig. 1), the basis functions N A,1 a 1 for a ∈ M A c are the same. In other words, a = 7, 11, and 15 are pointing to the same index a 1 . Moreover, a = 6, 10, and 14 are pointing to a common index a 1 , and a = 5, 9, and 13 are pointing to a common index a 1 . We group them as M A c1 and M A c0 . The position along the parametric line is expressed as We know that From Eqs. (24)- (27), we observe the apparent directional basis functions and control points along the parametric line to bẽ and Similar observations can be made for the adjacent element 2. For notational convenience, we map the control points of the cable as  Figure 7 shows the apparent basis functions. These basis functions can be identified by the knot vector Ξ Ξ Ξ = {0, 0, 0, 1, 2, 2, 2, 3, 3, 3}, with 0 ≤ Ξ ≤ 2 for the membrane and 2 ≤ Ξ ≤ 3 for the cable.

Remark 7
The linear transformation of the knot vector values does not change the basis functions. When we connect two patches, a linear transformation can be applied in each patch. For example, if we scale the knots for the cable by a factor 0.5, we get Ξ Ξ Ξ = {0, 0, 0, 1, 2, 2, 2, 2.5, 2.5, 2.5}. Although the parametric space would be different, this would not change the element-wise functions and would not influence the discretization. However, the choice of the scale factor influences how the cable and membrane are connected. In fact, the desired smoothness can only be achieved with the correct scale factor. We will explain this more in Remark 8, after explaining the procedure for connecting the cable and membrane. Figure 8 shows the functions after the first and second knot removals at Ξ = 2. At each knot removal, one basis function vanishes, which means that the corresponding control point is removed. It is assumed that the geometry will not change after the knot removals. Therefore, the following relationship holds in each element: where k is the element-wise index, and N k and N k are the directional basis functions before and after the knot removals. Using the Bézier extraction operators, we can rewrite Eq. (40) Because B l are linearly independent, we obtain the following relationship for l = 0, . . . , p: For notational convenience, we rearrange this as The continuity desired after the second knot removal can be expressed by the equations We note that A and B are the elements in the membrane and the cable, and those positions are calculated from Eq. (43) for the membrane and cable elements. Figure 10 shows an example of the conversions for both the membrane and cable elements.  where s is the ratio of the neighboring nonzero knot spans of the cable (ΔΞ L ) and membrane (ΔΞ M ):

Remark 8 As explained in
The derivation can be found in Appendix B. With Eq. (49), we select s. In the example of Fig. 1, that gives us s = 1.

Remark 9
To have the highest continuity in using B-splines with polynomial order p, we need to satisfy p conditions Fig. 10 An example of the context for Eq. (43) in control-point conversions for both the membrane and cable elements. The membrane (thick green curve) and the apparent control points (dark green circles). The cable (thick brown curve) and the control points (dark brown circles). The control points after the conversions are light green circles (membrane) and light brown circles (cable). If Eqs. (46) and (47) are satisfied after the conversions, the curve will have the desired continuity and the two knots can be removed without changing the geometry that are generalization of those in Eqs. (46) and (47). The conditions arẽ for k = 1, . . . , p.
Then, we remove the control points 22 and 21. The new positions of the control points 7, 11, 15, 6, 10, and 14 are obtained from Eq. (43). With the basis functions corresponding to those, we can form the basis functions in the cable element. This can be done in the same way we did in connecting the cable and membrane with C 0 continuity. The corresponding basis functions, indexed by The control points and Bézier extraction operators for the membrane elements also change. We provide the Bézier extraction operators in Appendix A.3, and Eq. (43) can be used for obtaining Figures 11, 12, 13 6,10,14,7,11,15, and 20, and we can see the first three as replacements for 22, and the second three as replacements for 21.

Remark 10
In the membrane, up to (2 p + 1)×( p − 1) elements need to be modified, (2 p+1) in the edge direction, and ( p − 1) in the other direction. In the cable, up to p elements neighboring the membrane need to be modified.

Remark 11
The control points 7, 11, and 15 are distinct because of the membrane.

Remark 12
The strain and rotational freedoms are transferred between the two structures. Even if they are not based on a shell model or bending-stabilized cable, this smoothness requirement would add some kind of bending effect to the solution, which will be smaller and smaller with mesh refinement. Fig. 17 The membrane and cable meshes before connecting them. The checkerboard pattern and alternating colors are for differentiating between the elements, and circles represent the control points

Test computations
We use four structure models: membrane-cable, membranebending-stabilized-cable, shell-cable, and shell-bendingstabilized-cable. We use meshes with C 0 and C 2 continuity between the membrane and cable. In total, we compute eight test cases. The membrane does not have bending stiffness; the shell does. The membrane includes a wrinkling model [85]. The shell model is from [66]. The cable includes a slacking model, which precludes compression. The bendingstabilized cable is based on the model in [79], where the bending stress is represented by the second area moment and the curvature. Figure 17 shows the mesh before connecting the membrane and cable. The material properties for the membrane, shell, cable, and bending-stabilized cable are shown in Tables 1 and 2. The meshes with C 0 and C 2 continuity between the membrane and cable will be called "C 0 Mesh" and "C 2 Mesh." Both are based on cubic B-splines. Figures 18 and 19 show the meshes. We hook the membrane and cables at the two corners and two ends shown in those figures, all at the same elevation. The gravity is 9.81 m/s 2 . We gradually shorten the

Remark 13
To obtain the settled solution in each test case computed, depending on the test case, we either use the steady-state formulation or the unsteady formulation with a relatively large time-step size. The unsteady formulation is used to avoid the matrix singularities associated with slacking. Figure 20 shows the settled solution for all eight cases. Looking at the cases with C 0 Mesh, we can clearly see that the rotational freedom is not transferred between the membrane and cables. In the cases with C 2 Mesh, we see smoothness along the 1D structure and its parametric-line continuation in the 2D structure, as expected, even if the 2D structure is membrane or the 1D structure is cable. Naturally, how local the smoothness is in computations with the C 2 Mesh, i.e., how small the radius of curvature is, depends on the element Therefore, the decision to seek smoothness or just continuity would depend on those element sizes.

Concluding remarks
We have presented a T-splines computational method and its implementation for structural analysis where structures with different parametric dimensions are connected with continuity and smoothness. We derived basis functions that give us the desired smoothness between structures with 2D and 1D parametric dimensions. The derivation involves proper selection of a scale factor for the knot vector of the 1D structure and results in new control-point locations. While the method description focused on C 0 and C 1 continuity, paths to higherorder continuity were marked where needed. In presenting the method and its implementation, we referred to the 2D structure as "membrane" and the 1D structure as "cable." The method and its implementation are of course applicable also to other 2D-1D cases, such as shell-cable and shell-beam structures. When the membrane and cable are connected with smoothness, the strain and rotational freedoms are transferred between the two structures. For easy and efficient implementation of the method, we introduced the Bézier extraction row operators used in obtaining the basis functions. We presented test computations for combinations for four structural models and two meshes. The structural models were membranecable, membrane-bending-stabilized-cable, shell-cable, and shell-bending-stabilized-cable. The meshes were with C 0 continuity and C 2 continuity. The computations clearly demonstrate how the method performs in the classes of structural mechanics problems targeted.

A.1 Before connecting the cable and membrane
The example in Fig. 1 has two quadratic B-spline patches, one for the membrane and one for the cable. The membrane knot vectors in the two directions are Ξ Ξ Ξ 1 M = {0, 0, 0, 1, 2, 2, 2} and Ξ Ξ Ξ 2 M = {0, 0, 0, 1, 2, 3, 3, 3}. The cable knot vector is Ξ Ξ Ξ 1 L = {0, 0, 0, 1, 1, 1}. Figure 21 shows a mesh with the same basis functions as those in the mesh in Fig. 1. Just the geometry is more complex. We provide in Table 3 the global basis function indices used in each element, and corresponding element-wise indices for the directional basis functions. In Table 4, we provide the Bézier extraction row operators corresponding to those element-wise indices.
A.2 After connecting the cable and membrane with C 0 continuity Figure 22 shows the mesh obtained by upgrading the mesh in Fig. 21 to C 0 continuity in connecting the membrane and cable. The new mesh is no longer represented by B-splines. We provide in Table 5 the global basis function indices used in each element and the corresponding element-wise indices for the directional basis functions. In Table 6, we provide the Bézier extraction row operators corresponding to those element-wise indices. Figure 23 shows the mesh obtained by upgrading the mesh in Fig. 21 to C 1 continuity in connecting the membrane and cable. Like the mesh in Fig. 22, it is no longer represented by B-splines. We provide in Table 7 the global basis function Fig. 21 The membrane and cable elements before connecting them. The basis functions are the same as those in Fig. 1  Fig. 22 The membrane and cable elements after connecting them with C 0 continuity. The basis functions are the same as those in Fig. 2 indices used in each element and the corresponding elementwise indices for the directional basis functions. In Table 8, we provide the Bézier extraction row operators corresponding to those element-wise indices. Table 3 The membrane and cable elements before connecting them. For each element e, the global indices a or b, and their element-wise indices k = a 1 and a 2 or k = b 1 . The corresponding Bézier extraction row operators are provided in Table 4   Table 4 The membrane and cable elements before connecting them. For element e and element-wise index k, the Bézier extraction row operator for the directional basis functions

B Derivation of the smoothness constraints given by Eqs. (48) and (49)
We do the derivation for a general case of the knot removals in connecting a membrane and cable. For p = 2, with the knot spans of {ΔΞ 1 , ΔΞ 2 , ΔΞ 3 }, where ΔΞ 2 is the knot span for the element in consideration, the Bézier extraction Table 5 The membrane and cable elements after connecting them with C 0 continuity. For each element e, the global indices a or b, and their element-wise indices k = a 1 and a 2 or k = b 1 . The corresponding Bézier extraction row operators are provided in Table 6 Fig. 23 The membrane and cable elements after connecting them with C 1 continuity. The basis functions are the same as those in Fig. 6 operator can be written as The membrane knot spans are {(ΔΞ M ) 1 , ΔΞ M , 0}. The last knot span is zero because it is at the end of the patch. Table 6 The membrane and cable elements after connecting them with C 0 continuity. For element e and element-wise index k, the Bézier extraction row operator for the directional basis functions Table 7 The membrane and cable elements after connecting them with C 1 continuity. For each element e, the global indices a or b, and their element-wise indices k = a 1 and a 2 or k = b 1 . The corresponding Bézier extraction row operators are provided in Table 8 With that, from Eq. (53), the Bézier extraction operator for the membrane is expressed as where c 1 = where s is as defined in Eq. (50). From that, we get and With that, we obtaiñ and We substitute them into Eqs. (46) and (47), and obtain the following equations: After some rearrangement, we get