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에 도식화되어 있다.
2.1 기하학적 예비 사항
쉘 운동학을 정의하기 위해 필요한 미분기하학적 개념을 정리한다. 본 절의 유도는 Niordson(1985)과 Bastek과 Kochmann(2023)을 따른다. 라틴 첨자 는 1, 2, 3의 값을, 그리스 첨자 𝛼, 𝛽, ...는 1, 2의 값을 갖는 것으로 한다. 중간면의 위치 벡터의 직교좌표 성분을 로 두면, 매개변수 영역 위의 곡선좌표 를 이용하여 같이 매개변수화할 수 있다.
여기서 는 에 대해 연속 미분 가능한 세 개의 독립적인 단일값 함수이다. 쉘의 국소 양을 기술하기 위해, 중간면 각 점에서 공변 국소 기저를 다음 식 (1)과 같이 정의한다.
여기서, “,”는 좌표 에 대한 편미분, “×”는 3차원 외적, “”는 벡터 의 노름을 나타낸다. 는 단위 법선 벡터이므로 자기 자신과 쌍대 관계()를 가진다. 공변 기저로부터 중간면의 제1 기본형식(메트릭 텐서)이 다음 식 (2)와 같이 유도된다.
여기서, “”는 3차원 내적이다. 공변과 반변 기저 사이의 관계
로부터, 반변 메트릭 는 다음 조건으로 정의된다(식 (4)).
여기서, 는 Kronecker 델타이다(Niordson, 1985). 메트릭 텐서 와 그 역 는 첨자의 올림과 내림에 사용된다. 예를 들어,
메트릭 텐서는 표면의 내재(intrinsic) 기하를 나타내며, 이를 이용하여 곡선의 길이, 곡선 사이의 각, 그리고 면적이 표현된다. 임의의 3차원 좌표계에서 중간면 위의 요소 면적 는 다음 식 (6)과 같이 표현된다.
여기서, 이다. 한편 내재 기하만으로는 3차원 공간상의 표면 형태가 유일하게 결정되지 않으며(Niordson, 1985), 표면이 접평면에 대해 어떻게 휘어 있는지는 외재(extrinsic) 기하를 나타내는 제2 기본형식 또는 곡률 텐서 로 표현된다.
3차원 공간상의 점은 중간면의 좌표 로부터 법선 방향으로 만큼 떨어진 위치에 있으며, 그 좌표를 로 표기한다(식 (8)).
여기서, 두께 는 중간면 위에서 일정하다고 가정한다. 식 (3), (5), (7)을 이용하면, 3차원 기저 벡터 와 중간면 기저 의 관계는 다음과 같이 표현된다.
곡선 표면 위 접선 벡터장 의 미분은 다음과 같다.
여기서, 첫째 항은 공변 성분 의 변화를, 둘째 항은 국소 기저 의 공간적 변화를 나타낸다. 기저 벡터가 위치에 의존하므로 단순 편미분은 좌표 변환에 대해 텐서적으로 변환되지 않는다. 텐서적 일관성을 유지하기 위해 기저 벡터의 변화를 Gauss-Weingarten 관계식으로 표현한다.
여기서, 는 Christoffel 기호로, 메트릭 텐서 또는 기저 벡터로부터 다음 식 (12)와 같이 결정된다.
식 (11)을 식 (10)에 대입하면 다음 식 (13)과 같다.
첫째 항은 공변 성분의 공변 미분을 정의하며, 표면 위에서 2차 텐서로 변환된다(식 (14)).
따라서 표면 위 접선 벡터장의 미분은 공변 미분을 사용하여 다음과 같이 표현된다.
첫째 항 는 표면 위 내재 변화를, 둘째 항 는 표면 곡률에 의한 외재 변화를 나타낸다. 마지막으로 중간면의 제3 기본형식은 다음 식 (16)과 같이 정의되며,
위의 세 가지 기본형식 는 모두 대칭 2차 텐서이다.
2.2 Naghdi 쉘 모형과 변형 측도
얇은 곡면 구조의 역학적 응답을 기술하기 위해 Cosserat 이론에 기반한 Naghdi 5-parameter 쉘 모형을 채택한다. 본 모형은 Reissner-Mindlin 운동학 가정, 즉 director(법선 섬유)가 변형 전후에 직선성과 길이를 보존한다는 가정에 기반한다. 또한 director가 무한히 얇다고 가정하므로 그 축을 중심으로 한 회전(drilling)은 무시된다(Bastek and Kochmann, 2023). 이에 따라 변위장은 중간면 변위와 두께 방향 좌표에 비례하는 회전 기여의 합으로 다음과 같이 표현된다.
여기서, 중간면 변위는 다음 식 (18)과 같이 표현된다.
여기서, 는 반변 기저 에 따른 접선(면내) 성분, 는 에 따른 법선(면외) 성분이다. 회전 벡터는 중간면에 접하며, no-drilling 가정에 따라 다음과 같이 표현된다.
결과적으로 미지 자유도는 의 다섯 개이다. 각 점의 변형은 중간면의 메트릭과 곡률의 기하학적 변화로 특성화되며, 이를 Green-Lagrange 변형률 텐서 로 나타낸다. 본 연구에서는 미소 변형 영역을 고려하므로 그 선형화 형태를 사용한다(식 (20)).
여기서, 와 는 각각 변형 및 비변형 형상의 3차원 메트릭 텐서이다. 변형률 텐서를 다섯 매개변수 로 표현하기 위해, 식 (15)와 (17)을 사용하여 변위장의 편미분을 다음 식 (21)과 같이 계산한다.
식 (7)과 식 (15)를 이용하면, 중간면 변위의 편미분은 다음 식 (22)와 같다.
회전의 편미분은 다음 식 (23)과 같다.
또한 두께 좌표 𝜉3에 대한 변위장의 편미분은 다음 식 (24)와 같다.
식 (3)의 기저 벡터의 쌍대 관계와 식 (9)의 3차원 기저 와 중간면 기저 의 관계를 이용하면, 변형률 텐서는 다음 식 (25)와 같이 표현된다.
쉘이 얇다는 가정 아래에서 두께 좌표 𝜉3에 대해 2차인 항은 무시된다. 결과적으로 Reissner-Mindlin 쉘 모형의 범위 안에서, 변형률은 다음 세 가지 변형 성분으로 분해된다.
여기서, 는 각각 막 변형률, 휨 변형률, 횡전단 변형률을 나타낸다.
2.3 Duhamel-Neumann 분해와 열 변형 측도
본 연구에서는 쉘의 두께 방향으로 기준 온도 대비 변화량이 선형적으로 분포한다고 가정하고 이를 외부 입력으로 부여한다.
여기서, 는 기준 온도로부터의 변화량, 는 중간면의 기준 온도 대비 변화량, 은 두께 방향 온도 구배를 나타낸다.
등방성 선형 열탄성에서 Duhamel-Neumann 관계식(Boley and Weiner, 2012)에 따르면, 임의의 재료 점에서 전체 변형률 은 탄성 변형률 과 열 변형률 의 합으로 분해되며, 등방성 열팽창에 의해 열 변형률은 다음과 같이 표현된다.
여기서, 는 열팽창 계수이다.
식 (29)의 두께 방향 선형 분포를 식 (30)에 대입하고, 2.2절의 막 변형률 및 휨 변형률과 동일한 방식으로 두께 좌표에 대한 영차 항 및 1차 항을 분리하면, 열 변형률 또한 막 성분과 휨 성분으로 분해된다(식 (31), (32), (33)).
여기서, 는 막 방향 열 변형률, 는 휨 방향 열 변형률을 나타내며, 횡전단 방향에는 열변형 기여가 발생하지 않는다. 결과적으로 쉘의 탄성 변형률은 식 (26), (27), (28)의 전체 변형률 에서 위의 열 성분을 뺀 것으로 표현된다.
2.4 탄성 에너지, 열 등가하중 및 변분 손실 함수
최소 포텐셜 에너지 원리에 따르면, 쉘의 평형 형상은 전체 포텐셜 에너지 범함수를 최소화하는 상태에 해당한다. 영률, 푸아송비, 쉘 두께를 각각 로 나타내면, 등방성, 균질, 선형 탄성 쉘의 탄성 에너지는 탄성 변형률(식 (34), (35), (36))에 기반하여 다음과 같이 표현된다.
여기서, 와 는 다음 식 (38), (39)와 같이 정의된 막/휨 및 전단 탄성 텐서이다.
여기서, 식 (38)의 탄성 텐서는 두께 방향 평면응력 조건 하에서 유도된 것이며, 𝜅는 두께 방향의 비균일 전단 응력 분포를 보정하는 전단 보정 계수이다(Bastek and Kochmann, 2023). 등방재료의 경우 일반적으로 가 사용되며, 본 연구에서도 이 값을 채택한다.
식 (34), (35), (36)의 탄성 변형률을 식 (37)에 대입하여 전개하면, 탄성 에너지는 전체 변형률에 대한 이차 형식 , 전체 변형률에 선형으로 의존하는 열 변형 일 , 그리고 변위장에 무관한 상수 항으로 분해된다(식 (40), (41), (42)).
여기서,
상수 항은 변위장의 변분에 영향을 주지 않으므로 본 연구에서는 무시한다. 식 (42)의 에서 막 방향 기여는 막 강성 과 중간면 온도 변화 의 곱으로, 휨 방향 기여는 휨 강성 와 두께 방향 온도 구배 의 곱으로 표현되며, 이러한 형태를 통칭하여 열 등가하중(thermal equivalent loads)이라 한다. 횡전단 방향에는 열 변형이 발생하지 않으므로() 열 등가하중 기여도 존재하지 않는다.
이로써 외력 일 를 포함한 쉘의 전체 포텐셜 에너지 𝛱는 식 (41)의 로부터 다음과 같이 표현된다.
본 연구에서 다루는 세 가지 검증 문제에는 외력이 작용하지 않으므로 이다. 따라서 PINN 학습에 사용되는 최종 손실 함수는 다음과 같이 단순화된다.
기존 약형 쉘 PINN의 손실 함수(Bastek and Kochmann, 2023)와 비교하면, 본 정식화는 식 (42)에 해당하는 두 개의 열 항만 추가로 포함한다. 이는 Duhamel-Neumann 분해가 변분 구조에 자연스럽게 결합되는 결과로, 신경망 구조나 자동미분 흐름의 어떠한 수정도 없이 열변형을 통합할 수 있음을 의미한다.
3. 물리정보 신경망 구현
본 장에서는 쉘의 열탄성 평형 형상을 근사하기 위해 사용된 물리정보 신경망(PINN) 프레임워크를 소개한다. 먼저 Deep Ritz Method의 개요를 간략히 정리한다. 이 방법에서는 피드포워드 신경망이 미지장을 표현하며, 콜로케이션 점에서 평가된 전체 포텐셜 에너지 범함수를 직접 최소화함으로써 학습된다. 다음으로 이 정식화를 Naghdi 쉘 모형에 구체화한다. 참조 중간면 위에 PINN을 구성하여 다섯 개의 일반화된 쉘 변수를 예측하고, 기하 정보 기반의 trial function을 통해 Dirichlet 경계조건을 강제하며, 쉘 변형률을 조립하기 전에 global Cartesian 좌표계의 변위를 국소 표면 좌표계로 변환한다. 자동미분을 사용하여 신경망 출력은 제2장에서 유도한 탄성 에너지 밀도와 열 등가하중 일 밀도로 사상되며, 이들의 합이 손실 함수를 정의한다. 마지막으로 네트워크 구조, 콜로케이션 점 샘플링 전략, 최적화 절차 등 주요 구현 세부 사항을 정리한다. 본 장에서 다루는 PINN 아키텍처의 전체 구조는 Fig. 2에 도식화되어 있다.
3.1 Deep Ritz Method
본 연구에서는 Deep Ritz Method(E and Yu, 2018)를 사용하여 문제의 변분 구조를 PINN에 내장한다. 미지장은 매개변수 𝜏를 가진 피드포워드 신경망 𝒩로 근사된다. 공간 좌표 가 주어지면, 신경망 출력 를 다음 식 (45)와 같이 정의한다.
이 출력은 참해(exact solution) 의 근사()를 제공한다.
다수의 역학 문제는 전체 포텐셜 에너지 범함수의 최소화로 정식화된다(식 (46)).
여기서, 𝓕는 에너지 밀도이다. Deep Ritz Method에서는 이 범함수를 신경망 근사에 대해 평가하여 직접 학습 목적함수로 사용한다(식 (47)).
수치적 구현을 위해 영역 적분은 영역 𝜔에서 균일하게 샘플링된 콜로케이션 점 에 대한 이산합으로 Monte Carlo 적분을 통해 근사된다(E and Yu, 2018). 신경망 출력과 그 미분으로 평가된 에너지 밀도를 로 표기하면, 손실은 다음 식 (48)과 같이 근사된다.
여기서, 은 콜로케이션 점의 개수이다. 𝜏에 대해 손실을 최소화하면 전체 포텐셜 에너지를 근사적으로 최소화하는 신경망 𝒩를 얻으며, 따라서 시스템의 평형 형상을 표현한다.
본 손실 함수는 강형 잔차의 norm이 아니라 전체 포텐셜 에너지 범함수 그 자체이므로, 학습 과정에서 음수 값을 가질 수 있다. 즉 loss의 절대적 크기가 아니라 최소화 여부가 수렴의 지표이다.
3.2 Naghdi 쉘 모형에의 적용
이제 일반적인 Deep Ritz 정식화를 Naghdi 쉘 모형에 구체화한다. 본 PINN 아키텍처는 Bastek과 Kochmann(2023)이 개발한 쉘 프레임워크에 기반하며, 본 연구에서는 이를 확장하여 열 등가하중 항을 동일한 에너지 기반 정식화에 추가함으로써 통합한다. 신경망은 참조 중간면 𝜔 위에 직접 구성된다. 각 곡선좌표 에 대해 PINN 𝒩는 다섯 개의 일반화된 쉘 운동학 변수를 반환한다.
여기서, 는 global Cartesian 좌표계에서의 중간면 변위 세 성분이고, 는 접선 기저 에 대한 회전 벡터의 두 성분이다.
Dirichlet 경계조건은 기하 정보 기반의 trial function을 사용하여 신경망 표현에 직접 통합된다. 𝜕𝜔의 부분 영역 를 변위 또는 회전이 지정된 경계라 하고, 를 그에 대응하는 Dirichlet 데이터라 하자. 에서 𝜓=0이고 𝜔의 내부에서 𝜓≠0인 매끄러운 스칼라 함수 를 도입한다. 그리고 에서 가 되도록 경계 데이터를 내부로 확장하는 particular field 를 정의한다. 에너지 범함수에 들어가는 일반화된 쉘 장은 다음 식 (50)과 같이 표현된다.
구성에 의해 위에서 가 성립하므로, 학습 가능한 매개변수 𝜏의 어떠한 선택에 대해서도 Dirichlet 경계조건이 정확히 만족된다. Homogeneous Dirichlet 경계조건(=0)의 경우 =0으로 두면 일반화된 장은 다음 식 (51)과 같이 단순화된다.
본 연구의 세 가지 검증 문제에 사용한 trial function 𝜓는 Table 1에 정리되어 있다.
Table 1.
Trial functions for hard boundary condition.
| Benchmark | Boundary type | Applied to | |
| Plate | Fully pinned (4 edges) | only | |
| Plate | Cantilever (left edge encastre) | All 5 components | |
| Cylinder | Clamped (both ends) | All 5 components |
식 (49)에서 신경망 출력의 변위 부분 는 global Cartesian 좌표계에서 정의되는 반면, Naghdi 쉘 운동학에 들어가는 변위는 기저 벡터 와 법선 로 펼쳐진 국소 면 좌표계의 성분 로 표현된다. 한편 회전 벡터 𝜃는 5-parameter 쉘 모형에서 쉘 법선 주위의 drilling 회전이 무시되므로 두 성분만 필요하며, 본 PINN에서는 이 제약을 강제하기 위해 출력 를 처음부터 접선 기저 에 대한 성분으로 매개화함으로써 쉘 법선 주위의 회전을 구조적으로 배제한다. 따라서 trial function을 통한 경계 조건 부여 후, 의 변위 부분은 쉘 변형률을 조립하기 전에 국소 좌표계로 변환되어야 한다. 다음과 같은 선형 사상을 도입한다(식 (52)).
여기서, 는 global Cartesian 좌표계의 변위장, 는 국소 면 좌표계의 대응 변위장이다. 는 공변 기저 벡터 (각각 global Cartesian 성분으로 표현된 행 벡터)를 세로로 쌓은 3 × 3 행렬이며, 는 의 각 기저 방향 성분을 추출한다. 변환 행렬 는 표면 위치에 의존하므로, 의 공변 미분을 평가할 때 이 위치 의존성이 고려된다.
원통 쉘의 둘레 방향 좌표 𝜉2는 양 끝 𝜉2=±0.5가 동일한 물리 점에 대응하는 주기 좌표이므로, 신경망 입력에서 𝜉2를 다음과 같이 변환하여 주기성을 hard constraint로 보장한다.
이로써 신경망의 입력 차원이 2에서 3으로 증가하며, trial function 자체는 변환되지 않은 원래 좌표를 사용한다.
가 주어지면, 막, 휨, 횡전단 변형률 은 에 대한 자동미분을 사용하여 Naghdi 쉘 운동학으로부터 계산된다. 이 변형률과 변위장, 그리고 쉘 기하 정보를 사용하여 식 (41)의 과 식 (42)의 를 평가한다. 본 연구의 세 가지 검증 문제에는 외력이 작용하지 않으므로 쉘의 전체 포텐셜 에너지는 로 단순화되며(식 (43), (44) 참조), 식 (46)에 따라 다음 식 (54)와 같이 쓸 수 있다.
여기서, 는 중간면 매핑의 Jacobian 를 포함한 단위 매개변수 면적당 전체 에너지 밀도이다. 명시적으로 의 한 점에서(식 (55)),
Deep Ritz 정식화에서 손실 함수는 전체 포텐셜 에너지의 근사값으로 정의된다. 즉 를 최소화하는 신경망 매개변수 𝜏를 구한다. 𝜔에 대한 적분을 콜로케이션 점 에 대한 이산합으로 근사하면 다음 식 (56)과 같다.
여기서, 는 참조 매개변수 영역의 전체 면적이다.
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는 연속성을 가져 매끄러운 함수의 근사와 고차 미분 계산에 적합하다. ReLU는 =0에서 1계 도함수가 불연속이며 2계 도함수가 거의 모든 곳에서 0이 되어 고차 미분에 부적합하다. 한편 hyperbolic tangent(tanh)는 깊은 신경망에서 vanishing gradient 문제를 일으킨다. 이러한 이유로 두 활성화 함수는 본 약형 PINN 구현에 적합하지 않다.
콜로케이션 점은 매개변수 영역 [-0.5, 0.5]2 위에서 Sobol 시퀀스 기반의 quasi-Monte Carlo(QMC) 샘플링으로 =4096개를 생성한다. 균일한 의사난수에 비해 Sobol 시퀀스는 낮은 불일치도(low discrepancy)를 가지므로 적분 수렴이 빠르다.
얇은 쉘의 열변형 문제에서 손실 함수의 절대값은 수준으로 매우 작아, 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은 =4096 콜로케이션 점에 대한 full-batch 평가로 구성된다. 수치 정확도를 위해 double precision (float64)을 사용한다.
4. 수치 검증
본 장에서는 제2장의 정식화와 제3장의 PINN 구현이 쉘의 열변형을 잘 포착하는지 확인하기 위해 세 가지 검증 문제의 결과를 제시한다. fully pinned 평판(두께 구배 부여), cantilever 평판( 부여), 양단 고정 원통 쉘(중간면 균일 가열 부여)의 결과를 Abaqus 유한요소 해석과 정량 비교한다.
4.1 검증 설정
세 가지 케이스 모두에서 다음의 무차원 재료 물성과 쉘 두께 를 사용한다. 각 케이스의 기하 및 열하중은 Table 2에 정리되어 있다.
Table 2.
Geometry and thermal loading of the three benchmark problems.
| Case | Geometry | Dirichlet BC | ||
| Fully pinned plate | 1 × 1 square | 4 edges, | 0 | 10 |
| Cantilever plate | 1 × 1 square | Left edge encastre (all 5 dof) | 0 | 10 |
| Cylinder | Both ends encastre (all 5 dof) | 10 | 0 |
FEM 참조 해는 Abaqus/Standard의 선형 정적 해석(NLGEOM OFF)으로 얻었으며, 본 PINN 구현과 동일한 재료 물성과 경계조건을 적용하였다. 모든 모델에서 열팽창의 기준 온도는 0으로 설정하였다. 평판의 두께 방향 온도 구배()는 Abaqus predefined temperature field의 ‘Gradient through shell section’ 옵션을 사용하여 reference magnitude 0, thickness gradient 10으로 부여하였고( = 0, = 10), 원통 쉘의 균일 가열()은 ‘Constant through region’ 옵션으로 magnitude 10을 부여하였다( = 10, = 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을 평가한 후 다음과 같이 상대 오차를 계산한다.
여기서, 는 각 변위 또는 회전 성분이다.
4.2 Fully pinned 평판(두께 온도 구배)
Fig. 3은 1 × 1 정사각형 평판에 두께 방향 선형 온도 변화 ( = 0, = 10)을 부여한 문제이다. 네 변에서 변위 세 성분이 0으로 고정되며 회전은 자유롭다.
에 의해 두께 방향으로 선형 분포하는 열 변형률이 발생하며, 이를 식 (42)에 따라 적분하면 변 단위 길이당 열 모멘트 이 부과된다 (는 휨 강성). 이 열 모멘트는 두께 양 면 사이의 열팽창 차이에서 비롯되며, Fig. 3과 같이 자유 상태에서는 평판을 균일하게 휘게 하는 효과를 갖는다. 네 변에서 = 0이라는 강제 조건이 이 자유 휨을 억제하므로, 변에서 0이고 중앙에서 최대인 bowl 형태의 변위 분포가 형성된다. 정사각형 영역에 방향이 모두 동일하게 구속되어 있으므로 해는 두 방향에 대칭이다.
Fig. 3(a)~(c)와 같이 변위 의 상대 오차는 2.27%로 3% 이내에서 일치하며, Fig. 3(j)와 같이 수렴이 잘 된 것을 확인할 수 있다. 하지만, Fig. 3(d)~(i)에서 볼 수 있듯이 회전 성분 의 오차는 각각 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 평판이다. 온도 조건으로는 = 10( = 0)이 부여된다.
Fully-pinned plate와 같은 형태의 균일 열 모멘트 이 부과되나, 좌측 한 변에서만 고정되고 나머지가 자유로운 비대칭 경계조건에 의해 응답 양상은 크게 다르다. 좌측 고정단에서 우측 자유단으로 갈수록 처짐 가 단조 감소하며, 자유단(우측 변)이 가장 크게 휘어 내려간다. 두 회전 성분 중 은 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)와 같이 의 오차는 2.52%로 3% 이내에 들어와 정확도가 높으며, Fig. 4(j)와 같이 빠르게 수렴된 것을 확인할 수 있다. 그러나 Fig. 4(g)~(i)에서 𝜃2의 오차는 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 방향으로 매끄러운 단조 변화를 가지므로 오차가 2.51%로 변위와 유사한 수준에 머문다.
4.4 Clamped 원통 쉘 (균일 가열)
Fig. 5는 반지름 = 1, 길이 = 2인 원통 쉘에 균일 가열 = 10( = 0)을 부여한 문제이다. 양단에서 5개 자유도 모두 0으로 고정된다.
Clamped 원통 쉘의 거동은 고전 쉘 이론에서 표준적으로 다루는 두 영역으로 나뉜다(Boley and Weiner, 2012; Timoshenko and Woinowsky-Krieger, 1959). 양단에서 멀리 떨어진 내부에서는 굽힘 변형이 작아 막 응답이 지배적이며, 은 평면응력의 자유 팽창에 해당하는 균일 값 에 도달한다. 한편 양단에서는 encastre 조건으로 양단의 모든 변위와 회전을 강제한다. 두 영역의 응답이 서로 불일치하므로, 양단 근방의 좁은 영역에서 굽힘 변형이 활성화되어 막 해와 양단 조건을 매끄럽게 잇는다. 이 transition 영역을 경계층(boundary layer)이라 부른다.
축대칭, 평면응력 가정 하에서 반경 방향 평형 방정식은 4차 상미분방정식 ) 형태로 환원된다. 이 방정식의 전체 해는 위에서 도입한 특수해 와 동질 해(homogeneous solution)의 합으로 표현되며, 동질 해는 양단으로부터의 거리 𝜁에 따라 감쇠 진동하는 형태 를 갖는다. 여기서 이며, 1/𝛽는 진폭이 1/e로 감쇠하는 감쇠 길이(decay length)이다. 본 케이스 ()에서 1/𝛽≈0.174로 길이 =2의 약 9% 수준에 해당하고, 실제 경계층 폭은 그 몇 배에 걸쳐 진폭이 무시할 수 있을 정도로 감쇠하는 영역이다. 본 케이스는 이므로 양단의 경계층은 서로 충분히 분리되어 있어 양단을 독립적으로 처리할 수 있으며, 한쪽 끝의 두 진폭 A, B는 그 끝의 두 BC( = 0, 회전 = 0)로부터 결정된다. 그 결과 은 양단에서 0에서 시작하여 경계층 안에서 로 단조 증가하는 것이 아니라, Fig. 5(a)~(b)와 같이 진동의 첫 peak에서 를 잠시 넘어선 후 감쇠 진동하며 에 수렴한다. 이 첫 peak가 경계층 overshoot이며, Fig. 5(a)~(b)의 heatmap에서 양단 근처의 진한 띠로 관찰된다.
Fig. 5(a)~(c)에서 의 오차는 0.66%, Fig. 5(d)~(f)에서 의 오차는 1.63%로 세 케이스 중 가장 높은 변위 정확도를 보인다. Fig. 5(c), (f)에서 두 변위 성분의 오차는 모두 양단 경계층 근처에 집중되어 나타나며, 이는 굽힘 변형이 활성화되는 transition 영역과 대응한다.
회전 성분으로는 둘레 방향 축에 대한 물리적 회전각 𝜑를 평가하였다. PINN 출력 𝜃1은 식 (19)의 회전 벡터 표현 에서 방향의 좌표 성분이며, 본 원통 매개변수화에서 이므로 둘레 방향 축에 대한 물리적 회전각(라디안)은 로 환산된다. FEM에서 동일한 양은 다음과 같이 얻어진다. Abaqus 출력 UR1, UR2는 각각 global X축, Y축에 대한 회전 성분을 의미하며, 원주 좌표 )에 위치한 점의 둘레 방향 단위 벡터 에 회전 벡터를 투영하면 를 얻는다.
Fig. 5(g)-(i)에서 𝜑의 오차는 4.17%로, 변위장의 0.66-1.63%보다는 다소 크나, fully pinned 평판의 회전 오차(4.16-4.27%)와 비슷한 수준이다. 본 케이스는 둘레 방향에 무관한 축대칭 해를 가지며, 단면에서 PINN이 예측한 의 둘레 방향 표준편차는 1.22 × 10-5로, 같은 위치 FEM 변위 peak 대비 약 10-3 수준에 불과해 축대칭이 잘 보존됨을 확인할 수 있다. 이는 식 (53)의 주기성 변환이 효과적으로 적용됨을 시사한다.
4.5 종합 논의
Table 3에 정리한 바와 같이, 세 가지 검증 문제 모두에서 변위장의 상대 오차는 3% 이내에 있고, 회전장의 상대 오차는 cantilever 평판의 𝜃2를 제외하면 2~4% 수준에 머물러 본 약형 PINN 프레임워크가 Naghdi 쉘의 열변형 거동을 잘 포착함을 확인할 수 있다. 특히 원통 쉘의 경계층 overshoot까지 잘 재현한 것은 본 PINN이 쉘의 국소 변동까지 포착할 수 있음을 보여준다.
Table 3.
Relative errors for the three benchmark problems.
| Case | error of | error of 𝜃1 | error of 𝜃2 |
| Fully pinned plate | 2.27% | 4.16% | 4.27% |
| Cantilever plate | 2.52% | 2.51% | 9.04% |
| Case | error of | error of | error of 𝜑 |
| Cylinder | 0.66% | 1.63% | 4.17% |
회전장의 정확도는 변위장에 비해 다소 낮으며, 특히 cantilever 평판의 𝜃2에서 9.04%의 가장 큰 오차가 나타났다. 이는 자유 변에서의 Neumann 경계조건이 약형 손실 안에서만 부과되는 본 정식화의 본질적 한계에서 비롯된다. 이 한계를 완화하기 위한 향후 방향으로는 다음 세 가지를 고려할 수 있다. 첫째, 자유 변 에서 응력합력과 모멘트가 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 유한요소 해석 결과와 비교하여 검증하였다. 모든 예제에서 변위장의 상대 오차는 3% 이내로 나타났으며, 회전장의 상대 오차도 대부분 2~4% 수준으로 나타났다. 다만 cantilever 평판의 회전 성분에서는 9.04%의 상대 오차가 관찰되었으며, 이는 자유 경계에서의 회전장 예측이 변위장에 비해 상대적으로 민감할 수 있음을 보여준다. 전반적으로 제안된 약형 PINN은 주어진 열하중을 받는 쉘 구조의 변위 및 회전 응답을 유한요소 해석 결과와 비교적 잘 재현하였다.
본 연구는 열전달 방정식까지 포함하는 완전한 열-구조 연성 해석이 아니라, 주어진 온도장에 의해 유도되는 쉘 구조의 열탄성 응답을 대상으로 하였다. 또한 검증은 선형 열탄성 범위의 대표 benchmark 문제에 한정되었다. 향후 연구에서는 비균일 온도장의 공간적 복잡성, 재료 및 기하 비선형성, 다양한 경계조건을 갖는 쉘 구조로 정식화를 확장할 필요가 있다. 나아가 PINN의 물리 손실과 데이터 적합 항을 함께 사용할 수 있는 장점을 활용하며, 측정된 변형장으로부터 열하중 또는 재료 물성을 추정하는 역문제와 데이터 융합 기반 쉘 해석으로의 확장도 가능할 것으로 기대된다.







