Fig. 1은 고유치 해석 대상 임의 형상 음향 공동을 보여준다. 일반적인 고유치 해석 이론에서 음향 공동의 내부 음압은 식 (1)과 같은 헤름홀츠 방정식(Helmholtz equation) 형태의 지배방정식을 만족한다⁽²⁰⁾.

Fig. 1 
Arbitrarily shaped acoustic cavity with rigid-wall boundary Γ, which is discretized with boundary nodes P₁, P₂, ... P_N

여기서 k는 파수(wavenumber)이며 각주파수 ω와 음속 c로 구성된 식 k = ω/c에 의해 구해진다. 그리고 p(r)은 Fig. 1에서 음향 공동 내부의 한 점 P에서의 음압(sound pressure)을 나타내며, r은 점 P에 대한 위치 벡터를 의미한다.

음향 공동의 경계 조건은 강체벽(rigid-wall) 경계조건과 개방(pressure-released) 경계 조건이 존재하며, 이 논문에서는 가장 일반적인 강체벽 경계 조건만을 고려한다. 강체벽 경계 조건은 식 (2)와 같이 1차 미분 방정식의 형태를 가진다.

여기서 p(r_Γ)은 음향 공동 경계 Γ위의 한 점에서의 음압을 나타내며, r_Γ는 Fig. 1에서와 같이 음향 공동 경계 상의 한 점에 대한 위치 벡터이며 n은 경계 상의 한 점에서의 법선 방향을 의미한다.

NDIF법을 적용하기 위하여 Fig. 1과 같이 음향 공동의 경계 Γ는 N개의 노드들 P₁, P₂, ... P_N 으로 이산화 된다. 다음으로, s번째 경계 노드 P_s에 대한 무차원 동영향 함수를 식 (3)과 같이 정의한다.

여기서 J₀(⋯)는 제1종 0차 베셀 함수를 나타내며, 벡터 r_s는 경계 Γ에 위치한 s번째 노드에 대한 위치 벡터를 의미한다. 식 (3)의 무차원 동영향 함수는 물리적으로 노드 P_s에서 단위 음압이 ω의 주파수로 발생되었을 때 음향 공동의 내부 점 P에서의 음압의 크기를 의미한다. 식 (3)은 지배방정식 식 (1)을 정확히 만족하며 음향 공동의 음압을 구하기 위한 기저 함수로 사용된다.

음향 공동 내부 음압 p(r)은 경계 노드에서 정의된 무차원 동영향 함수 식 (3)의 선형 결합으로 식 (4)와 같이 근사화 된다.

여기서 A_s는 s번째 노드에서 정의된 무차원 동영향 함수가 내부 음압에 미치는 기여도를 나타내는 계수이다. 참고로, 식 (4) 또한 지배방정식 식 (1)을 만족한다.

음향 공동의 경계 Γ에 주어진 강체벽 경계 조건 식 (2)는 N개의 경계 노드 P₁, P₂, ... P_N에 대해 식 (5)와 같이 이산화 된다.

여기서 r_i와 n_i는 Fig. 1과 같이 경계 노드 P_i에 대한 위치 벡터와 법선 벡터를 각각 의미한다.

근사화된 음압 식 (4)를 이산화된 경계조건 식 (5)에 대입하면 식 (6)을 얻을 수 있다.

식 (6)에서 법선 방향 n_i에 대한 미분을 수행하면, 식 (7)이 얻어진다.

여기서 J₁(⋯)은 제1종 1차 베셀 함수를 나타낸다.

식 (7)을 행렬식의 형태로 정리하면 식 (8)과 같이 된다.

여기서 SM(k)는 파수 k를 변수로 가지는 시스템 행렬이고, A는 기여도 벡터를 의미하며 벡터 A의 s번째 성분은 A_s이다. 시스템 행렬 SM(k)의 성분은 식 (9)에 의해 주어진다.

음향 공동의 고유치는 식 (8)에서 주어진 시스템 행렬 SM(k)의 판별식을 계산함으로써 추출될 수 있다. 음향 공동의 고유치들은 식 (10)을 만족하는 파수 값들과 일치한다.

하지만, 식 (10)에 의해 구해진 고유치들 속에는 음향 공동고유치 이외의 허위 고유치들(spurious eigenvalues)이 포함되어 있다. 이들 허위 고유치들은 음향 공동과 동일한 형상을 가진 멤브레인의 고유치들에 해당됨이 과거 연구^(2,3)에서 밝혀졌다. 허위 고유치가 발생하는 이유는 시스템 행렬 SM(k)이 멤브레인에 대한 시스템 행렬 SM_mem(k)와 식 (11)과 같이 종속되어 있기 때문이다.

여기서 행렬 SM_net(k)은 허위 고유치를 제외한 음향 공동 고유치만을 제공하는 시스템 행렬을 의미한다. 식 (11)의 양변에 대한 판별식 값을 구하면 식 (12)와 같이 된다.

식 (12)에서 det(SM_net(k)) = 0을 만족하는 파수 값들이 음향 공동의 고유치와 정확히 일치하며, 식 (12)를 이용하면 det(SM_net(k)) = 0은 식 (13)으로 대체 가능하다.

식 (13)으로 부터 허위고유치를 제외한 음향 공동의 고유치만을 구할 수 있다. 아울러, 음향 공동의 j번째 고유모드 형상은 j번째 고유치를 식 (8)에 대입하여 j번째 기여도 벡터를 구하고, 이 기여도 벡터의 성분 값들과 j번째 고유치를 식 (4)에 대입하는 방법에 의해 구해진다.

이전 연구에서는^(18,19) 음향 공동의 노드의 수를 증가시킬 때 발생하는 시스템 행렬의 저주파수 영역 발산 문제를 해결하기 위해 식 (13) 대신 식 (14)를 사용하여 시스템 행렬의 판별식을 구하는 방안을 제안하였다.

여기서 λic와 λim은 각각 시스템 행렬 SM과 SM_(mem)의 i번째 고유치, S는 시스템 행렬의 크기, k_b는 SM이 풀 랭크 값을 가지기 시작하는 파수를 의미하며, R_c(k)와 R_m(k)는 파수 k에 따라 변하는 SM과 SM_(mem)의 랭크 함수(rank function)를 각각 의미한다.

식 (14)를 이용하여 음향 공동의 고유치를 추출할 경우 저주파수 고유치는 추출이 되나, 관심 주파수 범위에 따라서 저주파수 고유치들의 정밀도가 변하는 문제점이 발생한다. 이를 개선하기 위해 이 논문에서는 식 (14)를 발전시켜 식 (15)를 제안한다.

식 (15)의 추가적인 장점은 파수 값 k_b를 기준으로 양쪽 구간에 대해 두개의 식으로 분리하여 판별식을 계산해야하는 번거로움을 극복한 것이다.

이 연구에서 제안된 식 (15)를 검증하기 위해, Fig. 2와 같은 임의 사각 형상 음향 공동이 고려되어진다. Fig. 2(a)와 Fig. 2(b)는 음향 공동의 경계가 16개와 24개의 노드로 이산화된 경우를 각각 보여준다.

Fig. 2 
Quadrilateral acoustic cavity of arbitrary shape, which is discretized with 16 and 24 nodes, respectively

관심 파수 범위의 최소값 k_min을 0으로 설정한 후, 기본 NDIF법에서 제안했던 식 (13)을 이용하여 판별식 곡선을 그리면 Fig. 3과 같다⁽²⁾. 두 개의 판별식 곡선 중에서, 16개 노드를 사용한 점선에서는 음향 공동의 6개의 고유치 C₁ ~ C₆가 성공적으로 구해지나, 24개의 노드를 사용한 실선에서는 저주파수 구간에서 시스템 행렬이 발산하여 저주파수 고유치 C₁ ~ C₃이 구해지지 않음을 확인할 수 있다. 참고로, 곡선에 표시된 M₂와 M₃는 음향공동과 같은 형상을 멤브레인의 두 번째와 세 번째 고유치를 의미한다.

Fig. 3 
Determinant curves for the quadrilateral acoustic cavity by the original NDIF method using Eq. (13) (dotted line: 16 nodes, solid line: 24 nodes)

위의 문제점을 해결하기 위해 최근의 연구에서 제안한 식 (14)를 이용하여 24개의 노드를 사용한 경우(Fig. 2(b))에 대한 판별식 곡선을 구하면 Fig. 4와 같다^(18,19). Fig. 4에서 확인할 수 있듯이, 관심 파수 범위의 최소값(k_min)이 0.0, 0.5, 1.0, 1.5로 변화됨에 따라 판별식 곡선의 개형이 변화되고, 저주파수 고유치(C₁ ~ C₃)들의 크기 변화가 발생함을 Table 1에서 확인할 수 있다. 특히, 관심 주파수 범위의 최소값 k_min을 0으로 한 경우에는 저주파수 고유치 C₁ ~ C₃가 구해지지 않음도 확인할 수 있다. 참고로, 고주파수 고유치 C₄ ~ C₆은 k_b보다 크기 때문에 k_min의 변화에 영향을 받지 않는다는 것을 식 (14)에서 확인할 수 있다.

Fig. 4 
Determinant curve for the quadrilateral acoustic cavity with 24 nodes by the recently developed NDIF method using Eq. (14) as k_min is changed

Fig. 4와 Table 1과 같이, 관심 주파수 범위가 변화됨에 따라 고유치 값의 변화가 발생하는 것은 식 (14)가 이론적으로 완성도가 떨어짐을 의미하므로, 이러한 문제점을 극복하고자 이 연구에서는 식 (15)가 제안되었다. 식 (15)를 이용하여 노드 24개를 사용한 임의 사각 형상 음향 공동에 대한 판별식 곡선을 그리면 Fig. 5와 같다. 추가적으로, 식 (15)를 계산할 때 사용되어지는 음향 공동과 멤브레인의 시스템 행렬 랭크 곡선은 Fig. 6에서 주어진다.

Fig. 5 
Determinant curve for the quadrilateral acoustic cavity with 24 nodes by the proposed NDIF method using Eq. (15) that is independent of k_min

Fig. 5의 저주파수 영역을 살펴 보면, 고유치 C₁ ~ C₃에 해당되는 골(trough)들이 곡선의 불연속 지점들과 뒤섞여 있음을 확인할 수 있다. 더구나 점선 타원으로 마크된 네 번째 고유치 C₄는 골(trough)이 아니라 마루(crest)의 형태로 나타남을 확인할 수 있다. 이는 Fig. 6에서 점선 타원으로 마크한 고유치 C₄에 해당하는 파수에서 랭크 R_c(k)가 불연속적으로 변하기 때문인 것으로 분석된다. 일반적으로 NDIF법에서 고유치는 판별식 곡선에서 골의 형태로 나타나는데, 위의 경우는 예외적인 상황으로 실제 고유치를 허위 고유치(spurious eigenvalue)로 오인할 상황으로 보인다. 그래서 Fig. 6의 랭크 값들의 불연속적인 변화때문에 발생하는 Fig. 5의 불연속 특성을 제거하기 위해 식 (16)이 제안된다.

Fig. 6 
Rank functions R_c(k) and R_m(k) of the quadrilateral acoustic cavity and membrane discretized with 24 nodes, respectively

여기서 ∏λaddck는 음향 공동 시스템 행렬(SM)의 랭크 변화에 의해 추가적으로 곱해지는 고유치들의 곱을 의미하며, ∏λaddmk는 멤브레인 시스템 행렬(SM_(mem))의 랭크 변화에 의해 추가적으로 곱해지는 고유치들의 곱을 의미한다.

식 (16)을 이용하여 구한 Fig. 7의 판별식 곡선을 살펴보면, Fig. 5에서 나타난 불연속성이 완전히 제거되었음을 확인할 수 있다. 이로 인해 고유치 C₄를 포함한 6개의 고유치 C₁ ~ C₆에 해당되는 골들이 명확히 나타남도 확인 가능하다. 이들 고유치들은 Table 2에 요약하였다. Table 2에서 961개의 많은 노드를 사용한 ANSYS결과와 단지 24개 노드를 사용한 제안된 방법의 결과를 비교해보면, 최대 오차가 1.16 %로 제안된 방법이 정확하고 타당한 결과를 제공한다고 말할 수 있다. 참고로, 여기서 고려된 임의 사각 형상 음향 공동의 경우 엄밀해가 존재하지 않기 때문에 제안된 방법과 ANSYS 중 어느 방법이 정확하다고 단정할 수는 없다.

Fig. 7 
Determinant curve for the quadrilateral acoustic cavity with 24 nodes by the proposed NDIF method using Eq. (16) to remove the discontinuity of the curve

이 연구에서 제안된 방법이 곡선 경계를 가진 임의 형상 음향 공동에서도 타당한지를 확인하기 위해 Fig. 8과 같은 음향 공동이 이 예제 검증에서 고려되어진다. Fig. 8(a)와 Fig. 8(b)는 음향 공동의 경계가 16개와 24개의 노드로 이산화된 경우를 각각 보여준다.

Fig. 8 
Arbitrarily shaped acoustic cavity discretized with 16 and 24 boundary nodes, respectively

관심 파수 범위의 최소값을 0으로 설정한 후, 기본 NDIF법에서⁽²⁾ 제안했던 식 (13)을 이용하여 판별식 곡선을 그리면 Fig. 9와 같다. 두 개의 판별식 곡선 중에서, 16개 노드를 사용한 점선에서는 음향 공동의 6개의 고유치 C₁ ~ C₆가 성공적으로 구해지나, 24개의 노드를 사용한 실선에서는 저주파수 구간에서 시스템 행렬이 발산하여 저주파수 고유치 C₁과 C₂가 구해지지 않음을 확인할 수 있다.

Fig. 9 
Determinant curves for the arbitrarily shaped acoustic cavity by the original NDIF method using Eq. (13) (dotted line: 16 nodes, solid line: 24 nodes)

위의 문제점을 해결하기 위해 최근의 연구^(18,19)에서 제안한 식 (14)를 이용하여 24개의 노드를 사용한 경우(Fig. 2(b))에 대해 판별식 곡선을 구하면 Fig. 10과 같다. Fig. 10에서 확인할 수 있듯이, 관심 파수 범위의 최소값(k_min)이 0.0, 0.5, 1.0, 1.5로 변화됨에 따라 판별식 곡선의 개형이 변화되고 저주파수 고유치(C₁과 C₂)들의 크기 변화가 발생함도 확인되었다. 특히, 관심 주파수 범위의 최소값 k_min을 0으로 한 경우에는 저주파수 고유치가 구해지지 않음도 확인할 수 있다.

Fig. 10 
Determinant curve for the arbitrarily shaped acoustic cavity with 24 nodes by the recently developed NDIF method using Eq. (14) as k_min is changed

위의 문제점을 극복하고자 이 연구에서 제안된 식 (15)를 이용하여 노드 24개를 사용한 임의 형상 음향 공동에 대한 판별식 곡선을 그리면 Fig. 11과 같다. 추가적으로, 식 (15)를 계산할 때 사용되어지는 임의 형상 음향 공동과 멤브레인의 시스템 행렬 랭크 곡선은 Fig. 12에서 주어진다.

Fig. 11 
Determinant curve for the arbitrarily shaped acoustic cavity with 24 nodes by the proposed NDIF method using Eq. (15) that is independent of k_min

Fig. 12 
Rank functions R_c(k) and R_m(k) of the arbitrarily shaped acoustic cavity and membrane discretized with 24 nodes, respectively

Fig. 11의 판별식 곡선에서 저주파수 영역을 살펴보면, 고유치 C₁ ~ C₃에 해당되는 골들이 곡선의 불연속 지점들과 뒤섞여 있음을 확인할 수 있다. 이러한 불연속 특성을 제거하기 위해 이 연구에서 제안된 식 (16)을 사용하여 판별식 곡선을 그리면 Fig. 13과 같다. Fig. 11에서 나타난 불연속성이 Fig. 13에서 완전히 제거되었음을 확인할 수 있으며, 고유치 C₁ ~ C₃을 포함한 6개의 고유치 C₁ ~ C₆에 해당되는 골(trough)들이 명확히 나타남도 확인 가능하다. 이들 고유치들이 요약된 Table 3을 살펴보면, 1321개의 많은 노드를 사용한 ANSYS 결과와 단지 24개 노드를 사용한 제안된 방법이 0.26 % 오차 이내로 거의 일치함이 확인된다.

Fig. 13 
Determinant curve for the arbitrarily shaped acoustic cavity with 24 nodes by the proposed NDIF method using Eq. (16) to remove the discontinuity of the curve