Phân tích bài toán nứt phẳng đàn hồi tuyến tính bằng phần tử hữu hạn mở rộng tứ giác nội suy kép (xcq4)

 15 Tạp chí Khoa học Lạc Hồng Số 04 
Journal of Science of Lac Hong University
Vol.4 (12/2015), pp. 15-20
Tạp chí Khoa học Lạc Hồng
Số 4 (12/2015), trang 15-20
PHÂN TÍCH BÀI TOÁN NỨT PHẲNG ĐÀN HỒI TUYẾN TÍNH BẰNG PHẦN 
TỬ HỮU HẠN MỞ RỘNG TỨ GIÁC NỘI SUY KÉP (XCQ4) 
Optimization of linear elastic fracture mechanics by using an extended 
consecutive – Interpolation quadrilateral element (XCQ4) method
Nguyễn Đình Dư1, Nguyễn Duy Phích2
1dinhdu85@gmail.com
Khoa Kỹ Thuật Công trình Trường Đại học Lạc Hồng, Đồng Nai, Việt Nam
Đến tòa soạn 10/1/2015; Chấp nhận đăng: 1/2/2015
Tóm tắt. Bài báo trình bày một phần tử hữu hạn mở rộng tứ giác (XCQ4) dựa trên thủ tục nội suy kép với ứng suất liên tục tại 
nút để mô phỏng trường ứng suất vùng lân cận đỉnh vết nứt hai chiều. Khác với phương pháp truyền thống, hàm xấp xỉ trong 
nghiên cứu này bao gồm giá trị tại nút và trung bình cộng giá trị đạo hàm của nó từ bước đầu tiên. Mục tiêu chính của bài 
viết này là nhằm giới thiệu một sự phát triển của phần tử CQ4 được công bố gần đây với kỹ thuật làm giàu nhằm tính chính 
xác hệ số cường độ ứng suất tại đỉnh nứt (SIFs). Sự chính xác của phần tử XCQ4 trong nghiên cứu này được chứng minh 
thông qua các ví dụ số có hình học từ đơn giản đến phức tạp. Độ chính xác cũng như tốc độ hội tụ thu được SIFS từ XCQ4 là 
cao hơn XQ4 truyền thống với điều kiện như nhau. Phần tử mới XCQ4 có thể được mở rộng và áp dụng cho các bài toán 
phức tạp hơn trong thực tế. 
Từ khoá: Phần tử hữu hạn mở rộng; Cơ học rạn nứt; Hệ số cường độ ứng suất; Kỹ thuật làm giàu; Nội suy kép 
Abstract. This work presents a novel extended 4-node quadrilateral finite element (XCQ4) method based on the consecutive-
interpolation procedure (CIP) with continuous nodal stress in oder toaccurately model singular stress fields near crack tips 
of two-dimensional (2D) elastic cracked solids. Incontrast with conventional FEM methods, the approximation functions 
constructed based on CIP involve both nodal values and averaged nodal gradients as interpolation conditions. As a 
pioneering extension of a recently developed CQ4 element associated with enrichment method, the proposed XCQ4 
extracts the stress intensity factors (SIFs) at the crack tipsprecisely. Accuracy and convergence of the SIFs results obtained 
by the proposed method are as high as those of the standard XQ4 solutions. Importantly, the proposed XCQ4 
elementmethod is highly promising for use in other complex engineering problems. 
Keywords: Fracture; Extended finite element method; Stress intensity factors; Smooth nodal stress; Enrichment; Consecutive-
interpolation 
1. GIỚI THIỆU 
Mô hình hóa chính xác trường chuyển vị và biến dạng 
đỉnh nứt vẫn đang là một vấn đề thách thức trong cơ học 
rạn nứt. Việc dự đoán chính xác trường ứng suất cạnh đỉnh 
nứt đóng một vai trò quan trọng trong việc bảo trì, dự đoán 
tuổi thọ và đánh giá sự an toàn của vật liệu cũng như kết 
cấu trong tương lai.
Trong một vài thập kỷ qua, phương pháp phần tử hữu 
hạn (FEM) được xem là một công cụ số hiệu quả giải quyết 
nhiều bài toán kỹ thuật phức tạp mà lời giải giải tích không 
đáp ứng được. Tuy nhiên, khi giải quyết bài toán nứt gặp 
nhiều khó khăn do sự liên kết cấu trúc phần tử phải được 
cập nhập một cách liên tục trong suốt quá trình tái chia 
lưới. Do đó, có nhiều nỗ lực nhằm cải tiến phương pháp 
FEM để phù hợp với yêu cầu kỹ thuật, một trong những 
phương pháp đó chính là phần tử hữu hạn mở rộng 
(XFEM).
Trên cơ sở đó, trong phạm vi nghiên cứu này, một nghiên 
cứu mới cho bài toán nứt của kết cấu đàn hồi hai chiều 
được sử dụng bởi phần tử hữu hạn mở rộng tứ giác nội suy 
kép (XCQ4) được trình bày. Kết quả số thu được từ phương 
pháp được so sánh với XFEM, các phương pháp số khác 
cũng như lời giải chính xác nhằm kiểm chứng hiệu suất của 
XCQ4.
2. XÂY DỰNG PHẦN TỬ XCQ4 PHÂN TÍCH NỨT 
Về cơ bản mà nói, XCQ4 được nghiên cứu trong bài báo 
này là một phiên bản cải tiến của chúng tôi về CQ4 được 
tích hợp các chức năng làm giàu để mô phỏng chính xác 
trường biến dạng của vết nứt.
2.1 Phần tử CQ4 và các thuộc tính 
Mô tả chi tiết phần tử CQ4 có thể tìm thấy trong [1]. Tuy 
nhiên, để tiện theo dõi, tác giả xin trình bày ngắn gọn CQ4 
được trình bày trong bài báo này. Một điểm cần nội suy có 
tọa độ x(x, y) trong phần tử tứ giác có bốn nút lần lượt i, j, 
k, m được minh họa trong Hình 1. Như mô tả trong hình vẽ, 
chúng ta có các miền phần tử Si, Sj, Sk và Sm là các miền 
chứa tất cả các phần tử có liên quan lần lượt đến các nút i, j, 
k và m. Như vậy, những nút hỗ trợ cho điểm x trong phần 
tử CQ4 bao gồm tất cả các nút của miền phần tử Si, Sj, Sk và 
Sm.
Nguyên Đình Dư, Nguyên Duy Phích
 16 Tạp chí Khoa học Lạc Hồng Số 04 
Hình 1. Hình minh họa CQ4 trong 2D
Theo đó, miền hỗ trợ CQ4 cho điểm x thì rộng hơn miền 
hỗ trợ Q4 trong FEM chuẩn. Lời giải gần đúng tại điểm x
có thể viết như sau:
(1)
Trong phương trình (1), hàm dạng nội suy kép được xác 
định như sau [1], [2]:
(2)
Trong đó, dl là chuyển vị tại nút, còn là hàm dạng 
của nút i và ns là tổng số nút hỗ trợ có liên quan đến điểm x.
Trong công thức nội suy, việc xây dựng đạo hàm trung bình 
tại nút i (những nút khác tương tự) được thực hiện như [1] 
[2], có thể viết như sau:
(3)
Trong phương trình (3), thành phần là đạo hàm 
của được tính trong phần tử e, và we là hàm trọng số
của phần tử e ÎSi, nó được định nghĩa như sau:
, với e ÎSi (4)
và De là diện tích của phần tử e.
Trong phương trình (2), các hàm số i, ix và iylà các đa 
thức cơ sở có liên quan đến nút i phải thỏa mãn điều kiện 
sau:
i (xl) = il, i,x (xl) = 0 , i,y (xl) = 0
ix (xl) = 0 , ix,x (xl) = il, ix,y (xl) = 0 (5) (13)
iy (xl) = 0 , iy,x (xl) = 0 , iy,y (xl) = il
Trong đó, l lần lượt là một trong các chỉ số i, j, k và m, và
(6)
Cần lưu ý rằng, các điều kiện trên cũng cần áp dụng 
tương tự cho cho các hàm số khác như j , jx , jy; k , kx , 
kyvà m , mx , my. Sau cùng, các hàm đa thức cơ sở trong 
CQ4 được định nghĩa như sau:
(7)
Trong phương trình (7), p = 1/2 và các hàm số j , jx , 
jy; k , kx , ky ; m , mx , my cũng được thực hiện tương tự
bằng cách xoay vòng các chỉ số i, j, k và m. Theo đó, Li, Lj, 
Lk và Lm là hệ trục tọa độ diện tích của điểm cần nội suy x 
trong phần tử tứ giác i, j, k và m.
Tất cả các tính chất của phần tử CQ4 có thể được tìm 
thấy trong [1] và một điểm quan trọng cần nhắc đến ở đây 
đó là hàm dạng của CQ4 thì trơn và mịn hơn Q4 chuẩn 
được thể hiện trong Hình 2.
Hình 2. Hàm dạng của Q4 (a), CQ4 (b) trong 2D.
2.2Mở rộng xấp xỉ CQ4 cho nứt 
Cũng giống như ý tưởng chủ đạo của XFEM [3] là sử
dụng một chuyển vị xấp xỉ có khả năng mô tả chính xác sự
bất liên tục và kỳ dị của vùng lân cận đỉnh nứt. 
Do đó, xấp xỉ mở rộng CQ4 của chuyển vị cho vết nứt có 
thể viết như sau:
(a) 
(b) 
Phân tích bài toán nứt phẳng đàn hồi tuyến tính bằng phần tử hữu hạn mở rộng tứ giác nội suy kép (XCQ4)
 17 Tạp chí Khoa học Lạc Hồng Số 04 
(8)
Trong đó, đại diện cho hàm dạng CQ4 liên quan đến 
nút i để nội suy trường chuyển vị tiêu chuẩn, Jcut là tập hợp 
các nút mở rộng thuộc phần tử bị vết nứt cắt qua, Ktip là tập 
hợp các nút mở rộng thuộc phần tử chứa đầu vết nứt, ui là 
vectơ chuyển vị của nút i. H(x) là hàm Heaviside có giá trị
+1 trên đường nứt và – 1 dưới đường nứt, aj là chuyển vị
bậc tự do mở rộng thuộc Jcut. Cuối cùng, là bậc tự do 
chuyển vị tại nút K thuộc Ktip phù hợp với hàm mở rộng 
tiệm cận .
Đối với bài toán đàn hồi tuyến tính đẳng hướng 2D, hàm 
được định nghĩa như sau:
(9)
Trong đó, (r, ) là hệ tọa độ cực địa phương có gốc tọa độ là 
đỉnh nứt.
2.3 Phương trình rời rạc 
Các hệ thống rời rạc của phương trình cân bằng tuyến 
tính trong điều kiện đàn hồi biến dạng nhỏ có thể được thể
hiện như sau:
Kd = F (10)
Trong đó, K là ma trận độ cứng tổng thể, d là vectơ chuyển 
vị nút tổng thể bao gồm nút tiêu chuẩn và nút mở rộng. Đối 
với phần tử được làm giàu, ma trận độ cứng phần tử thu 
được:
(11)
Trong khi đối với phần tử tiêu chuẩn thì 
(12)
với 
; r, s@u, a, b (13)
Trong đó, B ma trận tính biến dạng, là đạo hàm đối xứng 
rời rạc của ma trận hàm dạng mở rộng. Các thành phần của 
B được cho như sau:
; ; 
(14)
Trong phương trình (10), F là vectơ lực nút tác động từ
bên ngoài và bản thân của phần tử dẫn đến vectơ lực tổng 
thể phần tử được viết như sau:
(15)
Đối với phần tử tiêu chuẩn thì 
(16) 
Các thành khác là của phần tử mở rộng, chúng được định 
nghĩa chi tiết như sau:
(17)
(18)
(19)
với và lần lượt là lực bản thân trên một đơn vị thể tích 
và lực tác dụng lên biên.
3. TÍCH PHÂN TƯƠNG TÁC VÀ CÁCH TÍNH SIFS 
Để thu được SIFs, dạng miền của tích phân tương tác 
được chọn trong bài viết này [4]. Hai trạng thái của tấm 
chứa vết nứt được xem xét. Trạng thái #1 là 
trạng thái thực tương ứng với sự làm việc hiện tại, còn trạng 
thái #2 là trạng thái ảo, trạng thái xấp xỉ, là 
trạng thái được chọn như là trường tiệm cận khu vực đầu 
đỉnh nứt cho mode I và mode II. Mối quan hệ giữa tích 
phân tương tác và SIFs đa mode như sau:
(20)
Trong đó, M(1,2) là tích phân tương tác và được định nghĩa
như sau:
(21)
Trong (20), là tham số vật liệu và được định nghĩa như 
sau:
(22)
và trong (21), .
Nếu trạng thái #2 được giả định là trường tiệm cận của 
mode I với KI(2) = 1 và KII(2) = 0, sau đó SIF mode I cho 
trạng thái #1 thu được như sau:
(23)
Ứng suất phẳng
Biến dạng phẳng
Nguyên Đình Dư, Nguyên Duy Phích
 18 Tạp chí Khoa học Lạc Hồng Số 04 
Tương tự cách tính trên, SIF mode II cho trạng thái #1 
thu được như sau:
(24)
Tích phân M(1,2) trong (21) không phù hợp lắm với phần 
tử hữu hạn, gây rất nhiều khó khăn trong tính toán. Vì thế, 
nó được chuyển sang tích phân miền bằng cách nhân thêm 
một hàm trọng số biên q(x), nó được minh họa như Hình 3, 
giá trị của nó thay đổi từ 1 khi ở gần đỉnh nứt và bằng 0 
trên đường viền bên ngoài theo quy định. Cuối cùng, tích 
phân tương tác trong hình thể mới được xác định bởi:
(25)
Một vấn đề cần lưu ý thêm khi thực hiện nội suy kép 
trong hàm xấp xỉ, đó là miền hỗ trợ cho các phần tử thuộc 
J-domain trong XCQ4 là lớn hơn XQ4 chuẩn. J-domain 
được thể hiện trong Hình 4a cho XCQ4 và 4b cho XQ4. Do 
đó hệ số SIFS thu được từ XCQ4 có độ chính xác cao hơn 
XQ4 và được kiểm chứng trong các ví dụ số.
Hình 3. Hình ảnh minh họa hàm trọng số q(x) dùng trong tích 
phân tương tác: (a) mặt bằng, (b) mặt đứng
(a) (b)
Hình 4. Hình ảnh minh họa J-domain giữa Q4 (a) và XCQ4 (b)
4. VÍ DỤ SỐ 
Hai ví dụ về bài toán cơ học rạn nứt đàn hồi tuyến tính 
hai chiều bao gồm mode đơn và hỗn hợp được xem xét. Hệ
số cường độ ứng suất (SIFs) được tính bằng tích phân 
tương tác thu được từ XCQ4 sau đó so sánh với XQ4 tiêu 
chuẩn hoặc các phương pháp số hiện có khác. Các phân tích 
về tổng số bậc tự do, cách chia lưới quy tắc hay bất quy tắc, 
bán kính tích phân tương tác nhằm đánh giá mức độ hội tụ
cũng như độ chính xác của phương pháp. Trong các ví dụ
số, điều kiện biến dạng phẳng được áp dụng.
4.1 Bài toán 1 
Trong ví dụ này, một tấm hữu hạn chứa vết nứt cạnh chịu
lực cắt phân bố phân bố đều tác dụng vào cạnh biên trên 
như Hình 5.a. Tấm có chiều rộng W = 7, chiều cao L = 16 
và chiều dài vết nứt a = 3.5. Lực cắt có cường độ t =1 được 
áp vào biên trên của tấm. Dưới tác dụng của lực này, tấm sẽ
bị phá hoại theo hai trạng thái đó là mode 1 và mode 2. Lời 
giải chính xác [1] cho kết quả KI = 34.0, KII = 4.55, là điều 
kiện dùng để so sánh với XCQ4.
(a) (b)
Hình 5. Tấm hữu hạn chứa vết nứt cạnh chịu lực cắt (a) và lực 
kéo (b)
Bảng 1. Độ hội tụ của SIFs với mật độ lưới của tấm hữu hạn 
chịu lực cắt
R Phương 
Pháp
Mật độ lưới
9´19 19´39 29´59 39´79
2 K1 XQ4 31.3742 33.0055 33.4473 33.6480
XCQ4 31.7094 33.1934 33.6175 33.8123
K2 XQ4 4.4075 4.4632 4.4775 4.4839
XCQ4 4.4987 4.5301 4.5635 4.5553
2.5 K1 XQ4 31.1810 32.8235 33.2693 33.4731
XCQ4 31.5365 33.0080 33.4352 33.6361
K2 XQ4 4.3953 4.4522 4.4698 4.4774
XCQ4 4.4467 4.4764 4.4838 4.4911
3 K1 XQ4 31.2015 32.8516 33.3274 33.5318
XCQ4 31.3362 33.0001 33.4795 33.6808
K2 XQ4 4.3912 4.4477 4.4664 4.4737
XCQ4 4.4265 4.4563 4.4709 4.4989
Sự ảnh hưởng của tích phân J phụ thuộc vào bán kính R 
được mô tả như Hình 6. 
Trong nghiên cứu này, bán kính J-domain R được khảo 
sát ở ba giá trị khác nhua bao gồm R = 2, R = 2.5, R = 3. 
Mật độ chia lưới hay tổng số bậc tự do là một trong những 
yếu tố ảnh hưởng rất lớn đến kết quả tính toán. Theo đó, 
bốn cách chia lưới có quy tắc 9´19; 19´39; 29´59 và 
39´79 được xem xét. Kết quả tính toán của mỗi miến tích 
phân cùng với cách chia lưới được thể hiện trong Bảng 1. 
Kết quả số thu được là rất tốt, sự chính xác của hệ số SIFs 
càng tăng gần với lời giải chính xác khi số phần tử tăng lên. 
Đồng thời, giá trị trong bảng 1 cũng cho thấy rằng kết quả
SIFs thu được từ XCQ4 là tốt hơn XQ4 truyền thống. Tầm 
ảnh hưởng của bán kính J-domain khi tính hệ số SIFs ở chế
độ đa mode là nhẹ nhàng và không đáng kể. Tuy nhiên, một 
bán kính miền không quá lớn cũng như không quá nhỏ
được chọn nhằm đảm bảo tính chính xác của phương pháp. 
Do đó, R = 2; 2.5 hoặc 3 có thể được dùng.
Sai số tương đối của hệ số SIFs khi tính bằng XCQ4 và 
XQ4 ở mode I và mode II so với lời giải chính xác được 
phác họa ở Hình 6 cho các bán kính miền J khác nhau. Dễ
dàng nhận thấy tốc độ hội tụ của cả hai mode là tương tự
nhau. 
Điều đó có nghĩa rằng tốc độ hội tụ của XCQ4 cho các 
trường hợp R khác nhau thì tốt hơn XQ4 cho cả hai mode. 
Do đó, khi thực hiện với XCQ4 thì ít suy biến hơn XQ4 
truyền thống.
(a)
(b)
Phân tích bài toán nứt phẳng đàn hồi tuyến tính bằng phần tử hữu hạn mở rộng tứ giác nội suy kép (XCQ4)
 19 Tạp chí Khoa học Lạc Hồng Số 04 
Hình 6. Kết quả hội tụ của mode I và mode II 
thu được từ XCQ4 và XQ4 với R khác nhau
4.2 Bài toán 2 
Tấm hữu hạn chịu lực kéo như sơ đồ tính thể hiện trong 
Hình 5.b. 
Lực kéo có cường độ =1 đặt tại biên trên của tấm, với 
điều kiện của bài toán như vậy thì chỉ có mode I bị phá 
hoại. Tương tự như lực cắt ở trên, ta cũng đi xét sự ảnh 
hưởng của bán kính miền tích phân J nhưng chỉ hai giá trị R
= 2 và R = 2.5 được xem xét. Độ hội tụ của hệ số SIF vào 
mật độ lưới cũng được phân tích. Kết quả số hệ số SIF 
mode I được kiểm chứng với lời giải chính xác được cho 
bởi Ewalds và Wanhill [5].
(26)
Hình 7. Kết quả hội tụ của mode I thu được từ XCQ4 và XQ4 với 
hai giá trị R
Trong đó, hệ số điều chỉnh C được tính như sau:
(27)
Kết quả tính toán và giá trị lời giải chính xác được thể
hiện trong Bảng 2. Như mong đợi, dễ dàng quan sát được 
sự chính sát của kết quả số thu được khi tính hệ số SIF 
mode I bởi XCQ4 là tốt hơn XQ4 thông thường. Đường 
biểu diễn kết quả số SIF của nghiên cứu này tiệm cận với 
lời giải chính xác khi tăng số lượng phần tử được thể hiện ở
Hình 7. Cũng giống như ví dụ ở trước, sự ảnh hưởng của 
bán kính miền tích phân J cũng tìm thấy trong ví dụ này.
Bảng 2. Độ hội tụ của KI phụ thuộc vào mật độ lưới
R Phương
Pháp
Mật độ lưới Chính
xác
9´19 19´39 29´59 39´79
2 XQ4 8.7292 9.1050 9.2137 9.2649 9.3721
XCQ4 8.7808 9.1499 9.2591 9.3105
2.5 XQ4 8.6869 9.0602 9.1681 9.2192
XCQ4 8.7362 9.1012 9.2108 9.2632
Hơn nữa, sai số tương đối của hệ số SIF mode I được rời 
rạc hóa theo bán kính miền R của tích phân J được minh 
họa trong Hình 8 rõ ràng cho thấy thêm lần nữa là XCQ4 
tốt hơn XQ4.
Nguyên Đình Dư, Nguyên Duy Phích
 20 Tạp chí Khoa học Lạc Hồng Số 04 
Hình 8. Trường ứng suất von Mises
giữa XQ4 (a) và XCQ4 (b)
Để tiếp tục giải thích biểu hiện sự hiệu quả của phương 
pháp XCQ4 trong vấn đề phân bố ứng suất. Ứng suất von 
Mises trong ví dụ này với kích thước chia lưới có quy tắc 
39´79 được thể hiện trong Hình 8.
Kết quả cho thấy rằng việc sử dụng phần tử XCQ4 với 
thủ tục nội suy kép được đề xuất trong nghiên cứu này cho 
trường ứng suất mịn màng hơn XQ4 truyền thống.
5. KẾT LUẬN 
Một phương pháp số mở rộng bằng cách sử dụng các hàm 
số cơ bản của CQ4 để mô hình hóa chính xác các vết nứt đã 
được phát triển bởi kỹ thuật làm giàu. Nó bao gồm hàm bất 
liên tục Heaviside và một trường tiệm cận gần đỉnh nứt.
Nghiên cứu đã sử dụng các mở rộng của phần tử CQ4 
(XCQ4), nó đặc biệt phù hợp cho mô hình vết nứt, để phân 
tích hệ số cường độ ứng suất và một số vấn đề trong bài 
toán cơ rạn nứt đàn hồi tuyến tính 2D. Trong mỗi trường 
hợp nghiên cứu, các SIFs và độ hội tụ của nó được tìm thấy 
khớp với các phân tích tương ứng như XQ4 nói chung và 
các phương pháp số khác nói riêng. 
Nói chung, các kết quả số thu được từ XCQ4 cho thấy
đạt một hiệu suất cao vì nó mang lại sự chính xác cao hơn 
phần tử thông thường. Mặt khác, XCQ4 còn thể hiện được 
trường ứng suất mịn màng ở đỉnh nứt mà không thể dễ dàng 
thu được tương tự ở phần tử thông thường.
6. TÀI LIỆU THAM KHẢO 
[1] Q.T. Bui, Q.D.Vo, C.H. Zhang, D.D. Nguyen, “A 
consecutive-interpolation quadrilateral element (CQ4): 
Formulation and applications,” Finite Elem Anal Des, vol. 
84, pp. 14-31, February 2014.
[2] Y.T. Yang, X.H.Tang, H. Zheng, “A three-node triangular 
element with continuous nodal stress,” Comput Struct, vol. 
141, pp. 46-58, 2014.
[3] N. Moes, J. Dolbow, T. Belytschko, “A finite element 
method for crack growth without remeshing,” Int. J.
Numer. Meth. Engng, vol. 46, pp. 131-150, 1999.
[4] Q.T. Bui,C.H. Zhang, “Extended finite element simulation 
of stationary dynamic cracks in piezoelectric solids under 
impact loading,” Comput Mater Sci, vol. 62, pp. 243-257, 
2012.
[5] H. Ewalds andR. Wanhill, Fracture Mechanics, Edward 
Arnold: New York, 1993.
TIỂU SỬ TÁC GIẢ 
Nguyễn Đình Dư 
Năm sinh 1985, Tuy Phước, Bình Định. Tốt nghiệp Đại học và Thạc sĩ tại Trường Đại học Bách Khoa Thành 
phố Hồ Chí Minh năm 2009 và 2013. Hiện anh đang là Giảng viên cơ hữu tại Khoa Kỹ thuật Công trình. Lĩnh 
vực nghiên cứu: Phương pháp số, Phần tử hữu hạn, Cơ học vật rắn biến dạng.