Journal of the Computational Structural Engineering Institute of Korea. 2018. 127-132
https://doi.org/10.7734/COSEIK.2018.31.3.127

ABSTRACT


MAIN

1 서 론

추계론적 유한요소해석(stochastic finite element method) 은 구조시스템 내의 구조인수들의 공간적, 시간적 임의성 (randomness)에 의한 구조의 응답에서의 변화도 산정을 목적 으로 한다. 임의성을 가지는 인수로는 재료탄성계수, 포아송비 등의 재료인수의 임의성과 단면적, 두께 등의 기하학적 인수, 그리고 시간변화하중 및 무작위하중(random excitation) 등을 수 있다(Noh and Yoon, 2010). 이들 인수들은 임의성을 가 지는 것으로 가정되지만, 이들은 확률적으로 특정의 확률분포 (확률밀도함수)와 확률특성(평균, 표준편차)을 가지는 것으로 모델링하는 것이 일반적이며, 기능경사재료의 경우에도 유사한 모델링이 사용된다(Xu et al., 2016; Shaker et al., 2008). 기대할 수 있는 바와 같이 특성 인수의 임의성 정도는 인수 마다 다르며, 또한 재료는 물로 구조의 종류에 따라서도 다를 수 있다. 간단한 예로서, 콘크리트 재료탄성계수의 임의성 정도는 구조용 강재료의 그것에 비해 클 것으로 기대할 수 있다.

추계론적 해석을 통하여 얻은 응답변화도는 구조 거동의 신뢰도 해석 등에 적용할 수 있다. 그러나 응답변화도 자체만 으로는 신뢰도해석에 응용하는 데에 한계가 있으며, 원칙적으로는 응답의 확률밀도를 구하는 것이 필요하다. 그러나 수치적 해석 법(numerical method)에 속하는 몬테카를로 해석(Monte Carlo simulation: MCS)외의 해석적 방법(analytical methods)에서는 응답의 확률밀도함수를 제시하는 것은 제한 적이라고 할 수 있으며, 일부의 해석적 방법만이 응답의 확률 밀도함수를 제시하고 있다(Grigoriu, 2000).

추계론적 해석은 확정론적 해석(deterministic analysis)을 포함하는 해석으로 확정론적 해석에 비하여 그 복잡도가 상대 적으로 높다고 할 수 있다. 그런 이유로 정정해석에 대한 연구가 주를 이루어 왔으며, 동적거동에 대한 추계론적 해석은 무작위 하중에 의한 거동 해석 등 상대적으로 제한적으로 수행되어 왔다 (To, 1992).

본 연구에서는 구조의 고유치 및 고유벡터의 해석에서 구조 재료의 공간적 임의성을 고려한 추계론적 유한요소해석 정식 화를 제안하고자 한다. 특히 복합재료의 일종인 기능경사재료 보 (functionally graded material beam)의 고유벡터의 변화도 산정을 위한 해석 정식화를 섭동법에 기초하여 제시한다. 정식 화에서 고려한 불확실인수로는 기능경사재료 보의 중앙면에서 보의 축을 따르는 방향으로 공간적 불확실성을 가지는 재료탄 성계수를 택하였다. 몬테카를로 해석을 함께 수행하여 제안한 정식화가 합리적인 결과를 주는지 확인하였다.

2 본 론

이 장에서는 고유치 문제에서 구조강성의 공간적 임의성에 의한 고유벡터의 변화도를 산정하기 위한 정식화에 대해서 논 한다.

2.1 확률변수

본 연구에서 대상으로 하는 불확실인수(random parameter) 는 구조재료의 탄성계수이다. 탄성계수는 구조영역(D ) 내에서 위치에 따라 다른 값을 가질 수 있으며, 이는 다음의 식으로 표현할 수 있다.

(1)
E(x)=Eo[1+f(x)],xD

여기서, Eo는 탄성계수의 평균치를 나타내며, f(x) 는 탄성계 수의 공간적 임의성(spatial randomness)을 나타내는 추계 장 함수(stochastic field function)이다. 유한요소법을 적용한 해석의 경우, 확률변수(random variable)는 개별 유한요소 내에서 평균적인 의미를 가지는 것으로 가정할 수 있다. 중앙 점법(Mid-point rule)을 적용할 경우 확률변수 벡터는 r=riT,ri=r(xc),(i=1,2,,ne) 가 된다. 여기서 ne 는 구조 모델링에 사용된 유한요소의 개수이고, xc는 개별 유한요소 중앙점의 좌표이다. 국부평균법(local average method)을 적용한다면 개별 유한요소 내 np개의 위치에 대하여 ri=rkk(xk)/np,(k=1,2,,np)로 나타낼 수 있다(Zhu et al., 1992). 어느 경우이든, 식 (1)로부터 확률변수의 평균치벡터 ro는 0이 된다(2중 반복 아래첨자는 합을 의미하며, 이 규칙을 지속 적용한다.).

2.2 고유치의 선형 변화도

고유치문제는 다음의 식으로 표현된다.

(2)
(KΛM)X=0

여기서, K 는 강성행렬, Λ는 고유치행렬, M 은 질량행렬, 그리고 X 는 고유벡터행렬을 의미한다. 고유치행렬 Λ를 평균 확률변 수를 중심으로 선형 Taylor 전개하면 다음과 같다.

(3)
Λ=Λo+ΔriΛri+

여기서, r은 탄성계수에 대한 확률변수를 나타내며, Λ0 = Λ(ro)를 의미한다. 확률변수 r의 평균중심 변량은 Δr = r - ro 로 표기하였다.

식 (3)에서 고유치의 평균(mean)과 공분산(covariance)은 각각 다음과 같이 쓸 수 있다.(4a)

(4a)
E[Λ]=Λo

(4b)
Coυ[Λ]=E[(ΛΛo)(ΛΛo)T]

식 (4b)의 공분산은 식 (3)으로부터,

(5)
Coυ[Λ]=Λri(Λrj)TE[rirj]

식 (5)에는 고유치의 확률변수에 대한 1차 미분이 사용되고 있음을 알 수 있다. 이는 식 (2)에 대한 편미분을 통하여 다음과 같이 얻을 수 있으며, 이 결과를 이용하면, 고유치에 대한 공분산 식 (5)는 다음과 같은 최종식으로 귀결된다.

(6)
Coυ[Λ]=XoTKriXoXoT(Krj)TXoE[rirj]

2.3 고유벡터의 선형 변화도

고유벡터의 변화도를 얻기 위하여 식 (2)의 고유치문제 식의 모든 행렬에 대하여 식 (3)과 같이 선형 Taylor전개를 수행하고, 이를 원식 (2)에 대입하여 정리하면, 동일 미분차수에 대하여 다음의 결과를 얻을 수 있다.

(7)
0차:(KoΛoMo)Xo=0VoXo=0

(8)
1차:(KriΛriMo)Xo+VXri=0

여기서, Vo=KoΛoMo이다. 식 (7)에 사용된 평균항은 개별 행렬에 대한 Taylor 전개식의 첫 항들이며, 식 (8)로부터 고유벡터 첫 미분항 Xri 를 다음과 같이 얻을 수 있다.

(9)
Xri=Vo1(KriΛriMo)Xo

그러므로 고유벡터를 식 (10)과 같이 선형 Taylor 전개하면,(11)

(10)
X=Xo+ΔriXri+

(11)
E[X]=Xo

고유벡터의 평균값은 다음과 같다.

이는 고유벡터의 평균치가 확정론적 해석의 결과와 동일함을 의미한다.

고유벡터의 공분산 Cov[X]는 다음과 같다.

(12)
Coυ[X]=E[(XXo)(XXo)T]=E[ΔriΔrjXri(Xrj)T]=E[rirj]Xri(Xrj)T

고유벡터의 확률변수에 대한 첫 미분식인 식 (9)를 이용하면, 식 (12)로 표현되는 고유벡터의 공분산은 다음과 같이 재정리 된다.

(13)
Coυ[X]=Vo1(KriΛriMo)Xo·=XoT(KrjΛrjMo)TVoTE[rirj]

식 (13)에서 고유치의 확률변수에 대한 1차 미분은 식 (6)에 적용된 결과를 이용한다. 식 (6), (13) 등에 나타난 확률변수 곱에 대한 평균은 함수의 정의에 따라 자기상관함수로 변형되며, 본 연구에서는 식 (14)로 주어진 지수함수를 사용한다.

(14)
E[rirj]=E[(ririo)(rjrjo)]=σff2exp{ξ24d2},ξ=xjxi

자기상관함수에서 σff2 는 추계장(f)의 분산, d는 불확실장의 상관관계거리(correlation distance), 그리고 ξ는 두 점 xi, xj 사이의 상대거리벡터를 나타낸다.

2.4 응답변화도에 대한 수치실험

앞 절에서 제시한 정식화에 대한 검증을 위하여 몬테카를로 해석을 적용한 수치실험을 수행한다. 본 연구에서는 불확실인 수로 택한 탄성계수의 공간적 불확실성을 식 (14)의 자기상관 함수에 적합한 추계장 표본으로 생성한다. 추계장 표본의 수치 생성에는 스펙트럼모사법(Shinozuka and Deodatis, 1996) 을 적용하였고, 해석에 사용한 표본의 수는 5,000으로 하였다. 상관관계거리 d는 변수로 설정하여 상관관계거리별 추계장에 대한 고유벡터의 변화도를 제시하고자 한다.

상관관계거리는 구조 내 임의 두 점에서의 불확실인수에 대한 상호상관관계를 정의하게 되며, 몬테카를로 해석의 경우 상관 관계거리는 수치생성 된 추계장을 직접 적용하여 실험적으로 고려하며, 제안한 정식화에서는 식 (14)의 함수를 식 (12)에 적용하여 이론적으로 고려하게 된다.

2.5 수치예제 및 고찰

해석 대상구조는 Fig. 1과 같이 길이 L =20.0m, 두께 t =1.0m, 질량밀도는 ρ =2,780kg/m3인 단순지지 기능경사 보이다.

https://cdn.apub.kr/journalsite/sites/jcoseik/2018-031-03/2COSEIK-31-127/images/COSEIK-31-127_F1.jpg
Fig. 1.

Simply supported FGM beam

탄성계수는 두께 방향으로 식 (15)를 따라 변화하며, 중앙 에서의 탄성계수는 E =70.0GPa이다.

(15)
E(z)=Eeβ|z|

보의 상·하단(zt/2)에서의 탄성계수는 2E 로 가정한다. 탄성계수의 분산계수는 σf =0.1로 가정하였다. 탄성계수의 불확실 분포는 z = 0에서 x축을 따라 임의적으로 분포하는 것만을 고려하며, 두께 방향 탄성계수의 변화는 x에서의 값이 식 (15)를 따라 분포하는 것으로 가정한다.

고유모드형상은 질량행렬에 대하여 정규화된 형상을 사용한다.

Fig. 2는 탄성계수의 공간적 불확실성에 의한 고유치의 공분산(COV(λ))을 보여준다. 가로축은 상관관계거리 d의 보의 길이 L 에 대한 비율로 나타내었다. 그래프에서 보는 바와 같이 추계장의 상관관계 거리가 증가함에 따라 고유치 공분산은 초기에 가정한 추계장의 공분산값(0.1)으로 수렴함을 알 수 있다. 또한 제안한 정식화(PM: perturbation method)에 의한 결과는 몬테카를로 해석(MCS)과 매우 유사한 결과를 보이고 있다.

https://cdn.apub.kr/journalsite/sites/jcoseik/2018-031-03/2COSEIK-31-127/images/COSEIK-31-127_F2.jpg
Fig. 2.

COV of eigen-value in terms of correlation distance d

Fig. 3은 예제 구조의 첫 번째 고유모드에 대한 평균치 (mean)와 표준편차(SD) 그리고 분산계수(COV)를 나타낸 것으로, 각각 상관관계거리 d=0.1, 5, 100인 경우에 대한 결과를 제시하고 있다.

https://cdn.apub.kr/journalsite/sites/jcoseik/2018-031-03/2COSEIK-31-127/images/COSEIK-31-127_F3.jpg
Fig. 3.

그림에서 보는 바와 같이 고유모드의 평균형상은 추계장의 상관관계거리 d에는 무관하게 모두 같은 것을 볼 수 있다. 그러나 고유모드 형상의 표준편차는 상관관계거리에 따라 다르게 나타나며, 상관관계거리와 비례관계에 있지도 않은 것을 알 수 있다. 현재 예의 경우, d=5에서 다른 경우보다 큰 표준 편차가 나타나고 있으며, 특히 d가 큰 경우에는 표준편차 값이 상대적으로 적게 나타나는 것을 볼 수 있다. 그 이유는, 상관 관계거리 d가 클 경우 추계장은 f(x)= c의 상태로 근접하기 때문으로, 이 경우 모든 표본에서 유사한 형태의 고유모드형상 이 나타나게 되어 이들에 대한 분산정도인 표준편차는 매우 작은 값을 나타내게 된다.

특히, 표준편차의 형상은 평균모드형상과 다름을 알 수 있다. 모드형상의 경우 접선경사가 0인 점이 하나인 반면, 표준편차는 두 개의 점이 존재함을 볼 수 있다.

Fig. 3(c)에서 볼 수 있는 바와 같이, 모드형상의 분산계수는 구조계에 걸쳐 동일한 값을 가지지 않으며 불확실인수의 분산 계수보다 적은 값을 나타냄을 알 수 있다.

섭동법에 따른 해석결과는 5,000개의 표본을 수치 생성 하여 해석한 몬테카를로 해석과 거의 동일한 결과를 나타내었다. Table 1은 섭동법과 몬테카를로 해석의 결과 비교로서, 모드 형상의 표준편차가 다른 경우보다 크게 나타난 d=5.0에서의 결과를 비교하고 있다. Table 1에서 보는 바와 같이 섭동법은 몬테카를로 해석과 거의 유사한 결과를 보여주고 있다. 오차는 몬테카를로 해석에 대한 제안 방법의 비율로 산정하였다. 몬테 카를로 해석에 의한 모드형상의 평균, 표준편차, 그리고 분산 계수는 Fig. 3의 그림과 거의 동일하며, 이에 대한 별도의 그림은 제시하지 않았다.

Table 1

Comparison of standard deviation between perturbation and monte carlo analyses(1st mode, d=5.0)

xPM(×E5)MCS(×E5)Error(%)
11.9641.976-0.645
23.7323.757-0.664
35.1315.166-0.692
46.0266.070-0.730
56.3366.386-0.783
66.0396.091-0.860
75.1755.226-0.978
83.8513.897-1.181
92.2772.314-1.592
101.2271.243-1.309

Fig. 4는 예제 구조의 두 번째 고유모드에 대한 섭동법 해석 결과로서 평균모드(Fig 4a), 표준편차(Fig 4b), 그리고 공분 산(Fig 4c)을 상관관계거리(d)에 따라 보여주고 있다. 결과의 정성적 특성은 첫 번째 고유모드와 동일하게 설명할 수 있다. 두 번째 모드의 경우 접선경사가 0인 점이 평균모드에서는 2개, 표준편차에서는 이보다 한 개가 많은 3개의 점에서 나타나고 있음을 볼 수 있다. 경계를 제외한 점에서 표준편차가 최소가 되는 위치는 평균모드의 접선경사가 0인 점과 반드시 일치하는 것은 아니다.

https://cdn.apub.kr/journalsite/sites/jcoseik/2018-031-03/2COSEIK-31-127/images/COSEIK-31-127_F4.jpg
Fig. 4.

Fig. 4(c)는 1차 모드와 동일하게, 모드형상의 분산계수가 구조계에 걸쳐 동일한 값을 가지지 않으며 불확실인수의 분산 계수보다 적은 값을 가짐을 보여주고 있다. Fig. 4(c)에서 x=10에서의 분산계수 값은, Fig 4(a)에서 보는 바와 같이 평균 고유모드가 0이기 때문에 매우 큰 값을 가지에 되어, x=9와 x=11에서의 값을 평균하여 임의 표기하였다.

Table 2에 보인 바와 같이, 두 번째 고유모드의 경우 몬테 카를로 해석 대비 제안 해석법의 오류 비율은 첫 번째 고유모 드와 유사한 것으로 나타남을 볼 수 있으며, 두 해석이 서로 매우 유사함을 알 수 있다.

Table 2

Comparison of standard deviation between perturbation and monte carlo analyses (2nd mode, d=5.0)

xPM(×E5)MCS(×E5)Error(%)
14.8584.879-0.423
28.3898.426-0.438
39.6429.687-0.470
48.2338.278-0.540
54.4584.493-0.770
62.1592.174-0.697
77.8917.929-0.476
813.46413.545-0.601
917.34417.473-0.736
1018.73118.901-0.895

해석결과에서 고찰할 수 있는 또 하나의 특징은 고유모드의 차수가 증가할수록 표준편차의 차수도 증가한다는 것이다. Fig. 3과 Fig. 4의 비교에서 알 수 있듯이, 첫 번째 고유모드 에 대한 표준편차 값은 두 번째 고유모드에 대한 표준편차 값 보다 작다. 즉, 고유모드의 변화도는 차수가 증가하면서 증가 하는 경향을 나타낸다.

3 결 론

본 연구에서는 구조재료에 나타나는 재료상수의 공간적 불확 실성에 의한 구조 고유모드의 변화도 산정에 대하여 논하였다. 유한 요소해석을 이용한 정식화를 제시하였으며, 예제를 통하여 그 타 당성을 보여 주었다. 해석결과의 타당성을 검토하기 위하여 몬테 카를로 해석을 동시에 수행하여 제안한 정식화에 의한 결과와 비교하였으며, 두 해석결과가 서로 매우 유사함을 보여 주었다.

고유모드의 표준편차는 평균모드에 비하여 접선경사가 0인 점이 하나씩 더 많게 나타났으며, 표준편차가 최소인 점과 평균 모드의 최대변위점은 일반적으로 일치하지는 않는다는 것을 알 수 있었다. 또한 고유모드 표준편차의 크기는 모드의 차수가 증가할수록 함께 커지는 것으로 나타난다. 즉, 고차의 모드가 저차의 모드에 비하여 더 큰 불확실성을 나타낸다. 모드형상의 분산정도는 추계장의 상관관계거리에 따라 상이하게 나타나며, 비례관계에 있지 않는 것으로 고찰되었다. 또한 정적거동에 의한 변위에서의 분산과 다르게 고유모드의 분산은 불확실인 수의 분산계수보다 적은 값을 가지며, 구조계에 걸쳐 동일한 값을 가지지도 않는다.

본 연구에서 제안한 결과는 추후 확장 연구를 통하여 구조계 동적거동에서의 응답변화도 해석에 적용할 수 있을 것으로 기대된다.

감사의 글

본 연구는 국토교통부/국토교통과학기술진흥원의 지원으로 수행된 연구 과제(과제번호 13IFIP-C113546-01)입니다.

References

1
Grigoriu, M., 2000. Stochastic Mechanics.. Int. J. Solids Struct., 37, pp.197-214.
10.1016/S0020-7683(99)00088-8
2
Noh, H.C. and Yoon, Y.C., 2010. Effect of Random oisson ?(tm)s Ratio on the Response Variability of Composite Plates.. J. Comput. Struct. Eng. Inst. Korea, 23, pp.727-737.
3
Shaker, A., Adberkragnabm, W., Tawfik, M. and Sadek, E., 2008. Stochastic Finite Element Analysis of the Free Vibration of Functionally Graded Material Plates.. Comput. Mech., 41, pp.707-714.
10.1007/s00466-007-0226-2
4
Shinozuka, M. and Deodatis, G., 1996. Simulation of Multi-dimensional Gaussian Stochastic Fields by Spectral Representation.. Appl. Mech. Rev., 49, pp.29-53.
10.1115/1.3101883
5
To, C.W.S., 1992. A Stochastic Version of the Newmark Family of Algorithms for Discretized Dynamic Systems.. Comput. Struc., 3, pp.667-673.
10.1016/0045-7949(92)90399-K
6
Xu, Y., Qian, Y. and Song, G., 2016. Stochastic Finite Element Method for Free Vibration Characterisitcs of Random FGM Beams.. Appl. Math. Model., 40, pp.10238-10253.
10.1016/j.apm.2016.07.025
7
Zhu, W.Q., Ren, Y.J. and Wu, W.Q., 1992. Stochastic FEM based on Local Average of Random Vector Fields.. J. Eng. Mech., 118, pp.496-511.
10.1061/(ASCE)0733-9399(1992)118:3(496)
페이지 상단으로 이동하기