I. 서  론

대형 수평축 풍력터빈은 매우 저속으로 회전하고 그 날개-깃(blade)이 보통 3개 이하임에 따라 날개-깃 통과 주파수 조화소음(BPF 조화소음)은 저주파수(low frequency) /초저주파수(infrasound) 영역에 크게 기여한다. 최근 풍력터빈 저주파수/초저주파수 소음 관련 평가기법, 특성 및 영향에 대한 연구들이 수행되고 있으며, 이는 향 후 국가별 풍력터빈 인증 또는 평가 시 저주파수/초저주파수 소음 평가가 포함 또는 강화 될 수 있음을 암시한다.[1-3]

Ffowcs Williams와 Hawkings 식은 유동 내 움직이는 음원에 대한 일반적인 음장 예측식이다. 이를 이용한 회전자(rotor) 또는 팬(fan)의 방사소음 예측의 경우, 주어진 날개-깃 상 압력분포를 그대로 음원으로 적용하여 소음을 예측할 수 있다.[4] Lighthill[5]의 음향상사식을 바탕으로 Lowson[6] 은 회전하는 음원에 적용하기에 편리한 형태의 음장 예측식을 유도하였다. Lowson의 예측식은 음원을 크게 힘의 섭동항, 속도항 그리고 가속항으로 구분한다. 이 때, 힘 섭동항과 가속도항은 원음장(far field)에 기여하고, 속도항은 근음장(near field)을 대표하는 BPF 조화소음원들이다.

후류(wake)와 날개-깃 간 상호작용에 의한 영향을 무시 할 수 있으며 회전축 방향 유동이 일정하다면, 일반적으로 열린공간에서 회전자 또는 팬은 회전에 따른 날개-깃 주위의 유동환경(유속 및 받음각)이 고정적이라고 가정할 수 있다. 반면, 풍력터빈은 고도 별 풍속분포에 따라 회전각도 별 날개-깃 유동환경이 변한다. 그러므로 풍력터빈의 날개-깃 주위의 평균화된 유동장을 해석하더라도, 그 날개-깃 상 압력장이 시간에 따라 변한다. 즉, 풍력터빈 BPF 조화소음에 대한 힘 섭동항의 기여도를 기대 할 수 있다.

본 연구에서, 풍력터빈의 저주파수 예측을 위하여 복합 전산 공력음향학 방법(복합 CAA 방법)을 이용한다. 이 방법은 유동장에 대한 음압의 영향이 작다는 가정 하에 음원이 되는 유동장을 우선 계산하고 음원을 모델링 한 후 음장의 계산을 수행하는 방법이다.[7]

복합 CAA 방법 내 음장 계산은 Lowson 음향상사식(Lowson’s acoustic analogy), 준-실험 자가 광대역 소음 모델[8] 과 준-실험 유입 광대역 소음 모델[9]들을 통해 수행한다. 준-실험 자가 광대역 소음 모델과 준-실험 유입 광대역 소음 모델은 기존의 연구들에서 풍력터빈 소음을 예측하기 위해 적용되어 왔으며, 본 논문 내에선 단지 저주파수 영역에서 정량적인 예측결과 비교를 위해 고려한다.[10] 즉, Lowson의 음향상사식을 이용한 BPF 조화소음원에 대응하는 저주파수 터빈 소음에 대한 시간영역에서의 예측 및 분석이 본 연구의 주요 목적이다.

Lowson 음향상사식, 준-실험 자가/유입 광대역 소음 모델은 하나의 GUI(Graphic User Interface) 프로그램으로 제작하였으며, 각 음압 예측 시 입력 값 계산을 위해 XFOIL[11]을 이용하여 공력해석을 한다. 이 때, MPI(Message Passing Interface)를 이용하여 병렬계산을 한다. 개략적인 터빈 소음 예측 프로그램 구성도는 Fig. 1과 같다. GUI 프로그램은 총 4개의 모듈로 구성되었다. C# 을 이용하여 GUI 환경과 후처리 과정을 구현하였으며, Fortran를 바탕으로 XFOIL을 병렬해석 한다.

각 구성 모듈의 기능을 이해함으로서, 소음 예측의 일련의 과정을 설명할 수 있다.

우선 터빈 소음 예측 시 각 모듈 중 날개-깃 정보(blade data) 모듈로부터 날개-깃의 구성 에어포일, 비틀림각(twist angle), 길이(chord)의 회전자 반경방향 분포와 장착각(pitch angle), 장착축(pitch axis) 및 비틀림축(twist axis)등을 설정한다.

다음으로, 터빈 구동조건(operating condition) 모듈에서 고도 별 풍속 분포와 회전속도, 대기온도와 상대습도 등을 입력한다.

날개-깃 정보와 터빈 구동조건을 설정한 후 공력해석(aerodynamic analysis) 모듈에서, XFOIL은 자가 광대역 소음 예측 코드(BPM code)의 입력 값인 경계층 변위 두께(boundary layer displacement thickness)와 Lowson의 음향상사식의 입력 값인 날개-깃 상 압력계수(pressure coefficient)를 계산한다. 이 때, BEM 이론(blade element momentum theory)[12]에 따라 날개-깃을 여러 요소를 분할하고 회전 각도별 각 요소 주위 유동 조건을 데이터베이스화한다. 그리고 유동 조건 데이터베이스를 입력 값으로 하여 각 요소에 대한 공력해석은 XFOIL을 통해 수행한다.

그리고 대기 중 난류강도(turbulent intensity)와 난류크기(turbulent length scale)는 지표면 거칠기 가정에 따라 구해지고 Lowson의 유입 광대역 소음 예측식(Lowson’s formula)의 입력 값이 된다.[10]

회전각도 및 날개-깃 요소별 경계층 변위 두께, 압력계수, 난류강도 및 크기는 데이터베이스화 되어 소음 예측(noise prediction) 모듈에 입력된다. 풍력터빈에 관련한 두 광대역 소음 예측 모델의 자세한 적용 및 관련 수식은 Zhu의 논문으로 부터 참고 할 수 있으며, 2.2와 2.3 절에 관련 식을 개략적으로 설명하였다.[10]

본 논문에서는 터빈 저주파수 소음의 주된 분석은 Lowson 음향상사식을 이용한다. 이 때 음원을 섭동항, 속도항과 가속도항으로 분리하고, 각 주파수 기여도를 분석한다. 끝으로, 예측결과와 측정결과를 비교하여 정량적으로 풍력터빈 저주파수 소음을 분석한다. 정량적 비교 시 예측 결과는 준-실험 자가 광대역 소음 예측 모델[8]과 준-실험 유입 광대역 소음 예측 모델[9]을 포함한다.

II. 소음 예측식

2.1 고도 별 풍속 분포 모델링

날개-깃 주위의 유동환경을 결정하기 위하여 풍력터빈으로 입사되는 고도 별 풍속을 모델링하였다. Table 1에 IEC61400-11[13]에 적용 할 수 있는 다양한 지형조건에 따른 지표면 거칠기 값을 정리하여 나타내었다. 이 값과 IEC 61400-11 에서 정의하고 있는 고도 별 풍속 관련 식을 이용하여 날개-깃 요소 별 입사 풍속을 계산한다. 해당 식은 다음과 같다.

여기서, z0ref는 참조 지표면 거칠기 값으로 50 mm 이다. z0는 설치된 풍력터빈 주변을 가장 유사하게 묘사하는 지표면 거칠기 값에 대응하며, Table 1을 참조하여 50 mm로 설정하였다. H는 풍력터빈 회전자 중심높이를 의미하고, zref 는 참조높이로 10 m에 해당한다. Vs는 표준풍속(standardized wind speed)을 나타내며, Vz는 높이 z에 위치한 풍속계로부터 구한 풍속을 뜻한다. 본 연구에선, 식(1)을 이용하여 주어진 표준풍속으로부터 고도 z 에 따른 풍속 Vz의 분포를 역으로 구하였다. 해석 절차를 고려하면, 식(1)은 터빈 구동 조건 모듈에서 실행한다.

2.2 준-실험 자가 광대역 소음 예측 모델

Brooks 등은 NACA0012 에어포일에 대한 방대한 실험을 통해 소음 및 유동특성을 측정하였다.[8] 그들은 실험 결과들을 이용하여, Ffowcs Williams와 Hall[14] 이 제시한 뒷전에서의 음압 산란식을 기준으로 측정소음을 유동 특성값으로부터 함수화 하였다. 자가 광대역 소음 모델 중 난류경계층-뒷전(turbulent boundary layer-trailing edge, TBL-TE) 소음과 박리-실속(separation and stall, S-S) 소음을 풍력터빈 주요 소음원으로 고려한다.[15] 이 들 주요 소음 모델들에 대한 준-실험 자가 광대역 소음 예측식은 다음과 같다.[8]

식(2)에서 (4) 내, SPLp과 SPLs 는 각각 압면(pressure side)과 부압면(suction side)에 대응하는 TBL-TE 소음 모델이며, SPLS-S은 S-S 소음 모델에 대응한다. 각 모델에 따른 예측 결과는 1/3 옥타브밴드레벨의 음압이다. δ s* 날개 부압면에서의 경계층 변위 두께, δ p* 는 압면에서의 경계층 변위 두께, M 는 마하수, ΔL은 스팬의 길이를 나타낸다. 그리고 r은 에어포일 뒷전 중심에서 수음점간 거리이며, D는 3차원 상 수음점 위치에 따른 방사음의 방향 함수(directivity function)이다. 각 식에서 로그함수 형태의 우측 첫 항은 Ffowcs Williams와 Hall의 뒷전 산란식으로부터 유도된 것이고, 나머지 우측 항들은 실험결과와 산란식간의 주파수 별 차이를 보정하기 위한 식들이다.[8]

본 연구에선, 공력해석 모듈에서 각 날개-깃 요소에 대하여 XFOIL로부터 압면과 부압면의 경계층 변위 두께를 계산하고, 소음 예측 모듈에선 식(2)에서(4)를 이용하여 자가 광대역 소음을 예측한다.

2.3 준-실험 유입 광대역 소음 예측 모델

Lowson은 Amiet[16]의 유입 광대역 소음에 대한 예측 모델을 확장하여, 대기 중의 난류 크기와 난류 강도로부터 유입 광대역 소음을 예측하는 모델을 제시하였다. 이 때, 난류 크기와 난류 강도는 지형 조건 가정에 따른 실험적 모델로 구한다.[10] Lowson이 제안한 준-실험 유입 광대역 소음 예측 모델은 다음과 같다.

여기서, ρ0 는 대기밀도, L 는 난류크기, I 는 난류강도, K = π f c/u 는 파수, c는 에어포일 코드, u 는 평균풍속, r 는 수음점과의 거리, d 는 날개 너비, LFC 는 저주파수 보정항, S 는 압축성 시어함수, β 2 = 1 - M 2이다. SPLHINF 과 SPLLINF는 각각 고주파수와 저주파수에 대한 유입 광대역 소음의 1/3 옥타브밴드레벨이다. 본 논문에서는 저주파수에서 고주파수로 천이구간을 고려한 다음의 준-실험 유입 광대역 소음 예측식을 이용하였다.

식(9)의 입력 값인 난류 크기와 난류 강도는 각각 다음 식으로부터 계산한다.[10]

여기서, 

이며, Z 는 지면으로부터 높이, z0는 Table 1과 같이 묘사된 지표면 거칠기 이다. 이로서, 소음 예측 모듈에선 식(5)에서 (12)을 이용하여 유입 광대역 소음을 예측한다.

2.4 Lowson 음향상사식

 Lowson의 음향상사식은 일정한 속도로 병진운동하는 음원에 대한 다음의 음장 지배식으로부터 유도할 수 있다.[6]

이 때, xi, yi(i=1,2,3)은 각각 수음점과 음원에 대한 좌표계를 의미하며, Mr은 벡터 r 방향의 음원속도 성분에 대한 마하수, 벡터 r은 음원에서 수음점으로 향하는 벡터, r은 벡터 r의 크기, a0는 음속, Fi는 점-힘(point force)을 나타내고 ‘[ ]’은 지연시간(retarded time)에서 값을 나타낸다. 식(13)은 고정된 속도로 병진운동을 하는 점-힘이 시간에 대하여 섭동하지 않을 경우 음이 발생하지 않음을 나타낸다. 하지만 힘이 일정하더라도 회전하는 경우, 회전에서 생기는 수음점 대비 가속도의 영향으로 소음이 발생한다. 이러한 현상은 일반적인 회전기계에서 발생하는 BPF 조화소음의 원인이며, 관련 소음을 Gutin 소음이라고도 한다.[17]

날개-깃 상 그리고 주위 난류 유동장의 분포 소음원(distributed acoustic source)에 의한 음장은 다음 식으로부터 계산 할 수 있다.

여기서, g는 임의의 음원으로서 사극자(quadrupole), 쌍극자(dipole), 단극자(monopole)로 모델링 할 수 있다. 저마하수에서는 쌍극자 소음원의 기여도가 크기 때문에 식(14)에서 g = -∂Fi/∂xi로 표현 할 수 있다. 쌍극자 소음원의 움직임을 3 차원 Dirac δ 함수로 정의하면, 공간에서 움직이는 분포 점-힘들에 의한 음장은 다음과 같이 나타낼 수 있다.

Lowson은 식(14)에서 지연시간을 고려한 δ 함수의 공간/시간 미분을 고려하여, 임의로 움직이며 섭동하는 힘에 의한 원음장과 근음장에 대한 식을 아래와 같이 제시하였다.

원음장:

근음장:

식(16) 내 ∂Fi/∂t는 힘 섭동항, ∂Mr/∂t는 가속도항, 식(17)은 속도항을 대표한다. 이 때, 날개-깃 상 Mr, r, xi, yi는 BEM 을 이용하여 구한다.[12] 풍력터빈의 BPF 조화소음을 식(16)과 (17)을 이용하여 계산하며, 입력 값인 점-힘은 XFOIL로부터 구한 각 요소 별 압력분포로부터 구한다.[6] 이 때, 난류천이는 자유천이(free transition) 조건이다. 풍속 6 m/s 기준, XFOIL 계산 시 날개-깃 팁(tip) 상 레이놀즈 수와 마하수는 각각 2.8 × 106와 0.23이다.

III. 예측 결과

제시한 방법에 따라 풍력터빈의 저주파수 소음 예측을 수행하였다. 저주파수 소음의 그 주파수 범위에 대한 각 국가 간 정의는 상이하며, 본 논문에서는 그 범위를 20-160 Hz로 정의하였다. 예측 검증을 위한 소음 측정 시 위치는 IEC61400-11 2nd의 기준위치 1(reference position 1) 이다.[13] 대략적인 측정위치는 터빈으로부터 약 120 m 떨어진 하류 위치이다. 예측 시, 실제 운용중인 풍력터빈의 날개-깃 설계정보 및 그 구동조건을 입력 값으로 사용하였다.

Figs. 2와 3은 표준풍속 6 m/s에서 BPF 조화소음 내 힘 섭동항, 가속도항과 속도항의 기여도를 각각 시간과 주파수 영역에서 분석한 결과를 나타낸다. Fig. 2의 시간영역 예측 결과를 통하여 힘 섭동항의 전체음압에 대한 기여도를 확인 할 수 있다. 또한, 주어진 위치에서 원음장/근음장 예측식으로부터 구한 각 음압 간 위상차로 인하여 전체 음압의 최대 진폭이 근접음원 단독 기여 일 때 보다 작다는 것을 알 수 있다. Fig. 3의 주파수영역 예측 결과는 힘 섭동항이 고차 BPF에서 큰 기여도를 가지는 것을 나타낸다. 예측 결과에서 근음장 예측식 내 속도항의 기여도가 1st BPF에서 가장 크다. 그리고 시간영역 결과에서 확인할 수 있었던 수음위치에서 각 음원항간 위상차에 의한 간섭으로 인한 영향으로 1st BPF에서 전체음압이 속도항에 대응하는 음압 보다 작게 나오는 것을 알 수 있다.

Figs. 4와 5는 각각 표준풍속이 6 m/s 그리고 9 m/s 일 때, 풍력터빈 측정 소음 대비 세 종류의 복합 CAA방법에 따른 예측결과들 간 비교를 나타낸다.[18] 이 들 결과는 본 연구 내 각 BPF 소음원 항의 주파수 기여도 분석과 연계하여 이해한다. 각 비교는 IEC 61400-11 기준으로 소음을 측정한 결과(Measurement), 식(16)과 (17)로부터 예측한 결과(Lowson code), Brooks, Pope와 Marcolini의 준-실험 자가 광대역 소음모델과 Lowson의 준-실험 유입 광대역 소음모델에 의한 예측 간 합(BPM-L code), 나머지는 모든 예측 결과의 합(Total)을 포함한다.[18]

Fig. 4에서 풍속이 6 m/s인 경우 저주파수 영역에서 광대역 소음 음압크기가 BPF 조화소음 성분의 그 크기보다 높다는 것을 알 수 있다. 그러나 Fig. 5와 같이 풍속이 9 m/s로 증가하면, 저주파수(80 Hz 이하) 영역에서 광대역소음 성분보다 BPF 소음의 기여도가 높아진다. 이는 시간, 주파수 분석에서 확인할 수 있었던 고차 BPF 소음에 크게 기여하는 힘 섭동항의 영향 때문이다. 또한, 정량적으로 실측결과와 예측결과가 잘 일치하는 것을 확인할 수 있다. 특히 준-경험식 광대역 예측모델에 추가적으로 저주파 소음의 예측 모델을 포함함으로써, 풍속 증가에 따른 80 Hz 이하에서의 방사음압의 증가 경향을 잘 반영하고 있음을 알 수 있다.

IV. 결  론

Lowson 의 음향상사식과 XFOIL 코드를 이용한 복합 CAA 방법으로 부터 풍력터빈의 주요 저주파 소음 성분인 BPF 조화소음을 예측하였다. BPF 조화소음을 각각 힘 섭동항, 가속도항과 속도항으로 분리하여 그 주파수 기여도를 분석하였다. 그 결과, 속도항과 가속도항이 저차 BPF에서 크게 기여하는 반면 힘 섭동항은 고차 BPF에서 주로 기여함을 확인 하였다. 픙속이 6 m/s와 9 m/s 일 경우에 대한, 방사소음 측정과 예측 결과 모두 80 Hz 이하에서 풍속 변화에 따른 민감한 음압크기 변동 양상을 보였다. 이 급격한 변화는 고차 BPF 에 대응하는 Lowson 식 내 힘 섭동항 때문인 것으로 확인하였다. 또한 정량적 의미에서 복합 CAA 방법에 기초한 BPF 조화소음 예측은 광대역 소음원만에 의한 예측결과 보다 저주파수에서 예측결과 정확성을 향상시킬 수 있다는 것을 확인하였다.

본 연구에서, 고도 별 풍속 분포에 관련한 힘 섭동항이 저주파수 영역에서 큰 기여도를 가짐을 확인하였다. 이는 풍력터빈 평가 장소와 설치 장소간의 고도 별 풍속 기울기 분포가 상이 할 경우, 저주파수 영역에서 소음 평가 결과는 실제 설치지역에서의 풍력터빈 저주파수 소음 특성을 대표하기 힘들 수 있음을 암시한다. 이와 같은 문제를 해결하기 위해 본 예측 방식을 응용 할 수 있을 것으로 기대한다.