Ứng dụng phân tích đẳng hình học trong bài toán biến dạng phẳng của lý thuyết đàn hồi

THE INTERNATIONAL CONFERENCE ON MARINE SCIENCE AND TECHNOLOGY 2016 
HỘI NGHỊ QUỐC TẾ KHOA HỌC CÔNG NGHỆ HÀNG HẢI 2016 271 
Ứng dụng phân tích đẳng hình học trong bài toán biến dạng phẳng 
của lý thuyết đàn hồi 
Application of isogeometric analysis in plane strain 
Phạm Quốc Hoàn 
Trường Đại học Hàng hải Việt Nam, 
hoanpq.ctt@vimaru.edu.vn 
Tóm tắt 
Bài báo này trình bày tổng quan về phương pháp phân tích đẳng hình học và ứng dụng của 
phương pháp này trong tính toán bài toán phẳng của lý thuyết đàn hồi. 
Từ khóa: Phân tích đẳng hình học, Nurb, hàm dạng, B - Spline, biến dạng phẳng, bài toán phẳng. 
Abstract 
In this paper, we present an overview of Isogeometric analysis method and its application in 
plane strain of two-dimensional elasticity. 
Keywords: Isogeometric, Nurb, basic function, B - Spline, plane strain, two-dimensional elasticity. 
1. Đặt vấn đề 
Phương pháp số hay còn gọi là giải tích số là môn khoa học chuyên nghiên cứu cách giải 
gần đúng, đa phần là các phương trình, các bài toán xấp xỉ hàm số và các bài toán tối ưu. Ngày này, 
các phương pháp tính số đã và đang phát triển rất mạnh mẽ và trở thành công cụ để giải toán các 
bài toán khoa học kỹ thuật như phương pháp sai phân hữu hạn, phương pháp phần tử hữu hạn, 
phương pháp phần tử biên, phương pháp phân tích đẳng hình học,... Trong đó phương pháp phân 
tích đẳng hình học là một phương pháp mới được phát triển với nhiều ưu điểm và có khả năng ứng 
dụng cao trong các bài toán liên quan đến kết cấu. 
2. Tổng quan về IGA và hình học Nurb[2, 3] 
Phương pháp phân tích đẳng hình học “Isogeometric analysis” hay được viết tắt là IGA 
được giáo sư T.Huges và các cộng sự đưa ra vào năm 2005 dựa trên ý tưởng sử dụng các hàm mô 
phỏng trong các phần mềm CAD cho việc chia lưới và chia phần tử tính toán trong phương pháp 
phần tử hữu hạn. Khác với phần tử hữu hạn sau khi xây dựng được mô hình trên các phần mềm 
CAD cần có bước chia lại lưới sau đó chuyển qua quá trình phân tích và tối ưu hóa kết cấu, IGA sử 
dụng luôn mô hình kết cấu từ CAD trong bước phân tích và tối ưu hóa kết cấu. IGA sử dụng trực 
tiếp mô hình từ CAD, thay vì sử dụng đa thức Lagrange cho việc xấp xỉ hình học và chuyển vị sẽ 
sử dụng B - Spline và Nurbs (Non uniform rational B-Spline). 
2.1. Vector nút (knot vector) và hàm cơ sở 
B-Spline là hàm cơ sở dùng để biểu diễn đường cong từ năm 1972 còn Nurbs là dạng tổng 
quát hóa của đường cong B - Spline với khả năng biểu diễn chính xác các đường “conic” (đường 
tròn, đường elip). Nurbs bắt đầu được sử dụng trong thiết kế kỹ thuật từ năm 1972. Ban đầu Nurbs 
chỉ được ưu tiên trong lĩnh vực xe hơi, hàng không, ngày nay Nurbss đã có mặt trong tất cả các gói 
CAD chuẩn. 
Véctơ Knot (hay còn gọi là véctơ nút) 1 2 1, ,..., n p   + += là một chuỗi các giá trị tham số 
không giảm, I < i+1, i = 1 ÷ n + p + 1; trong đó i R là knot (nút) thứ i, i được gọi là chỉ số nút, 
p là bậc của hàm cơ sở còn n là số hàm cơ sở được sử dụng để xây dựng đường cong B - Spline. 
Trong trường hợp tổng quát, một đường cong B - Spline sẽ không đi qua hai điểm điều 
khiển đầu và cuối, nó chỉ đi qua nếu điểm knot đầu và cuối có bội là p+1 và một vector knot như 
vậy được gọi là vector knot mở. Hiện nay, thông thường khi thiết kế một đường cong, chỉ yêu cầu 
định rõ điểm đầu và điểm cuối nên vector knot trong các chương trình CAD là mở. 
Nếu cho một véctơ Knot thì các hàm cơ sở được định nghĩa đệ quy bắt đầu với p = 0: 
THE INTERNATIONAL CONFERENCE ON MARINE SCIENCE AND TECHNOLOGY 2016 
HỘI NGHỊ QUỐC TẾ KHOA HỌC CÔNG NGHỆ HÀNG HẢI 2016 272 
1 1
( ),0
0 
khi i i
Ni
  

 +
= 
 (1) 
 Với p > 1 
1
( ) ( ) ( ), , 1 1, 1
1 1
i piN N Ni p i p i p
i p i p i
  
  
   
-- + +
= +- + -
- -+ + + +
 (2) 
2.2. Đường cong và bề mặt B - Spline 
Đường cong B - Spline được xác định bằng tổ hợp tuyến tính của các điểm điều khiển và 
các hàm cơ sở tương ứng 
( ) ( ),
1
n
Ni p i
i
 = 
=
C P (3) 
Pi: điểm điều khiển (control point) thứ i; 
, ( )i pN  : là hàm cơ sở B - Spline thứ i có bậc là p. 
Bề mặt B - Spline được định dựng bằng tích ten xơ của các hàm cơ sở B - Spline một chiều, 
hai vector Knot và một lưới các điểm điều khiển Pi,j hai chiều n×m: 
, , ,
1
( , ) ( ) ( )   
=
= 
n
i p J p i J
i
S N M P (4) 
M và N là các hàm cơ sở B - Spline theo phương  và  tương ứng. 
Bề mặt B - Spline còn được định nghĩa như sau: 
( , ) ( , )   = 
nxm
A A
A
S N P (5) 
NA(,): là hàm dạng liên quan đến nút A. 
2.3. Hình học Nurb 
Hình học Nurbs là sự tổng quát hóa của hình học B - Spline. Tương tự như đường S - pline, 
đường Nurbs bậc p được định nghĩa bởi: 
R = 
=
( ) ( ),
1
n
i p i
i
C P (6) 
N
R
N
 

 
=

=
( ),
( )
,
( ),
1
i p i
i p n
i p i
i
i: trọng số; 
Ri,p: hàm hữu tỷ một chiều. 
Trọng số là các đại lượng vô hướng lớn hơn không, các trọng số không nhất thiết phải bằng 
nhau, nếu tất cả các trọng số bằng nhau thì Nurbs trở thành B - Spline. 
Còn đối với các mặt Nurbs, hàm cơ sở hữu tỉ hai chiều được biểu diễn như sau: 
N
R
N M
  
 
  
=
 
= =
( )M ( ), , ,,
( , )
,j
( ) ( ), j,q ,
1 1
i p j q i jp q
i n m
i p i j
i j (7) 
2.4. Các phương pháp làm mịn lưới hình học 
Ta có 3 phương pháp làm mịn lưới sau: 
THE INTERNATIONAL CONFERENCE ON MARINE SCIENCE AND TECHNOLOGY 2016 
HỘI NGHỊ QUỐC TẾ KHOA HỌC CÔNG NGHỆ HÀNG HẢI 2016 273 
- Phương pháp làm mịn lưới h: được thực hiện thông qua kỹ thuật chèn knot mà vẫn giữ 
nguyên bậc của hàm cơ sở. Khi thực hiện theo phương pháp làm mịn h vector knot được biến đổi 
tương ứng. 
Ví dụ chèn thêm 1 nút vào 0,0,0,1,1,1=
ta được vecter knot mới 0,0,0,0.51,1,1= , khi đó 
các hàm cơ sở cũng thay đổi tương ứng: 
Hình 1. Hàm cơ sở được làm mịn sau khi chèn thêm knot 
- Phương pháp làm mịn lưới p: được thực hiện thông qua kỹ thuật tăng bậc của các hàm cơ 
sở dùng trong biểu diễn hình học, cứ mỗi quá trình tăng bậc thì bội của knot được tăng lên trong khi 
không thêm bất kỳ giá trị knot nào. 
Hình 2. Hàm cơ sở được làm mịn sau tăng bậc 
- Phương pháp làm mịn lưới k: khác với phương pháp làm mịn h và và p, phương pháp làm 
mịn k là phương pháp mới so với phương pháp phần tử hữu hạn. Phương pháp làm mịn k tăng độ 
mịn của lưới bằng cách vừa tăng bậc của hàm cơ sở vừa chèn thêm nút để đảm bảo tính liên tục của 
hàm cơ sở không giảm. 
3. Giải bài toán phẳng của lý thuyết đàn hồi theo phương pháp phân tích đẳng hình học 
3.1. Dạng yếu và dạng mạnh [4] 
Dạng mạnh của bài toán phẳng được phát biểu như sau: cho g và các điều kiện biên 
“Drichlet” hay còn gọi là điều kiện biên chính hoặc điều kiện biên hình học: u = u trên biên và 
điều kiện biên “Neumann” hay còn gọi là điều kiện biên tự nhiên hoặc điều kiện biên lực trên biên 
q = q yêu cầu tìm u sao cho: 
   . 0 + =g (8) 
D us=  (9) 
Dạng yếu được viết như sau: cho g và các điều kiện biên “Drichlet” hay còn gọi là điều kiện 
biên chính hoặc điều kiện biên hình học: u = u trên biên và điều kiện biên “Neumann” hay còn gọi 
THE INTERNATIONAL CONFERENCE ON MARINE SCIENCE AND TECHNOLOGY 2016 
HỘI NGHỊ QUỐC TẾ KHOA HỌC CÔNG NGHỆ HÀNG HẢI 2016 274 
là điều kiện biên tự nhiên hoặc điều kiện biên lực trên biên q = q yêu cầu tìm u sao cho với mọi 
hàm thử (hàm trọng số) u: 
( ) T TTu D ud u qd u gds s t
    = +     (10) 
Hình 3. Phương pháp làm mịn k 
3.2. Phương trình ma trận, ma trận độ cứng và vector lực cục bộ 
Phương trình ma trận trong phương pháp phân tích đẳng hình học và phương pháp phần tử 
hữu hạn cũng có dạng như nhau: 
Kd = F (11) 
Trong đó: 
K: ma trận độ cứng; 
d: vector chuyển vị; 
F: vector lực nút; 
Ma trận độ cứng và vector lực cục bộ như sau: 
( ( )) ( )
( ) ( )
e e Tk k B x DB x da bab e
e ef f N x gd N x qda a ae t
 = =   
 = =  +    
 (12,13) 
Trong đó: 
0
0
Na
x
NaBa
y
N Na a
y x
 
 
 
= 
 
  
  
 (14) 
Để đưa vào lập trình: 
e Tk B DBd
e
=   (15) 
THE INTERNATIONAL CONFERENCE ON MARINE SCIENCE AND TECHNOLOGY 2016 
HỘI NGHỊ QUỐC TẾ KHOA HỌC CÔNG NGHỆ HÀNG HẢI 2016 275 
1 2 ... 0 0 ... 0
0 0 ... ...
1 2 1 2... ...
N N Nnen
x x x
B
N N N N N Nnen nen
y y y x x x
   
   
= 
      
       
 (16) 
e T Tef N gd N gd
e t
= +    (17) 
Trong đó N là ma trận hàm dạng được định nghĩa cho bài toán hai chiều như sau: 
0... 0 0 01 2
0 0 0 0 ...1 2
N N Nnen
N
NN N nen
= 
 (18) 
nen: số lượng hàm dạng của mỗi phần tử. 
3.3. Áp dụng phương pháp đẳng hình học trong bài toán ống dày chịu áp lực phân bố đều 
(bài toán Lame) 
Xét ứng suất trong ống dày có đường kính ngoài b và đường kính trong a có mặt trong và 
mặt ngoài chịu áp lực phân bố đều phía trong và bên n goài là Pa và Pb. 
Giá trị ứng suất trong ống được tính toán theo giải tích như sau [1]: 
 (19) 
 (20) 
Ứng dụng tính toán cho đường ống dẫn chất lỏng thành dày chịu áp lực bên trong về bên 
ngoài thành ống cụ thể trong trường hợp đồng ống của dự án khí điện đạm Cà Mau như sau: 
Bán kính trong a = 21,11 cm; 
Bán kính ngoài b = 22,86 cm; 
Áp lực bên ngoài thành ống Pb = 50 T/m2 = 500 kN/m2 = 0.5 kN/cm2; 
Áp lực trong thành ống Pa = 14,76 MPa = 1,476 kN/cm2; 
Mô đun đàn hồi E = 20700 kN/cm2, hệ số nở hông là 0,3. 
Đường ống có chiều dài rất lớn so với kích thước của tiết diện nên có thể tính toán theo lý 
thuyết bài toán biến dạng phẳng của lý thuyết đàn hồi, do tính chất đối xứng ta chỉ cần tính toán cho 
¼ mô hình. Sử dụng IGA với số lượng phần tử theo 2 phương là 2x2. 
Hình 4. Mô hình hóa 
theo phương pháp IGA 
Hình 5. Mô hình hóa bằng Sap2000 
2 2 2 2
2 2 2 2 2
( )
( )
a b a b
r
P a P b b a P P
b a r b a

- -
= +
- -
2 2 2 2
2 2 2 2 2
( )
( )
a b a bP a P b b a P P
b a r b a

- -
= -
- -
THE INTERNATIONAL CONFERENCE ON MARINE SCIENCE AND TECHNOLOGY 2016 
HỘI NGHỊ QUỐC TẾ KHOA HỌC CÔNG NGHỆ HÀNG HẢI 2016 276 
Để so sánh với phương pháp phần tử hữu hạn ta tiến hành tính toán cho ¼ mô hình bằng phần 
mềm Sap 2000 số lượng phần tử là 4x2 đồng thời chia lưới ảo 2x2. 
Kết quả tính toán tại bán kính trung bình so sánh với lý thuyết được cho ở bảng 1. 
Kết quả do Sap 2000 sai khác so với lý thuyết khá lớn do số lượng phần tử nhỏ, đây là hạn 
chế của phương pháp phần tử hữu hạn khi tính toán cho trường hợp kết cấu có hình dạng bậc cao. 
4. Kết luận 
Phương pháp phân tích đẳng hình học trong các bài toán phẳng của lý thuyết đàn hồi cho kết 
quả tương đối chính xác so với lý thuyết đàn hồi trong khi số lượng phần tử nhỏ hơn rất nhiều so 
với tính toán bằng phương pháp phần tử hữu hạn. So với phương pháp phần tử hữu hạn do việc sử 
dụng hàm cơ sở B - Spline và Nurb đặc biệt có hiệu quả khi tính toán cho các vật thể có hình dạng 
bậc cao. 
Bảng 1. Kết quả đo chiều cao sóng 
Góc mở 
LTĐH IGA Sap2000 
xx yy xx yy xx yy 
0 -0.95 11.26 -0.95 11.26 0.14 18.85 
 /8 0.83 9.47 0.84 9.48 1.91 14.55 
 /4 5.15 5.15 5.16 5.15 8.31 8.31 
3 /8 9.47 0.83 9.48 0.83 14.55 1.91 
 /2 11.26 -0.95 11.26 -0.95 18.85 0.14 
Tài liệu tham khảo 
[1]. Nguyễn Văn Ngọc. Bài giảng đại học Cơ học trong môi trường liên tục. Đại học Hàng hải 
Việt Nam. 2014 
[2]. Đỗ Văn Hiến, Châu Nguyên Khánh, Nguyễn Xuân Hùng. Isogeometric analysis of plane - 
curved beams. presented at the Hội nghị Cơ học kỹ thuật toàn quốc. Đà Nẵng. 2015. 
[3]. Nguyễn Xuân Hùng. Isogeometric analysis, from theory to application. VGU. 2014. 
[4]. Nguyễn Xuân Hùng. Phân tích đẳng hình học cầu nối hợp nhất giữa mô hình mô phỏng và 
thiết kế. NXB Đại học Quốc gia TP. Hồ Chí Minh. 2015.