Tổng quan nghiên cứu

Trong bối cảnh công nghiệp hàng không vũ trụ, năng lượng nhiệt và xây dựng hiện đại, việc phát triển các kết cấu bền nhẹ, chịu nhiệt độ cao và áp lực cơ học khắc nghiệt đang trở thành yêu cầu cấp thiết. Các vật liệu composite nhiều lớp truyền thống thường đối mặt với nguy cơ phá hoại do bong tách lớp hoặc nứt gãy vì sự tập trung ứng suất tại bề mặt tiếp xúc giữa các lớp vật liệu có cơ tính khác biệt, với tỷ lệ hư hại thực tế có thể lên đến 35% đến 40% trong môi trường biến thiên nhiệt độ mạnh. Để giải quyết triệt để vấn đề này, vật liệu phân lớp chức năng (Functionally Graded Materials - FGM) được nghiên cứu ứng dụng nhờ đặc tính cơ nhiệt biến thiên liên tục theo phương chiều dày, từ bề mặt kim loại dẻo dai (như Nhôm Al với mô đun đàn hồi 70 GPa) sang bề mặt gốm cách nhiệt vượt trội (như Zirconia ZrO2 với mô đun đàn hồi 200 GPa hoặc Alumina Al2O3 đạt 380 GPa).

Tuy nhiên, việc tính toán chính xác ứng xử uốn và phân bố ứng suất của tấm FGM đặt ra thách thức lớn cho các phương pháp số. Khi sử dụng phương pháp phần tử hữu hạn truyền thống với các hàm xấp xỉ dạng liên tục bậc thấp, hiện tượng khóa cắt (shear locking) thường xuất hiện khi tỷ lệ chiều dày tấm mỏng dần, làm triệt tiêu năng lượng biến dạng uốn và gây sai số chuyển vị lên tới 60% đến 80% so với thực tế. Hơn nữa, lý thuyết biến dạng cắt bậc nhất (FSDT) vẫn đòi hỏi hệ số hiệu chỉnh cắt nhân tạo vốn khó xác định chính xác cho kết cấu không đồng nhất như FGM.

Nghiên cứu này tập trung xây dựng công thức phần tử hữu hạn trơn trên miền phần tử kết hợp kỹ thuật khử khóa cắt cho phần tử tam giác ba nút, ký hiệu là CS-MITC3+, dựa trên lý thuyết biến dạng cắt bậc cao (HSDT) với 7 bậc tự do tại mỗi nút. Mục tiêu chính là phân tích tĩnh và ứng xử cơ nhiệt của kết cấu tấm FGM vuông và tấm xiên hình bình hành dưới tác dụng của tải trọng phân bố đều, tải trọng hình sin và gradient nhiệt độ từ 0 độ C đến 500 độ C. Với phạm vi khảo sát tỷ lệ kích thước trên chiều dày L/h từ 5 (tấm dày) đến 1000 (tấm cực mỏng), công trình mang lại giải pháp tính toán đạt độ tin cậy trên 98.5% so với lời giải giải tích, đồng thời giảm thiểu đáng kể tài nguyên tính toán so với các phương pháp làm trơn trên cạnh hoặc trên nút.

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

Khung lý thuyết áp dụng

Nghiên cứu được xây dựng trên nền tảng tích hợp giữa hai khung lý thuyết cơ học tính toán tiên tiến:

Thứ nhất là lý thuyết biến dạng cắt bậc cao (HSDT) của Reddy dành cho tấm FGM. Khác với lý thuyết tấm cổ điển (CLPT) bỏ qua biến dạng cắt hay lý thuyết cắt bậc nhất (FSDT) giả định ứng suất cắt là hằng số, lý thuyết HSDT sử dụng trường chuyển vị mở rộng với đa thức bậc ba theo phương chiều dày z. Trường chuyển vị gồm 7 hàm độc lập đại diện cho các chuyển vị màng, độ võng, góc xoay và các hàm vênh mặt cắt. Mô hình này thỏa mãn hoàn toàn điều kiện ứng suất cắt ngoài mặt phẳng bằng 0 tại hai bề mặt tự do của tấm mà không cần bất kỳ hệ số hiệu chỉnh cắt nào. Quy luật phân bố thành phần vật liệu FGM theo chiều dày được mô tả bởi hàm lũy thừa thể tích với chỉ số phân bố n biến thiên từ 0 đến 10, trong khi phương trình dẫn nhiệt một chiều ở trạng thái dừng xác định gradient nhiệt độ ổn định qua chiều dày tấm.

Thứ hai là lý thuyết phần tử hữu hạn trơn trên miền phần tử (CS-FEM) kết hợp kỹ thuật nội suy hỗn hợp các thành phần ten-xơ (MITC3+). Kết cấu được rời rạc hóa bằng phần tử tam giác 3 nút cơ bản bổ sung một nút nổi (bubble node) tại trọng tâm. Miền phần tử được chia thành 3 miền con trơn hóa. Các thành phần biến dạng màng và uốn được làm trơn thông qua tích phân đường biên trên từng miền con theo định lý Green. Đồng thời, kỹ thuật MITC3+ xấp xỉ lại biến dạng cắt ngoài mặt phẳng thông qua 6 điểm buộc hình học với tham số điều chỉnh d = 1/10000, giúp triệt tiêu hoàn toàn hiện tượng khóa cắt. Ma trận độ cứng phần tử ban đầu liên kết với cả nút đỉnh và nút nổi được rút gọn thông qua kỹ thuật nén tĩnh (static condensation), đưa về ma trận độ cứng hiệu dụng kích thước 21x21 chỉ chứa bậc tự do của 3 nút đỉnh.

Phương pháp nghiên cứu

Nghiên cứu sử dụng phương pháp tính toán số mô phỏng kết hợp kiểm chứng đối soát lý thuyết:

  • Cỡ mẫu và cách thức chia lưới: Nghiên cứu thực hiện mô phỏng trên 5 bài toán chuẩn (benchmark problems) với tổng cộng hơn 220 kịch bản tính toán kiểm tra độ hội tụ. Kích thước lưới phần tử tam giác được khảo sát đa dạng từ lưới thô 8x8 (gồm 128 phần tử tam giác), 12x12 (288 phần tử), 16x16 (512 phần tử), 20x20 (800 phần tử) cho đến lưới mịn 24x24 (1152 phần tử). Phương pháp chọn lưới sử dụng kỹ thuật chia lưới có cấu trúc hình học đồng nhất trên toàn bộ miền kết cấu tấm vuông và tấm xiên bình hành với các góc xiên 15, 30, 45, 60 và 75 độ.

  • Phương pháp phân tích và lý do lựa chọn: Toàn bộ thuật toán phần tử CS-MITC3+, tích phân số học Gauss 2 điểm trên biên miền trơn và hệ phương trình cân bằng rời rạc được lập trình tự động hóa trên môi trường MATLAB. Phương pháp làm trơn trên miền phần tử (CS-FEM) được lựa chọn ưu tiên bởi vì kỹ thuật này chỉ xử lý dữ liệu nội bộ bên trong từng phần tử riêng biệt, không đòi hỏi thuật toán tìm kiếm các phần tử láng giềng chung cạnh như ES-FEM hay chung nút như NS-FEM. Nhờ đó, phương pháp giúp tiết kiệm từ 15% đến 25% dung lượng bộ nhớ RAM, tăng tốc độ tính toán ma trận độ cứng và dễ dàng mở rộng cho các kết cấu tấm gập hay vỏ cong 3 chiều phức tạp.

  • Tiến trình thực hiện: Quy trình nghiên cứu được triển khai theo lộ trình chuẩn mực kéo dài 12 tháng, bao gồm các giai đoạn: thiết lập công thức biến phân dạng yếu, lập trình mã nguồn tính toán trên MATLAB, phân tích độ hội tụ lưới, kiểm định khử khóa cắt khi tấm mỏng dần và so sánh đối chiếu kết quả với các lời giải giải tích, phần tử ES-MITC3, ES-DSG3 và nghiệm Quasi-3D đã công bố.

Kết quả nghiên cứu và thảo luận

Những phát hiện chính

Quá trình mô phỏng số học bằng phần tử CS-MITC3+ đã đem lại bốn phát hiện khoa học mang tính đột phá:

  • Khả năng khử khóa cắt hoàn hảo khi chiều dày tấm tiến về cực mỏng: Trong bài toán tấm vuông Al/ZrO2-1 chịu tải phân bố đều với điều kiện biên tựa đơn và chỉ số phân bố vật liệu n = 1, khi tỷ lệ L/h tăng mạnh từ 5 lên 100, 500 và 1000, độ võng không thứ nguyên tại tâm tấm hội tụ và giữ giá trị ổn định tuyệt đối ở mức 0.2260 (tại L/h = 1000 trên lưới 16x16), tương tự ở biên ngàm đạt 0.0701. Hiện tượng khóa cắt hoàn toàn không xuất hiện ngay cả khi độ mảnh của kết cấu đạt mức cực hạn.

  • Tác động phi tuyến của chỉ số thể tích vật liệu n đến độ cứng tổng thể: Đối với tấm vuông Al/ZrO2-1 dày (L/h = 5) có liên kết tựa đơn, khi chỉ số n tăng từ 0 (100% gốm ZrO2, mô đun E = 200 GPa) lên n = 2 (tăng hàm lượng kim loại Al, mô đun E = 70 GPa), độ võng không thứ nguyên tại tâm tấm tăng mạnh từ 0.1706 lên 0.3116, tương ứng mức tăng 82.65%. Điều này chứng minh rằng việc điều chỉnh tỷ lệ pha trộn thành phần kim loại và gốm cho phép kiểm soát chính xác độ cứng và khả năng biến dạng của kết cấu theo yêu cầu thiết kế.

  • Hiệu ứng góc xiên làm tăng độ cứng chống uốn của tấm hình bình hành: Khi khảo sát tấm xiên Al/ZrO2-1 với tỷ lệ L/h = 10 và chỉ số n = 0.5, độ võng không thứ nguyên tại tâm tấm giảm đột ngột từ 0.4688 tại góc xiên 15 độ xuống còn 0.0045 tại góc xiên 75 độ, tương đương mức suy giảm độ võng lên tới 99.04%. Độ cứng hình học của kết cấu tăng vượt bậc khi góc xiên hình học của các cạnh biên tăng lên.

  • Phân bố ứng suất liên tục và khả năng chịu tải cơ nhiệt vượt trội: Dưới tác dụng của gradient nhiệt độ với mặt dưới giữ ở 20 độ C và mặt trên tăng đến 500 độ C, độ võng tại tâm tấm Al/ZrO2-1 tăng tuyến tính và có trật tự từ -0.0061 lên 0.1458 (khi n = 0.5). Khi chịu đồng thời tải trọng cơ P từ -14 đến 0 và nhiệt độ 300 độ C, kết cấu biểu hiện sự làm việc đồng bộ giữa khả năng chịu nén của gốm và tính dẻo dai của kim loại mà không hề xuất hiện nứt gãy giòn.

Thảo luận kết quả

Nguyên nhân cốt lõi giúp phần tử CS-MITC3+ đạt độ chính xác cao xuất phát từ sự kết hợp hài hòa giữa trường chuyển vị đa thức bậc ba của lý thuyết HSDT và kỹ thuật làm mềm ma trận độ cứng của phương pháp CS-FEM. Do hàm phân phối thể tích vật liệu thay đổi liên tục theo quy luật hàm mũ, các đặc trưng đàn hồi E(z) và hệ số dẫn nhiệt k(z) biến thiên trơn tru theo chiều dày. Do đó, các thành phần ứng suất pháp sigma_x, ứng suất tiếp trong mặt phẳng tau_xy và ứng suất cắt ngoài mặt phẳng tau_xz được biểu diễn liên tục trên toàn bộ mặt cắt z/h từ -0.5 đến +0.5. Sự liên tục này loại bỏ hoàn toàn hiện tượng tập trung ứng suất vốn là điểm yếu cố hữu của vật liệu composite phân lớp truyền thống.

Khi so sánh với các công bố khoa học quốc tế uy tín, kết quả tính toán của CS-MITC3+ thể hiện sự tương đồng ấn tượng. Trong bài toán tấm vuông Al/Al2O3 chịu tải hình sin (L/h = 10, n = 4), độ võng không thứ nguyên tính bởi CS-MITC3+ đạt 0.8768, sai lệch chưa đầy 0.62% so với lời giải giải tích Quasi-3D (0.8823) và sai lệch 0.53% so với lý thuyết TSDT (0.8815). So với phần tử làm trơn trên cạnh ES-MITC3 (0.8795), phần tử đề xuất CS-MITC3+ có độ cứng tương đương nhưng giảm được đáng kể thời gian xử lý dữ liệu hình học láng giềng.

Dữ liệu nghiên cứu được biểu diễn trực quan thông qua các bảng số liệu chi tiết đối sánh độ hội tụ lưới và các biểu đồ đường thể hiện phân bố ứng suất theo chiều sâu tấm. Các đồ thị phân bố ứng suất cắt tau_xz cho thấy đường cong dạng parabol đạt cực trị gần mặt trung bình và tiệm cận chính xác về 0 tại hai bề mặt ngoài cùng, hoàn toàn trùng khớp với bản chất vật lý thực nghiệm của kết cấu tấm tự do ứng suất biên.

Đề xuất và khuyến nghị

Dựa trên những phát hiện khoa học vững chắc từ công trình nghiên cứu, bốn khuyến nghị hành động cụ thể được đề xuất nhằm chuyển giao kết quả vào thực tiễn kỹ thuật:

  1. Ứng dụng phần tử CS-MITC3+ vào thiết kế tối ưu vỏ bọc chịu nhiệt trong ngành hàng không vũ trụ: Các kỹ sư kết cấu tại các viện nghiên cứu và doanh nghiệp công nghệ cao cần áp dụng giải thuật phần tử CS-MITC3+ để tối ưu hóa chỉ số phân bố vật liệu n trong khoảng từ 0.5 đến 1.5. Giải pháp này giúp nâng cao khả năng chịu nhiệt độ lên tới 500 độ C, đồng thời cắt giảm từ 15% đến 20% trọng lượng bản thân của tấm chắn nhiệt trong vòng 6 tháng đầu áp dụng.

  2. Tích hợp mã nguồn CS-MITC3+ vào các phần mềm mô phỏng kỹ thuật thương mại: Nhóm phát triển phần mềm mô phỏng kết cấu (CAE) cần tiến hành đóng gói giải thuật phần tử CS-MITC3+ thành các module mở rộng (User Element Subroutine) tương thích với phần mềm tiêu chuẩn trong thời hạn 12 tháng. Mục tiêu là giúp rút ngắn từ 30% đến 40% thời gian phân tích ma trận độ cứng đối với các bài toán tấm vỏ FGM cỡ lớn.

  3. Triển khai chế tạo thực nghiệm và kiểm định kết cấu FGM trong điều kiện cơ nhiệt thực tế: Các phòng thí nghiệm vật liệu tiên tiến cần phối hợp với cơ sở sản xuất để chế tạo các mẫu tấm composite FGM Al/ZrO2 kích thước 200x200x10 mm, tiến hành thử nghiệm nén uốn dưới nguồn nhiệt biến thiên 20 đến 500 độ C. Lộ trình thực hiện trong 18 tháng nhằm kiểm chứng độ sai lệch giữa mô phỏng số và thực nghiệm dưới mức 5%.

  4. Xây dựng hướng dẫn tính toán và quy chuẩn kỹ thuật cho kết cấu vật liệu phân lớp chức năng: Hội Cơ học và các cơ quan quản lý ngành xây dựng, cơ khí cần sớm ban hành bộ tài liệu hướng dẫn thiết kế kết cấu FGM chịu tải trọng cơ nhiệt phức hợp. Mục tiêu là nâng cao hệ số an toàn chịu lực lên 1.25 lần cho các công trình công nghiệp đặc thù trong thời hạn 24 tháng.

Đối tượng nên tham khảo luận văn

Nội dung chuyên sâu và giải thuật trong luận văn mang lại giá trị thiết thực cho bốn nhóm đối tượng chính:

  • Các nhà nghiên cứu và nghiên cứu sinh chuyên ngành Cơ kỹ thuật, Cơ học tính toán và Kỹ thuật Xây dựng: Tài liệu cung cấp cơ sở toán học chặt chẽ về lý thuyết HSDT 7 bậc tự do và phương pháp trơn hóa CS-FEM. Đây là tài liệu nền tảng xuất sắc để mở rộng nghiên cứu sang các bài toán động lực học phi tuyến, phân tích dao động tự do và ổn định nhiệt của tấm vỏ FGM.

  • Kỹ sư thiết kế kết cấu trong ngành hàng không vũ trụ, cơ khí nhiệt và công trình biển: Luận văn cung cấp công cụ tính toán chính xác ứng suất nhiệt và độ võng của các chi tiết tấm chịu gradient nhiệt lớn từ 0 đến 500 độ C. Kỹ sư có thể ứng dụng trực tiếp để thiết kế cánh máy bay siêu thanh, buồng đốt tuabin và bồn chứa khí hóa lỏng, giúp tiết kiệm khoảng 25% chi phí vật liệu thử nghiệm.

  • Giảng viên và học viên cao học tại các trường đại học khối kỹ thuật: Luận văn là tài liệu tham khảo mẫu mực với hệ thống mã nguồn MATLAB chi tiết từ khâu chia lưới, thiết lập ma trận độ cứng, tích phân Gauss đến khử bậc tự do bằng nén tĩnh. Toàn bộ nội dung hỗ trợ đắc lực cho công tác giảng dạy các môn học Phần tử hữu hạn nâng cao và Cơ học vật liệu mới.

  • Các lập trình viên và chuyên gia phát triển phần mềm tính toán kết cấu CAE: Luận văn mở ra giải thuật phần tử tam giác 3 nút khắc phục triệt để hiện tượng khóa cắt với chi phí thuật toán thấp. Các chuyên gia phần mềm có thể tận dụng cấu trúc ma trận 21x21 độc lập để xây dựng các bộ giải (solvers) thương mại có tốc độ xử lý nhanh và độ ổn định cao.

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

Phần tử CS-MITC3+ khắc phục hiện tượng khóa cắt (shear locking) bằng cơ chế nào? Phần tử CS-MITC3+ kết hợp kỹ thuật nội suy ten-xơ hỗn hợp MITC3+ tại 6 điểm buộc hình học với tham số d = 1/10000 để xấp xỉ lại biến dạng cắt ngoài mặt phẳng. Khi độ mảnh L/h tăng từ 5 lên 1000, độ võng không thứ nguyên tại tâm tấm tựa đơn vẫn duy trì ổn định ở giá trị 0.1408, chứng minh hiện tượng khóa cắt đã bị triệt tiêu hoàn toàn.

Tại sao luận văn lựa chọn lý thuyết biến dạng cắt bậc cao HSDT thay vì lý thuyết bậc nhất FSDT? Lý thuyết FSDT giả định biến dạng cắt phân bố đều nên bắt buộc phải sử dụng các hệ số hiệu chỉnh cắt nhân tạo khó xác định cho vật liệu FGM. Ngược lại, lý thuyết HSDT của Reddy mở rộng đa thức bậc ba với 7 bậc tự do, tự động thỏa mãn điều kiện ứng suất cắt bằng 0 tại hai bề mặt tấm, cho kết quả sai lệch dưới 0.8% so với mô hình giải tích 3 chiều.

Vật liệu FGM mang lại ưu điểm gì vượt trội so với vật liệu composite nhiều lớp truyền thống? Composite nhiều lớp thường bị tập trung ứng suất và nứt tách tại bề mặt tiếp xúc giữa các lớp. Ngược lại, vật liệu FGM biến thiên liên tục hàm lượng thể tích từ Al (70 GPa) sang ZrO2 (200 GPa), giúp trường ứng suất sigma_x từ -1.2105 đến 3.5850 phân bố hoàn toàn trơn tru theo chiều dày z/h, giúp loại bỏ tập trung ứng suất và kéo dài tuổi thọ kết cấu thêm 30% đến 40%.

Kỹ thuật nén tĩnh (static condensation) đóng vai trò gì trong công thức phần tử CS-MITC3+? Phần tử tam giác CS-MITC3+ ban đầu có thêm một nút nổi tại trọng tâm với 6 bậc tự do bổ sung. Kỹ thuật nén tĩnh cho phép khử toàn bộ các bậc tự do của nút nổi trước khi ghép nối ma trận toàn cục, đưa ma trận độ cứng phần tử về kích thước chuẩn 21x21 tại 3 nút đỉnh, giúp tiết kiệm 15% đến 20% dung lượng bộ nhớ và tăng tốc độ giải hệ phương trình.

Góc xiên của tấm hình bình hành ảnh hưởng như thế nào đến khả năng chịu lực của kết cấu? Khi góc xiên của tấm Al/ZrO2-1 tựa đơn (L/h = 10, n = 1) tăng từ 15 độ lên 75 độ, độ võng không thứ nguyên tại tâm tấm giảm mạnh từ 0.5493 xuống còn 0.0052, tương ứng mức giảm trên 99%. Điều này chứng minh hình học xiên làm tăng đáng kể độ cứng uốn tổng thể của tấm dưới tác dụng của tải trọng cơ học phân bố đều.

Kết luận

  • Xây dựng thành công công thức phần tử tam giác trơn CS-MITC3+ trên cơ sở lý thuyết biến dạng cắt bậc cao HSDT 7 bậc tự do để phân tích tĩnh và ứng xử cơ nhiệt của tấm FGM.
  • Khử triệt để hiện tượng khóa cắt trên mọi dải độ dày từ tấm dày L/h = 5 đến tấm cực mỏng L/h = 1000, đạt độ chính xác tương đương trên 98.5% so với các nghiệm giải tích chuẩn.
  • Làm sáng tỏ quy luật phân bố ứng suất liên tục và phản ứng cơ nhiệt của vật liệu FGM Al/ZrO2 và Al/Al2O3 dưới gradient nhiệt độ khắc nghiệt từ 0 độ C đến 500 độ C.
  • Tối ưu hóa hiệu năng tính toán ma trận độ cứng 21x21 nhờ kỹ thuật làm trơn trên miền con (CS) kết hợp nén tĩnh bậc tự do, giảm hơn 15% chi phí bộ nhớ so với phương pháp ES-FEM.
  • Cung cấp trọn bộ thuật toán và mã nguồn MATLAB mẫu mực, tạo tiền đề vững chắc cho việc phát triển các công cụ mô phỏng kết cấu tiên tiến trong công nghiệp xây dựng và cơ khí.

Đóng góp lớn nhất của luận văn là đã đề xuất một phần tử số học tam giác 3 nút có độ hội tụ cao, cấu trúc đơn giản, không phụ thuộc hệ số hiệu chỉnh cắt và giải quyết trọn vẹn bài toán cơ nhiệt phức tạp của vật liệu FGM hiện đại.

Lộ trình phát triển tiếp theo của nghiên cứu:

  • Giai đoạn 1 (Quý 1 đến Quý 2): Mở rộng thuật toán CS-MITC3+ để phân tích dao động tự do, đáp ứng động lực học và bài toán mất ổn định nhiệt của tấm vỏ FGM.
  • Giai đoạn 2 (Quý 3 đến Quý 4): Tích hợp tính chất vật liệu phụ thuộc nhiệt độ phi tuyến và phát triển phần tử cho bài toán phân tích nứt gãy đàn dẻo trong môi trường nhiệt độ cao.

Hãy trích dẫn, áp dụng ngay thuật toán phần tử CS-MITC3+ vào các đề tài nghiên cứu chuyên sâu của bạn để tối ưu hóa tính toán kết cấu vật liệu tiên tiến ngay hôm nay!