2.1 극한 기온 정의

극한기후현상은 기온 및 강수량 등이 통계적 평년값 대비 비정상적으로 큰 편차를 보이는 현상으로, 사전 정의된 기준값(상대적·절대적 임계치)를 초과하거나 미달하는 기상학적 현상을 의미한다. 국가별 기준이 상이하나, 국제적으로는 주로 상·하위 퍼센타일 기반의 통계적 경계값을 기준으로 정의된다. 예를 들어, 한국 기상청은 평년 대비 현저히 높은 또는 낮은 값을 보이는 이상 현상을 극한기후로 분류하며, 이상고온은 일최고·최저기온이 90퍼센타일 초과, 이상저온은 10퍼센타일 미만일 경우로 정의하고 있다.

2.2 일반화 파레토 분포(GPD)

일반화 파레토 분포(GPD)는 특정 임계치(threshold)를 넘는 극단치들의 극한분포로 정의된다. 서로 동일하고 독립적인 극단치 $X_{1},\: X_{2},\: .... ,\: X_{n}$이 특정 분포 $F$에서 추출되었다고 할 때, $y=X-u$는 특정 임계치($u$)를 넘는 값으로 정의한다면, $y$의 조건부 확률분포는 식 (1)과 같이 나타낼 수 있다. 식 (1)에 나타난 $F_{u}(y)$는 Balkema-de Hann-Pickands 정리에 따라 식 (2)와 같은 GPD 분포함수로 수렴하게 된다. 여기서 $\beta$는 규모 매개변수(scale parameter), $\xi$는 형상 매개변수(shape parameter)이다. 만약, $\xi >0$이면 Type Ⅰ파레토분포(Pareto, Cauchy, t-분포 등), 이면 Type Ⅱ파레토분포(Uniform, Beta 분포 등), $\xi =0$이면 지수분포 형태를 가진다^[17].

GPD는 특정 임계치를 초과하는 극단값의 확률분포이므로 임계치 설정에 따라 초과치 표본 수가 결정되며, 이는 모수 추정 결과에 직접적인 영향을 미친다. 본 연구에서는 90퍼센타일 값을 임계치로 선정하였으며, 이는 IPCC 및 한국 기상청의 극한 기온 현상 정의에 기반한 것이다. 일반적으로 임계치는 MEF-그래프 및 Hill-그래프 등을 통해 최적값을 탐색하여 결정하나, 본 연구에서는 기상학적 정의에 따라 90퍼센타일을 사용하고 사후적으로 Kolmogorov-Smirnov 검정과 AIC 기준을 통해 모델 적합도를 검증하였다. GPD에서 각 모수는 임계치 를 초과하는 자료들로부터 최대우도추정 방법(Maximum Likelihood, ML)을 적용하여 추정하였다. ML 추정에 필요한 로그우도함수는 식 (3)과 같다^[17].

2.3 대표 기상관측소 선정 및 데이터베이스 구축

본 연구는 RTS GMLC 계통을 분석 대상으로 선정함에 따라, 각 Area의 대표 위·경도 지점을 기준으로 대표 기상관측소를 매핑하였다. 기온 정보는 미국 국립해양대기청(NOAA)에서 제공하는 Integrated Surface Database(ISD)의 Global Hourly 관측 자료를 활용하였으며, 다음 두 가지 기준을 적용하여 대표 기상관측소를 선정하였다. 첫째, 각 Area의 대표 지점으로부터 공간적 거리가 가장 인접한 관측소를 선정하였다. 둘째, 극한 기온 통계 추정의 신뢰성 확보를 위해 1990-2024년 기간 중 최소 30년 이상의 연속적 관측 데이터가 확보된 관측소를 선별하였다. 각 지역별 최종 선정된 대표 관측소의 세부 정보는 표 2에 정리하였다.

선정된 관측소의 원시 데이터는 관측소별 상이한 관측 시각을 가지므로, 00분 기준 1시간 간격 시계열로 재구성하였으며, 결측 구간에 대해서는 선형 보간법을 적용하였다.

이후 POT(Peak Over Threshold) 방법의 핵심 가정인 초과치 간의 통계적 독립성을 충족하기 위해, 전처리된 시간별 데이터로부터 일별 최고기온을 산출하고, 임계값을 연속적으로 초과하는 구간을 하나의 클러스터(cluster)로 정의한 후, 각 클러스터 내에서 최대값만을 추출하도록 하였다. 이를 통해 동일한 기상 이벤트로 인한 연속적인 고온 현상을 하나의 독립적인 극값 사건으로 처리하여 자기상관에 의한 종속성 문제를 제거하고, GPD 기반 통계 분석의 전제 조건인 초과치 간의 독립성을 강화하고자 하였다. 그 결과, Area 1의 초과값 개수는 427개에서 153개(35.8%)로, Area 2는 318개에서 125개(39.3%)로, Area 3은 379개에서 136개(35.9%)로 축소되어 통계적으로 독립적인 극값 집합이 도출되었다.

2.4 GPD 기반 극한 기온 시나리오 생성 및 검증

일반화 파레토 분포(GPD)를 기반으로 극한 기온 시나리오를 개발하고 그 적합성을 통계적으로 검증하였다. 본 연구에서는 2.3절의 전처리 과정을 바탕으로 여름철(6-9월) 기간에 대해 90 퍼센타일 임계치 기반의 POT 기법을 적용하였다. 여름철을 분석 대상으로 선정한 근거는 냉방 수요 급증에 따른 전력 부하 증가 및 고온에 의한 전력 설비의 열적 취약성이 동시에 발생하는 시기이기 때문이다. 실제로 2003년 8월 14일 북미 대정전 사태(U.S.-Canadian blackout)와 1996년 7월 2일 미국 남부 아이다호 및 유타주 과부하 사례 등 주요 정전 사례들이 여름철에 집중되어 발생한 바 있다.

모델의 신뢰성 검증을 위해 실측 누적분포함수(Empirical CDF) 및 추정된 GPD 누적분포함수 간 비교 분석을 수행한 결과는 그림 1에 제시하였다. 분포 적합도의 통계적 유의성 평가를 위한 Kolmogorov-Smirnov 검정 결과, Area 1의 p-value는 0.732, Area 2는 0.612, Area 3은 0.725로 산정되었다. 모든 지역에서 유의수준 $\alpha =0.05$ 기준을 상회하는 p-value를 나타내어, 추정된 GPD 모델이 실측 극한 기온 데이터의 분포 특성을 통계적으로 유의하게 재현함을 확인하였다. 추가적으로 지수분포, 감마분포, Webihull 분포 등 다른 확률분포와의 AIC(Akaike Information Criterion) 기준 비교 분석 결과, GPD 분포가 가장 낮은 수치를 보여 극한값 모델링에 최적임을 확인하였다.

그림 1. 지역별 실측 CDF 및 GPD의 CDF 비교 결과

Fig. 1. Comparison of empirical and GPD cumulative distribution functions by area

3.1 온도 변화에 따른 부하 및 역률 변화 모델링

전력 시스템에서 외기온도는 부하 패턴에 직접적인 영향을 미치는 주요 요인 중 하나이다. 본 연구에서는 2020년 1월부터 12월까지의 시간별 부하 데이터 및 기온 데이터를 활용하여 각 지역별 정규화 부하-기온 간의 상관관계를 정량적으로 분석하였다. 온도에 따른 부하 변화는 식 (4)와 같은 3차 다항식 회귀모델로 모델링하였다. 정규화된 부하는 온도 변화에 따른 상대적인 부하 민감도를 정량적으로 비교하기 위해 사용되었으며, 각 시점의 부하를 해당 지역의 평균 부하로 나누어 표준화한 값이다. 외기온도가 변화할 경우, 부하 버스 $i$의 유효전력 $P_{i}$은 식 (5)와 같이 초기 온도에서의 기준 부하 $P_{0,\: i}$에 온도 변화에 따른 부하 증가율을 곱하여 온도 변화 시의 실제 유효전력을 모델링하였다.

고온에서 냉방 부하가 증가하며, 냉방 부하는 일반적인 부하에 비해 무효전력 소비가 크고 역률이 낮은 특성을 갖는다. 본 연구에서는 온도 상승에 따른 역률 변화는 선형적 감소 관계로 가정하여 식 (6)과 같이 모델링하였으며, 이에 따른 무효전력은 식 (7)과 같이 가정하였다. 이때, $pf_{i}$는 버스 $i$의 역률, $pf_{0,\: i}$는 초기 온도에서의 역률, $k_{i}^{pf}$는 온도에 따른 역률 변화 계수를 의미한다.

3.2 온도 변화에 따른 동적 송전용량(DLR) 산정

송전선로의 허용전류는 대기온도를 비롯한 기상조건에 따라 차등 적용이 가능하다. 본 연구에서는 온도 조건을 고려한 송전선로 허용전류 산정을 통해 온도 변화에 따른 계통 구성요소의 취약도를 모델링하고자 하였다. 도체 주변의 기상요소 모니터링을 통해 기온, 풍속, 풍향 및 일사량을 측정하고, IEEE Std 738에 기반한 열평형 방정식을 적용하여 허용전류를 산정할 수 있다. 복잡한 대기온도 예측기술을 통해 송전용량을 정밀하게 차등화하려는 다양한 연구가 진행되고 있으나, 본 연구에서는 온도 변화에 따른 송전용량 변화를 선로 온도 증가 시 DLR이 선형적으로 감소하는 관계로 단순화하여 식 (8)과 같이 모델링하였다^[14]. 이때, $V_{ij}^{rated}$ 및 $V_{ij}$는 송전선로 정격전압(kV) 및 per unit 전압이며, $k_{ij}$는 도체 종류에 따라 결정되는 기울기 계수이다. 단순화를 위해 모든 송전선로에 동일한 $k_{ij}=0.02k A/℃$를 적용하였으며, $c_{ij}$는 식 (9)과 같이 산정된다. 송전선로 $l:i→j$의 초기 정격용량은 $T_{0,\: ij}=(T_{0,\: i}+T_{0,\: j})/2$을 기준으로 결정된다고 가정하였다.

3.3 시뮬레이션 초기 조건 및 계통 상태 설정

전력계통 연쇄사고 시뮬레이션에서는 여름철 평균 온도 조건을 가정하여 유효전력, 무효전력, 역률, 송전용량 등의 상태변수에 대한 초기 조건을 설정하였다. 극한 기온 시나리오 생성 시에는 각 지역의 GPD 분포로부터 상위 10% 이상의 확률 구간에서 온도 증가분을 샘플링하고, 무작위로 선택된 부하 버스를 중심으로 반경 50km 내의 모든 버스에 동일한 온도 증가를 적용하여 Temperature Disturbance Area를 정의하였다. 각 지역의 GPD 분포를 통해 산출된 극한 기온 분위값 및 초기 온도 대비 온도 증가분에 대한 정량적 분석 결과는 표 3 및 표 4에 제시하였다. 반경 50km는 실제 기상 현상의 국지적 영향 범위를 고려하여 설정되었으며, 이에 대한 적용 예시는 그림 2에 나타냈다. 정의된 영역 내 모든 버스에 대해 온도 증가로 인한 냉방 부하 증가, 송전용량(DLR) 저하, 역률 감소에 따른 무효전력 증가, 설비 탈락 확률 증가 등의 영향을 종합적으로 반영하였다.

그림 2. RTS-GMLC 계통에서 무작위로 선택된 부하 모선을 중심으로 한 극한 기온 발생 지역(반경 50km)

Fig. 2. Extreme temperature impact regions based on random load bus centroids within the RTS-GMLC system

4.1 연쇄사고 시뮬레이션 알고리즘

연쇄사고 시뮬레이션 알고리즘은 그림 3과 같다. Monte Carlo 방법을 통해 지역별 극한 기온 분포로부터 샘플링된 기온 상승분($\triangle T$)을 특정 지역의 부하 모선에 적용하고, 초기 시점에서 온도 상승에 따른 유효전력 부하 증가, 무효전력 증가, 역률 저하 및 동적송전용량 감소를 동시에 반영하였다. 선로 및 발전기 탈락 확률은 각각 과부하율 및 무효전력 초과 시간 누적량에 기반한 수학적 확률모델을 통해 계산하였으며, 보호 및 제어 조치로는 저전압 부하 차단 및 발전기 재급전을 적용하였다. 해당 로직은 ^[14]의 연구를 바탕으로 구현하였다.

그림 3. 연쇄사고 시뮬레이션 flowchart

Fig. 3. Flowchart of the cascading failure simulation

4.2 선로 탈락 확률 모델링

본 연구에서 선로 탈락은 (1) 과부하율에 따라 정의된 확률함수 기반으로 Uniform 난수와 비교하여 발생하도록 하였으며, (2) 한계 용량 초과 시 즉시 차단되는 방식이 아니라 동적 송전용량을 초과하여 운전될 때의 과부하 상태를 시간에 따라 누적하여 열적 한계에 도달하면 탈락하도록 하였다. 먼저, 송전선로 $l:i→j$에 대하여 과부하율 $R_{ij(T)}= F_{ij}(T_{s_{L}})/\overline{F_{ij}^{d}}(T_{ij})$라고 할 때, 과부하율에 따른 송전선로 탈락 확률 함수 $P_{ij}^{trip}$는 식 (10)과 같이 표현된다.

여기서 $f_{t}$는 식 (11)과 같이 구간별로 정의된다. 이때, $b_{1}=(\ln p_{1}-\ln p_{2})/(-ε)$, $a_{1}= p_{1}/e^{b_{1}}$, $b_{2}=(\ln p_{2}-\ln p_{3})/(1+ε-K)$, $a_{2}= p_{2}/e^{b_{2}(1+ε)}$이다. 또한, $ε=0.01$, $K=1.5$, $p_{1}=0.001$, $p_{2}=0.3$, $p_{3}=1$로 가정한다.

다음으로, 선로가 동적 송전용량을 초과하여 운전될 때 시뮬레이션 반복 $k+1$에서 선로 $l:i → j$의 과부하에 따른 탈락 시간 $Δt_{ij,\: k+1}$은 식 (12)와 같이 계산된다. 여기서 $\overline{o_{ij}}$는 열적 한계로서 정격용량의 150% 과부하를 20초간 지속할 수 있는 한계값을 의미하며, $\sum_{m=0}^{k}\triangle o_{ij,\: m}$는 누적된 부하량, $F_{ij}(t_{k})$는 $k$번째 반복에서의 선로 조류, $\overline{F_{ij}^{d}}(t_{k})$는 동적 송전용량을 의미한다. 해당 식은 과부하 상태에서 열적 한계에 도달하기까지 남은 시간을 계산하는 것으로, 누적 과부하 지속시간이 열적 한계 $\overline{o_{ij}}$에 도달하면 선로가 열적 손상으로 인해 탈락하게 된다고 가정한다.

4.3 발전기 탈락 확률 모델링

본 연구에서 발전기 탈락은 (1) 무효전력 출력에 따라 정의된 확률함수 기반으로 발생하는 확률적 탈락과 (2) 무효전력이 허용 한계를 20\% 초과한 상태가 누적 30분 이상 지속될 때 발생하는 과여자 탈락으로 구분된다. 모든 발전기는 AVR을 통해 단자 전압을 제어하나, 무효전력($Q$)이 허용 한계를 초과하거나 임계값에 근접할 경우, 발전기 탈락 확률이 증가한다. 먼저, 발전기 $i$에 대하여 무효전력 출력 $Q_{i}$가 허용 범위 $[\underline{Q_{i}},\: \overline{Q_{i}}]$를 벗어날 때의 탈락 확률 함수 $P_{i}^{trip}$은 식 (13)과 같이 표현된다.

여기서 $f_{g}$는 식 (14)와 같이 구간별로 정의된다. 이때, $b_{3}=(\ln p_{6}-\ln p_{5})/(\underline{K_{Q_{i}}}-\underline{Q_{i}}+\underline{ε})$, $a_{3}=p_{6}/e^{b_{3}\underline{K_{Q_{i}}}}$, $b_{4}=(\ln p_{5}-\ln p_{4})/(-\underline{ε})$, $a_{4}=p_{5}/e^{b_{4}(\underline{Q_{i}}-\underline{\epsilon})}$, $b_{5}=(\ln p_{4}-\ln p_{5})/(-\overline{ε})$, $a_{5}=p_{4}/e^{b_{5}\overline{Q_{i}}}$, $b_{6}=(\ln p_{5}-\ln p_{6})/(\overline{Q_{i}}+\overline{\epsilon}-\overline{K_{Q_{i}}})$, $a_{6}=p_{5}/e^{b_{6}(\overline{Q_{i}}+\overline{\epsilon})}$ 이다. 또한, $\underline{K_{Q_{i}}}=1.5\underline{Q_{i}}(Q_{i}<0)$, $\underline{K_{Q_{i}}}=-0.5(\underline{Q_{i}}=0)$, $\overline{K_{Q_{i}}}=1.5\overline{Q_{i}}$, $p_{4}=0.001$, $p_{5}=0.3$, $p_{6}=1$로 가정한다.

다음으로, 발전기가 무효전력 허용 한계를 초과하여 운전될 때 시뮬레이션 반복 $k+1$에서 발전기 $i$의 과여자에 따른 탈락 시간 $Δt_{g,\: k+1}$은 식 (15)와 같이 계산된다. 여기서 $\overline{o_{g}}$는 무효전력 허용 한계를 20% 초과한 상태를 30분간 지속할 수 있는 한계값을 의미한다. 해당 식은 무효전력 한계 초과 상태에서 무효전력 운전 한계에 도달하기까지 남은 시간을 계산하는 것으로, 누적 과여자 지속시간이 운전 한계 $\overline{o_{g}}$에 도달하면 과여자 보호장치 동작에 의해 발전기가 탈락하게 된다고 가정한다.

4.4 제어 및 보호 조치(Remedial Action)

사고대응을 위한 수정조치로는 저전압 부하차단(Undervoltage Load Shedding, UVLS) 및 발전기 재급전(Generator Redispatch)을 모델링하였다.