Research Paper

Journal of the Computational Structural Engineering Institute of Korea. 31 October 2026. 345-358
https://doi.org/10.7734/COSEIK.2026.39.5.345

ABSTRACT


MAIN

  • 1. 서 론

  • 2. Naghdi 쉘 정식화와 열변형 도입

  •   2.1 기하학적 예비 사항

  •   2.2 Naghdi 쉘 모형과 변형 측도

  •   2.3 Duhamel-Neumann 분해와 열 변형 측도

  •   2.4 탄성 에너지, 열 등가하중 및 변분 손실 함수

  • 3. 물리정보 신경망 구현

  •   3.1 Deep Ritz Method

  •   3.2 Naghdi 쉘 모형에의 적용

  •   3.3 구현 세부 사항

  • 4. 수치 검증

  •   4.1 검증 설정

  •   4.2 Fully pinned 평판(두께 온도 구배)

  •   4.3 Cantilever 평판(두께 온도 구배)

  •   4.4 Clamped 원통 쉘 (균일 가열)

  •   4.5 종합 논의

  • 5. 결 론

1. 서 론

쉘 구조는 항공우주, 자동차, 압력용기, 열교환기 등 다양한 산업 분야에서 핵심적인 구성 요소로 사용되며, 실제 운용 환경에서 기계 하중과 함께 열하중을 받는 경우가 많다. 예를 들어 가스 터빈의 연소실 라이너, 우주 발사체의 추진제 탱크, 발전소의 압력 용기, 전자 패키지의 박막 구조 등에서는 비균일한 온도 분포에 의해 발생하는 열변형이 전체 구조 거동에 직접적인 영향을 미친다. 쉘 이론은 3차원 탄성 문제를 중간면에 정의된 2차원 모델로 환원함으로써 해석의 효율성을 크게 향상시키며(Koiter and Simmonds, 1973; Krätzig, 1980; Niordson, 1985), 그중 Naghdi 5-parameter 쉘 이론은 횡전단 변형을 포함하므로 비교적 두꺼운 쉘 구조까지 기술할 수 있다. 그러나 일반적인 형상, 경계조건, 하중 조건을 갖는 쉘 방정식의 해석해는 제한적인 경우에만 알려져 있으며, 실제 공학 문제에서는 대부분 수치해석에 의존한다.

쉘 구조의 수치해석에는 유한요소법(Finite Element Method, FEM)이 표준적으로 사용되어 왔으며, 다양한 쉘 요소와 안정화 기법이 개발되어 산업용 해석 프로그램에서 폭넓게 활용되고 있다. 국내에서도 degenerated 쉘 유한요소의 개선(Choi and Yoo, 1990), 향상된 쉘 요소를 이용한 유한요소해석(Heo et al., 2000), 그리고 기하학적으로 엄밀한 쉘 요소의 등기하해석(Kim et al., 2024) 등 쉘 구조 전산해석을 위한 방법론 연구가 지속적으로 수행되어 왔다. 특히 정적 순방향 해석의 경우 FEM은 계산 효율성과 정확도 측면에서 매우 성숙한 방법론이다. 한편 최근 쉘 해석의 관심은 단일 순방향 해석을 넘어 열-기계 결합 문제, 측정된 변형으로부터 미지의 열원이나 재료 물성을 추정하는 역문제, 그리고 실험 데이터와 물리 모델을 함께 활용하는 데이터 융합 문제로 확장되고 있다. 이러한 문제들은 FEM 기반 접근으로도 다룰 수 있으나, 목적함수의 민감도 계산, 반복적인 매개변수 해석, 또는 외부 최적화 및 데이터 동화 알고리즘과의 결합이 추가적으로 요구되는 경우가 많다.

이러한 배경에서 물리정보 신경망(Physics-Informed Neural Network, PINN; Lagaris et al., 1998; Raissi et al., 2019)이 계산역학 문제를 해석하기 위한 하나의 대안적 접근으로 주목받고 있다. PINN은 자동미분(automatic differentiation)을 이용하여 임의의 위치에서 미분 가능한 연속 근사해를 제공하며, 물리 법칙과 측정 데이터를 단일 손실함수 안에 함께 포함할 수 있다. 특히 변분 원리에 기반한 약형(weak-form) 또는 Deep Ritz 정식화(E and Yu, 2018)는 강형 잔차 기반 PINN에 비해 요구되는 미분 차수가 낮고, 전체 에너지 또는 약형 잔차를 손실함수로 사용한다는 점에서 구조역학 문제와 자연스럽게 연결된다. Bastek과 Kochmann(2023)은 이러한 변분 기반 접근을 Naghdi 쉘 이론에 적용하여 막 변형, 휨 변형, 횡전단 변형의 에너지와 외력 일을 하나의 손실함수로 구성하는 약형 쉘 PINN 정식화를 제시하였다.

평판 및 쉘 구조에 PINN을 적용한 연구는 최근 활발히 이루어지고 있다. Li 등(2021)과 Zhuang 등(2021)은 Kirchhoff 평판 휨 문제에 PINN을 적용하였으며, Yan 등(2022)은 PINN과 Extreme Learning Machine을 결합하여 복합재 평판 및 얕은 원통 쉘의 기계적 거동을 해석하였다. 쉘 구조에 대해서는 Bastek과 Kochmann(2023)이 Naghdi 5-parameter 쉘 이론에 기반한 약형 PINN 프레임워크를 제시하고, 다양한 곡면 쉘의 기계 하중 문제에 대해 FEM과 비교 가능한 정확도를 보였다. 그러나 기존의 쉘 PINN 연구는 주로 기계 하중에 의한 변형을 대상으로 하였으며, 열하중에 의한 쉘 구조의 열변형은 충분히 다루어지지 않았다. 한편 PINN을 열탄성 문제에 적용한 연구도 진행되어 왔다. Harandi 등(2024)은 3차원 솔리드의 열-기계 결합 문제에 대해 변위, 응력, 온도, 열유속을 함께 예측하는 혼합 정식화(mixed formulation) 기반 PINN을 제시하였다. 그러나 이러한 열탄성 PINN 연구의 대상은 대부분 3차원 솔리드 또는 1차원 구조에 집중되어 있으며, 평판 및 곡면 쉘의 열변형을 약형 쉘 이론 안에서 직접 다룬 연구는 아직 제한적이다.

본 연구의 목적은 FEM을 대체하는 범용 해석 도구를 제시하는 데 있지 않다. 정적 순방향 해석에 한정하면 FEM은 이미 충분히 검증된 방법이며, 계산 속도와 절대 정확도 측면에서 PINN이 FEM보다 우수하다고 일반화하기는 어렵다. 그럼에도 쉘 구조의 열변형 해석에 약형 PINN을 도입하는 것은 다음과 같은 의미를 갖는다. 첫째, 변분 기반 PINN의 손실함수는 전체 포텐셜 에너지의 형태로 구성되므로, 열변형과 같은 추가적인 물리 효과를 에너지 항의 형태로 자연스럽게 포함할 수 있다. 둘째, PINN은 물리 법칙과 측정 데이터를 동일한 학습 목적함수 안에서 결합할 수 있으므로, 향후 측정된 변형장으로부터 열하중이나 재료 물성을 추정하는 역문제 및 데이터 융합 문제로 확장될 가능성이 있다. 따라서 주어진 온도장에 의한 쉘 구조의 열변형을 약형 PINN으로 정식화하고 검증하는 것은, 열-기계 쉘 문제를 위한 데이터 기반 변분 해석의 기초 단계로 볼 수 있다.

본 연구에서는 Naghdi 5-parameter 쉘 이론에 기반하여, 주어진 열하중에 의한 쉘 구조의 열변형을 해석하기 위한 약형 PINN 프레임워크를 제시한다. 열변형은 Duhamel-Neumann 관계식에 따라 전체 변형률을 탄성 변형률과 열 변형률로 분해하여 도입한다. 이 분해를 표준 쉘 에너지 범함수에 대입하면, 기존 약형 PINN의 손실함수에 열에 의한 막 및 굽힘 등가하중 항이 추가되는 형태로 정식화된다. 따라서 신경망 구조나 학습 방식의 근본적인 변경없이 열하중 효과를 변분 손실함수에 포함할 수 있다. 본 연구의 주요 기여는 다음과 같다. 첫째, Duhamel-Neumann 분해를 이용하여 Naghdi 쉘의 약형 PINN 손실함수에 열변형을 일관되게 포함하는 정식화를 제시한다. 둘째, pinned 평판, cantilever 평판, 양단 고정 원통 쉘의 세 가지 열하중 benchmark 문제를 통해 제안된 방법을Abaqus 유한요소 해석 결과와 정량적으로 비교한다. 셋째, 이러한 검증을 통해 주어진 열하중을 받는 쉘 구조의 열탄성 응답을 약형 PINN으로 재현할 수 있음을 보이고, 향후 다물리 결합 및 데이터 융합 쉘 해석으로의 확장 가능성을 논의한다.

본 논문의 구성은 다음과 같다. 제2장에서는 Naghdi 쉘 이론과 열변형을 포함하는 변분 정식화를 제시하고, 제3장에서는 약형 손실함수에 기반한 PINN 구현 방법을 설명한다. 제4장에서는 세 가지 검증 문제에 대한 해석 결과를 Abaqus 유한요소 해석 결과와 비교하고, 제5장에서 결론과 향후 연구 방향을 제시한다.

2. Naghdi 쉘 정식화와 열변형 도입

본 장에서는 PINN 학습의 손실 함수로 사용될 전체 포텐셜 에너지를 구성하기 위한 정식화를 제시한다. 먼저 쉘 운동학을 기술하는 데 필요한 미분기하학적 예비 사항을 정리하고, Naghdi 5-parameter 쉘 이론에 기반한 변형 측도(strain measure)와 탄성 에너지를 유도한다. 이어서 Duhamel-Neumann 열탄성 관계식에 따라 열 등가하중을 도입함으로써, 기존 약형 쉘 PINN의 변분 손실 함수에 자연스럽게 결합되도록 한다. 본 정식화는 Bastek과 Kochmann(2023)의 약형 쉘 PINN 정식화를 기반으로 하며, 본 연구의 기여는 열 항의 도입에 있다. 본 장에서 다루는 쉘 운동학의 전체 구도는 Fig. 1에 도식화되어 있다.

https://cdn.apub.kr/journalsite/sites/jcoseik/2026-039-05/N0040390508/images/Figure_jcoseik_39_05_08_F1.jpg
Fig. 1.

Overview of the Naghdi shell kinematics

2.1 기하학적 예비 사항

쉘 운동학을 정의하기 위해 필요한 미분기하학적 개념을 정리한다. 본 절의 유도는 Niordson(1985)과 Bastek과 Kochmann(2023)을 따른다. 라틴 첨자 i,j,…는 1, 2, 3의 값을, 그리스 첨자 𝛼, 𝛽, ...는 1, 2의 값을 갖는 것으로 한다. 중간면의 위치 벡터의 직교좌표 성분을 xi로 두면, 매개변수 영역 ω⊂ℝ2 위의 곡선좌표 (ξ1,ξ2)를 이용하여 xi=fi(ξ1,ξ2),(ξ1,ξ2)∈ω같이 매개변수화할 수 있다.

여기서 fi는 ξα에 대해 연속 미분 가능한 세 개의 독립적인 단일값 함수이다. 쉘의 국소 양을 기술하기 위해, 중간면 각 점에서 공변 국소 기저를 다음 식 (1)과 같이 정의한다.

(1)
aα=f,αi≡∂fi∂ξα,a3=a3=a1×a2a1×a2

여기서, “,”는 좌표 ξα에 대한 편미분, “×”는 3차원 외적, “v”는 벡터 v의 l2 노름을 나타낸다. a3는 단위 법선 벡터이므로 자기 자신과 쌍대 관계(a3=a3)를 가진다. 공변 기저로부터 중간면의 제1 기본형식(메트릭 텐서)이 다음 식 (2)와 같이 유도된다.

(2)
aαβ=aβα=aα·aβ=f,αi·f,βi

여기서, “·”는 3차원 내적이다. 공변과 반변 기저 사이의 관계

(3)
aα·aβ=δαβ

로부터, 반변 메트릭 aαβ는 다음 조건으로 정의된다(식 (4)).

(4)
aαγaγβ=δβα

여기서, δβα는 Kronecker 델타이다(Niordson, 1985). 메트릭 텐서 aαβ와 그 역 aαβ는 첨자의 올림과 내림에 사용된다. 예를 들어,

(5)
vα=aαβvβ,vα=aαβvβ

메트릭 텐서는 표면의 내재(intrinsic) 기하를 나타내며, 이를 이용하여 곡선의 길이, 곡선 사이의 각, 그리고 면적이 표현된다. 임의의 3차원 좌표계에서 중간면 위의 요소 면적 dA는 다음 식 (6)과 같이 표현된다.

(6)
dA=∥a1×a2∥dξ1dξ2=det(aαβ)dξ1dξ2=a1/2dξ1dξ2

여기서, a=det(aαβ)이다. 한편 내재 기하만으로는 3차원 공간상의 표면 형태가 유일하게 결정되지 않으며(Niordson, 1985), 표면이 접평면에 대해 어떻게 휘어 있는지는 외재(extrinsic) 기하를 나타내는 제2 기본형식 또는 곡률 텐서 bαβ로 표현된다.

(7)
bαβ=bβα=a3·aα,β=-a3,β·aα

3차원 공간상의 점은 중간면의 좌표 xi로부터 법선 방향으로 ξ3∈[-h2,h2] 만큼 떨어진 위치에 있으며, 그 좌표를 Xi로 표기한다(식 (8)).

(8)
Xi=Fi(ξ1,ξ2,ξ3)=xi+ξ3a3=fi(ξ1,ξ2)+ξ3a3

여기서, 두께 h는 중간면 위에서 일정하다고 가정한다. 식 (3), (5), (7)을 이용하면, 3차원 기저 벡터 gi와 중간면 기저 ai의 관계는 다음과 같이 표현된다.

(9)
gα=∂Fi∂ξα=aα+ξ3a3,α=aα-ξ3bαβaβg3=a3

곡선 표면 위 접선 벡터장 v=vαaα의 미분은 다음과 같다.

(10)
∂v∂ξβ=vα,βaα+vαa,βα

여기서, 첫째 항은 공변 성분 vα의 변화를, 둘째 항은 국소 기저 aα의 공간적 변화를 나타낸다. 기저 벡터가 위치에 의존하므로 단순 편미분은 좌표 변환에 대해 텐서적으로 변환되지 않는다. 텐서적 일관성을 유지하기 위해 기저 벡터의 변화를 Gauss-Weingarten 관계식으로 표현한다.

(11)
a,βα=-Γγβαaγ+bβαa3

여기서, Γγβα는 Christoffel 기호로, 메트릭 텐서 또는 기저 벡터로부터 다음 식 (12)와 같이 결정된다.

(12)
Γγβα=12aαδ(∂βaγδ+∂γaβδ-∂δaβγ)=aα·aγ,β

식 (11)을 식 (10)에 대입하면 다음 식 (13)과 같다.

(13)
∂v∂ξβ=(vα,β-Γαβγvγ)aα+vαbβαa3

첫째 항은 공변 성분의 공변 미분을 정의하며, 표면 위에서 2차 텐서로 변환된다(식 (14)).

(14)
vα|β=vα,β-Γαβγvγ

따라서 표면 위 접선 벡터장의 미분은 공변 미분을 사용하여 다음과 같이 표현된다.

(15)
∂v∂ξβ=vα|βaα+vαbβαa3

첫째 항 vα|βaα는 표면 위 내재 변화를, 둘째 항 vαbβαa3는 표면 곡률에 의한 외재 변화를 나타낸다. 마지막으로 중간면의 제3 기본형식은 다음 식 (16)과 같이 정의되며,

(16)
cαβ=bαλbλβ

위의 세 가지 기본형식 aαβ,bαβ,cαβ는 모두 대칭 2차 텐서이다.

2.2 Naghdi 쉘 모형과 변형 측도

얇은 곡면 구조의 역학적 응답을 기술하기 위해 Cosserat 이론에 기반한 Naghdi 5-parameter 쉘 모형을 채택한다. 본 모형은 Reissner-Mindlin 운동학 가정, 즉 director(법선 섬유)가 변형 전후에 직선성과 길이를 보존한다는 가정에 기반한다. 또한 director가 무한히 얇다고 가정하므로 그 축을 중심으로 한 회전(drilling)은 무시된다(Bastek and Kochmann, 2023). 이에 따라 변위장은 중간면 변위와 두께 방향 좌표에 비례하는 회전 기여의 합으로 다음과 같이 표현된다.

(17)
U(ξ1,ξ2,ξ3)=u(ξ1,ξ2)+ξ3θ(ξ1,ξ2)

여기서, 중간면 변위는 다음 식 (18)과 같이 표현된다.

(18)
u(ξ1,ξ2)=uα(ξ1,ξ2)aα(ξ1,ξ2)+u3(ξ1,ξ2)a3(ξ1,ξ2)

여기서, uα는 반변 기저 aα에 따른 접선(면내) 성분, u3는 a3에 따른 법선(면외) 성분이다. 회전 벡터는 중간면에 접하며, no-drilling 가정에 따라 다음과 같이 표현된다.

(19)
θ(ξ1,ξ2)=θα(ξ1,ξ2)aα(ξ1,ξ2)

결과적으로 미지 자유도는 u1,u2,u3,θ1,θ2의 다섯 개이다. 각 점의 변형은 중간면의 메트릭과 곡률의 기하학적 변화로 특성화되며, 이를 Green-Lagrange 변형률 텐서 εij로 나타낸다. 본 연구에서는 미소 변형 영역을 고려하므로 그 선형화 형태를 사용한다(식 (20)).

(20)
εij=12(gij*-gij)≈12(gi·U,j+gj·U,i)

여기서, gij*=gi*·gj*와 gij=gi·gj는 각각 변형 및 비변형 형상의 3차원 메트릭 텐서이다. 변형률 텐서를 다섯 매개변수 (u1,u2,u3,θ1,θ2)로 표현하기 위해, 식 (15)와 (17)을 사용하여 변위장의 편미분을 다음 식 (21)과 같이 계산한다.

(21)
∂U∂ξα=∂u∂ξα+ξ3∂θ∂ξα

식 (7)과 식 (15)를 이용하면, 중간면 변위의 편미분은 다음 식 (22)와 같다.

(22)
∂u∂ξα=∂∂ξα(uλaλ+u3a3)=uλ,αaλ+uλa,αλ+u3,αa3+u3a,α3=uλ|αaλ+bαλuλa3+u3,αa3-u3bλαaλ

회전의 편미분은 다음 식 (23)과 같다.

(23)
∂θ∂ξα=∂∂ξα(θλaλ)=θλ,αaλ+θλa,αλ=θλ|αaλ+bαλθλa3

또한 두께 좌표 𝜉3에 대한 변위장의 편미분은 다음 식 (24)와 같다.

(24)
∂U∂ξ3=θλaλ

식 (3)의 기저 벡터의 쌍대 관계와 식 (9)의 3차원 기저 gi와 중간면 기저 ai의 관계를 이용하면, 변형률 텐서는 다음 식 (25)와 같이 표현된다.

(25)
εαβ=12(uα|β+uβ|α)-bαβu3+ξ312(θα|β+θβ|α-bβλuλ|α-bαλuλ|β)+cαβu3+(ξ3)212(bβλθλ|α+bαλθλ|β)εα3=12(θα+u3,α+bαλuλ)ε33=0

쉘이 얇다는 가정 아래에서 두께 좌표 𝜉3에 대해 2차인 항은 무시된다. 결과적으로 Reissner-Mindlin 쉘 모형의 범위 안에서, 변형률은 다음 세 가지 변형 성분으로 분해된다.

(26)
eαβ=12(uα|β+uβ|α)-bαβu3
(27)
kαβ=12(θα|β+θβ|α-bβλuλ|α-bαλuλ|β)+cαβu3
(28)
γα=θα+u3,α+bαλuλ

여기서, eαβ,kαβ,γα는 각각 막 변형률, 휨 변형률, 횡전단 변형률을 나타낸다.

2.3 Duhamel-Neumann 분해와 열 변형 측도

본 연구에서는 쉘의 두께 방향으로 기준 온도 대비 변화량이 선형적으로 분포한다고 가정하고 이를 외부 입력으로 부여한다.

(29)
∆T(ξ1,ξ2,ξ3)=∆T0(ξ1,ξ2)+ξ3∆T1(ξ1,ξ2)

여기서, ∆T는 기준 온도로부터의 변화량, ∆T0는 중간면의 기준 온도 대비 변화량, ∆T1은 두께 방향 온도 구배를 나타낸다.

등방성 선형 열탄성에서 Duhamel-Neumann 관계식(Boley and Weiner, 2012)에 따르면, 임의의 재료 점에서 전체 변형률 εtotal은 탄성 변형률 εel과 열 변형률 εth의 합으로 분해되며, 등방성 열팽창에 의해 열 변형률은 다음과 같이 표현된다.

(30)
εαβtotal=εαβel+εαβth,εαβth=αT∆Taαβ

여기서, αT는 열팽창 계수이다.

식 (29)의 두께 방향 선형 분포를 식 (30)에 대입하고, 2.2절의 막 변형률 및 휨 변형률과 동일한 방식으로 두께 좌표에 대한 영차 항 및 1차 항을 분리하면, 열 변형률 또한 막 성분과 휨 성분으로 분해된다(식 (31), (32), (33)).

(31)
eαβth=αT∆T0aαβ
(32)
kαβth=αT∆T1aαβ
(33)
γαth=0

여기서, eαβth는 막 방향 열 변형률, kαβth는 휨 방향 열 변형률을 나타내며, 횡전단 방향에는 열변형 기여가 발생하지 않는다. 결과적으로 쉘의 탄성 변형률은 식 (26), (27), (28)의 전체 변형률 eαβtotal,kαβtotal,γαtotal에서 위의 열 성분을 뺀 것으로 표현된다.

(34)
eαβel=eαβtotal-αT∆T0aαβ
(35)
kαβel=kαβtotal-αT∆T1aαβ
(36)
γαel=γαtotal

2.4 탄성 에너지, 열 등가하중 및 변분 손실 함수

최소 포텐셜 에너지 원리에 따르면, 쉘의 평형 형상은 전체 포텐셜 에너지 범함수를 최소화하는 상태에 해당한다. 영률, 푸아송비, 쉘 두께를 각각 E,ν,h로 나타내면, 등방성, 균질, 선형 탄성 쉘의 탄성 에너지는 탄성 변형률(식 (34), (35), (36))에 기반하여 다음과 같이 표현된다.

(37)
Πel=12∫ωhCαβληeαβeleληel+h312Cαβληkαβelkληel+hκDαβγαelγβela1/2dξ1dξ2

여기서, Cαβλη와 Dαβ는 다음 식 (38), (39)와 같이 정의된 막/휨 및 전단 탄성 텐서이다.

(38)
Cαβλη=E2(1+ν)aαλaβη+aαηaβλ+2ν1-νaαβaλη
(39)
Dαβ=E2(1+ν)aαβ

여기서, 식 (38)의 탄성 텐서는 두께 방향 평면응력 조건 하에서 유도된 것이며, 𝜅는 두께 방향의 비균일 전단 응력 분포를 보정하는 전단 보정 계수이다(Bastek and Kochmann, 2023). 등방재료의 경우 일반적으로 κ=56가 사용되며, 본 연구에서도 이 값을 채택한다.

식 (34), (35), (36)의 탄성 변형률을 식 (37)에 대입하여 전개하면, 탄성 에너지는 전체 변형률에 대한 이차 형식 Πtotal, 전체 변형률에 선형으로 의존하는 열 변형 일 Wth, 그리고 변위장에 무관한 상수 항으로 분해된다(식 (40), (41), (42)).

(40)
Piel=Pitotal-Wth+(변위장무관상수)

여기서,

(41)
Πtotal =12∫ωhCαβληeαβtotal eληtotal +h312Cαβληkαβtotal kληtotal +hκDαβγαtotal γβtotal a1/2dξ1dξ2
(42)
Wth =∫ωCαβληhαTΔT0aληeαβtotal +h312αTΔT1aληkαβtotal a1/2dξ1dξ2

상수 항은 변위장의 변분에 영향을 주지 않으므로 본 연구에서는 무시한다. 식 (42)의 Wth 에서 막 방향 기여는 막 강성 Hαβλη=hCαβλη과 중간면 온도 변화 ΔT0의 곱으로, 휨 방향 기여는 휨 강성 Dαβλη=h312Cαβλη와 두께 방향 온도 구배 ΔT1의 곱으로 표현되며, 이러한 형태를 통칭하여 열 등가하중(thermal equivalent loads)이라 한다. 횡전단 방향에는 열 변형이 발생하지 않으므로(γαth=0) 열 등가하중 기여도 존재하지 않는다.

이로써 외력 일 Wext를 포함한 쉘의 전체 포텐셜 에너지 𝛱는 식 (41)의 Πtotal로부터 다음과 같이 표현된다.

(43)
Π[u,θ]=Πtotal -Wth -Wext 

본 연구에서 다루는 세 가지 검증 문제에는 외력이 작용하지 않으므로 Wext =0이다. 따라서 PINN 학습에 사용되는 최종 손실 함수는 다음과 같이 단순화된다.

(44)
𝓛=Πtotal [u,θ]-Wth [u,θ]

기존 약형 쉘 PINN의 손실 함수(Bastek and Kochmann, 2023)와 비교하면, 본 정식화는 식 (42)에 해당하는 두 개의 열 항만 추가로 포함한다. 이는 Duhamel-Neumann 분해가 변분 구조에 자연스럽게 결합되는 결과로, 신경망 구조나 자동미분 흐름의 어떠한 수정도 없이 열변형을 통합할 수 있음을 의미한다.

3. 물리정보 신경망 구현

본 장에서는 쉘의 열탄성 평형 형상을 근사하기 위해 사용된 물리정보 신경망(PINN) 프레임워크를 소개한다. 먼저 Deep Ritz Method의 개요를 간략히 정리한다. 이 방법에서는 피드포워드 신경망이 미지장을 표현하며, 콜로케이션 점에서 평가된 전체 포텐셜 에너지 범함수를 직접 최소화함으로써 학습된다. 다음으로 이 정식화를 Naghdi 쉘 모형에 구체화한다. 참조 중간면 위에 PINN을 구성하여 다섯 개의 일반화된 쉘 변수를 예측하고, 기하 정보 기반의 trial function을 통해 Dirichlet 경계조건을 강제하며, 쉘 변형률을 조립하기 전에 global Cartesian 좌표계의 변위를 국소 표면 좌표계로 변환한다. 자동미분을 사용하여 신경망 출력은 제2장에서 유도한 탄성 에너지 밀도와 열 등가하중 일 밀도로 사상되며, 이들의 합이 손실 함수를 정의한다. 마지막으로 네트워크 구조, 콜로케이션 점 샘플링 전략, 최적화 절차 등 주요 구현 세부 사항을 정리한다. 본 장에서 다루는 PINN 아키텍처의 전체 구조는 Fig. 2에 도식화되어 있다.

https://cdn.apub.kr/journalsite/sites/jcoseik/2026-039-05/N0040390508/images/Figure_jcoseik_39_05_08_F2.jpg
Fig. 2.

PINN architecture used in this study. The fully-connected network on the left takes the parametric coordinates (ξ1,ξ2) as input and outputs the five generalized shell variables (u1^,u2^,u3^,θ1^,θ2^).

3.1 Deep Ritz Method

본 연구에서는 Deep Ritz Method(E and Yu, 2018)를 사용하여 문제의 변분 구조를 PINN에 내장한다. 미지장은 매개변수 𝜏를 가진 피드포워드 신경망 𝒩로 근사된다. 공간 좌표 ξ∈ω⊂ℝn가 주어지면, 신경망 출력 U^τ(ξ)를 다음 식 (45)와 같이 정의한다.

(45)
U^τ(ξ):=𝓝(ξ)

이 출력은 참해(exact solution) U(ξ)의 근사(U≈Uτ^)를 제공한다.

다수의 역학 문제는 전체 포텐셜 에너지 범함수의 최소화로 정식화된다(식 (46)).

(46)
ΠU^τ=∫ω𝓕ξ,U^τ(ξ),∇U^τ(ξ)dξ

여기서, 𝓕는 에너지 밀도이다. Deep Ritz Method에서는 이 범함수를 신경망 근사에 대해 평가하여 직접 학습 목적함수로 사용한다(식 (47)).

(47)
𝓛(τ):=ΠU^τ

수치적 구현을 위해 영역 적분은 영역 𝜔에서 균일하게 샘플링된 콜로케이션 점 {ξi}i=1N에 대한 이산합으로 Monte Carlo 적분을 통해 근사된다(E and Yu, 2018). 신경망 출력과 그 미분으로 평가된 에너지 밀도를 𝓕(ξi;τ)로 표기하면, 손실은 다음 식 (48)과 같이 근사된다.

(48)
𝓛(τ)≈|ω|N∑i=1N𝓕ξi;τ

여기서, N은 콜로케이션 점의 개수이다. 𝜏에 대해 손실을 최소화하면 전체 포텐셜 에너지를 근사적으로 최소화하는 신경망 𝒩를 얻으며, 따라서 시스템의 평형 형상을 표현한다.

본 손실 함수는 강형 잔차의 norm이 아니라 전체 포텐셜 에너지 범함수 그 자체이므로, 학습 과정에서 음수 값을 가질 수 있다. 즉 loss의 절대적 크기가 아니라 최소화 여부가 수렴의 지표이다.

3.2 Naghdi 쉘 모형에의 적용

이제 일반적인 Deep Ritz 정식화를 Naghdi 쉘 모형에 구체화한다. 본 PINN 아키텍처는 Bastek과 Kochmann(2023)이 개발한 쉘 프레임워크에 기반하며, 본 연구에서는 이를 확장하여 열 등가하중 항을 동일한 에너지 기반 정식화에 추가함으로써 통합한다. 신경망은 참조 중간면 𝜔 위에 직접 구성된다. 각 곡선좌표 ξ=ξ1,ξ2∈ω에 대해 PINN 𝒩는 다섯 개의 일반화된 쉘 운동학 변수를 반환한다.

(49)
U^τ(ξ)=u1^(ξ),u2^(ξ),u3^(ξ),θ1^(ξ),θ2^(ξ)

여기서, u1^,u2^,u3^는 global Cartesian 좌표계에서의 중간면 변위 세 성분이고, θ1^,θ2^는 접선 기저 a1,a2에 대한 회전 벡터의 두 성분이다.

Dirichlet 경계조건은 기하 정보 기반의 trial function을 사용하여 신경망 표현에 직접 통합된다. 𝜕𝜔의 부분 영역 ΓD를 변위 또는 회전이 지정된 경계라 하고, UD를 그에 대응하는 Dirichlet 데이터라 하자. ΓD에서 𝜓=0이고 𝜔의 내부에서 𝜓≠0인 매끄러운 스칼라 함수 ψ:ω→ℝ를 도입한다. 그리고 ΓD에서 Up=UD가 되도록 경계 데이터를 내부로 확장하는 particular field Up:ω→ℝ5를 정의한다. 에너지 범함수에 들어가는 일반화된 쉘 장은 다음 식 (50)과 같이 표현된다.

(50)
Uτ(ξ)=Up(ξ)+ψ(ξ)U^τ(ξ)

구성에 의해 ΓD 위에서 Uτ(ξ)=Up(ξ)가 성립하므로, 학습 가능한 매개변수 𝜏의 어떠한 선택에 대해서도 Dirichlet 경계조건이 정확히 만족된다. Homogeneous Dirichlet 경계조건(UD=0)의 경우 Up=0으로 두면 일반화된 장은 다음 식 (51)과 같이 단순화된다.

(51)
Uτ(ξ)=ψ(ξ)U^τ(ξ)

본 연구의 세 가지 검증 문제에 사용한 trial function 𝜓는 Table 1에 정리되어 있다.

Table 1.

Trial functions for hard boundary condition.

Benchmark Boundary type ψ(ξ1,ξ2) Applied to
Plate Fully pinned (4 edges) ((ξ1)2-0.25)((ξ2)2-0.25)u1,u2,u3 only
Plate Cantilever (left edge encastre) ξ1+0.5 All 5 components
Cylinder Clamped (both ends) (ξ1+0.5)(ξ1-0.5) All 5 components

식 (49)에서 신경망 출력의 변위 부분 u¯=(u1^,u2^,u3^)는 global Cartesian 좌표계에서 정의되는 반면, Naghdi 쉘 운동학에 들어가는 변위는 기저 벡터 a1,a2와 법선 a3로 펼쳐진 국소 면 좌표계의 성분 u로 표현된다. 한편 회전 벡터 𝜃는 5-parameter 쉘 모형에서 쉘 법선 주위의 drilling 회전이 무시되므로 두 성분만 필요하며, 본 PINN에서는 이 제약을 강제하기 위해 출력 (θ1^,θ2^)를 처음부터 접선 기저 a1,a2에 대한 성분으로 매개화함으로써 쉘 법선 주위의 회전을 구조적으로 배제한다. 따라서 trial function을 통한 경계 조건 부여 후, Uτ의 변위 부분은 쉘 변형률을 조립하기 전에 국소 좌표계로 변환되어야 한다. 다음과 같은 선형 사상을 도입한다(식 (52)).

(52)
u(ξ)=S(ξ)u¯(ξ),S=-a1--a2--a3-

여기서, u¯는 global Cartesian 좌표계의 변위장, u는 국소 면 좌표계의 대응 변위장이다. S는 공변 기저 벡터 a1,a2,a3(각각 global Cartesian 성분으로 표현된 행 벡터)를 세로로 쌓은 3 × 3 행렬이며, Su¯는 u¯의 각 기저 방향 성분을 추출한다. 변환 행렬 S는 표면 위치(ξ1,ξ2)에 의존하므로, u의 공변 미분을 평가할 때 이 위치 의존성이 고려된다.

원통 쉘의 둘레 방향 좌표 𝜉2는 양 끝 𝜉2=±0.5가 동일한 물리 점에 대응하는 주기 좌표이므로, 신경망 입력에서 𝜉2를 다음과 같이 변환하여 주기성을 hard constraint로 보장한다.

(53)
ξ2↦(sin(2π(ξ2+0.5)),cos(2π(ξ2+0.5)))

이로써 신경망의 입력 차원이 2에서 3으로 증가하며, trial function 자체는 변환되지 않은 원래 좌표를 사용한다.

Uτ가 주어지면, 막, 휨, 횡전단 변형률 (eαβtotal,kαβtotal,γαtotal)은 (ξ1,ξ2)에 대한 자동미분을 사용하여 Naghdi 쉘 운동학으로부터 계산된다. 이 변형률과 변위장, 그리고 쉘 기하 정보를 사용하여 식 (41)의 Πtotal과 식 (42)의 Wth를 평가한다. 본 연구의 세 가지 검증 문제에는 외력이 작용하지 않으므로 쉘의 전체 포텐셜 에너지는 Π=Πtotal-Wth로 단순화되며(식 (43), (44) 참조), 식 (46)에 따라 다음 식 (54)와 같이 쓸 수 있다.

(54)
ΠUτ=∫ω𝓕(ξ;τ)dξ1dξ2

여기서, 𝓕(ξ;τ)는 중간면 매핑의 Jacobian a1/2=det(aαβ)를 포함한 단위 매개변수 면적당 전체 에너지 밀도이다. 명시적으로 ξ∈ω의 한 점에서(식 (55)),

(55)
𝓕(ξ;τ)=12hCαβληeαβtotal eληtotal +h312Cαβληkαβtotal kληtotal +hκDαβγαtotal γβtotal -hCαβληαTΔT0aληeαβtotal +h312CαβληαTΔT1aληkαβtotal ·a1/2

Deep Ritz 정식화에서 손실 함수는 전체 포텐셜 에너지의 근사값으로 정의된다. 즉 Π[Uτ]를 최소화하는 신경망 매개변수 𝜏를 구한다. 𝜔에 대한 적분을 콜로케이션 점 {ξi}i=1N⊂ω에 대한 이산합으로 근사하면 다음 식 (56)과 같다.

(56)
𝓛(τ):=ΠUτ≈|ω|N∑i=1N𝓕ξi;τ

여기서, |ω|=∫ωdξ1dξ2는 참조 매개변수 영역의 전체 면적이다.

3.3 구현 세부 사항

본 PINN 구현은 Bastek과 Kochmann(2023)이 제안한 미소 변형 Naghdi 쉘 프레임워크에 기반한다. 본 연구에서는 제2장에서 제시한 정식화에 따라 세 가지 검증 문제(fully pinned 평판, cantilever 평판, 양단 고정 원통 쉘)를 다루며, 각 문제에 대한 Dirichlet 경계조건은 Table 1의 trial function을 통해 강제된다.

다른 딥러닝 기반 솔버와 마찬가지로 다양한 하이퍼파라미터가 지정되어야 한다. 여기에는 네트워크 구조(예: 층의 수, 층당 뉴런 수, 활성화 함수)와 학습 매개변수(예: 콜로케이션 점의 수와 분포, 가중치 초기화, 다양한 손실 항의 상대 가중치)가 포함된다. 옵티마이저의 선택과 그 설정 또한 중요한 요소이다. 이러한 선택이 PINN 성능에 미치는 영향은 여러 최근 연구(Jagtap et al., 2020; Raissi et al., 2019; Zhuang et al., 2021)에서 분석되었다. 본 연구의 관찰 결과는 이러한 연구들에서 보고된 정성적 경향과 일치하므로, 여기서는 선택된 하이퍼파라미터에 대한 간략한 설명으로 한정하고, 보다 상세한 분석은 위 연구들을 참조할 수 있다.

본 연구에서 제시된 모든 결과는 표준 피드포워드 다층 퍼셉트론(MLP)을 사용하였다. 네트워크는 3개의 은닉층, 각 50개의 뉴런으로 구성된다. 활성화 함수로는 Gaussian Error Linear Unit(GELU; Hendrycks and Gimpel, 2016)를 채택하였다. GELU는 C∞ 연속성을 가져 매끄러운 함수의 근사와 고차 미분 계산에 적합하다. ReLU는 x=0에서 1계 도함수가 불연속이며 2계 도함수가 거의 모든 곳에서 0이 되어 고차 미분에 부적합하다. 한편 hyperbolic tangent(tanh)는 깊은 신경망에서 vanishing gradient 문제를 일으킨다. 이러한 이유로 두 활성화 함수는 본 약형 PINN 구현에 적합하지 않다.

콜로케이션 점은 매개변수 영역 [-0.5, 0.5]2 위에서 Sobol 시퀀스 기반의 quasi-Monte Carlo(QMC) 샘플링으로 N=4096개를 생성한다. 균일한 의사난수에 비해 Sobol 시퀀스는 낮은 불일치도(low discrepancy)를 가지므로 적분 수렴이 빠르다.

얇은 쉘의 열변형 문제에서 손실 함수의 절대값은 O(10-9) 수준으로 매우 작아, L-BFGS 옵티마이저의 Strong-Wolfe line search에서 수치 정확도 문제를 일으킬 수 있다. 본 구현에서는 손실 함수에 상수 인자 109를 곱하여 그래디언트를 안정적인 수치 범위로 끌어올린다. 이는 손실의 최소점 위치를 변경하지 않으며 옵티마이저의 line search 안정성을 향상시키는 것으로 확인되었다.

신경망은 Deep Ritz 정식화(E and Yu, 2018)를 사용하여 전체 포텐셜 에너지 범함수를 최소화함으로써 학습한다. 본 구현은 PyTorch(Paszke et al., 2019) 기반이다. 신경망 매개변수는 PyTorch의 limited-memory BFGS(L-BFGS) 솔버(Liu and Nocedal, 1989)와 Strong-Wolfe line search를 사용하여 최적화한다. PyTorch의 기본 가중치 초기화를 사용하며, 평판 문제는 500epoch, 원통 쉘은 100epoch을 학습한다. 매 epoch은 N=4096 콜로케이션 점에 대한 full-batch 평가로 구성된다. 수치 정확도를 위해 double precision (float64)을 사용한다.

4. 수치 검증

본 장에서는 제2장의 정식화와 제3장의 PINN 구현이 쉘의 열변형을 잘 포착하는지 확인하기 위해 세 가지 검증 문제의 결과를 제시한다. fully pinned 평판(두께 구배 ΔT1 부여), cantilever 평판(ΔT1 부여), 양단 고정 원통 쉘(중간면 균일 가열 ΔT0 부여)의 결과를 Abaqus 유한요소 해석과 정량 비교한다.

4.1 검증 설정

세 가지 케이스 모두에서 다음의 무차원 재료 물성과 쉘 두께 E=1,ν=0.3,αT=10-3,h=0.05,κ=5/6를 사용한다. 각 케이스의 기하 및 열하중은 Table 2에 정리되어 있다.

Table 2.

Geometry and thermal loading of the three benchmark problems.

Case Geometry Dirichlet BC ΔT0ΔT1
Fully pinned plate 1 × 1 square 4 edges, u1=u2=u3=0 0 10
Cantilever plate 1 × 1 square Left edge encastre (all 5 dof) 0 10
Cylinder R=1,L=2 Both ends encastre (all 5 dof) 10 0

FEM 참조 해는 Abaqus/Standard의 선형 정적 해석(NLGEOM OFF)으로 얻었으며, 본 PINN 구현과 동일한 재료 물성과 경계조건을 적용하였다. 모든 모델에서 열팽창의 기준 온도는 0으로 설정하였다. 평판의 두께 방향 온도 구배(ΔT1)는 Abaqus predefined temperature field의 ‘Gradient through shell section’ 옵션을 사용하여 reference magnitude 0, thickness gradient 10으로 부여하였고(ΔT0 = 0, ΔT1 = 10), 원통 쉘의 균일 가열(ΔT0)은 ‘Constant through region’ 옵션으로 magnitude 10을 부여하였다(ΔT0 = 10, ΔT1 = 0). 세 케이스 모두 4절점 감소 적분 쉘 요소(S4R)를 사용하였으며, 메쉬는 평판은 100 × 100 요소(101 × 101 노드), 원통은 둘레 100 × 축 31 요소(100 × 32 노드)로 구성하였다. 평판의 경우 메쉬를 절반으로 세분화하여도 최대 변위가 변하지 않았다. 원통 쉘의 경우 최대 반경 방향 변위가 기준 메쉬에서 1.381 × 10-2, 1/2 세분화 시 1.375 × 10-2, 1/5 세분화 시 1.374 × 10-2로, 세분화에 따른 변화량이 각각 0.43%, 0.07%로 감소하여 사용한 메쉬가 수렴 영역에 있음을 확인하였다. 정확도 정량을 위해 FEM 노드 위치에서 PINN을 평가한 후 다음과 같이 상대 L2 오차를 계산한다.

(57)
L2fIINN,fFEM=fIINN-fFEM2fFEM2

여기서, f는 각 변위 또는 회전 성분이다.

4.2 Fully pinned 평판(두께 온도 구배)

Fig. 3은 1 × 1 정사각형 평판에 두께 방향 선형 온도 변화 ΔT(ξ3)=ξ3ΔT1 (ΔT0 = 0, ΔT1 = 10)을 부여한 문제이다. 네 변에서 변위 세 성분이 0으로 고정되며 회전은 자유롭다.

https://cdn.apub.kr/journalsite/sites/jcoseik/2026-039-05/N0040390508/images/Figure_jcoseik_39_05_08_F3.jpg
Fig. 3.

Comparison between PINN and Abaqus FEM for the fully pinned plate under the thickness gradient ΔT1 = 10

ΔT1에 의해 두께 방향으로 선형 분포하는 열 변형률이 발생하며, 이를 식 (42)에 따라 적분하면 변 단위 길이당 열 모멘트 MT∝DαTΔT1이 부과된다 (D=Eh3/[12(1-ν2)]는 휨 강성). 이 열 모멘트는 두께 양 면 사이의 열팽창 차이에서 비롯되며, Fig. 3과 같이 자유 상태에서는 평판을 균일하게 휘게 하는 효과를 갖는다. 네 변에서 uz = 0이라는 강제 조건이 이 자유 휨을 억제하므로, 변에서 0이고 중앙에서 최대인 bowl 형태의 변위 분포가 형성된다. 정사각형 영역에 ξ1,ξ2 방향이 모두 동일하게 구속되어 있으므로 해는 두 방향에 대칭이다.

Fig. 3(a)~(c)와 같이 변위 uz의 상대 L2 오차는 2.27%로 3% 이내에서 일치하며, Fig. 3(j)와 같이 수렴이 잘 된 것을 확인할 수 있다. 하지만, Fig. 3(d)~(i)에서 볼 수 있듯이 회전 성분 θ1,θ2의 오차는 각각 4.16%, 4.27%로 변위보다 다소 크게 나타난다. 본 케이스의 trial function은 Table 1에 정리된 바와 같이 변위에만 곱해지며, 따라서 변위는 네 변에서 hard constraint로 정확히 강제되는 반면 회전은 trial function 없이 신경망 출력을 그대로 사용한다. 그 결과 변 위에서의 모멘트가 0이라는 Neumann 경계조건이 약형 손실의 stationarity 조건을 통해서만 간접적으로 부과된다. cantilever(좌측 변에서 회전 강제) 및 cylinder(양단에서 회전 강제) 케이스에서는 trial function이 5개 성분 모두에 곱해져 회전이 일부 변에서 hard constraint로 강제되지만, 본 케이스에서는 회전이 어느 변에서도 강제되지 않는 것이 특수하다. 이러한 차이가 변위(약 2%)에 비해 회전(약 4%)이 다소 큰 오차를 보이는 주된 원인이다. 실제로 Fig. 3(f), (i)에서 회전 오차는 네 변 경계를 따라 집중되어 나타나며, 이는 회전이 경계에서 hard constraint로 강제되지 않는다는 위 설명과 일치한다.

4.3 Cantilever 평판(두께 온도 구배)

Fig. 4는 1 × 1 정사각형 평판의 좌측 변에서 5개 자유도 (변위 3개와 회전 2개)가 모두 0으로 고정(encastre)되고, 나머지 세 변은 자유로운 cantilever 평판이다. 온도 조건으로는 ΔT1 = 10(ΔT0 = 0)이 부여된다.

https://cdn.apub.kr/journalsite/sites/jcoseik/2026-039-05/N0040390508/images/Figure_jcoseik_39_05_08_F4.jpg
Fig. 4.

Comparison between PINN and Abaqus FEM for the cantilever plate under the thickness gradient ΔT1 = 10

Fully-pinned plate와 같은 형태의 균일 열 모멘트 MT∝DαTΔT1이 부과되나, 좌측 한 변에서만 고정되고 나머지가 자유로운 비대칭 경계조건에 의해 응답 양상은 크게 다르다. 좌측 고정단에서 우측 자유단으로 갈수록 처짐 uz가 단조 감소하며, 자유단(우측 변)이 가장 크게 휘어 내려간다. 두 회전 성분 중 θ1은 director가 𝜉1 방향으로 기울어지는 양으로, 𝜉1 방향 미분이 𝜉1 방향 굽힘 곡률 𝜅11에 해당하며, 이는 cantilever의 주 굽힘 방향이다. 반면 𝜃2는 director가 𝜉2 방향으로 기울어지는 양으로, 𝜉2 방향 미분이 𝜉2 방향 굽힘 곡률 𝜅22에 해당한다. 본 케이스의 trial function 𝜓=𝜉1+0.5는 Table 1과 같이 5개 자유도 모두에 곱해지므로 좌측 고정단 𝜉1 = -0.5에서는 변위와 회전 모두 hard constraint로 정확히 만족된다. 그러나 자유 변(우측 𝜉10 = 0.5, 상하 𝜉2 = ±0.5)에서는 trial function에 의한 추가 강제가 없고, 모멘트와 응력합력이 0이어야 한다는 Neumann 경계조건이 약형 손실의 stationarity 조건을 통해 간접적으로만 만족된다. 특히 𝜉2 방향에는 좌측 고정단과 같은 hard constraint 출발점이 없으므로, 𝜃2의 분포는 상하 자유 측면의 약형 Neumann BC에 의해 결정된다.

Fig. 4(a)~(c)와 같이 uz의 L2 오차는 2.52%로 3% 이내에 들어와 정확도가 높으며, Fig. 4(j)와 같이 빠르게 수렴된 것을 확인할 수 있다. 그러나 Fig. 4(g)~(i)에서 𝜃2의 L2 오차는 9.04%로 본 연구의 세 가지 케이스 중 가장 큰 오차를 보인다. 이는 위 분석과 같이 𝜃2 분포가 자유 측면의 Neumann BC에 의해 직접 결정되는 반면, 그 BC가 본 PINN의 trial function 구조에서는 hard constraint로 부과되지 않기 때문이다. 이에 따라 Fig. 4(i)에서 𝜃2 오차는 우측 자유단과 상하 자유 변 근처에 주로 분포한다. 반면 𝜃1은 Fig. 4(d)~(f)에서 볼 수 있듯이 좌측 고정단에서 hard constraint로 0이 강제되어 𝜉1 = -0.5에서의 값이 정확히 고정되며, cantilever의 주 굽힘 방향에 해당하여 𝜉1 방향으로 매끄러운 단조 변화를 가지므로 L2 오차가 2.51%로 변위와 유사한 수준에 머문다.

4.4 Clamped 원통 쉘 (균일 가열)

Fig. 5는 반지름 R = 1, 길이 L = 2인 원통 쉘에 균일 가열 ΔT0 = 10(ΔT1 = 0)을 부여한 문제이다. 양단에서 5개 자유도 모두 0으로 고정된다.

https://cdn.apub.kr/journalsite/sites/jcoseik/2026-039-05/N0040390508/images/Figure_jcoseik_39_05_08_F5.jpg
Fig. 5.

Comparison between PINN and Abaqus FEM for the clamped- end cylindrical shell under the uniform midsurface heating ΔT0 = 10

Clamped 원통 쉘의 거동은 고전 쉘 이론에서 표준적으로 다루는 두 영역으로 나뉜다(Boley and Weiner, 2012; Timoshenko and Woinowsky-Krieger, 1959). 양단에서 멀리 떨어진 내부에서는 굽힘 변형이 작아 막 응답이 지배적이며, uradial은 평면응력의 자유 팽창에 해당하는 균일 값 wp=RαTΔT0(1+ν)에 도달한다. 한편 양단에서는 encastre 조건으로 양단의 모든 변위와 회전을 강제한다. 두 영역의 응답이 서로 불일치하므로, 양단 근방의 좁은 영역에서 굽힘 변형이 활성화되어 막 해와 양단 조건을 매끄럽게 잇는다. 이 transition 영역을 경계층(boundary layer)이라 부른다.

축대칭, 평면응력 가정 하에서 반경 방향 평형 방정식은 4차 상미분방정식 Dw''''+Eh/R2w=(Eh/R)αTΔT0(1+ν)D=Eh3/121-ν2,w=uradial ) 형태로 환원된다. 이 방정식의 전체 해는 위에서 도입한 특수해 wp와 동질 해(homogeneous solution)의 합으로 표현되며, 동질 해는 양단으로부터의 거리 𝜁에 따라 감쇠 진동하는 형태 e-βζ·[Acos(βζ)+Bsin(βζ)]를 갖는다. 여기서 β=[3(1-ν2)/(R2h2)]1/4이며, 1/𝛽는 진폭이 1/e로 감쇠하는 감쇠 길이(decay length)이다. 본 케이스 (R=1,h=0.05,ν=0.3)에서 1/𝛽≈0.174로 길이 L=2의 약 9% 수준에 해당하고, 실제 경계층 폭은 그 몇 배에 걸쳐 진폭이 무시할 수 있을 정도로 감쇠하는 영역이다. 본 케이스는 1/β≪L이므로 양단의 경계층은 서로 충분히 분리되어 있어 양단을 독립적으로 처리할 수 있으며, 한쪽 끝의 두 진폭 A, B는 그 끝의 두 BC(uradial = 0, 회전 = 0)로부터 결정된다. 그 결과 uradial은 양단에서 0에서 시작하여 경계층 안에서 wp로 단조 증가하는 것이 아니라, Fig. 5(a)~(b)와 같이 진동의 첫 peak에서 wp를 잠시 넘어선 후 감쇠 진동하며 wp에 수렴한다. 이 첫 peak가 경계층 overshoot이며, Fig. 5(a)~(b)의 uradial heatmap에서 양단 근처의 진한 띠로 관찰된다.

Fig. 5(a)~(c)에서 uradial의 L2 오차는 0.66%, Fig. 5(d)~(f)에서 uradial의 L2 오차는 1.63%로 세 케이스 중 가장 높은 변위 정확도를 보인다. Fig. 5(c), (f)에서 두 변위 성분의 오차는 모두 양단 경계층 근처에 집중되어 나타나며, 이는 굽힘 변형이 활성화되는 transition 영역과 대응한다.

회전 성분으로는 둘레 방향 축에 대한 물리적 회전각 𝜑를 평가하였다. PINN 출력 𝜃1은 식 (19)의 회전 벡터 표현 θ=θαaα에서 a1 방향의 좌표 성분이며, 본 원통 매개변수화에서 |a1|=1/L이므로 둘레 방향 축에 대한 물리적 회전각(라디안)은 φ=θ1/L로 환산된다. FEM에서 동일한 양은 다음과 같이 얻어진다. Abaqus 출력 UR1, UR2는 각각 global X축, Y축에 대한 회전 성분을 의미하며, 원주 좌표 ϑ=(2π(ξ2+0.5))에 위치한 점의 둘레 방향 단위 벡터 e。=(-sinϑ,cosϑ,0)에 회전 벡터를 투영하면 φFEM=-UR1·sinϑ+UR2·cosϑ를 얻는다.

Fig. 5(g)-(i)에서 𝜑의 L2 오차는 4.17%로, 변위장의 0.66-1.63%보다는 다소 크나, fully pinned 평판의 회전 오차(4.16-4.27%)와 비슷한 수준이다. 본 케이스는 둘레 방향에 무관한 축대칭 해를 가지며, z=L/2 단면에서 PINN이 예측한 uradial의 둘레 방향 표준편차는 1.22 × 10-5로, 같은 위치 FEM 변위 peak 대비 약 10-3 수준에 불과해 축대칭이 잘 보존됨을 확인할 수 있다. 이는 식 (53)의 주기성 변환이 효과적으로 적용됨을 시사한다.

4.5 종합 논의

Table 3에 정리한 바와 같이, 세 가지 검증 문제 모두에서 변위장의 상대 L2 오차는 3% 이내에 있고, 회전장의 상대 L2 오차는 cantilever 평판의 𝜃2를 제외하면 2~4% 수준에 머물러 본 약형 PINN 프레임워크가 Naghdi 쉘의 열변형 거동을 잘 포착함을 확인할 수 있다. 특히 원통 쉘의 경계층 overshoot까지 잘 재현한 것은 본 PINN이 쉘의 국소 변동까지 포착할 수 있음을 보여준다.

Table 3.

Relative L2 errors for the three benchmark problems.

Case L2 error of uzL2 error of 𝜃1L2 error of 𝜃2
Fully pinned plate 2.27% 4.16% 4.27%
Cantilever plate 2.52% 2.51% 9.04%
Case L2 error of uradialL2 error of uaxialL2 error of 𝜑
Cylinder 0.66% 1.63% 4.17%

회전장의 정확도는 변위장에 비해 다소 낮으며, 특히 cantilever 평판의 𝜃2에서 9.04%의 가장 큰 오차가 나타났다. 이는 자유 변에서의 Neumann 경계조건이 약형 손실 안에서만 부과되는 본 정식화의 본질적 한계에서 비롯된다. 이 한계를 완화하기 위한 향후 방향으로는 다음 세 가지를 고려할 수 있다. 첫째, 자유 변 ∂ωN에서 응력합력과 모멘트가 0이 되는 조건을 1차원 경계 적분 항으로 손실에 명시 부과하는 방법이다. 둘째, 회전장의 표현력을 강화하기 위해 신경망 구조를 개선하는 방법으로, 예를 들어 신경망의 capacity를 키우거나 adaptive activation function(Jagtap et al., 2020)을 도입할 수 있다. 셋째, 응력합력과 모멘트를 별도의 출력장으로 함께 예측하여 평형 조건을 강형으로 부과하는 mixed formulation(Harandi et al., 2024)을 적용하는 방법이다.

5. 결 론

본 연구에서는 주어진 열하중에 의한 쉘 구조의 열변형 해석을 위해 Naghdi 5-parameter 쉘 이론에 기반한 약형 물리정보 신경망 프레임워크를 제시하였다. 전체 변형률을 탄성 변형률과 열 변형률로 분해하는 Duhamel-Neumann 관계식을 쉘 에너지 범함수에 적용함으로써, 열하중의 효과를 막 및 굽힘 등가하중 항으로 약형 손실함수에 포함하였다. 이 정식화는 기존 약형 쉘 PINN의 신경망 구조를 변경하지 않고 열변형 효과를 도입할 수 있다는 점에서 간결한 확장 형태를 갖는다.

제안된 방법은 pinned 평판, cantilever 평판, 양단 고정 원통 쉘의 세 가지 문제에 적용되었으며, Abaqus 유한요소 해석 결과와 비교하여 검증하였다. 모든 예제에서 변위장의 상대 L2 오차는 3% 이내로 나타났으며, 회전장의 상대 L2 오차도 대부분 2~4% 수준으로 나타났다. 다만 cantilever 평판의 회전 성분에서는 9.04%의 상대 오차가 관찰되었으며, 이는 자유 경계에서의 회전장 예측이 변위장에 비해 상대적으로 민감할 수 있음을 보여준다. 전반적으로 제안된 약형 PINN은 주어진 열하중을 받는 쉘 구조의 변위 및 회전 응답을 유한요소 해석 결과와 비교적 잘 재현하였다.

본 연구는 열전달 방정식까지 포함하는 완전한 열-구조 연성 해석이 아니라, 주어진 온도장에 의해 유도되는 쉘 구조의 열탄성 응답을 대상으로 하였다. 또한 검증은 선형 열탄성 범위의 대표 benchmark 문제에 한정되었다. 향후 연구에서는 비균일 온도장의 공간적 복잡성, 재료 및 기하 비선형성, 다양한 경계조건을 갖는 쉘 구조로 정식화를 확장할 필요가 있다. 나아가 PINN의 물리 손실과 데이터 적합 항을 함께 사용할 수 있는 장점을 활용하며, 측정된 변형장으로부터 열하중 또는 재료 물성을 추정하는 역문제와 데이터 융합 기반 쉘 해석으로의 확장도 가능할 것으로 기대된다.

References

1

Bastek, J.H., Kochmann, D.M. (2023) Physics-Informed Neural Networks for Shell Structures, Eur. J. Mech. A/Solids, 97, 104849.

10.1016/j.euromechsol.2022.104849
2

Boley, B.A., Weiner, J.H. (2012) Theory of Thermal Stresses, John Wiley &Sons, New York, p.608.

3

Choi, C.K., Yoo, S.W. (1990) Improved Degenerated Shell Finite Elements for Analysis of Shell Structures, J. Comput. Struct. Eng. Inst. Korea, 3(1), pp.97~107.

4

E, W., Yu, B. (2018) The Deep Ritz Method: A Deep Learning-Based Numerical Algorithm for Solving Variational Problems, Commun. Math. Stat., 6(1), pp.1~12.

10.1007/s40304-018-0127-z
5

Harandi, A., Moeineddin, A., Kaliske, M., Reese, S., Rezaei, S. (2024) Mixed Formulation of Physics-Informed Neural Networks for Thermo-Mechanically Coupled Systems and Heterogeneous Domains, Int. J. Numer. Methods Eng., 125, e7388.

10.1002/nme.7388
6

Hendrycks, D., Gimpel, K. (2016) Gaussian Error Linear Units (GELUs), arXiv:1606.08415.

7

Heo, J.-J., Kim, H.-K., Kim, J.-S. (2000) The Finite Element Analysis of Shell Structures Using Improved Shell Element, J. Comput. Struct. Eng. Inst. Korea, 13(4), pp.449~459.

8

Jagtap, A.D., Kawaguchi, K., Karniadakis, G.E. (2020) Adaptive Activation Functions Accelerate Convergence in Deep and Physics-Informed Neural Networks, J. Comput. Phys., 404, 109136.

10.1016/j.jcp.2019.109136
9

Kim, M.-G., Song, Y., Lee, H., Kim, J. (2024) Computational Analysis of Geometrically Exact Shell Elements Using Multipatch IsoGeometric Analysis, J. Comput. Struct. Eng. Inst. Korea, 37(5), pp.345~352.

10.7734/COSEIK.2024.37.5.345
10

Koiter, W.T., Simmonds, J.G. (1973) Foundations of Shell Theory, in: Becker, E., Mikhailov, G.K. (Eds.), Theoretical and Applied Mechanics, Springer, Berlin, pp.150~176.

10.1007/978-3-642-65590-6_11
11

Krätzig, W.B. (1980) Introduction to General Shell Theory, in: Olszak, W. (Ed.), Thin Shell Theory: New Trends and Applications, Springer, Vienna, p.61.

10.1007/978-3-7091-2442-0_1
12

Lagaris, I.E., Likas, A., Fotiadis, D.I. (1998) Artificial Neural Networks for Solving Ordinary and Partial Differential Equations, IEEE Trans. Neural Netw., 9(5), pp.987~1000.

10.1109/72.712178
13

Li, W., Bazant, M.Z., Zhu, J. (2021) A Physics-Guided Neural Network Framework for Elastic Plates: Comparison of Governing Equations-Based and Energy-Based Approaches, Comput. Methods Appl. Mech. Eng., 383, 113933.

10.1016/j.cma.2021.113933
14

Liu, D.C., Nocedal, J. (1989) On the Limited Memory BFGS Method for Large Scale Optimization, Math. Program., 45, pp.503~528.

10.1007/BF01589116
15

Niordson, F.I. (1985) Shell Theory, Vol. 29, North-Holland Series in Applied Mathematics and Mechanics, Elsevier Science, Amsterdam, p.408.

16

Paszke, A., Gross, S., Massa, F., Lerer, A., Bradbury, J., Chanan, G., Killeen, T., Lin, Z., Gimelshein, N., Antiga, L., Desmaison, A., Köpf, A., Yang, E., DeVito, Z., Raison, M., Tejani, A., Chilamkurthy, S., Steiner, B., Fang, L., Bai, J., Chintala, S. (2019) PyTorch: An Imperative Style, High-Performance Deep Learning Library, Adv. Neural Inf. Process. Syst., 32, pp.8024~8035.

17

Raissi, M., Perdikaris, P., Karniadakis, G.E. (2019) Physics-Informed Neural Networks: A Deep Learning Framework for Solving Forward and Inverse Problems Involving Nonlinear Partial Differential Equations, J. Comput. Phys., 378, pp.686~707.

10.1016/j.jcp.2018.10.045
18

Timoshenko, S.P., Woinowsky-Krieger, S. (1959) Theory of Plates and Shells, 2nd ed., McGraw-Hill, New York, p.580.

19

Yan, C.A., Vescovini, R., Dozio, L. (2022) A Framework Based on Physics-Informed Neural Networks and Extreme Learning for the Analysis of Composite Structures, Comput. Struct., 265, 106761.

10.1016/j.compstruc.2022.106761
20

Zhuang, X., Guo, H., Alajlan, N., Zhu, H., Rabczuk, T. (2021) Deep Autoencoder Based Energy Method for the Bending, Vibration, and Buckling Analysis of Kirchhoff Plates with Transfer Learning, Eur. J. Mech. A/Solids, 87, 104225.

10.1016/j.euromechsol.2021.104225
페이지 상단으로 이동하기