퇴화 쉘요소를 사용한 원통형 복합재 패널의 비선형 좌굴후 현상해석

1. 서 론

항공 우주산업의 발달과 더불어 가볍고 강성이 강한 재료의 개발과 설계에 대한 연구가 집중적으로 이루어졌다. 항공기 및 우주구조물은 그 사용 용도의 특성상 가볍고, 열에 대한 저항력이 크고, 강성이 높은 특성을 요구한다. 이러한 구조적 요구사항을 가장 적절하게 만족시킬 수 있는 재료가 복합재료이다.

복합재는 비강성, 비강도가 높고 열, 크랙 및 기타 구조적으로도 안정하기 때문에 거의 대부분의 항공우주 구조물은 복합재로 제작된다. 항공기의 날개 및 기타 구조물을 제작함에 있어서, 원통형의 얇은 패널의 사용은 필수적이다. 뿐만 아니라 압력용기, 보일러, 저장탱크 등 다양한 구조에 적용되는 원통형 구조물의 좌굴 설계는 오랜시간동안 연구의 주요 주제가 되었다^[1,2].

원통형 구조물은 그 구조의 형상적 특성상 좌굴현상이 발생하기 쉽다. 원통의 좌굴현상에 대한 초기연구에서 Batdorf^[1], Sobel^[2], Hoff^[3] 등은 이론적 접근방법을 시도하였다. 이들은 좌굴하중, 변형율, 응력 등을 계산하는데 이론적인 수식에 의한 접근법을 사용하였다.

컴퓨터의 발달과 더불어 수치해석을 통한 좌굴해석에 대하여 많은 발전이 있었다. 근래 더 높은 차원의 비선형 좌굴후 거동 해석에 대한 연구가 활발히 이루어 졌다^[5-7]. Koiter^[4]는 초기의 좌굴후 거동을 해석하기 위하여 섭동 접근법을 제시하였다. 이 방법은 나중에 다른 연구자에 의해서 좌굴후 비선형성이 강한 거동의 해석에 많이 적용되었다^[8,9].

Snap-through와 snap-back 현상이 발생하는 좌굴해석에서는 하중제어와 변위제어를 동시에 사용하는 Rike’s 방법을 사용하거나, arc-length 방법이 오랜시간동안 활발히 적용되었다^[10-12]. 좌굴`은 매우 불안정한 상태의 특이점에서 해를 구하는 것이기 때문에 해의 발산이 빈번하게 발생하여 수치해를 수렴시키는 것은 매우 어렵다. 특히 복합재의 경우에는 그 거동 현상이 등방성 재료에 비해서 더욱 불규칙적이고 비선형성이 강하기 때문에 비선형해석에서 가장 어려운 문제라고 할 수 있다.

본 연구에서는 이러한 비선형성이 강한 좌굴현상을 해석하는데 많이 사용되는 arc-length 방법을 사용하여 좌굴 임계점 분기 및 경로 해석을 하였다. 복합재 원통형 패널의 새로운 좌굴거동 경로를 찾기 위하여 arc-length 방법을 branch switching^[13] 수치해석 방법과 결합하였다.

등방성재료 해석을 통해서 새로운 좌굴 경로를 찾기 위하여 몇 개의 경로 분기 섭동점(bifurcation point)을 찾았고 이 해는 기존의 결과들과^[13,14] 매우 잘 일치 하였다.

이러한 섭동점을 찾기 위하여 유한요소 정식화에 모드형상 벡터(eigen vector)를 사용하였다. 일단 좌굴 경로에서 섭동점이 찾아지면, 2차 좌굴 경로(secondary path)를 따라가기 위하여 주 평형 경로(primary equilibrium path)로 부터 경로가 교체 된다. 이 거동 경로의 교체를 위하여 branch-switching 방법이 사용되었다. 좌굴이 발생할 때의 모드 값(eigen value)은 무한대로 발산한다. 접선 강성행렬과 일치하는 제로 모드 값(zero eigenvalue)의 모드 벡터는 2차 경로를 가리킨다. 이 새로운 2차 경로는 주 경로에서 분기한 모드 벡터의 섭동이다.

2. 유한요소 해석 정식화

2.1 퇴화 쉘요소 정식화

퇴화 쉘요소의 유한요소 정식화를 위한 벡터의 정의를 Fig. 1에 나타내었다. Fig. 1에서 단위벡터 e_r, e_s, e_t는 서로 직교하지 않으나, $$ e_{\overline{r}},e_{\overline{s}},e_t$$는 서로 직교하는 벡터이다. 자연좌표계를 r, s, t, 전역 좌표계를 x, y, z로 정의하였고 q개의 노드를 가지는 요소에 대하여 정식화하였다.

Fig. 1

Vector definitions and FEA element models of degenerated shell element

복합재 3D 쉘요소를 정의하기 위하여 두께방향으로 적층된 복합재 모델을 사용하고 강성행렬을 구성할 때 각 층에 정의된 재료의 강성행렬을 계산하여 모든 층에 대한 강성행렬을 구성하여 합산하는 방식으로 구성하였다. 쉘요소는 각 절점에 5개의 자유도를 가지고 있다.

쉘요소에서 변위를 정의하기 위해서는 형상함수를 N_k로 정의하고, α_k, β_k를 절점에서의 방향각으로 정의를 하면 Eq. (1)과 같이 나타낼 수 있다.

$\[ \begin{array}{l}u\left(r,s,t\right)=\sum _{k=1}^q{N}_k{u}_k+\frac{t}{2}\sum _{k=1}^q{a}_k{N}_k\left(-{}^0{V}_{2x}^k{\alpha }_k\times {}^0{V}_{1x}^k{\beta }_k\right)\\ v\left(r,s,t\right)=\sum _{k=1}^q{N}_k{v}_k+\frac{t}{2}\sum _{k=1}^q{a}_k{N}_k\left(-{}^0{V}_{2y}^k{\alpha }_k\times {}^0{V}_{1y}^k{\beta }_k\right)\\ w\left(r,s,t\right)=\sum _{k=1}^q{N}_k{w}_k+\frac{t}{2}\sum _{k=1}^q{a}_k{N}_k\left(-{}^0{V}_{2z}^k{\alpha }_k\times {}^0{V}_{1z}^k{\beta }_k\right)\end{array} \]$

(1)

여기서 각 절점의 방향을 나타내는 벡터 V^k는 아래의 Eq. (2)와 같이 단위벡터와 e와 직교벡터 V_n의 조합으로 계산된다.

$\[ {{}^{0}V}_1^k=\frac{e_y\times {{}^{0}V}_n^k}{‖e_y\times {{}^{0}V}_n^k‖} \]$

(2)

$\[ {\mathrm{V}}_{\mathrm{n}}^{\mathrm{k}}=-{{}^{0}V}_2^{\mathrm{k}}{\alpha }_k+{{}^{0}V}_1^{\mathrm{k}}{\beta }_k \]$

(3)

V_n은 쉘요소의 면에서 수직인 벡터로 Eq. (3)과 같이 정의된다. 자연좌표계와 전역좌표계의 상관성을 정의하기 위한 관계식을 정의하기 위하여 전역좌표계의 변위를 자연좌표계로 미분을 하면 Eq. (4)와 같다.

$\[ \left[\begin{array}{c}\frac{\partial u}{\partial r}\\ \frac{\partial u}{\partial s}\\ \frac{\partial u}{\partial t}\end{array}\right]=\sum _{k=1}^q\left[\begin{array}{c}\frac{\partial N_k}{\partial r}\left[\begin{array}{ccc}1& {tg}_{1x}^k& {tg}_{2x}^k\end{array}\right]\\ \frac{\partial N_k}{\partial s}\left[\begin{array}{ccc}1& {tg}_{1x}^k& {tg}_{2x}^k\end{array}\right]\\ N_k\left[\begin{array}{ccc}0& g_{1x}^k& g_{2x}^k\end{array}\right]\end{array}\right]\left[\begin{array}{c}{u}_k\\ {\alpha }_k\\ {\beta }_k\end{array}\right] \]$

(4)

여기서 $$ g_1^k=-\frac{1}{2}{a}_k {}^0{V}_2^k, { g}_2^k=\frac{1}{2}{a}_k {}^0{V}_1^k $$ 는 자연좌표계와 전역좌표계로의 변환을 위한 좌표회전 벡터의 의미를 가진다.

유한요소에서는 전역좌표계의 실제 면적을 자연좌표계의 면적으로 변환하는데 이때 사용되는 것이 자코비안(J)이다. 자코비안은 Eq. (5)와 같다.

$\[ \frac{\partial }{\partial \mathbf{x}}={\mathbf{J}}^{-1}\frac{\partial }{\partial \mathbf{r}} \]$

(5)

$\[ \left[\begin{array}{c}\frac{\partial u}{\partial x}\\ \frac{\partial u}{\partial y}\\ \frac{\partial u}{\partial z}\end{array}\right]=\sum _{k=1}^q\left[\begin{array}{c}\begin{array}{ccc}\frac{\partial N_k}{\partial x}& g_{1x}^k{G}_x^k& g_{2x}^k{G}_x^k\end{array}\\ \begin{array}{ccc}\frac{\partial N_k}{\partial y}& g_{1x}^k{G}_y^k& g_{2x}^k{G}_y^k\end{array}\\ \begin{array}{ccc}\frac{\partial N_k}{\partial z}& g_{1x}^k{G}_z^k& g_{2x}^k{G}_z^k\end{array}\end{array}\right]\left[\begin{array}{c}{u}_k\\ {\alpha }_k\\ {\beta }_k\end{array}\right] \]$

(6)

Eq. (5)의 자코비안을 행렬의 형태로 나타내면 Eq. (6)과 같다. Eq. (6)에 포함된 형상함수에 대한 전역 좌표계 x에 대한 미분은 다음과 같다.

$\[ \begin{array}{c}\frac{\partial N_k}{\partial x}=J_{11}^{-1}\frac{\partial N_k}{\partial r}+J_{12}^{-1}\frac{\partial N_k}{\partial s}\\ G_x^k=t\left(J_{11}^{-1}\frac{\partial N_k}{\partial r}+J_{12}^{-1}\frac{\partial N_k}{\partial s}\right)+J_{13}^{-1}{N}_k\end{array} \]$

(7)

여기서 J_i^-_j¹는 J^-1 행렬의 한 요소이다.

3. 원통형 패널 좌굴 해석

3.1 형상모델

원통형 패널의 형상에 대하여 등방성재료 및 직교이방성 복합재에 대하여 해석을 수행하였다. Fig. 2에 나타낸 바와 같이 곡률반경이 r = 2,540 mm, 두께 t = 6.35 mm, 길이 l = 508 mm인 등방성재료의 패널에 대해서 해석을 하였다.

Fig. 2

Geometry models of buckling analysis (isotropic mat.)

패널의 양끝단에서는 단순지지되는 경계조건이 주어졌다. 패널의 중앙에는 수직하중이 작용하고 있다. 본 모델에 대해서 E = 3,102.75 MPa, ν = 0.3의 값을 가지는 등방성 재료에 대해서 좌굴해석이 수행되었다. Fig. 3은 직교이방성 재료의 복합재의 물성치를 가지는 패널을 나타내고 있다. 등방성 재료와 같은 치수, 형상을 가지고 있다. 재료의 물성치는 EL = 3,300 MPa, ET = 1,100 MPa, G_LT = 660 MPa, ν_LT = 0.3 이다. 원통형 쉘 패널은 중앙에 하중을 가하면 좌굴이 발생하는 전형적인 구조이며 많은 연구가 수행되었다.

Fig. 3

Geometry models of buckling analysis (orthotropic mat.)

Fig. 2, 3의 모델에 대하여 좌굴거동의 주 경로, 2차경로 등을 해석하고 검토하였다.

3.2 등방성 재료 해석

좌굴거동은 모드형상과 매우 밀접한 관계가 있다. 외력에 의한 외란으로 강성이 가장 약한 부분으로 변위가 형성 되는 것을 나타내는 것이 모드이다. 때문에 좌굴거동은 이러한 모드의 형상과 일치하게 형성된다. 좌굴 하중은 이 모드에 대하여 계산된다.

등방성재료에 대해서 좌굴거동을 해석하기 위하여 모달 해석을 먼저 수행하였다.

Fig. 2에 나타낸 형상에 대하여 저차부터 6개의 모드 및 모드형상을 계산하였으며 결과를 Fig. 4에 나타내었다. 좌굴이 각각의 모드에 상응하는 형상으로 나타나기 때문에 각각의 모드에 대한 좌굴 거동 경로를 본 연구에서 계산하였다. Fig. 2에 나타낸 원통형 패널의 A, B, C, D 네 점에서 하중 변위 관계를 나타내면 Fig. 5와 같다. A점 에서의 거동은 널리 알려진 바와 같이 snap-through와 snap-back이 발생하며 C점에서는 더욱 복잡한 snap-back 현상이 나타난다. 한 구조물에서 좌굴은 다양한 형태로, 이론적으로 모드 수 만큼의 좌굴형상 개수가 존재한다. Fig. 6은 2.2절에서 언급한 한계점, 섭동점, 임계점 등의 정식화에 의한 좌굴의 2차 분기 경로를 나타낸다. 즉 패널의 중앙에 하중을 가하였을 때 주 좌굴거동 경로는 Fig. 5(a)와 같다. 2차 좌굴 경로는 이 주 좌굴 경로를 벗어나 새로운 경로를 따라서 좌굴 거동이 일어나는 것을 말한다. 2차 좌굴 경로는 주 좌굴 모드를 기반으로 2차, 3차, 4차, 5차 등 다양한 모드의 좌굴이 합성되는 것을 의미한다.

Fig. 4

Modal analysis results of cylindrical panel for isotropic material

Fig. 5

Force vs. displacement curves at location A, B, C, D in Fig. 2

Fig. 6

Nonlinear bifurcation buckling paths of hinged cylindrical isotropic shell, (a) deformation of secondary path_1, (b) deformation of secondary path_2, (c) deformation of secondary path_3

좌굴 경로의 새로운 분기점의 발생이 의미하는 것은 임의의 지점에서 좌굴의 주 모드가 아닌 다른 모드가 발생한다는 것을 의미한다. Fig. 5의 하중-변위 그래프에서 (a)-(b), (a)-(c), (a)-(d)의 좌굴 모드가 합성되는 것을 나타낸다. 즉 Fig. 6(a)는 패널의 중앙점에서 모드 2의 좌굴이 발생할 때 하중-변위 곡선을 나타낸다. 이 결과는 기존의 Zhou’s^[14]의 결과들과 잘 일치한다.

3.3 직교이방성 재료 해석

직교이방성 복합재로 구성된 구조물의 좌굴현상은 더욱 복잡하다. 방향성에 따라서 강성이 상이하기 때문에 좌굴의 거동은 등방성 재료에 비해서 불규칙적인 현상을 나타낸다. Fig. 7은 직교이방성 물성치를 가지는 단순지지된 원통형 패널의 하중-변위 곡선을 나타낸다. Fig. 3의 A점과 C점의 거동이 등방성 재료와 비교하였을 때 비선형성이 더욱 증가 되었다.

Fig. 7

Force vs. displacement curves at location A, B, C, D in Fig. 3

복합재의 하중-변위 곡선에서 곡선의 기울기가 급격하게 변하는 것은 전형적인 특징이다. 등방성 재료의 경우와 마찬가지로 snap-through와 snap-back 현상이 명확하게 나타난다. 직교이방성 복합재 패널의 경우에도 2차 좌굴 경로를 계산하였으며 Fig. 8에 나타내었다. Fig. 8(a)의 2차 경로_1의 경우에는 등방성 재료와 비슷한 경향을 나타낸다. 2차 경로_2의 경우 등방성과 유사하나 더욱 과도하게 변위가 형성된 것을 나타낸다. 경로_3의 경우는 하중-변위 곡선에서 특이점이 다수 존재한다.

Fig. 8

Nonlinear bifurcation buckling paths of hinged cylindrical isotropic shell, (a) deformation of secondary path_1, (b) deformation of secondary path_2, (c) deformation of secondary path_3

이와 같이 좌굴거동에서 2차 경로가 중요한 이유는 외란이 불규칙적으로 가해졌을 때 구조체의 좌굴 형태는 복잡하게 발생하며 구조체의 좌굴 경로가 다양한 형태로 진행될 수 있다는 것이다. 다양한 2차 경로에서 어느 경로로 거동을 나타낼 것인가는 외란의 특성에 기인하는 것이다. 다만 이러한 좌굴의 다양한 경로에 대한 구조체의 응력, 변형율, 파괴 등에 대한 자료를 미리 획득함으로써 구조체의 설계를 보다 신뢰성 있게 진행할 수 있다. 다양한 좌굴의 모드에 대한 이론적 비선형 수치해석으로 경로를 찾았다는 것이 중요한 의미를 가진다.

복합재 하중-변위 곡선에서 snap-back과 snap-through 현상이 다수 존재하며 이것은 특이점 즉 접선강성행렬 K가 singularity 행렬이 된다는 것을 의미한다. 본 연구결과를 바탕으로 다양한 적층복합재에 대한 2차 좌굴 경로의 계산도 신뢰성 있게 계산이 가능하며, 이것은 좌굴 연구에 있어서 큰 진전이라 할 수 있다.

4. 결 론

원통형 쉘 패널의 등방성 재료와 직교이방성재료에 대하여 좌굴 및 좌굴후 현상을 해석하였다. 기존의 연구에서 많이 수행되지 않았던 2차 좌굴 경로를 찾는 방법에 대하여 제시하였다. 본 연구에서는 3D 퇴화 쉘요소를 사용하여 복합재 좌굴 해석을 하였다.

쉘요소의 유한요소 해석 정식화를 제시하였고 이 요소를 사용하여 arc-length 방법과 branch-switching 방법을 결합하여 성공적으로 좌굴의 2차 경로를 찾았다.

비선형성이 심한 복합재 패널에 대하여 2차 좌굴경로를 찾는 것은 수렴성에 문제가 있기 때문에 가장 어려운 문제로 간주된다. 본 연구에서는 3D 복합재 쉘요소의 체계적인 유한요소해석 정식화를 통하여 좌굴 2차 경로를 성공적으로 신뢰성있게 계산할 수 있는 체계를 제시하였다. 추후 원통형 좌굴 거동 구조물의 최적화 설계를^[15] 위해서도 2차경로 해석은 필수적이라 할 수 있다.

Acknowledgments

본 연구는 대구연구개발 특구 기술사업화과제(2018-DG-RD-0003-01-101)의 연구수행으로 진행되었습니다.

References

• Batdorf, S. B., (1947), A Simplified Method of Elastic-Stability Analysis for Thin Cylindrical Shells, 2-Modified Equilibrium Equation, NACA Technical Note 1342.
• Sobel, L. H., (1964), Effects of Boundary Conditions on the Stability of Cylinders Subject to Lateral and Axial Pressures, AIAA Journal, 2(8), p1437-1440. [https://doi.org/10.2514/3.2572]
• Hoff, N. J., Soong, T. C., (1965), Buckling of Circular Cylindrical Shells in Axial Compression, International Journal of Mechanical Sciences, 7(7), p489-520. [https://doi.org/10.1016/0020-7403(65)90050-0]
• Koiter, W. T., (1970), The Stability of Elastic Equilibrium, Technical Report AFFDL-TR-70-25.
• Wagner, W., Wriggers, P., (1988), A Simple Method for the Calculation of Postcritical Branches, Engineering Computations, 5(2), p103-109. [https://doi.org/10.1108/eb023727]
• Wardle, B. L., Lagace, P. A., Tudela, M. A., (2004), Buckling of Transversely Loaded Composite Shells, Part 2: Numerical Analysis, AIAA Journal, 42(7), p1464-1464.
• Bayat, Y., EkhteraeiToussi, H., (2017), Exact Solution of Thermal Buckling and Post-Buckling of Composite and SMA Hybrid Composite Beam by Layerwise Theory, Aerospace Science and Technology, 67, p484-494. [https://doi.org/10.1016/j.ast.2017.04.029]
• Kordkheili, S. A. H., Soltani, Z., (2018), A Layerwise Finite Element for Geometrically Nonlinear Analysis of Composite Shells, Composite Structures, 196, p355-364. [https://doi.org/10.1016/j.compstruct.2017.12.022]
• Battini, J. M., Pacoste, C., Eriksson, A., (2003), Improved Minimal Augmentation Procedure for the Direct Computation of Critical Points, Computer Methods in Applied Mechanics and Engineering, 192, p2169-2185. [https://doi.org/10.1016/S0045-7825(03)00254-8]
• Hu, H. T., Chen, H. C., (2018), Buckling Optimization of Laminated Truncated Conical Shells Subjected to External Hydrostatic Compression, Composite Prat B, 135, p95-109. [https://doi.org/10.1016/j.compositesb.2017.09.065]
• Sayyad, A. S., Ghugal, Y. M., (2017), Bending, Buckling and Free Vibration of Laminated Composite and Sandwich Beams: A Critical Review of Literature, Composite Structures, 171, p486-504. [https://doi.org/10.1016/j.compstruct.2017.03.053]
• Sze, K. Y., Liu, X. H., Lo, S. H., (2017), Buckling, Postbuckling and Progressive Failure Analyses of Composite Laminated Plates Under Compressive Loading, Composite Part B, 120, p143-151. [https://doi.org/10.1016/j.compositesb.2017.03.066]
• Ou, X., Zhang, X., Zhang, R., Yao, X., Han, Q., (2018), Weak Form Quadrature Element Analysis on Nonlinear Bifurcation and Post-Buckling of Cylindrical Composite Laminates, Composite Structures, 188, p266-277. [https://doi.org/10.1016/j.compstruct.2018.01.007]
• Zhou, Y., Stanciulescu, I., Eason, T., Spottswood, M., (2015), Nonlinear Elastic Buckling and Postbuckling Analysis of Cylindrical Panels, Finite Elements in Analysis and Design, 96, p41-50. [https://doi.org/10.1016/j.finel.2014.12.001]
• Kim, Y. H., Han, S. Y., (2018), Topology Shape Optimization Scheme for Nonlinear Structures Based on Artificial Bee Colony Algorithm, Journal of the Korean Society of Manufacturing Technology Engineers, 27(4), p329-338.