Research Paper

Journal of the Computational Structural Engineering Institute of Korea. 31 December 2023. 355-364
https://doi.org/10.7734/COSEIK.2023.36.6.355

ABSTRACT


MAIN

  • 1. 서 론

  • 2. 폭발 해석 방법론

  •   2.1 폭발 하중의 분류

  •   2.2 유한요소 모델링 기법

  •   2.3 Kingery and Bulmash(K-B) 곡선

  • 3. 수치해석 모델

  •   3.1 유한요소 모델

  •   3.2 재료모델

  • 4. 해석 결과 분석

  •   4.1 폭발 수치해석 모델의 검증

  •   4.2 폭발 충격파 전파 양상

  •   4.3 폭발압력과 격자망 크기에 대한 관계 분석

  •   4.4 TNT 중량에 따른 최적 격자망 크기 예측

  • 5. 결 론

1. 서 론

지난 2022년 발생한 러시아-우크라이나 전쟁은 폭발로 인한 비극적인 인명 손실과 심각한 사회 혼란을 극명하게 보여주었다. 휴전 국가인 대한민국은 북한의 여러 차례에 걸친 미사일 발사를 포함한 지속적인 군사적 도발의 위험에 처해있으며, 국가 기반 시설 및 군사 시설과 민간 지역이 잠재적인 폭발 위험에 노출되어 있다. 따라서 이러한 위험으로부터 구조물을 보호하는 방호 관련 연구에 대한 수요가 높아지고 있다.

폭발은 물질이 급격한 화학변화나 물리적인 변화를 일으켜 높은 온도와 에너지를 발생시키는 현상이다. 폭발의 유형 중 특히 폭굉은 음속보다 빠른 속도로 일어나는 화염의 전파를 의미한다(Choi and Sin, 2014). 대기 중에서는 폭굉을 통하여 큰 파괴력을 지니는 충격파 즉, 폭발 과압(overpressure)이 발생하는데 이는 폭발 발생 시 직접적인 인명사고를 동반할 수 있는 주된 위험요소이다(Bransford and Gayle, 1965, Larcher, 2007).

대기 중에서 발생하는 폭발 충격파의 효과 및 폭발압력을 파악하기 위해 많은 실험이 수행되었으며, 계측된 실험 결과를 기반으로 한 폭발압력-시간 이력에 관한 경험식이 다수 제시되었다(Brode, 1955; Henrych and Major, 1980; Kingery and Bulmash, 1984; Kinney and Graham, 2013; Mills, 1988). 여러 실험 데이터를 근거로 개발된 미 국방부 문서 UFC 3-340-02에는 폭발압력-시간 이력 곡선과 방폭설계에 대한 방호 가이드라인이 상세하게 제시되어 있다(US Department of Defense, 2008). 하지만, 이러한 폭발실험에는 안전성 및 경제성 측면의 제약이 따르며, 실제 실험 현장에서 매우 짧은 시간에 일어나는 화학반응을 측정하기에는 어려움이 존재한다.

이러한 실험의 단점을 보완하면서 폭발 현상을 자세하게 분석하기 위한 방법으로 유한요소 프로그램을 활용한 수치해석 연구가 진행되고 있다(Gang et al., 2012). Filice 등(2022)은 구조물에 대한 폭발하중을 조사하기 위해 근거리 폭발에 대한 작은 중량의 폭발 파라미터를 결정하였다. Xu 등(2020)은 유한요소 프로그램과 이론적 해석을 바탕으로 폭발물 이동속도에 따른 폭발압력을 조사하였다. Jung과 Hong(2021)은 1995년 테러 공격을 받은 Alfred Murrah 건물을 대상으로 폭발 수치해석을 수행하였고, 폭발 하중에 따른 부재의 저항력을 조사하였다.

보다 정확한 폭발 수치해석을 위해 전산유체역학 및 고체역학 시뮬레이션에서 주로 사용되는 ALE(Arbitrary Lagrangian- Eulerian) 기법이 고려될 수 있다. ALE 기법은 격자망(Mesh) 왜곡을 효과적으로 처리하고 계산 효율성을 개선할 수 있어 큰 변형이나 움직이는 경계와 같은 문제를 해석하는 데 유용하다. Cheng 등(2013)은 근거리 폭발로 인한 충격파를 예측하는 데 ALE 기법이 효과적임을 나타냈다. Wu 등(2022)은 폭발실험과 ALE 기법을 사용한 수치해석을 통해 철근 콘크리트 패널에 미치는 효과를 조사하였다.

하지만 ALE 기법을 사용한 폭발 수치해석은 대기를 구성하는 격자망의 크기에 따라 해석 결과의 정확도에 차이가 크게 발생한다고 알려져 있다. 또한 요소 개수가 많아지면 계산시간이 많이 소요되기 때문에 원거리 폭발 해석에서는 적용하기가 어렵고 근거리 폭발에서만 적용이 가능하다는 단점이 있다(Kim, 2007). 더불어, 폭발과 같은 충격파의 전파 속도가 빠른 동적 유한요소해석에서 ALE 기법을 사용한다면 격자망 크기가 일정 수준 이하로 줄어들 때 폭발압력이 수렴하지 않거나(Kwasniewski et al., 2010, Wu et al., 2022) 부정확한 결과가 나타날 수 있다고 알려져 있다(Nassiri and Kinsey, 2016, Nassiri et al., 2015, Wang et al., 2021).

이전까지 수행되었던 ALE 기반 폭발 연구는 폭발물의 중량과 측정지점 중 하나를 고정시켜 최적 격자망 크기를 산정한 후, 폭발압력을 측정하는 과정을 거쳤다. 하지만 폭발물의 중량에 대해 근거리부터 원거리까지 해당하는 넓은 측정 영역에서 높은 정확도를 보이는 대기 영역의 격자망 크기에 관한 연구는 부족한 실정이다.

본 연구에서는 대기 중 폭발물의 충격파 경향을 넓은 범위에서 예측할 수 있는 최적 격자망 크기에 관한 연구를 수행한다. 상용 소프트웨어인 LS-DYNA의 3차원 ALE 기법을 이용하여 대기 중 폭발 현상을 구현한다. 대기의 격자망 크기를 11개의 경우로 나눠 실험값과 비교가 어려운 근접거리를 제외한 영역들에 대해 최대 폭발압력과 충격파의 양상을 조사한다. 각 격자망 크기에 따른 입사압력을 UFC 3-340-02에서 제안한 실험값과 비교하여 수치해석 결과의 정확도를 계산함에 따라서 각 TNT 중량에 대한 최적의 대기 격자망 크기를 도출하고, 이와 함께 TNT 중량에 따른 최적 격자망 크기의 상관관계식을 제안한다. 본 연구결과는 대기 중의 폭발 현상을 이해하기 위한 향상된 수치 모델의 개발에 기여할 수 있을 것으로 사료된다.

2. 폭발 해석 방법론

2.1 폭발 하중의 분류

US Department of Defense(2008)에 따르면 외부에서 일어나는 폭발은 발생 위치와 구조물 및 지면의 관계에 따라 자유 대기 중 폭발(Free-air blast), 대기 중 폭발(Air blast), 지표면 폭발(Surface blast)로 분류된다. 이 중 자유 대기 중 폭발은 충격파가 대기 중에서 전파되어 구조물이나 지면 등의 영향을 받지 않고 구조물에 직접 입사하는 경우이다. 자유 대기 중에서 폭발이 일어나게 되면 폭발물 중심에서 급격하게 전파되는 구형의 충격파가 발생하고, 매우 빠른 속도로 대기압보다 큰 압력으로 충격전단이 생성된다. 이후 충격전단을 뒤따르는 기류의 전파로 충격전단의 압력은 감소하고 발생하는 충격량은 누적되며 대기압보다 낮은 압력을 나타낸다. 이러한 연쇄적인 폭발현상은 Fig. 1과 같이 폭발 충격파의 압력-시간 이력 곡선으로 나타난다. 폭발압력이 대기압보다 큰 상태를 정압단계라 하며, 작은 단계를 부압단계라 한다. 부압단계의 경우 최대압력이 정압단계와 비교해서 무시할 정도로 작게 나타나므로 폭발 해석에서는 주로 정압단계만을 고려한다(Børvik et al., 2009).

https://static.apub.kr/journalsite/sites/jcoseik/2023-036-06/N0040360601/images/Figure_jcoseik_36_06_01_F1.jpg
Fig. 1.

Pressure-time history curve of an explosion (Choi et al., 2013)

2.2 유한요소 모델링 기법

연속체의 움직임을 다루는 고전적인 방법으로는 주로 라그랑지안(Lagrangian) 기법과 오일러리안(Eulerian) 기법이 사용된다. Fig. 2는 라그랑지안, 오릴러리안 기법을 사용해 설계된 격자망이 변형에 대응하는 방법을 각각 나타낸 것이다. 라그랑지안 기법에서는 연속체가 움직일 때 유한요소의 절점이 물체와 함께 움직이거나 변형되어(Fig. 2a), 고체 구조물과 같이 물질 좌표계 중심의 모델링이 적합한 경우에 주로 사용된다. 반면에 오일러리안 기법은 절점이 물체와 함께 움직이지 않는 공간 좌표계가 중심이 되므로(Fig. 2b), 유체와 같은 물질을 모델링할 때 사용하는 것이 유리하다.

https://static.apub.kr/journalsite/sites/jcoseik/2023-036-06/N0040360601/images/Figure_jcoseik_36_06_01_F2.jpg
Fig. 2.

Mesh deformation: (a) Lagrangian description and (b) Eulerian description (Jang, 2007)

라그랑지안 기법은 격자망이 물체를 따라가기 때문에 작은 변형에 대하여 빠른 계산이 가능하지만, 폭발과 같은 대변형의 문제에서는 격자망의 과도한 변형을 유발하여 해석 시간 단위(time step)를 줄어들게 하고, 결과적으로 전체 해석 시간을 기하급수적으로 증가시키거나 해석이 중단되게 한다(Jang, 2007). 반면에 오일러리안 기법은 공간에 고정된 점에서 물체를 바라보기 때문에 고정된 요소들을 통하여 유동 해석을 수행할 수 있다. 그러나 요소들을 통한 유동량이 많아지면 에너지 손실이 발생하여 해석의 신뢰도가 떨어지고, 두 개 이상의 유체를 구현하는 경우에는 경계면의 추적이 어렵다는 단점이 있다.

ALE 기법은 전술한 두 가지 기법의 장점은 적용하고 단점은 최소화하기 위해 개발되었다. 해당 기법은 유동적이고 가변적인 경계조건에서도 높은 정확도의 해석을 수행할 수 있다. 또한, 독립적인 참조 좌표계가 물체에 설정되어 격자망이 물체의 움직임과는 별개로 이동하기 때문에 비틀리지 않고 재정렬된다. 이러한 특징으로 ALE 기법은 대변형에서도 격자망의 왜곡을 줄일 수 있으며 계산의 정확도를 높일 수 있다.

Fig. 3은 ALE 기법의 계산과정을 설명한다. ALE 기법에서는 분해 연산기법을 이용하여 라그랑지안 단계와 오일러리안 이류 단계로 나누어 계산한다. 요소의 변형이 크게 발생하면 완화 알고리즘을 이용하여 절점들을 이동시키는 라그랑지안 단계를 통해 요소를 수정한다. 이후 오일러리안 이류 단계에서는 요소의 물성치를 이동 전 절점들에서 이동 후 수정된 절점들에 대해 변환시킨다. ALE 기법은 라그랑지안 기법에서 허용되는 것보다 더 큰 연속체의 왜곡을 처리할 수 있으며, 오일러리안 기법보다 더 높은 해상도를 가질 수 있어 압력의 전파 및 대변형 해석에 유리하다(Donea et al., 2004).

https://static.apub.kr/journalsite/sites/jcoseik/2023-036-06/N0040360601/images/Figure_jcoseik_36_06_01_F3.jpg
Fig. 3.

Advection calculation process in ALE method (Gröger et al., 2022)

2.3 Kingery and Bulmash(K-B) 곡선

자유 대기 중 폭발에서는 폭발압력이 구조물 및 지면의 간섭 없이 입사압력의 형태로 나타나기 때문에, 폭발 수치해석결과와 실험을 통한 입사압력의 결과를 비교하여 정밀성을 검증할 수 있다. Fig. 4는 Kingery와 Bulmash에 의해 개발되어 UFC 3-340-02에서 제공하는 K-B 곡선을 나타내며, 환산 거리(Z)에 따른 최대 입사압력(Pso), 충격량 및 폭발압력 도달시간 등 다양한 폭발 매개변수 값들을 제공한다. 환산 거리(Z)는 폭발물의 중량(W)과 폭발물의 중심으로부터 떨어진 이격 거리(R)을 식 (1)에 대입하여 구할 수 있다.

(1)
Z=RW3

https://static.apub.kr/journalsite/sites/jcoseik/2023-036-06/N0040360601/images/Figure_jcoseik_36_06_01_F4.jpg
Fig. 4.

K-B curve in spherical TNT explosion (US Department of Defense, 2008)

폭발물의 중량은 TNT를 기준으로 계산한다. TNT가 아닌 폭발물에 대해서는 반응열, 최대압력, 충격량 등을 기준으로 하는 TNT 등가계수로 환산하여 사용한다(US Department of Defense, 2008).

K-B 곡선이 제시하는 환산 거리 중 폭발물과 근접한 거리(Z ≤ 0.8m/kg1/3)에서는 폭발 매개변수의 측정이 매우 어렵고, 측정하더라도 정밀하게 검토되지 못한다(Shin and Kim, 2016). 따라서 본 연구에서는 신뢰성 있는 검증을 위해 환산 거리 0.8 m/kg1/3 이후의 측정지점에서 측정된 최대 입사압력을 K-B 곡선과 비교한다.

3. 수치해석 모델

3.1 유한요소 모델

3차원 ALE 기법이 적용된 폭발 모델의 구축 및 해석을 위하여, 상용 유한요소 프로그램 LS-DYNA가 사용된다. Fig. 5는 본 해석에서 사용된 수치해석 모델을 보여준다.

https://static.apub.kr/journalsite/sites/jcoseik/2023-036-06/N0040360601/images/Figure_jcoseik_36_06_01_F5.jpg
Fig. 5.

Finite element geometry model

해석 모델에서 사용된 TNT와 대기는 ALE 요소로 모델링되며, 본 해석은 TNT의 중량이 25, 50, 75, 100kg인 경우를 가정하여 수행한다. 폭발물 중량이 커질수록 충격파의 전파 거리도 증가하기 때문에, 같은 환산 거리에 대해 분석하기 위하여 대기 영역의 크기도 증가시키며 해석을 수행해야 한다. TNT의 형상은 구형으로 설정하며, 중량에 따른 TNT의 반지름과 이에 맞추어 확장된 대기 영역의 크기를 Table 1에 나타냈다.

Table 1

TNT radius and air domain size

TNT weight (kg) TNT radius (m) Air domain size (m)
(x × y × z)
25 0.154 9.0 × 1.588 × 0.794
50 0.194 11.0 × 2.0 × 1.0
75 0.222 12.0 × 2.290 × 1.145
100 0.245 13.0 × 2.520 × 1.260

본 모델에서 적용된 경계조건은 Fig. 6에 나타내었다. 자유 대기 중 폭발을 구현하기 위해서는 경계면에서 폭발 충격파가 반사되거나 중첩되지 않고 최초의 입사파 형태로 전파되어야 한다. 이를 구현하기 위해서 NON_REFLECTING 키워드를 이용하여 반사파의 생성을 제거함으로써 자유 대기 평면을 모사한다.

https://static.apub.kr/journalsite/sites/jcoseik/2023-036-06/N0040360601/images/Figure_jcoseik_36_06_01_F6.jpg
Fig. 6.

Boundary condition of the numerical model

또한, ALE 기법은 연속체를 유한개의 요소로 나누어 해석하여 계산시간이 비교적 크기 때문에, 계산효율을 최적화하기 위한 목적으로 대칭경계조건을 설정한다. z축을 기준으로 절반의 영역을 모델링하며, 모델이 대칭을 이루는 x-y면에 대해서는 대칭 경계조건을 부여한다. Fig. 7은 폭발 수치해석 모델에서 해석결과의 분석을 위해 설정한 측정지점을 보여준다. 해석결과의 계측을 위한 측정지점은 TNT의 중심으로부터 x축으로 0.2m의 간격을 두어 모델링 된 대기 영역의 끝까지 설정한다.

https://static.apub.kr/journalsite/sites/jcoseik/2023-036-06/N0040360601/images/Figure_jcoseik_36_06_01_F7.jpg
Fig. 7.

Pressure measurement location for the numerical model

ALE 기법을 사용하여 폭발 수치해석을 하기 위해 ALE Multi- Material Group(AMMG) 키워드를 적용할 수 있다. AMMG는 라그랑지안 요소가 오일러리안 요소와 같이 이동하게 하는 조직화 절차로 연속체의 움직임을 계산한다. 대기와 TNT를 하나의 그룹으로 설정하여 ALE 기법을 이용할 수 있으며, AMMG로 정의된 재료들의 영역을 할당하는 키워드인 INITIAL_ VOLUME_FRACTION_GEOMETRY를 통해 TNT의 크기와 좌표를 입력한다.

3.2 재료모델

자유 대기 중 폭발 시뮬레이션을 구축하기 위해서 대기와 TNT의 재료적 특성을 나타낼 수 있는 물성값을 지정해야 한다. 또한, ALE 기법의 적용을 위해서 재료의 부피 및 밀도에 따른 압력의 변화를 설명할 수 있는 적절한 상태방정식이 필요하다. LS-DYNA의 MAT_NULL 재료모델을 사용하여 대기의 특성을 정의할 수 있으며(LSTC, 2013), 대기의 상태방정식은 선형 다항식 형태로 식 (2)과 같은 방정식을 사용한다.

(2)
P=C0+C1μ+C2μ2+C3μ3+(C4+C5μ+C6μ2)E0

여기 식 (2)에서 E0는 단위 체적당 초기에너지, 𝜇는 역학 점도이다.

(3a)
C0+C1+C2+C3+C6=0
(3b)
C4=C5=γ-1=0.4

대기를 이상기체로 가정했을 때, 식 (3a)식 (3b)에 적용된 재료모델과 상태방정식 매개변수는 Table 2에 정리하였으며, 선형 다항식의 매개변수인 C4E0를 조정하여 0.1MPa의 표준 대기압을 적용한다.

Table 2

Air material parameters (Sonntag et al., 2003)

MAT_NULL
RO (kg/m3) PC MU
1.25 0 0
EOS_LINEAR_POLYNOMIAL
C0C1C2C3C4C5C6E0 (Pa) V0
0 0 0 0 0.4 0.4 0 2.5e+6 1

폭발물인 TNT를 시뮬레이션에서 모사하기 위한 재료적 특성으로 LS-DYNA의 MAT_HIGH_EXPLOSIVE_BURN 재료모델을 사용하여 밀도와 초기속도를 부여한다. 또한, 폭발물의 거동을 설명하기 위해 Johnson-Wilkinson and Lee 상태방정식을 사용하며(식 (4)), 폭발물의 매개변수는 Table 3에 정리하였다.

(4)
P=A1-ωR1Ve-R1V+B1-ωR2Ve-R2V+ωE0V

여기서 P는 폭발압력, V는 상대 부피로 폭발물의 폭발 전・후의 밀도 비로 계산한다.

Table 3

TNT material parameters (Dobratz and Crawford, 1985)

MAT_HIGH_EXPLOSIVE_BURN
RO (kg/m3) D (m/s) PCJ (GPa)
1,630 6,930 21
EOS_JWL
A (GPa) B (GPa) R1R2 𝜔 E0 (Pa) V0
3.71 0.032 4.15 0.95 0.4 0.07 1

Table 3에서 RO는 TNT의 밀도, D는 폭굉 속도, PCJ는 압력 Chapman-Jouget 압력을 나타내며, A, B,R1,R2,𝜔는 TNT 상태방정식 계수로서 본 논문에서는 Dobratz와 Crawford가 실험을 통해 결정된 수치를 참고한다(Dobratz and Crawford, 1985).

4. 해석 결과 분석

3장의 내용을 바탕으로 폭발 수치해석모델을 구축하고 대기의 격자망 크기를 변경하며 해석을 진행한다. 이를 위해 고려된 4개의 중량에 대한 폭발 해석을 수행함과 동시에 대기의 격자망 크기를 0.015m부터 0.045m까지 11단계로 세분화하여 측정되는 입사압력을 조사한다.

4.1 폭발 수치해석 모델의 검증

실제와 같은 폭발 거동을 모사하기 위해 사용한 폭발 수치해석 모델에 대해 검증을 수행한다. 유한요소해석 프로그램LS-DYNA의 LOAD_BLAST_ENHANCED(LBE)는 2.3절에 언급한 K-B 곡선을 전산화하여 나타내는 키워드이다. LBE 기능을 이용하면 Kingery와 Bulmash의 폭발실험으로 얻어진 폭발압력-시간 이력곡선을 폭발 하중으로 적용하여 시뮬레이션을 할 수 있다(Shin et al., 2022). 본 연구에서는 3차원 ALE 기법을 이용한 수치해석 모델의 검증을 위해 특정한 측정지점에서의 폭발압력-시간 이력곡선을 LBE 기능으로 얻어진 결과와 비교한다. Fig. 8은 연구에서 사용한 ALE 기법과 LBE 기능으로 나타난 폭발압력-시간 이력곡선을 비교한 그림이다. 계측 지점은 TNT의 중심으로부터 4.4m의 이격 거리를 가진 곳에서 동일하게 측정되었으며, Table 4에는 ALE 기법과 LBE 기능에서 나타난 각 중량의 최대 입사압력과 LBE 기능의 결과에 대한 ALE 기법의 오차가 계산되어 있다.

https://static.apub.kr/journalsite/sites/jcoseik/2023-036-06/N0040360601/images/Figure_jcoseik_36_06_01_F8.jpg
Fig. 8.

Comparison of the pressure-time histories between LBE and ALE at a distance of 4.4m

Table 4

Maximum incident pressure and errors between LBE and ALE methods

Mass of TNT
(kg)
Maximum incident pressure (MPa) Errors
(%)
LBE ALE
25 0.367 0.368 -0.366
50 0.623 0.639 -2.627
75 0.839 0.860 -2.531
100 1.050 1.046 0.417

LBE와 ALE의 거리에 따른 폭발압력-시간 이력 곡선을 비교한 결과, 충격파 도달시간 및 그래프의 경향이 서로 일치하며, TNT의 중량이 커질수록 폭발압력은 빠르게 도달되고 입사압력은 더욱 크게 측정된다. 이를 통하여 본 연구에서 구축한 3차원 ALE 기법을 사용한 해석 모델은 높은 정확도로 폭발 충격파의 전파양상을 예측할 수 있다.

4.2 폭발 충격파 전파 양상

폭발의 영향을 받는 대기 영역이 모두 요소로 이루어져 있는 3차원 ALE 기법의 장점을 활용하여, 시간에 따라 대기로 전파되는 충격파의 압력을 관찰한다. Fig. 9는 수행된 수치해석의 결과로, 시간에 따라 폭발 충격파가 대기를 통해파 전파되는 형상을 등압선으로 나타낸 것이다. 구형의 TNT가 폭발하면 폭발물의 중심에서 발생하는 충격파의 형태 또한 구형으로 전파된다. 대기의 경계면에서는 자유 대기 평면을 모사한 키워드의 설정으로 충격파가 반사되지 않고 무한대기 평면이 사실적으로 구현된 것을 확인할 수 있다. 또한, 폭발물의 중심으로부터 멀리 떨어질수록 충격파 전단의 압력이 감소하고, TNT의 중량이 커질수록 같은 이격 거리를 가진 측정지점에서 폭발 충격파가 비교적 빠르게 도달된다.

https://static.apub.kr/journalsite/sites/jcoseik/2023-036-06/N0040360601/images/Figure_jcoseik_36_06_01_F9.jpg
Fig. 9.

Blast shock wave propagation by TNT weight: (a) 25kg;(b) 50kg;(c) 75kg; and(d) 100kg

4.3 폭발압력과 격자망 크기에 대한 관계 분석

폭발 수치해석을 통해 각 측정지점의 최대 입사압력 값을 기록하고 측정된 최대 입사압력 값을 실험값인 UFC 3-340-02의 K-B 곡선과 비교한다(US Department of Defense, 2008). 이를 위해 전체 측정 구간에서의 평균적인 예측 정확도를 알기 위해 다음과 같은 평균 제곱 오차(MSE, Mean Square Error)를 활용한다(식 (5)).

(5)
MSE=1ni=1n(Pi-Pi^)2

평균 제곱 오차는 측정된 입사압력(Pi^)을 K-B 곡선에서 얻어진 입사압력(Pi)과의 오차 제곱을 표본 개수로 나눈 것으로, 계산된 수치가 0에 가까울수록 해석 모델의 최대 입사압력 예측 성능이 우수하다.

폭발 수치해석에서 측정된 최대 입사압력과 K-B 곡선의 입사압력 간 계산된 평균 제곱 오차를 Table 5에 나타내었으며, 폭발의 정확한 예측이 어려운 근접거리와 K-B 곡선에서 입사압력이 대기압 수준까지 떨어지는 원거리를 제외한 환산 거리가 0.8m/kg1/3 < Z < 2.6m/kg1/3인 측정지점에서의 최대 입사압력을 계산한다.

Table 5

Mean square error (×10-4) of the incident pressure by mass of TNT

Element
size (m)
Mass of TNT (kg)
25 50 75 100
0.0150 5.3404 13.3225 130.3371 20.2176
0.0175 1.7149 8.5318 13.0532 15.9182
0.0200 0.9076 4.3900 8.9996 14.2494
0.0225 5.4327 1.7177 4.7272 9.1734
0.0250 8.7972 1.5885 2.6444 3.5422
0.0275 17.2375 3.2402 1.2184 2.5279
0.0300 32.4375 2.2126 1.6559 1.6361
0.0325 44.5457 12.2165 4.7376 1.2473
0.0350 62.5001 19.0535 7.3188 2.7043
0.0375 77.7399 28.3884 14.7746 5.6941
0.0400 85.6785 33.4823 20.9505 7.0607
0.0425 125.6843 53.2763 30.0126 15.0190
0.0450 141.4677 73.2802 38.9541 21.5930

다양한 크기의 대기 격자망을 사용한 수치해석 결과는 TNT 중량별 가장 낮은 평균 제곱 오차를 보여주는 격자망의 크기가 다르다는 것을 보여주며, 이를 통하여 TNT의 중량을 고려한 입사압력을 높은 정확도로 예측하기 위해서는 서로 다른 대기의 격자망 크기를 사용하여야 함을 알 수 있다.

Fig. 10Table 4에서 나타난 각 수치해석의 최대 제곱 오차를 격자망의 크기에 따라 그래프로 표현한 것이다. 중량별 평균 제곱 오차는 격자망의 크기가 커질수록 특정 지점에서 0에 수렴하다가 다시 증가하는 형태를 보인다. 이는 격자망의 크기가 작아질수록, 근접거리와 인접한 측정점의 압력이 과도하게 크게 측정되어 평균 제곱 오차가 증가한 것으로 사료된다. 따라서 격자망의 크기와 평균 제곱 오차는 단순 비례하는 것이 아니며 폭발물 중량별로 정밀한 결과를 보장할 수 있는 격자망의 크기가 특정 구간에 존재하는 것을 알 수 있다. 이는 폭발 수치해석에서 대기 영역의 격자망 크기가 최대 입사압력 예측에 중요한 변수로 작용함을 의미한다. 표현된 그래프에서 평균 제곱 오차의 최저점을 표시하여 해당 점에서의 격자망 크기를 최적 격자망 크기로 선택하였으며. 수치해석 결과를 통하여 선정된 TNT의 중량별 최적 격자망 크기를 Table 6에 나타냈다.

https://static.apub.kr/journalsite/sites/jcoseik/2023-036-06/N0040360601/images/Figure_jcoseik_36_06_01_F10.jpg
Fig. 10.

Mean square error by mesh size

Table 6

Optimal mesh size of air domain by mass of TNT

TNT mass (kg) 25 50 75 100
Mesh size (m) 0.0200 0.0250 0.0275 0.0325

Fig. 11은 TNT 중량 50kg의 해석에서 가장 작은 격자망 크기(0.012m), 가장 큰 격자망 크기(0.045m)와 얻어진 최적 격자망 크기(0.025m)에서 측정된 최대 입사압력을 K-B 곡선의 최대입사압력과 비교하여 계산한 상대오차를 환산 거리에 따라 나타낸 것이다. 같은 환산 거리를 기준으로 오차를 비교하면 격자망의 크기가 가장 큰 모델은 환산 거리 2.25m/kg1/3 이전까지 K-B 곡선보다 입사압력을 과소평가하여 측정하는 것을 볼 수 있다. 대기 영역의 격자망 크기가 작은 모델에서는 모든 측정 영역에서 오차가 음의 값으로 계산된다. 이는 격자망의 크기가 줄어들수록 최대 입사압력 값을 과대평가하여 예측하는 것을 의미한다. 그에 비해, 최적 격자망 크기에서 계산된 오차는 대부분의 환산 거리 영역에서 0을 기준으로 가장 적게 벗어난 결과를 보여주며, 넓은 영역의 환산 거리에서 입사압력을 높은 정확도로 예측할 수 있다.

https://static.apub.kr/journalsite/sites/jcoseik/2023-036-06/N0040360601/images/Figure_jcoseik_36_06_01_F11.jpg
Fig. 11.

Error with K-B curve of incident pressure example (50kg of TNT)

4.4 TNT 중량에 따른 최적 격자망 크기 예측

ALE 기법을 사용한 폭발해석은 대기의 격자망 크기에 따라 측정되는 폭발압력이 달라지기 때문에, 대기의 격자망 크기를 적절하게 선정하는 것이 중요하다. Fig. 12는 본 연구에서 제안하는 대기의 최적 격자망 크기를 25kg, 50kg, 75kg, 100kg의 TNT 중량에 따라 나타낸 것이다. 또한, 연구결과를 바탕으로 선형회귀를 통해 대기의 최적 격자망 크기를 예측할 수 있는 추세선을 표시한다. TNT 중량과 대기 영역 격자망 크기의 상관관계를 분석한 결과로, 다양한 중량의 폭발 해석에서 효율적인 결과를 나타내는 대기 영역의 최적 격자망 크기를 얻을 수 있다.

https://static.apub.kr/journalsite/sites/jcoseik/2023-036-06/N0040360601/images/Figure_jcoseik_36_06_01_F12.jpg
Fig. 12.

Linear regression result of optimal mesh size according to TNT weight

본 연구에서 제안한 네 가지의 중량에 따른 최적의 격자망 크기를 바탕으로 추세선을 그려보았을 때 TNT 중량과 최적 격자망 크기는 식 (6)과 같은 선형적인 관계를 보인다.

(6)
y=0.0002x+0.0163

여기서 x는 TNT중량(kg), y는 격자망 크기(m)이다. 선형회귀를 통해 나타낸 추세선의 결정계수(R-squared)는 0.9846으로, 본 논문에서 제안된 격자망 크기의 경향을 신뢰성 있게 설명한다. 식 (6)에 대한 검증을 위하여, TNT 150kg, 200kg 중량에 대해 제안된 식으로부터 도출한 격자망 크기로 수치해석을 수행하였다. 중량에 대해 도출된 대기의 최적 격자망 크기는 각각 0.0463m, 0.0563m이고, 해석의 영역은 근접거리 이후부터 측정하였다. Fig. 13는 환산 거리에 따른 입사압력을 나타낸 것이며, K-B 곡선과 수치해석의 결과를 비교하면 유사한 수준으로 입사압력을 예측함을 보인다.

https://static.apub.kr/journalsite/sites/jcoseik/2023-036-06/N0040360601/images/Figure_jcoseik_36_06_01_F13.jpg
Fig. 13.

Incident pressure comparison between K-B curve and proposed model

5. 결 론

폭발은 큰 위험을 동반하기 때문에 폭발 현상에 대한 분석과 이를 사실적으로 구현할 수 있는 수치해석 방법의 개발이 중요하다. 대부분의 기존 연구는 특정 환산 거리를 가진 단일 측정지점에서 폭발압력을 사실적으로 측정하기 위한 격자망 크기를 조사하였다. 그러나 넗은 환산 거리 영역을 가진 다수의 측정지점에서 정확도 높은 폭발압력을 측정하기 위한 연구는 부족한 실정이다.

본 연구는 대기의 격자망 크기가 폭발압력에 미치는 영향을 조사하고자 상용 유한요소 해석 프로그램인 LS-DYNA를 사용하여 3차원 ALE 기법으로 구현된 자유 대기 중 폭발 수치해석 모델을 구현하였다. 사용된 폭발해석모델의 검증을 위해 LBE 기능을 사용한 폭발압력-시간이력 곡선과 ALE 기법의 해석결과를 비교하였으며, 검증된 수치해석 모델로 다양한 범위의 격자망 크기에서 4가지의 중량을 가진 TNT의 폭발상황을 시뮬레이션하였다. 이를 통하여 TNT 중량에 따라 가장 높은 정확도를 보이는 대기의 최적 격자망 크기를 선정하였다.

폭발 해석에서 측정된 최대 입사압력은 대기의 격자망 크기에 따라 다르게 나타났다. 가장 높은 정확도를 보이는 최적의 격자망 크기는 TNT 중량에 따라 다르다는 것이 확인되었다. 최적의 격자망 크기를 대기에 적용하면, 넓은 환산 거리 영역에서 K-B 곡선의 입사압력과 높은 일치성을 보였다. 연구결과를 바탕으로 TNT 중량과 대기의 최적 격자망 크기가 나타내는 관계를 조사한 결과, 두 요인 사이에 선형적인 관계가 있음을 도출하였다. 결론적으로 폭발상황에서 TNT의 중량이 증가할 경우, 폭발압력의 전파경로인 대기의 격자망 크기도 중량에 비례하게 증가시키는 것이 정확한 해석을 위한 필수사항임을 알 수 있었다.

수치해석의 결과가 넓은 범위의 환산거리에서 실험값과 높은 일치성을 보임에 따라 본 연구에서 제안한 TNT 중량과 최적 격자망 크기의 관계식이 신뢰성 있는 폭발 해석 모델을 구축하는 데 활용될 수 있을 것으로 기대된다. 이런 정밀한 수치해석 모델은 폭발 상황에서 구조물의 설계 및 안전성 평가에 범용적으로 사용될 수 있다.

본 논문에서는 입사압력을 기준으로 3차원 ALE 기법을 사용하여 높은 정확도를 보이는 최적 격자망 크기를 제안하였다. 하지만, 대부분의 폭발 현상에서는 충격파가 지면이나 주변 구조물에 닿으면서 발생하는 반사압력이 더 중요하게 고려되어 진다. 따라서 추후 연구에서는 지면이나 구조물을 모사하여 발생하는 반사압력과 구조물 사이의 응답에 관한 해석을 진행할 예정이다.

Acknowledgements

본 연구는 원자력안전위원회의 재원으로 한국원자력안전재단의 지원을 받아 수행한 원자력안전연구사업의 연구결과입니다(No. 00242257).

References

1
Børvik, T., Hanssen, A.G., Langseth, M., Olovsson, L. (2009) Response of Structures to Planar Blast Loads-A Finite Element Engineering Approach, Comput. & Struct., 87(9-10), pp.507~ 520. 10.1016/j.compstruc.2009.02.005
2
Bransford, J.W., Gayle, J.B. (1965) Size and Duration of Fireballs from Propellant Explosions, No.NASA-TM-X-53314, National Space and Aeronautics Administration, USA.
3
Brode, H.L. (1955) Numerical Solutions of Spherical Blast Waves, J. Appl. Phys., 26(6), pp.766~775. 10.1063/1.1722085
4
Cheng, D.S., Hung, C.W., Pi, S.J. (2013) Numerical Simulation of Near-Field Explosion, J. Appl. Sci. & Eng., 16(1), pp.61~ 67.
5
Choi, H., Shin, D.H., Kim, M.S., Kim, D.J., Kim, H., Lee, Y.H. (2013) A Study on the Proposal of Residual Strength Equation of Reinforced Concrete Columns under Blast Load, J. Korean Soc. Hazard Mitig., 13(1), pp.17~24. 10.9798/KOSHAM.2013.13.1.017
6
Choi, M., Sin, P. (2014) An Analysis on the Main Causes of Explosion in Industrial Settings and a Study on its Prevention, J. Korean Soc. Hazard Mitig., 14(2), pp.209~216. 10.9798/KOSHAM.2014.14.2.209
7
Cook, R.D. (2007) Concepts and Applications of Finite Element Analysis, John Wiley & Sons.
8
Dobratz, B.M., Crawford, P.C. (1985) LLNL ExplosivesHandbook: Properties of Chemical Explosives and Explosive Simulants, Report No. UCRL-52997-Chg.2, Lawrence Livermore National Laboratory, USA.
9
Donea, J., Huerta, A., Ponthot, J.P., Rodríguez‐Ferran, A. (2004) Arbitrary L Agrangian-E Ulerian Methods, Encyclopedia of Computational Mechanics.
10
Filice, A., Mynarz, M., Zinno, R. (2022) Experimental and Empirical Study for Prediction of Blast Loads, Appl. Sci., 12(5), p.2691. 10.3390/app12052691
11
Gang, H.G., Jeong, Y.U., Hong, J.U. (2012) Overview of Blast Load Analysis, Comput. Struct. Eng., 25(4), pp.57~62.
12
Gröger, B., Wang, J., Bätzel, T., Hornig, A., Gude, M. (2022) Modelling and Simulation Strategies for Fluid-Structure- Interactions of Highly Viscous Thermoplastic Melt and Single Fibres-A Numerical Study, Mater., 15(20), p.7241. 10.3390/ma1520724136295308PMC9609902
13
Henrych, J., Abrahamson, G.R. (1980) The Dynamics of Explosion and Its Use, J. Appl. Mech., 47(1), p.218. 10.1115/1.3153619
14
Jang, I.H. (2007) Sloshing Response Analysis of LNG Carrier Tank Using Fluid Structure Interaction Analysis Technique of LS-DYNA 3D, Doctoral Dissertation, Korea Maritime and Ocean University.
15
Jung, J.W., Hong, J.W. (2021) Blast Response Simulation of the Alfred Murrah Building Reinforced by Use of HPFRCC, Int. J. Concr. Struct. & Mater., 15(1), pp.1~16. 10.1186/s40069-020-00442-9
16
Kim, J.H. (2007) Shock Response Analysis under Underwater Explosion for Underwater Ship using ALE Technique, J. Korean Soc. Mar. Environ. & Energy,10(4), pp.218~226.
17
Kingery, C.N., Bulmash, G. (1984) Airblast Parameters from TNT Spherical Air Burst and Hemispherical Surface Burst, US Army Armament and Development Center, Ballistic Research Laboratory.
18
Kinney, G.F., Graham, K.J. (2013) Explosive Shocks in Air, Springer Science & Business Media. p.269.
19
Kwasniewski, L., Balcerzak, M., Wojciechowski, J. (2010) A Feasibility Study on Modeling Blast Loading using ALE Formulation, In Cost Action C, 26, pp.127~132.
20
Larcher, M. (2007) Simulation of the Effects of an Air Blast Wave, Luxemburg Office of Official Publications of the European Communities.
21
LSTC, L.D. (2013) Keyword User's Manual Volume II-Material Models, Livermore Software, 887.
22
Nassiri, A., Chini, G., Vivek, A., Daehn, G., Kinsey, B. (2015) Arbitrary Lagrangian-Eulerian Finite Element Simulation and Experimental Investigation of Wavy Interfacial Morphology during High Velocity Impact Welding, Mater. & Des., 88, pp.345~358. 10.1016/j.matdes.2015.09.005
23
Nassiri, A., Kinsey, B. (2016) Numerical Studies on High-Velocity Impact Welding: Smoothed Particle Hydrodynamics(SPH) and Arbitrary Lagrangian-Eulerian (ALE), J. Manuf. Processes, 24, pp.376~381. 10.1016/j.jmapro.2016.06.017
24
Mills, C.A. (1988) The Design of Concrete Structures to Resist Explosions and Weapon Effects, FIP Notes, 2, pp.11~15.
25
Shin, H.S., Kim, S.W., Moon, J.H. (2022) Applicability Analysis of the FE Analysis Method Based on the Empirical Equation for Near-field Explosions, J. Comput, Struct, Eng. Inst. Korea, 35(6), pp.333~342. 10.7734/COSEIK.2022.35.6.333
26
Shin, J., Kim, J. (2016) Effects of Near-Field Explosion of High Explosives on Blast Overpressures and Impulses, J. Korean Soc. Hazard Mitig., 16(6), pp.21~28. 10.9798/KOSHAM.2016.16.6.21
27
Sonntag, R.E., Borgnakke, C., Van Wylen, G.J., Van Wyk, S. (2003) Fundamentals of Thermodynamics, 6thed., New York: Wiley.
28
US Department of Defense (2008) Structures to Resist the Effects of Accidental Explosions, UFC 3-340-02.
29
Wang, S., Islam, H., Soares, C.G. (2021) Uncertainty due to Discretization on the ALE Algorithm for Predicting Water Slamming Loads, Mar. Struct., 80, p.103086. 10.1016/j.marstruc.2021.103086
30
Wu, J., Liu, Z., Yu, J., Xu, S. (2022) Experimental and Numerical Investigation of Normal Reinforced Concrete Panel Strengthened with Polyurea under Near-Field Explosion, J. Build. Eng., 46, p.103763. 10.1016/j.jobe.2021.103763
31
Xu, Q.P., Su, J.J., Li, Z.R., Liu, Y., Jiang, H.Y., Huang, F.L. (2020) Air Blast Pressure Characteristics of Moving Charge, J. Phys.: Conf. Ser., 1507(3), p.032052. IOP Publishing. 10.1088/1742-6596/1507/3/032052
페이지 상단으로 이동하기