Phân tích kết cấu tấm bằng phần tử mitc3+ được làm trơn trên cạnh (es-mitc3+)

PHÂN TÍCH KẾT CẤU TẤM BẰNG PHẦN TỬ MITC3+  
ĐƢỢC LÀM TRƠN TRÊN CẠNH (ES-MITC3+)  
SUCCESSFUL DESIGNS OF MITC3 + SECTION  
MADE IN THE EYES (ES-MITC3 +)  
Trần Văn Chơn1, Châu Đình Thành2  
1 Trường đại học Sư phạm Kỹ thuật TP.HCM  
2 Trường đại học Sư phạm Kỹ thuật TP.HCM  
T M T T  
Bài báo này nhằm mục đích nghiên cứu kỹ thuật làm trơn trên cạnh để phát triển cho  
phần tử MITC3+, để tạo nên một phần tử mới được gọi là phần tử ES-MITC3+ có khả năng phân  
tích tĩnh cho tấm đồng nhất và đẳng hướng. Do sử dụng kỹ thuật MITC3+ để khử hiện tượng  
khóa cắt nên phần tử này có thể được dùng để giải quyết cho các bài toán tấm dày và tấm mỏng.  
Phần tử ES-MITC3+ được sử dụng để tính toán cho một số bài toán tấm điển hình nhằm  
đánh giá mức độ chính xác và hội tụ. Thông qua các kết quả ví dụ số minh họa đạt được, ta thấy  
phần tử ES-MITC3+ có khả năng phân tích tĩnh cho cả tấm mỏng và tấm dày với độ chính xác  
tương đương các phần tử tấm tam giác 3 nút cùng loại.  
T   h    
ABSTRACT:  
In this thesis, the edge-based smoothed strain method (ES) has been developed for MITC3+ plate  
elements. The proposed element called ES-MITC3+ is able to statically analyze homogeneous  
and isotropic plates. The shear-locking phenomenon is removed by using the MITC3+ technique.  
Therefore, the ES-MITC3+ plate element can be used to solve both thin and thick plates.  
The accuracy and convergence of the ES-MITC3+ are evaluated by analyzing the  
defection and moment of some benchmark plate problems. The numerical results show that the  
ES-MITC3+ is capable of computing behaviors of both thin and thick plates with accuracy  
similar to other three-node triangular plate elements.  
Keywords:  
được sử dụng làm khung chịu lực cho các  
công trình khác.  
1. ĐẶT VẤN ĐỀ  
Sự phát triển của khoa học kỹ thuật tác  
động ngành xây dựng, đặc biệt là kết cấu tấm.  
Vì tấm mang lại vẽ thẫm mĩ, khả năng chịu lực  
tốt, đồng thời đảm bảo các yêu cầu kiến trúc  
và thỏa mãn kết cấu của công trình xây dựng  
hiện nay.  
2. TỔNG QUAN NGHIÊN CỨU  
2.1. Các nghiên cứu trong nƣớc  
Những năm gần đây, phương pháp  
phần tử hữu hạn trơn [14] đã được phát triển  
và áp dụng thành công trong việc phân tích  
nhiều bài toán kết cấu khác nhau. Kết quả tính  
toán số cho thấy trong rất nhiều bài toán, PP  
PTHH trơn có độ chính xác và tốc độ hội tụ tốt  
hơn PP PTHH truyền thống. Theo phương  
pháp PTHH trơn, trường biến dạng được trung  
bình trên các miền làm trơn được định nghĩa  
trên miền con của phần tử (CS), hoặc trên  
Kết cấu tấm là một trong những kết cấu  
phổ biến hiện nay, do đặc trưng mỏng, nhẹ,  
khả năng chịu uốn, vượt nhịp lớn, nên tấm đã  
được ứng dụng vào các bộ phận công trình  
như: sàn, tường, vách, làm kết cấu mái, bao  
che cho các công trình; ngoài ra còn rất nhiều  
các ứng dụng khác của tấm phục vụ cho thi  
công xây dựng côngtrình….Thậm chí tấm còn  
64  
miền giữa 2 phần tử chung cạnh (ES) hoặc trên thì kỹ thuật làm trơn trên cạnh thường có ưu  
miền giữa các phần tử chung nút (NS).  
điểm hơn. Thông qua khảo sát ví dụ số, do tính  
không phụ thuộc vào thứ tự đánh số nút của  
phần tử, kỹ thuật khử khóa cắt MITC3 khi kết  
hợp với kỹ thuật làm trơn cho kết quả tốt phần  
tử làm trơn sử dụng kỹ thuật DGS3 hoặc  
MIN3. Ngoài ra, kỹ thuật khử khóa cắt  
MITC3+ mới được phát triển gần đây, sử dụng  
thêm nút nổi nên có khả năng tính toán tốt hơn  
phần tử sử dụng kỹ thuật khử khóa cắt MITC3.  
Đối với phần tử tam giác 3 nút dùng  
phân tích tấm đồng nhất đẳng hướng, phương  
pháp PTHH trơn đã kết hợp với các kỹ thuật  
khử khóa cắt khác nhau như MIN3, DGS3,  
MITC3 và MITC3+ trong các nghiên cứu sau:  
Nguyen-Thoi và cộng sự [15] phát  
triển phương pháp làm trơn trên miền phần tử  
(CS) với kỹ thuật khử khóa cắt MIN3 tạo  
thành phần tử CS-MIN3 dùng phân tích tĩnh  
và tần số của tấm.  
2
Nguyen-Xuan và cộng sự [16] nghiên  
cứu phần tử ES-DGS3 kết hợp giữa làm trơn  
trên cạnh (ES) và khử khóa cắt DSG3 để phân  
tích tĩnh và dao động tấm Reissner-Mindlin.  
4
3
1
Nguyen-Xuan và cộng sự [17] xây  
dựng phần tử NS-DSG3 dùng phương pháp  
làm trơn trên nút (NS) và khử khóa cắt DSG3  
để phân tích tĩnh và dao động tấm Reissner-  
Mindlin.  
Hình 1: Phần tử tam giác 3 nút với  
1 nút nổi (bubble) ở trọng tâm  
MITC3+  
Từ những nhận xét trên và tình hình  
nghiên cứu trong và ngoài nước, luận văn này  
tập trung nghiên cứu sự kết hợp giữa phương  
pháp làm trơn trên cạnh (ES) và kỹ thuật khử  
khóa cắt MITC3+ để xây dựng công thức  
PTHH tấm ES-MITC3+ dùng phân tích tấm  
đồng nhất đẳng hướng.  
Nguyen-Thoi và cộng sự [18] phân  
tích tĩnh và dao động của tấm đồng nhất bằng  
phần tử CS-DSG3 sử dụng phương pháp làm  
trơn trên miền phần tử (CS) và kỹ thuật DSG3.  
Nguyễn Duy Quang [19] trong luận  
văn thạc sĩ đã nghiên cứu phương pháp làm  
trơn trên cạnh (ES) kết hợp với kỹ thuật khử  
khóa cắt MITC3 tạo thành phần tử ES-MITC3  
phân tích tĩnh và dao động của tấm [20].  
2.2. Các nghiên cứu ngoài nƣớc  
Lý thuyết tấm đồng nhất đẳng hướng  
được chia thành lý thuyết tấm mỏng Kirchhoff  
và lý thuyết tấm dày Mindlin. Lý thuyết tấm  
mỏng Kirchhoff giả thuyết bỏ qua biến dạng  
cắt ngoài mặt phẳng, phù hợp với ứng xử thực  
tế của các loại tấm có tỉ số chiều dài / chiều  
dày tấm lớn. Lý thuyết tấm dày Mindlin kể  
đến biến dạng cắt ngoài mặt phẳng và phù hợp  
với ứng xử của các loại tấm có tỉ số chiều dài /  
chiều dày tấm nhỏ [1].  
Nguyễn Văn Dũng [21] trong luận  
văn thạc sĩ đã xây dựng phần tử NS-MITC3  
bằng cách kết hợp phương pháp làm trơn trên  
nút (NS) và kỹ thuật MITC3 để phân tích tĩnh  
và dao động tấm đồng nhất đẳng hướng và tấm  
composite nhiều lớp [22].  
Võ Ngọc Tuyển [23] đã nghiên cứu  
sự kết hợp của phương pháp làm trơn trên  
miền phần tử (CS) và kỹ thuật khử khóa cắt  
MITC3 thì thấy tương tương với việc không  
làm trơn do tất cả trường biến dạng đều là  
hằng số trên miền phần tử. Do đó, tác giả đã sử  
dụng kỹ thuật khử khóa cắt MITC3+ với kỹ  
thuật làm trơn trên miền con của phần tử để  
xây dựng phần tử CS-MITC3+ phân tích tĩnh  
kết cấu tấm [24].  
Để phân tích, tính toán ứng xử tấm  
chịu uốn, các phương pháp giải tích và phương  
pháp số khác nhau đã được đề xuất. Trong đó,  
phương pháp tính toán số được sử dụng nhiều  
nhất và ưu thế nhất là phương pháp phần tử  
hữu hạn (PP PTHH). Phương pháp phần tử  
hữu hạn được bắt nguồn từ những yêu cầu giải  
các bài toán phức tạp về lý thuyết đàn hồi,  
phân tích kết cấu trong xây dựng và kỹ thuật  
hàng không. Hrennikoff, McHenry và Richard  
Trong các nghiên cứu trên thì nếu  
cùng sử dụng một loại kỹ thuật khử khóa cắt  
65  
Courant là những người khởi đầu sự phát triển sau những năm 1950 trong việc phân tích kết  
của phương pháp phần tử hữu hạn vào những cấu khung máy bay và công trình xây dựng đã  
năm 1940 [2]. Sự phát triển chính thức của thu được nhiều kết quả ở Berkeley trong  
phương pháp phần tử hữu hạn được bắt đầu những năm 1960 trong ngành xây dựng.  
vào nửa sau những năm 1950 khi việc phân  
Gần đây, Liu đã đề nghị một phương pháp  
tích kết cấu khung máy bay và công trình xây  
phần tử hữu hạn trơn dựa trên biên (ES-FEM)  
dựng đã thu được nhiều kết quả ở Berkeley  
cho phân tích tĩnh, tự do và dao động của các  
trong những năm 1960. Trong ngành xây  
tấm 2D. Kết quả số có độ chính xác cao đã  
dựng, phương pháp phần tử hữu hạn được sử  
chứng minh rằng ES-FEM có cho ra kết quả  
dụng nhiều trong phân tích tấm vì những ưu  
tốt: (1) mô hình ES-FEM thường được tìm  
điểm nổi trội như kết quả tính toán có độ chính  
thấy độ hội tụ cao và thậm chí nhiều chính xác  
xác cao, mạnh mẽ khi phân tích các tấm phức  
hơn FEM sử dụng phần tử 4 cạnh (FEM-Q4)  
tạp, kết hợp tốt các lý thuyết tấm.  
với cùng số lượng nút; không có kiểu năng  
Xấp xỉ PTHH cho trường chuyển vị lượng ảo nào và do vậy phương pháp này làm  
của tấm mỏng theo lý thuyết Kirchhoff đòi hỏi việc tốt và ổn định theo thời gian hơn cho phân  
hàm liên tục dạng C1, nghĩa là đạo hàm cấp 1 tích động và khi làm theo phương pháp này thì  
của xấp xỉ chuyển vị vẫn liên tục. Trong khi đơn giản hơn, và hiệu quả tính toán tốt hơn  
đó, trường chuyển vị của tấm dày theo lý FEM khi sử dụng cùng số nút.  
thuyết Mindlin chỉ cần dùng hàm xấp xỉ dạng  
3. PHƢƠNG PHÁP NGHIÊN CỨU  
C0. Việc xây dựng hàm xấp xỉ PTHH dạng C0,  
Từ các nghiên cứu trong, ngoài nước  
và đề xuất của một số tác giả trước thì phương  
pháp nghiên cứu theo các bước sau được thực  
hiện:  
đặc biệt là đối với các phần tử đẳng tham số,  
dễ hơn rất nhiều so với xây dựng hàm xấp xỉ  
dạng C1. Tuy nhiên, nếu công thức PTHH của  
tấm dùng hàm xấp xỉ dạng C0 thuần túy thì sẽ  
không tính toán đúng ứng xử khi chiều dày  
tấm mỏng do xấp xỉ dạng C0 không loại bỏ  
được biến dạng cắt ngoài mặt phẳng khi chiều  
dày tấm mỏng. Thậm chí, xấp xỉ PTHH dạng  
C0 còn làm cho biến dạng cắt ngoài mặt phẳng  
tăng khi chiều dày tấm giảm. Điều đó dẫn kết  
kết quả chuyển vị giảm khi chiều dày tấm  
giảm, ngược hoàn toàn với ứng xử thực tế của  
tấm. Hiện tượng này gọi là hiện tượng khóa  
cắt. Vì vậy, nhiều nghiên cứu khác nhau đã đề  
xuất các phương pháp khử khóa cắt cho công  
thức PTHH tấm sử dụng xấp xỉ dạng C0 nhằm  
loại bỏ hiện tượng khóa cắt khi sử dụng các  
phần tử này phân tích tấm mỏng. Các phương  
pháp khử khóa cắt phổ biến như phương pháp  
tích phân giảm (RI) [3] hay tích phân lựa chọn  
(SI) [4], phương pháp giả sử biến dạng tự  
nhiên (ANS) [5], phương pháp giả sử biến  
dạng nâng cao (EAS) [6], phương pháp cải  
tiến biến dạng cắt tấm Mindlin (MIN3) [7],  
phương pháp rời rạc chênh lệch cắt (DSG3)  
[8], hoặc phương pháp nội suy các thành phần  
ten-xơ hỗn hợp (MITC3, MITC3+, MITC4,  
MITC4+, MITC7, MITC9, …) [9]–[13].  
3.1 Nghiên cứu lý thuyết.  
Lý thuyết tấm biến dạng cắt, hay còn  
được gọi là lý thuyết tấm dày chịu uốn, hoặc  
còn gọi là lý thuyết tấm Reissner-Mindlin do  
Levy-Reissner-Hencky-Mindlin thực hiện, dựa  
trên giả thiết của 2 tác giả Reissner-Mindlin,  
giống giả thiết lý thuyết dầm Timoshenko chịu  
cắt. Trong đó, tùy theo tỉ số h/a giữa bề dày  
tấm h và kích thước nhỏ nhất của mặt phẳng  
tấm a mà ta có thể chia tấm thành hai loại sau:  
Tấm dày: h/a > h/5 hoặc tỉ lệ h/a 1/20  
Tấm mỏng: 1/20 h/a h/5 và độ võng lớn  
nhất wmax h/4  
Trong tấm mỏng, các ứng suất màng  
là rất nhỏ so với các ứng suất gây ra bởi sự  
uốn tấm do tải trọng vuông góc tấm gây ra khi  
tấm có độ võng nhỏ. Tuy nhiên nếu tấm mỏng  
có độ võng lớn wmax h/4 các ứng suất do uốn  
này bị ảnh hưởng nhiều bởi các ứng suất  
màng. Khi đó phải tính toán với lý thuyết tấm  
có biến dạng lớn.  
3.1.1. Trƣờng chuyển vị  
Alexander Hrennikoff và Richard  
Courant là những người khởi đầu sự phát triển  
của phương pháp PTHH, được bắt đầu vào nửa  
Theo giả thuyết tấm Mindlin-  
Reissner, khi tấm phẳng chịu tải trọng phân bố  
66  
q thì đoạn thẳng vuông góc mặt trung bình tấm  
sau khi biến dạng sẽ chuyển vị võng xuống 1  
đoạn thẳng w theo phương z (trong đó w là 1  
hàm phụ thuộc x, y) và 2 chuyển vị góc xoay  
3.1.3. Trƣờng ứng suất  
β , β (là 2 góc xoay của đoạn thẳng vuông  
x
y
góc với mặt trung bình của tấm lần lượt quay  
quanh trục y, x) như Hình 2.  
Hình 3: Các thành phần ứng suất của  
Cả 3 thành phần này w, β , β được  
x
y
tấm  
biểu diễn như sau:  
Dựa vào định luật Hooke, theo lý  
thuyết ứng suất phẳng, từ trường biến dạng  
(0.2) cho ra trường ứng suất;  
E
x   
y   
(x vy )  
(vx y )  
1v2  
E
w
w
w
1v2  
(0.3)  
z 0  
E
xy  
xz  
yz  
xy  
2(1v)  
E
Hình 2: Trường chuyển vị trong  
tấm  
xz  
yz  
2(1v)  
E
u zx x, y  
2(1v)  
(0.1)  
v zy x, y  
Trong đó: x, y, z là các ứng suất pháp theo  
phương x, y, z; xy, xz, yz là các ứng suất tiếp  
trên mặt có vectơ pháp tuyến là trục x, y, z; E  
là mô-đun đàn hồi khi kéo nén của vật liệu; ν:  
là hệ số Poisson.  
w w x, y  
Trong đó: u, v, w lần lượt là chuyển vị thẳng  
theo phương x, y, z; β , β lần lượt là góc xoay  
x y  
của đoạn thẳng vuông góc với mặt trung bình  
quanh trục y và trục x với chiều dương qui ước  
được định nghĩa như Hình 2.  
Trường ứng suất (0.3) có thể được  
viết lại dưới dạng ma trận như sau:  
3.1.2. Trƣờng biến dạng  
1 v  
v 1  
0
0
Từ trường chuyển vị (0.1), ta có  
trường biến dạng:  
x   
x   
(0.4)  
E
y  
y  
1v2  
u  
x  
v  
y  
x  
x  
y  
xy  
1v xy  
x   
y   
z  
z  
0 0  
2   
xz   
(0.2)  
E
y  
xz   
(0.5)  
z 0  
yz  
yz  
2(1v)  
y  
u v  
y x  
x  
y  
xy  
z  
Thế quan hệ giữa biến dạng và  
chuyển vị cho bởi (0.2) vào ứng suất (0.4) và  
x  
u w  
z x  
w  
x  
w  
xz  
yz  
x   
y   
(0.5), ta có:  
v w  
x  
x  
y  
x  
x  
y  
z
z
z y  
y  
1
v
0
0
1
v
0
0
x   
E
1v2  
Ez  
1v2  
y  
v
1
   
v
1
Trong đó: x, y, z là các biến dạng dài theo  
phương x, y, z; xy, xz, yz lần lượt là các biến  
dạng trượt trong mặt phẳng Oxy, Oxz, Oyz.  
y  
y  
xy  
1v  
2
1v  
2
0
0
0
0
  
y  
x  
x  
y  
z
y   
   
y  
x  
x  
67  
Ez   
1
v
0
0
M
1 v  
0
0
(0.6)  
(0.7)  
x   
(0.11)  
σ   
v
1
κ
1v2  
M M y D v  
1
κ
1v  
2   
Mxy  
1v  
2   
0
0
0 0  
w  
Trong đó, M = [Mx My Mxy]T và D là độ cứng  
x   
y   
xz   
E
x  
w  
y  
E
τ   
γ
trụ của tấm với  
yz  
2(1v)  
2(1v)  
Eh3  
12(1v2 )  
D   
(0.12)  
Trong đó: = [x y xy]T;  
T  
y  
  
x  
x  
Lực cắt trên 1 đơn vị chiều dài:  
y   
κ   
(0.8)  
(0.9)  
h/2  
x y y  
x  
Q dA 1dz  
xz   
x
xz  
Và  
= [xz yz]T và  
A
h/2  
h/2  
(0.13)  
T  
w  
x  
w  
y  
γ x   
y   
Q dA 1dz  
yz   
y
yz  
A
h/2  
Trong đó, Qx, Qy là lực cắt trên 1 đơn vị dài  
lần lượt nằm trong mặt phẳng Oxz và Oyz.  
3.1.4. Nội lực  
Tách 1 đơn vị dài trên chiều dày h của  
tấm. Khi đó, hợp lực các ứng suất phân bố trên  
Thế ứng suất cắt (0.7) vào lực cắt  
đơn vị này theo chiều dày h được gọi là nội (0.13), ta có véc-tơ lực cắt;  
lực, hay còn gọi là các thành phần ứng lực của  
tấm.  
Q
   
x
1 0  
0 1  
Q   
kGh  
γ
(0.14)  
   
Qy  
   
Trong đó, G = E/[2(1+)] là mô-đun đàn hồi  
trượt, k = 5/6 là hệ số hiệu chỉnh cắt dùng để  
hiệu chỉnh sự phân bố theo chiều dày tấm của  
các ứng suất tiếp xz, yz từ dạng không phải  
hàm bậc nhất sang dạng hàm bậc nhất.  
4.1 Lập trình tính toán.  
Hình 4: Các thành phần nội lực  
4.1.1. Công thức PTHH tấm tam giác  
của tấm  
Moment trên 1 đơn vị chiều dài:  
h/2  
zx 1dz  
M zdA   
x
x
A
h/2  
(0.10)  
h/2  
zy 1dz  
M zdA   
y
y
A
h/2  
h/2  
Mxy  
zxy 1dz  
zdA   
xy  
A
h/2  
Hình 5: Phần tử tấm tam giác 3 nút với  
Trong đó, Mx, My lần lượt là moment uốn trên  
1 đơn vị dài mặt cắt quay quanh trục y, x có  
chiều dương như Hình 4; Mxy là moment xoắn  
trên 1 đơn vị dài mặt cắt quay quanh trục trục  
x.  
chiều dương bậc tự do tại nút  
Rời rạc mặt trung bình của tấm bằng  
các phần tử tam giác 3 nút. Trường chuyển vị  
ở (0.1) được xấp xỉ thông qua chuyển vị nút  
của phần tử và nút nổi [10]  
Thế ma trận ứng suất (0.6) vào  
moment (0.10), ta có:  
68  
3
h 2  
h 2  
εTσdz dAk  
γT τdz dA qwdA  
w N w  
i
i
   
   
i1  
A
h 2  
A
h 2  
A
e   
e   
e
4
(0.15)  
(0.21)  
N   
x
i
yi  
i1  
Dựa vào mối quan hệ giữa ứng suất  
4
và biến dạng cho bởi công thức (0.6) và (0.7),  
nguyên lý công ảo (0.21) có thể viết lại:  
  N   
y
i
xi  
i1  
Trong đó, wi, xi, yi lần lượt là độ võng, góc  
xoay quay quanh trục x và trục y của nút i.  
Chiều dương của các chuyển vị nút được định  
nghĩa ở Hình 5. Ni là các hàm dạng định nghĩa  
trong hệ tọa độ tự nhiên , như sau:  
1
v
0
0
h
2
h
2
Ez  
1v2  
E
zκT  
v
1
κdz dAk  
γT  
γdz dA qwdA  
   
   
2(1v)  
A
h  
2
A
h  
2
A
e
e
e
1v  
2
0
0
κT DbκdAk γT DsγdA qwdA  
N 19 1  
1
A
A
A
e
e
e
(0.16)  
N 9 1  
2
Trong đó,  
N 9 1  
3
N 27 1  
4
1
v
0
0
1
v
0
0
1
v
0
0
h
2
Ez  
Eh3  
Db   
z
v
1
dz   
v
1
D v  
1
1v2  
12(1v2 )  
Thế (0.15) vào (0.8) và (0.9), ta được;  
h  
2
1v  
2
1v  
2
1v  
2
0
0
0
0
0
0
i   
Ni  
x  
N  
x  
yi  
0
0
0
(0.23)  
w
i   
4
4
4
Ni  
y  
Ni  
y  
h 2  
κ   
xi  
0
xi  
Bbd  
1 0  
0 1  
1 0  
0 1  
  
  
i
i
Ds k  
G
dz kGh  
i1  
i1  
i1  
(0.24)  
yi  
Ni  
y  
Ni  
x  
Ni N  
h 2  
i   
yi  
xi  
0
x  
y  
Thế (0.15), (0.17) và (0.18) vào  
(0.22), ta có phương trình cân bằng rời rạc  
theo phương pháp phần tử hữu hạn;  
(0.17)  
Ni  
x  
N  
x  
Ni  
i
w
Niyi  
w
0
Ni  
i   
i
4
4
4
γ   
xi  
Bsd  
i
i
4
4
4
4
3
dT N 0 0 TdA  
Ni  
i1  
i1  
i1  
(
Bbd )T Db ( Bbd )dA  
(
Bsd )T Ds ( Bsd )dA   
q  
Nixi  
w
Ni  
0
yi  
i
i
i
i
i
i
i
i
i
i
i
i1  
i1  
i1  
i1  
i1  
x  
 x  
A
e
A
e
A
e
kb ks d f k d f  
(0.18)  
(0.25)  
e   
e
e
e
e
e
e
Trong đó, di = [wi xi yi]T với i = 1..3 và d4 =  
k kb ks  
kb , ks  
lần lượt là  
[0 x4 y4]T và  
Trong đó,  
với  
e
e
e
e
e
ma trận độ cứng uốn và cắt của phần tử e có;  
i   
N  
x  
bT  
i
b
b
j
0
0
b
e,ij  
(0.26)  
k B D B dA  
(0.19)  
(0.20)  
A
e
Ni  
y  
Bbi 0  
0
s
e,ij  
sT  
i
s
s
j
(0.27)  
k B D B dA  
Ni N  
i   
A
e
0
x  
y  
fe là véc-tơ lực của phần tử có;  
N  
i
f q N 0 0 T dA  
0
Ni  
x  
Bs   
(0.28)  
e,i  
i
i
A
Ni  
 x  
e
Ni  
0
T
T
T
Và de = [d1 d2 d3 0 x4 y4]T là véc-tơ  
chuyển vị của phần tử.  
Theo nguyên lý công ảo giữa năng  
lượng biến dạng đàn hồi của tấm chịu uốn và  
ngoại lực q phân bố tác dụng vuông góc với  
mặt trung bình của tấm thì [25]  
4.1.2. Công thức PTHH tấm tam giác  
MITC3+  
Do các thành phần chuyển vị nút  
được xấp xỉ bằng các hàm dạng bậc thấp nên  
đạo hàm bậc 1 của xấp xỉ chuyển vị là hằng số.  
69  
Vì vậy, biến dạng cắt ngoài mặt phẳng ở  
công thức (0.9) không thể bằng 0. Tức là, khi  
tấm mỏng biến dạng cắt ngoài mặt phẳng vẫn  
tồn tại. Điều này không phù hợp với ứng xử  
thực tế của tấm mỏng và dẫn đến kết quả, tấm  
càng mỏng thì độ võng càng lớn. Đây chính là  
hiện tượng khóa cắt xảy ra với các phần tử  
dùng hàm xấp xỉ chuyển vị bậc thấp.  
Hình 6: Vị trí điểm buộc của phương pháp  
Hiện tượng khóa cắt đã được nhiều  
tác giả nghiên cứu và đề xuất các phương pháp  
giải quyết khác nhau. Trong luận văn này,  
phương pháp nội suy các thành phần ten-xơ  
hỗn hợp MITC3+ [10] được sử dụng để loại  
bỏ hiện tượng khóa cắt của phần tử tấm tam  
giác 3 nút có thêm 1 nút nổi ở trọng tâm phần  
tử.  
MITC3+ [10]  
Bảng 1: Tọa độ các điểm buộc của phương  
pháp khử khóa cắt MITC3+ với d = 1/10000  
[10]  
Điểm buộc  
1/6  
2/3  
1/6  
2/3  
1/6  
A
B
C
D
E
F
1/6  
Ý tưởng cơ bản của phương pháp  
MITC3+ xấp xỉ lại biến dạng cắt ngoài mặt  
phẳng thông qua các giá trị của biến dạng cho  
bởi (0.18) tính tại các điểm buộc có vị trí như  
Hình 6 và tọa độ cho trong Bảng 1. Xấp xỉ của  
biến dạng cắt ngoài mặt phẳng là;  
1/3 + d  
1/3 – 2d  
1/3 + d  
1/3 – 2d  
1/3 + d  
1/3 + d  
Do đó, ma trận độ cứng cắt tính theo  
phương pháp MITC3+ có dạng;  
2
1
1
3
1
3
1
B   
B
s
sT  
i
s
s
j
ˆ
(0.31)  
ˆ
xz xz   
C C c(31)  
ˆ
ˆ
ˆ
yz   
yz   
xz  
ke,ij  
B D B dA  
3
2
2
1
3
1
A
e
A   
A
C C c(31)  
ˆ
ˆ
yz xz   
Và ma trận độ cứng phần tử MITC3+  
yz   
yz   
xz  
3
2
3
được tính bởi;  
(0.29)  
b
s
bT  
i
b
b
j
sT  
i
s
s
j
ˆ
ˆ
B D B dAB D B dA  
ˆ  
ke,ij ke,ij ke,ij  
c     F E  
F
xz  
D
A
A
Với  
. Trong đó,  
ˆ
e
e
xz
  
yz
  
yz  
(0.32)  
xIz , yIz là giá trị biến dạng cắt ngoài mặt  
phẳng tại các điểm buộc I = A, B, C, D, E, F.  
4.1.3. Công thức phần tử tấm tam giác ES-  
MITC3+  
Thế các giá trị biến dạng cắt  
ngoài mặt phẳng tại các điểm buộc được tính  
theo chuyển vị nút phần tử di cho bởi (0.18)  
vào xấp xỉ biến dạng cắt ngoài mặt phẳng  
trong công thức (0.29), ta có thể tìm được  
quan hệ giữa biến dạng cắt ngoài mặt phẳng  
xấp xỉ theo kỹ thuật MITC3+ và chuyển vị nút  
phần tử di như sau:  
Trong luận văn này, công thức phần  
tử tấm tam giác ES-MITC3+ được thiết lập  
bằng cách làm trơn trường biến dạng uốn trên  
miền làm trơn (k) nối trọng tâm 2 phần tử  
chung cạnh với 2 nút của cạnh chung này như  
Hình 7 và giữ nguyên trường biến dạng cắt  
ngoài mặt phẳng xấp xỉ theo phương pháp  
MITC3+.  
4
T  
s
i
ˆ
(0.30)  
ˆ
ˆ
ˆ
γ xz yz  
B d  
Theo phương pháp làm trơn trên cạnh  
(ES) [14], trường biến dạng uốn xấp xỉ như  
sau:  
i
i1  
1
A
(0.33)  
κd  
κ
(k )  
Trong đó, A(k) là diện tích miền làm trơn Ω(k).  
70  
N(k )  
1
A
(0.37)  
κ
B
i1  
Vậy, ma trận độ cứng uốn sau khi  
được làm trơn trên miền chung cạnh của 2  
phần tử có dạng  
Hình 7: Phân chia miền làm trơn trên cạnh phần  
A
k
(0.38)  
tử  
Thế (0.17) vào (0.33), ta được;  
Để cải thiện độ chính xác, ma trận Ds  
được hiệu chỉnh bằng cách đưa vào hệ số ổn  
định như sau [26]  
Nn(k )  
Nn(k )  
  
1
A
1
A(k)  
(0.34)  
κ
Bbd d   
Bbdd  
i
i
i
i
(k )  
(k )  
i1  
i1  
(k)  
kGh3  
1 0  
0 1  
Trong đó, Nn là số nút của miền làm trơn  
s
ˆ
(0.39)  
D   
Ω(k), Nn = 3 đối với miền làm trơn của cạnh  
(k)  
h2 le2  
làm trơn nằm trên biên và Nn(k) = 4 đối với các  
Trong đó, le là chiều dài lớn nhất của cạnh  
miền làm trơn còn lại.  
phần tử.  
Áp dụng định lý Green, ta có:  
Và ma trận độ cứng phần tử ES-  
N  
x  
N  
y  
i d  Ninxd;  
i d  Ninyd  
MITC3+ được tính bởi  
(k )  
(k )  
(k )  
(k )  
1
A
s
j
(0.40)  
dA  
kij k
(0.35)  
A
e
Với nx, ny lần lượt là hình chiếu theo phương x  
và y của véc-tơ n pháp tuyến với biên (k) của  
miền (k) như Hình 8. Với các hàm dạng bậc  
thấp, các tích phân (0.35) được tính chính xác  
bằng cách sử dụng 2 điểm Gauss trên mỗi  
cạnh của biên miền làm trơn (k).  
Hiện tại, độ cứng phần tử ES-  
MITC3+ ke liên quan bậc tự do phần tử de =  
T
T
T
[d1 d2 d3 0 x4 y4]T. Để loại bỏ bậc tự do w4  
= 0, các hàng và cột của ma trận ke được loại  
bỏ (hàng 10 và cột 10). Các chuyển vị x4 y4  
được tính theo chuyển vị nút phần tử d1, d2, d3  
theo phương pháp nén tĩnh như sau:  
11  
ke k1e2 de  
   
   
fe  
k d f   
   
e
e
e
   
ke21 ke22 d4  
0
   
(0.41)  
   
12  
11  
ke de ke d4 fe (*)  
ke21de ke22d4 0 (**)  
1  
Hình 8: Véc-tơ pháp tuyến  
với biên của miền làm trơn  
Từ (**) suy ra d   k22 ke21de (0.42)  
4
e
Và thế vào (*), ta có;  
Dùng (0.35) và công thức Bib cho bởi  
(0.19), ta có thể viết  
1  
k1e1 k1e2 k22 ke21 de fe kede fe  
e
(0.43)  
0
0
Ninxd  
0
( k )  
1  
(0.36)  
11  
12  
ke ke ke k22 ke21  
Trong đó,  
(0.44)  
b
e
Bi d  0  
Ninyd  
Ninxd  
B
  
( k )  
( k )  
Với k11,k22 ,k12 k21 T lần lượt là các thành  
0
N n d  
  
i
y
e
e
e
e
( k )  
( k )  
phần của ma trận độ cứng ke liên quan bậc tự  
do tại các nút đỉnh phần tử, nút nổi và tương  
tác giữa bậc tự do nút ở đỉnh và nút nổi;  
Thế (0.36) vào (0.34), ta có quan hệ  
giữa biến dạng uốn trơn và chuyển vị nút phần  
tử  
T  
de w x1 y1 w2 x2 y2 w3 x3 y3  
là  
1
các bậc tự do tại nút đỉnh của phần tử và  
71  
T
toàn giống nhau, phần tử ES-MITC3+ đã vượt  
qua điều kiện patch test vì có khả năng biểu  
diễn chính xác trường chuyển vị và biến dạng,  
ứng suất hoặc nội lực.  
(0.45)  
f q N 0 0 N 0 0 N 0 0 dA  
e
1
2
3
A
e
Như vậy, ma trận độ cứng phần tử  
ES-MITC3+ có kích thước 9x9 và liên quan  
các bậc tự do tại nút đỉnh của phần tử. Kết quả,  
ma trận độ cứng phần tử ES-MITC3+ được lập  
ghép vào ma trận độ cứng kết cấu thông qua  
các bậc tự do tại nút đỉnh phần tử.  
Bảng 3: Kết quả chuyển vị và mô-men tại nút 5  
của bài toán patch test  
Phƣơng pháp  
w5  
Mx5  
My5  
Mxy5  
x5  
y5  
ES-MITC3+ 0.6422 1.13 -0.64 -0.0111-0.0111-0.0033  
5.1 So sánh kết quả tính toán với các  
nghiên cứu khác.  
Lời giải giải tích 0.6422 1.13 -0.64 -0.0111-0.0111-0.0033  
5.1.2. Tấm hình vuông ngàm 4 cạnh chịu  
tải phân bố đều  
Trong bào báo này, một số bài toán tấm  
đồng nhất đẳng hướng được giải bằng phần tử  
đề xuất ES-MITC3+ để xác định độ võng và  
mô-men uốn. Kết quả tính toán được so sánh  
với các phần tử tam giác 3 nút cùng loại như  
ES-DSG3 [27], MITC3 [19], ES-MITC3 [19]  
hoặc phần tử tứ giác 4 nút MITC4 [27]. Trong  
tất cả các kết quả số của phần tử ES-MITC3+,  
hệ số ổn định = 0.1 được chọn giống như hệ  
số ổn định đã được dùng cho các phần tử ES-  
DSG3 và ES-MITC3.  
Cho tấm vuông cạnh a, dày h, ngàm 4  
cạnh biên (w = θ = θ = 0 tại biên) và chịu tải  
x
y
trọng đều q = 1 như Hình 10. Vật liệu làm tấm  
đồng nhất, đẳng hướng, và có E = 1092000 và  
= 0.3. Độ võng wc và mô-men Mc tại tâm  
tấm được xác định trong 2 trường hợp: tấm  
mỏng có tỉ số h/a = 0.001 và tấm dày có h/a =  
0.1. Để so sánh với kết quả tham khảo, độ  
võng và mô-men tại tâm tấm được chuẩn hóa  
như sau:  
5.1.1. Bài toán patch test  
100D  
qL4  
10  
qL2  
Để kiểm tra khả năng xấp xỉ trường  
biến dạng và ứng suất của phần tử ES-  
MITC3+, xét tấm chữ nhật dày h = 0.01 có tọa  
độ nút và lưới phần tử như Hình 9. Tấm làm  
bằng vật liệu có mô-đun đàn hồi E = 105 và hệ  
số Poisson = 0.25. Tấm chịu độ võng cưỡng  
bức w = (1 + x + 2y + x2 + xy + y2 ) / 2.  
wc wc  
; Mc Mc  
y
q=1  
q=1  
x
1 (0.2; 0.12)  
4 (0.44; 0.12)  
a
a
Hình 10: Tấm vuông ngàm 4 cạnh  
5 (0.3; 0.6)  
chịu tải phân bố đều  
Tấm được chia lưới 2NN phần tử  
tam giác 3 nút đều như Hình 11. Trong đó, N =  
4, 8, 10, 12 và 16 là số phần tử trên mỗi cạnh  
của tấm. Kết quả tính toán độ võng và mô-men  
chuẩn hóa bằng phần tử ES-MITC3+ ứng với  
các loại lưới khác nhau cho bài toán tấm mỏng  
h/a = 0.001 và tấm dày h/a = 0.1 được thể hiện  
độ hội tụ của độ võng và mô-men chuẩn hóa  
tại tâm tấm được biểu diễn trong biểu đồ Hình  
12 và Hình 13.  
2 (0.2; 0)  
3 (0.44; 0)  
Hình 9: Tọa độ nút phần tử  
của bài toán patch test  
Theo lý thuyết tấm mỏng, bỏ qua biến  
dạng cắt ngoài mặt phẳng, từ độ võng w tính  
được góc xoay quanh trục x, y tương ứng  
x  w y, y  w x , và các thành phần  
biến dạng, ứng suất hoặc mô-men. Kết quả  
tính toán chuyển vị và nội lực tại nút 5 có tọa  
độ x = 0.3 và y = 0.6 bằng lời giải lý thuyết  
tấm mỏng và phần tử ES-MITC3+ được cho  
trong Bảng 2  
Với kết quả chuyển vị và mô-men của  
phần tử ES-MITC3+ và lời giải lý thuyết hoàn  
72  
khảo, độ võng và mô-men tại tâm tấm được  
chuẩn hóa như sau:  
100D  
qL4  
10  
qL2  
wc wc  
; Mc Mc  
(a)  
(b)  
y
Hình 11: Mô phỏng tấm vuông  
ngàm 4 cạnh bằng (a) lưới N = 4 và  
(b) lưới N = 16  
Kết quả được thể hiện qua hình ảnh sau:  
q=1  
q=1  
x
a
a
Hình 14: Mô hình tấm vuông tựa đơn 4 cạnh  
Tấm được chia lưới 2NN phần tử  
tam giác 3 nút như HìnhHình15. Trong đó, N  
= 4, 8, 10, 12 và 16 là số phần tử trên mỗi cạnh  
của tấm. Kết quả tính toán độ võng và mô-men  
chuẩn hóa bằng phần tử ES-MITC3+ ứng với  
các loại lưới khác nhau cho bài toán tấm mỏng  
h/a = 0.001 và tấm dày h/a = 0.1. Hội tụ của  
độ võng và mô-men chuẩn hóa tại tâm tấm  
được biểu diễn trong Hình16 và Hình 17.  
(a)  
(b)  
Hình 12: Độ hội tụ của độ võng chuẩn hóa tại tâm  
tấm khi(a) h/a = 0.001 và (b) h/a = 0.1  
(a)  
(b)  
Hình 13: Độ hội tụ của mô-men chuẩn hóa tại tâm  
tấm khi (a) h/a = 0.001 và (b) h/a = 0.1  
Từ biểu đồ Hình 12 cho thấy đối với  
tấm mỏng và tấm dày phần tử ES-MITC3+  
cho ra giá trị độ võng hội tụ đến lời giải giải  
tích khi lưới thay đổi từ thô cho đến mịn  
nhưng chưa tốt bằng kết quả cho bởi các phần  
tử ES-MITC3, MITC4, MITC3, ES-DSG3.  
Tuy nhiên, biểu đồ Hình 13 cho thấy khi so  
sánh với lời giải giải tích, giá trị mô-men tại  
tâm tấm cho bởi phần tử ES-MITC3+ tốt hơn  
kết quả của các phần tử ES-MITC3, MITC4,  
MITC3, ES-DSG3.  
(a)  
(b)  
Hình15: Mô phỏng tấm vuông tựa đơn 4  
cạnh bằng (a) lưới N = 4 và (b) lưới N =  
16  
Kết quả được thể hiện qua hình ảnh sau:  
5.1.3. Tấm hình vuông tựa đơn 4 cạnh chịu  
tải phân bố đều  
Cho tấm vuông cạnh a, dày h, tựa đơn  
(a)  
(b)  
4 cạnh biên (w = θ = θ = 0 tại biên) và chịu  
x
y
Hình16: Độ hội tụ của độ võng chuẩn hóa tại tâm tấm  
tải trọng phân bố đều q = 1 như Hình 14. Vật  
liệu làm tấm đồng nhất, đẳng hướng, và có E =  
1092000 và = 0.3. Độ võng wc và mô-men  
Mc tại tâm tấm được xác định trong 2 trường  
hợp: tấm mỏng có tỉ số h/a = 0.001 và tấm dày  
có h/a = 0.1. Để so sánh với kết quả tham  
khi(a) h/a = 0.001 và (b) h/a = 0.1  
73  
(a)  
(b)  
Hình 17: Độ hội tụ của mô-men chuẩn hóa tại tâm tấm  
khi(a) h/a = 0.001 và (b) h/a = 0.1  
Từ biểu đồ Hình16 cho thấy, đối với  
tấm mỏng và tấm dày phần tử ES-MITC3+  
cho ra kết quả giá trị độ võng tại tâm tấm hội  
tụ đến lời giải giải tích tốt hơn kết quả cho bởi  
các phần tử ES-MITC3, MITC3, MITC3+,  
ES-DSG3. Tuy nhiên, giá trị độ võng cho bởi  
phần tử ES-MITC3+ không tốt bằng phần tử  
tứ giác 4 nút MITC4.  
(a)  
(b)  
Hình 19: Mô phỏng tấm hình thoi Razzaque  
tựa đơn ở cạnh trên và cạnh dưới, 2 cạnh  
bên tự do bằng (a) lưới N = 4 và (b) lưới N  
= 16  
Kết quả được thể hiện qua hình ảnh sau:  
Đối với mô-men tại tâm tấm vuông  
tựa đơn 4 cạnh, Hình 17 cho thấy, phần tử ES-  
MITC3+ cho kết quả tốt nhất so với các phần  
tử khác trong trường hợp tấm mỏng. Nhưng  
đối với tấm dày phần tử ES-MITC3+ cho kết  
quả không tốt bằng phần tử ES-MITC3 nhưng  
vẫn tốt hơn các phần tử còn lại.  
(a)  
(b)  
Hình 20: Độ hội tụ của độ võng (a) và mô-men (b)  
chuẩn hóa tại tâm tấm Razzaque tựa đơn ở cạnh trên  
và cạnh dưới, 2 cạnh bên tự do khi h/a = 0.001  
5.1.4. Tấm hình thoi tựa đơn ở cạnh trên  
và cạnh dƣới, 2 cạnh bên tự do và chịu tải  
phân bố đều (tấm Razzaque)  
Từ biểu đồ Hình 20 cho thấy trong  
trường hợp này khi tấm hình thoi phải chia  
lưới bằng các phần tử không phải tam giác  
vuông cân như ở ví dụ 4.2 và 4.3, độ võng và  
mô-men chuẩn hóa cho bởi phần tử ES-  
MITC3+ hội tụ đến lời giải bằng phương pháp  
sai phân hữu hạn [28] với sai số ít hơn kết quả  
của các phần tử MITC4, MITC3, ES-MITC3.  
Đặc biệt, phần tử ES-MITC3+ đã cải thiện  
đáng kể độ chính xác của mô-men tại tâm tấm  
so với phần tử MITC3.  
Xét tấm hình thoi Razzaque cạnh a,  
dày h với h/a = 0.001. Tấm có cạnh trên, cạnh  
dưới tựa đơn và 2 cạnh bên tự do nghiêng 1  
góc = 60º. Tấm chịu tải trọng phân bố đều q  
= 1 như Hình18. Vật liệu làm tấm đồng nhất,  
đẳng hướng và có E = 1092000, = 0.3. Để  
so sánh với kết quả tham khảo, độ võng wc và  
mô-men My tại tâm tấm được chuẩn hóa như  
sau:  
100D  
qL4  
10  
qL2  
5.1.5. Tấm hình thoi tựa đơn 4 cạnh chịu  
tải phân bố đều (tấmMorley)  
wc wc  
; My My  
y
y
q=1  
q=1  
30°  
x
q=1  
q=1  
a
a
Hình 21: Mô hình tấm hình thoi Morley  
x
Xét tấm hình thoi Morley cạnh a, dày  
h với h/a = 0.01 như Hình 21. Tấm tựa đơn 4  
cạnh và chịu tải trọng phân bố đều q = 1. Vật  
liệu làm tấm đồng nhất, đẳng hướng và có E =  
1092000, = 0.3. Độ võng wc và mô-men  
chính nhỏ nhất McMIN tại tâm tấm được được  
chuẩn hóa để so sánh với tài liệu tham khảo  
như sau:  
a
Hình18: Mô hình tấm hình thoi Razzaque  
Tấm được chia lưới 2NN phần tử  
tam giác 3 nút đều như Hình Hình 19. Trong  
đó, N = 2, 4, 6, 8, 12 và 16 là số phần tử trên  
mỗi cạnh của tấm. Kết quả tính toán độ võng  
và mô-men chuẩn hóa bằng phần tử ES-  
MITC3+. Sự hội tụ của độ võng và mô-men  
chuẩn hóa tại tâm tấm khi N tăng dần từ N = 2  
đến N = 16 được biểu diễn trong Hình 20.  
74  
0.001. Vật liệu làm tấm đồng nhất, đẳng  
hướng và có E = 1092000, = 0.3. Độ võng w  
tại điểm A và B của tấm (xem HìnhHình 24)  
được chuẩn hóa bởi w/(qL4/Eh3) để tiện so  
sánh với các kết quả tham khảo.  
1000D  
100  
wc wc  
; Mc McMIN  
4
2
qL  
qL  
Tấm được chia lưới 2NN phần tử  
tam giác 3 nút đều hoặc không đều như Hình  
22. Trong đó, N = 4, 8, 10, 12 và 16 là số phần  
tử trên mỗi cạnh của tấm. Kết quả tính toán độ  
Tấm được chia lưới 2NN phần tử  
võng và mô-men chuẩn hóa bằng phần tử ES- tam giác 3 nút đều Hình 25. Trong đó, N = 2,  
MITC3+. Độ hội tụ của độ võng và mô-men  
4, 8 và 16 là số phần tử trên mỗi cạnh của tấm.  
chuẩn hóa tại tâm tấm được biểu diễn trong  
Hình 23.  
Kết quả tính toán độ võng chuẩn hóa tại điểm  
A và điểm B bằng phần tử ES-MITC3+. Độ  
hội tụ của độ võng chuẩn hóa tại điểm A và  
điểm B của tấm được biểu diễn trong Hình 26,  
Hình 27 và Hình 28.  
(a)  
(b)  
Hình 22: Mô phỏng tấm hình thoi Morley  
tựa đơn chịu tải phân bố đều bằng(a) lưới  
đều với N = 4 và (b) lưới không đều với N  
= 16  
Kết quả được thể hiện qua hình ảnh sau:  
(a)  
(b)  
Hình 25: Mô phỏng tấm hình thoi ngàm 1  
cạnh, chịu tải phân bố đều bằng (a) lưới N =  
4 và (b) lưới N = 16  
Kết quả được thể hiện qua hình ảnh sau:  
(a)  
(b)  
Hình 23: Độ hội tụ của (a) độ võng và (b) mô-men  
chính nhỏ nhất chuẩn hóa tại tâm tấm khi h/a =  
0.01  
Từ Hình 23 cho thấy, so với lời giải  
của Morley [29] độ chính xác của độ võng và  
mô-men chính nhỏ nhất chuẩn hóa tại tâm tấm  
cho bởi phần tử ES-MITC3+ chỉ tốt hơn phần  
tử MITC3.  
(a)  
(b)  
Hình 26: Độ hội tụ của độ võng chuẩn hóa tại (a)  
điểm A và (b) tại điểm B với = 20° khi h/a =  
0.001  
5.1.6. Tấm hình thoi ngàm 1 cạnh chịu tải  
phân bố đều  
y
A
B
q=1  
q=1  
x
a
a
Hình 24: Mô hình tấm hình thoi ngàm 1  
(a)  
(b)  
cạnh, chịu tải phân bố đều  
Hình 27: Độ hội tụ của độ võng chuẩn hóa tại (a)  
điểm A và (b) tại điểm B với = 40° khi h/a =  
0.001  
Cho tấm hình thoi ngàm 1 cạnh, chịu  
tải phân bố đều như Hình 24. Tấm dày h, cạnh  
a và có 2 cạnh bên nghiêng với phương thẳng  
đứng góc β = 20º, 40º hoặc 60º. Tỉ số h/a =  
75  
Hình 30: Mô hình chia lưới 1/4 tấm  
tròn sử dụng (a) 6, (b) 24, (c) 54 và (d)  
96 phần tử tam giác  
(a)  
(b)  
Kết quả được thể hiện qua hình ảnh sau:  
Hình 28: Độ hội tụ của độ võng chuẩn hóa tại (a)  
điểm A và (b) tại điểm B với = 60° khi h/a =  
0.001  
Từ biểu đồ Hình 26, Hình 27 và Hình  
28 cho thấy, độ võng chuẩn hóa tại điểm A và  
B cho bởi phần tử ES-MITC3+ tương tự với  
các phần tử THS, QHS, MITC3, ES-MITC3  
và có kết quả sai số tương đối lớn so với lời  
giải của Reissner và Stein [29], đặc biệt trong  
trường hợp tại điểm B với = 60o.  
(a)  
(b)  
Hình 31: Độ hội tụ của độ võng chuẩn hóa tại  
tâm tấm tròn ngàm khi (a) h/R = 0.02 và (b) h/R  
= 0.2  
5.1.7. Tấm hình tròn ngàm chịu tải phân  
bố đều  
(a)  
(b)  
Cho tấm tròn bán kính R = 5, dày h  
với h/R = 0.02 hoặc 0.2. Tấm ngàm theo chu vi  
và chịu tải trọng phân bố đều q = 1 như Hình  
Hình 29. Vật liệu làm tấm đồng nhất, đẳng  
hướng và có E = 1092000, = 0.3.  
Hình 32: Độ hội tụ của mô-men chuẩn hóa tại tâm  
tấm tròn ngàm khi (a) h/R = 0.02 và (b) h/R = 0.2  
Từ biểu đồ Hình 31, cho thấy, so sánh  
với lời giải giải tích [1], cả 2 trường hợp tấm  
mỏng (h/R = 0.02) và tấm dày (h/R = 0.2) từ  
lưới thô đến lưới mịn phần tử ES-MITC3+ cho  
ra kết quả giá trị độ võng và mô-men rất tốt.  
Kết quả độ võng và mô-men tại tâm tấm của  
phần tử ES-MITC3+ tốt hơn các phần tử tam  
giá 3 nút ES-DSG3, MITC3 và ES-MITC3,  
ngoại trừ trường hợp mô-men tại tâm tấm khi  
h/R = 0.2 giải bằng phần tử ES-MITC3.  
y
q=1  
q=1  
x
R
Hình 29: Mô hình tấm tròn với điều  
kiện biên ngàm  
Sử dụng tính chất đối xứng, 1/4 tấm  
được chia lưới tam giác với 6, 24, 54 hoặc 96  
phần tử như Hình 30. Kết quả tính toán độ  
võng và mô-men tại tâm tấm bằng phần tử ES-  
MITC3+ ứng với các loại lưới khác nhau. Độ  
hội tụ của độ võng và mô-men tại tâm tấm  
được biểu diễn trong Hình 31.  
6.1 Nhận xét kết quả và kết luận.  
Bài báo đã nói về cách phát triển công  
thức phần tử ES-MITC3+ bằng cách làm trơn  
biến dạng uốn bằng phương pháp PTHH trơn  
trên cạnh (ES) và xấp xỉ biến dạng cắt bằng kỹ  
thuật MITC3+ để khử hiện tượng khóa cắt.  
Phần tử ES-MITC3+ có sử dụng nút nổi để  
76  
tăng bậc xấp xỉ của góc xoay nhưng vẫn là techniques in the finite element analysis of  
phần tử tam giác 3 nút khi kết nối phần tử plates,” Nuclear Engineering and Design, vol.  
trong ma trận độ cứng tổng thể do các góc 46, no. 1, pp. 203–222, Mar. 1978.  
xoay tại nút nổi đã được tính theo chuyển vị  
[5] G. M. Kulikov and S. V. Plotnikova, “A  
tại các nút đỉnh nhờ kỹ thuật nén tĩnh. Biến  
family of ANS four-node exact geometry shell  
dạng uốn được trung bình trên miền làm trơn  
elements in general convected curvilinear  
định nghĩa bởi các đoạn thẳng nối 2 trọng tâm  
coordinates,” Int. J. Numer. Meth. Engng., vol.  
với 2 đỉnh của cạnh chung giữa 2 phần tử. Do  
83, no. 10, pp. 1376–1406, Sep. 2010.  
đó, ma trận độ cứng uốn của phần tử có thể  
[6] U. Andelfinger and E. Ramm, “EAS-  
elements for two-dimensional, three-  
dimensional, plate and shell structures and  
their equivalence to HR-elements,” Int. J.  
Numer. Meth. Engng., vol. 36, no. 8, pp.  
1311–1337, Apr. 1993.  
chuyển từ tích phân mặt sang tích phân đường  
trên biên của miền làm trơn. Phần tử ES-  
MITC3+ có khả năng biểu diễn chính xác  
trường chuyển vị và biến dạng, ứng suất hoặc  
nội lực khi giải bài toán patch test. Tuy nhiên,  
do có sử dụng thêm nút nổi, tức là ma trận độ  
cứng không phải là hằng số trên phần tử, phần  
tử ES-MITC3+ có thời gian tính toán lâu hơn  
các phần tử tam giác 3 nút khác có ma trận độ  
cứng là hằng số được viết ở dạng tường minh.  
[7] A. Tessler and T. J. R. Hughes, “A three-  
node Mindlin plate element with improved  
transverse shear,” Computer Methods in  
Applied Mechanics and Engineering, vol. 50,  
no. 1, pp. 71–101, Jul. 1985.  
Phần tử đề xuất ES-MITC3+ được sử  
dụng để phân tích một số bài toán tấm hình  
vuông, hình tròn, hình thoi chịu tải trọng phân  
bố đều. Tương tự các phần tử tam giác 3 nút  
có làm trơn cùng loại như ES-DSG3, ES-  
MITC3, kết quả cho bởi phần tử ES-MITC3+  
hội tụ đến lời giải giải tích, ngoại trừ trường  
hợp tấm hình thoi ngàm 1 cạnh có độ hội tụ  
chưa tốt. So với các công thức PTHH ES-  
DSG3, ES-MITC3, phần tử ES-MITC3+  
thường cho kết quả tính toán mô-men tốt hơn  
vì phần tử có sử dụng nút nổi trong xấp xỉ góc  
xoay. Trong một số trường hợp, phần tử đề  
xuất có thể cho kết quả tốt hơn phần tử tứ giác  
4 nút MITC4. Phần tử ES-MITC3+ thường  
cho kết quả lời giải tấm mỏng tốt hơn lời giải  
tấm dày.  
[8] K.-U. Bletzinger, M. Bischoff, and E.  
Ramm, “A unified approach for shear-locking-  
free triangular and rectangular shell finite  
elements,” Computers & Structures, vol. 75,  
no. 3, pp. 321–334, Apr. 2000.  
[9] P.-S. Lee and K.-J. Bathe, “Development  
of MITC isotropic triangular shell finite  
elements,” Computers & Structures, vol. 82,  
no. 11–12, pp. 945–962, May 2004.  
[10] Y. Lee, P.-S. Lee, and K.-J. Bathe,  
“The MITC3+ shell element and its  
performance,” Computers & Structures, vol.  
138, pp. 12–23, Jul. 2014.  
[11] K.-J. Bathe and E. N. Dvorkin, “A  
four-node plate bending element based on  
Mindlin/Reissner plate theory and a mixed  
interpolation,” Int. J. Numer. Meth. Engng.,  
vol. 21, no. 2, pp. 367–383, Feb. 1985.  
TÀI LIỆU THAM KHẢO  
[1] S. P. Timoshenko and S. Woinowsky-  
Krieger, Theory of Plates and Shells, Second  
Edition. McGraw-Hill, 1959.  
[12] Y. Ko, P.-S. Lee, and K.-J. Bathe, “The  
MITC4+ shell element and its performance,”  
Computers & Structures, vol. 169, pp. 57–68,  
Jun. 2016.  
[2] D. L. Logan, A first course in the Finite  
Element Method, Fourth. Thomson.  
[3] O. C. Zienkiewicz, R. L. Taylor, and J. M.  
Too, “Reduced integration technique in  
general analysis of plates and shells,” Int. J.  
Numer. Meth. Engng., vol. 3, no. 2, pp. 275–  
290, Apr. 1971.  
[13] K.-J. Bathe, F. Brezzi, and S. W. Cho,  
“The MITC7 and MITC9 Plate bending  
elements,” Computers & Structures, vol. 32,  
no. 3, pp. 797–814, Jan. 1989.  
[14] G. R. Liu and T. Nguyen-Thoi,  
Smoothed Finite Element Methods. CRC  
Press, 2010.  
[4] T. J. R. Hughes, M. Cohen, and M.  
Haroun, “Reduced and selective integration  
77  
[15] T. Nguyen-Thoi, P. Phung-Van, H. composite nhiều lớp bằng phần tử MITC3 có  
Luong-Van, H. Nguyen-Van, and H. Nguyen- biến dạng được trung bình trên miền nút phần  
Xuan, “A cell-based smoothed three-node tử (NS-MITC3),” in Tuyển tập công trình Hội  
Mindlin plate element (CS-MIN3) for static nghị  ho  học toàn quốc “Vật liệu và Kết cấu  
and free vibration analyses of plates,” Comput Composite  Cơ học, Công nghệ và Ứng dụng,”  
Mech, vol. 51, no. 1, pp. 65–81, Apr. 2012.  
Nha Trang, Vietnam, 2016, pp. 613–620.  
[16] H. Nguyen-Xuan, G. R. Liu, C. Thai- [23] T. Võ-Ngọc, “Phân tích kết cấu tấm  
Hoang, and T. Nguyen-Thoi, “An edge-based bằng phần tử biến dạng trơn CS-MITC3,”  
smoothed finite element method (ES-FEM) Luận văn Thạc sĩ, Trường ĐH Sư phạm Kỹ  
with stabilized discrete shear gap technique for thuật Tp.HCM, 2017.  
analysis  
of  
Reissner–Mindlin  
plates,”  
[24] T. Chau-Dinh, T. Vo-Ngoc, and P.  
Nguyen-Hoang, “A cell-based smoothed three-  
node plate finite element with a bubble node  
for static analyses of both thin and thick  
Computer Methods in Applied Mechanics and  
Engineering, vol. 199, no. 9–12, pp. 471–489,  
Jan. 2010.  
[17] H. Nguyen-Xuan, T. Rabczuk, N. plates,” Vietnam Journal of Mechanics, vol.  
Nguyen-Thanh, T. Nguyen-Thoi, and S. 39, no. 3, pp. 229–243, Sep. 2017.  
Bordas, “A node-based smoothed finite  
[25] K.-J.  
Procedures. Prentice Hall International, Inc.,  
1996.  
Bathe,  
Finite  
Element  
element method with stabilized discrete shear  
gap technique for analysis of Reissner–  
Mindlin plates,” Comput Mech, vol. 46, no. 5,  
pp. 679–701, Oct. 2010.  
[26] M. Lyly, R. Stenberg, and T. Vihinen,  
“A stable bilinear element for the Reissner-  
Mindlin plate model,” Computer Methods in  
Applied Mechanics and Engineering, vol. 110,  
no. 3–4, pp. 343–357, Dec. 1993.  
[18] T. Nguyen-Thoi, P. Phung-Van, C.  
Thai-Hoang, and H. Nguyen-Xuan, “A cell-  
based smoothed discrete shear gap method  
(CS-DSG3) using triangular elements for static  
and free vibration analyses of shell structures,”  
International Journal of Mechanical Sciences,  
vol. 74, pp. 32–45, Sep. 2013.  
[27] T. Nguyen-Thoi, P. Phung-Van, H.  
Nguyen-Xuan, and C. Thai-Hoang, “A cell-  
based smoothed discrete shear gap method  
using triangular elements for static and free  
vibration analyses of Reissner–Mindlin  
plates,” Int. J. Numer. Meth. Engng, vol. 91,  
no. 7, pp. 705–741, Aug. 2012.  
[19] Q. Nguyễn-Duy, “Phân tích kết cấu  
tấm bằng phần tử biến dạng trơn ES-MITC3,”  
Luận văn Thạc sĩ, Trường ĐH Sư phạm Kỹ  
thuật Tp.HCM, 2016.  
[28] A. Razzaque, “Program for triangular  
bending elements with derivative smoothing,”  
Int. J. Numer. Meth. Engng., vol. 6, no. 3, pp.  
333–343, Jan. 1973.  
[20] T. Chau-Dinh, Q. Nguyen-Duy, and H.  
Nguyen-Xuan, “Improvement on MITC3 plate  
finite element using edge-based strain  
smoothing enhancement for plate analysis,”  
Acta Mech, vol. 228, no. 6, pp. 2141–2163,  
Jun. 2017.  
[29] L. S. D. Morley, Skew Plates and  
Structures. Pergamon Press, Oxford, 1963.  
[30] M. Rezaiee-Pajand and M. Karkon,  
“Hybrid stress and analytical functions for  
analysis of thin plates bending,” Latin  
American Journal of Solids and Structures,  
vol. 11, no. 4, pp. 556–579, Aug. 2014.  
[21] D. Nguyễn-Văn, “Phân tích kết cấu  
tấm bằng phần tử biến dạng trơn NS-MITC3,”  
Luận văn Thạc sĩ, Trường ĐH Sư phạm Kỹ  
thuật Tp.HCM, 2016.  
[22] T. Châu Đình and D. Nguyễn Văn,  
“Phân tích tĩnh và dao động riêng tấm  
Tác giả chịu trách nhiệm bài viết:  
Họ tên: Trần Văn Chơn  
78  
Đơn vị: Khoa Xây dựng, Đại học Sư phạm Kỹ thuật TP.HCM, Việt Nam  
Điện thoại: 0939686791  
79  
BÀI BÁO KHOA HỌC  
THỰC HIỆN CÔNG BỐ THEO QUY CHẾ ĐÀO TẠO THẠC SỸ  
Bài báo khoa học của học viên  
có xác nhận và đề xuất cho đăng của Giảng viên hướng dẫn  
Bản tiếng Việt ©, TRƯỜNG ĐẠI HỌC SƯ PHẠM KỸ THUẬT TP. HỒ CHÍ MINH và TÁC GIẢ  
Bản quyền tác phẩm đã được bảo hộ bởi Luật xuất bản và Luật Sở hữu trí tuệ Việt Nam.  
Nghiêm cấm mọi hình thức xuất bản, sao chụp, phát tán nội dung khi chưa có sự đồng ý  
của tác giả và Trường Đại học Sư phạm Kỹ thuật TP. Hồ Chí Minh.  
ĐỂ CÓ BÀI BÁO KHOA HỌC TỐT, CẦN CHUNG TAY BẢO VỆ TÁC QUYỀN!  
Thực hiện theo MTCL & KHTHMTCL Năm học 2018-2019 của Thư viện Trường Đại học  
Sư phạm Kỹ thuật Tp. Hồ Chí Minh.  
pdf 17 trang Mãnh Khiết 12/01/2024 16300 Free
Bạn đang xem tài liệu "Phân tích kết cấu tấm bằng phần tử mitc3+ được làm trơn trên cạnh (es-mitc3+)", để tải tài liệu gốc về máy hãy click vào nút Download ở trên

File đính kèm:

  • pdfphan_tich_ket_cau_tam_bang_phan_tu_mitc3_duoc_lam_tron_tren.pdf