Acta Polytechnica CTU Proceedings doi:10.14311/APP.2017.13.0039 Acta Polytechnica CTU Proceedings 13:39–42, 2017 © Czech Technical University in Prague, 2017 available online at http://ojs.cvut.cz/ojs/index.php/app ON DEVELOPMENT OF THE INTERACTIVE DESIGN TOOL BASED ON ISOGEOMETRIC ANALYSIS Edita Dvořáková∗, Bořek Patzák Czech Technical University in Prague, Faculty of Civil Engineering, Thákurova 7, 166 29 Prague 6, Czech Republic ∗ corresponding author: edita.dvorakova@fsv.cvut.cz Abstract. Isogeometric analysis is a new concept of Finite element method which has been proposed to bridge the gap between the CAD systems and the FEM solvers. In Isogeometric analysis, the same basis functions, typically splines or NURBS, are used for geometry description as well as for approximation of unknowns and thus the same model can be shared between CAD and IGA systems. This results in a higher accuracy and overall efficiency of the analysis. Many isogeometric elements have been already proposed and implemented into the existing finite element solvers, but the automatic connection with a CAD system is usually missing. Our goal is to develop such a connection and to provide a tool which would interactively run the analysis when the model changes. This approach can enormously enhance the design process as it can provide the basic knowledge about the structural behavior already in conceptual design phase. Keywords: curved beam element, discrete shear gap, interactive design, isogeometric analysis. 1. Introduction Finite element solvers are the integral part of a nowa- days structural design. Once the architectural model is developed in CAD (Computer-Aided Design) sys- tem, the structure is further analyzed in FEM (Finite Element Method) solver. There are two independent models (CAD model and FEM model) and thus when one of them changes during the design process the other has to be correspondingly (and usually manu- ally) adjusted, leading to the higher time-consumption. Moreover, the inaccuracies in geometry description are usually introduced as the CAD and FEM models use different functions for geometry representation. The standard finite element method usually uses polynomial basis functions, while the CAD systems are based on splines. The use of spline functions for geometry and unknown approximations in finite element analysis enables to share one model between architectural design and structural analysis. This idea has been proposed by Huges et al. [1] and it is called Isogeometric analysis (IGA). Since 2005, when Isogeometric analysis has been introduced, many researchers have been focusing on its development and the results proved the advantages of IGA over standard FEM in many fields. Our focus is placed on its application for the analysis of curved beams as the curved geometries can especially benefit from the isogeometric approach. We use the formu- lation of curved beam element presented by Bouclier et al. [2] with the locking treatment based on Dis- crete Shear Gap (DSG) method [3, 4]. The element has been implemented into OOFEM finite element solver [5] and our main goal is to develop a connection between CAD system and the solver to provide a tool encapsulating the analysis within the architectural design environment. 2. Curved beam element The presented two-dimensional curved beam element is based on Timoshenko beam theory. It has three inde- pendent unknowns, axial displacement ut, transverse displacement un and rotation θ. A strain-displacement matrix B defined as ε = Br, (1) where ε = {εm,γs,χb}T and r = {ut,un,θ}T , is de- rived from the relations for membrane strain εm, trans- verse shear strain γs, and bending strain χb εm = u′t − un R , γs = ut R + u′n −θ, χb = θ′, (2) where s ∈ 〈0,L〉 runs along the midline of a beam and prime indicates a derivative with respect to s. By omitting the terms which are divided by radius of curvature R the straight beam formulation is obtained. A material matrix D results from N = EAεm, Q = kGAγs, M = EIχb, (3) where N,Q and M are axial force, transverse shear force and bending moment, respectively. Young’s modulus E, shear modulus G, shear coefficient k, area A and moment of inertia I are the material and 39 http://dx.doi.org/10.14311/APP.2017.13.0039 http://ojs.cvut.cz/ojs/index.php/app Edita Dvořáková, Bořek Patzák Acta Polytechnica CTU Proceedings cross-section characteristics. A stiffness matrix is then calculated as K = ∫ L 0 BT DB ds. (4) To evaluate the stiffness matrix the Gaussian quadra- ture involving p + 1 Gauss points over each knotspan is used, where p is a degree of approximation. 2.1. NURBS-based formulation In case of IGA, the Langrange polynomial basis func- tions used in FEM are replaced by splines. The pre- sented element is formulated using Non-Uniform Ra- tional B-Splines (NURBS) [6], which is the most widespread and developed technology in CAD indus- try. A NURBS curve is defined by its degree, control points with weights, and knotvector. NURBS basis functions are generated by weighting the B-splines functions, which are piecewise polynomial functions and are a special subset of NURBS functions. NURBS curve is given as the linear combination of basis func- tions and control points coordinates. The knotvec- tor is a set of non-decreasing parametric coordinates (knots) which influence the mapping. The knots divide a parametric space into knotspans (usually referred to as “elements” in IGA) and each knotspan is influenced only by limited number of control points. Therefore the knotspan defines which control points affect partic- ular part of the domain. See Figure 1 for an example of a NURBS curve geometry. Very important aspect of Isogeometric analysis is a non-interpolatory nature of the basis functions, as the NURBS do not satisfy Kronecker-delta property. Another distinct feature of the NURBS basis functions is the higher inter-element continuity. While only C0 continuity is provided by traditional polynomial functions, in case of IGA up to Cp−1 continuity can be achieved. By repeating knots in a knotvector, the continuity at that knot can be artificially reduced. 2.2. Locking treatment Shear locking is a well-known problem of standard Timoshenko beam elements and this phenomena per- sists also in isogeometric formulation. This prob- lem occurs because the displacements and rotation are treated independently and approximated by the functions of the same order. From Eq. 2 for shear and bending strains it can be seen that formula for bending strain results in lower order term than for- mula for shear strain, but actually this should be vice-versa. Also the field inconsistency within the term for shear strain causes that zero shear strain cannot be satisfied along entire patch when the same order interpolation of unknowns is used. Several lock- ing removal techniques have been proposed to unlock isogeometric beam element including reduced inte- gration, B̄-method, DSG method, see [2, 3] for the references. In this study, the DSG-based formulation is used. While this method can increase computational cost (as well as B̄-method), it has a big advantage over reduced integration: for the reduced integration the recovered strains can be evaluated only at Gauss points, DSG method enables to evaluate them along the entire patch. Discrete shear gap method has been originally de- veloped for the standard finite elements [3] but Echter and Bischoff [4] have extended its use also for isoge- ometric analysis. The approach can be divided into several steps yielding the modified B matrix. The idea is not to satisfy the equation γs = ut R + u′n −θ (5) pointwise, but in integral sense. The shear contri- bution uγ hi n (called “shear gaps”) to the deflection is obtained by integration of Eq. 5 as uγ hi n = ∫ si 0 γhs ds = ∫ si 0 ut R + u′n −θ ds = B DSGr, (6) at collocation points si calculated as Greville abscissa of the control points [6]. Modified transverse displace- ments are interpolated using NURBS basis functions Ni uγ modh n = n∑ i=1 Niũ γhi n . (7) Please note, that discrete shear gaps ũγ hi n are non- interpolatory in isogeometric analysis. To expressed the discrete shear gaps by values at the control points the transformation matrix A is derived{ uγ h n } = A { ũγ h n } , Aij = Nj(si), (8) where { uγ h n } are interpolatory values of shear gaps at control points and Nj(si) is the jth-basis function evaluated at ith-collocation point. The modified shear strain is then given as γmod h s = n∑ i=1 N′iũ γhi n . (9) Combining Eqs. 6-9 results in the modified strain- displacement matrix B B = N′A−1BDSG, (10) which is used to evaluate stiffness matrix K (Eq. 4). 3. Interactive design tool Thanks to the isogeometric formulation, the presented element can be easily connected to CAD system. The possibility of seamless connection between CAD and FEA is the major benefit of IGA. This goal can be achieved even without use of isogeometric approach (see e.g. [7]), nevertheless the costly transformations 40 vol. 13/2017 On Development of the IGA based Interactive Design Tool knots control points s = L s = 0 s ξ ξ1 = 0 ξ2 = 0.5 ξ3 = 1 knot span patch y x Figure 1. Description of NURBS geometry. between CAD and FEM are necessary in such a work- flow. Firstly, the geometrical model is completed by analysis data, such as loads, boundary conditions and material characteristics, then the model is trans- formed into suitable computational model and the analysis mesh is generated, and finally the analysis is performed in FEM solver. To visualize the results within the CAD system the model has to be trans- formed back to the original model. With the use of Isogeometric analysis, once the model is completed with analysis data, the analysis can be performed di- rectly [8] and no transformation is needed for both the analysis and the visualization. For our purposes, the use of CAD system Rhinoceros (Rhino) [9] has been chosen. The ge- ometry representation of Rhino is based on NURBS, moreover Rhino enables plug-in development and use of built-in tool Grasshopper [10]. Grasshopper is a vi- sual programming interface within Rhino which can directly access NURBS geometry. In addition, vi- sual programming is intelligible even to a user with no programming background, so the user can define additional analysis data at this level. The workflow of the tool is illustrated in Figure 2. Firstly, the NURBS geometry specified by a user in Rhino is passed to the Grasshopper, where the model is completed with boundary conditions and material and crosssection characteristics. Currently, these data are provided in the text format similar to one used in OOFEM, but in the next phase of the development, the special plug-ins can be developed. Once all the analysis data are collected, the python script within Grasshopper together with other avail- able Grasshopper tools directly passes these data to OOFEM and run the analysis. When the analysis is finished, the output file is uploaded back to the Grasshoper, where the visualization of the results is carried out. Schematic illustration of Grasshopper environment with results visualized in Rhino can be seen in Figure 3. Very important feature of the developed tool is, that the Grasshopper (and in turn also OOFEM) re- computes the results immediately when something changes. No matter whether the user adjusts the GEOMETRY MODEL RESULTS ANALYSIS DATA GEOMETRY MATERIAL, CROSSSECTION BOUNDARY CONDITIONS ANALYSIS RESULTS VISUALISATION RHINOCEROS 3D GRASSHOPPER GRASSHOPPER OOFEM Figure 2. Interactive design tool. geometry or changes the boundary conditions, the results are interactively updated. This provides the immediate knowledge about behaviour of the designed structure and eases the process of finding optimal de- sign. In the future work the designer will receive graphical indication of sections where load-bearing capacity or stability limits are achieved. Such an ap- proach would allow to establish viable conceptual design from both architectural and structural point of view without the necessity of understanding details of structure behaviour. 4. Conclusions The formulation of the curved beam element with locking treatment has been presented and the element has been implemented into OOFEM finite element code. Thanks to the isogeometric formulation, the automatic connection of CAD system and the analysis has been provided or modified. The developed tool enables to display the results in Rhino immediately 41 Edita Dvořáková, Bořek Patzák Acta Polytechnica CTU Proceedings Figure 3. Schematic illustration of Grasshopper environment with results visualized in Rhino. The example corresponds to a cantilever circular arc beam (black) subjected to horizontal and vertical load at its tip. The deformed geometry is in red and the bending moment is in blue. after the geometry and analysis data are provided. This gives the user instant knowledge about behaviour of the designed structure and enables to see differences of alternative designs. Acknowledgements The financial support of this research by the Grant Agency of the Czech Technical University in Prague (SGS project No. SGS17/043/OHK1/1T/11) is gratefully acknowl- edged. References [1] T. J. R. Hughes, J. A. Cottrell, Y. Bazilevs. Isogeometric analysis: CAD, finite elements, NURBS, exact geometry and mesh refinement. Comput Methods Appl Mech Engrg 194:4135–4195, 2005. [2] R. Bouclier, T. Elguedj. Locking free isogeometric formulations of curved thick beams. Comput Methods Appl Mech Engrg 245-246:144–162, 2012. [3] K. Bletzinger, M. Bischoff, E. Ramm. A unified approach for shear-locking-free triangular and rectangular shell finite elements. Computers and Structures 75:321–334, 2000. [4] R. Echter, M. Bischoff. Numerical effiency, locking and unlocking of NURBS finite elements. Comput Methods Appl Mech Engrg 199:374–382, 2010. [5] B. Patzák. OOFEM project home page, 2017. http://www.oofem.org. [6] L. Piegl, W. Tiller. The NURBS Book. Springer-Verlag Berlin Heidelberg, New York, 1997. [7] L. Svoboda, J. Novák, L. Kurilla, J. Zeman. A framework for integrated design of algorithmic architectural forms. Advances in Engineering Software 72:109–118, 2014. [8] M. Hsu, C. Wang, A. Herrema, et al. An interactive geometry modeling and parametric desigh platform for isogeometric analysis. Computers and Mathematics with Applications 70:1481–1500, 2015. [9] Rhinoceros, 2017. http://www.rhino3d.com/. [10] Grasshopper, 2017. http://www.grasshopper3d.com/. 42 http://www.oofem.org http://www.rhino3d.com/ http://www.grasshopper3d.com/ Acta Polytechnica CTU Proceedings 13:40–43, 2017 1 Introduction 2 Curved beam element 2.1 NURBS-based formulation 2.2 Locking treatment 3 Interactive design tool 4 Conclusions Acknowledgements References