PHÁT TRIỂN PHẦN TỬ HỮU HẠN MỚI DỰA TRÊN NỘI SUY LIÊN TIẾP CHO BÀI TOÁN CƠ NHIỆT 2D VÀ 3D

Luận án Tiến sĩ Cơ Kỹ thuật: Phát triển phần tử hữu hạn mới dựa trên nội suy liên tiếp cho bài toán cơ nhiệt 2D/3D. Nghiên cứu chuyên sâu về kỹ thuật mới.

• Tác giả: Nguyễn Ngọc Minh
• Người hướng dẫn: Trương Tích Thiện
• Trường học: university of technology
• Chuyên ngành: Engineering Mechanics (Cơ kỹ thuật)
• Đề tài: Development Of New Finite Elements Based On Consecutive-interpolation For 2D And 3D Thermal-mechanical Problems
• Loại tài liệu: Phd Thesis
• Năm xuất bản: 2020
• Địa điểm: Hồ Chí Minh

Mục lục chi tiết

• LỜI CAM ĐOAN / DECLARATION
• LIST OF FIGURES
• LIST OF TABLES
• 1. CHAPTER 1: INTRODUCTION
• 1.1. Heat transfer and thermo-mechanical problems
• 1.2. Finite element method (FEM) and its issues
• 1.3. Trends in development of numerical methods
• 1.4. Original contributions of the dissertation
• 1.5. Scientific and practical meaning of the contributions by the dissertation
• 1.6. Scope of the dissertation
• 1.7. Outline of the dissertation
• 2. CHAPTER 2: LINEAR THERMO-ELASTIC PROBLEMS
• 2.1. Time integration scheme
• 2.1.1. Backward Euler scheme
• 3. CHAPTER 3: CONSECUTIVE-INTERPOLATION PROCEDURE FOR 1D AND 2D PROBLEMS
• 3.1. Issue of non-physically discontinuous nodal gradient in Finite Element Method (FEM): An example of two-node bar element (L2 element)
• 3.2. The consecutive-interpolation procedure (CIP) for two-node bar element: CL2 element
• 3.2.1. Calculation of CIP-based shape functions
• 3.2.2. First order derivative of CIP-based shape functions
• 3.2.3. Modification to retain the C0-continuity
• 3.3. The consecutive-interpolation procedure (CIP) for three-node triangular element (CT3) and four-node quadrilateral element (CQ4)
• 3.3.1. CIP-enhanced formulation for 2D domain
• 3.3.2. The CT3 element
• 3.3.3. The CQ4 element
• 4. CHAPTER 4: TWO-DIMENSIONAL DYNAMIC AND QUASI-STATIC THERMOELASTIC FRACTURE PROBLEMS IN ISOTROPIC MATERIALS
• 4.1. Numerical modeling of cracks
• 4.2. The extended consecutive-interpolation four-node quadrilateral element (XCQ4)
• 4.2.1. Enriched formulation for displacement field
• 4.2.2. Enriched formulation for temperature field
• 4.3. Computation of (dynamic) stress intensity factors (DSIFs) for thermo-elastic fracture problems
• 4.4. Crack growth modeling
• 4.5. Numerical results and discussion
• 4.5.1. Edge crack under constant flux (mode-I)
• 4.5.2. Static SIFs analysis: Square plate with a center crack
• 4.5.3. Static SIFs analysis: Rectangular plate with a slant center crack (mixed - mode)
• 4.5.4. Quasi-static crack propagation simulation of a slant edge crack in a cruciform panel
• 4.5.5. Dynamic SIFs analysis: Edge crack under quasi-static thermal shock
• 4.5.6. Dynamic SIFs analysis: Center crack under quasi-static thermal shock
• 4.5.7. Dynamic SIFs analysis: Curved crack under dynamic thermal shock
• 5. CHAPTER 5: TWO-DIMENSIONAL QUASI-STATIC THERMO-ELASTIC FRACTURE PROBLEMS IN ORTHOTROPIC MATERIALS
• 5.1. Formulation of XCQ4 element for linear thermo-elastic fracture problems in orthotropic media
• 5.1.1. The characteristic equation of an arbitrary orthotropic material
• 5.1.2. Enriched formulation for displacement
• 5.1.3. Enriched formulation for temperature
• 5.2. Evaluation of SIFs by Interaction integral
• 5.3. Crack growth modeling
• 5.3.1. Static SIFs analysis: Single edge notched specimen under mechanical tensile load
• 5.3.2. Static SIFs analysis: Rectangular epoxy/glass plate with a horizontal edge crack under constant flux
• 5.3.3. Static SIFs analysis: Rectangular epoxy/glass plate with a slant center crack
• 5.3.4. Static SIFs analysis: An anisotropic square plate with two parallel isothermal cracks
• 5.3.5. Quasi static crack propagation in an anisotropic cracked disc
• 5.3.6. Quasi static crack propagation of an edge crack in a rectangular plate under constant flux
• 5.3.7. Quasi static crack propagation in a perforated panel with a circular hole under constant heat flux
• 6. CHAPTER 6: DEVELOPMENT OF 3D CIP-BASED FINITE ELEMENTS
• 6.1. Generalized formulation to determine auxiliary functions
• 6.2. CIP-enhanced FEM for 3D linear heat transfer problems
• 6.2.1. Steady-state heat convection in a 3D complicated domain
• 6.2.2. Transient heat transfer in a plate with cylindrical hole
• 6.3. CIP-enhanced FEM for 3D linear elastic problems
• 6.3.1. Static analysis of a cantilever beam with T-shaped cross-section
• 6.3.2. Free vibration analysis of a hollow cylinder
• 6.3.3. Free vibration analysis of a composite sandwich beam
• 7. CHAPTER 7: CONCLUSIONS AND OUTLOOKS
• 7.1. XCQ4 element for analysis of linear thermo-elastic fracture problems
• 7.2. General formulation of CIP-enhanced elements

Tóm tắt

I. Tổng Quan Về Bài Toán Cơ Nhiệt và Phương Pháp FEM 55 ký tự

Bài toán cơ nhiệt là một lĩnh vực quan trọng trong kỹ thuật, xuất hiện trong nhiều ứng dụng thực tế như hàng không vũ trụ, kỹ thuật điện, cơ khí và xây dựng. Các bài toán này liên quan đến việc phân tích sự tương tác giữa trường nhiệt và trường ứng suất-biến dạng trong vật liệu. Sự thay đổi nhiệt độ có thể gây ra ứng suất nhiệt, dẫn đến biến dạng và thậm chí phá hủy cấu trúc. Do đó, việc mô phỏng và phân tích chính xác các bài toán cơ nhiệt là rất quan trọng để đảm bảo an toàn và độ tin cậy của các công trình và thiết bị. Phương pháp phần tử hữu hạn (FEM) là một công cụ mạnh mẽ để giải quyết các bài toán này. FEM chia miền giải thành các phần tử nhỏ hơn, xấp xỉ nghiệm trên mỗi phần tử và sau đó kết hợp các nghiệm này để có được nghiệm tổng thể. Phương pháp này có thể xử lý các hình dạng phức tạp, điều kiện biên khác nhau và các thuộc tính vật liệu không đồng nhất. Tuy nhiên, FEM truyền thống cũng có những hạn chế nhất định, đặc biệt là trong việc mô tả chính xác gradient trường (như gradient nhiệt độ hoặc ứng suất) tại các nút.

1.1. Giới thiệu bài toán truyền nhiệt và ứng dụng thực tế

Bài toán truyền nhiệt đóng vai trò then chốt trong vô số ứng dụng kỹ thuật, từ hệ thống điều hòa không khí phức tạp đến các thiết bị điện tử hàng ngày. Quá trình truyền nhiệt có thể xảy ra thông qua dẫn nhiệt, đối lưu, bức xạ hoặc kết hợp của các cơ chế này. Phân tích chính xác truyền nhiệt cho phép kỹ sư tối ưu hóa hiệu suất thiết bị, tăng cường độ bền và đảm bảo vận hành an toàn. Ví dụ, trong thiết kế động cơ đốt trong, hiểu rõ sự phân bố nhiệt giúp ngăn ngừa quá nhiệt cục bộ, cải thiện hiệu suất đốt cháy và kéo dài tuổi thọ động cơ. Tương tự, trong kỹ thuật điện, phân tích truyền nhiệt là yếu tố quyết định trong việc thiết kế hệ thống làm mát hiệu quả cho các linh kiện điện tử, tránh quá nhiệt và đảm bảo độ tin cậy của thiết bị. Theo [4], quá trình truyền nhiệt tuân theo định luật bảo toàn năng lượng.

1.2. Tổng quan phương pháp phần tử hữu hạn FEM và ưu nhược điểm

Phương pháp phần tử hữu hạn (FEM) là một kỹ thuật số mạnh mẽ, cho phép kỹ sư mô phỏng và phân tích hành vi của các cấu trúc và hệ thống phức tạp. Ưu điểm chính của FEM nằm ở khả năng xử lý các hình học phức tạp, vật liệu không đồng nhất và điều kiện biên đa dạng. FEM chia miền giải quyết thành nhiều phần tử hữu hạn nhỏ kết nối với nhau tại các nút. Phương pháp này xấp xỉ nghiệm trên từng phần tử bằng các hàm nội suy, sau đó lắp ráp các phương trình phần tử để tạo thành hệ phương trình toàn cục. Mặc dù rất hiệu quả, FEM truyền thống cũng có những hạn chế nhất định. Một trong những hạn chế chính là việc biểu diễn không liên tục của gradient trường (ví dụ: gradient nhiệt độ hoặc ứng suất) tại các nút, điều này có thể ảnh hưởng đến độ chính xác của kết quả, đặc biệt trong các bài toán có sự tập trung ứng suất cao.

II. Thách Thức Với FEM Độ Chính Xác và Tính Liên Tục 59 ký tự

Một trong những thách thức chính của phương pháp phần tử hữu hạn (FEM) truyền thống là tính không liên tục của trường đạo hàm (gradient) tại các nút. Điều này không phù hợp về mặt vật lý, vì các trường vật lý thực tế thường có tính liên tục. Sự không liên tục này có thể dẫn đến sai số trong kết quả, đặc biệt là trong các bài toán mà gradient trường đóng vai trò quan trọng, ví dụ như bài toán truyền nhiệt hoặc bài toán có sự tập trung ứng suất. Để khắc phục vấn đề này, nhiều phương pháp đã được phát triển, bao gồm các phương pháp có độ liên tục cao (high-order continuity) như phương pháp Đẳng hình học (Isogeometric Analysis) và phương pháp không lưới (Meshfree methods). Tuy nhiên, các phương pháp này thường phức tạp và tốn kém về mặt tính toán. Luận án này tập trung vào việc phát triển một phần tử hữu hạn mới có khả năng khắc phục vấn đề này mà vẫn duy trì tính đơn giản và hiệu quả của FEM truyền thống.

2.1. Vấn đề không liên tục gradient trong FEM truyền thống

Trong phương pháp phần tử hữu hạn truyền thống, các hàm nội suy thường chỉ đảm bảo tính liên tục của trường giá trị (ví dụ: nhiệt độ, chuyển vị) giữa các phần tử, nhưng không đảm bảo tính liên tục của trường đạo hàm (ví dụ: gradient nhiệt độ, ứng suất). Điều này có nghĩa là giá trị đạo hàm có thể khác nhau ở hai phía của một nút chung giữa hai phần tử, dẫn đến sự không liên tục tại nút. Sự không liên tục này không phù hợp với các hiện tượng vật lý thực tế, nơi mà các trường đạo hàm thường biến đổi mượt mà. Sự không liên tục gradient có thể gây ra sai số trong kết quả, đặc biệt khi phân tích các bài toán có sự tập trung gradient cao, chẳng hạn như bài toán truyền nhiệt có nguồn nhiệt cục bộ hoặc bài toán cơ học có vết nứt.

2.2. Các phương pháp khắc phục Đẳng hình học và phương pháp không lưới

Để giải quyết vấn đề không liên tục gradient trong FEM, các phương pháp tiên tiến hơn như phương pháp Đẳng hình học (Isogeometric Analysis) và phương pháp không lưới (Meshfree methods) đã được phát triển. Phương pháp Đẳng hình học sử dụng các hàm cơ sở NURBS (Non-Uniform Rational B-Splines) để mô tả cả hình học và trường nghiệm, đảm bảo tính liên tục cao hơn. Phương pháp không lưới loại bỏ hoàn toàn sự phụ thuộc vào lưới phần tử, sử dụng các hàm gần đúng dựa trên các điểm rời rạc, cho phép mô phỏng các bài toán có biến dạng lớn và tái cấu trúc. Tuy nhiên, cả hai phương pháp này đều có độ phức tạp cao hơn và đòi hỏi chi phí tính toán lớn hơn so với FEM truyền thống, gây khó khăn cho việc áp dụng rộng rãi trong thực tế.

III. Kỹ Thuật Nội Suy Liên Tiếp Giải Pháp FEM Mới 58 ký tự

Kỹ thuật nội suy liên tiếp (Consecutive-Interpolation Procedure - CIP) là một phương pháp mới được đề xuất để cải thiện độ chính xác và tính liên tục của phương pháp phần tử hữu hạn mà không làm tăng đáng kể độ phức tạp tính toán. CIP tích hợp thông tin về gradient trường tại các nút vào quá trình nội suy, giúp tạo ra trường gradient liên tục hơn so với FEM truyền thống. Phương pháp này sử dụng không chỉ giá trị tại nút mà còn giá trị đạo hàm trung bình tại nút, từ đó cải thiện độ chính xác của trường đạo hàm. Ưu điểm chính của CIP là nó vẫn duy trì thuộc tính Kronecker quan trọng trong tính toán số, sử dụng cùng một lưới phần tử với FEM truyền thống và không làm tăng số lượng bậc tự do. Luận án này tập trung vào việc hệ thống hóa và phát triển CIP để tạo ra một nhóm phần tử hữu hạn mới phù hợp với nhiều miền bài toán, từ một chiều đến ba chiều, và được áp dụng để phân tích bài toán tương tác cơ-nhiệt.

3.1. Nguyên lý hoạt động của kỹ thuật nội suy liên tiếp CIP

Kỹ thuật nội suy liên tiếp (CIP) dựa trên ý tưởng sử dụng thông tin về gradient trường tại các nút để cải thiện quá trình nội suy. Trong FEM truyền thống, các hàm nội suy chỉ sử dụng giá trị trường (ví dụ: nhiệt độ, chuyển vị) tại các nút. CIP bổ sung thêm thông tin về gradient trường tại các nút, thường được tính bằng cách trung bình gradient từ các phần tử lân cận. Thông tin gradient này được sử dụng để xây dựng các hàm nội suy mới, đảm bảo tính liên tục của gradient trường giữa các phần tử. Do đó, CIP có thể tạo ra trường gradient liên tục và chính xác hơn so với FEM truyền thống.

3.2. Ưu điểm của CIP so với các phương pháp FEM truyền thống

So với FEM truyền thống, CIP mang lại nhiều ưu điểm đáng kể. Thứ nhất, nó cải thiện độ chính xác của kết quả bằng cách giảm thiểu sự không liên tục của trường đạo hàm tại các nút. Thứ hai, CIP vẫn duy trì thuộc tính Kronecker, cho phép đơn giản hóa quá trình lắp ráp và giải hệ phương trình. Thứ ba, CIP sử dụng cùng một lưới phần tử với FEM truyền thống, giúp dễ dàng tích hợp vào các phần mềm FEM hiện có. Thứ tư, CIP không làm tăng số lượng bậc tự do, giúp giảm chi phí tính toán so với các phương pháp có độ liên tục cao khác.

3.3. Duy trì thuộc tính Kronecker và không tăng bậc tự do

Một trong những lợi thế quan trọng của kỹ thuật nội suy liên tiếp là nó duy trì thuộc tính Kronecker. Điều này có nghĩa là các hàm nội suy có giá trị bằng 1 tại nút tương ứng và bằng 0 tại các nút khác. Thuộc tính Kronecker giúp đơn giản hóa việc xây dựng ma trận độ cứng và vectơ tải trọng, đồng thời đảm bảo rằng các điều kiện biên được áp dụng một cách chính xác. Ngoài ra, CIP không làm tăng số lượng bậc tự do so với FEM truyền thống. Điều này là do CIP sử dụng thông tin gradient hiện có tại các nút, thay vì thêm các nút mới hoặc các biến số mới.

IV. Ứng Dụng XCQ4 Bài Toán Cơ Nhiệt 2D Vật Liệu Dị Hướng 57 ký tự

Luận án này trình bày việc phát triển và ứng dụng phần tử hữu hạn mới dựa trên CIP, được gọi là XCQ4, để giải quyết bài toán tương tác cơ nhiệt trong vật liệu dị hướng. XCQ4 là một phần tử tứ giác 4 nút được tăng cường bằng kỹ thuật CIP, cho phép mô phỏng chính xác các trường gradient và sự tập trung ứng suất trong vật liệu dị hướng. Phương pháp này được áp dụng để phân tích các bài toán có vết nứt, trong đó sự tập trung ứng suất và gradient nhiệt độ cao đóng vai trò quan trọng. Kết quả cho thấy rằng XCQ4 có độ chính xác cao hơn so với các phần tử FEM truyền thống trong việc dự đoán các yếu tố cường độ ứng suất và đường lan truyền vết nứt trong vật liệu dị hướng chịu tải trọng cơ nhiệt.

4.1. Phát triển phần tử XCQ4 cho bài toán cơ nhiệt 2D

Phần tử XCQ4 là một phần tử tứ giác 4 nút được phát triển dựa trên kỹ thuật nội suy liên tiếp (CIP), đặc biệt được thiết kế để giải quyết các bài toán tương tác cơ nhiệt trong không gian hai chiều. Phần tử này được tăng cường bằng các hàm nội suy bổ sung, cho phép mô tả chính xác hơn trường gradient và sự tập trung ứng suất. Công thức XCQ4 bao gồm cả trường chuyển vị và trường nhiệt độ, cho phép mô phỏng sự tương tác giữa hai trường này một cách hiệu quả.

4.2. Ứng dụng phân tích bài toán có vết nứt trong vật liệu dị hướng

Một trong những ứng dụng quan trọng của phần tử XCQ4 là phân tích các bài toán có vết nứt trong vật liệu dị hướng. Trong các bài toán này, sự tập trung ứng suất và gradient nhiệt độ cao xảy ra xung quanh đầu vết nứt. XCQ4 có khả năng mô phỏng chính xác các trường này, cho phép dự đoán chính xác các yếu tố cường độ ứng suất (SIF) và đường lan truyền vết nứt. Vật liệu dị hướng, chẳng hạn như vật liệu composite, có các thuộc tính cơ học và nhiệt khác nhau theo các hướng khác nhau, làm cho việc phân tích trở nên phức tạp hơn.

4.3. So sánh hiệu quả XCQ4 với phần tử hữu hạn truyền thống

So sánh với các phần tử hữu hạn truyền thống, XCQ4 cho thấy độ chính xác và độ tin cậy cao hơn trong việc dự đoán các yếu tố cường độ ứng suất và đường lan truyền vết nứt trong vật liệu dị hướng chịu tải trọng cơ nhiệt. Các thử nghiệm số đã chứng minh rằng XCQ4 có khả năng hội tụ nhanh hơn và cung cấp kết quả chính xác hơn, đặc biệt trong các bài toán có sự tập trung ứng suất cao. Điều này làm cho XCQ4 trở thành một công cụ hiệu quả để phân tích và thiết kế các cấu trúc làm từ vật liệu dị hướng.

V. Phát Triển Phần Tử Hữu Hạn 3D Dựa Trên Kỹ Thuật CIP 56 ký tự

Luận án tiếp tục mở rộng kỹ thuật nội suy liên tiếp (CIP) để phát triển các phần tử hữu hạn 3D mới, phù hợp với các bài toán cơ nhiệt phức tạp trong không gian ba chiều. Các phần tử 3D này được xây dựng dựa trên công thức tổng quát để xác định các hàm phụ trợ, cho phép tích hợp CIP vào các loại phần tử khác nhau, chẳng hạn như phần tử tứ diện và phần tử hexahedron. Các phần tử 3D CIP được áp dụng để phân tích các bài toán truyền nhiệt và đàn hồi tuyến tính trong các miền phức tạp, cho thấy độ chính xác và hiệu quả cao hơn so với các phần tử FEM truyền thống.

5.1. Công thức tổng quát xác định hàm phụ trợ trong CIP 3D

Việc phát triển các phần tử hữu hạn 3D dựa trên kỹ thuật nội suy liên tiếp (CIP) đòi hỏi một công thức tổng quát để xác định các hàm phụ trợ. Công thức này phải có khả năng tích hợp CIP vào các loại phần tử khác nhau, chẳng hạn như phần tử tứ diện và phần tử hexahedron. Các hàm phụ trợ này được sử dụng để kết hợp thông tin gradient tại các nút vào quá trình nội suy, đảm bảo tính liên tục của trường đạo hàm.

5.2. Ứng dụng CIP 3D cho bài toán truyền nhiệt và đàn hồi tuyến tính

Các phần tử 3D CIP được áp dụng để phân tích các bài toán truyền nhiệt và đàn hồi tuyến tính trong các miền phức tạp. Trong bài toán truyền nhiệt, CIP giúp mô phỏng chính xác sự phân bố nhiệt độ và gradient nhiệt độ trong các cấu trúc 3D. Trong bài toán đàn hồi tuyến tính, CIP giúp dự đoán chính xác ứng suất và biến dạng trong các cấu trúc chịu tải trọng khác nhau. Các thử nghiệm số đã chứng minh rằng các phần tử 3D CIP có độ chính xác và hiệu quả cao hơn so với các phần tử FEM truyền thống.

VI. Kết Luận Hướng Phát Triển Cho Phần Tử Hữu Hạn Mới 59 ký tự

Luận án này đã trình bày việc phát triển một nhóm phần tử hữu hạn mới dựa trên kỹ thuật nội suy liên tiếp (CIP), có khả năng giải quyết các bài toán cơ nhiệt 2D và 3D một cách hiệu quả. Phương pháp này cải thiện độ chính xác và tính liên tục của FEM truyền thống mà không làm tăng đáng kể độ phức tạp tính toán. Các phần tử mới, đặc biệt là XCQ4 và các phần tử 3D CIP, đã được chứng minh là có độ chính xác cao hơn trong việc dự đoán các yếu tố cường độ ứng suất và đường lan truyền vết nứt trong vật liệu dị hướng. Trong tương lai, nghiên cứu có thể được mở rộng để áp dụng CIP cho các bài toán phi tuyến, các bài toán động lực học và các bài toán có sự tương tác giữa nhiều trường vật lý.

6.1. Tóm tắt kết quả nghiên cứu và ưu điểm của phương pháp CIP

Nghiên cứu này đã thành công trong việc phát triển một phương pháp mới, kỹ thuật nội suy liên tiếp (CIP), để cải thiện độ chính xác và tính liên tục của phương pháp phần tử hữu hạn (FEM) trong việc giải quyết các bài toán cơ nhiệt. CIP đã được chứng minh là có nhiều ưu điểm so với FEM truyền thống, bao gồm độ chính xác cao hơn, tính liên tục tốt hơn, khả năng hội tụ nhanh hơn và không làm tăng đáng kể độ phức tạp tính toán.

6.2. Hướng phát triển tiếp theo cho phần tử hữu hạn dựa trên CIP

Trong tương lai, nghiên cứu có thể được mở rộng theo nhiều hướng khác nhau. Một hướng đi là áp dụng kỹ thuật nội suy liên tiếp (CIP) cho các bài toán phi tuyến, chẳng hạn như các bài toán có độ dẻo, độ nhớt và sự biến dạng lớn. Một hướng khác là phát triển các phần tử CIP cho các bài toán động lực học, cho phép mô phỏng các hiện tượng như va chạm và rung động. Ngoài ra, CIP có thể được tích hợp vào các phần tử có khả năng mô phỏng sự tương tác giữa nhiều trường vật lý, chẳng hạn như sự tương tác giữa cơ học, nhiệt học và điện từ học. Cuối cùng, việc tối ưu hóa các thuật toán CIP và phát triển các phần mềm chuyên dụng sẽ giúp CIP trở nên dễ tiếp cận hơn và được ứng dụng rộng rãi trong thực tế.

Nội dung chính

Tổng quan luận án

Hiện tượng truyền nhiệt và tương tác cơ nhiệt là các quá trình vật lý xuất hiện phổ biến trong nhiều lĩnh vực kỹ thuật như hàng không, cơ khí, xây dựng và năng lượng hạt nhân. Sự thay đổi nhiệt độ dẫn đến hiện tượng giãn nở nhiệt và sinh ra ứng suất nhiệt trong kết cấu. Khi ứng suất nhiệt vượt quá giới hạn cho phép, sự phá hủy cấu trúc có thể xảy ra, điển hình như nứt mặt đường asphalt hoặc nứt do sốc nhiệt trong các bình phản ứng áp lực của nhà máy điện hạt nhân.

Để giải quyết các bài toán cơ nhiệt mô tả bởi hệ phương trình đạo hàm riêng, Phương pháp Phần tử Hữu hạn (Finite Element Method - FEM) truyền thống đã được sử dụng rộng rãi. Tuy nhiên, FEM truyền thống bộc lộ một số hạn chế nội tại:

1. Trường đạo hàm thu được (như trường ứng suất, biến dạng hoặc gradient nhiệt độ) bị gián đoạn phi vật lý tại các nút phần tử do tính liên tục $C^0$ của các hàm dạng Lagrange.
2. Phương pháp gặp các hiện tượng khóa (locking phenomena) như khóa thể tích (volume locking) đối với vật liệu gần như không chịu nén và khóa trượt (shear locking) trong các bài toán uốn dầm, tấm, vỏ.
3. Độ chính xác của trường ứng suất và thông số phá hủy phụ thuộc lớn vào mức độ mịn của lưới; khi mô phỏng nứt phát triển hoặc biên chuyển động, việc tái chia lưới (remeshing) phức tạp và tốn kém tài nguyên tính toán.

Mục tiêu nghiên cứu của luận án được xác định rõ:

• Xây dựng và hệ thống hóa nhóm phần tử hữu hạn mới tích hợp Kỹ thuật Nội suy Liên tiếp (Consecutive-Interpolation Procedure - CIP) cho các bài toán truyền nhiệt và cơ học đàn hồi tuyến tính từ miền 1D, 2D đến 3D.
• Phát triển phần tử tứ giác bốn nút mở rộng dựa trên CIP (XCQ4) để phân tích cơ học phá hủy cơ nhiệt tuyến tính (LTEFM) ở trạng thái tĩnh, tựa tĩnh và động đối với vật liệu đẳng hướng và trực hướng.
• Đề xuất một công thức tổng quát để xác định các hàm phụ trợ cho kỹ thuật CIP, loại bỏ sự phụ thuộc vào từng loại phần tử riêng biệt và cho phép tích hợp trực tiếp vào các mã tính toán FEM sẵn có.

Đối tượng và phạm vi nghiên cứu:

• Đối tượng nghiên cứu: Các kết cấu vật rắn đàn hồi nguyên vẹn và kết cấu chứa vết nứt chịu tác dụng của tải trọng cơ học và trường nhiệt độ.
• Phạm vi vật liệu và hình học: Vật liệu đàn hồi tuyến tính tuân theo định luật Hooke, giả thiết biến dạng nhỏ và chuyển vị bé; vật liệu đẳng hướng và vật liệu trực hướng có xét đến góc định hướng trục vật liệu; các mô hình hình học từ một chiều (1D), hai chiều (2D: biến dạng phẳng, ứng suất phẳng) đến ba chiều (3D). Đối với bài toán phá hủy nứt, phạm vi tập trung vào không gian 2D.

Tổng quan tài liệu và vị trí của luận án

Để khắc phục các nhược điểm của FEM truyền thống, các nhà nghiên cứu trên thế giới đã phát triển nhiều phương pháp số khác nhau:

• Phương pháp Phân tích Đẳng hình học (Isogeometric Analysis - IGA) do Hughes và cộng sự [16] đề xuất, sử dụng các hàm B-splines và NURBS làm hàm dạng. Phương pháp này cho phép kiểm soát bậc và tính liên tục bậc cao của hàm dạng, đã được áp dụng trong phân tích tấm/vỏ [17, 18, 19], dòng thấm chưa bão hòa [20, 21] và bài toán đàn dẻo [22, 23]. Tuy nhiên, IGA có chi phí tính toán phức tạp và không thỏa mãn thuộc tính delta Kronecker, dẫn đến khó khăn khi áp đặt các điều kiện biên Dirichlet (cần dùng phương pháp nhân tử Lagrange hoặc hàm phạt).
• Các phương pháp không lưới (Meshless methods) như Element Free Galerkin (EFG) [24, 25], Reproducing Kernel Particle Method (RKPM) [26, 27], và Radial Point Interpolation Method (RPIM) [28, 29]. Nhóm phương pháp này loại bỏ lưới phần tử nhưng việc xây dựng hàm dạng phụ thuộc vào các tham số do người dùng chọn mà không có tiêu chuẩn cố định, đồng thời đa phần không thỏa mãn thuộc tính delta Kronecker.
• Phương pháp Phần tử Biên (Boundary Element Method - BEM) [30, 31], chỉ cần rời rạc hóa trên biên miền tính toán nhưng đòi hỏi phải xác định nghiệm cơ bản (fundamental solution) và gặp khó khăn khi trích xuất dữ liệu bên trong miền.
• Phương pháp Phần tử Hữu hạn Làm trơn (Smoothed FEM - S-FEM) [32] với bốn biến thể: cell-based [33, 34], node-based [35, 36], edge-based [37, 38] và face-based [39, 40]. S-FEM làm trơn trường đạo hàm nhưng đòi hỏi thêm các điểm phụ trợ không phải là nút trường.
• Kỹ thuật Nội suy Liên tiếp (CIP): Được giới thiệu ban đầu cho phần tử tam giác 3 nút (CT3) [41] và phần tử tứ giác 4 nút (CQ4) [42] trong bài toán đàn hồi 2D; mở rộng cho cơ học phá hủy đàn hồi tuyến tính [43, 44] và phần tử tấm/vỏ theo lý thuyết biến dạng trượt bậc nhất [45].

Khoảng trống nghiên cứu được luận án giải quyết: Trước thời điểm nghiên cứu của luận án, kỹ thuật CIP chỉ mới được xây dựng riêng lẻ cho từng phần tử 2D; chưa có công thức tổng quát cho miền 3D; việc áp dụng CIP vào bài toán liên kết cơ nhiệt và bài toán phá hủy do nứt nhiệt (đặc biệt là vật liệu trực hướng có xét đến góc định hướng) chưa được khảo sát đầy đủ. Luận án đã lấp đầy khoảng trống này bằng việc xây dựng hệ thống phần tử hoàn chỉnh (CL2, CT3, CQ4, XCQ4, CTH4, CHH8) và thiết lập công thức tổng quát hóa cho kỹ thuật CIP.

Cơ sở lý thuyết và phương pháp nghiên cứu

Khung lý thuyết và phương trình chi phối

Hệ phương trình chi phối bài toán cơ nhiệt tuyến tính bao gồm phương trình cân bằng động lượng và phương trình truyền nhiệt quá độ:

• Phương trình cân bằng động lượng: $\nabla \cdot \boldsymbol{\sigma} + \mathbf{b} = \rho \ddot{\mathbf{u}}$
• Định luật Hooke mở rộng có giãn nở nhiệt: $\boldsymbol{\sigma} = \mathbf{C} : (\boldsymbol{\varepsilon} - \boldsymbol{\varepsilon}^T)$, với biến dạng nhiệt $\boldsymbol{\varepsilon}^T = \boldsymbol{\alpha} \Delta T \mathbf{I}$, trong đó $\Delta T = T - T_{ref}$
• Quan hệ biến dạng - chuyển vị: $\boldsymbol{\varepsilon} = \frac{1}{2}(\nabla \mathbf{u} + \nabla \mathbf{u}^T)$
• Phương trình truyền nhiệt: $-\nabla \cdot \mathbf{q} + Q = \rho c \dot{T}$, với định luật dẫn nhiệt Fourier $\mathbf{q} = -k \nabla T$
• Điều kiện biên cơ học: Chuyển vị cho trước $\mathbf{u} = \bar{\mathbf{u}}$ trên $\Gamma_u$, lực mặt $\boldsymbol{\sigma} \cdot \mathbf{n} = \bar{\mathbf{t}}$ trên $\Gamma_t$
• Điều kiện biên nhiệt: Nhiệt độ cho trước (Dirichlet) trên $\Gamma_1$, dòng nhiệt cho trước (Neumann) trên $\Gamma_2$, đoạn nhiệt trên $\Gamma_3$ ($\mathbf{q} \cdot \mathbf{n} = 0$), và đối lưu nhiệt (Robin) trên $\Gamma_4$ ($\mathbf{q} \cdot \mathbf{n} = h(T - T_a)$)

Dạng yếu và hệ phương trình rời rạc

Áp dụng nguyên lý biến phân và định lý Green, hệ phương trình rời rạc theo không gian có dạng: $$\mathbf{M}{TT} \dot{\hat{\mathbf{T}}} + \mathbf{K}{TT} \hat{\mathbf{T}} = \mathbf{F}T$$ $$\mathbf{M}{uu} \ddot{\hat{\mathbf{u}}} + \mathbf{K}{uT} \hat{\mathbf{T}} + \mathbf{K}{uu} \hat{\mathbf{u}} = \mathbf{F}_u$$

Trong đó:

• $\mathbf{M}{TT}$, $\mathbf{K}{TT}$, $\mathbf{F}_T$ lần lượt là ma trận dung nhiệt, ma trận dẫn nhiệt và vector tải nhiệt.
• $\mathbf{M}{uu}$, $\mathbf{K}{uu}$, $\mathbf{K}_{uT}$, $\mathbf{F}_u$ lần lượt là ma trận khối lượng, ma trận độ cứng kết cấu, ma trận liên kết nhiệt - cơ và vector tải cơ học.

Sơ đồ tích phân thời gian

• Trường đạo hàm thời gian bậc nhất của nhiệt độ $\dot{T}$ được tích phân bằng sơ đồ Euler ẩn (Backward Euler scheme): $z(t + \Delta t) \approx z(t) + \Delta t \dot{z}(t + \Delta t)$, nhằm loại bỏ các dao động giả (spurious oscillations).
• Trường gia tốc $\ddot{\mathbf{u}}$ và vận tốc $\dot{\mathbf{u}}$ được tích phân bằng sơ đồ Newmark không điều kiện dừng với các tham số ổn định $\beta = 0.25$ và $\gamma = 0.5$.

Kỹ thuật nội suy liên tiếp (CIP)

Trong CIP, hàm trường $u(\mathbf{x})$ được xấp xỉ không chỉ bằng giá trị nút mà còn kết hợp với gradient trung bình có trọng số tại nút: $$u(\mathbf{x}) = \sum_{I=1}^n \left( \phi_I u^{[I]} + \phi_{Ix} u_{,x}^{[I]} + \phi_{Iy} u_{,y}^{[I]} \right) = \sum_{I=1}^n R_I(\mathbf{x}) u_I = \mathbf{R} \mathbf{u}$$

Gradient trung bình tại nút $I$ được tính toán trên tập hợp các phần tử $S_I$ kết nối với nút $I$: $$u_{,x}^{[I]} = \sum_{e \in S_I} w_e u_{,x}^{[e]}(\mathbf{x}_I)$$ Hệ số trọng số $w_e$ tương ứng là tỷ lệ chiều dài (1D), tỷ lệ diện tích (2D) hoặc tỷ lệ thể tích (3D) của phần tử $e$ so với tổng phần tử trong $S_I$.

Các hàm dạng CIP ($R_I$) giữ nguyên thuộc tính delta Kronecker ($R_I(\mathbf{x}J) = \delta{IJ}$), có miền hỗ trợ mở rộng sang các phần tử lân cận và không làm tăng tổng số bậc tự do của hệ thống trên cùng một lưới phần tử.

Nội dung chính theo từng chương

Chương 1: Giới thiệu tổng quan (Introduction)

Chương 1 trình bày tổng quan về các bài toán truyền nhiệt và tương tác cơ nhiệt trong kỹ thuật, nêu bật các ví dụ thực tế về hư hỏng do nhiệt. Tác giả phân tích nguồn gốc và cấu trúc của phương pháp phần tử hữu hạn truyền thống, đồng thời chỉ ra 4 nhược điểm cố hữu của FEM (phụ thuộc lưới, méo phần tử, gián đoạn gradient tại nút, và hiện tượng khóa). Xu hướng phát triển các phương pháp số hiện đại (IGA, Meshless, BEM, S-FEM, CIP) được tổng thuật chi tiết. Chương này xác định 4 đóng góp mới của luận án, nêu rõ ý nghĩa khoa học, phương pháp luận, phạm vi nghiên cứu và cấu trúc tổng thể của tài liệu.

Chương 2: Bài toán cơ nhiệt tuyến tính (Linear Thermo-elastic Problems)

Chương 2 thiết lập các phương trình vi phân chi phối toàn diện cho bài toán cơ nhiệt liên kết tuyến tính dựa trên các định luật nhiệt động lực học và đàn hồi học. Tác giả xây dựng dạng biến phân (weak-form) thông qua việc áp dụng tích phân từng phần và định lý Green cho các trường biến thiên. Các ma trận rời rạc $(\mathbf{M}{TT}, \mathbf{K}{TT}, \mathbf{M}{uu}, \mathbf{K}{uu}, \mathbf{K}_{uT})$ và vector tải $(\mathbf{F}_T, \mathbf{F}_u)$ được dẫn xuất chi tiết cho bài toán 2D (ứng suất phẳng và biến dạng phẳng). Sơ đồ tích phân thời gian Backward Euler và sơ đồ Newmark được tích hợp để giải quyết hệ phương trình vi phân cấp hai theo thời gian. Năng lượng biến dạng đàn hồi $\mathcal{U}(u)$ và năng lượng nhiệt tương đương $\mathcal{U}(T)$ được định nghĩa làm chuẩn đánh giá sai số.

Chương 3: Kỹ thuật nội suy liên tiếp cho bài toán 1D và 2D (Consecutive-Interpolation Procedure for 1D and 2D Problems)

Chương 3 đi sâu vào bản chất toán học của CIP:

• Phân tích hiện tượng gián đoạn đạo hàm tại nút thông qua ví dụ phần tử thanh hai nút (L2).
• Xây dựng phần tử thanh hai nút nội suy liên tiếp (CL2), thiết lập các hàm phụ trợ $\phi_i, \phi_{ix}$ và chứng minh tính trơn của hàm dạng $R_I$ cũng như đạo hàm bậc nhất của nó qua các nút.
• Đề xuất kỹ thuật khôi phục tính liên tục $C^0$ tại các mặt tiếp xúc vật liệu hoặc biên hình học bằng cách thay thế gradient trung bình bằng gradient nút cục bộ của phần tử.
• Mở rộng CIP cho phần tử tam giác ba nút (CT3) và phần tử tứ giác bốn nút (CQ4) trên miền 2D.
• Khảo sát các ví dụ số:
  1. Thanh 1D chịu lực khối: Tốc độ hội tụ của CL2 nằm giữa phần tử tuyến tính và bậc hai truyền thống.
  2. Dầm console chịu tải trọng cắt parabol: So sánh 4 loại phần tử (T3, Q4, CT3, CQ4). Kết quả cho thấy CQ4 và CT3 cho trường ứng suất pháp $\sigma_{xx}$ liên tục và trơn.
  3. Ống trụ dày chịu ứng suất nhiệt: Mô phỏng chính xác phân bố nhiệt độ và ứng suất vòng (hoop stress).

Bảng tổng hợp so sánh thời gian tính toán và độ chính xác chuẩn năng lượng giữa phần tử Q4 và CQ4 cho bài toán dầm console:

Cỡ lưới phần tử	Thời gian Q4 (s)	Sai số năng lượng Q4 (%)	Thời gian CQ4 (s)	Sai số năng lượng CQ4 (%)
$12 \times 3$	9,031	18,25	11,080	3,12
$24 \times 6$	10,120	9,84	12,450	1,45
$48 \times 12$	12,560	5,12	15,230	0,72
$96 \times 24$	18,340	2,60	22,100	0,35
$192 \times 48$	28,252	1,31	35,420	0,18

Ghi chú số liệu: Để đạt mức sai số xấp xỉ 3%, phần tử CQ4 chỉ cần lưới thô $12 \times 3$ với thời gian 11,080 giây, trong khi phần tử Q4 truyền thống phải cần đến lưới $192 \times 48$ với thời gian 28,252 giây.

Chương 4: Bài toán nứt cơ nhiệt động và tựa tĩnh 2D trong vật liệu đẳng hướng (Two-Dimensional Dynamic and Quasi-Static Thermoelastic Fracture Problems in Isotropic Materials)

Chương 4 phát triển phần tử tứ giác bốn nút mở rộng dựa trên CIP (XCQ4):

• Tích hợp hàm bước Heaviside để mô phỏng sự gián đoạn qua thân vết nứt và các hàm nhánh tiệm cận đầu vết nứt (branch functions) để mô phỏng kỳ dị ứng suất tại đỉnh nứt cho cả trường chuyển vị và trường nhiệt độ.
• Đánh giá hệ số cường độ ứng suất tĩnh ($K_I, K_{II}$) và động (DSIFs) bằng tích phân tương tác ($M$-integral).
• Mô phỏng quá trình lan truyền nứt tựa tĩnh dựa trên tiêu chuẩn ứng suất tiếp tuyến cực đại (Maximum Tangential Stress - MTS).
• Khảo sát các bài toán chuẩn:
  1. Vết nứt biên chịu dòng nhiệt không đổi (Mode-I).
  2. Tấm vuông có vết nứt tâm ở hai điều kiện biên nứt: đoạn nhiệt (adiabatic) và đẳng nhiệt (isothermal).
  3. Tấm chữ nhật có vết nứt nghiêng ở tâm chịu tải hỗn hợp (mixed-mode).
  4. Mô phỏng lan truyền nứt trong tấm hình chữ thập (cruciform panel) dưới 4 trường hợp tải trọng cơ nhiệt.
  5. Vết nứt biên và vết nứt tâm chịu sốc nhiệt tựa tĩnh.
  6. Vết nứt cong trong lớp vật liệu Bismuth chịu sốc nhiệt động, mô phỏng sóng đàn hồi lan truyền tại các thời điểm $t = 0,5,\mu\text{s}$ đến $t = 4,\mu\text{s}$.

Chương 5: Bài toán nứt cơ nhiệt tựa tĩnh 2D trong vật liệu trực hướng (Two-Dimensional Quasi-Static Thermo-Elastic Fracture Problems in Orthotropic Materials)

Chương 5 mở rộng phần tử XCQ4 cho môi trường trực hướng:

• Thiết lập phương trình đặc trưng của vật liệu trực hướng bất kỳ; tích hợp góc định hướng trục chính vật liệu $\beta$ vào các hàm làm giàu chuyển vị và nhiệt độ.
• Đánh giá hệ số cường độ ứng suất bằng tích phân tương tác mở rộng cho vật liệu dị hướng.
• Mô phỏng các bài toán nứt trong vật liệu composite và gỗ:
  1. Mẫu thử có rãnh khía biên làm từ gỗ vân sam (Picea abies) chịu kéo cơ học.
  2. Tấm composite epoxy/thủy tinh có vết nứt biên chịu dòng nhiệt không đổi.
  3. Tấm composite epoxy/thủy tinh có vết nứt nghiêng ở tâm.
  4. Tấm vuông dị hướng chứa hai vết nứt song song đẳng nhiệt.
  5. Lan truyền nứt tựa tĩnh trong đĩa nứt Brazilian dị hướng với góc nghiêng vết nứt ban đầu $\theta = 45^\circ$.
  6. Khảo sát ảnh hưởng của tỷ số mô đun đàn hồi $E_1/E_2$ và góc định hướng $\beta$ ($60^\circ$) đến quỹ đạo nứt.
  7. Tấm có lỗ tròn chịu dòng nhiệt không đổi, so sánh quỹ đạo nứt giữa vật liệu đẳng hướng và trực hướng tại các góc $\beta$ khác nhau ($0^\circ, 30^\circ, 45^\circ, 60^\circ, 90^\circ$).

Chương 6: Phát triển các phần tử hữu hạn 3D dựa trên nội suy liên tiếp (Development of 3D CIP-Based Finite Elements)

Chương 6 xây dựng công thức tổng quát cho kỹ thuật CIP trên không gian 3D:

• Thiết lập công thức tổng quát xác định các hàm phụ trợ, áp dụng cho phần tử tứ diện 4 nút (CTH4) và phần tử hình hộp lục diện 8 nút (CHH8).
• Ứng dụng trong bài toán truyền nhiệt 3D:
  1. Dẫn nhiệt - đối lưu dừng trong miền hình học 3D phức tạp (so sánh với nghiệm lưới mịn 71.540 phần tử HH8).
  2. Truyền nhiệt quá độ trong tấm có lỗ hình trụ.
• Ứng dụng trong bài toán đàn hồi tuyến tính 3D:
  1. Phân tích tĩnh dầm console có tiết diện chữ T chịu lực cắt (lưới 480 phần tử lục diện, đối chiếu nghiệm chuẩn 8.800 phần tử HH8).
  2. Phân tích dao động tự do của kết cấu vỏ trụ rỗng: so sánh 4 tần số dao động riêng không thứ nguyên đầu tiên với nghiệm giải tích.
  3. Phân tích dao động tự do của dầm sandwich composite 3 lớp (vật liệu carbon/epoxy và lõi polyurethane): xác định 5 mode dao động đầu tiên trong mặt phẳng $L-t$ và đối chứng với dữ liệu thực nghiệm.

Chương 7: Kết luận và hướng phát triển (Conclusions and Outlooks)

Chương 7 tổng kết toàn bộ các kết quả đạt được về mặt lý thuyết và tính toán số: khẳng định tính khả thi của phần tử XCQ4 trong phân tích nứt cơ nhiệt và công thức tổng quát hóa của CIP từ 1D đến 3D. Đồng thời, chương này đưa ra các khuyến nghị về hướng nghiên cứu tiếp theo.

Kết quả và những đóng góp mới

Luận án tự xác định các đóng góp mới về mặt khoa học và kỹ thuật tính toán bao gồm:

1. Phát triển phần tử XCQ4 cho bài toán nứt cơ nhiệt 2D trong vật liệu đẳng hướng: Tích hợp thành công kỹ thuật nội suy liên tiếp vào phương pháp phần tử hữu hạn mở rộng (XFEM), tạo ra phần tử XCQ4 có khả năng biểu diễn trơn trường gradient nhiệt độ và ứng suất, nâng cao độ chính xác khi xác định các hệ số cường độ ứng suất tĩnh ($K_I, K_{II}$) và động (DSIFs) theo thời gian.
2. Mở rộng mô hình cho vật liệu trực hướng: Xây dựng công thức phần tử mở rộng xét đến ảnh hưởng của góc định hướng trục chính vật liệu ($\beta$) và mức độ dị hướng ($E_1/E_2$), giải thích cơ chế lan truyền nứt nhiệt dưới tải trọng hỗn hợp trong vật liệu composite và gỗ.
3. Phát triển nhóm phần tử hữu hạn 3D dựa trên CIP: Xây dựng thành công các phần tử khối 3D gồm tứ diện 4 nút (CTH4) và lục diện 8 nút (CHH8), ứng dụng hiệu quả cho bài toán truyền nhiệt dừng, truyền nhiệt quá độ, phân tích tĩnh và dao động riêng của kết cấu composite 3D.
4. Đề xuất công thức tổng quát hóa cho kỹ thuật CIP: Thiết lập quy trình giải tích tổng quát để xác định các hàm phụ trợ cho mọi loại phần tử từ 1D đến 3D. Đóng góp này giải quyết điểm nghẽn phụ thuộc vào loại phần tử của CIP trước đây, cho phép lập trình cài đặt CIP như một mô-đun bổ sung (add-on) trong các phần mềm FEM hiện có mà không làm tăng bậc tự do hệ thống.

Hạn chế và hướng nghiên cứu tiếp

Dựa trên phạm vi và nội dung được tác giả trình bày trong văn bản:

• Hạn chế về phạm vi hình học phá hủy: Mô hình bài toán cơ học phá hủy và lan truyền nứt nhiệt mới chỉ được thiết lập và kiểm chứng trên miền hai chiều (2D), chưa mở rộng cho bài toán nứt không gian ba chiều (3D).
• Hạn chế về mô hình vật liệu: Các phân tích chỉ giới hạn trong phạm vi đàn hồi tuyến tính (định luật Hooke) và giả thiết biến dạng bé; chưa xét đến ứng xử phi tuyến vật liệu (đàn dẻo, từ biến khi nhiệt độ cao) hoặc biến dạng hình học lớn.
• Hướng nghiên cứu tiếp: Mở rộng phương pháp CIP và XCQ4 cho bài toán nứt 3D dưới tác dụng của tải trọng nhiệt - cơ; phát triển phần tử CIP cho các bài toán phi tuyến vật liệu và bài toán tấm/vỏ phức hợp; nghiên cứu tối ưu hóa thuật toán tính toán ma trận để giảm chi phí thời gian khi tính các hàm phụ trợ trên lưới quy mô lớn.

Giá trị tham khảo

Luận án là tài liệu tham khảo cho các đối tượng:

• Nghiên cứu sinh và học viên cao học chuyên ngành Cơ kỹ thuật, Cơ học tính toán: Cung cấp cơ sở lý thuyết chi tiết, các bước biến phân, dẫn xuất ma trận phần tử CIP (1D, 2D, 3D), kỹ thuật làm giàu XFEM và thuật toán tích phân tương tác đánh giá $K_I, K_{II}$.
• Giảng viên và nhà nghiên cứu phương pháp số: Cung cấp hệ thống nghiệm chuẩn giải tích, số liệu đối chuẩn thực nghiệm (như bài toán dầm sandwich composite, tấm nứt chịu tải cơ nhiệt) và phương pháp xây dựng hàm dạng thỏa mãn thuộc tính delta Kronecker.
• Kỹ sư phát triển phần mềm mô phỏng: Cung cấp thuật toán tổng quát hóa CIP để tích hợp trực tiếp vào các mã nguồn FEM thương mại hoặc mã nguồn mở (như FEAP, OOFEM) nhằm nâng cao độ chính xác của trường ứng suất mà không cần thay đổi cấu trúc bậc tự do.

Câu hỏi thường gặp

1. Kỹ thuật Nội suy Liên tiếp (CIP) khác gì so với phương pháp phần tử hữu hạn (FEM) truyền thống?
Trong FEM truyền thống, xấp xỉ trường chỉ sử dụng giá trị hàm tại nút, dẫn đến đạo hàm bậc nhất (ứng suất, biến dạng, dòng nhiệt) bị gián đoạn tại nút phần tử ($C^0$). Trong CIP, xấp xỉ trường kết hợp cả giá trị nút và gradient trung bình có trọng số tại nút. Nhờ đó, trường đạo hàm thu được liên tục và trơn qua các nút trên cùng một lưới phần tử mà không làm tăng số lượng bậc tự do.

2. Điểm khác biệt cơ bản giữa CIP và các phương pháp bậc cao khác như IGA hay phương pháp không lưới là gì?
Phương pháp IGA và hầu hết các phương pháp không lưới không thỏa mãn thuộc tính delta Kronecker ($R_I(\mathbf{x}J) \neq \delta{IJ}$), đòi hỏi các kỹ thuật phức tạp (như nhân tử Lagrange hoặc hàm phạt) để áp đặt điều kiện biên Dirichlet. Ngược lại, hàm dạng CIP vẫn bảo toàn nguyên vẹn thuộc tính delta Kronecker, giúp việc áp đặt điều kiện biên chuyển vị và nhiệt độ diễn ra trực tiếp như trong FEM truyền thống.

3. Tại sao phần tử CQ4 lại đạt hiệu quả thời gian tính toán tốt hơn Q4 khi muốn đạt cùng mức độ chính xác?
Mặc dù trên cùng một lưới phần tử, CQ4 tốn nhiều thời gian tính toán hơn Q4 cho một bước lặp do phải tính các hàm phụ trợ, nhưng tốc độ hội tụ của CQ4 cao hơn nhiều. Để đạt cùng mức sai số năng lượng xấp xỉ 3% trong bài toán dầm chịu uốn, CQ4 chỉ cần lưới $12 \times 3$ (mất 11,080 giây), trong khi Q4 phải chia lưới mịn tới $192 \times 48$ (mất 28,252 giây).

4. Kỹ thuật CIP xử lý các vị trí có biên hình học hoặc tiếp giáp giữa hai vật liệu khác nhau (nơi gradient bị gián đoạn vật lý) như thế nào?
Để khôi phục tính liên tục $C^0$ tại các vị trí này, công thức CIP được điều chỉnh bằng cách thay thế "gradient trung bình tại nút" bằng chính "gradient nút cục bộ của phần tử". Khi áp dụng điều chỉnh này cho toàn bộ miền tính toán, phần tử CIP sẽ suy biến hoàn toàn về phần tử FEM truyền thống.

5. Luận án đã sử dụng những sơ đồ tích phân thời gian nào để giải bài toán cơ nhiệt quá độ và động lực học?
Luận án sử dụng sơ đồ Euler ẩn (Backward Euler scheme) cho phương trình truyền nhiệt để tính đạo hàm cấp một của nhiệt độ theo thời gian $\dot{T}$ (nhằm tránh dao động giả), và sử dụng sơ đồ Newmark không điều kiện dừng (với $\beta = 0,25$ và $\gamma = 0,5$) cho phương trình chuyển động để tính đạo hàm cấp hai của chuyển vị $\ddot{\mathbf{u}}$.

Kết luận

Luận án của tác giả Nguyễn Ngọc Minh đã giải quyết trọn vẹn việc phát triển và hệ thống hóa nhóm phần tử hữu hạn mới dựa trên kỹ thuật nội suy liên tiếp (CIP) áp dụng cho các bài toán cơ nhiệt từ 1D, 2D đến 3D. Các mô hình phần tử mới (CL2, CT3, CQ4, XCQ4, CTH4, CHH8) đã khắc phục được nhược điểm gián đoạn phi vật lý của trường đạo hàm trong FEM truyền thống, duy trì thuộc tính delta Kronecker và nâng cao độ chính xác cũng như hiệu quả tính toán. Công thức tổng quát hóa được đề xuất mở ra khả năng tích hợp kỹ thuật CIP vào các hệ thống phần mềm tính toán kỹ thuật hiện hành.

Trích đoạn nội dung tài liệu

VIETNAM NATIONAL UNIVERSITY HO CHI MINH CITY UNIVERSITY OF TECHNOLOGY NGUYỄN NGỌC MINH DEVELOPMENT OF NEW FINITE ELEMENTS BASED ON CONSECUTIVE-INTERPOLATION FOR 2D AND 3D THERMAL-MECHANICAL PROBLEMS PHD THESIS IN ENGINEERING HO CHI MINH CITY - 2020 VIETNAM NATIONAL UNIVERSITY HO CHI MINH CITY UNIVERSITY OF TECHNOLOGY NGUYỄN NGỌC MINH DEVELOPMENT OF NEW FINITE ELEMENTS BASED ON CONSECUTIVE-INTERPOLATION FOR 2D AND 3D THERMAL-MECHANICAL PROBLEMS Major: Engineering Mechanics (Cơ kỹ thuật) Codes: 62520101 Independent Examiner 1: Associate Prof. NGUYỄN VĂN HIẾU Independent Examiner 2: Associate Prof. NGUYỄN QUỐC HƯNG Examiner 1: Dr. TRỊNH ANH NGỌC Examiner 2: Associate Prof.

CHÂU ĐÌNH THÀNH Examiner 3: Associate Prof. LƯƠNG VĂN HẢI SCIENTIFIC SUPERVISORS: 1. TRƯƠNG TÍCH THIỆN 2. BÙI QUỐC TÍNH LỜI CAM ĐOAN / DECLARATION Tác giả xin cam đoan đây là công trình nghiên cứu của bản thân tác giả.

Các kết quả nghiên cứu và các kết luận trong luận án này là trung thực, và không sao chép từ bất kỳ một nguồn nào và dưới bất kỳ hình thức nào. Việc tham khảo các nguồn tài liệu (nếu có) đã được thực hiện trích dẫn và ghi nguồn tài liệu tham khảo đúng quy định. This doctoral dissertation is the outcome of my original research, conducted at the Ho Chi Minh city University of Technology, VNU – HCM, Viet Nam. I declare that this document is my own work and has not submitted for any other degree or qualification except as specified.

All the references are cited in the document. Tác giả luận án / Author Chữ ký / Signature i TÓM TẮT LUẬN ÁN Luận án trình bày một nhóm phần tử hữu hạn mới dựa trên sự tích hợp kỹ thuật nội suy liên tiếp vào phương pháp phần tử hữu hạn truyền thống. Với kỹ thuật này, không chỉ giá trị tại nút mà cả giá trị đạo hàm trung bình tại nút cũng được sử dụng trong quá trình xấp xỉ. Nhờ đó trường đạo hàm thu được từ nhóm phần tử hữu hạn mới là một trường liên tục, thay vì bất liên tục tại nút (không hợp lý về mặt vật lý) như phương pháp truyền thống.

Sự cải thiện về tính liên tục cũng mang đến độ chính xác cao hơn của lời giải xấp xỉ. Tuy nhiên khác với các phương pháp có độ liên tục bậc cao hiện hành như phương pháp Đẳng hình học hay phương pháp không lưới, phương pháp phần tử hữu hạn nội suy liên tiếp vẫn duy trì thuộc tính Kronecker quan trọng trong tính toán số. Thêm vào đó, phương pháp đề xuất sử dụng cùng một lưới phần tử với phương pháp phần tử hữu hạn truyền thống và không làm tăng số lượng bậc tự do. Kỹ thuật nội suy liên tiếp ban đầu được giới thiệu riêng lẻ cho phần tử tam giác 3 nút và phần tử tứ giác 4 nút và áp dụng với bài toán đàn hồi tuyến tính hai chiều.

Trong luận án này, phương pháp được hệ thống hóa và phát triển nâng cao để tạo ra một nhóm phần tử hữu hạn mới phù hợp với nhiều miền bài toán từ một chiều đến ba chiều, và được áp dụng để phân tích bài toán tương tác cơ-nhiệt. Phương pháp đề xuất tiếp tục được mở rộng để khảo sát ứng xử của miền chứa dạng bất liên tục như vết nứt, với vật liệu đẳng hướng và vật liệu trực hướng. ii ABSTRACT This dissertation presents a new group of finite elements based on the integration of consecutive-interpolation procedure into the traditional Finite Element Method. With this technique, not only the nodal values but also the averaged nodal gradients are included in approximation process.

As a result, the gradient fields obtained by the new group of finite elements are smooth, unlike the fields which are (non-physically) discontinuous at nodes delivered by traditional Finite Element Method. The improvement on continuity results in higher accuracy of approximated solution as well. On the other hand, unlike the other higher-order methods such as the Isogeometric Analysis and the Meshfree methods, the novel Consecutive-interpolation Finite Element Method possesses the important Kronecker-delta property. Furthermore, the proposed method employs the same discretization mesh with the traditional Finite Element Method and does not increase the number of degrees of freedom.

The consecutive-interpolation procedure was initially introduced separately for the 3- node triangular element and the 4-node quadrilateral element to be used in analysis of two-dimensional linear elastic problems. In this dissertation, the method is further developed to form a new class of finite elements which is suitable for domains from 1D to 3D. The new group of finite elements (being integrated with consecutive-interpolation procedure) is applied to analyze the thermo-mechanical problems. The proposed method is also extended to study behaviors of bodies containing discontinuities such as cracks, for both isotropic and orthotropic materials.

iii CONTENTS LỜI CAM ĐOAN / DECLARATION .iv LIST OF FIGURES. viii LIST OF TABLES. xiii CHAPTER 1 INTRODUCTION .1 Heat transfer and thermo-mechanical problems .1 Finite element method (FEM) and its issues. 2 Trends in development of numerical methods .5 Original contributions of the dissertation.

7 Scientific and practical meaning of the contributions by the dissertation .8 Scope of the dissertation .8 Outline of the dissertation .9 CHAPTER 2 LINEAR THERMO-ELASTIC PROBLEMS .13 Time integration scheme .1 Backward Euler scheme. 16 CHAPTER 3 CONSECUTIVE-INTERPOLATION PROCEDURE FOR 1D AND 2D PROBLEMS .17 Issue of non-physically discontinuous nodal gradient in Finite Element Method (FEM): An example of two-node bar element (L2 element). 17 The consecutive-interpolation procedure (CIP) for two-node bar element: CL2 element.1 Calculation of CIP-based shape functions .2 First order derivative of CIP-based shape functions .3 Modification to retain the C0-continuity. 25 The consecutive-interpolation procedure (CIP) for three-node triangular element (CT3) and four-node quadrilateral element (CQ4) .1 CIP-enhanced formulation for 2D domain .2 The CT3 element .3 The CQ4 element .37 CHAPTER 4 TWO-DIMENSIONAL DYNAMIC AND QUASI-STATIC THERMOELASTIC FRACTURE PROBLEMS IN ISOTROPIC MATERIALS .38 Numerical modeling of cracks .39 The extended consecutive-interpolation four-node quadrilateral element (XCQ4) .1 Enriched formulation for displacement field .2 Enriched formulation for temperature field .43 Computation of (dynamic) stress intensity factors (DSIFs) for thermo-elastic fracture problems .45 Crack growth modeling.

48 Numerical results and discussion .1 Edge crack under constant flux (mode-I) .2 Static SIFs analysis: Square plate with a center crack .3 Static SIFs analysis: Rectangular plate with a slant center crack (mixed - mode) .4 Quasi-static crack propagation simulation of a slant edge crack in a cruciform panel .5 Dynamic SIFs analysis: Edge crack under quasi-static thermal shock .6 Dynamic SIFs analysis: Center crack under quasi-static thermal shock .7 Dynamic SIFs analysis: Curved crack under dynamic thermal shock .76 CHAPTER 5 TWO-DIMENSIONAL QUASI-STATIC THERMO-ELASTIC FRACTURE PROBLEMS IN ORTHOTROPIC MATERIALS .78 v Formulation of XCQ4 element for linear thermo-elastic fracture problems in orthotropic media .1 The characteristic equation of an arbitrary orthotropic material .2 Enriched formulation for displacement .3 Enriched formulation for temperature. 83 Evaluation of SIFs by Interaction integral. 83 Crack growth modeling.1 Static SIFs analysis: Single edge notched specimen under mechanical tensile load .2 Static SIFs analysis: Rectangular epoxy/glass plate with a horizontal edge crack under constant flux.3 Static SIFs analysis: Rectangular epoxy/glass plate with a slant center crack .4 Static SIFs analysis: An anisotropic square plate with two parallel isothermal cracks .5 Quasi static crack propagation in an anisotropic cracked disc .6 Quasi static crack propagation of an edge crack in a rectangular plate under constant flux .7 Quasi static crack propagation in a perforated panel with a circular hole under constant heat flux.106 CHAPTER 6 DEVELOPMENT OF 3D CIP-BASED FINITE ELEMENTS .107 Generalized formulation to determine auxiliary functions .107 CIP-enhanced FEM for 3D linear heat transfer problems .1 Steady-state heat convection in a 3D complicated domain.2 Transient heat transfer in a plate with cylindrical hole .113 CIP-enhanced FEM for 3D linear elastic problems .1 Static analysis of a cantilever beam with T-shaped cross-section .2 Free vibration analysis of a hollow cylinder .3 Free vibration analysis of a composite sandwich beam .124 CHAPTER 7 CONCLUSIONS AND OUTLOOKS .1 XCQ4 element for analysis of linear thermo-elastic fracture problems .2 General formulation of CIP-enhanced elements .144 vii LIST OF FIGURES Figure 1. Thermal crack in an asphalt pavement [7].

Illustration of a finite element mesh for a two-dimensional heat transfer problem. Convergence of thermal energy corresponding to mesh size. Illustration of a two-node bar element (L2) in global system of coordinates (left) and in natural coordinates (right). A one-dimensional domain being discretized by two L2 elements.

Shape function associated with node 2 (global) of the Example provided in Figure 3.2, computed by traditional FEM and by CIP-enhanced FEM. CIP-based shape functions R1, R2 and R3 associated with the nodes given in Figure 3. First-order derivative of the shape functions associated with node 2 of the Example given in Figure 3.2, computed by traditional FEM and CIP-enhanced FEM .4: Sketch of 1D bar subjected to body load. The bar is uniformly discretized by 5 eleements .4: Comparison of stress calculated by standard FEM and the proposed CIP-enhanced FEM.4: Convergence rate with respect to number of degrees of freedom (DOFs) of the standard linear FEM, quadratic FEM and CIP approach.

Visualization of the shape function for a CQ4 element. Visualization of the first order derivative of shape function for a CQ4 element. Sketch of the CIP-enhanced four-node quadrilateral element (CQ4) in a 2D finite element mesh [42] .1: Cantilever beam being subject to parabolic shear load31 Figure 3.1: (a) The triangular mesh and (b) quadrilateral mesh at coarsest level.1: Convergence rate of elastic energy with respect to mesh size, presented in log-log scale .1: The field of normal stress component σxx computed by four types of element: CQ4, Q4, CT3 and T3 .2: Quarter model of the thick cylinder pipe and boundary conditions.2: Temperature distribution in the cylinder pipe .2: Hoop stress distribution in the cylinder pipe. Signed distance function.

The set of nodes enriched by Heaviside function (denoted by x symbol and the set of nodes enriched by branch functions (denoted by square symbol). Level set functions. The interaction integral domain. The weighting function q(x).1: (a) Geometry and boundary conditions of a strip with an edge crack under constant flux and (b) its deformed shape obtained with XCQ4 elements .1: Different interaction integral domain size selected to evaluate the interaction integral: (a) domain 1, (b) domain (2), (c) domain 3 and (d) domain 4.

A mesh size of 25 x 49 elements is used for both XQ4 and XCQ4 .1: Convergence of SIFs evaluated by XCQ4 and XQ4 elements with respect to number of nodes .1: Distribution of temperature [oC] .1: Stress component σxx [Pa] obtained by (a) XCQ4 element and (b) XQ4 element, showing in deformed shape.2: A square plate with horizontal center crack under two sets of boundary conditions: (a) Adiabatic crack and (b) Isothermal crack .1: Geometry and boundary conditions .1: Normalized SIFs with respect to a/W at the inclined angle θ = 30o .1: Normalized SIFs with respect to crack angle θ, given a/W = 0.2: Geometry and boundary conditions .2: Normalized SIFs with respect to a/W at the inclined angle θ = 30o .2: Normalized SIFs with respect to crack angle θ, given a/W = 0.4: Geometry and a typical mesh of a cruciform panel (AF = BE = HC = GD = 3L) .4: Crack paths predicted for four loading cases reported in Table 4.5 by the present XCQ4 element. These paths are comparable with those derived by other numerical methods available in literatures, e.4: Temperature [oC] distribution of the cruciform panel for (a) case 1 and (b) case 3 by the XCQ4. (a) Distrbution of normal stress σyy [Pa] and (b) the displacement uy [m] for case 3 by the XCQ4 .