I. 서 론

음향 포물선 방정식 법은 해양에서 단일 음원에서 방사된 중저주파 음파 전달을 모의하는데 가장 많이 사용되는 방법이다.^[1] 포물선 방정식법으로 구현된 모델로 전 세계적으로 RAnge-dependent Model(RAM)이 널리 알려져 있다.^[2],[3]

수중음파전달에서는 먼 거리에서 관찰한 값이 중요한 경우가 많기 때문에, 단일 음원을 사용하여 모델링을 해도 대부분은 소나방정식을 통해 음파전달 예측을 하는데 크게 불편함은 없다. 하지만 센서 배열의 성능을 시험하거나, 근거리 음파전달, 음파 굴절에 의한 영향 등 을 모의하고 싶을 때는 방향성이 있는 음원에 대한 포물선 방정식 모델이 필요할 것이다.

본 연구에서는 방향성이 있는 음원을 모의하기 위한 두 가지 방법을 제안한다. 첫 번째 방법은 전 방향성 음원을 모사하는 함수인 Delta 함수를 필터링 하는 것이다. 필터링 하는 구간에 따라 구현하고자하는 방향성 음원의 특성이 달라질 것이다. 본 연구에서는 수학적으로 방향성 음원의 빔 크기와 필터링 간격 간의 관계를 보였다.

두 번째는 단일 음원의 self-starter 해^[4]에 필터를 적용하는 것이다. 필터를 통해서 self-starter 해의 원하는 부분의 파수 성분만 남겨놓는 것이다. 참고로 본 연구와 연구 목적은 다르지만 위와 같은 방법은 이미 1997년에 포물선 방정식 self-starter 해를 안정화시키기 위한 방법으로 Cederberg et al.^[5],[6]에 의해 제안되었다. 그렇지만 방향성 있는 음원 연구에 적용한 것은 본 연구가 처음이다.

II장에서는 포물선 방정식법 및 제안된 방법에 대한 이론적인 설명을 기술하였으며, III장에서는 Pekeris 환경 및 심해 환경에 대한 수치실험 결과를 기술했다. 마지막으로 V장은 결론이다.

II. 이 론

2.1 Split-Step Padè 포물선 방정식 모델

2차원 포물선 방정식의 Split-step 해는 다음과 같이 표현된다.^[2]

여기서 $r$은 수평거리를 나타내며, $k_0$는 참조 파수이고, $∆r$은 거리방향의 증분이다. X 연산자는 다음과 같다.

위 식에서 $k_0$는 참조 파수이고, $\rho$는 매질의 밀도, $k$는 매질의 파수이다. $\alpha$는 에너지 보존을 만족시키기 위한 계수이다. $\partial /\partial z$는 깊이방향에 대한 편미분을 나타낸다.

Eq. (1)의 포물선 방정식을 풀기 위해서는 초기해가 필요하다. 일반적으로 포물선 방정식에서는 초기 해는 self-starter를 이용하여 계산을 한다. Collins에 의해 2차원 원통형 좌표계에서 제안된 초기해는 다음과 같다.^[5],[6]

위 식에서 $\delta$(-)는 Delta 함수를 의미하며, $z_s$는 음원이 위치한 수심이다. 참고로 유한차분법에서 수치적으로 Delta 함수를 적용할 때, 음원이 위치한 노드에서는 1/$∆z$을 놓고 그 외에는 0을 놓는다. 여기서 $∆z$는 깊이 격자의 간격이다.

2.2 Delta 함수의 필터링

수학적으로 Delta 함수는 푸리에 변환을 이용하여 다음과 같이 표현할 수 있다.

Eq. (4)의 우측의 스펙트럼에서 $[k_l,k_u]$의 구간에 해당되는 파수 성분만 뽑아낸다고 하면, 근사된 Delta 함수 ${\delta }_k$는 아래와 같이 정리할 수 있다.

여기서 $∆=(k_u-k_l)/2$, $k_c=(k_u+k_l)/2$이다.

Eq. (5)를 보면, 밴드패스 필터링 된 Delta 함수는 해석적으로 sinc 함수의 형태로 나타나는 것을 확인할 수 있다. 위 식에서 k는 물리적으로 음파의 수직 파수의 의미를 갖는다. 음파의 전파 각도를 스침각 $\theta$라고 할 때, 수직 파수는 각각 $k_l=ksin{\theta }_l$, $k_u=ksin{\theta }_u$로 치환을 할 수 있다. 이때 ${∆}_{\theta }={\theta }_u-{\theta }_l$이 방향성 있는 음원의 빔 폭이 된다.

이론상으로는 Eq. (5)에서 $k_l$과 $k_u$의 값을 조정하면 가느다란 빔을 만들 수 있으니, 현실적으로 수치 영역이 유한하기 때문에 한계가 있다. 보다 자세히 설명하면, Eq. (5)의 sinc 함수의 주엽의 null-to-null 크기는 다음과 같이 정리할 수 있다.

여기서 $c_0$는 음원 위치에서의 수중 음속이며, $f$는 주파수이다.

만약에 10도의 빔 폭을 갖는 음원을 만들려고 하고 수중음속이 1500 m/s라 할 때, ${∆}_z\approx 9000m/f$가 된다. 100 Hz의 음원을 사용한다고 할 때, 주엽의 크기는 90 m가 된다. Eq. (5)의 sinc 함수에서 주요 에너지 성분을 담을 수 있는 영역은 경험적으로 주엽 크기의 최소 5배 이상이어야 한다. 그렇기 때문에 위의 조건에서 수심은 최소 450 m 이상이어야 한다. 실제 한국 해 근처의 해역은 위 조건을 항상 만족하지 못한다. 그렇기 때문에 위의 필터링 된 Delta 함수는 고주파수의 음원(1 kHz 이상)에서 주로 사용가능할 것이다.

필터링 된 Delta 함수를 이용한 방향성 음원의 self-starter는 Eq. (3)의 포물선 방정식의 self-starter의 $\delta (z-z_s)$ 항을 Eq. (5)로 대체하는 것으로 쉽게 제작할 수 있다.

2.3 유리 함수를 이용한 self-starter 필터링

두 번째로 제안하는 방법은 첫 번째 방법과 달리 Eq. (3)의 self-starter 해에 유리 함수 필터를 직접 적용하는 방법이다.

주 내용을 기술하기에 앞서 Eq. (2)의 물리적인 의미를 고찰한다. 등 음속과 등 밀도의 해양환경을 가정할 때, Eq. (2)는 $x=(1/k_0^2)[{\partial }^2/\partial z^2+(k^2-k_0^2\left)\right]$로 단순화 된다. 평면파를 가정할 때, ${\partial }^2/\partial z^2=(jk\sin \theta {)}^2$이 되며, $X=(k^2/k_0^2){\cos }^2\theta -1$로 정리할 수 있다. 이 식에서 ‘–1’에 해당하는 성분만 제외하면, X 연산자는 결국 음파의 수평파수의 제곱이라는 것을 짐작할 수 있다. 그렇기 때문에 X 연산자를 필터링하는 임의의 함수 $R\left(X\right)$를 설계할 수 있다면, 수평파수에 대한 필터링이 가능할 것이다.

Collins는 안정된 해를 구하기 위해 Eq. (3)을 다음과 같이 2단계로 나누어서 계산을 했다.

우리는 Eq. (7)의 두 번째 줄에 아래와 같이 필터링 함수 $R\left(X\right)$를 도입했다.

즉, $R\left(X\right)$는 $p_{tmp}\left(X\right)$를 필터링 해주는 역할을 하는 임의의 필터이다. 본 연구에서는 Cederberg et al.^[5],[6]가 제안한 사각 필터를 사용했다.

여기서 $X_l$와 $X_u$는 필터에서 통과시키는 대역의 하한값과 상한값을 말한다.

Eq. (9)를 포물선 방정식 모델에 적용하기 위해서 아래와 같이 $R\left(X\right)$를 유리함수로 근사해야 한다.

위의 식에서 $a_i$와 $b_i$는 최소자승법을 적용하여 계산하였다.

Eq. (10)은 유리 함수의 형태이기 때문에 Eq. (8)의 포물선 방정식 알고리즘에 쉽게 적용할 수 있다. 하지만 Eq. (9)에서 대역 통과 구간의 크기가 작을수록 Eq. (10)에서 유리함수 근사의 차수가 높아질 것이다. 유리 함수 차수가 높아지면 Eq. (8)의 계산시간이 증가하나 보통 1 s 이내이다.

또한, Eq. (9)의 필터링 함수는 물리적으로 ‘수평파수의 제곱(X)’을 필터링하도록 설계 된 것이다. 만약에 음원에서 아랫방향으로 $\theta$의 각도로 전파하는 음파와 윗 방향으로 $-\theta$로 전파하는 음파가 있을 때, 이 두 개의 음파의 X값은 서로 동일할 것이다. 그렇기 때문에, 위의 필터링 방법으로는 깊이 방향에 대해 위로 가는 음파와 아래로 가는 음파를 구분할 수는 없다.

추가적으로 본 연구에서는 사각 필터를 유리함수를 이용한 self-starter에 적용하였으나, 다른 여러 형태의 필터를 적용하는 것도 가능하다. 예를 들어 가우시안 형태의 필터를 적용한다면, Eq. (9)는 아래와 같이 바꿀 수 있다.

여기서 $X_a=(X_u-X_l)/2$이고, $\sigma$는 사용자에 의해 결정된다. 하한값과 상한값을 위의 필터가 0.5 이하의 크기 손실이 발생하는 영역의 경계값으로 정의하면, $\sigma =(X_u-X_l)/\left(2\sqrt{2\ln 2}\right)$이 된다.

III. 수치 실험

본 논문에서 제안된 방법을 확인하기 위해 크게 천해와 심해 환경에 대해 수치실험을 수행했다.

천해 환경은 Pekeris 도파관으로 가정하였으며, 심해 환경은 동해 환경을 가정했다. Fig. 1은 계산에 사용된 해양환경 및 물성 값을 보여준다.

Fig. 1.

Ocean environments for numerical experiments. (a) Pekeris waveguide; (b) deep environment.

3.1 필터링 된 Delta 함수 예제

방향성 있는 음원을 구현하기 위해 ${\theta }_l=7.5°$, ${\theta }_u=12.5°$로 설정했다. 이때 빔 폭은 5도가 나오며 빔의 중심각 ${\theta }_0=5°$이다. 음원의 주파수는 1 kHz로 설정했으며 전체 수치 영역은 깊이 400 m이고, 깊이 방향으로 0.5 m의 간격으로 분할했다. Padè계수는 8차를 사용했다. Fig. 2는 음원의 깊이가 100 m일 때 필터링 된 Delta 함수를 사용한 경우와 사용하지 않은 경우를 보여준다.

Fig. 2.

Sound transmission loss in the Perkeris waveguide with (a) an omni-directional source of Delta function and (b) a directional source obtained by Eq. (5) at the frequency of 1 kHz. The unit of colorbar is dB.

Fig. 2(a)는 필터링을 하지 않은 경우의 음파전달손실을 보여준다. 음파가 전 방향으로 다중반사를 하며 진행하는 것을 알 수 있다. Fig. 2(b)는 위에서 서술한 빔 폭에 대해 방향성을 준 음원의 결과이다. 명확하게 음원의 방향성이 구현되는 것을 확인할 수 있다. 다만 수치영역에 부엽이 발생한다.

이것은 수치 영역이 유한하기 때문에 발생하는 Aliasing으로 볼 수 있다.

Fig. 3은 동일한 환경에서 주파수만 100 Hz로 낮추었을 때의 결과를 보여준다. Delta 함수의 유효한 영역에 비해 도파관의 깊이가 작기 때문에 정확하게 방향성 있는 음원이 구현되지 못하고 뭉그러지는 것을 확인할 수 있다.

Fig. 3.

Sound transmission loss in the Perkeris waveguide with a directional source obtained by Eq. (5) at the frequency of 100 Hz. The unit of colorbar is dB.

Fig. 4는 심해 환경에서 빔의 조향각에 따른 전달 손실의 결과를 보여준다. 빔 폭은 5도로 동일하며 음원의 주파수는 2 kHz를 사용했다. Fig. 4(a)는 비교를 위해 전 방향의 음원을 사용했을 때의 전달손실의 결과를 나타낸다. 빔의 중심각을 ${\theta }_0$라고 할 때 Fig. 4(b)는 ${\theta }_0=3°$, Fig. 4(c)는 ${\theta }_0=10°$로 설정했다. 여기서 +값은 음파의 방향이 아랫방향인 것을 말한다.

Fig. 4.

Sound transmission loss in the deep sea with (a) an omni-directional source, (b) a directional source with the center angle of ${\theta }_0=3.0°$, and (c) a directional source with the center angle of ${\theta }_0=10.0°$ at the frequency of 2 kHz. The directional sources are calculated with the filtered Delta function of Eq. (5). The unit of colorbar is dB.

Fig. 4(b)와 4(c)를 살펴보면, 심해 환경에서 나타나는 굴절파와 해저면 반사파가 Delta 함수의 필터링을 통해 명확하게 분리가 되는 것을 확인할 수 있다.

3.2 유리함수 필터링 예제

본 절에서는 유리함수 필터링을 이용한 방향성 있는 음원에 대한 수치실험을 진행했다. 우선 앞 절의 심해 환경에 대해 각각 ${\theta }_0=2.5°$, ${\theta }_0=17.5°$이고, 빔 폭은 5도인 방향성 있는 음원을 생각했다. 유리 함수 근사의 차수는 6차를 사용했다.

Fig. 5는 각각의 결과를 보여준다. Fig. 5(a)와 5(b)는 각각 ${\theta }_0=2.5°$이고 ${\theta }_0=17.5°$일 때이다 (실제는 유리함수 필터링의 특성에 따라 ${\theta }_0=\pm 2.5°$와 ${\theta }_0=\pm 17.5°$의 빔이 발생한다). 해당 영역에 대해서 필터링이 되는 것을 볼 수 있다. 특히 Eq. (5)의 필터링 된 Delta 함수를 사용하는 것에 비해 유리 함수 필터링이 보다 더 깨끗한 결과를 얻는 것을 확인할 수 있다. 이는 필터링 된 Delta 함수는 수치 영역의 유한성과 불연속성으로 인해 정확하게 구현하는 것에 어려움이 있으나, 유리 함수 필터를 사용한 방법은 음장의 파수를 필터링 시켜주는 방법이므로 보다 정확하다.

Fig. 5.

Sound transmission loss in the deep sea with directional sources (a) with the center angle of ${\theta }_0=2.5°$ and (b) with the center angle of ${\theta }_0=17.5°$ at the frequency of 2 kHz. The directional sources are calculated by the rational filtering of the Delta function of Eq. (8). The unit of colorbar is dB.

한편, Fig. 5(b)의 5 km 이내에서의 음장을 관찰하면 약하게 위 방향의 음파와 아랫방향의 음파가 동시에 발생하는 것을 확인할 수 있다.

이를 명확하게 보이기 위해 r = 10 m에서 유리함수 필터를 이용한 음장에 대해 수직파수 분석을 수행했다. 수직 파수 분석 결과는 r = 10 m에서 수직방향의 음장을 푸리에 변환하여 얻을 수 있다. Fig. 6은 그 결과이다. 전 방향의 음원을 사용하는 경우에 대략 $-1<k_z/k_0<1$에서 분포하는 음장이 유리함수 필터링에 의해 각각 필터링 되는 것을 명확하게 확인할 수 있다.

Fig. 6.

Spectrum of the pressure field at the range of 10 m. Zmax indicates the maximum depth used for the spectum analysis.

IV. 결 론

본 연구에서는 포물선 방정식 모델에 방향성 있는 음원을 구현할 수 있는 기법을 개발했다. 첫 번째 기법은 Delta 함수를 수학적으로 필터링하는 방법으로 간단하게 sinc 함수 형태로 얻어진다. 이 방법은 포물선 방정식 모델의 입력부에 구현하기는 매우 단순하나 실제 수치영역은 유한하며, 깊이 방향에 대해 불연속성이 존재하므로 실제 도파관 환경에서 구현을 해보면 부엽이 나타나는 것을 확인할 수 있다. 이것은 실제 수중에서 수직 선배열을 이용하여 방향성 있는 음장을 구현할 때 발생하는 문제들과 매우 유사하다.

두 번째 방법은 유리함수 필터링을 이용한 방법이다. 유리함수 필터는 포물선 방정식에 나타나는 연산자인 X 연산자를 필터링 하는 방법으로, 본 연구에서는 방향성 있는 음원을 만들기 위해 위 방법을 적용하였다. 위의 방법은 음장 자체를 필터링하는 것이기 때문에 보다 명료한 음장이 얻어진다는 장점이 있으나, 깊이방향에 대해서 위로 진행하는 음파와 아래로 진행하는 음파를 구분하지 못하는 단점이 존재한다. 위의 두 방법을 서로 조합하면 보다 정확도 높은 필터링이 가능하리라 예상되나 응용적인 연구에 속하므로 본 연구에서는 다루지 않을 것이다.