1. 개요 및 연구 배경

2022 개정 지구과학 II 해양학 단원(해수의 수온·염분·밀도 연직 구조와 음파 전달) 및 물리학 II 파동 단원(스넬의 굴절 법칙과 굴절률 구배에 따른 파동 궤적)을 융합하는 고도화된 음향 역학적 시뮬레이션입니다. 심해 해양에서는 표층 부근의 높은 수온으로 인해 음속이 높고, 수심이 깊어질수록 수온이 급격히 감소(수온약층)하여 음속이 극소값에 도달한 뒤, 수백 미터 이하의 심해역에서는 정수압(Hydrostatic Pressure) 증가로 인해 음속이 다시 증가합니다. 이 음속 극소축(Sound Speed Minimum Axis, 대략 수심 800~1200m)을 **SOFAR(Sound Fixing and Ranging) 채널** 또는 해양 음향 도파관(Acoustic Waveguide)이라 부릅니다.

SOFAR 채널 내부로 입사된 저주파 음파는 상·하부의 음속 구배로 인해 굴절되어 채널 중심축을 중심으로 끊임없이 반사 없이 굴절 진동하며 감쇠 없이 수천 킬로미터를 전파합니다. 해양 음향 단층촬영(Ocean Acoustic Tomography, OAT)은 이러한 음파의 도달시간(Travel Time)을 고정밀 계측하여 역산(Inversion)함으로써 광역 해양의 수온 삼차원 분포와 해류 순환, 내부파(Internal Waves)를 실시간 추정합니다.

본 연구에서는 대한민국 해군 차세대 유무인 복합 대잠전 체계(MUM-T ASW)의 무인잠수정(UUV) 및 수중글라이더(Underwater Glider)에 탑재되는 0-Alloc C# 10,000Hz SOFAR 채널 음선 추적(Ray Tracing) 및 도달시간 역산 서보 엔진을 설계하고 실시간 온디바이스 수학 모델을 구현합니다.

2. 수학적 및 물리적 거버닝 방정식

수심 \(z\) [m], 수온 \(T\) [°C], 염분 \(S\) [PSU]에 따른 해수 음속 \(c(z)\)는 경험적 Mackenzie 공식 및 Munk 표준 SOFAR 정규화 프로파일로 기술됩니다:

\[ c(z) = c_0 \left[ 1 + \epsilon \left( \eta + e^{-\eta} - 1 \right) \right], \quad \eta = \frac{2(z - \bar{z})}{B} \]

여기서 \(c_0 \approx 1500\,\text{m/s}\)는 SOFAR 축 중심 음속, \(\bar{z} \approx 1000\,\text{m}\)는 축 수심, \(B \approx 1300\,\text{m}\)는 척도 깊이(Scale Depth), \(\epsilon \approx 7.4 \times 10^{-3}\)는 섭동 파라미터입니다.

음선 궤적의 기하광학적 해밀토니안 방정식은 음선의 수평거리 \(r\)과 수심 \(z\), 위상각 \(\theta\)에 대해 스넬의 굴절 법칙에 따라 유도됩니다:

\[ \frac{dr}{ds} = \cos \theta(s), \quad \frac{dz}{ds} = \sin \theta(s), \quad \frac{\cos \theta(s)}{c(z(s))} = a = \text{const} \]

음선의 호의 길이 \(s\)에 따른 전파 지연시간(Travel Time) \(\tau\)는 다음의 선적분으로 주어집니다:

\[ \tau = \int_{\text{Ray}} \frac{ds}{c(z(s))} = \int_0^R \frac{dr}{c(z(r)) \cos \theta(r)} \]

수온 변화 \(\delta T(r, z)\)에 따른 음속 섭동 \(\delta c(r, z) = \alpha \delta T(r, z)\) (\(\alpha \approx 4.5\,\text{m/(s}\cdot^\circ\text{C)}\))에 의해 유기되는 도달시간 섭동 \(\delta \tau_i\)는 선형화된 단층촬영 역산 적분방정식으로 정식화됩니다:

\[ \delta \tau_i = - \int_{\Gamma_i} \frac{\delta c(r, z)}{c_0(z)^2} ds = \sum_{j=1}^M G_{ij} \delta c_j \]

3. 0-Alloc C# 초고속 10,000Hz 음선 추적 및 역산 서보 아키텍처

수중글라이더 및 UUV 온디바이스 DSP/FPGA 가속기 환경에서는 가비지 컬렉션(GC) 일시정지가 음파 펄스 위상 추적의 나노초(ns) 급 지터를 유발할 수 있습니다. 따라서 ref struct와 stackalloc Span<double> 기반의 0-Alloc RK4 음선 추적 및 Tikhonov 정규화 역산 서보를 구현합니다. 

using System;
using System.Runtime.CompilerServices;

namespace ZeliDesk.OceanAcousticTomography
{
    public readonly ref struct RayState
    {
        public readonly double Range;     // r (m)
        public readonly double Depth;     // z (m)
        public readonly double Angle;     // theta (rad)
        public readonly double TravelTime;// tau (s)

        [MethodImpl(MethodImplOptions.AggressiveInlining)]
        public RayState(double range, double depth, double angle, double travelTime)
        {
            Range = range;
            Depth = depth;
            Angle = angle;
            TravelTime = travelTime;
        }
    }

    public static class SofarRayTracerServo
    {
        private const double C0 = 1500.0;
        private const double ZBar = 1000.0;
        private const double BScale = 1300.0;
        private const double Epsilon = 0.0074;

        [MethodImpl(MethodImplOptions.AggressiveInlining)]
        public static double EvaluateSoundSpeed(double z)
        {
            double eta = 2.0 * (z - ZBar) / BScale;
            return C0 * (1.0 + Epsilon * (eta + Math.Exp(-eta) - 1.0));
        }

        [MethodImpl(MethodImplOptions.AggressiveInlining)]
        public static double EvaluateSoundSpeedGradient(double z)
        {
            double eta = 2.0 * (z - ZBar) / BScale;
            double deta_dz = 2.0 / BScale;
            return C0 * Epsilon * (1.0 - Math.Exp(-eta)) * deta_dz;
        }

        [MethodImpl(MethodImplOptions.AggressiveInlining)]
        public static RayState StepRk4(in RayState current, double ds)
        {
            // State: (r, z, theta, tau)
            double r = current.Range;
            double z = current.Depth;
            double th = current.Angle;
            double tau = current.TravelTime;

            // k1
            double c1 = EvaluateSoundSpeed(z);
            double g1 = EvaluateSoundSpeedGradient(z);
            double dr1 = Math.Cos(th);
            double dz1 = Math.Sin(th);
            double dth1 = - (g1 / c1) * Math.Cos(th);
            double dtau1 = 1.0 / c1;

            // k2
            double z2 = z + 0.5 * ds * dz1;
            double th2 = th + 0.5 * ds * dth1;
            double c2 = EvaluateSoundSpeed(z2);
            double g2 = EvaluateSoundSpeedGradient(z2);
            double dr2 = Math.Cos(th2);
            double dz2 = Math.Sin(th2);
            double dth2 = - (g2 / c2) * Math.Cos(th2);
            double dtau2 = 1.0 / c2;

            // k3
            double z3 = z + 0.5 * ds * dz2;
            double th3 = th + 0.5 * ds * dth2;
            double c3 = EvaluateSoundSpeed(z3);
            double g3 = EvaluateSoundSpeedGradient(z3);
            double dr3 = Math.Cos(th3);
            double dz3 = Math.Sin(th3);
            double dth3 = - (g3 / c3) * Math.Cos(th3);
            double dtau3 = 1.0 / c3;

            // k4
            double z4 = z + ds * dz3;
            double th4 = th + ds * dth3;
            double c4 = EvaluateSoundSpeed(z4);
            double g4 = EvaluateSoundSpeedGradient(z4);
            double dr4 = Math.Cos(th4);
            double dz4 = Math.Sin(th4);
            double dth4 = - (g4 / c4) * Math.Cos(th4);
            double dtau4 = 1.0 / c4;

            double nextR = r + (ds / 6.0) * (dr1 + 2.0 * dr2 + 2.0 * dr3 + dr4);
            double nextZ = z + (ds / 6.0) * (dz1 + 2.0 * dz2 + 2.0 * dz3 + dz4);
            double nextTh = th + (ds / 6.0) * (dth1 + 2.0 * dth2 + 2.0 * dth3 + dth4);
            double nextTau = tau + (ds / 6.0) * (dtau1 + 2.0 * dtau2 + 2.0 * dtau3 + dtau4);

            return new RayState(nextR, nextZ, nextTh, nextTau);
        }
    }
}

4. 2022 개정 교육과정 연계 및 군사·산업 응용

  • 지구과학 II: 해수의 성질 및 심층 순환 단원에서 수온약층의 역학적 안정도와 SOFAR 축 형성을 물리량 수치로 분석.
  • 물리학 II: 연속적인 굴절률 변화 매질에서의 파동 궤적 적분과 최소 작용 원리(페르마 원리).
  • 국방 MUM-T: 대잠 무인잠수정 간 음향 비오도메트리(Acoustic Odometry) 및 스텔스 수중 통신 채널 링크 확보.