구조/음향 연성을 고려한 에너지흐름유한요소해석을 위해 Fig. 1과 같이 음향요소와 구조요소간 생성되는 면경계(area junction)에 대한 파동전달해석(wave transmission analysis)이 요구된다. 식 (1)~(3)은 각각 유체1에서 θ_i로 평판에 입사하는 음압 입사파(p_i)와 θ_r로 반사하는 음압 반사파(p_r), 평판의 굽힘파(ξ)를 보여준다.

Fig. 1 
Wave transmission analysis of semi infinite air-plate-air area junction

여기서 ω는 파동의 각주파수(angular frequency, rad/s)이고 k_1x=k₁cosθ_i와 k_1y=k₁sinθ_r은 음파의 파수(wavenumber, rad/m)의 각각 x, y방향 성분을 나타낸다. A, B, L은 각각 음파의 입사파, 반사파, 평판 굽힘파의 복소 계수(complex coefficients)고, 경계면 방향 파수 매칭에 의해 k_1y=k_y이다.

유체 매질1 측 평판 표면에서 음압(p)과 매질 입자속도(particle velocity, v) 간의 관계식을 얻기 위해 식 (4)의 오일러공식(Euler formula)을 이용하면 식 (5)와 같은 각 파동의 크기 사이의 관계식을 유도할 수 있다.

여기서 ρ와 ρ₁은 각각 유체 밀도(kg/m³)와 유체1의 밀도이며 v는 jωξ이다.

유체 매질2에서 식 (6)으로 표현된 투과파(transmitted wave) 압력(p_t)의 크기(C)와 유체 표면(x=0) 입자 수직 변위(ξ)와 식 (4)를 이용하면 식 (7)이 유도된다.

여기서 k_2x=k₂cosθ_t와 k_2y=k₂sinθ_t는 음파의 파수(k₂)의 각각 x, y 방향 성분을 나타낸다. 또한 경계면 방향 파수 매칭에 의해 k_2y=k_y이고 ρ₂는 유체2의 밀도이다.

양쪽 유체 매질의 읍압에 의해 가진되는 평판의 x방향 수직변위의 운동방정식은 식 (8)과 같이 표현할 수 있다.

여기서 D_s는 평판의 굽힘강성(flexible rigidity)이고 ρ_s 와 h는 평판의 밀도와 두께이며, p^-는 유체 매질1에서 식 (1), (2)의 p_i와 p_r의 합을 의미하고 p⁺는 유체 매질2에서 투과파 p_t를 의미한다.

식 (1)~(7)을 이용하여 각 파동의 계수를 구하고 입사각에 대한 입사 면적을 고려한 파워투과반사계수를 구하면 식 (9)~(11)과 같다.

여기서 τ_aa12는 음향매질 1에서 2로의 파워투과계수를 의미하고 τ_ap는 음향1에서 평판 굽힘파로 파워투과계수, γ_aa11는 음향공간 1에서 1로의 파워반사계수를 의미하며, 각 인텐시티(intensity, W/m)는 식 (12)~(15)로 표현할 수 있다.

여기서 c₁과 c₂는 각각 유체1,2에서의 음파의 위상속도이다.

분산장(diffuse field)의 파워투과반사계수를 얻기 위해 입사각에 대한 평균된 식 (16)을 사용한다.

두께 0.01m의 철(E=2.1x10¹¹Pa, ρ=7800kg/m³, ν=0.28) 재질의 무한 평판에 입사파 쪽 유체는 공기, 투과파 쪽 유체는 해수인 모델에 대해 음향-구조-음향 면연성 파동전달 수치해석 결과는 Fig. 2와 같다. 해석 결과는 음향 파워반사계수(γ_aa11)가 가장 큰 값을 가져 대부분 입사 파워가 반사되며 반대편 해수 음장으로의 파워투과계수(τ_aa12)는 주파수가 높아짐에 따라 작아지는 것을 확인할 수 있다. 유체1에서 유체2로의 주파수 증가에 따른 파워투과계수가 작아지면 격벽 투과손실(transmission loss, TL)은 반대로 커진다는 의미이므로 TL의 주파수에 따른 일반적인 특성을 잘 반영하고 있다. 해당 파동전달해석 모델에서 전체 파동의 파워는 보존되기 때문에 반대로 평판으로의 투과계수(τ_ap)와 공기중으로 반사계수(γ_aa11)는 높아지는 경향을 보이는 것을 확인할 수 있다.

Fig. 2 
Wave transmission analysis of semi infinite fluid-plate-fluid area junction(acoustic excitation)

유체/구조 연성해석을 통해 동시에 실내소음과 평판의 진동응답 예측을 위해서는 Fig. 3에서와 같이 구조 가진의 경우 접한 유체로의 파동전달해석도 수행되어야 한다. 본 연구에서는 기존 편측 유체 접수 경우 파동전달해석 (Bitsie, 1996)을 비상관관계(uncorrelated)인 양측 유체 접수 모델로 새롭게 확장하였다.

Fig. 3 
Power transfer relationship from structure to acoustic

두께가 일정한 평판에서 음향으로 전달된 파워는 식 (17)로 표현할 수 있다.

여기서 η_rad,1과 η_rad,2는 각각 음향 공간1,2로의 방사손실계수(radiation loss factor)이고 e는 평판의 굽힘 에너지밀도(flexural energy density)이다.

에너지중첩원리(principle of energy superposition)에 의해 평판에서 인텐시티는 산란 인텐시티(I_scat)와 입사 인텐시티(I_inc)의 합으로 표현된다.

여기서 c_gf는 평판 굽힘파의 에너지전달속도이며 평판에서 평판으로 파워반사계수(γ_ppff)는 산란 인텐시티와 입사 인텐시티의 비로 표현할 수 있다.

식 (17)과 식 (20)을 이용하면 평판에서 음향공간1,2로의 파워투과계수와 평판으로 파워반사계수는 식 (21)~(23)과 같이 얻어진다.

여기서 σ_rad,1과 σ_rad,2는 각각 음향 공간1,2로의 방사효율(radiation efficiency)이며, 음향임피던스 비(acoustic impedance ratio, β)는 식 (24), (25)와 같다.

공기(c_air=343 m/s, ρ_air=1.205 kg/m³)와 해수(c_sea=1500 m/s, ρ_sea=1024.7 kg/m³)로 접한 철(E=2.1x10¹¹Pa, ρ=7800 kg/m³, ν=0.28) 평판(1 m×1m×0.01m)모델에서 평판으로부터 굽힘파가 입사할 경우 파동전달해석이 수행되었다. 방사효율은 밀도가 높은 유체의 경우 부가질량효과가 반영되지 않지만 밀도 차이로 인한 효과를 고려할 수 있고 모든 주파수 대역에서 값을 수식으로 얻기 위해 Leppington의 식을 사용하였다 (Leppington, 1996). Fig. 4(a)는 유한 평판과 무한 평판에서 공기와 해수로의 방사효율을 보여주는데 해수의 임계주파수(23kHz)는 공기의 임계주파수(1.2kHz)보다 높기 때문에 해수의 임계주파수 아래에서는 공기로의 방사효율이 높다. 식 (21), (22)를 통해 10kHz 이내에서는 평판에서 입사된 파워는 공기 음장보다 해수 음장으로 더 많이 투과되는 것을 확인할 수 있다. 또한 Fig. 4(c)에서 파워 보존에 의해 모든 주파수 영역에서 파워투과반사계수의 합이 1이 되는 것을 확인할 수 있다.

Fig. 4 
Wave transmission analysis of fluid-plate-fluid area junction(plate excitation)

일반적으로 해수와 같은 밀도 높은 유체가 접수되면 방사가 잘 안되는 임계주파수 아래에서 부가질량효과로 인해 구조물의 진동특성이 진공 상태나 공기 중에서와 크게 달라진다. 하지만 기존 연구에서는 접수로 인한 구조물의 굽힘파의 파수 변화와 이에 따른 에너지전달속도(energy transfer velocity, m/s)의 변화를 고려하지 않았다 (Kim, 2007; Kwon, 2009; Kwon et al., 2012).

본 연구에서는 해수와 같은 고밀도 유체(heavy fluid)가 접수했을 경우에도 구조진동 에너지흐름해석의 신뢰성을 높이기 위해 임의의 유체가 양측접수 상태에서 평판의 굽힘파의 에너지전달속도를 유도하였다.

접수효과를 반영한 굽힘파의 에너지전달속도 유도를 위해 파동가진(wave excitation)을 받는 일차원 평판의 운동방정식과 2차원 음향방정식(helmholtz equation)은 다음과 같다.

여기서 u는 평판의 수직변위(m), κ는 가진파동의 파수(rad/m)이다. f는 가진 파동의 크기이며 k는 음파의 파수이다.

식 (27)의 임계주파수 아래에서 음압 파동해를 식 (3)의 오일러공식을 이용하면 다음과 같이 표현할 수 있다.

여기서 v는 음향입자속도의 크기(m/s), ρ₀는 음향매질의 밀도이다.

식 (26)에 음장에 의한 유체하중 식 (28)를 적용하면 다음과 같다.

식 (29)에서 f=0과 편측 유체장을 양측 유체장으로 확대하고 두 유체장은 서로 독립적(uncorrelated)이라고 가정하면 식 (30)과 같은 분산관계(dispersion relation)을 얻을 수 있다.

여기서 k₁, k₂는 각각 유체음장 1, 2에서 음파의 파수를 의미하고 ρ₁, ρ₂는 각 유체의 밀도를 나타낸다.

식 (30)을 만족하는 κ를 얻으면 밀도가 높은 유체장에 접수된 평판의 보정 굽힘 파수(corrected flexural wavenumber, κ_p,c)를 얻을 수 있고 이를 이용해 굽힘파의 보정 에너지전달속도(c_gf,c)를 유도하면 식 (31)과 같다.

유도된 보정 굽힘파수와 보정 굽힘 에너지전달속도의 검증을 위한 수치해석을 양측에 해수가 접수되어 있을 경우와 해수와 공기가 각각 접수해 있는 0.01 m의 두께를 가진 steel 평판에 대해 수행하였다. 뉴턴-랩슨법(Newton-Raphson method)을 이용해 분산관계 식 (30)의 근을 구해 보정 굽힘 파수를 계산하고 식 (31)에 대입해 보정 에너지전달속도를 예측하였다. Fig. 5는 평판 양쪽 모두 해수가 접한 경우와 공기와 해수가 접한 경우에 대해 보정 굽힘 파수와 보정 에너지전달속도를 보여준다.

Fig. 5 
Fluid loading effects in flexural wavenumber and group velocity of immersed plate.

일반적으로 접수로 인한 부가질량효과는 음향 방사가 잘 발생하지 않는 임계주파수 아래에서 극대화 되는데, 이에 따라 해당모델의 해수의 임계주파수 23 kHz 아래에서는 질량 증가효과 극대화로 보정 굽힘 파수가 진공 상태의 굽힘 파수 대비 커지고 음향 방사가 잘 발생하는 임계주파수 위 주파수 대역에서는 부가질량효과가 미미해져 진공상태의 굽힘 파수와 거의 일치하게 된다. Fig. 5(a)와 5(c)에서 모두 밀도가 높은 해수 접수효과가 반영된 보정 굽힘 파수는 임계주파수보다 높은 주파수영역에서는 진공 상태의 굽힘 파수와 같은 값을 갖지만 임계주파수보다 낮은 주파수대역에서는 접수 효과가 커서 파수의 차이가 있음을 알 수 있다. Fig. 5(b)와 5(d)에서의 보정 에너지전달속도에서도 보정 굽힘 파수와 동일한 현상이 나타난다. 또한 공기의 접수효과는 임계주파수(1.2 kHz)가 상대적으로 해수보다 낮고 무거운 유체가 더 크기 때문에 공기-해수가 접한 평판보다 해수-해수가 접한 평판의 보정된 굽힘 파수와 에너지전달속도가 진공상태의 값과 차이가 큰 것을 확인할 수 있다. Gu and Fuller(1992)는 임계주파수 아래에서 부가질량효과와 감쇠 효과로 인해 해수에 접한 평판의 고유진동수와 응답이 진공 상태의 값보다 고유진동수는 낮아지고 응답도 작아지는 결과를 도출하였다. 이런 결과로 유추해 볼 때, 진공상태의 굽힘 파수(κ=(ρ_shω²/D_s )^0.25)가 질량이 증가할수록 커지기 때문에 본 연구에서 얻어진 해수에 접한 보정 굽힘 파수가 진공 상태의 굽힘 파수보다 커지는 경향과 해수에 접한 보정 에너지전달속도가 진공 상태의 값보다 작아지는 경향과도 일치함을 알 수 있다.

3차원 유체 공간의 점음원(acoustic point source)이 있을 경우 평면파에 대한 에너지지배방정식은 다음과 같다 (Cho, 1993).

여기서, ∇²는 3차원 공간에 대한 라플라시안(Laplacian)이며, c_ga는 음향 파동의 에너지전달속도(group speed, m/s), η는 음향공간의 내부손실계수(internal loss factor), ω는 가진 중심주파수(rad/s), <e_a>는 음향에너지밀도(J/m³), Π_in,a는 입력 음향파워, (x₀,y₀,z₀)은 점음원의 위치를 의미한다.

일반적으로 에너지흐름유한요소법의 유한요소정식화를 위해 식 (33)과 같은 변분식(variational formulation)을 얻고 Galerkin 근사화를 이용한다 (Park, 2019).

여기서 벡터 n은 영역 V의 경계면 S에 바깥쪽으로의 수직 단위벡터를 나타내고, ψ는 시험함수(trial function)로 Lagrange 형상함수(shape function)를 사용하여 근사화한다.

본 연구에서는 기존 연구에서 직육면체 형상의 음향공간만 유한요소해석 (Kim, 2007)을 할 수 있었던 한계를 극복하기 위해 임의의 형상의 실내음장을 에너지흐름유한요소해석할 수 있도록, 사면체(tetrahedron)요소를 정식화하였다 (Yoon, 2020). Fig. 7의 3차원 육면체, 사면체 음향요소를 위한 마스터요소(master element)와 형상함수(shape function)는 각각 식 (34), (35)와 같다.

여기서 N과 L_i는 각각 육면체요소와 사면체요소의 1차 형상함수를 나타내며, 사면체요소는 면적 좌표계(area coordinates)를 사용하고 V_total과 V_pijk는 각각 사면체 전체면적과 내부 점 P와 꼭지점 i,j,k로 이루어진 사면체의 면적을 의미한다.

3차원 에너지흐름유한요소 정식화를 검증하기 위해 길이 1 m의 정육면체 음향공간의 중심((x₀,y₀,z₀) = (0m,0m,0m)에 3000 Hz의 1 W의 음향 파워 음원(acoustic power source)를 입력하고 모든 외곽 면에서 음향 파워 0이라는 경계조건을 부여하였다. 유한요소해석결과를 검증하기 위한 엄밀해(exact solution)는 식 (32)의 푸리에급수해(Fourier series solution)로 얻어진다 (Yoon, 2020).

여기서 계수 E_mnl는 다음과 같이 계산할 수 있다.

Fig. 6에서 보여진 대로 육면체요소(27000개)와 사면체요소(135000개)를 사용한 에너지흐름유한요소해석결과와 푸리에급수해(60개 모드)와 에너지밀도 레벨과 공간 분포가 잘 일치한다.

Fig. 6 
Validation of energy flow finite element analysis using tetrahedron element (f=3000 Hz)

Fig. 7 
Master elements for hexahedron and tetrahedron

3차원 실내 음장과 2차원 평판구조물로 이뤄진 복합시스템의 유체/구조 연성을 고려한 에너지흐름유한요소해석을 위해 식 (33)을 확장하여 표현할 수 있다.

여기서 <e>는 에너지밀도, ψ는 시험함수(trial function)이고 Ω_e는 평판구조물에서는 진동응답의 2차원 영역, 음향공간에서는 소음 응답의 3차원 영역을 의미한다. 또한 Γ_e는 평판구조물 사이의 1차원 선경계(line junction)와 음향요소와 구조요소 사이의 2차원 면경계(area junction)을 나타낸다.

식 (38)을 행렬 형태로 표현하면 다음과 같은 식으로 나타낼 수 있다 (Yoon, 2020).

식 (39)에서 아래 첨자 f,l,s,a는 각각 평판의 굽힘파, 면내 종파, 면내 전단파, 음파를 의미하고 [J]_S-S는 평판의 선경계에서 연결요소행렬, [J]_S-A는 평판 요소와 음향 요소의 면경계에서 연결요소행렬, [J]_A-A는 음향 요소와 음향 요소 사이의 면경계에서 연결요소행렬을 의미하며, 각 행렬식은 식 (40)~(42)로 얻어진다.

여기서 i와 j는 에너지 자유도(degree of freedom)를 의미하고 Ki,je은 요소의 강성행렬, Ji,je은 연결요소행렬(joint element matrix)이며 N은 형상함수 행렬, [I]는 단위행렬, [P]는 파워투과반사계수 행렬, [C]는 에너지전달속도 행렬이다.

Fig. 8은 두 개의 방으로 이루어진 직육면체 형상의 해수에 잠긴 몰수체를 보여준다. 2개의 음향공간(1 m*1 m*1 m)과 0.02 m의 두께를 가진 steel 평판(1 m*1 m) 11개로 구성된다. 평판의 이력감쇠계수(hysteresis damping)는 0.001이고 가진 굽힘 파워는 (x₀,y₀,z₀)=(0.3 m, 0.2 m, 0 m)에 10 W로 가정하였다.

Fig. 8 
Complicated plate structure fully immersed in sea-water or air

Fig. 9는 공기 중에 놓여있는 상태의 진동 소음 에너지밀도를 보여준다. 복합평판구조물의 굽힘 진동에너지는 가진 평판에서 가장 큰 값을 갖고 경계에서 파동 변환에 따라 감소하는 경향을 잘 보여주며, 평판 면내 진동에너지의 경우 가진 평판과 직각으로 연결된 곳에서 가장 큰 레벨을 갖고 굽힘 진동에너지와 마찬가지로 거리에 따라 감소하는 경향을 잘 보여준다.

Fig. 9 
Energy flow finite element analysis of complicated plate structure in air(2000Hz)

음향 에너지밀도의 경우는 고체음으로 인해 가진 평판이 있는 음향공간에서 레벨이 높고 중간 격벽으로 인한 투과손실(transmission loss)로 인접 음향공간의 음향 에너지밀도는 저하되는 것을 확인할 수 있다. Fig. 9(d)는 음향 인텐시티 레벨의 분포로서 고체음의 특성대로 외벽에서 중심쪽으로 음향파워가 전파되는 것을 보여준다.

Fig. 10은 해수에 잠긴 구조물의 진동 소음 에너지밀도의 결과이다. 평판이 해수에 잠기게 되면 Fig. 5에서의 결과와 같이 임계주파수 이내에서는 접수효과로 인해 굽힘 파수가 커지게 되고 에너지전달속도는 작아지게 된다. 따라서 이는 공간적인 감쇠 요소(ηω/c_g)가 커지는 효과로 나타난다. 또한 방사효율의 효과보다 해수의 임피던스(ρ₀c₀=1.53X10⁶ Rayls)가 공기에 비해 약 3700배 크기 때문에 방사손실계수(ρ₀c₀σ_rad/ωρ_ph_p)도 공기에 비해 크다. 이는 Fig. 4(b)의 평판에서 유체로의 파워전달계수에서도 확인할 수 있다. 따라서 해수에 접한 평판은 공기 중의 평판에 비해 에너지감쇠가 크게 나타난다. Fig. 10(a)에서와 같이 진동파워는 가진점으로부터 공급되지만 굽힘 에너지밀도는 유일하게 해수와 접하지 않고 공기와 접한 중심 격벽에서 가장 높은 것을 확인할 수 있다.

Fig. 10 
Energy flow finite element analysis of complicated plate structure immersed in seawater (2000 Hz)

Fig. 11은 평판과 음향공간의 크기는 Fig. 8에서의 모델과 동일하지만 전체 구조물의 외부가 모두 해수로 둘러싸여 있지 않고 측면 중심부 0.5m까지 높이까지만 해수에 잠겨있다고 가정한다.

Fig. 11 
Complicated plate structure partially immersed in seawater

Fig. 12(e)에서 굽힘 진동 인텐시티가 가진점에서 주변으로 퍼져나감에도 불구하고 Fig. 12(a)를 보면 해수에 절반 잠긴 구조물의 경우 해수에 접하지 않은 천장 판넬(ceiling panel)의 굽힘 에너지밀도가 가장 높은 것을 확인할 수 있다. 이는 앞에서 확인한 대로 해수 접수로 인한 판넬 내부감쇠(internal damping)가 높아졌기 때문이다.

Fig. 12 
Energy flow finite element analysis of complicated plate structure partially immersed in seawater (2000 Hz)

Fig. 12(b)는 면내 진동 에너지밀도의 분포를 나타내는데 접수효과로 인해서는 면내 진동에너지 감쇠에는 영향을 미치지 않기 때문에 경계에서 발생하는 면내 진동에너지 특성으로 인해 굽힘 가진 평판에서 높은 면내 에너지밀도를 보인다.

Fig. 12(d)는 전면부 판넬(1 m≦x≦2 m, y=0 m, 0 m≦z≦1 m)의 굽힘 진동 에너지밀도의 상하 분포를 나타낸 그림이다. 굽힘 가진점이 바닥 평판임에도 불구하고 바닥으로부터 높이 0.5 m까지 접수로 인해 감쇠가 커진 효과를 확인할 수 있다.