Phân tích tấm reissner-mindlin bằng phần tử es+ns-mitc3
PHÂN TÍCH TẤM REISSNER-MINDLIN BẰNG PHẦN TỬ
ES+NS-MITC3
Châu Đình Thành a,*, Hồ Thị Còn b, và Lê Phương Bình a
a Khoa Xây Dựng, Trường ĐH Sư Phạm Kỹ Thuật TP.HCM,
Quận Thủ Đức, TP. Hồ Chí Minh
b Trung Tâm Giáo Dục Nghề Nghiệp – Giáo Dục Thường Xuyên,
Huyện An Phú, Tỉnh An Giang
Tóm tắt
Trong nghiên cứu này, phương pháp phần tử hữu hạn biến dạng trơn (S-FEM) dựa trên cạnh
và trên nút phần tử (ES-FEM và NS-FEM) được kết hợp để phát triển cho phần tử tấm tam
giác 3 nút sử dụng kỹ thuật MITC3 để khử hiện tượng khóa cắt. Cách tiếp cận này dựa trên
trên phương pháp -FEM với tham số là hệ số hiệu chỉnh tỉ lệ miền được làm trơn trên
cạnh và trên nút phần tử. Phần tử đề xuất ES+NS-MITC3 có các trường biến dạng uốn và cắt
được làm trơn một phần trên cạnh (ES) và phần còn lại trên nút (NS) phần tử. Phần tử
ES+NS-MITC3 được áp dụng để phân tích một số loại tấm Reissner-Mindlin. Kết quả số cho
thấy, phần tử ES+NS-MITC3 cho kết quả tốt hơn các loại phần tử chỉ sử dụng phương pháp
ES-FEM hoặc NS-FEM, đặc biệt là trong các trường hợp tấm mỏng.
Từ khóa: Tấm Reissner-Mindlin, MITC3, ES-FEM, NS-FEM.
Abstract
In this thesis, the edge and node-based smoothed strain finite element method (ES-FEM and
NS-FEM) were combined to develop for 3-node using the MITC3 technique to attenuate the
shear-locking. This approach is based on the -FEM in which the parameter is used to tune
the contribution ratio of the edge-based and node-based. The ES+NS-MITC3 element in
which the bending and shear strain fields are smoothed on a part of the edge-based domains
and the other on the node-based domains. The ES+NS-MITC3 element phenomenon is
employed some of the similar Reissner-Mindlin plates. Numerical results show that, the
ES+NS-MITC3 element can give results better than those of the similar elements using the
ES-FEM or NS-FEM only, particular in the thin plates.
Keywords: Plates Reissner-Mindlin, MITC3, ES-FEM, NS-FEM.
1
1. Giới thiệu
cạnh và trên nút ES-FEM và NS-FEM đã
được tác giả Chen và LI [18] đề xuất cho
các vấn đề về tấm Reissner – Mindlin, là
phương pháp làm trơn dựa trên phiên bản
Kết cấu tấm được sử dụng phổ biến trong
xây dựng công trình hiện nay. 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ì thế mà lý thuyết tấm Kirchoff [1]
và lý thuyết tấm của Mindlin [2] được sử
dụng phổ biến. Với ý tưởng thiết lập công
thức PTHH [3], cho loại phần tử đơn giản
có thể áp dụng tính toán cho cả tấm dày,
tấm mỏng với độ chính xác cao, việc xây
dựng hàm xấp xỉ chuyển vị C0 đơn giản
hơn so với xây dựng hàm xấp xỉ C1 nhưng
khi xấp xỉ chuyển vị C0 sẽ dẫn đến hiện
tượng khóa cắt “shear locking ”. Đến năm
1985, Bathe và Dvorkin [4 -5] đã tìm ra
phương pháp nội suy hỗn hợp các thành
phần ten-xơ (MITC - Mixed Interpolation
Tensorial Components) và phương pháp
MITC nhanh chóng giải quyết được vấn đề
“khoá cắt”. Sau đó, phương pháp MITC đã
được phát triển rất thành công, cho các loại
phần tử tứ giác 4 nút như MITC4, 8 nút
MITC8 [6], 9 nút MITC9 [7], 16 nút
MITC16 [8],… và MITC4+ [9]. Nối tiếp có
rất nhiều các nghiên cứu, phương pháp
PTHH trơn đã được phát triển cho các phần
tử tấm tam giác 3 nút sử dụng kỹ thuật khử
khóa cắt DGS3 [10], MIN3 [11], MITC3
[12] và MITC3+ [13] như các phần tử, NS-
DSG3 [14], ES-DSG3 [15], ES-MITC3
[16], NS-MITC3 [17],…
(
beta PTHH FEM).
Kỹ thuật làm trơn trên cạnh (ES) và trên
nút (NS) phần tử của phương pháp PTHH
làm trơn đã được kết hợp, để có thể tận
dụng lợi thế về tính chất độ cứng cho bởi cả
hai kỹ thuật làm trơn này. Cách tiếp cận
này còn gọi là beta PTHH (-FEM).
Gần đây, tác giả Li and He [19] phát triển
và tối ưu hóa vị trí của các điểm tích phân
Gauss dựa trên phần tử tam giác trên cạnh
được sử dụng trong ES-FEM kết hợp với
NS-FEM để xây dựng FEM được làm
trơn tổng quát cho các bài toán tấm.
Trong nghiên cứu này, sự kết hợp của làm
trơn trên cạnh (ES) và trên nút (NS) của
phần tử sẽ được phát triển cho phần tử tấm
tam giác 3 nút sử dụng kỹ thuật khử khóa
cắt MITC3, để phân tích kết cấu tấm đồng
nhất đẳng hướng, qua đó đánh giá kết quả
với ý tưởng kết hợp giữa 2 phương pháp
làm trơn khác nhau trên 1 phần tử sẽ cho
kết quả tính toán chính xác hơn và tốc độ
hội tụ nhanh hơn so với chỉ sử dụng
phương pháp làm trơn đơn lẻ trên 1 phần tử
như ES-MITC3 [ 16], NS-MITC3 [ 17] qua
so sánh với các phương pháp hiện tại.
Ngoài ra phương pháp PTHH trơn dựa trên
1
chiều dương qui ước được định nghĩa như
2. Lý thuyết tấm biến dạng cắt
Hình 2.1.
bậc nhất
2.2. Trường biến dạng.
Trong bài báo này, lý thuyết biến dạng cắt
bậc nhất [2] được áp dụng cho tấm đồng
nhất đẳng hướng.
Trường biến dạng được xác định từ trường
chuyển vị (2.1)
2.1. Trường chuyển vị.
u
x
v
y
x
x
x
z
z
Theo giả thuyết tấm Mindlin-Reissner, khi
tấm phẳng chịu tải trọng phân bố 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ị
y
y
y
0
z
(2.2)
y
u v
x
y
xy
z
y x
x
u w
w
x
xz
z x
x
w
y
v w
z y
góc xoay β , β (là 2 góc xoay của đoạn
x
y
yz
y
thẳng vuông góc với mặt trung bình của
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.
tấm lần lượt quay quanh trục y, x) như
Hình 2.1. Cả 3 thành phần này w, β , β
x
y
được biểu diễn như sau:
2.3. Trường ứng suất.
Xét kết cấu tấm được phân tích theo trạng
thái ứng suất phẳng Hình 2.2:
Hình 2.1: Trường chuyển vị trong tấm
u z x, y
x
Hình 2.2: Các thành phần ứng suất phân
bố theo bề dày tấm
v z x, y
y
(2.1)
Theo định luật Hooke mối quan hệ giữa
w w x, y
ứng suất và biến dạng ở Hình 2.2 được
biểu diễn như sau:
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à
x y
góc xoay 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
3
E
x
y
x
x
y
2
z
z
1
x
E
1
v
0
0
x
y
E
y
x
y
2
v
1
1
2
y
(2.3)
0
1 v
1 v
2
Z
0
0
xy
y
E
x
z
xy
xy
y
x
2(1
)
)
E
xz
yz
xz
yz
x
2(1
E
2(1
x
1
v
0
0
Ez
)
y
v
1
2
y
1 v
1 v
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.
0
0
y
2
x
y
x
1 v
0
Ez
1 v2
(2.6)
σ
v 1
0
κ
1 v
Theo giả thuyết, bỏ ứng suất pháp theo
phương z (z = 0). Khi đó trường ứng suất
(2.3) có thể được viết lại dưới dạng ma
trận như sau
0 0
2
w
x
E
x
w
y
xz
yz
2(1 v)
y
(2.7)
Ứng suất do biến dạng uốn:
E
τ
γ
x
1 v
v 1
0
0
x
2(1 v)
E
(2.4)
(2.5)
y
y
Trong đó: = [x y xy]T; κ là vectơ biến
1 v2
xy
1 v xy
0 0
dạng uốn, được xác định bởi
2
x
x
y
xz
E
xz
yz
yz
2(1 v)
(2.8)
κ
Thế quan hệ giữa biến dạng và chuyển vị
cho bởi (2.2) vào ứng suất (2.4) và
y
x
y
y
x
(2.5), ta có:
Ứng suất do biến dạng cắt
= [xz yz]T
và γ là vectơ biến dạng cắt, được xác định
như sau
4
T
w
x
w
y
M
γ x
y
(2.9)
1 v
0
x
(2.11)
M M D v 1 0 κ
y
2.4. Nội lực trong tấm.
1v
M
0 0
xy
2
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 đơn vị này được gọi là nội lực, hay
còn gọi là các thành phần ứng lực của tấm.
Trong đó, M = [Mx My Mxy]T và
D
là độ
cứng trụ của tấm, đặc trưng về hình học và
vật liệu của tấm khi chịu uốn:
Eh3
12(1 v2 )
(2.12)
(2.13)
D
Lực cắt trên 1 đơn vị chiều dài:
h/2
Q
xzdA
yzdA
1dz
1dz
xz
x
A
h/2
h/2
Hình 2.3: Các thành phần nội lực của tấm
Để đưa bài toán ba chiều về bài toán hai
Q
yz
h/2
y
A
chiều, lấy tích phân biểu thức (2.4) theo
chiều dày tấm, ta nhận được các thành
phần nội lực như sau:
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.
h/2
z 1dz
M z dA
x
Thế ứng suất cắt (2.7) vào lực cắt (2.13), ta
x
x
A
h/2
có véc-tơ lực cắt
h/2
z 1dz
(2.10)
M z dA
y
y
y
A
h/2
Q
x
1 0
0 1
Q
kGh
γ
(2.14)
h/2
Qy
z 1dz
M z dA
xy
xy
xy
A
h/2
E
Với
là mô-đun đàn hồi trượt;
G
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 2.3; Mxy là
moment xoắn trên 1 đơn vị dài mặt cắt
quay quanh trục trục x.
2(1v)
k
= 5/6 là hệ số hiệu chỉnh cắt.
3. Công thức phần tử ES+NS-
MITC3
Thế ma trận ứng suất do biến dạng uốn
(2.7) vào trường moment (2.10), ta có:
5
3.1. Công thức PTHH tam giác 3 nút
cho tấm Reissner-Mindlin.
N1 1 ; N2 ; N3
wi là chuyển vị thẳng theo phương
nút thứ i.
z
của
θxi, θyi lần lượt là góc xoay của đoạn thẳng
vuông góc với mặt trung bình quay quanh
trục x và trục y có chiều như Hình 3.1 của
nút thứ i.
Công thức (3.1) biểu diễn dưới dạng ma
trận:
Hình 3.1: Phần tử tam giác 3 nút.
w
x
y
N
0
0
0
w
Xét trên mặt trung bình của tấm, khi ta
thực hiện rời rạc bề mặt này thành các
phần tử tam giác 3 nút để tính toán thì mỗi
nút phần tử sẽ được xác định gồm 3 thành
phần: chuyển vị w và 2 góc xoay (βx, βy).
Mỗi phần tử tam giác tính toán sẽ gồm 3
nút. Mỗi nút sẽ có 3 bậc
i
3
(3.3)
0
0
Ni x
i1
y
Ni
0
i
N
i
Từ công thức (2.8), vectơ biến dạng uốn
được viết lại dưới dạng ma trận:
x
x
y
tự do. Giả thiết với vật liệu đồng nhất và
đẳng hướng có module đàn hồi E và hệ số
poisson .
y
x
y
y
x
(3.4)
Xấp xỉ trường chuyển vị:
i
N
x
0
0
w
3
Ni
y
Ni N
0
0
0
x
y
i1
i
i
x
y
Bbi
Hình 3.2: Hệ tọa độ tổng thể.
3
Trong đó, Bbi là ma trận quan hệ giữa biến
dạng uốn và chuyển vị và được xác định
như sau:
w N w
i
i
i1
3
(3.1)
N
x
i
yi
i1
B B Bb2 Bb3
b
b1
3
N
y
i
xi
0
b c
0
0
0
0
c
0
0
d
c
0
0
0
b
0
0
a
i1
(3.5)
1
0
0
0
d a
2A
Trong đó: Ni là các hàm dạng trong hệ tọa
độ tự nhiên (,).
e
d a b c
d
a
b
7
Với Ae là diện tích của phần tử tam giác có
tọa độ nút (x1,y1), (x2,y2), (x3,y3) như Hình
3.3 và
Với
1
0
0
Eh3
12(1 2
Cb
1
(3.11)
)
a x2 x1; b y2 y1;
(3.6)
1
0 0
c y3 y1; d x3 x1
2
1 0
0 1
Eh
Cs
(3.12)
2(1
)
3.2. Kỹ thuật MITC3.
Ý tưởng cơ bản của kỹ thuật MITC3 dùng
để khử hiện tượng khóa cắt cho phần tử
tấm tam giác 3 nút là xấp xỉ lại các biến
dạng cắt ngoài mặt phẳng sao cho các biến
dạng cắt này là hằng số trên cạnh phần tử.
Lee và Bathe [20] đã xây dựng các hàm
xấp xỉ lại biến dạng cắt ngoài mặt phẳng
thông qua giá trị biến dạng cắt tính từ xấp
xỉ chuyển vị (3.1) tính tại các điểm buộc.
Từ đó, biến dạng cắt ngoài mặt phẳng xấp
xỉ theo kỹ thuật MITC3 có thể được biểu
diễn theo chuyển vị nút phần tử.
Hình 3.3: Tọa độ phần tử tam giác 3 nút
Tương tự, từ công thức (2.9), vectơ biến
dạng cắt được viết lại dưới dạng ma trận:
w
x
w
i
x
y
y
(3.7)
N
x
N
y
w
0
Ni
3
x
i1
i
Ni
0
y
i
B
si
Theo phương pháp PTHH truyền thống, từ
các quan hệ giữa biến dạng và chuyển vị,
hệ phương trình rời rạc của phần tử tấm
tam giác 3 nút dùng phân tích tĩnh tấm
Reissner–Mindlin được viết dưới dạng:
Ku = F (3.8)
Để có thể dễ dàng áp dụng các kỹ thuật
làm trơn cho biến dạng cắt ngoài mặt
phẳng xấp xỉ theo kỹ thuật MITC3, Chau-
Dinh và cộng sự [21] đã tính toán ma trận
quan hệ giữa biến dạng cắt ngoài mặt
phẳng theo kỹ thuật MITC3 và chuyển vị
nút phần tử theo tọa độ nút phần tử bằng
cách sử dụng 1 điểm Gauss.
Trong đó, u là véc-tơ chuyển vị nút kết
cấu, F là véc-tơ tải trọng tại các nút của kết
cấu và được lắp ghép từ các véc-tơ tải
trọng tại nút phần tử
Do đó, công thức PTHH tấm tam giác 3
nút MITC3 có biến dạng cắt (3.7) được
T
fe N 0 0 N 0 0 N 0 0 qdA (3.9)
1
2
3
w
A
e
3
thay bằng MITC3
BMITC3 x (3.13)
ke BTb CbBbdA BTs CsBsdA
si
i1
y
A
A
i
e
e
8
Với:
(b-c)(bc)
6
(d-a)(bc)
(b-c)
A
e
1
2A
MITC3
6
B
B
s1
(b-c)(ad)
(d-a)(ad)
e
(d-a) -A -
-
e
6
6
bc c (bc) ac
d
6
c
(bc)
1
2
2
3
2
MITC3
s2
2A
bd
2
c
6
ad
d
6
e
d
(ad)
(ad)
Hình 3.4: Chia lưới phần tử tam giác làm
trơn trên cạnh và trên nút.
Đối với miền làm trơn kết hợp giữa miền
2
bc b
bd
2
a
6
a
6
b
(bc)
(bc)
1
2A
MITC3
s3
2
6
B
ac
2
b
6
ad
2
e
a
(ad)
(ad)
trơn trên nút và miền làm trơn trên cạnh,
các phần tử tam giác 3 nút trong miền
chuẩn PTHH được chia thành miền bao
quanh một nút để tạo thành miền làm trơn
(3.14)
Và công thức ma trận độ cứng phần tử
(3.10) được viết lại cho phần tử MITC3:
T
MITC3
e
T
b
MITC3
s
MITC3
s
k
ˆ
k
B
C B dA
B
C
B
dA
(3.15)
trên nút NS (miền chấm chấm trên Hình
b
b
s
A
A
e
e
3.5) và phần còn lại chung cạnh giữa 2
Trong đó, BMs ITC3 Bs1
B
B
s3
MITC3
MITC3
s2
MITC3
ES
phần tử là miền làm trơn trên cạnh
3.3. Công thức phần tử ES+NS-MITC3.
(miền gạch chéo trên Hình 3.5).
Xét 1 tấm được chia thành các phần tử tam
giác 3 nút. Theo phương pháp PTHH trơn,
các trường biến dạng được trung bình trên
miền làm trơn giữa các phần tử chung nút
k
kES
hoặc giữa 2 phần tử chung cạnh
NS
như Hình 3.4. Trong đó, miền làm trơn
trên nút được xác định bằng cách nối
Hình 3.5: Miền làm trơn dựa trên lưới tam
giác cho mô hình ES+NS-MITC3
kNS
Để điều chỉnh tỉ lệ giữa 2 miền làm trơn,
trọng tâm với trung điểm các cạnh của các
phần tử chung nút (Hình 3.4a). Miền làm
hệ số hiệu chỉnh được định nghĩa nhằm
xác định tỉ lệ chiều dài trên cạnh phần tử
tham gia vào miền làm trơn trên nút (NS)
và trên cạnh (ES) như sau
kES
trơn trên cạnh
được định nghĩa bằng
cách nối trọng tâm của 2 phần tử chung
cạnh với 2 đỉnh của cạnh chung (Hình
L
AEkS
ANkS
lần lượt là diện tích
L L 1
; L2
L
(3.16)
3.4b). Gọi
và
1
3
2
k
của miền làm trơn trên cạnh
và trên
ES
Trong đó, L là chiều dài cạnh phần tử, L1,
L3 là chiều dài của cạnh phần tử tham gia
vào miền làm trơn trên nút, và L2 là chiều
k
nút
.
NS
9
dài của cạnh phần tử tham gia vào miền
làm trơn trên cạnh.
Và kỹ thuật làm trơn trên nút
e
N
k
1
1 i
A B
e bi
ˆ
B
i1
;
Khi đó, miền làm trơn trên cạnh ES có
bi
2
k
3
1
A
NS
(3.20)
diện tích
e
k
N
1
2
MITC3
1 i MITC3
A B
ˆ
B
1
si
e si
k
AES
AEkS
3
1
A
i1
(3.17)
NS
3
Ai
Trong đó,
là diện tích của phần tử thứ i
k
e
ˆ
Và miền làm trơn trên nút NS có diện
Nke
tham gia vào miền làm trơn và
là số
tích
phần tử trong miền làm trơn k.
1
k
ANS 1 2 Ak
ˆ
(3.18)
NS
Kết quả, ma trận độ cứng của phần tử
3
ES+NS-MITC3 được xác định như sau
Làm trơn các biến dạng uốn và cắt ngoài
keESNS kee
(3.21)
MITC3
mặt phẳng của phần tử tấm tam giác 3 nút
MITC3 trên các miền làm trơn
và
ES
Với:
T
k
MITC3
T
b
MITC3
MITC3
ˆ
k
C
(3.22)
(3.23)
NS bằng cách trung bình biến dạng giữa
e
b
b
s
s
s
s
A
A
e
e
các
T
MITC3
e
T
b
MITC3
s
MITC3
s
ˆ
ˆ
ˆ
ˆ
ˆ
k
B
C B dA
B
B
dA
b
b
A
A
phần tử nằm trong miền làm trơn với diện
tích miền làm trơn, các ma trận quan hệ
giữa biến dạng và chuyển vị nút phần tử
được cho bởi kỹ thuật làm trơn trên cạnh
e
e
Như vậy, từ công thức (3.9) ma trận độ
cứng phần tử ES+NS-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ử.
Thế (3.21) vào phương trình cân bằng
(3.8), sử dụng ngôn ngữ lập trình Matlab
để lập trình tính toán chuyển vị và nội lực
của tấm.
e
N
k
1
1 i
A B
e bi
B
i1
;
bi
2 k
3
A
ES
(3.19)
e
k
N
1
MITC3
si
1 i MITC3
A B
B
2 k
A
ES
e si
3
i1
ES-MITC3 [16], NS-MITC3 [17], ES-
DSG3 [15], QHS [22]. Tham số hiệu chỉnh
tỉ lệ tham gia của miền làm trơn trên nút và
trên cạnh được chọn = 0.6.
4. Ví dụ số
Bài báo thực hiện việc phân tích tĩnh một
số bài toán tấm đồng nhất đẳng hướng, sau
đó so sánh kết quả chuyển vị và mô-men
cho bởi phần tử ES+NS-MITC3 với kết
quả tính toán của các phần tử tấm tam giác
3 nút tham khảo gồm phần tử MITC3 [12],
10
tỉ số h/L = 0.1 hoặc 0.001, và chịu tải trọng
phân bố đều q = 1. Tấm được làm bằng vật
liệu đồng nhất đẳng hướng, có các thông
số vật liệu: mô-đun đàn hồi E = 1092000
và hệ số Poisson = 0.3.
4.1. Bài toán kiểm tra.
1 (0.2; 0.12)
4 (0.44; 0.12)
5 (0.3; 0.6)
Tấm được mô hình thành bằng các phần tử
tam giác 3 nút với lưới chia đều N×N phần
tử như Hình 4.3. Trong đó, N là số phần tử
chia trên mỗi cạnh của tấm. Để khảo sát sự
hội tụ của phần tử đề xuất, các lưới phần
tử khác nhau với N = 4, 8, 10, 12, 16 được
khảo sát. Nhằm mục đích so sánh với lời
giải giải tích và lời giải tham khảo [23],
kết quả độ võng và mô-men không thứ
nguyên tại tâm tấm được tính theo công
thức sau:
2 (0.2; 0)
3 (0.44; 0)
Hình 4.1: Bài toán kiểm tra
Bảng 4.1: Kết quả bài toán kiểm tra
Phần tử
ES+NS-
MITC3
Exact
w5
θx5
θy5
mx5
my5
mxy5
0,6422
0,6422
1,13 – 0,64 – 0,0111 – 0,0111 – 0,0033
1,13 – 0,64 – 0,0111 – 0,0111 – 0,0033
Trong bài toán này, phần tử ES+NS-
MITC3 được kiểm tra khả năng tìm lại kết
quả chuyển vị và nội lực của tấm chữ nhật
có chiều dày h = 0.01 và có kích thước
như Hình 4.1. Tấm làm bằng vật liệu đồng
nhất đẳng hướng có mô-đun đàn hồi E =
105 và hệ số Poisson = 0.3.
4
;
w w / qL /100D
c
c
2
Tấm được rời rạc thành các phần tử tam
giác 3 nút có tọa độ như Hình 4.1. Kết quả
độ võng, góc xoay và mô-men tại nút 5 tọa
độ x = 0.3 và y = 0.6 cho bởi phần tử
ES+NS-MITC3 và lời giải chính xác được
cho trong Bảng 4.1 cho thấy phần tử
NS+ES-MITC3 có khả năng tìm lại được
chính xác giá trị chuyển vị và nội lực của
bài toán kiểm tra.
M M / qL /10
c
c
Trong đó, wc và Mc lần lượt là độ võng và
Eh3
12(1 v2 )
mô-men tại tâm tấm, và
D
là độ cứng trụ của tấm.
4.2. Tấm hình vuông ngàm và tấm tựa
đơn 4 cạnh chịu tải phân bố điều.
Xét ví dụ tấm hình vuông ở 2 trường hợp
ngàm 4 cạnh như Hình 4.2(a) hoặc tựa đơn
4 cạnh như Hình 4.2(b). Cả 2 trường hợp
tấm vuông cạnh chiều dài L, độ dày h với
11
Hình 4.5: Kết quả hội tụ chuẩn hóa độ
võng tại tâm tấm dày vuông ngàm với tỉ lệ
h/L = 0.1
Hình 4.2: (a) Tấm vuông đồng nhất đẳng
hướng biên ngàm 4 cạnh; (b) Tấm vuông
đồng nhất đẳng hướng biên tựa đơn 4
cạnh.
Phần tử ES+NS-MITC3 cho kết quả độ
võng ở tấm mỏng tốt hơn tấm dày. Tuy
nhiên ở trường họp này chỉ vượt trôi hơn
phần tử NS-MITC3 nhưng lại kém hơn
phần tử ES-MITC3, ES–DSG3, MITC3.
Hình 4.3: Chia lưới phần tử hữu hạn cho
tấm.
Kết quả độ võng và mô-men không thứ
nguyên tại tâm tấm ứng với các loại lưới
khác nhau cho bởi phần tử đề xuất và các
phần tử ES-MITC3, NS-MITC3, ES-
DSG3, MITC3
Hình 4.6: Kết quả hội tụ chuẩn hóa
moment tại tâm tấm mỏng vuông ngàm
với tỉ lệ h/L = 0.001
Hình 4.7: Kết quả hội tụ chuẩn hóa
moment tại tâm tấm dày vuông ngàm với tỉ
lệ h/L = 0.1
Hình 4.4: Kết quả hội tụ chuẩn hóa độ
võng tại tâm tấm mỏng vuông ngàm với tỉ
lệ h/L = 0.001
Đối với cả tấm dày và mỏng từ lưới 8×8
đến 16×16 phần tử ES+NS-MITC3 cho kết
quả mô-men tại tâm tấm khá tốt và hội tụ
nhanh đến lời giải chính xác khi chia lưới
12
càng mịn. Ở tấm dày Hình 4.7, phần tử
ES+NS-MITC3 cho kết quả mô-men hội
tụ tương đương phần tử NS-MITC3 và
vượt qua phần tử ES-MITC3, ES-DSG3,
MITC3. Khi tấm mỏng Hình 4.6, phần tử
ES+NS-MITC3 cho kết quả mô men hội tụ
tương đương phần tử ES-MITC3 và vượt
qua phần tử NS-MITC3, ES-DSG3,
MITC3. Tóm lại, phần tử ES+NS-MITC3
cho kết quả mô-men rất tốt khi giải tấm
dày hoặc tấm mỏng.
Hình 4.10: Kết quả hội tụ chuẩn hóa
moment tại tâm tấm dày vuông tựa đơn
với tỉ lệ h /L = 0.001
Khi sử dụng phần tử ES+NS-MITC3 giải
bài toán tấm vuông đồng nhất đẳng hướng
có biên tựa đơn 4 cạnh chịu tải trọng phần
bố đều trong trường hợp tấm mỏng lẫn tấm
dày nhìn chung đều cho ra kết quả độ võng
khá tốt, cho ra kết quả chính xác hơn các
phần tử ES-MITC3, NS-MITC3, ES-
DSG3, và MITC3.
Hình 4.11: Kết quả hội tụ chuẩn hóa
moment tại tâm tấm mỏng vuông tựa đơn
với tỉ lệ h /L = 0.1
Hình 4.8: Kết quả hội tụ chuẩn hóa độ
võng tại tâm tấm mỏng vuông tựa đơn
với tỉ lệ h/L = 0.001
Bên cạnh đó ở Hình 4.10- 4.11 phần tử
ES+NS-MITC3 cho kết quả mô-ment tại
tâm tấm với 2 loại lưới thô và lưới mịn hội
tụ như nhau và gần trùng với lời giải chính
xác và là phần tử cho kết quả tốt nhất so
với các phần tử ES-MITC3, NS-MITC3,
ES-DSG3, và MITC3.
Hình 4.9: Kết quả hội tụ chuẩn hóa độ
võng tại tâm tấm dày vuông tựa đơnvới tỉ
lệ h /L = 0.1
4.3. Tấm hình thoi tựa đơn ở cạnh trên
và dưới, 2 cạnh bên tự do và chịu tải
phân bố đều (tấm Razzaque).
Ở ví dụ này, ta xét tấm hình thoi Razzaque
[24] đồng nhất đẳng hướng nghiên một
góc = 600 so với trục x. Tấm chịu tác
dụng của tải phân bố đều q = 1 như Hình
4.12. Tấm có cạnh kích thước L = 100,
13
chiều dày h = 0.1, tựa đơn cạnh trên và
cạnh dưới, và 2 cạnh bên tự do. Vật liệu
làm tấm có mô-đun đàn hồi E = 10.92 và
hệ số Poisson = 0.3.
phần tử ES+NS-MITC3 cho kết quả tốt
hơn phần tử ES-MITC3, NS-MITC3,
MITC3, và kết quả độ hội tụ gần như trùng
với lời giải Razzaque.
Tấm cũng được mô hình thành các phần tử
tam giác, với các lưới 4×4, 8×8, 12×12,
16×16, và 24×24 như Hình 4.13. Kết quả
4
độ võng chuẩn hóa
w w 100D / qL
c
c
và
mô-men
chuẩn
hóa
Hình 4.15: Kết quả hội tụ chuẩn hóa mô-
men tại tâm tấm hình thoi Razzaque với tỉ
lệ h /L = 0.001
2
tại tâm tấm cho bởi
M M / 10/ qL
y
y
phần tử ES+NS-MITC3 và tài liệu tham
Qua Hình 4.14 và Hình 4.15 ta thấy kết
khảo [25]
quả hội tụ chuẩn hóa độ võng và mô-men
tại tâm tấm hình thoi phần tử ES+NS-
MITC3 cho kết quả chuẩn hóa độ võng và
mô-men gần trùng với lời giải Razzaque.
Phương pháp ES+NS-MITC3 cải thiện
đáng kể độ chính xác của độ võng và
moment so với các phương pháp ES-
MITC3, NS-MITC3, và MITC3.
Hình 4.12: Mô hình tấm hình thoi
Razzaque
4.4. Tấm hình thoi tựa đơn 4 cạnh chịu
tải phân bố đều (tấm Morley).
Hình 4.13: Các cách chia lưới N×N của
tấm hình thoi Razzaque thành các phần tử
tam giác 3 nút (với N = 4, 8, 12, 16, 24).
Tiếp theo, ta xét tấm đồng nhất đẳng
hướng hình thoi Morley tựa đơn chịu tải
phân bố đều như Hình 4.16. Lý thuyết tấm
này được nghiên cứu bởi Morley [26].
Tấm hình thoi nghiêng một góc = 300 so
với trục x, cạnh L = 100, chiều dày h = 0.1,
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 có mô-đun đàn
hồi E = 10.92 và hệ số Poisson = 0.3.
Hình 4.14: Kết quả hội tụ chuẩn hóa độ
võng tại tâm tấm hình thoi Razzaque
với tỉ lệ h/L = 0.001
Từ Hình 4.14 cho thấy giá trị độ võng đạt
được của tấm hình thoi Razzaque bằng
14
Tương tự như ví dụ 4.2, tấm được chia đều
thành các phần tử tam giác 3 nút, với các
loại lưới 4×4, 8×8, 10×10, 12×12, và
16×16.
MITC3. Điều này cho thấy, phần tử
ES+NS-MITC3 có hiệu quả rất tốt khi giải
bài toán tấm Morley.
4.5. Tấm hình tròn ngàm chịu tải phân bố
đều.
Cho tấm tròn bán kính R = 5, chiều dày h
với h/R = 0.2 hoặc 0.02, ngàm trên biên và
chịu tải phân bố đều q = 1 như Hình 4.18.
Tấm được làm bằng vật liệu đồng nhất
đẳng hướng có mô-đun đàn hồi E = 10.92
và hệ số Poisson = 0.3.
Hình 4.16: Mô hình tấm hình thoi Morley
Sử dụng tính chất đối xứng, 1/4 hình tấm
được mô hình bằng các phần tử tam giác 3
nút với các loại lưới như Hình 4.20.
Hình 4.17: Kết quả hội tụ chuẩn hóa độ
võng tại tâm tấm hình thoi Morley tựa đơn
với tỉ lệ h/L = 0.001
Hình 4.19: Tấm tròn với điều kiện biên
ngàm
Hình 4.18: Kết quả hội tụ chuẩn hóa
moment tại tâm tấm hình thoi Morley tựa
đơn với tỉ lệ h/L = 0.001
Từ Hình 4.17 và Hình 4.18 ta thấy độ võng
và moment tấm hình thoi Morley tựa đơn
sử dụng phần tử ES+NS-MITC3 cho kết
quả độ võng hội tụ nhanh hơn phần tử ES-
DSG3, NS-MITC3, MITC3. Tuy nhiên,
kết quả của phần tử đề xuất kém hơn phần
tử ES-MITC3. Mô-men chuẩn hóa tại tâm
tấm của phần tử ES+NS-MITC3 lại hội tụ
nhanh hơn giá trị mô-men cho bởi phần tử
ES-MITC3, MITC3, ES-DSG3, và NS-
Hình 4.20: Mô hình chia lưới 1/4 tấm tròn
có sử dụng phần tử tam giác 3 nút
Kết quả khảo sát độ võng và mô-men tại
tâm tấm khi tỉ số h/R thay đổi với các loại
lưới khác nhau, kết quả so sánh độ võng và
mô-men cho bởi phần tử ES+NS-MITC3
và các kết quả tham khảo [27] được thể
hiện qua Hình 4.21-4.24.
15
Hình 4.23: Kết quả hội tụ mô-men
của tấm tròn ngàm (h/R = 0.02),
chịu tải phân bố đều.
Hình 4.21: Kết quả hội tụ của độ võng của
tấm tròn ngàm (h/R = 0.02), chịu tải phân
bố đều.
Hình 4.24: Kết quả hội tụ mô-men của tấm
tròn ngàm (h/R = 0.2),
chịu tải phân bố đều.
Từ Hình 4.21-4.24 cho thấy tấm tròn đồng
nhất đẳng hướng có biên ngàm giải bằng
phần tử ES+NS-MITC3 cho ra kết quả rất
tốt đối với tấm mỏng và tấm dày. Trong
trường hợp lưới thô hoặc lưới mịn, độ
võng và mô-men đều hội tụ gần với lời
giải chính xác và hội tụ nhanh hơn các
phần tử ES-DSG3, MITC3, NS-MITC3,
và ES-MITC3.
Hình 4.22: Kết quả hội tụ của độ võng của
tấm tròn ngàm (h/R = 0.2),
chịu tải phân bố đều.
5. Kết luận
Từ những kết quả tính toán độ võng và
mô-men các bài toán tấm đồng nhất đẳng
hướng với hình dạng và điều kiện biên
khác nhau và so sánh với các phần tử tấm
tam giác 3 nút cùng loại, có thể rút ra một
số kết luận như sau: Phần tử ES+NS-
MITC3 cho kết quả tương đối tốt, khá gần
với lời giải chính xác đối với một số loại
tấm khác nhau như tấm hình vuông, tấm
hình thoi Razzaque, Morley và tấm tròn có
điều kiện biên tựa đơn hoặc ngàm. Đặc
biệt, kết quả độ võng và mô-men cho bởi
phần tử ES+NS-MITC3 rất tốt đối với tấm
tròn, tấm hình thoi Morley, Razzaque, và
tấm vuông tựa đơn chịu tải phân bố đều.
Ngoại trừ trường hợp tấm vuông ngàm cho
kết quả hội tụ chưa tốt nhưng so với phần
tử NS-MITC3, ES-DSG3, MITC3, QHS
thì phần tử ES+NS-MITC3 cải thiện đáng
kể độ hội tụ và chính xác về độ võng và
16
mô-men tại điểm khảo sát. Tuy nhiên, so
với các loại phần tử chỉ sử dụng 1 kỹ thuật
làm trơn thì phần tử ES+NS-MITC3 tốn
nhiều thời gian hơn để tính toán ma trận độ
cứng phần tử do phải xác định miền làm
trơn trên cạnh (ES) và trên nút (NS) phần
tử.
Tài liệu tham khảo
[1] K. T. S. Oldham, The doctrine of
description: Gustav Kirchhoff, classical
physics, and the “purpose of all science”
in 19th-century Germany. University of
California, Berkeley, 2008.
[6]
K.-J. Bathe and E. N. Dvorkin, “A
formulation of general shell elements—
the use of mixed interpolation of
tensorial components,” International
Journal for Numerical Methods in
Engineering, vol. 22, no. 3, pp. 697–722,
Mar. 1986.
[2] R. D. Mindlin, ‘Influence of Rotary
Inertia and Shear on Flexural Motions of
Isotropic, Elastic Plates,’ Journal of
Applied Mechanics, Vol. 18, No. 1, 1951,
[7] K.-J. Bathe, P.-S. Lee, and J.-F. Hiller,
“Towards improving the MITC9
pp.
31-
shell element,” Computers
&
38. - References - Scientific Research
Publish.”
Structures, vol. 81, no. 8, pp. 477–
489, May 2003.
[3] G. R. Liu and N. T. Trung, Smoothed
Finite Element Methods, 1 edition.
Boca Raton: CRC Press, 2010.
[8] M. L. Bucalem and K.-J. Bathe,
“Higher-order MITC general shell
elements,” International Journal for
Numerical Methods in Engineering,
vol. 36, no. 21, pp. 3729–3754, Nov.
1993.
[4] E. Dvorkin and K.-J. Bathe, “A
Continuum mechanics based four-
node shell element for general
nonlinear analysis,” Engineering
Computations, vol. 1, pp. 77–88,
Dec. 1984.
[9] Y. Ko, P.-S. Lee, and K.-J. Bathe, “A
new MITC4+ shell element,”
Computers & Structures, vol. 182,
pp. 404–418, Apr. 2017.
[5] K.-J. Bathe and E. N. Dvorkin, “A
four-node plate bending element
based on Mindlin/Reissner plate
theory and a mixed interpolation,”
International Journal for Numerical
Methods in Engineering, vol. 21, no.
2, pp. 367–383. 1984.
[10] 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.
17
[11] 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.
[16] Q. Nguyen-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.
[17] D. Nguyen-Van ''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.
[12] H.-M. Jeon, P.-S. Lee, and K.-J.
Bathe, “The MITC3 shell finite
element enriched by interpolation
covers,” Computers & Structures,
vol. 134, pp. 128–142, Apr. 2014.
[18] F. Wu, W. Zeng, L. Y. Yao, M. Hu,
Y. J. Chen, and M. S. Li, “Smoothing
Technique Based Beta FEM (βFEM)
for Static and Free Vibration
[13] Y. Lee, P.-S. Lee, and K.-J. Bathe,
"The MITC3+ shellelement and its
performance
Computers
&
Analyses
of
Reissner–Mindlin
Structures'', vol. 138, pp. 12-23, Jul.
2014.
Plates,” Int. J. Comput. Methods, p.
1845006, Mar. 2018.
[14] Nguyen – Thoi T, Liu GR, Nguyen –
Xuan H. Additional properties of the
node – based smoothed finite element
method (NS – FEM) for solid
mechanics problems. International
Journal of Computational Methods,
6(4): 633 – 666, 2009.
[19] E. Li and Z. He ''Development of a
perfect match system in the
improvement of eigenfrequencies of
free
vibration,”
Applied
Mathematical Modelling, vol. 44, pp.
614–639, Apr. 2017.
[20] P.-S. Lee and K.-J. Bathe,
“Development of MITC isotropic
triangular shell finite elements,”
Computers & Structures, vol. 82, no.
11, pp. 945–962, May 2004.
[15] P. Phung-Van, Chien H. Thai, T.
Nguyen-Thoi, H. Nguyen-Xuan,
''Static and free vibration analyses of
composite and sandwich plates by an
edge-based smoothed discrete shear
gap method (ES-DSG3) using
[21] T. Chau-Dinh, Q. Nguyen-Duy, and
H. Nguyen-Xuan, “Improvement on
MITC3 plate finite element using
triangular elements
based on
layerwise theory'', Composites Part
B: Engineering, Volume 60, Pages
227–238, April 2014.
edge-based
strain
smoothing
enhancement for plate analysis,” Acta
Mechanica, vol. 228, Feb. 2017.
2
[22] Chen Wanji and Y. K. Cheung.
Refined 9-Dof triangular Mindlin
plate elements. International journal
for Static and Free Vibration
Analyses of Reissner–Mindlin
Plates,'' Int. J. Comput. Methods, p.
1845006, Mar. 2018.
for
engineering; 51:1259–1281. 2001.
[23] S. P. Timoshenko and
Woinowsky-Krieger,
numerical
methods
in
[26] L. S. D. Morley, Skew plates and
structures. Oxford; New York:
Pergamon Press; distributed in the
Western Hemisphere by Macmillan,
New York, 1963.
S.
Theory of
Plates and Shells, Second Edition.
McGraw-Hill, 1959.
[24] A.
Razzaque,
“Program
for
[27] R. Ayad, A. Rigolot, and N. Talbi,
“An improved three-node hybrid-
mixed element for Mindlin/Reissner
plates,” International Journal for
Numerical Methods in Engineering,
vol. 51, no. 8, pp. 919–942, 2001.
triangular bending elements with
derivative smoothing,” International
Journal for Numerical Methods in
Engineering, vol. 6, no. 3, pp. 333–
343, Jan. 1973.
[25] J. Chen, and M. S. Li, ''Smoothing
Technique Based Beta FEM (βFEM)
Tác giả chịu trách nhiệm bài viết:
Họ tên: Hồ Thị Còn
Đơn vị: Trung Tâm Giáo Dục Nghề Nghiệp – Giáo Dục Thường Xuyên,
Huyện An Phú, Tỉnh An Giang
Điện thoại: 0939550094
Email: [email protected]
3
Bạn đang xem tài liệu "Phân tích tấm reissner-mindlin bằng phần tử es+ns-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:
phan_tich_tam_reissner_mindlin_bang_phan_tu_esns_mitc3.pdf

