BỘ GIÁO DỤC VÀ ĐÀO TẠO
TRƯỜNG ĐẠI HỌC CÔNG NGHỆ TP. HCM
---------------------------
DƯƠNG GIẢ DỰNG
TỐI ƯU HÓA ĐA MỤC TIÊU KẾT CẤU DÀN
DỰA TRÊN ĐỘ TIN CẬY
SỬ DỤNG PHƯƠNG PHÁP NSGA-II
LUẬN VĂN THẠC SĨ
Chuyên ngành: Kỹ thuật xây dựng công trình dân dụng và công nghiệp
Mã số ngành: 60580208
TP. HỒ CHÍ MINH, tháng 10 năm 2017
BỘ GIÁO DỤC VÀ ĐÀO TẠO
TRƯỜNG ĐẠI HỌC CÔNG NGHỆ TP. HCM
---------------------------
DƯƠNG GIẢ DỰNG
TỐI ƯU HÓA ĐA MỤC TIÊU KẾT CẤU DÀN
DỰA TRÊN ĐỘ TIN CẬY
SỬ DỤNG PHƯƠNG PHÁP NSGA-II
LUẬN VĂN THẠC SĨ
Chuyên ngành: Kỹ thuật xây dựng công trình dân dụng và công nghiệp
Mã số ngành: 60580208
CÁN BỘ HƯỚNG DẪN KHOA HỌC: PGS. TS. NGUYỄN THỜI TRUNG
TP. HỒ CHÍ MINH, tháng 10 năm 2017
CÔNG TRÌNH ĐƯỢC HOÀN THÀNH TẠI
TRƯỜNG ĐẠI HỌC CÔNG NGHỆ TP. HCM
Cán bộ hướng dẫn khoa học: PGS. TS. NGUYỄN THỜI TRUNG
Luận văn Thạc sĩ được bảo vệ tại Trường Đại học Công nghệ TP. HCM
ngày 04 tháng 10 năm 2017.
Thành phần Hội đồng đánh giá Luận văn Thạc sĩ gồm:
TT Họ và tên Chức danh Hội đồng
1 TS. Khổng Trọng Toàn Chủ tịch
2 TS. Trần Tuấn Nam Phản biện 1
3 TS. Đào Đình Nhân Phản biện 2
4 TS. Nguyễn Văn Giang Ủy viên
5 TS. Nguyễn Hồng Ân Ủy viên, Thư ký
Xác nhận của Chủ tịch Hội đồng đánh giá Luận sau khi Luận văn đã được
sửa chữa.
Chủ tịch Hội đồng đánh giá luận văn
TRƯỜNG ĐH CÔNG NGHỆ TP. HCM CỘNG HÒA XÃ HỘI CHỦ NGHĨA VIỆT NAM
PHÒNG QLKH – ĐTSĐH Độc lập – Tự do – Hạnh phúc
TP. HCM, ngày 10 tháng 10 năm 2017
NHIỆM VỤ LUẬN VĂN THẠC SĨ
Họ tên học viên: DƯƠNG GIẢ DỰNG Giới tính: Nam
Ngày, tháng, năm sinh: 22/03/1989 Nơi sinh: Bình Định
Chuyên ngành: Kỹ thuật xây dựng công trình dân dụng và công nghiệp
MSHV: 1541870003
I- Tên đề tài:
Tối ưu hóa đa mục tiêu kết cấu dàn dựa trên độ tin cậy sử dụng phương pháp NSGA-II.
II- Nhiệm vụ và nội dung:
- Thiết lập bài toán tối ưu hoá đa mục tiêu tiền định.
- Thiết lập bài toán tối ưu hoá đa mục tiêu dựa trên độ tin cậy.
- Nghiên cứu giải thuật NSGA-II.
- Nghiên cứu bài toán phân tích độ tin cậy cho kết cấu và giải thuật tối ưu hóa dựa trên
độ tin cậy một vòng lặp đơn xác định.
- Áp dụng thuật toán NSGA-II để tìm lời giải tối ưu bài toán kết cấu dàn 2D và dàn
3D đã được thiết lập:
III- Ngày giao nhiệm vụ: 26/09/2016
IV- Ngày hoàn thành nhiệm vụ: 31/03/2017
V- Cán bộ hướng dẫn: PGS. TS. NGUYỄN THỜI TRUNG
CÁN BỘ HƯỚNG DẪN KHOA QUẢN LÝ CHUYÊN NGÀNH
i
LỜI CAM ĐOAN
Tôi xin cam đoan Luận văn “Tối ưu hóa đa mục tiêu kết cấu dàn dựa trên độ tin
cậy sử dụng phương pháp NSGA-II” dưới sự hướng dẫn của PGS. TS. Nguyễn Thời
Trung là công trình nghiên cứu của riêng tôi. Các số liệu, kết quả nêu trong Luận văn
là trung thực và chưa từng được ai công bố trong bất kỳ công trình nào khác. Tôi xin
cam đoan rằng mọi sự giúp đỡ cho việc thực hiện Luận văn này đã được cảm ơn và
các thông tin trích dẫn trong Luận văn đã được chỉ rõ nguồn gốc.
Học viên thực hiện Luận văn
Dương Giả Dựng
ii
LỜI CẢM ƠN
Đầu tiên tôi xin gửi lời cảm ơn chân thành đến thầy PGS. TS. Nguyễn Thời
Trung, thầy là người đã giảng dạy cho tôi môn “Phương pháp phần tử hữu hạn” và
môn “Tối ưu hóa kết cấu dựa trên độ tin cậy” trong quá trình học cao học. Thầy cũng
là người hướng dẫn tôi hoàn thành luận văn này với một lòng nhiệt huyết, kiến thức
chuyên sâu và lời giảng dạy xúc tích chuyên nghiệp. Thầy cũng là người tạo động lực
và truyền cảm hứng cho tôi để tôi quyết định thực hiện luận văn chuyên ngành kết
cấu, cụ thể là lĩnh vực đánh giá độ tin cậy kết cấu.
Kế đến tôi xin gửi lời cám ơn đến các thầy cô đã giảng dạy tôi trong quá trình
học cao học tại Trường đại học Công nghệ TP.HCM. Những kiến thức giá trị của các
thầy cô đã giúp tôi rất nhiều cho việc hoàn thành luận văn này.
Tôi cũng xin cảm ơn ThS. Hồ Hữu Vịnh, ThS. Võ Duy Trung là người đã
hướng dẫn cho tôi sử dụng ngôn ngữ lập trình Matlab, đã hỗ trợ tận tình cho tôi trong
suốt quá trình làm luận văn. Tôi cũng xin cám ơn các anh chị em đang làm việc tại
Viện Khoa học Tính toán của Trường đại học Tôn Đức Thắng, là một môi trường cởi
mở, thân thiện. Nhờ sự giúp đỡ nhiệt tình của mọi người, tôi đã hoàn thành tốt luận
văn này.
Cuối cùng, tôi xin gửi lời cảm ơn đến cơ quan, đồng nghiệp và gia đình đã luôn
động viên, hỗ trợ và tạo điều kiện thuận lợi nhất cho tôi trong suốt khóa học vừa qua.
TP. HCM, ngày 10 tháng 10 năm 2017
Dương Giả Dựng
iii
TÓM TẮT
Đề tài nghiên cứu tính toán tối ưu hoá đa mục tiêu cho kết cấu dàn 2D và 3D
dựa trên độ tin cậy. Bài toán tối ưu hóa đa mục tiêu được thành lập với hàm mục tiêu
là cực tiểu khối lượng và cực tiểu chuyển vị của toàn bộ hệ kết cấu. Biến thiết kế là
diện tích mặt cắt ngang của các thanh dàn. Hàm ràng buộc là các ràng buộc ứng xử
của kết cấu dàn như: ràng buộc chuyển vị nút, ràng buộc ứng suất, ràng buộc về tần
số dao động. Ứng xử của kết cấu dàn được phân tích bằng phương pháp phần tử hữu
hạn (PP-PTHH) – phần tử thanh hai nút tuyến tính. Giải thuật di truyền
Nondominated Sorting Genetic Algorithm II (NSGA-II) được sử dụng để giải bài
toán tối ưu sau khi được thành lập. Độ tin cậy cũng như tính hiệu quả của cách tiếp
cận trong luận văn này được đánh giá và kiểm chứng bằng cách so sánh các kết quả
hiện tại với các kết quả đã công bố trước đó.
Từ khóa: tối ưu hóa đa mục tiêu, tối ưu hóa đa mục tiêu dựa trên độ tin cậy,
NSGA-II (nondominated sorting genetic algorithm II), kết cấu dàn, chuyển vị, ứng
suất, tần số dao động
iv
ABSTRACT
The research aims to formulate and solve multi-objective reliability-based
design optimization problems for 2D and 3D frame structures. The objective
functions are to minimize the weight and the maximum displacement of the entire
structures. The cross-section areas of bars are considered as design variables. The
constraint functions include the behavior constraints of structures such as node
displacements, element stresses and natural frequencies. The behavior of structures
is analyzed by the Finite Element Method (FEM) – using two-node linear truss
elements. Then, Nondominated Sorting Genetic Algorithm II (NSGA-II) is applied
to solve the multi-objective reliability-based design optimization problems. The
reliability and the effectiveness of the proposed approach in this study are evaluated
and verified by comparing the achieved results with those available in the literature.
Keywords: Multi-objective optimization, Multi-objective optimization based
on reliability, NSGA-II (nondominated sorting genetic algorithm II), truss structures,
displacements, stresses, frequency constraint.
v
MỤC LỤC
1.1. Giới thiệu và đặt vấn đề ................................................................................... 1
1.2. Tình hình nghiên cứu thế giới .......................................................................... 6
1.3. Tình hình nghiên cứu trong nước ..................................................................... 7
1.4. Nội dung nghiên cứu ........................................................................................ 9
1.5. Phạm vi nghiên cứu ........................................................................................ 10
1.6. Phương Pháp nghiên cứu ............................................................................... 10
Thành lập ma trận độ cứng, ma trận khối lượng và véc-tơ lực ............... 11
Phương trình cân bằng tĩnh và động ........................................................ 20
2.1. Phương pháp phần tử hữu hạn cho kết cấu dàn ............................................. 11
Tối ưu hóa đa mục tiêu ............................................................................ 22
Bài toán tối ưu hóa đa mục tiêu dựa trên độ tin cậy ................................ 26
Tối ưu hóa đa mục tiêu dựa trên độ tin cậy cho kết cấu dàn ................... 27
2.2. Tối ưu hóa đa mục tiêu dựa trên độ tin cậy cho bài toán kết cấu ................... 22
Chuyển đổi không gian thiết kế xác suất sang không gian thiết kế tiền định
2.3. Giải thuật tối ưu hóa dựa trên độ tin cậy một vòng lặp đơn xác định ............ 28
Thành lập ràng buộc tiền định ................................................................. 32
............................................................................................................... 31
vi
Một số định nghĩa và thuật toán sử dụng trong phương pháp NSGA–II 33
Vòng lặp chính của phương pháp NSGA-II ............................................ 39
2.4. Giải thuật di truyền phân loại không trội NSGA-II ....................................... 33
2.5. Giải thuật một vòng lặp đơn RBDO sử dụng phương pháp NSGA-II ........... 41
Ví dụ 1 – Phân tích chuyển vị và ứng suất của dàn 11 thanh .................. 44
Ví dụ 2 – Phân tích tần số dao động của dàn 10 thanh ............................ 45
Ví dụ 3 – Kiểm chứng code lập trình giải thuật NSGA-II ...................... 46
Ví dụ 4 – Kiểm chứng code lập trình RBDO-NSGA-II .......................... 48
3.1. Kiểm chứng code lập trình Matlab ................................................................. 43
Bài toán kết cấu dàn phẳng 10 thanh ....................................................... 50
Bài toán kết cấu dàn không gian 72 thanh ............................................... 61
Bài toán kết cấu dàn phẳng 200 thanh ..................................................... 75
3.2. Bài toán tối ưu hóa đa mục tiêu kết cấu dàn .................................................. 50
4.1. Kết luận .......................................................................................................... 91
4.2. Hướng phát triển ............................................................................................ 94
vii
DANH MỤC CÁC TỪ VIẾT TẮT
Phần Tử Hữu Hạn PTHH
Finite Element Method FEM
MOEA/D Multi-Objective Evolutionary Algorithm Based on Decomposition
NSGA-II Non-Dominated Sorting Genetic Algorithm -II
Improved Differential Evolution IDE
Double Loop Method DLM
Differential Evolution DE
First Order Reliability Method FORM
Genetic Algorithm GA
Single Loop Deterministic Method SLDM
Reliability Based Design Optimization RBDO
Dertemined Optimization DO
viii
DANH MỤC CÁC BẢNG
Bảng 2.1: Bảng so sánh ưu nhược điểm giữa tối ưu hóa đơn mục miêu và đa mục tiêu.
................................................................................................................. 23
Bảng 3.1: Kết quả chuyển vị trong các nút. .............................................................. 44
Bảng 3.2: Kết quả ứng suất trong các thanh. ............................................................ 45
Bảng 3.3: Kết quả tần số dao động trong hệ kết cấu. ................................................ 46
Bảng 3.4: Kết quả tối ưu Pareto. ............................................................................... 48
Bảng 3.5: Thông số bài toán kết cấu dàn phẳng 10 thanh. ....................................... 50
Bảng 3.6: Kết quả tối ưu hóa DO kết cấu dàn 10 thanh. .......................................... 54
Bảng 3.7: Kết quả tối ưu hóa RBDO kết cấu dàn 10 thanh. ..................................... 58
Bảng 3.8: Kết quả thiết kế RBDO của bài toán kết cấu 10 thanh ở các mức độ an toàn
khác nhau. ................................................................................................ 60
Bảng 3.9: Thông số bài toán kết cấu dàn không gian 72 thanh. ............................... 62
Bảng 3.10: Hai trường hợp tải tác dụng lên kết cấu dàn không gian 72 thanh. ........ 62
Bảng 3.11: Kết quả tối ưu hóa DO kết cấu dàn 72 thanh. ........................................ 68
Bảng 3.12: Kết quả tối ưu hóa RBDO kết cấu dàn 72 thanh. ................................... 73
Bảng 3.13: Kết quả thiết kế RBDO của bài toán kết cấu 72 thanh ở các mức độ an
toàn khác nhau. ........................................................................................ 74
Bảng 3.14: Nhóm thanh diện tích trong kết cấu dàn 200 thanh. ............................... 77
Bảng 3.15: Thông số bài toán kết cấu dàn phẳng 200 thanh. ................................... 78
Bảng 3.16: Kết quả tối ưu hóa DO kết cấu dàn 200 thanh. ...................................... 82
Bảng 3.17: Kết quả tối ưu hóa RBDO kết cấu dàn 200 thanh. ................................. 85
Bảng 3.18: Kết quả thiết kế RBDO của bài toán kết cấu 200 thanh ở các mức độ an
toàn khác nhau. ........................................................................................ 88
ix
DANH MỤC CÁC BIỂU ĐỒ, ĐỒ THỊ, SƠ ĐỒ, HÌNH ẢNH
Hình 1.1: Tòa nhà TTTM. ........................................................................................... 1
Hình 1.2: Bãi đậu xe tại TP. Hà Nội. .......................................................................... 1
Hình 1.3: SVĐ Điện Biên. .......................................................................................... 2
Hình 1.4: Nhà ga Milano Centrale, Ý. ........................................................................ 2
Hình 1.5: Dàn khoan Hibernia, Canada. ..................................................................... 2
Hình 1.6: Cột truyền tải điện. ...................................................................................... 2
Hình 2.1: Kết cấu dàn trong hệ tọa độ địa phương và tổng thể. ............................... 12
Hình 2.2: Kết cấu dàn 2D trong hệ tọa độ địa phương và tổng thể. ......................... 12
Hình 2.3: Chuyển vị của phần tử dàn trong hệ tọa độ địa phương và tổng thể. ....... 13
Hình 2.4: Phần tử dàn Ωe với nút 1 và nút 2 mỗi đầu và tải trọng tác dụng lên hai nút.
................................................................................................................. 14
Hình 2.5: Mô tả kết quả của một bài toán tối ưu hóa, với x* là lời giải tối ưu. ......... 25
Hình 2.6: Mô tả sự ảnh hưởng của các yếu tố ngẫu nhiên lên kết quả tối ưu hoá. ... 26
Hình 2.7: Mô tả kết quả của một bài toán tối ưu hóa dựa trên độ tin cậy, với x* là lời
giải tối ưu. ................................................................................................ 27
Hình 2.8: Thiết lập bài toán tối ưu hóa đa mục tiêu kết cấu dàn. ............................. 28
Hình 2.9: Sơ đồ giải thuật vòng lặp đôi. ................................................................... 29
Hình 2.10: Sơ đồ giải thuật vòng lặp tuần tự. ........................................................... 30
Hình 2.11: Sơ đồ giải thuật một vòng lặp đơn xác định. .......................................... 31
Hình 2.12: Minh họa không gian thiết kế của bài toán. ............................................ 32
Hình 2.13: Thuật toán xác định vùng an toàn (crowding-distance-assignment). ..... 37
Hình 2.14: Sơ đồ giải thuật NSGA-II. ...................................................................... 40
Hình 2.15: Các bước của vòng lặp chính giải thuật NSGA-II. ................................. 41
Hình 2.16: Sơ đồ giải thuật một vòng lặp đơn sử dụng giải thuật NSGA-II ............ 42
Hình 3.1: Mô hình bài toán dàn phẳng 11 thanh. ...................................................... 44
Hình 3.2: Mô hình bài toán dàn phẳng 10 thanh. ...................................................... 46
Hình 3.3: Kết quả bài toán tối ưu. ............................................................................. 47
x
Hình 3.4: Kết quả bài toán tối ưu. ............................................................................. 49
Hình 3.5: Mô hình kết cấu dàn phẳng 10 thanh. ....................................................... 50
Hình 3.6: Mô hình bài toán tối ưu DO kết cấu dàn 10 thanh. ................................... 52
Hình 3.7: Kết quả tối ưu hóa DO kết cấu dàn 10 thanh. ........................................... 53
Hình 3.8: Mô hình bài toán tối ưu RBDO kết cấu dàn 10 thanh. ............................. 55
Hình 3.9: Kết quả tối ưu hóa RBDO kết cấu dàn 10 thanh....................................... 59
Hình 3.10: Kết quả thiết kế RBDO của bài toán kết cấu 10 thanh ở các mức độ an
toàn khác nhau. ........................................................................................ 61
Hình 3.11: Mô hình kết cấu dàn không gian 72 thanh. ............................................. 63
Hình 3.12: Mô hình tối ưu DO bài toán kết cấu dàn 72 thanh. ................................. 64
Hình 3.13: Kết quả tối ưu hóa DO kết cấu dàn 72 thanh. ......................................... 66
Hình 3.14: Mô hình bài toán tối ưu RBDO kết cấu dàn 72 thanh. ........................... 69
Hình 3.15: Kết quả tối ưu hóa RBDO kết cấu dàn 72 thanh. .................................... 70
Hình 3.16: Kết quả thiết kế RBDO của bài toán kết cấu 72 thanh ở các mức độ an
toàn khác nhau. ........................................................................................ 75
Hình 3.17: Mô hình kết cấu dàn phẳng 200 thanh .................................................... 76
Hình 3.18: Mô hình bài toán tối ưu DO kết cấu dàn 200 thanh. ............................... 79
Hình 3.19: Kết quả tối ưu hóa DO kết cấu dàn 200 thanh. ....................................... 80
Hình 3.20: Mô hình bài toán tối ưu RBDO kết cấu dàn 200 thanh. ......................... 84
Hình 3.21: Kết quả tối ưu RBDO-NSGA-II kết cấu dàn 200 thanh. ........................ 87
Hình 3.22: Kết quả thiết kế RBDO của bài toán kết cấu 200 thanh ở các mức độ an
toàn khác nhau. ........................................................................................ 90
xi
DANH MỤC KÝ HIỆU
Ma trận độ cứng tổng thể kết cấu K
Ma trận khối lượng tổng thể kết cấu M
Véc-tơ ngoại lực tác động lên các nút của kết cấu F
Véc-tơ diện tích các thanh A
Véc-tơ chuyển vị u
ω Tần số góc
f Tần số dao động tự nhiên
βt Chỉ số độ tin cậy yêu cầu
Khối lượng cộng thêm tại các nút W0
P Tải trọng tác dụng lên hệ kết cấu
E Mô-đun đàn hồi
ρ Khối lượng riêng vật liệu
umin; umax Giới hạn chuyển vị
σmin; σmax Giới hạn ứng suất
Giới hạn tần số dao động
1
TỔNG QUAN
1.1. Giới thiệu và đặt vấn đề
Kết cấu dàn là loại kết cấu phổ biến và ngày càng được sử dụng rộng rãi trong
nhiều công trình xây dựng đặc biệt là các công trình xây dựng dân dụng và công
nghiệp. Bởi vì kết cấu dàn có rất nhiều ưu điểm so với các loại kết cấu khác (bê tông
cốt thép, kết cấu gỗ, kết cấu gạch - đá v.v.) như: khẩu độ vượt lớn, nhẹ, kinh tế và đặc
biệt về phương diện kiến trúc có thể tạo được nhiều hình dáng khác nhau có thể kể
đến như: vòm cầu, vòm trụ, vòm yên ngựa v.v. mà hiện nay đang được sử dụng phổ
biến ở nhiều công trình trên thế giới. Một số công trình tiêu biểu sử dụng kết cấu dàn
có thể kể đến như: công trình dân dụng cao tầng (Hình 1.1), bãi đậu xe (Hình 1.2),
các cột truyền tải điện (Hình 1.6), cột truyền thông, dàn khoan (Hình 1.5), mái che
cho các công trình sân vận động (Hình 1.3), nhà thi đấu, cung thể thao, trung tâm
thương mại, nhà ga (Hình 1.4), xưởng sửa chữa bảo dưỡng máy bay v.v.
Hình 1.1: Tòa nhà TTTM. (Nguồn: http://batdongsanviet.com.vn) Hình 1.2: Bãi đậu xe tại TP. Hà Nội. (Nguồn: https://www.baotintuc.vn)
2
Hình 1.3: SVĐ Điện Biên. (Nguồn:http://www.phongnhajsc.com) Hình 1.4: Nhà ga Milano Centrale, Ý. (Nguồn: https://diaoconline.vn)
Hình 1.5: Dàn khoan Hibernia, Canada. (Nguồn: http://kienthuc.net.vn) Hình 1.6: Cột truyền tải điện. (Nguồn:https://pxhere.com)
Trong thực tế, để việc thiết kế, chế tạo kết cấu dàn vừa đảm bảo độ bền, độ ổn
định, cũng như đảm bảo các yếu tố cạnh tranh như tiết kiệm vật liệu, v.v., việc thiết
lập và giải các bài toán tối ưu thiết kế cho kết cấu dàn luôn là một vấn đề quan trọng
và nhận được sự quan tâm lớn của nhiều nhà nghiên cứu trên thế giới và trong nước.
Trong việc thiết kế kết cấu dàn làm việc trong môi trường chịu tải trọng, có
hai mục tiêu quan trọng nhất cần phải được thỏa mãn. Thứ nhất, khối lượng kết cấu
nên càng nhẹ càng tốt để giảm chi phí xây dựng. Thứ hai, kết cấu phải đáp ứng được
yêu cầu làm việc trong cả môi trường tĩnh và động. Cụ thể, đối với kết cấu làm việc
trong môi trường tĩnh, kết cấu sẽ bị phá hủy nếu tải trọng bên ngoài tác động vào quá
3
lớn làm ứng xử kết cấu vượt giới hạn chuyển vị hay vượt giới hạn ứng suất của kết
cấu. Trong bài toán thiết kế tối ưu để tránh các tình huống này xảy ra, các ràng buộc
về chuyển vị và ứng suất cần được áp đặt. Tương tự đối với kết cấu làm việc trong
môi trường động, nếu tần số dao động của tải trọng động lực học cộng hưởng tại một
trong những tần số riêng của kết cấu, thì kết cấu sẽ bị khuếch đại dao động và dẫn
đến sự phá hoại của kết cấu. Trong bài toán thiết kế tối ưu, để tránh các tình huống
này xảy ra, các ràng buộc động học cũng cần được áp đặt.
Thông thường, để giải các bài toán tối ưu dàn, các nhà nghiên cứu thường tập
trung giải bài toán đơn mục tiêu liên quan đến cực tiểu chuyển vị, hoặc cực đại tần
số, hoặc cực tiểu khối lượng. Điều này sẽ dẫn đến những bất lợi sau: (1) chỉ tìm được
kết quả thiết kế tối ưu thỏa một mục tiêu duy nhất, trong khi đó không xét ảnh hưởng
của nó đến các mục tiêu khác; (2) chỉ cho được một giá trị tối ưu duy nhất. Trong khi
đó, người sử dụng thực tế luôn có nhu cầu tối ưu nhiều mục tiêu cùng lúc và muốn
có nhiều kết quả thiết kế tối ưu để lựa chọn tùy theo khả năng và yêu cầu của họ.
Chính vì vậy, bài toán tối ưu đa mục tiêu được phát triển và nghiên cứu cho kết cấu
dàn. Trong bài toán tối ưu đa mục tiêu, các mục tiêu thường mâu thuẫn với nhau (ví
dụ như chi phí sản xuất và độ bền kết cấu) và các mục tiêu này được tối ưu cùng lúc
với nhau. Lời giải của bài toán tối ưu đa mục tiêu sẽ là một tập hợp của các giải pháp
tối ưu có quan hệ đối nghịch giữa các hàm mục tiêu. Do đó không có một nghiệm
duy nhất thỏa mãn các mục tiêu này, mà thay vào đó là một tập hợp các nghiệm thỏa
mãn các mục tiêu này, và tập hợp nghiệm này được gọi là tập nghiệm Pareto [1].
Để giải bài toán tối ưu hóa đa mục tiêu, một số phương pháp tối ưu cổ điển đã
được đề xuất như: phương pháp tổng trọng số (Weight Sum Method); phương pháp
ràng buộc ε (ε - Constraint Method); v.v. Những phương pháp này tuy có ưu điểm là
đơn giản và đảm bảo tìm được nghiệm trên đường tối ưu Pareto lồi, nhưng cũng tồn
tại những hạn chế như: chỉ tìm được tối đa một nghiệm trong một lần giải, gặp khó
khăn khi giải quyết bài toán có miền mục tiêu phi tuyến, bị phụ thuộc nhiều vào trọng
số (wi) và các ràng buộc (εi). Để giải quyết các khó khăn trên, năm 2002 Deb và cộng
sự [2] đã đề xuất một phương pháp mới là Elitist Non-Dominated Sorting Genetic
4
Algorithm (NSGA-II). Đây là một phương pháp tối ưu hóa đa mục tiêu thuộc nhóm
tìm kiếm trực tiếp dựa trên nguyên lý tiến hóa (Evolutionary Multi-objective
Optimization: EMO). Khác với các phương pháp tối ưu hóa cổ điển, khi chỉ tìm được
một nghiệm duy nhất trên đường Pareto, phương pháp EMO này tìm được cả một tập
hợp nghiệm tối ưu trên đường Pareto. Hơn nữa, so với các phương pháp khác như:
giải thuật di truyền GA (Genetic Algorithm); giải thuật tiến hóa DE (Differential
Evolution) [3]; giải thuật đàn kiến ACO (Ant Colony Optimization) [4]; giải thuật
bầy đàn PSO (Particle Swarm Optimization) [5]; giải thuật MOEA/D (A Multi-
Objective Evolutionary Algorithm Based on Decomposition) [6], v.v., phương pháp
này còn có thời gian tính toán khá nhanh.
Một vấn đề nữa cũng luôn nhận được sự quan tâm của các nhà khoa học cũng
như các nhà thiết kế đó là kết quả của bài toán tối ưu hoá luôn nằm trên các ràng buộc
của bài toán. Tức là kết quả tối ưu sẽ nằm ngay ranh giới cho phép của bài toán. Tuy
nhiên dưới tác động của các yếu tố ngẫu nhiên (điều kiện khí hậu, mưa, gió, lũ lụt,
hạn hán, động đất, v.v), hoặc yếu tố con người (sự dao động của ngoại lực tác động
trong quá trình sử dụng, sai số trong tính toán, sai số do mô hình, sai sót do thi công,
sản xuất, v.v) làm cho dữ liệu đầu vào trở thành các biến ngẫu nhiên phân bố quanh
giá trị thiết kế ban đầu và tuân theo một qui luật phân phối xác suất nhất định, thì ứng
xử kết cấu thực tế sẽ phân bố ngẫu nhiên quanh nghiệm thiết kế tối ưu nằm trên ranh
giới của miền an toàn và miền phá hủy. Khi đó xác suất phá hủy của kết cấu sẽ tăng
cao và xấp xỉ khoảng 50%. Trong trường hợp này kết cấu sẽ không đảm bảo an toàn
như mong muốn của người thiết kế. Do đó việc tính toán xác suất phá hủy của kết
cấu (hay phân tích độ tin cậy kết cấu) trong bài toán thiết kế tối ưu là rất cần thiết khi
ta xem xét yếu tố phân bố ngẫu nhiên của các biến thiết kế. Để đánh giá độ tin cậy
của kết cấu, một số phương pháp phổ biến thường được sử dụng như: phương pháp
đánh giá độ tin cậy bậc nhất FORM (First Order Reliability Method), độ tin cậy bậc
hai SORM (Second Order Reliability Method) hay phương pháp mô phỏng Monte
Carlo. Như vậy để đảm bảo độ tin cậy nhất định cho kết quả thiết kế tối ưu khi xem
xét đến yếu tố ngẫu nhiên, việc kết hợp giữa giải thuật tối ưu hóa với các giải thuật
5
phân tích độ tin cậy là rất cần thiết. Khi đó bài toán được gọi là bài toán tối ưu hóa
dựa trên độ tin cậy RBDO (Reliability Based Design Optimization).
Giải một bài toán tối ưu hóa đa mục tiêu dựa trên độ tin cậy (RBDO), cần có
hai vòng lặp tính toán gồm một vòng lặp tối ưu hóa và một vòng lặp đánh giá độ tin
cậy, nên chi phí tính toán để giải một bài toán tối ưu dựa trên độ tin cậy là tương đối
lớn. Để giải bài toán tối ưu dựa trên độ tin cậy, Li và cộng sự [7] đưa ra 3 cách tiếp
cận chính bao gồm: Giải thuật vòng lặp đôi DLM (Double Loop Method); Giải thuật
vòng lặp tách rời hay còn gọi là giải thuật lặp tuần tự DDLM (Decouple Double Loop
Method); và giải thuật một vòng lặp đơn xác định SLDM (Single Loop Deterministic
Method). Trong 3 giải thuật trên, giải thuật một vòng lặp đơn xác định được xem là
một trong những giải thuật hiệu quả nhất cho đến thời điểm hiện tại [8], vì khả năng
cân đối tốt giữa chi phí tính toán và độ chính xác của nghiệm tối ưu.
Từ những lập luận trên, có thể thấy rằng tối ưu hóa đa mục tiêu cho kết cấu
dàn dựa trên độ tin cậy có ý nghĩa quan trọng và việc đưa ra một phương pháp giải
tối ưu phù hợp giúp đảm bảo sự vận hành an toàn của kết cấu khi đưa vào sử dụng và
giảm chi phí xây dựng là yêu cầu cần thiết hiện nay. Do đó, luận văn này được thực
hiện với các mục đích sau:
+ Nghiên cứu áp dụng giải thuật NSGA-II trong việc giải các bài toán tối ưu
đa mục tiêu kết cấu dàn nhằm khắc phục các nhược điểm của các nghiên cứu trước
đây. Nghiên cứu xem xét đồng thời các ràng buộc tĩnh học và động học tác động lên
kết cấu;
+ Áp dụng giải thuật NSGA-II để giải bài toán tối ưu hóa dựa trên độ tin cậy
bằng cách kết hợp nó với giải thuật tối ưu hóa dựa trên độ tin cậy một vòng lặp đơn
xác định với chỉ số độ tin cậy cho trước nhằm đảm bảo độ an toàn cho phép kết cấu
dàn khi đưa vào vận hành sử dụng trong các điều kiện ngẫu nhiên cho trước.
Hiệu quả và độ tin cậy của phương pháp đề xuất trong luận văn sẽ được đánh
giá khi so sánh kết quả của phương pháp hiện tại với kết quả của các phương pháp
khác đã được nghiên cứu trước đây.
6
1.2. Tình hình nghiên cứu thế giới
Trong những năm gần đây, việc áp dụng các giải thuật tối ưu cho các bài toán
tối ưu trong lĩnh vực xây dựng được phát triển rất mạnh mẽ. Nhiều cách tiếp cận cũng
như nhiều phương pháp tối ưu hóa khác nhau đã được đề xuất để giải bài toán tối ưu
hóa cho kết cấu dàn. Một số nghiên cứu ứng dụng có thể tìm thấy ở các tài liệu tham
khảo:
- Nghiên cứu tối ưu hóa đơn mục tiêu chưa xét đến độ tin cậy có thể kể đến
các công trình nghiên cứu của Farshchin và công sự (2016) [9]; A.Kaveh và Ghazaan
(2015) [10]; Mohsen Khatibinia và Naseralavi (2014) [11]; Pezeshk và cộng sự
(2000) [12]; Camp và cộng sự (2005) [13]; Degertekin (2008) [14]; Saka (2009) [15];
Kaveh và Talatahari (2010) [16]; Kaveh và Abbasgholiha (2011) [17]; Dogan và
SaKa (2012) [18]; Kaveh và cộng sự (2009) [19]; Li và cộng sự (2007) [20]; Lee và
cộng sự (2004) [21]; Camp và cộng sự (1998) [22] v.v. Nhiều giải thuật tìm kiếm đã
được sử dụng trong các nghiên cứu này bao gồm: giải thuật tiến hóa DE (Differential
Evolution); giải thuật di truyền GA (Genetic Algorirthm); giải thuật đàn kiến ACO
(Ant Colony Optimization); giải thuật tìm kiếm hòa hợp lấy cảm hứng từ âm nhạc
HS (Hamony Search algorithm); giải thuật tối ưu hóa bầy đàn PSO (Particle Swarm
Optimization); và giải thuật về vụ nổ lớn BB-BC (Big Bang – Big Crunch).
- Các nghiên cứu về tối ưu hóa đa mục tiêu có thể kể đến như: Duong-Gia. D
và cộng sự (2017) [23] đã giải bài toán tối ưu hóa khối lượng và tần số kết cấu dầm
Composite sử dụng giải thuật NSGA-II; Coello và Christiansen (2000) [24] nghiên
cứu tối ưu hóa khối lượng, chuyển vị và ứng suất của các kết cấu dàn sử dụng giải
thuật di truyền GA (Genetic Algorirthm); Omkar và cộng sự (2008) [25] nghiên cứu
tối ưu hoá khối lượng và tổng chi phí của tấm composite sử dụng giải thuật VEPSO
(Vector evaluated particle swarm optimization); Gurugubelli và Kallepalli (2014)
[26] nghiên cứu tối ưu hoá khối lượng và chuyển vị của dầm consol sử dụng giải thuật
NSGA–II (Non-Dominated sorting Genetic Algorithm). Ngoài ra còn có các nghiên
cứu của Coelho và cộng sự (2011) [27]; Poldee và cộng sự (2012) [1]; Poldee và cộng
sự (2013) [28]; Angelo và cộng sự (2015) [29] v.v.
7
- Trong lĩnh vực tối ưu hóa dựa trên độ tin cậy, việc kết hợp các thuật toán tối
ưu hóa tìm kiếm trực tiếp với các phương pháp đánh giá độ tin cậy cũng nhận được
nhiều sự quan tâm. Một số nghiên cứu gần đây có thể kể đến như: Dimou và
Koumousis (2009) [30] nghiên cứu tối ưu hóa dựa trên độ tin cậy cho kết cấu dàn sử
dụng giải thuật tối ưu hóa PSO; Ghasemi và Yousefi (2011) [31] nghiên cứu tối ưu
dựa trên độ tin cậy cho khung thép 2D, 3D sử dụng giải thuật di truyền sửa đổi kết
hợp với phương pháp phân tích độ tin cậy bậc hai SORM; Shayanfar và cộng sự
(2014) [32] nghiên cứu tối ưu dựa trên độ tin cậy cho khung thép 2D sử dụng phương
pháp GA kết hợp với phương pháp phân tích độ tin cậy bậc nhất FORM. Ngoài ra
còn có nghiên cứu của Lee và cộng sự (2002) [33] v.v.
Gần đây Yu Su và cộng sự (2016) [34] đã áp dụng giải thuật DEMO để giải
bài toán tối ưu hóa đa mục tiêu kết cấu dàn với hai mục tiêu gồm cực tiểu khối lượng
và cực tiểu mô-men. Tuy nhiên nghiên cứu này chỉ xét đến ràng buộc về chuyển vị
và ứng suất chứ chưa xét đến ràng buộc tần số dao động. Và mới nhất Lu-Kai Song
và cộng sự (2017) [35] đã áp dụng giải thuật DMOPSO để giải tối ưu hóa đa mục tiêu
dựa trên độ tin cậy, với đối tượng của nghiên cứu là Turbine.
Theo dữ liệu thu thập của tác giả, cho đến thời điểm hiện tại, việc phát triển
và áp dụng giải thuật NSGA-II vào các bài toán tối ưu hóa dựa trên độ tin cậy cho
các bài toán kết cấu dàn vẫn chưa được thực hiện. Đặc biệt chưa có nghiên cứu về sự
kết hợp của NSGA-II với giải thuật tối ưu hóa dựa trên độ tin cậy một vòng lặp đơn
xác định chịu các ràng buộc tĩnh học và động học.
1.3. Tình hình nghiên cứu trong nước
Theo thu thập của tác giả cho đến thời điểm hiện tại, việc nghiên cứu và áp
dụng tối ưu hóa trong lĩnh vực xây dựng nói chung và kết cấu dàn nói riêng vẫn còn
hạn chế. Các nghiên cứu chủ yếu được công bố ở các hội nghị chuyên ngành trong
nước, và ít có các công bố được đăng trên các tạp chí uy tín trên thế giới. Đặc biệt là
bài toán tối ưu hóa dựa trên độ tin cậy, nó vẫn còn khá mới mẻ với nhiều nhà nghiên
cứu ở Việt Nam. Có thể điểm qua một số nghiên cứu sau đây:
8
- Luận văn thạc sĩ của Trần Văn Dần [36] (Đại học Bách Khoa TP.HCM, năm
2013): tối ưu hóa dựa trên độ tin cậy tấm composite laminate bằng giải thuật di truyền.
Trong nghiên cứu này, hàm mục tiêu là cực tiểu khối lượng và tổng năng lượng biến
dạng của tấm chịu các ràng buộc về chiều dày, góc hướng sợi, chuyển vị, và ứng suất.
- Luận văn thạc sĩ của Thái Bình Quốc [37] (ĐH Công Nghệ TP.HCM, năm
2015): tối ưu hóa kết cấu dàn với biến diện tích rời rạc sử dụng phương pháp DE cải
tiến. Trong nghiên cứu này, hàm mục tiêu là cực tiểu khối lượng của kết cấu chịu các
ràng buộc về chuyển vị và ứng suất.
- Luận văn thạc sĩ của Nguyễn Thụy Đoan Nhi [38] (Đại học Tôn Đức Thắng,
năm 2014): tối ưu hóa đa mục tiêu kết cấu tấm gấp composite dùng phương pháp
NSGA-II. Trong nghiên cứu này, hàm mục tiêu là cực tiểu khối lượng và cực tiểu
năng lượng biến dạng của tấm chịu các ràng buộc về chuyển vị và ứng suất.
- Luận văn thạc sĩ của Hồ Hữu Vịnh [39] (Đại học Bách Khoa TP.HCM, năm
2016): phát triển giải thuật tiến hóa DE cho bài toán tối ưu hóa dựa trên độ tin cậy
với biến thiết kế rời rạc.
- Tối ưu hóa đa mục tiêu kết cấu dầm Composites nhiều lớp chịu ràng buộc về
tần số của Dương Giả Dựng, Nguyễn Thanh Liêm, Võ Duy Trung, Hồ Hữu Vịnh và
Nguyễn Thời Trung [40] (Hội nghị Khoa học toàn quốc Vật liệu và Kết cấu
Composites Cơ học, Công nghệ và ứng dụng, năm 2016). Trong nghiên cứu này, các
tác giả đã sử dụng giải thuật NSGA-II để tính toán tối ưu cho dầm composite nhiều
lớp với biến thiết kế là chiều dày (t) và tỷ lệ phân bố sợi trong tấm (rf).
Gần đây luận văn thạc sĩ của Bùi Trần Vĩnh Thái (Đại học Mở, năm 2016)
[41]: đã áp dụng giải thuật NSGA-II kết hợp với giải thuật tối ưu hóa dựa trên độ tin
cậy một vòng lặp đơn xác định cho kết cấu. Tuy nhiên luận văn chỉ xét đến ràng buộc
về chuyển vị chứ chưa xét đến ràng buộc tần số dao động, và đối tượng nghiên cứu
là khung thép. Và mới nhất là luận văn thạc sĩ của Lâm Thị Phúc Hạnh (ĐH Công
Nghệ TP.HCM, năm 2016) [42]: tối ưu hóa dựa trên độ tin cậy kết cấu khung thép sử
dụng giải thuật một vòng lặp đơn xác định. Luận văn cũng đã sử dụng giải thuật tối
ưu hóa dựa trên độ tin cậy một vòng lặp đơn xác định. Tuy nhiên luận văn chỉ mới
9
tối ưu đơn mục tiêu cho kết cấu, và ràng buộc của bài toán chỉ xét đến độ bền của kết
cấu và chuyển vị, chứ chưa xét đến ràng buộc về tần số dao động, và đối tượng nghiên
cứu là khung thép.
Kết luận
Hầu hết các nghiên cứu được liệt kê bên trên đều đạt được những kết quả phù
hợp theo yêu cầu. Tuy nhiên các nghiên cứu trên chỉ thực hiện giải các bài toán tối
ưu hoá mà không xét đến độ tin cậy hoặc nếu có xét đến độ tin cậy thì chỉ với bài toán
tối ưu hoá đơn mục tiêu, chứ chưa có hoặc có rất ít các công trình nghiên cứu về tối
ưu hoá đa mục tiêu dựa trên độ tin cậy cho kết cấu dàn. Ngoài ra chi phí tính toán đối
với bài toán tối ưu hóa này vẫn còn lớn vì các tác giả còn sử dụng các phương pháp
cũ, chưa có nhiều cải tiến. Chính vì vậy nhằm giúp các nhà thiết kế có thêm các
phương án thiết kế tối ưu cho các sản phẩm thiết kế của mình mà vẫn đảm bảo an
toàn cho kết cấu vận hành trong các điều kiện ngẫu nhiên, luận văn này được thực
hiện nhằm thiết lập và giải bài toán tối ưu hoá đa mục tiêu dựa trên độ tin cậy cho kết
cấu dàn sử dụng giải thuật một vòng lặp đơn xác định và dùng phương pháp NSGA-
II để tính toán. Ngoài ra việc tính toán tối ưu hóa cho kết cấu dàn còn có ý nghĩa giúp
phát triển và mở rộng ứng dụng các giải thuật tối ưu cho các bài toán xây dựng. Từ
đó giúp thúc đẩy sự phát triển của hướng nghiên cứu tối ưu hoá ứng dụng nhằm đáp
ứng cùng lúc nhiều mục tiêu khác nhau như giảm chi phí và tiết kiệm thời gian của
các doanh nghiệp xây dựng.
1.4. Nội dung nghiên cứu
Nghiên cứu nhằm thực hiện những mục tiêu sau:
- Thiết lập bài toán tối ưu hoá đa mục tiêu tiền định.
- Thiết lập bài toán tối ưu hoá đa mục tiêu dựa trên độ tin cậy.
- Tìm hiểu phương pháp PTHH dùng cho kết cấu dàn.
- Tìm hiểu giải thuật NSGA-II (Non-Dominated sorting Genetic Algorithm).
- Tìm hiểu bài toán phân tích độ tin cậy cho kết cấu và giải thuật tối ưu hóa dựa
trên độ tin cậy một vòng lặp đơn xác định.
- Áp dụng thuật toán NSGA-II để tìm lời giải tối ưu bài toán đã được thiết lập:
10
+ Tối ưu hoá đa mục tiêu kết cấu dàn phẳng 2D và dàn không gian 3D khi không
xét đến độ tin cậy.
+ Tối ưu hoá đa mục tiêu kết cấu dàn phẳng 2D và dàn không gian 3D có xét
đến độ tin cậy.
1.5. Phạm vi nghiên cứu
Nghiên cứu được thực hiện trong phạm vi sau: kết cấu dàn phẳng 2D và dàn
không gian 3D. Bài toán tối ưu hóa được thành lập với hai mục tiêu là cực tiểu khối
lượng của kết cấu và cực tiểu chuyển vị của kết cấu dựa trên độ tin cậy chịu ràng buộc
về chuyển vị, ứng suất và tần số dao động.
1.6. Phương Pháp nghiên cứu
Sử dụng phương pháp phần tử hữu hạn để tính tính toán phân tích ứng xử kết
cấu dàn (tính chuyển vị, ứng suất, khối lượng, tần số dao động) trên nền phần mềm
tính toán Matlab.
Thành lập bài toán tối ưu đa mục tiêu tiền định và tối ưu đa mục tiêu dựa trên
độ tin cậy với mức độ tin cậy được xác định trước.
Dựa trên bài toán tối ưu đa mục tiêu tiền định và tối ưu đa mục tiêu dựa trên
độ tin cậy, thành lập các hàm mục tiêu và ràng buộc.
Chuyển đổi các ràng buộc xác suất thành các ràng buộc tiền định bằng giải
thuật vòng lặp đơn (Single Loop Method).
Áp dụng giải thuật NSGA-II để tìm kết quả tối ưu cho bài toán.
11
CƠ SỞ LÝ THUYẾT
Trong chương này, luận văn trình bày lý thuyết tổng quát về phương pháp
phần tử hữu hạn (PP-PTHH) cho kết cấu dàn, lý thuyết tối ưu hoá đa mục tiêu và độ
tin cậy để thiết lập bài toán tối ưu hoá đa mục tiêu dựa trên độ tin cậy cho kết cấu dàn
2D và 3D. Để giải bài toán đã được thiết lập, lý thuyết về giải thuật một vòng lặp đơn
xác định Single-Loop Deterministic Method (SLDM) và giải thuật di truyền phân loại
không trội NSGA-II cũng được trình bày ở phần tiếp theo trong Chương này.
2.1. Phương pháp phần tử hữu hạn cho kết cấu dàn
Để xác định các ràng buộc cho bài toán tối ưu hóa như ràng buộc ứng suất,
chuyển vị, và tần số dao động, v.v., phương pháp phần tử hữu hạn sử dụng phần tử
thanh hai nút tuyến tính được sử dụng. Chi tiết cho việc thành lập và giải bài toán kết
cấu dàn sẽ được trình bày ở phần này [43].
Các phần tử dàn này sẽ được ứng dụng để phân tích chuyển vị, ứng suất và tần
số dao động các hệ kết cấu dàn cả trong không gian hai chiều (dàn phẳng 2D) và ba
chiều (dàn không gian 3D). Trong hệ dàn phẳng 2D, chuyển vị có hai thành phần
trong các hướng x và y. Còn trong hệ dàn không gian 3D, chuyển vị có 3 thành phần
trong các hướng x, y, và z. Trong hệ kết cấu dàn, các phần tử dàn được liên kết với
nhau bởi các khớp ở hai đầu, vì vậy chỉ có lực dọc (không có mô-men) được truyền
giữa các thanh. Để tiện cho việc trình bày, phần này chỉ đề cập các phần tử dàn có
diện tích mặt cắt ngang giống nhau.
Thành lập ma trận độ cứng, ma trận khối lượng và véc-tơ lực
Dàn phẳng 2D
Xét kết cấu dàn nằm trong mặt phẳng OXYZ như trong Hình 2.1, với OXYZ là
hệ tọa độ tổng thể và O’x’y’z’ là hệ tọa độ địa phương.
Y
O’
O
X
Z
12
Hình 2.1: Kết cấu dàn trong hệ tọa độ địa phương và tổng thể.
trong đó hệ dàn phẳng 2D, hệ tọa độ tổng thể chỉ còn là OXY và không có thành phần
Y
Hệ tọa độ tổng thể
X
O
Z như trong Hình 2.2.
Hình 2.2: Kết cấu dàn 2D trong hệ tọa độ địa phương và tổng thể.
Chuyển vị của kết cấu dàn 2D được biểu diễn theo hệ trục tọa độ địa phương
và theo hệ trục tọa độ tổng thể như trong Hình 2.3, trong đó là chuyển vị trong
hệ tọa độ địa phương; là chuyển vị trong hệ trục tọa độ tổng thể. Mối
liên hệ giữa chuyển vị trong hai hệ tọa độ này được biểu diễn như sau:
Chuyển vị trong hệ tọa độ địa phương:Equation Chapter 2 Section 1
(2.1)
và chuyển vị trong hệ tọa độ tổng thể:
13
(2.2)
Phương trình chuyển đổi tọa độ tổng thể sang tọa độ địa phương:
(2.3)
và:
(2.4)
Điều này có thể viết dưới dạng véc-tơ như sau:
(2.5)
với:
e
e
v2
u2
e
sinθ
u2
e
cosθ
u1
e
v1
θ
e
u1
(2.6)
Hình 2.3: Chuyển vị của phần tử dàn trong hệ tọa độ địa phương và tổng thể.
Đầu tiên, ma trận độ cứng, véc-tơ lực và ma trận chuyển vị của dàn 2D được
thành lập trong hệ trục địa phương. Sau đó, sử dụng quan hệ (2.5), các ma trận trong
hệ trục tọa độ tổng thể sẽ được thành lập.
Chi tiết việc thành lập ma trận độ cứng, ma trận khối lượng, véc-tơ lực trong
hệ trục tọa độ địa phương sẽ được trình bày ở mục tiếp theo.
14
2.1.1.1.1 Hệ tọa độ địa phương
a) Hàm chuyển vị và các hàm dạng phần tử
Xét một phần tử dàn Ωe với nút 1 và 2 ở hai đầu như chỉ trong Hình 2.4. Chiều
dài của phần tử là le. Trục địa phương x được lấy theo hướng dọc trục thanh với điểm
gốc O là điểm 1. Trong hệ tọa độ địa phương, chỉ có 1 bậc tự do (BTD) ở mỗi nút của
phần tử. Vì vậy có tổng cộng 2 BTD cho phần tử.
Hình 2.4: Phần tử dàn Ωe với nút 1 và nút 2 mỗi đầu và tải trọng tác dụng lên hai
nút.
Đặt là hàm chuyển vị (hàm nghiệm xấp xỉ) của phần tử Ωe trong hệ tọa
độ địa phương. Hàm chuyển vị này được xấp xỉ dựa trên ma trận của các hàm dạng
và véc-tơ chuyển vị tại nút:
(2.7)
trong đó là véc-tơ chuyển vị nút trong hệ tọa độ địa phương và được sắp xếp theo
thứ tự các nút như sau:
(2.8) Chuyển vị tại nút 1 của phần tử Ωe Chuyển vị tại nút 2 của phần tử Ωe
và là ma trận hàm dạng của phần tử. Cho phần tử dàn, rút gọn lại là
một véc-tơ và được viết theo thứ tự của các nút như sau:
Nút thứ 1 Nút thứ 2
(2.9)
trong đó , (I=1, 2) là hai hàm dạng tương ứng với hai nút của phần tử dàn và
được viết dưới dạng hiện như sau
15
(2.10)
b) Véc-tơ biến dạng phần tử và ma trận biến dạng - chuyển vị của phần
tử
Véc-tơ biến dạng phần tử trong hệ tọa độ địa phương được tính từ véc-tơ
chuyển vị nút của phần tử như sau:
(2.11)
Đặt là ma trận biến dạng – chuyển vị của phần tử Ωe, và
thay (2.11) ta được:
(2.12)
cho phần tử dàn được viết cụ thể lại là một véc-tơ như sau:
(2.13)
Thay hàm dạng , (I=1, 2) trong công thức (2.10) vào công thức (2.13),
véc-tơ có dạng:
(2.14)
c) Ma trận độ cứng phần tử
Ma trận độ cứng phần tử trong hệ tọa độ địa phương có dạng sau:
(2.15)
trong đó Ae là diện tích mặt cắt ngang của phần tử dàn; và ma trận các hằng số vật
liệu D được rút gọn lại thành mô-đun đàn hồi E cho trường hợp phần tử dàn.
d) Ma trận khối lượng phần tử
Ma trận độ cứng phần tử trong hệ tọa độ địa phương có dạng:
16
(2.16)
Xét trường hợp miền bài toán có tỉ trọng khối lượng ρ là hằng số, và thay
trong công thức (2.9) vào công thức (2.16), ta được:
(2.17)
e) Véc-tơ tải phần tử
Xét phần tử dàn chịu tải trọng tập trung tại nút như Hình 2.4.
Véc-tơ tải của phần tử dàn gồm hai thành phần là trọng lượng bản thân của
thanh dàn và lực tập trung tác dụng tại hai nút, được biểu diễn dưới dạng sau:
(2.18)
Trong trường hợp phần tử dàn, biên có thể là một hoặc cả hai điểm cuối 1
và 2 (tức x=0 và x=le) của phần tử và được viết cụ thể lại như sau:
(2.19)
trong đó t(0) và t(le) lần lượt là hai lực tập trung đặt tại 2 điểm cuối 1 và 2 của phần
tử. Tại các điểm cuối 1 và 2 này, ta có giá trị các hàm dạng như sau:
(2.20)
Thừa nhận rằng, phần tử dàn chịu tải trọng phân bố đều dọc theo trục
x, và 2 lực tập trung và tương ứng đặt tại 2 điểm cuối 1 và 2
17
của phần tử, khi đó thay bởi công thức (2.10), và bởi công
thức (2.20) vào véc-tơ trong công thức (2.19) ta được:
(2.21)
2.1.1.1.2 Hệ trục tọa độ tổng thể
Các ma trận và véc-tơ phần tử trong các công thức (2.15), (2.17) và (2.21)
được thành lập dựa trên hệ tọa độ địa phương, trong đó trục x trùng khớp với trục
chính tâm của thanh 1-2 như trong Hình 2.2. Trong thực tế, các thanh dàn được phân
bố theo nhiều hướng và vị trí khác nhau. Để có thể lắp ghép các ma trận và véc-tơ
phần tử để hình thành các ma trận và véc-tơ tổng thể, ta cần thực hiện phép biến đổi
tọa độ để hình thành các ma trận và véc-tơ phần tử trong hệ tọa độ tổng thể.
a) Ma trận biến đổi tọa độ
Ma trận biến đổi tọa độ cho phần tử dàn Te được cho bởi:
(2.22)
và có tính chất trực giao như sau:
(2.23)
Chiều dài của phần tử dàn có thể được tính toán dùng tọa độ tổng thể của
hai nút I và J của phần tử như sau:
(2.24)
b) Các ma trận và véc-tơ phần tử
Véc-tơ tải tại nút fe của phần tử Ωe trong hệ tọa độ tổng thể XY có dạng:
18
(2.25)
Ma trận độ cứng Ke của phần tử Ωe trong hệ tọa độ tổng thể XY có dạng:
(2.26)
Ma trận khối lượng Me của phần tử Ωe trong hệ tọa độ tổng thể XY có dạng:
(2.27)
Dàn không gian 3D
Tương tự như dàn phẳng 2D, tất cả việc thành lập của phép biến đổi tọa độ của
hệ dàn phẳng đều có thể áp dụng lại cho hệ dàn không gian 3D và chỉ cần thêm các
hàng và cột tương ứng với trục Z. Một cách ngắn gọn ta trình bày các kết quả của
phép biến đổi tọa độ cho hệ dàn không gian 3D như sau:
a) Ma trận biến đổi tọa độ
Véc-tơ chuyển vị nút trong hệ tọa độ địa phương liên hệ với véc-tơ chuyển
vị nút de trong hệ tọa độ tổng thể XYZ bởi phép biến đổi tọa độ như sau:
(2.28)
trong đó Te là ma trận biến đổi tọa độ cho phần tử dàn và được cho bởi:
19
(2.29)
với
là các góc cô-sin chỉ phương của trục của phần tử dàn. Ta dễ dàng có kết quả sau:
(2.30)
mà chỉ rằng ma trận Te là một ma trận trực giao. Chiều dài le của phần tử dàn có thể
được tính toán dùng tọa độ tổng thể của hai nút I và J của phần tử như sau:
(2.31)
Tương tự như véc-tơ chuyển vị tại nút de, véc-tơ tải tại nút fecủa phần tử Ωe
trong hệ tọa độ tổng thể XYZ có dạng:
(2.32)
và véc-tơ tải của phần tử trong hệ tọa độ địa phương cũng liên hệ với véc-tơ tải
của phần tử fe trong hệ tọa độ tổng thể thông qua phép biến đổi tọa độ như sau:
(2.33)
b) Các ma trận và véc-tơ phần tử
Từ các công thức (2.21) và (2.33), véc-tơ tải tại fe của phần tử Ωe trong hệ tọa
độ tổng thể XYZ có dạng cụ thể như sau:
20
(2.34)
Ma trận độ cứng Ke của phần tử Ωe trong hệ tọa độ tổng thể XYZ:
(2.35)
và ma trận khối lượng Me của phần tử Ωe trong hệ tọa độ tổng thể XYZ có dạng:
(2.36)
Phương trình cân bằng tĩnh và động
Phương trình cân bằng trong hệ toạ độ tổng thể
(2.37)
trong đó là véc-tơ của các thành phần gia tốc chuyển vị cho tất cả các nút trong
miền bài toán.
21
Phương trình tĩnh
Quá trình thành lập hệ phương trình rời rạc cho bài toán tĩnh học dựa trên công
thức (2.37) tuy nhiên đại lượng , khi đó phương trình cân bằng (2.37) trở
thành:
(2.38)
Khi đó véc-tơ chuyển vị tại nút của phần tử Ωe.
(2.39)
Phương trình dao động tự do
Tương tự như thành lập phương trình cân bằng tĩnh, quá trình thành lập hệ
phương trình rời rạc cho bài toán phân tích dao động tự do dựa trên công thức (2.37)
bằng cách loại bỏ đi véc-tơ ngoại lực fe, ta được hệ phương trình cho bài toán phân
tích dao động tự do (tự nhiên):
(2.40)
Để tìm tần số dao động tự do của hệ, ta thừa nhận nghiệm chuyển vị de có dạng
sau:
(2.41)
trong đó φe là biên độ của chuyển vị nút; ωe là tần số góc, t là thời gian và i là biến
phức sao cho i2 = -1. Thay công thức (2.41) vào (2.40), ta được:
(2.42)
hay
(2.43)
trong đó
(2.44)
Ta thấy, phương trình (2.42) hay (2.43) có một nghiệm là φe = 0 tương ứng
với trường hợp hệ đứng yên, không dao động. Vì vậy để tìm các nghiệm φe ≠ 0 tương
22
ứng với các trường hợp hệ có dao động thật sự, thì định thức của ma trận sau phải
bằng 0:
(2.45)
Khai triển phương trình (2.45) ta được một phương trình đa thức theo có
. bậc Nn. Giải phương trình đa thức này ta sẽ có Nn nghiệm
Tần số góc và tần số dao động tự nhiên liên hệ với nhau qua công thức sau:
(2.46)
hay
(2.47)
Thay các giá trị vào công thức (2.47) ta được các tần số dao động tự nhiên
tương ứng.
(2.48)
2.1.3 Hậu xử lý - tính ứng suất của dàn
Ứng suất Ne của phần tử trong hệ tọa độ tổng thể có dạng:
(2.49)
trong đó Be là ma trận biến dạng – chuyển vị của phần tử Ωe:
(2.50)
2.2. Tối ưu hóa đa mục tiêu dựa trên độ tin cậy cho bài toán kết cấu
Tối ưu hóa đa mục tiêu
Tính toán tối ưu là một lĩnh vực quan trọng của toán học có ảnh hưởng đến
hầu hết các lĩnh vực như kinh tế, kết cấu, giao thông, v.v. Một bài toán tối ưu có thể
gồm một hay nhiều hàm mục tiêu và có hoặc không có các điều kiện ràng buộc.
Nhưng trong thực tế hầu hết các trường hợp ra quyết định cho một vấn đề luôn đòi
hỏi phải có sự hòa hợp giữa hai hay nhiều mục tiêu cùng lúc.
23
Tối ưu hóa gồm hai loại: tối ưu hóa đơn mục tiêu và tối ưu hóa đa mục tiêu.
Các ưu, nhược điểm của hai loại tối ưu được trình bày chi tiết trong Bảng 2.1.
Bảng 2.1: Bảng so sánh ưu nhược điểm giữa tối ưu hóa đơn mục miêu và đa mục
tiêu.
Loại tối ưu Ưu điểm Nhược điểm
- Bài toán tối ưu hoá đơn mục tiêu
- Trong thực tế, một bài toán thường có nhiều hơn một mục là bài toán có một hàm mục tiêu
tiêu để tính toán và thường các cần tối ưu. Lời giải của bài toán tối
ưu đơn mục tiêu là một lời giải duy
mục tiêu này có sự tương quan với nhau, do đó bài toán tối ưu nhất. Do đó chi phí tính toán, thời Tối ưu hóa đơn mục tiêu không đáp ứng gian tính toán thấp. đơn mục tiêu - Thiết lập bài toán đơn giản, đễ
thực hiện. được điều này. - Khi dữ liệu đầu ra thay đổi thì phải thiết lập và giải lại bài toán - Các nhà khoa học nghiên cứu và từ đầu. phát triển sớm, có nhiều thuật toán
áp dụng.
- Bài toán tối ưu hoá đa mục tiêu là
bài toán sẽ có từ hai hàm mục tiêu
trở lên cần tối ưu. Lời giải của bài
toán tối ưu đa mục tiêu là một tập
hợp của các giải pháp tối ưu có
quan hệ đối nghịch giữa các hàm
- Chi phí tính toán lớn, thời gian tính toán lâu. - Thiết lập bài toán phức tạp hơn tối ưu đơn mục tiêu. - Ít giải thuật áp dụng. - Ra đời và phát triển chậm hơn tối ưu đơn mục tiêu. Tối ưu hóa đa mục tiêu. Do đó lời giải tối ưu của
mục tiêu bài toán tối ưu đa mục tiêu không
phải là một lời giải duy nhất như ở
bài toán tối ưu đơn mục tiêu mà sẽ
là một tập hợp của nhiều lời giải tối
ưu.
- Có sự tương quan giữa các mục
tiêu, giúp người thiết kế có cái nhìn
24
Loại tối ưu Ưu điểm Nhược điểm
tổng quát hơn, đúng đắn hơn về bài
toán.
- Dựa vào dữ liệu đầu ra, kỹ sư
thiết kế có thể đưa ra dạng tổng
quát của bài toán, từ đó chạy ra một
tập hợp các kết quả. Từ đó dựa vào
từng bài toán cụ thể mà lựa chọn
một nghiệm tối ưu thích hợp trong
tập hợp nghiệm.
Bảng 2.1 cho thấy, nhược điểm lớn nhất của tối ưu hóa đa mục tiêu so với tối
ưu đơn mục tiêu chính là chi phí tính toán lớn, thời gian tính toán lâu. Tuy nhiên với
sự phát triển không ngừng của khoa học máy tính, nhiều phiên bản máy tính mới có
tốc độ xử lý cao đã ra đời đặc biệt là các hệ siêu máy tính đã làm cho việc giải quyết
các bài toán lớn có chi phí tính toán cao trở nên dễ dàng hơn. Vì vậy, với những ưu
điểm nổi bật của mình, và các nhược điểm dần được khắc phục, tối ưu hóa đa mục
tiêu ngày càng được sự quan tâm của rất nhiều nhà khoa học trên thế giới và được
phát triển rất mạnh. Do đó, việc áp dụng tối ưu hóa đa mục tiêu để tính toán cho kết
cấu là thiết thực và mang lại nhiều lợi ích.
Như đã trình bày ở trên thì bài toán tối ưu hoá đa mục tiêu là bài toán sẽ có từ
hai hàm mục tiêu trở lên cần tối ưu. Lời giải của bài toán tối ưu đa mục tiêu sẽ là một
tập hợp của các giải pháp tối ưu có quan hệ đối nghịch giữa các hàm mục tiêu, do đó
công thức tổng quát của một bài toán tối ưu hóa đa mục tiêu biến tiền định như sau:
Chịu ràng buộc (2.51)
trong đó x là véc-tơ chứa các biến thiết kế; lần lượt là các hàm
ràng buộc dạng bất đẳng thức và đẳng thức liên quan đến ứng xử của kết cấu hay các
ràng buộc thiết kế; M là số lượng hàm mục tiêu; J là số lượng ràng buộc bất đẳng
25
thức; K là số lượng ràng buộc đẳng thức; là cận dưới và cận trên của biến
là các hàm mục tiêu. Khi M = 1 thì bài thiết kế xi; n là số lượng biến thiết kế;
toán trở thành tối ưu hóa đơn một tiêu.
Mục tiêu của bài toán là tìm kiếm các giá trị của biến thiết kế trong không gian
thiết kế sao cho cực tiểu hoặc cực đại hóa các hàm mục tiêu nhưng vẫn thỏa
mãn các điều kiện ràng buộc và Một điều cần lưu ý đối với loại bài toán
này là các giá trị thiết kế cũng như các tham số của mô hình bài toán là những giá trị
cố định và không thay đổi trong quá trình thiết kế hay sử dụng. Kết quả một nghiệm
tối ưu bất kỳ trên đường Pareto của bài toán tối ưu hóa đa mục tiêu tiền định (2.51)
có thể được biểu diễn bởi như Hình 2.5.
Hình 2.5: Mô tả kết quả của một bài toán tối ưu hóa, với x* là lời giải tối ưu.
Kết quả bài toán tối ưu hóa đa mục tiêu tiền định là tập hợp một chuỗi nghiệm,
trong đó tất cả các nghiệm này đều phải đảm bảo an toàn cho kết cấu hay thỏa mãn
các điều kiện ràng buộc khi các thông số cũng như các giá trị kết quả tối ưu của mô
hình bài toán là cố định. Vì vậy, khi có sự thay đổi ngẫu nhiên của các thông số trong
mô hình bài toán, khả năng xảy ra phá hủy của kết cấu có thể xấp xỉ gần 50% như
được mô tả ở Hình 2.6.
26
Hình 2.6: Mô tả sự ảnh hưởng của các yếu tố ngẫu nhiên lên kết quả tối ưu hoá.
Bài toán tối ưu hóa đa mục tiêu dựa trên độ tin cậy
Khác với bài toán tối ưu hóa đa mục tiêu tiền định (bài toán tối ưu không xét
đến ảnh hưởng của các yếu tố ngẫu nhiên), bài toán tối ưu hóa đa mục tiêu dựa trên
độ tin cậy RBDO (Reliability Based-Design Optimization) có kể đến sự ảnh hưởng
của các yếu tố ngẫu nhiên trong mô hình bài toán và được biểu diễn như sau:
(2.52) Biến thiết kế:
Chịu ràng buộc:
trong đó là hàm mục tiêu của bài toán; M là số lượng hàm mục tiêu của
bài toán; L là số lượng hàm ràng buộc cả bài toán; d là véc-tơ chứa các biến thiết kế
không phải là biến ngẫu nhiên; x là véc-tơ chứa các giá trị vừa là biến thiết kế vừa là
biến ngẫu nhiên; p là véc-tơ chứa các tham số ngẫu nhiên của mô hình bài toán;
là các ràng buộc độ tin cậy hay xác suất an toàn của kết cấu khi
chịu ảnh hưởng của các yếu tố ngẫu nhiên; là phần trăm độ tin cậy yêu
cầu và là chỉ số độ tin cậy yêu cầu.
Tương tự như bài toán tối ưu hóa đa mục tiêu tiền định, mục tiêu của bài toán
tối ưu hóa đa mục tiêu dựa trên độ tin cậy cũng đi tìm các giá trị của biến thiết kế
trong không gian thiết kế sao cho cực tiểu hoặc cực đại các hàm mục tiêu
và thỏa mãn các điều kiện ràng buộc về độ tin cậy, ràng buộc đẳng thức và ràng buộc
27
bất đẳng thức và Kết quả tối ưu của bài toán tối ưu hóa đa mục tiêu dựa
trên độ tin cậy (2.52) cũng là một hợp các nghiệm nằm trên đường Pareto. Có thể
biểu diễn một nghiệm bất kỳ trên đường Pareto như Hình 2.7.
Hình 2.7: Mô tả kết quả của một bài toán tối ưu hóa dựa trên độ tin cậy, với x* là
lời giải tối ưu.
Kết quả Hình 2.7 cho thấy ràng buộc xác suất sẽ giúp đưa giá trị thiết kế vào
vùng an toàn với phần trăm độ tin cậy cho trước ứng với phần diện tích màu xanh bên
trên hàm ràng buộc , và phần diện tích màu tím bên dưới ràng buộc
là xác suất phá hủy của kết cấu khi chịu ảnh hưởng ngẫu nhiên của hai
biến thế kế x1 và x2.
Tối ưu hóa đa mục tiêu dựa trên độ tin cậy cho kết cấu dàn
Ứng dụng các lý thuyết trên, bài toán tối ưu hóa trong luận văn này được thành
lập với các thành phần sau đây: hàm mục tiêu là cực tiểu khối lượng và cực tiểu
chuyển vị toàn hệ kết cấu dàn. Mô hình bài toán được xây dựng như Hình 2.8.
28
Bài toán tối ưu hoá tiền định Bài toán tối ưu hoá dựa trên độ tin cậy
Hình 2.8: Thiết lập bài toán tối ưu hóa đa mục tiêu kết cấu dàn.
trong đó weight(A) là khối lượng của toàn bộ kết cấu dàn; li là chiều dài tương ứng
của mỗi thanh; ρ là khối lượng riêng; e, n, nc lần lượt là tổng số thanh, tổng số nút và
b là ràng buộc ổn định của thanh thứ i khi nó chịu nén; P
tổng số thanh chịu nén trong hệ kết cấu; ui; σi lần lượt là chuyển vị nút và ứng suất
phần tử trong hệ kết cấu; σi
là tải trọng tác dụng lên hệ kết cấu; E là mô-đun đàn hồi; A là diện tích các thanh;
disp(A) là chuyển vị; f1, f2, f3 lần lượt là tần số dao động thứ nhất, thứ 2 và thứ 3 của
toàn bộ hệ kế cấu; K là ma trận độ cứng tổng thể kết cấu; là phần trăm độ tin
cậy yêu cầu và là chỉ số độ tin cậy yêu cầu; W0 là khối lượng cộng thêm tại các
nút của hệ kết cấu; F là véc-tơ ngoại lực tác động lên các nút của kết cấu.
2.3. Giải thuật tối ưu hóa dựa trên độ tin cậy một vòng lặp đơn xác định
Vì quá trình giải bài toán tối ưu hóa dựa trên độ tin cậy liên quan đến hai vòng
lặp gồm một vòng lặp tối ưu hóa và một vòng lặp đánh giá độ tin cậy, nên chi phí tính
toán cho việc giải một bài toán tối ưu dựa trên độ tin cậy là tương đối lớn so với bài
toán tối ưu hóa không xét đến độ tin cậy. Thông thường, có 3 cách tiếp cận chính để
giải một bài toán tối ưu hóa dựa trên độ tin cậy bao gồm [7]:
a) Giải thuật vòng lặp đôi (Double Loop Method), ở giải thuật này, vòng lặp
đánh giá độ tin cậy được lồng vào trong vòng lặp tối ưu hóa như được minh họa ở
Hình 2.9. Do việc kiểm tra độ tin cậy được thực hiện ở mỗi bước lặp trong quá trình
29
tìm kiếm lời giải tối ưu, nên phương pháp này thường cho kết quả chính xác. Tuy
nhiên, do việc phân tích độ tin cậy lặp lại nhiều lần trong quá trình tối ưu hóa nên chi
phí tính toán sẽ rất cao đặc biệt là khi bài toán RBDO có nhiều ràng buộc độ tin cậy.
Vì vậy, nó không hiệu quả khi áp dụng cho các bài toán lớn trong kỹ thuật. Chi tiết
của giải thuật này có thể tìm thấy ở tài liệu tham khảo [8].
Hình 2.9: Sơ đồ giải thuật vòng lặp đôi.
b) Giải thuật vòng lặp tách rời hay còn gọi là giải thuật lặp tuần tự (Decouple
Double Loop Method). Ở giải thuật này, vòng lặp tối ưu hóa và vòng lặp đánh giá độ
tin cậy được tách rời nhau như được minh họa ở Hình 2.10. Do vòng lặp tối ưu và
vòng lặp phân tích độ tin cậy thực hiện một cách độc lập nên chi phí tính toán sẽ giảm
đi đáng kể so với giải thuật vòng lặp đôi [44]. Tuy nhiên, vì độ tin cậy của bài toán
chỉ được kiểm tra sau khi nghiệm bài toán tối ưu tìm được nên độ chính xác của bài
toán sẽ thấp hơn so với giải thuật vòng lặp đôi. Chi tiết của giải thuật này có thể được
tìm thấy ở tài liệu tham khảo [44].
30
Hình 2.10: Sơ đồ giải thuật vòng lặp tuần tự.
c) Giải thuật vòng lặp đơn (Single Loop Method), ở giải thuật này, bài toán tối
ưu hóa dựa trên độ tin cậy sẽ được chuyển thành bài toán tối ưu hóa tiền định thông
qua việc chuyển đổi các ràng buộc xác suất của bài toán RBDO thành các ràng buộc
tiền định trước khi thực hiện việc giải bài toán tối ưu hóa như được minh hoạ ở Hình
2.11. Do đó, sau khi thực hiện việc chuyển các ràng buộc xác suất thành các ràng
buộc tiền định, việc giải bài toán RBDO lúc này chỉ là giải bài toán tối ưu hóa tiền
định nên chi phí tính toán sẽ giảm đi đáng kể so với hai cách tiếp cận trước đó. Tuy
nhiên, do việc chuyển ràng buộc xác suất sang ràng buộc tiền định là một bước xấp
xỉ nên kết quả tối ưu đạt được của bài toán RBDO cũng có độ chính xác không cao
so với cách tiếp cận giải thuật vòng lặp đôi. Mặc dù, kết quả có độ chính xác thấp hơn
giải thuật vòng lặp đôi, nhưng xét về chi phí tính toán và khả năng áp dụng để giải
các bài toán RBDO có chi phí lớn trong thực tế thì giải thuật này được xem là cách
tiếp cận hiệu quả nhất so với các giải thuật hiện tại [8]. Trong luận văn này, giải thuật
SLDM thuộc nhóm này được tác giả sử dụng. Chi tiết của giải thuật này sẽ được trình
bày ở các mục tiếp theo.
31
Hình 2.11: Sơ đồ giải thuật một vòng lặp đơn xác định.
Chuyển đổi không gian thiết kế xác suất sang không gian thiết kế tiền định
Ở bước này, không gian thiết kế của bài toán tối ưu được tạo thành từ các ràng
buộc xác suất sẽ được chuyển sang không gian thiết kế tiền định bởi việc chuyển đổi
các ràng buộc xác suất sang các ràng buộc tiền định. Việc chuyển đổi được thực hiện
bằng cách dịch chuyển các hàm trạng thái giới hạn trong các ràng buộc xác suất một
khoảng cách như ở Hình 2.12, ở đó là chỉ số độ tin cậy được xác định từ phần
trăm xác suất an toàn cho trước của bài toán thông qua công thức , với
là xác suất an toàn. Hình 2.12, đường cong nét liền màu xanh đại diện cho hàm
trạng thái giới hạn ban đầu trong ràng buộc xác suất , trong đó u là véc-tơ
chứa các tham số ngẫu nhiên và biến thiết kế ngẫu nhiên trong không gian chuẩn hóa;
đường cong nét đứt màu đỏ đại diện cho biên của ràng buộc tiền định đã được chuyển
đổi ; và miền được gạch bởi các đường nét liền màu xanh là miền khả
thi xác định. Việc chuyển đổi này nhằm đảm bảo khoảng cách nhỏ nhất giữa hai điểm
bất kỳ trên đường cong nét đứt màu đỏ và đường cong nét liền màu xanh luôn là .
Điều này cũng đảm bảo rằng kết quả tối ưu tìm được trên miền khả thi tiền định sẽ
thỏa mãn được các ràng buộc xác suất ban đầu của bài toán.
32
Hình 2.12: Minh họa không gian thiết kế của bài toán.
Thành lập ràng buộc tiền định
Đặt lần lượt là hàm trạng thái giới hạn và ràng buộc xác
và suất tiền định. Giả sử là một điểm bất kỳ thuộc , khi đó điểm có mật độ
xác suất cao nhất (MPP-Most Probable Point) trên hàm tương ứng với
điểm có thể đạt được bằng cách dịch chuyển điểm đến hàm một
khoảng (dịch chuyển vuông góc). Mối quan hệ giữa và điểm tương ứng
của nó được mô tả qua công thức sau:
(2.53)
trong đó n là véc-tơ đạo hàm chuẩn hóa được đánh giá ở trên hàm được
xác định như sau:
(2.54)
Theo Li và cộng sự [7], miền khả thi của bài toán có thể được viết lại như sau:
(2.55)
33
Như vậy, sau khi ràng buộc tiền định của bài toán được thành lập, bài toán
RBDO lúc này sẽ trở thành bài toán tối ưu hóa tiền định. Bài toán này được giải bằng
phương pháp NSGA-II.
2.4. Giải thuật di truyền phân loại không trội NSGA-II
Giải thuật NSGA-II được hình thành dựa trên khái niệm về tiến hóa trong tự
nhiên. Tính tối ưu của quá trình tiến hoá thể hiện ở chỗ thế hệ sau sẽ chứa nhiều đặc
tính tốt hơn so với thế hệ trước (phát triển hơn, hoàn thiện hơn và phù hợp với môi
trường hơn). Trong quá trình tiến hoá, các thế hệ mới được sinh ra để bổ sung và thay
thế cho thế hệ cũ. Cá thể nào phát triển hơn và thích ứng hơn với môi trường sẽ tồn
tại, cá thể nào kém thích ứng hơn sẽ bị đào thải.
Các cá thể mới sinh ra trong quá trình tiến hoá nhờ sự lai ghép ở thế hệ cha
mẹ. Một cá thể mới có thể mang những thuộc tính của cha mẹ (di truyền), cũng có
thể mang thuộc tính hoàn toán mới (đột biến). Nhìn chung các thuật toán tiến hoá,
tuy có những điểm khác biệt, nhưng đều mô phỏng bốn quá trình tiến hoá cơ bản: lai
ghép, đột biến, sinh sản, lựa chọn.
Trên cơ sở và ý tưởng từ GA (Genetic Agorithm – Giải thuật di truyền),
Srinivas và Deb đã đề xuất phương pháp NSGA (Non-Dominated Sorting Genetic
Algorithm) vào năm 1994 [45], phương pháp này chỉ khác phương pháp GA ở bước
lựa chọn Selection). Tuy nhiên, phương pháp NSGA vẫn tồn tại những hạn chế như:
thời gian tính toán chậm, phụ thuộc tham số điều khiển, v.v. Để khắc phục những hạn
chế của NSGA, Deb và cộng sự đã đề xuất phương pháp NSGA-II vào năm 2002 [2].
Phương pháp này không những khắc phục được những hạn chế của NSGA mà còn
đảm bảo sự đa dạng và duy trì được các cá thể tốt qua các thế hệ. Do đó, luận văn này
đề xuất lựa chọn phương pháp NSGA-II để tìm nghiệm tối ưu cho các bài toán cần
khảo sát.
Một số định nghĩa và thuật toán sử dụng trong phương pháp NSGA–II
Để thuận lợi hơn khi tiếp cận phương pháp NSGA – II, ta cần biết một số khái
niệm và định nghĩa sau:
34
Khái niệm về nghiệm Pareto
Bài toán tối ưu đa mục tiêu có nghiệm là một chuỗi nghiệm và tập hợp nghiệm
này gọi là nghiệm Pareto.
Khái niệm về sự trội (domination)
Hầu hết các thuật toán tối ưu đa mục tiêu đều sử dụng khái niệm về sự trội. nó
dùng để so sánh hai nghiệm (cá thể) với nhau.
Định nghĩa: Một nghiệm x(1) được gọi là trội so với nghiệm x(2), nếu thỏa điều
kiện (a) và (b):
a) Nghiệm x(1) không xấu hơn nghiệm x(2) ở tất cả các giá trị hàm mục tiêu, kí
hiệu là: với j = 1, 2,..., M.
b) Nghiệm x(1) phải tốt hơn nghiệm x(2) ở ít nhất một mục tiêu, kí hiệu là:
với ít nhất một .
Nếu bất kì các điều kiện ở trên bị vi phạm, nghiệm x(1) không trội so với nghiệm
x(2).
x(2)): Một số cách diễn đạt nếu x(1) trội so với nghiệm x(2) (kí hiệu x(1)
• x(2) bị trội bởi x(1); • x(2) không bị xấu hơn x(1).
Định nghĩa tập không trội (Non-dominated set)
Tập không trội P’ là tập con của tập P và nó chứa các cá thể (nghiệm) không
bị trội bởi bất kỳ cá thể nào của tập P.
Phương pháp Naive và Slow
Phương pháp này được ứng dụng để tìm tập (biên) không bị trội (non-
dominated set). Trong đó, mỗi cá thể i được so sánh với tất cả các cá thể khác trong
quần thể. Nếu cá thể i bị trội bởi một cá thể bất kỳ trong quần thể (nghĩa là, có một
cá thể nào đó tốt hơn cá thể i) thì tồn tại ít nhất một cá thể trong quần thể tốt hơn i ở
tất cả các hàm mục tiêu. Vậy, cá thể i nằm ngoài tập không bị trội. Ngược lại, nếu
không có cá thể nào trong quần thể trội hơn (tốt hơn) i thì i thuộc về tập không bị trội.
35
Áp dụng phương pháp này cho các nghiệm khác trong bộ dân số chúng ta sẽ tìm được tập không trội kí hiệu là P’ trong quần thể P có thước N
Thuật toán tìm tập không bị trội (non-dominated set)
Bước 1: Khởi tạo biến đếm của cá thể i =1 và tập chứa các cá thể không bị trội
Bước 2: Xét một cá thể bất kì j ϵ P (với j ≠ i), kiểm tra tính trội, nếu cá thể j
trội hơn cá thể i. Đến bước 4.
Bước 3: Nếu có nhiều cá thể ở trong quần thể P ta tăng biến j lên 1, j = j +1
trở lại bước 2; ngược lại thêm i vào tập không bị trội P’, .
Bước 4: Tăng biến i thêm 1, i =i +1. Khi i ≤ N, trở lại bước 2; ngược lại, dừng
và xuất ra tập không bị trội P’.
Phương pháp sắp xếp tập không bị trội của quần thể (Non-dominated
Sorting)
Hầu hết các giải thuật tiến hóa để giải quyết các bài toán tối ưu đa mục tiêu
chỉ yêu cầu tìm được biên không bị trội tốt nhất trong quần thể (biên không bị trội
thứ nhất). Các thuật toán này phân loại quần thể thành tập chứa các cá thể không bị
trội và các cá thể bị trội. Tuy nhiên, cũng có một số giải thuật yêu cầu sắp xếp toàn
bộ các cá thể của quần thể vào các tập không bị trội.
Thuật toán phân loại không bị trội (Non-dominated Sorting)
Bước 1: Tạo tập Pj (j = 1,2,...) để chứa tất cả các tập không bị trội và Pj là rỗng.
Khởi tạo biến đếm của tập không bị trội j = 1.
Bước 2: Sử dụng phương pháp Naive and Slow tìm tập không bị trội P’ của
quần thể P.
Bước 3: Cập nhật Pj = P’ và loại P’ ra khỏi tập P, P=P \ P’.
Bước 4: Nếu P ≠ tăng biến j thêm 1 tức là j = j +1, và quay lại bước 2.
Ngược lại, dừng vòng lặp và xuất ra các tập chứa các cá thể không bị trội Pi, (i=1,
2,..., j).
36
Thuật toán phân loại không bị trội nhanh
Các bước của thuật toán phân loại không bị trội nhanh (Fast Non-dominated
Sort):
. Với mọi j ≠ i, j ϵ P, thực hiện bước Bước 1: Với mỗi p ϵ P, np = 0 và Sp =
2 và bước 3.
Bước 2: Nếu thì . Ngược lại, nếu thì .
Bước 3: Nếu thì biên không bị trội thứ nhất , .
Bước 4: Trong khi thì thực hiện các bước tiếp theo.
Bước 5: Khởi tạo , sắp xếp tiếp các nghiệm không bị trội. Với mỗi
và .
Bước 5a: tăng .
Bước 5b: nếu thì .
Bước 6: Cập nhật và . Trở lại bước 4.
trong đó Sp là tập hợp các nghiệm mà nghiệm p trội hơn; np là số lượng các nghiệm
trội hơn nghiệm p; nq là số lượng các nghiệm trội hơn nghiệm q; F1 là biên không bị
trội thứ nhất; Fi là biên không bị trội hiện hành; H là biên không bị trội hiện hành tiếp
theo nếu nq = 0.
Thuật toán xác định giá trị vùng an toàn
Để ước lượng mật độ cá thể quanh một cá thể, lấy trung bình khoảng cách hai
điểm gần nhất dọc theo mỗi mục tiêu như Hình 2.13. Đại lượng idistance dùng để tính
kích thước hình chữ nhật bao quanh điểm i. Hình chữ nhật này gọi là vùng an toàn
(Crowding Distance).
37
f2 1
vùng an toàn
i-1
i i+1 f1 l
Hình 2.13: Thuật toán xác định vùng an toàn (crowding-distance-assignment).
Bước 1: Gọi số cá thể trong mặt I là l = |I|, với mỗi cá thể i trong tập hợp, khởi
tạo I [i]distance = 0.
Bước 2: Với mỗi mục tiêu m xếp hạng giá trị của hàm mục tiêu từ thấp đến
cao: I = sort(I, m).
Bước 3: Gán khoảng cách cho hai cá thể ở biên (cá thể 1 và cá thể l ) I [1]distance
= I [l]distance = ∞, các khoảng cách của các cá thể từ cá thể i = 2 đến (l - 1) được tính
theo công thức .
trong đó : I [i]distance là giá trị khoảng cách của cá thể thứ i trong tập hợp I; I [i].m là
giá trị hàm mục tiêu thứ m của cá thể thứ i trong tập hợp I; là giá trị lớn nhất
quần thể của hàm mục tiêu thứ m; là giá trị nhỏ nhất quần thể của hàm mục tiêu
thứ m.
Thuật toán so sánh kích thước vùng an toàn ( ≥ n)
Thuật toán này cần có hai thuộc tính:
Hạng không bị trội trong quần thể (rank: irank);
Vùng an toàn trong quần thể (Crowding Distance: idistance).
Cá thể i trội hơn (tốt hơn, chiến thắng) cá thể j (ký hiệu i ≥ n j) nếu thỏa hai
điều kiện sau:
Khi hạng của cá thể i nhỏ hơn hạng của cá thể j (irank < jrank);
Khi hạng của cá thể i bằng j thì khoảng cách an toàn của cá thể i phải lớn
hơn khoảng cách an toàn của cá thể j, (irank = jrank) và (idistance > jdistance).
38
Chọn lọc
Quá trình chọn lọc trực tiếp giữa hai cá thể (Crowded Binary Tournament
Selection) được thực hiện bằng thuật toán so sánh kích thước vùng an toàn.
Lai ghép
Sau khi lựa chọn được hai cá thể cha mẹ có đặc tính tốt, áp dụng phép lai trung
gian để tạo ra hai cá thể con. Mục đích của quá trình này là duy trì những đặc tính tốt
của thế hệ cha mẹ.
trong đó child1 là cá thể con 1; child2 là cá thể con 2; parent1 là cha mẹ 1; parent2
là cha mẹ 2; rand là hàm lấy một giá trị ngẫu nhiên; ratio là một số vô hướng. Nếu
ratio được lấy ngẫu nhiên trong khoảng [0, 1] thì cá thể con được tạo ra nằm trong
khoảng giữa của cha mẹ. Nếu muốn giải thuật kết thúc sớm (premature) ta nên chọn
ratio lớn hơn 1.
Đột biến
Đột biến là quá trình tạo ra ở cá thể con những tính trạng không có ở thế hệ
cha mẹ. Mục đích của quá trình này là tăng sự đa dạng của các cá thể con. Giải thuật
đột biến phổ biến nhất là cộng thêm vào các biến một số ngẫu nhiên phân bố đều (đột
biến Gaussian - Gaussian mutation).
trong đó child là cá thể con; parent là cha mẹ; randn là hàm lấy một giá trị ngẫu
nhiên; ub là cận trên của biến thiết kế; lb là cận dưới của biến thiết kế; scale là tham
số xác định độ lệch chuẩn của các số ngẫu nhiên; shrink là một số vô hướng, thuộc
khoảng [0, 1]; currGen là thế hệ hiện hành; maxGen là thế hệ lớn nhất.
Xử lý ràng buộc
Ràng buộc được xử lý bằng quá trình chọn lọc trực tiếp (Crowded Binary
Tournament Selection) giữa hai cá thể của quần thể, cá thể (nghiệm) nào tốt hơn sẽ
được chọn. Đối với bài toán có ràng buộc, nghiệm bài toán có thể nằm trong vùng
39
khả thi hoặc không khả thi. Như vậy, khi so sánh sẽ xảy ra ba trường hợp: (i) cả hai
nghiệm đều nằm trong vùng khả thi, (ii) một nghiệm khả thi một nghiệm không, (iii)
cả hai nghiệm không nằm trong vùng khả thi.
Nếu là bài toán đơn mục tiêu, chọn như sau cho mỗi trường hợp:
Trường hợp (i): chọn nghiệm có giá trị hàm mục tiêu tốt hơn;
Trường hợp (ii): chọn nghiệm khả thi;
Trường hợp (iii): chọn nghiệm có tổng vi phạm ràng buộc nhỏ hơn.
Bài toán tối ưu đa mục tiêu, trường hợp (ii) và (iii) có thể chọn dễ dàng, nhưng
gặp khó khăn ở trường hợp (i) vì có nhiều hàm mục tiêu. Nghĩa là, khi hai nghiệm
đều nằm trong vùng khả thi, chọn nghiệm nằm trên biên không trội có thứ hạng tốt
hơn. Trường hợp nằm cùng một biên, để duy trì các cá thể tốt trên biên này, sử dụng
khoảng cách an toàn để chọn ra những cá thể có giá trị khoảng cách an toàn lớn hơn.
Thuật toán
n x(j)), khi thỏa một trong các điều kiện sau:
Một nghiệm x(i) được gọi là ‘trội – ràng buộc’ so với nghiệm x(j) (kí hiệu x(i) ≥
• Nghiệm x(i) là nghiệm khả thi còn nghiệm x(j) thì không khả thi.
• Nghiệm x(i) và x(j) đều không khả thi, nhưng nghiệm x(i) có tổng vi phạm
ràng buộc nhỏ hơn.
• Nghiệm x(i) và x(j) đều khả thi, nghiệm x(i) không xấu hơn x(j) ở tất cả các giá
trị hàm mục tiêu và nghiệm x(i) phải tốt hơn nghiệm x(j) ở ít nhất một mục tiêu.
Vòng lặp chính của phương pháp NSGA-II
Quần thể ban đầu Pt có N cá thể và được chọn ngẫu nhiên. Quần thể Qt có N
cá thể và được tạo từ Pt. Để đảm bảo sự chọn lọc tự nhiên (Elitism) quần thể Rt được
tổ hợp từ quần thể Pt và Qt và có 2N cá thể. Phân loại nhanh các cá thể của Rt theo tiêu chí không bị trội (Fast Non-dominated Sort). Lớp F1 (là tập cá thể không bị trội
tốt nhất – Best Non-dominated) chứa các cá thể tốt nhất trong quần thể. Tiếp tục xếp
hạng các cá thể còn lại của Rt theo tiêu chí không bị trội được lớp F2, cứ thế tiếp tục được F3, v.v.
40
Nếu số lượng cá thể trong F1 bé hơn N cá thể. Thì thế hệ Pt+1 chứa tất cả các
cá thể của F1. Tuy nhiên, quần Pt+1 cần phải có đủ N cá thể. Nên số lượng cá thể còn
thiếu sẽ được bổ sung từ lớp F2, F3, v.v. bằng cách ứng dụng thuật toán so sánh kích
thước vùng an toàn (Crowded Comparison Operator: ≥ n) để sắp xếp các giá trị vùng
an toàn (Crowding Distance Value).
Thế hệ con Qt+1 được tạo ra từ Pt+1 bằng cách thực hiện bước chọn lọc trực
tiếp giữa hai cá thể (Crowded Binary Tournament Selection), lai chéo, đột biến như
trong Hình 2.14. Trong đó, áp dụng thuật toán so sánh kích thước vùng an toàn để
thực hiện bước chọn lọc trực tiếp giữa hai cá thể, lai chéo và đột biến.
Nếu số lượng cá thể trong F1 là N cá thể. Thì thế hệ Pt+1 là tất cả các cá thể của
F1, không bổ sung thêm cá thể từ lớp F2, F3, v.v. nữa. Như vậy số lần tính toán sẽ
giảm đáng kể.
Hình 2.14: Sơ đồ giải thuật NSGA-II.
Các bước của vòng lặp chính
Bước 1: Quần thể ban đầu Pt được chọn ngẫu nhiên, quần thể con Qt được tạo
ra từ quần thể Pt, kết hợp quần thể con và quần thể cha mẹ, Rt = Pt Qt. Tiến hành
phân loại không bị trội nhanh Rt để tạo được các biên không bị trội Fi với i = 1, 2,...
Bước này được thực hiện bằng thuật toán phân loại không bị trội nhanh.
Bước 2: Khởi tạo một quần thể mới Pt+1 = , cài đặt biến đếm i = 1. Khi
|Pt+1|+|Fi| 41 Bước 3: Dùng giá trị vùng an toàn để so sánh kích thước của chúng cho từng biên không bị trội và chọn ra N cá thể đầu tiên, Pt+1 = Pt+1 [0:N]. Bước này được thực hiện bằng hai thuật toán là xác định giá trị vùng an toàn và so sánh kích thước vùng an toàn. Bước 4: Tạo thế hệ con Qt+1 từ quần thể Pt+1 thông qua quá trình chọn lọc trực tiếp giữa hai cá thể, quá trình lai chéo và quá trình đột biến. Các bước của vòng lặp chính có thể tham khảo Hình 2.15. Quần thể ban đầu
Qt và Pt Phân loại không
trội nhanh 1
-
n Quần thể mới p
ặ
l g
n
ò
V Xác định vùng an toàn
So sánh kích thước vùng an
toàn Vùng an toàn Chọn lọc
Lai chéo
Đột biến Thế hệ con
Qt+1 và Pt+1 Vòng lặp n Xuất ra kết quả Hình 2.15: Các bước của vòng lặp chính giải thuật NSGA-II. 2.5. Giải thuật một vòng lặp đơn RBDO sử dụng phương pháp NSGA-II Qui trình để giải một bài toán tối ưu hóa dựa trên độ tin cậy sử dụng giải thuật một vòng lặp đơn xác định kết hợp với phương pháp NSGA-II (RBDO-NSGA-II) được mô tả như Hình 2.16. - Bước 1: Phân tích kết cấu. 42 - Bước 2: Dựa vào các nhu cầu thiết kế, thành lập bài toán tối ưu dựa trên độ tin cậy với mức độ tin cậy được xác định trước. - Bước 3: Chuyển đổi các ràng buộc xác suất thành các ràng buộc tiền định. - Bước 4: Áp dụng giải thuật NSGA-II để tìm kết quả tối ưu cho bài toán. Sơ đồ giải thuật được trình bày như sau: Hình 2.16: Sơ đồ giải thuật một vòng lặp đơn sử dụng giải thuật NSGA-II (RBDO-NSGA-II). 43 Trong chương này, các ví dụ số sẽ được thực hiện nhằm làm rõ ưu điểm của bài toán tối ưu hóa đa mục tiêu, những lợi ích khi xét đến độ tin cậy trong thiết kế kết cấu, và hiệu quả của các phương pháp đề xuất NSGA-II và RBDO-NSGA-II trong luận văn. Các ví dụ số được chia làm hai phần: Phần thứ nhất là xác minh tính chính xác và độ tin cậy của code lập trình Matlab, code lập trình giải thuật NSGA-II và code lập trình giải thuật RBDO-NSGA- II. Phần thứ hai là áp dụng phương pháp NSGA-II và RBDO-NSGA-II giải các bài toán tối ưu hóa đa mục tiêu kết cấu dàn. Các bước thực hiện đối với phần thứ hai như sau: Thiết lập và giải bài toán tối ưu hóa kết cấu dàn với hàm mục tiêu gồm cực tiểu khối lượng và cực tiểu chuyển vị chịu ràng buộc chuyển vị, ứng suất và tần số dao động. + Thiết lập và giải bài toán tối ưu hóa kết cấu dàn với hàm mục tiêu gồm cực tiểu khối lượng và cực tiểu chuyển vị chịu ràng buộc chuyển vị, ứng suất và tần số dao động có xét thêm ràng buộc về độ tin cậy. Hiệu quả và độ tin cậy của các phương pháp đề xuất sẽ được đánh giá khi so sánh kết quả hiện tại với kết quả của các phương pháp khác đã được nghiên cứu trước đây. Một số kết quả tối ưu thu được từ NSGA-II và RBDO-NSGA-II sẽ được đánh giá độ tin cậy bởi phương pháp đánh giá độ tin cậy bậc nhất FORM (First Order Reliability Method). 3.1. Kiểm chứng code lập trình Matlab Trong mục này, phương pháp phần tử hữu hạn được áp dụng để phân tích ứng xử cho kết cấu dàn (chuyển vị, ứng suất và tần số dao động). Kết quả phân tích đạt được của phương pháp phần tử hữu hạn FEM (Finite Element Method), code lập trình 44 NSGA-II và code lập trình RBDO-NSGA-II được so sánh với kết quả của các nghiên cứu trước đây. Ví dụ 1 – Phân tích chuyển vị và ứng suất của dàn 11 thanh Trong ví dụ này, phương pháp phần tử hữu hạn sẽ được áp dụng để phân tính chuyển vị cho các nút và ứng suất trong các thanh cho kết cấu dàn phẳng 11 thanh như trong Hình 3.1. Các thông số bài toán: E = 70 GPa, A =3*10-4 m2. Bài toán này được thực hiện nhằm mục đích xác minh tính chính xác và độ tin cậy của phương pháp phần tử hữu hạn FEM (Finite Elemesnt Method) tính toán chuyển vị và ứng của kết cấu dàn. Kết quả phân tích của phương pháp phần tử hữu hạn FEM (Finite Elemesnt Method) được so sánh với kết quả nghiên cứu của Ferreira [46]. Hình 3.1: Mô hình bài toán dàn phẳng 11 thanh. Kết quả phân tích của phương pháp phần tử hữu hạn FEM (Finite Elemesnt Method) và kết quả nghiên cứu của Ferreira [46] được trình bày ở Bảng 3.1 và Bảng 3.2. Bảng 3.1: Kết quả chuyển vị trong các nút. Ferreira [46] FEM Sai số (%) Nút 1
2
3
4
5
6 Chuyển vị
phương x
0
7.1429
5.2471
5.2471
10.4942
3.3513 Chuyển vị
phương y
0
-9.0386
-16.2965
-20.0881
0
-9.0386 Chuyển vị
phương x
0
7.1429
5.2471
5.2471
10.4942
3.3513 Chuyển vị
phương y
0
-9.0386
-16.2965
-20.0881
0
-9.0386 Chuyển vị
phương x
0.00
0.00
0.00
0.00
0.00
0.00 Chuyển vị
phương y
0.00
0.00
0.00
0.00
0.00
0.00 45 Bảng 3.2: Kết quả ứng suất trong các thanh. Ứng suất (σ) Thanh 1 Ferreira [46]
-210.9015 FEM
-210.9015 Sai số
(%)
0.00 2 122.4318 122.4318 0.00 3 62.5575 62.5575 0.00 4 -44.2349 -44.2349 0.00 5 -173.1447 -173.1447 0.00 6 -88.4697 -88.4697 0.00 7 62.5575 62.5575 0.00 8 -173.1447 -173.1447 0.00 9 -44.2349 -44.2349 0.00 10 122.4318 122.4318 0.00 11 -210.9015 -210.9015 0.00 Nhận xét: Bảng 3.1 và Bảng 3.2 cho thấy, kết quả chuyển vị và ứng suất được tính bằng phương pháp số giống kết quả nghiên cứu của Ferreira [46]. Điều đó có nghĩa là ứng suất trong các thanh dàn và chuyển vị trong các nút tính toán bởi Matlab bằng với kết quả nghiên cứu của Ferreira [46]. Do đó code phần tử hữu hạn để tính toán chuyển vị và ứng suất cho kết cấu dàn sử dụng trong luận văn là chính xác và có độ tin cậy cao. Ví dụ 2 – Phân tích tần số dao động của dàn 10 thanh Trong ví dụ này, phương pháp phần tử hữu hạn sẽ được áp dụng để xác định tần số dao động của hệ dàn phẳng 10 thanh như trong Hình 3.2. Kết cấu có khối lượng riêng là 2770 kg/m3, mô-đun đàn hồi là 6.89.1010 N/m2, khối lượng cộng thêm từ nút 1 đến nút 4 là 454 kg, tiết diện các thanh lần lượt là A1 = 35.15 cm2, A2 = 14.76 cm2, A3 = 35.11 cm2, A4 = 14.72 cm2, A5 = 0.65 cm2, A6 = 4.55 cm2, A7 = 23.71 cm2, A8 = 23.63 cm2, A9 = 12.38 cm2, A10 = 12.45 cm2. Bài toán này được thực hiện nhằm mục đích xác minh tính chính xác và độ tin cậy của phương pháp phần tử hữu hạn FEM (Finite Elemesnt Method) tính toán tần số dao động của kết cấu dàn. Kết quả phân 46 tích của phương pháp phần tử hữu hạn FEM (Finite Elemesnt Method) được so sánh với kết quả nghiên cứu của Ho-Huu.V và cộng sự [47]. Hình 3.2: Mô hình bài toán dàn phẳng 10 thanh. Kết quả phân tích của phương pháp phần tử hữu hạn FEM (Finite Element Method) và kết quả nghiên cứu của Ho-Huu.V và cộng sự [47] được trình bày trong Bảng 3.3. Bảng 3.3: Kết quả tần số dao động trong hệ kết cấu. f 1
2
3
4
5 Tần số
Ho-Huu.V [47]
7.00
16.19
20.00
20.00
28.55 FEM
7.00
16.19
20.00
20.00
28.55 Sai số
%
0.00
0.00
0.00
0.00
0.00 Nhận xét: Bảng 3.3 cho thấy, kết quả tần số dao động được tính toán bằng phương pháp số giống với kết quả nghiên cứu của Ho-Huu.V và cộng sự [47]. Điều đó có nghĩa là tần số dao động được tính toán bởi Matlab bằng với kết quả nghiên cứu của Ho-Huu.V và cộng sự [47]. Do đó code phần tử hữu hạn để tính toán tần số dao động cho kết cấu dàn sử dụng trong luận văn là chính xác và có độ tin cậy cao. Ví dụ 3 – Kiểm chứng code lập trình giải thuật NSGA-II Trong ví dụ này, giải thuật NSGA-II được áp dụng để giải bài toán tối ưu hóa sau: 47 Hàm mục tiêu: Với Và x = (x1,…,xn)Tϵ [0,1]n với n=10. Số lượng cá thể: 50. Số vòng lặp: 200. Kết quả tối ưu của phương pháp NSGA-II và kết quả nghiên cứu của Zhang và cộng sự [6] được trình bày trong Hình 3.3. Để làm rõ khái niệm về sự trội (domination) trong Mục 2.4.1.2 và thuật toán phân loại không trội, 4 nghiệm nằm trên đường Pareto (Hình 3.3) được sử dụng để so sánh. Kết quả so sánh được thể hiện như trong Bảng 3.4. Luận văn Zhang và cộng sự [6] Hình 3.3: Kết quả bài toán tối ưu. Nhận xét: Hình 3.3 cho thấy, kết quả tối ưu bằng phương pháp NSGA-II giống với kết quả nghiên cứu của Zhang và cộng sự [6]. Điều này cho thấy code lập trình giải thuật NSGA-II để giải bài toán tối ưu hóa đa mục tiêu sử dụng trong luận văn là chính xác và có độ tin cậy cao. 48 Kết quả Bảng 3.4 cho thấy, trong số bốn lời giải thì lời giải A có (f1(min), f1(max)) và lời giải C có (f1(max), f1(min)). Điều này có nghĩa là lời giải A có f1 nhỏ hơn nhưng f2 lớn hơn lời giải C. Do đó, không có lời giải nào trong hai lời giải này có thể được cho là tốt hơn liên quan đến cả hai mục tiêu. Khi điều này xảy ra giữa hai lời giải như thế thì chúng được gọi là các lời giải không trội. Tương tự cho hai cặp lời giải khác là A- C và B-C. Trong ba lời giải (A đến C) thì bất kỳ cặp lời giải nào cũng có thể so sánh đối với cả hai mục tiêu. Bảng 3.4: Kết quả tối ưu Pareto. NSGA-II Biến thiết kế Điểm A Điểm B Điểm C Điểm D 1 0.66 0 0.68 x1 0 0.00 0 0.00 x2 0 0.00 0 0.00 x3 0 0.00 0 0.00 x4 0 0.00 0 0.00 x5 0 0.00 0 0.00 x6 0 0.00 0 0.00 x7 0 0.00 0 0.00 x8 0 0.00 0 0.01 x9 0 0.00 0 0.00 x10 1 0.66 0 0.68 f1 0 0.18 1 0.19 f2 Ngoài ra Bảng 3.4 cũng cho thấy, khi so sánh lời giải D với lời giải B thì nhận thấy rằng lời giải B là tốt hơn trong cả hai mục tiêu. Khi điều này xảy ra trong việc so sánh hai lời giải thì lời giải B được gọi là trội so với lời giải D hoặc lời giải D bị trội bởi lời giải B. Tiếp tục so sánh lời giải D và C thì lời giải D là tốt hơn ở mục tiêu thứ hai nhưng lại kém hơn ở mục tiêu thứ nhất. Như vậy, trong trường hợp không có các lời giải A, B và bất kỳ lời giải không trội khác thì lời giải D sẽ trong cùng một nhóm với lời giải C. Tuy nhiên, do thực tế có lời giải B đã chỉ ra rằng các lời giải B 49 và C là không bị trội đối với nhau, trong khi lời giải D là một lời giải bị trội. Vì thế, lời giải D là một lời giải tối ưu phụ và không thật sự tối ưu để sử dụng. Điều này đúng với khái niệm của nghiệm tối ưu đa mục tiêu. Ví dụ 4 – Kiểm chứng code lập trình RBDO-NSGA-II Trong ví dụng này, giải thuật RBDO-NSGA-II được áp dụng để giải bài toán tối ưu hóa có xét đến độ tin cậy với hai biến thiết kế ngẫu nhiên x = (x1, x2) với độ lệch chuẩn σ = 0.03. Thông số bài toán như sau: Hàm mục tiêu: Số lượng cá thể: 100. Số vòng lặp: 200. Chịu ràng buộc: Kết quả tối ưu của phương pháp RBDO-NSGA-II và kết quả nghiên cứu của Deb và cộng sự [48] được trình bày trong Hình 3.4. Luận văn Deb 2009 [48] Hình 3.4: Kết quả bài toán tối ưu. Nhận xét: Hình 3.4 cho thấy, kết quả tối ưu phương pháp RBDO-NSGA-II trong các trường hợp β = 1.28 (tương ứng độ an toàn 89.97%), β = 2 (tương ứng độ 50 an toàn 97.72%), β = 3 (tương ứng độ an toàn 99.87%) giống với kết quả nghiên cứu của Deb và cộng sự [48]. Do đó code lập trình giải thuật RBDO-NSGA-II để giải bài toán tối ưu hóa đa mục tiêu có kể đến độ tin cậy sử dụng trong luận văn có độ tin cậy cao. 3.2. Bài toán tối ưu hóa đa mục tiêu kết cấu dàn Bài toán kết cấu dàn phẳng 10 thanh Mô hình kết cấu dàn phẳng 10 thanh được cho như Hình 3.5. Bài toán thiết kế DO (Dertemined Optimization) và RBDO (Reliability Based-Design Optimization) nhằm cực tiểu khối lượng và cực tiểu chuyển vị của toàn hệ kết cấu. Bảng 3.5: Thông số bài toán kết cấu dàn phẳng 10 thanh. Thông số Giá trị Mô-đun đàn hồi E 6.895x1010 N/m2 khối lượng riêng ρ 2767 kg/m3 454 kg Khối lượng bổ sung W0 Chiều dài l 9.144 m 444.82x103 N Lực P = 7 Hz, = 15 Hz, = 20 Hz Giới hạn tần số [9]; [11]; [49] Giới hạn chuyển vị [19]–[22] 5.08 cm Giới hạn ứng suất [19]–[22] 172.375x106 N/m2 Giới hạn biến thiết kế 0.6452 cm2 ≤ Ai ≤ 225.80 cm2 Note: 1 in.2 = 6.452 cm2; 4.45 lb = 1 N; 2.2046 lb = 1 kg. Hình 3.5: Mô hình kết cấu dàn phẳng 10 thanh. 51 Bài toán được tối ưu lần lượt qua 2 trường hợp: Không xét đến độ tin cậy (DO - Dertemined Optimization) Trong trường hợp này bài toán thiết kế DO nhằm cực tiểu khối lượng và cực tiểu chuyển vị của toàn hệ kết cấu. Mô hình toán học của bài toán DO được cho như Hình 3.6. Trong đó diện tích mặt cắt ngang của các thanh dàn là biến thiết kế. Giá trị cận dưới và cận trên của biến thiết kế, mô-đun đàn hồi E, lực tác dụng P, khối lượng riêng của vật liệu ρ, khối lượng bổ sung tại các nút W0 (từ nút 1 đến nút 4), giới hạn chuyển vị, ứng suất và tần số dao động được cho như trong Bảng 3.5. Bài toán này đã được nghiên cứu trước đó bởi nhiều tác giả như Kaveh và cộng sự [19]; Li và cộng sự [20]; Lee và cộng sự [21]; Camp và cộng sự [22] v.v. Các bài toán này đều tối ưu hóa đơn mục tiêu, với mục tiêu là cực tiểu khối lượng chịu các ràng buộc về chuyển vị và ứng suất. Tuy nhiên, cho đến nay các bài toán tối ưu hoá đa mục tiêu cực tiểu khối lượng, cực tiểu chuyển vị chịu ràng buộc về chuyển vị, ứng suất và tần số dao động vẫn chưa được thực hiện. Trong trường hợp này bài toán được tối ưu với 3 trường hợp ràng buộc khác nhau: + Ràng buộc chuyển vị và ứng suất (tĩnh học) với giới hạn chuyển vị của hệ là 5.08 cm. + Ràng buộc chuyển vị, ứng suất (tĩnh học) và tần số dao động (động học) với giới hạn chuyển vị của hệ là 5.08 cm. + Ràng buộc chuyển vị, ứng suất (tĩnh học) và tần số dao động (động học) với giới hạn chuyển vị của hệ là 10.16 cm. 52 a) Không xét đến ràng buộc tần số dao động b) Có xét đến ràng buộc tần số dao động Hình 3.6: Mô hình bài toán tối ưu DO kết cấu dàn 10 thanh. Kết quả tối ưu của bài toán trong trường hợp này được thể hiện qua Hình 3.7. Kết quả là 3 đường Pareto tương ứng với 3 trường hợp ràng buộc khác nhau với trục đứng (trục tung) thể hiện chuyển vị và trục ngang (trục hoành) thể hiện khối lượng của hệ dàn. Để đánh giá hiệu quả của phương pháp NSGA-II, một số nghiệm trên 3 đường Pareto (Hình 3.7) được sử dụng để so sánh. Kết quả so sánh được thể hiện như trong Bảng 3.6. Ngoài ra Bảng 3.6 còn thể hiện kết quả nghiên cứu của Kaveh và cộng sự [19]; Li và cộng sự [20]; Lee và cộng sự [21]; Camp và cộng sự [22]. Bảng 3.6 cho thấy, kết quả chia làm hai phần: + Phần thứ nhất: bài toán tối ưu hoá tiền định không xét ràng buộc tần số dao động, tại vị trí cực tiểu khối lượng, phương pháp NSGA-II cho kết quả tốt với giá trị là 2300 kg (điểm A). Kết quả xấp xỉ kết quả nghiên cứu của Lee và cộng sự [21] 2294.24 kg; Li và cộng sự [20] 2295.62 kg; Kaveh và cộng sự [19] 2293.64 kg và tốt hơn kết quả của Camp và cộng sự [22] 2302.60 kg. Điều này cho thấy giải thuật NSGA-II là giải thuật tối ưu hiệu quả và có độ tin cậy cao. + Phần thứ hai: bài toán tối ưu hoá tiền định có xét ràng buộc tần số dao động, tại vị trí cực tiểu khối lượng, phương pháp NSGA-II cũng cho kết quả tốt với giá trị là 2314 kg (điểm B). Kết quả này không tốt bằng với kết quả nghiên cứu của Lee và cộng sự [21]; Li và cộng sự [20]; Kaveh và cộng sự [19]. Tuy nhiên, khi xét đến tần số dao động tại mode 1, mode 2 và mode 3 thì kết quả các nghiên cứu tham khảo không thỏa f1 ≥ 7, f2 ≥ 15 và f3 ≥ 20. Do đó có thể thấy rằng tần số dao động ảnh hưởng đến kết quả tối ưu và nghiên cứu này đã khắc phục được khuyết điểm nghiên 53 cứu của [19]–[22] là chỉ xét ràng buộc tĩnh và nghiên cứu của [9], [11], [49] là chỉ xét ràng buộc động. Hình 3.7: Kết quả tối ưu hóa DO kết cấu dàn 10 thanh. Ngoài ra Bảng 3.6 cho thấy, khi tăng giá trị ràng buộc chuyển vị từ 5.08 cm lên 10.16 cm thì tại vị trí cực tiểu khối lượng (điểm C), khối lượng của kết cấu cũng thay đổi theo nhưng kết quả vẫn thỏa ràng buộc ứng suất và tần số dao động. Hình 3.7 cho thấy, kết quả bài toán tối ưu không phải là một lời giải duy nhất như ở bài toán tối ưu đơn mục tiêu mà sẽ là một tập hợp của nhiều lời giải tối ưu nằm trên đường Pareto. Trong đó nếu lựa chọn kết cấu có khối lượng nhỏ thì chuyển vị của kết cấu lớn và ngược lại. Ngoài ra Hình 3.7 cho thấy, đường Pareto P2 tiệm cận đường Pareto P1 có nghĩa là tại cũng một vị trí chuyển vị thì khối lượng tối ưu khi không xét ràng buộc tần số và khi có xét đến ràng buộc tần số dao động là gần bằng nhau, tuy nhiên các nghiệm trên đường Pareto P2 đã khắc phục được nhược điểm của đường Pareto P1 không chỉ xét đến ràng buộc về chuyển vị, ứng suất (tĩnh học) mà còn xét đến ràng buộc tần số dao động (động học). Đường Pareto P2 trùng với đường Pareto P3 cho thấy giới hạn của các ràng buộc trong bài toán có thể thay đổi tuy nhiên nghiệm tối ưu thì chỉ có một đường Pareto duy nhất, và độ dài hay ngắn của đường Pareto này phụ thuộc vào giới hạn của các ràng buộc đó. 54 Bảng 3.6: Kết quả tối ưu hóa DO kết cấu dàn 10 thanh. Luận văn
DO-NSGA-II Kaveh và
cộng sự
[19] Biến thiết kế
(cm2) Lee và
cộng sự
[21]
HS
194.53
0.66
146.52
98.52
0.66
3.51
48.65
139.11
138.40
0.65 Li và
cộng sự
[20]
PSO
198.10
0.65
149.47
97.96
0.65
3.56
48.13
135.35
138.77
0.65 A1
A2
A3
A4
A5
A6
A7
A8
A9
A10
Khối lượng (kg)
Chuyển vị (cm) 195.17
0.65
141.11
105.24
0.65
4.23
49.20
139.89
136.83
0.645
2300
5.0797
f1 =5.92
f2 =10.78 106.47
5.27
69.81
44.68
0.65
5.30
38.76
66.47
68.88
1.650
1216
10.159
f1 =9.2
f2 =15.1 Điểm D(***)
211.85
5.95
155.63
92.80
0.65
10.17
52.23
130.54
135.33
0.645
2348
5.058
f1 =11.89
f2 =15.76 2294.24 2295.62 HPSACO Điểm A(*) Điểm B(**) Điểm C(***)
194.40
5.33
150.01
97.15
0.65
6.10
47.65
133.55
143.97
0.645
2314
5.0797
f1 =11.2
f2 =15.02 Camp và
cộng sự
[22]
GA
186.59
0.65
155.30
90.07
0.65
3.61
49.62
141.62
142.52
0.65
2302.60
5.08
f1 =5.94
f2 =10.4 5.08
f1 = 5.96
f2 = 10.3 5.08
f1 =5.93
f2 =10.4 195.54
0.65
151.20
100.04
0.65
3.38
47.98
136.00
136.97
0.65
2293.64
5.08
f1 =5.93
f2 =10.3 Tần số f3 =19.14 f3 =20.23 f3 =20 f3 =21.29 f3 =18.6 f3 = 18.5 f3 =18.5 f3 =18.3 (**): Giới hạn chuyển vị 5.08 cm, có xét ràng buộc tần số. (***): Giới hạn chuyển vị 10.16 cm, có xét ràng buộc tần số. Chú thích: (*): Giới hạn chuyển vị 5.08 cm, không xét ràng buộc tần số. 55 Có xét đến độ tin cậy RBDO (Reliability Based-Design Optimization) Trong trường hợp này bài toán thiết kế RBDO nhằm cực tiểu khối lượng, cực tiểu chuyển vị của toàn hệ kết cấu được phát biểu như trong Hình 3.8. Mô-đun đàn hồi E, lực tác dụng P, khối lượng bổ sung tại các nút W0 (từ nút 1 đến nút 4) và khối lượng riêng của vật liệu ρ được xem là những biến ngẫu nhiên với giá trị trung bình được cho như trong Bảng 3.5. Tất cả các biến ngẫu nhiên tuân theo luật phân phối Gauss (normal distribution) và được thống kê độc lập với hệ số thay đổi là 0.05. Giá trị giới hạn chuyển vị, ứng suất và tần số dao động được cho như trong Bảng 3.5. Chỉ số độ tin cậy yêu cầu là βt = 3, tương ứng với xác suất an toàn là 99.87%. Diện tích mặt cắt ngang của các thanh dàn vừa là biến thiết kế vừa là biến ngẫu nhiên. Trong trường hợp này bài toán được tối ưu với ràng buộc chuyển vị, ứng suất (tĩnh học) và tần số dao động (động học) với giới hạn chuyển vị của hệ là 5.08 cm. Hình 3.8: Mô hình bài toán tối ưu RBDO kết cấu dàn 10 thanh. Bài toán này đã được nghiên cứu trước đó bởi nhiều tác giả như Shayanfar và cộng sự [32]; Lee và cộng sự [33] v.v. Các bài toán này đều là tối ưu hóa đơn mục tiêu, với mục tiêu là cực tiểu khối lượng chịu ràng buộc về chuyển vị và ứng suất và có xét đến độ tin cậy với chỉ số chỉ số độ tin cậy yêu cầu là βt = 3. Tuy nhiên, cho đến nay các bài toán tối ưu hoá đa mục tiêu gồm cực tiểu khối lượng, cực tiểu chuyển vị chịu ràng buộc về chuyển vị, ứng suất và tần số dao động có xét đến độ tin cậy vẫn chưa được thực hiện. 56 Tương tự như kết quả tối ưu DO, kết quả tối ưu bài toán này cũng là một đường Pareto như trong Hình 3.9 với trục đứng (trục tung) thể hiện chuyển vị và trục ngang (trục hoành) thể hiện khối lượng của hệ dàn. Để đánh giá hiệu quả của bài toán trong phần này, kết quả tối ưu DO-NSGA-II trong Mục 3.2.1.1 ứng với ràng buộc chuyển vị, ứng suất (tĩnh học) và tần số dao động (động học) với giới hạn chuyển vị của hệ là 5.08 cm cũng được xem xét. Để đánh giá thêm sự khác biệt giữa bài toán DO và bài toán RBDO, một số nghiệm nằm trên 2 đường Pareto P1 và Pareto P2 (Hình 3.9) sẽ được sử dụng để so sánh. Kết quả so sánh được thể hiện như trong Bảng 3.7. Ngoài ra, Bảng 3.7 còn thể hiện kết quả nghiên cứu của Shayanfar và cộng sự [32]; Lee và cộng sự [33]. Bảng 3.7 cho thấy tại vị trí cực tiểu khối lượng, kết quả tối ưu bài toán RDBO (điểm C) 3184 kg lớn kết quả tối ưu bài toán DO-NSGA-II (điểm A) 2314 kg. Tuy nhiên độ an toàn của kết cấu trong bài toán DO-NSGA-II chịu ảnh hưởng của các yếu tố ngẫu nhiên là rất thấp. Cụ thể, sử dụng phương pháp đánh giá độ tin cậy bậc nhất FORM để tính chỉ số độ tin cậy tại điểm A thì giá trị β lần lượt là β1= 0.04 tương ứng 51.74%, β2= 0.09 tương ứng 53.46%, β4= 0.03 tương ứng 51.36% và β5= 0.34 tương ứng 63.36%. Điều đó có nghĩa là mức an toàn của kết cấu chỉ đạt 51.36%. Tại điểm C, giá trị β lần lượt là β1= 3.00 tương ứng 99.87%, β2= 3.04 tương ứng 99.88%, β3= 13.41 tương ứng 100%, β4= 3.11 tương ứng với 99.91% và β5= 3.00 tương ứng với 99.87%. Điều đó có nghĩa kết cấu luôn đảm bảo mức độ an toàn cho trước là 99.87%. Bảng 3.7 cũng cho thấy rằng, khi xét tại cũng một vị trí chuyển vị và khối lượng (điểm B và điểm C). Tại điểm B, trong 5 giá trị β chỉ có 3 giá trị β ≥ 3 là β1 (ràng buộc ). Trong khi đó tại điểm C kết về chuyển vị), β3 (ràng buộc về ), β4 (ràng buộc về cấu luôn đảm bảo mức độ an toàn cho trước là 99.87%. Điều đó cho thấy hai bài toán DO và RBDO cho kết quả tối ưu giống nhau, tuy nhiên dưới tác động của các yếu tố ngẫu nhiên thì nghiệm tối ưu RBDO luôn đảm bảo cho kết cấu an toàn ở mức độ cho trước, còn nghiệm tối ưu DO thì không. Do đó bài toán DO chỉ có thể áp dụng để giải các bài toán trong điều kiện lý tưởng (bài toán lý thuyết), còn bài toán tối ưu RBDO là 57 một bài toán thực tế, kết quả bài toán luôn đảm bảo độ an toàn nhất định cho kết cấu vận hành trong những điều kiện thiết kế ngẫu nhiên cụ thể. Ngoài ra Bảng 3.7 cũng cho thấy, kết quả tối ưu trong luận văn lớn hơn so với kết quả nghiên cứu của Shayanfar và cộng sự [32]; nghiên cứu của Lee và cộng sự [33]. Điều này là do các nghiên cứu tham khảo không xét đến ràng buộc tần số dao động mà chỉ xét đến ràng buộc về chuyển vị và ứng suất. Qua đó cho thấy nghiên cứu của Shayanfar và cộng sự [32]; Lee và cộng sự [33] còn gặp hạn chế khi giải bài toán tối ưu mà không đưa ràng buộc tần số dao động vào bài toán. Trong khi tần số dao động sẽ ảnh hưởng đáng kể đến kết quả tối ưu. Do đó việc giải bài toán này trong luận văn đã giải quyết được hạn chế trong nghiên cứu của Shayanfar [32] và Lee [33] khi không xét đến đến ràng buộc của tần số dao động. 58 Bảng 3.7: Kết quả tối ưu hóa RBDO kết cấu dàn 10 thanh. Luận văn Shayanfar và
cộng sự [32] Lee và cộng
sự [33] Biến thiết kế
(cm2) DO-NSGA-II DLM-GA TPB - β A1
A2
A3
A4
A5
A6
A7
A8
A9
A10
Khối lượng (kg)
Chuyển vị (cm) Tần số 221.64
0.65
191.51
169.53
0.65
0.65
21.53
182.94
168.64
0.65
2817.43
-
f1 = 5.12
f2 = 6.68
f3 = 13.12 248.47
0.65
173.43
123.04
0.65
0.65
39.87
176.46
175.04
0.65
2787.17
-
f1 = 5.12
f2 = 6.7
f3 = 16.43 RBDO-NSGA-II
Điểm C
225.43
7.04
215.34
137.69
0.65
8.95
66.91
197.47
204.04
0.645
3184
3.72
f1 = 12.68
f2 = 17.10
f3 = 22.31 Điểm A
194.40
5.33
150.01
97.15
0.65
6.10
47.65
133.55
143.97
0.645
2314
5.0797
f1 =11.2
f2 =15.02
f3 =20.23 Điểm B
225.81
7.32
220.77
123.36
0.65
26.36
51.83
207.26
202.88
0.645
3184
3.72
f1 = 13.38
f2 = 17.17
f3 = 20.70
β1= 0.04 (51.74%) β1= 2.96 (99.84%) β1= 3.00 (99.87%) β2= 0.09 (53.46%) β2= 0.06 (52.26%) β2= 3.04 (99.88%) βFORM 3.060 2.999 β3= 11.46 (100%) β3= 14.13 (100%) β3= 13.41 (100%) β4= 0.03 (51.36%) β4= 3.26 (99.94%) β4= 3.11 (99.91%) β5= 0.34 (63.36%) β5= 0.94(82.68%) β5= 3.00 (99.87%) 59 Hình 3.9 cho thấy, kết quả bài toán tối ưu không phải là một lời giải duy nhất như ở bài toán tối ưu đơn mục tiêu mà sẽ là một tập hợp của nhiều lời giải tối ưu nằm trên đường Pareto. Ví dụ nếu lựa chọn kết cấu có khối lượng nhỏ thì chuyển vị của kết cấu lớn và ngược lại. Hình 3.9: Kết quả tối ưu hóa RBDO kết cấu dàn 10 thanh. Ngoài ra, Hình 3.9 cũng cho thấy đường Pareto P2 nằm trùng với đường Pareto P1, tuy nhiên đường P2 ngắn hơn đường P1. Điều này có nghĩa tập hợp các nghiệm của P1 là nghiệm tối ưu nhất, trong đó có cả nghiệm thỏa giá trị β và nghiệm không thỏa giá trị β. Còn khi xét đến độ tin cậy thì tất cả các nghiệm của P2 là nghiệm tối ưu nhất, và đều thỏa giá trị β cho trước. Kết quả đạt được cho thấy tính logic hợp lý của bài toán, kết quả tối ưu dựa trên độ tin cậy RBDO có giá trị lớn hơn so với kết quả tối ưu không xét đến độ tin cậy DO. Sự gia tăng chi phí này là hiển nhiên nhằm đảm bảo độ an toàn nhất định cho kết cấu vận hành trong những điều kiện thiết kế ngẫu nhiên cụ thể. Kết quả tối ưu nhằm khảo sát sự ảnh hưởng của các mức độ tin cậy (an toàn) khác nhau lên kết cấu được thể hiện trên Hình 3.10. Kết quả đạt được là 4 đường Pareto tương ứng với 4 giá trị β khác nhau (β =2; β =2.5; β =3; và β =3.5) với trục đứng (trục tung) 60 thể hiện chuyển vị và trục ngang (trục hoành) thể hiện khối lượng của hệ dàn. Kết quả tại vị trí cực tiểu khối lượng và cực đại chuyển vị của hệ dàn được trình bày trong Bảng 3.8. Kết quả Bảng 3.8 cho thấy khi mức độ tin cậy yêu cầu càng tăng (hay chỉ số độ tin cậy tăng) thì khối lượng của kết cấu cũng tăng theo. Điều này là do khả năng chịu lực của kết cấu sẽ càng cao (hay độ an toàn của kết cấu sẽ càng cao) khi đó diện tích của các thanh dàn càng lớn và dẫn đến khối lượng của kết cấu sẽ càng tăng. Vì vậy, tuỳ thuộc vào tầm quan trọng của kết cấu khi thiết kế mà người thiết kế nên đưa ra những mức độ an toàn phù hợp nhằm giảm chi phí và tránh sự lãng phí không cần thiết. Bảng 3.8: Kết quả thiết kế RBDO của bài toán kết cấu 10 thanh ở các mức độ an toàn khác nhau. Biến thiết kế
diện tích (cm2)
A1
A2
A3
A4
A5
A6
A7
A8
A9
A10
Khối lượng (kg)
Chuyển vị βFORM β = 2
(Điểm A)
225.79
6.81
184.52
116.97
0.65
7.12
62.36
173.34
181.40
0.645
2866
4.11
2.00
2.30
12.71
2.72
2.43 β = 2.5
(Điểm B)
225.80
6.92
201.79
131.26
0.65
6.96
65.99
187.13
184.29
0.645
3018
3.91
2.50
2.91
12.84
2.92
2.52 β = 3
(Điểm C)
225.43
7.04
215.34
137.69
0.65
8.95
66.91
197.47
204.04
0.645
3184
3.72
3.00
3.04
13.41
3.11
3.00 β = 3.5
(Điểm D)
225.78
7.45
224.42
153.49
0.65
8.51
73.57
214.05
212.93
0.645
3363
3.55
3.50
3.98
13.47
3.70
3.50 β1
β2
β3
β4
β 5 61 Hình 3.10, minh họa kết quả của 4 đường Pareto tương ứng với 4 giá trị β khác nhau. Cả 4 đường Pareto này nằm chồng lên nhau, tuy nhiên giá trị β càng lớn thì đường Pareto càng ngắn lại. Điều này cho thấy tính logic hợp lý của bài toán RBDO, khi mức độ tin cậy yêu cầu càng tăng (hay chỉ số độ tin cậy tăng) thì khối lượng của kết cấu cũng tăng theo, và chuyển vị của kết cấu càng nhỏ lại. Hình 3.10: Kết quả thiết kế RBDO của bài toán kết cấu 10 thanh ở các mức độ an toàn khác nhau. Bài toán kết cấu dàn không gian 72 thanh Bài toán này được thực hiện với mục đích: + Đánh giá hiệu quả của việc áp dụng giải thuật NSGA-II và RBDO-NSGA-II vào bài toán kết cấu dàn không gian. + Đánh giá hiệu quả của việc giải bài toán tối ưu khi xét thêm ràng buộc về tần số dao động. + Đánh giá ảnh hưởng của việc thay đổi các giới hạn ràng buộc đến kết quả tối ưu. Xét kết cấu dàn không gian 72 thanh như Hình 3.11. Bài toán thiết kế DO (Dertemined Optimization) và RBDO (Reliability Based-Design Optimization) nhằm 62 cực tiểu khối lượng và cực tiểu chuyển vị của toàn hệ kết cấu. Kết cấu có 72 thanh và được chia làm 16 nhóm tương ứng với 16 biến thiết thế bao gồm: (1) A1-A4, (2) A5-A12, (3) A13-A16, (4) A17-A18, (5) A19-A22, (6) A23-A30, (7) A31-A34, (8) A35-A36, (9) A37-A40, (10) A41-A48, (11) A49-A52, (12) A53-A54, (13) A55-A58, (14) A59-A66 (15), A67-A70, và (16) A71- A72. Bảng 3.9: Thông số bài toán kết cấu dàn không gian 72 thanh. Thông số Giá trị Mô-đun đàn hồi E 6.895x1010 N/m2 Khối lượng riêng ρ 2767 kg/m3 2268 kg Khối lượng bổ sung W0 Chiều dài l 1.524 m Giới hạn tần số [9]–[11]; [49] = 4 Hz, = 6 Hz Giới hạn chuyển vị [50]; [21]; [20]; [51] 0.635 cm Giới hạn ứng suất [50]; [21]; [20]; [51] 172.375x106 N/m2 Giới hạn biến thiết kế 0.6452 cm2 ≤ Ai ≤ 20.65 cm2 Bảng 3.10: Hai trường hợp tải tác dụng lên kết cấu dàn không gian 72 thanh. Trường hợp 2 Nút 17
18
19
20 Trường hợp 1
Py (kN)
22.25
0.0
0.0
0.0 Px (kN)
22.25
0.0
0.0
0.0 Pz (kN)
-22.25
0.0
0.0
0.0 Px (kN) 0.0
0.0
0.0
0.0 Py (kN) Pz (kN)
-22.25
-22.25
-22.25
-22.25 0.0
0.0
0.0
0.0 Thông số bài toán được cho như Bảng 3.9. Kết cấu được thiết kế cho hai trường hợp tải khác nhau như Bảng 3.10 ([50] [21] [20] [51]). 63 Hình 3.11: Mô hình kết cấu dàn không gian 72 thanh. Tương tự như bài toán dàn phẳng 10 thanh, bài toán này cũng được tối ưu lần lượt qua 2 trường hợp: Không xét đến độ tin cậy (DO - Dertemined Optimization) Trong trường hợp này bài toán thiết kế DO nhằm cực tiểu khối lượng và cực tiểu chuyển vị của toàn hệ kết cấu. Mô hình toán học của bài toán DO được cho như Hình 3.12. Trong đó diện tích mặt cắt ngang của các thanh dàn là biến thiết kế. Giá trị cận dưới và cận trên của biến thiết kế, mô-đun đàn hồi E, khối lượng riêng của vật liệu ρ, khối lượng bổ sung tại các nút W0 (từ nút 17 đến nút 20), giới hạn chuyển vị, ứng suất và tần số dao động được cho như trong Bảng 3.9. 64 a) Không xét đến ràng buộc tần số dao động b) Có xét đến ràng buộc tần số dao động Hình 3.12: Mô hình tối ưu DO bài toán kết cấu dàn 72 thanh. Bài toán này đã được nghiên cứu trước đó bởi nhiều tác giả như Erbatur và cộng sự [50]; Lee và cộng sự [21]; Li và cộng sự [20]; Perez và Behdinan [51] v.v. Các bài toán này đều là tối ưu hóa đơn mục tiêu, với mục tiêu là cực tiểu khối lượng chịu ràng buộc về chuyển vị và ứng suất. Tương tự như bài toán tối ưu kết cấu dàn 10 thanh, cho đến nay các bài toán tối ưu hóa đa mục tiêu kết cấu dàn 72 thanh gồm cực tiểu khối lượng, cực tiểu chuyển vị chịu ràng buộc về chuyển vị, ứng suất và tần số dao động vẫn chưa được thực hiện. Giống như bài toán tối ưu dàn phẳng 10 thanh, trong trường hợp này bài toán cũng được tối ưu với 3 trường hợp ràng buộc khác nhau: + Ràng buộc chuyển vị và ứng suất (tĩnh học) với giới hạn chuyển vị của hệ là 0.635 cm. + Ràng buộc chuyển vị, ứng suất (tĩnh học) và tần số dao động (động học) với giới hạn chuyển vị của hệ là 0.635 cm. + Ràng buộc chuyển vị, ứng suất (tĩnh học) và tần số dao động (động học) với giới hạn chuyển vị của hệ là 1.27 cm. Kết quả tối ưu của bài toán trong trường hợp này được thể hiện qua Hình 3.13. Kết quả cho thấy 3 đường Pareto tương ứng với 3 trường hợp ràng buộc khác nhau với trục đứng (trục tung) thể hiện chuyển vị và trục ngang (trục hoành) thể hiện khối lượng 65 của hệ dàn. Kết quả tại vị trí cực tiểu khối lượng được trình bày trong Bảng 3.11 với điểm A là kết quả trường hợp kết cấu chịu ràng buộc tĩnh học với giới hạn chuyển vị là 0.635 cm, điểm B là kết quả trường hợp kết cấu chịu ràng buộc tĩnh học và động học với giới hạn chuyển vị là 0.635 cm, điểm C là kết quả trường hợp kết cấu chịu ràng buộc tĩnh học và động học với giới hạn chuyển vị là 1.27 cm, điểm D kết quả tối ưu trường hợp kết cấu không chịu ràng buộc tĩnh học tại vị trí chuyển vị là 0.508 cm. Ngoài ra các kết quả trên cũng được so sánh với kết quả nghiên cứu của Erbatur và cộng sự [50]; Lee và cộng sự [21]; Li và cộng sự [20]; Perez và Behdinan [51] để đánh giá hiệu quả của bài toán. Kết quả Bảng 3.11 cho thấy rằng: + Đối với bài toán tối ưu hoá DO không chịu ràng buộc tần số dao động, tại vị trí cực tiểu khối lượng, phương pháp NSGA-II cho kết quả tốt với giá trị là 171.22 kg (điểm A). Kết quả này tốt hơn kết quả của Erbatur và cộng sự [50], Lee và cộng sự [21], Perez và Behdinan [51]. Điều này một lần nữa khẳng định rằng NSGA-II là giải thuật tối ưu hiệu quả và có độ tin cậy cao, không những có thể áp dụng cho những bài kết cấu dàn phẳng đơn giản mà còn có thể áp dụng cho những bài toán kết cấu dàn không gian với số lượng thanh lớn. + Đối với bài toán tối ưu hoá DO có xét đến ràng buộc tần số dao động, kết quả cho thấy tại vị trí cực tiểu khối lượng, phương pháp DO-NSGA-II cho kết quả tốt với giá trị là 325.86 kg (điểm B). Kết quả này không tốt bằng kết quả nghiên cứu của Erbatur và cộng sự [50]; Lee và cộng sự [21]; Li và cộng sự [20]; Perez và Behdinan [51]. Tuy nhiên khi xét về tần số dao động tại mode 1 và mode 3 thì các nghiên cứu tham khảo trên không thỏa f1 = 4 và f3 ≥ 6, do đó ta có thể khẳng định lần nữa rằng tần số dao động ảnh hưởng lớn đến kết cấu dàn và nghiên cứu này đã khắc phục được những khuyết điểm của những nghiên cứu trên khi không xét đến ràng buộc dao động. + Tại cùng một vị trí chuyển vị 0.5 cm, kết quả bài toán tối ưu không xét đến ràng buộc tần số có giá trị là 224.57 kg (điểm D), kết quả này tốt hơn kết quả tối ưu xét đến 66 ràng buộc tần số là 325.86 kg (điểm B và điểm C). Kết quả đạt được cho thấy tính logic hợp lý của bài toán, kết quả tối ưu khi xét nhiều ràng buộc có giá trị lớn hơn so với kết quả tối ưu xét ít ràng buộc. Sự gia tăng chi phí này là hiển nhiên nhằm đảm bảo độ an toàn nhất định cho kết cấu làm việc trong cả môi trường tĩnh và động. Ngoài ra Bảng 3.11 cũng cho thấy, khi tăng giá trị ràng buộc chuyển vị (từ 0.635 cm lên 1.27 cm) thì giá trị khối lượng và chuyển vị của hệ không thay đổi (điểm B và điểm C). Điều đó có nghĩa rằng việc tăng giá trị chuyển vị không đồng nghĩa với việc khối lượng của kết cấu nhỏ đi một cách tương ứng vì trong bài toán này, ta không những xét giới hạn chuyển vị mà còn xét ràng buộc về ứng suất và tần số dao động. Do đó kết quả này cho thấy một bài toán tối ưu khi xét đến nhiều ràng buộc sẽ cho kết quả càng chính xác và các ràng buộc này có ảnh hưởng qua lại với nhau. Hình 3.13: Kết quả tối ưu hóa DO kết cấu dàn 72 thanh. Kết quả Hình 3.13 một lần nữa cho thấy, kết quả tối ưu không phải là một lời giải duy nhất như ở bài toán tối ưu đơn mục tiêu mà sẽ là một tập hợp của nhiều lời giải tối ưu. Trong đó nếu lựa chọn kết cấu có khối lượng nhỏ thì chuyển vị của kết cấu lớn và ngược lại. 67 Hình 3.13 cũng cho thấy kết quả tối ưu là 3 đường Pareto P1, Pareto P2 và Pareto P3. Kết quả này cho thấy: + Vùng I (giá trị khối lượng biến thiên từ 171.22 kg (điểm A) đến 353 kg (điểm F) tương ứng chuyển vị giảm từ 0.643 cm xuống 0.318 cm), Pareto P2 và Pareto P3 trùng lên nhau, đường Pareto P1 nằm bên trái Pareto P2 và Pareto P3. Cả ba đường Pareto này đều hội tụ tại điểm F. Điều này có nghĩa là tại cùng một vị trí chuyển vị, kết quả khối lượng bài toán tối ưu không xét đến ràng buộc tần số luôn nhỏ hơn kết quả bài toán tối ưu xét đến ràng buộc tần số. Do đó, có thể khẳng định lần nữa rằng tần số dao động có ảnh hưởng lớn đến bài toán tối ưu do đó việc xét đến ràng buộc tần số dao động là cần thiết và nên áp dụng cho các bài toán kết cấu dàn chịu tải trọng động. Việc xét đến ràng buộc tần số dao động tương tự như việc chúng ta đưa thêm độ tin cậy vào trong kết cấu, làm cho kết cấu bền vững hơn. + Khi khối lượng kết cấu tăng dần, chuyển vị kết cấu càng nhỏ, kết cấu càng bền vững thì ba đường Pareto sẽ hội tụ lại một điểm (Điểm F). Sau điểm hội tụ cả 3 đường Pareto này trùng lên nhau và tất các các nghiệm nằm trên ba đường Pareto sau điểm hội tụ là các nghiệm tối ưu duy nhất (vùng II - giá trị khối lượng biến thiên từ 353 kg (điểm F) đến 977 kg (điểm G) tương ứng chuyển vị giảm từ 0.318 cm xuống 0.148 cm). Điều này có nghĩa là khi khối lượng kết cấu càng lớn, chuyển vị kết cấu càng nhỏ, kết cấu càng bền vững và các nghiệm tối ưu nằm trên đường Pareto là duy nhất, kể cả khi không xét đến ràng buộc tần số dao động. Ngoài ra Hình 3.13 cũng cho thấy sự tương quan giữa chuyển vị và khối lượng kết cấu, sự ảnh hưởng của các ràng buộc lên bài toán tối ưu. Từ đó chúng ta có thể thiết lập và giải một bài toán tối ưu hóa đa mục tiêu cho kết cấu, và tùy theo các điều kiện đầu ra mà lựa chọn một nghiệm thích hợp nằm trên đường Pareto. Tất cả các nghiệm nằm trên đường Pareto đều phải thỏa tất cả các điều kiện ràng buộc. 68 Bảng 3.11: Kết quả tối ưu hóa DO kết cấu dàn 72 thanh. Luận văn
DO-NSGA-II Biến thiết kế
(diện tích cm2) Điểm A(*) Điểm D(*) Điểm B(**) Điểm C(***) A1
A2
A3
A4
A5
A6
A7
A8
A9
A10
A11
A12
A13
A14
A15
A16 Khối lượng (kg)
Chuyển vị (cm) Tần số Perez và
Behdinan
[51]
1.04
3.29
3.20
3.63
3.32
3.53
0.65
0.71
8.44
3.33
0.65
0.65
11.24
3.35
0.65
0.65
173.23
-
f1 = 1.8
f3 = 3.9 11.65
3.15
0.65
0.65
8.90
3.37
0.65
0.65
3.94
3.21
0.65
0.65
1.02
3.27
3.18
3.40
171.22
0.643
f1 = 2.8
f3 = 3.8 Erbatur
và cộng
sự [50]
1.00
3.45
3.10
3.36
2.97
3.42
0.77
1.06
7.45
3.77
0.65
0.65
11.32
3.26
0.68
0.68
174.98
-
f1 = 1.7
f3 = 4.0 Lee và
cộng sự
[21]
11.55
3.36
0.65
0.65
7.93
3.37
0.65
0.65
3.34
3.25
0.65
0.65
1.01
3.53
2.85
4.45
172.04
-
f1 = 2.8
f3 = 3.9 Li và
cộng sự
[20]
11.97
3.25
0.65
0.65
8.08
3.26
0.65
0.65
3.21
3.28
0.65
0.65
0.65
3.39
2.54
3.45
167.67
-
f1 = 2.7
f3 = 3.9 15.77
4.74
1.49
1.26
9.57
3.55
1.02
0.67
5.04
3.28
0.65
0.65
1.18
4.81
3.59
7.82
224.57
0.508
f1 = 3.12
f3 = 4.25 17.11
7.69
1.14
0.65
12.77
7.76
0.65
0.65
7.81
7.84
0.65
0.65
3.43
8.50
0.65
0.65
325.86
0.503
f1 = 4.0
f3 = 6.0 17.11
7.83
1.39
0.65
12.25
7.70
0.65
0.65
8.46
8.09
0.65
0.65
3.31
8.12
0.65
0.65
325.86
0.503
f1 = 4.0
f3 = 6.0 (**): Giới hạn chuyển vị 0.635 cm, có xét ràng buộc tần số. (***): Giới hạn chuyển vị 1.27 cm, có xét ràng buộc tần số. Chú thích: (*): Giới hạn chuyển vị 0.635 cm, không xét ràng buộc tần số. 69 Có xét đến độ tin cậy RBDO (Reliability Based-Design Optimization) Trong trường hợp này bài toán thiết kế RBDO nhằm cực tiểu khối lượng, cực tiểu chuyển vị của toàn hệ kết cấu như được phát biểu trong Hình 3.14. Mô-đun đàn hồi E, lực tác dụng P, khối lượng bổ sung tại các nút W0 (từ nút 17 đến nút 20) và khối lượng riêng của vật liệu ρ được xem là những biến ngẫu nhiên với giá trị trung bình được cho như trong Bảng 3.9. Tất cả các biến ngẫu nhiên tuân theo luật phân phối Gauss (normal distribution) và được thống kê độc lập với hệ số thay đổi là 0.05. Giá trị giới hạn chuyển vị, ứng suất và tần số dao động được cho như trong Bảng 3.9. Chỉ số độ tin cậy yêu cầu là βt = 3, tương ứng với xác suất an toàn là 99.87%. Diện tích mặt cắt ngang của các thanh dàn vừa là biến thiết kế vừa là biến ngẫu nhiên. Trong trường hợp này bài toán được tối ưu với ràng buộc chuyển vị, ứng suất (tĩnh học) và tần số dao động (động học) với giới hạn chuyển vị của hệ là 0.635 cm. Hình 3.14: Mô hình bài toán tối ưu RBDO kết cấu dàn 72 thanh. Dựa trên thông tin được thu thập cho đến nay của tác giả, bài toán thiết kế RBDO cho kết cấu này cũng chưa được thực hiện. Giống như kết quả tối ưu DO, kết quả tối ưu bài toán này cũng là một đường Pareto như trong Hình 3.15 với trục đứng (trục tung) thể hiện chuyển vị và trục ngang (trục hoành) thể hiện khối lượng của hệ dàn. Để đánh giá hiệu quả của bài toán trong phần này, kết quả tối ưu DO-NSGA-II trong Mục 3.2.1.1 ứng với ràng buộc chuyển vị, ràng buộc ứng suất (tĩnh học) và ràng buộc tần số dao động (động học) với giới hạn chuyển vị của hệ là 0.635 cm cũng được xem xét. Để đánh giá thêm sự khác biệt giữa 70 bài toán DO và bài toán RBDO, một số nghiệm nằm trên 2 đường Pareto P1 và Pareto P2 (Hình 3.15) sẽ được sử dụng để so sánh. Kết quả so sánh được thể hiện như trong Bảng 3.12, trong đó điểm A và điểm B là hai nghiệm thuộc Pareto P1 và điểm C và điểm D là 2 nghiệm thuộc Pareto P2. Hình 3.15: Kết quả tối ưu hóa RBDO kết cấu dàn 72 thanh. Kết quả Bảng 3.12 cho thấy chỉ số độ tin cậy ảnh hưởng lớn đến kết quả tối ưu. Cụ thể tại vị trí cực tiểu khối lượng, kết quả tối ưu bài toán RDBO (điểm C) 412.9 kg lớn kết quả tối ưu bài toán DO-NSGA-II (điểm A) 326 kg. Tuy nhiên khi sử dụng phương pháp đánh giá độ tin cậy bậc nhất FORM để tính chỉ số độ tin cậy tại điểm A thì giá trị β lần lượt là β1= 2.25 tương ứng 98.77%, β2= 11.88 tương ứng 100%, β3= 0.00 tương ứng 42.91% và β4= 0.00 tương ứng 43.11%. Điều đó có nghĩa là mức an toàn của kết cấu chỉ đạt 42.91%. Tại điểm C, giá trị β lần lượt là β1= 4.10 tương ứng 99.99%, β2= 14.19 tương ứng 100%, β3= 3.00 tương ứng 99.87% và β4= 3.00 tương ứng với 99.87%. Điều đó có nghĩa kết cấu luôn đảm bảo mức độ an toàn cho trước là 99.87%. Bảng 3.12 cũng cho thấy, khi xét tại một vị trí chuyển vị (điểm B và điểm C) thì khối lượng của kết cấu khi xét đến chỉ số độ tin cậy lớn hơn so với khi không xét đến chỉ số độ tin cậy, trong trường hợp này khối lượng tăng xấp xỉ 78.9 kg. Tuy nhiên khi xét 71 đến độ an toàn của kết cấu thì tại điểm B giá trị β lần lượt là β1= 4.2 tương ứng 100%, β2= 11.96 tương ứng 100%, β3= 0.00 tương ứng 48.26% và β4= 0.004 tương ứng với 50.16%. Điều đó có nghĩa là kết cấu chỉ đảm bảo an toàn ở mức rất thấp là 48.26%. Trong khi đó tại điểm C kết cấu luôn đảm bảo mức độ an toàn cho trước là 99.87%. Do đó, kết quả đạt được cho thấy tính logic hợp lý của bài toán, kết quả tối ưu dựa trên độ tin cậy RBDO có giá trị lớn hơn so với kết quả tối ưu không xét đến độ tin cậy DO. Sự gia tăng chi phí này là hiển nhiên nhằm đảm bảo độ an toàn nhất định cho kết cấu vận hành trong những điều kiện thiết kế ngẫu nhiên cụ thể. Ngoài ra Bảng 3.12 cũng cho thấy rằng, khi xét tại một vị trí chuyển vị, cùng vị trí khối lượng (điểm D và điểm E) khi sử dụng phương pháp đánh giá độ tin cậy bậc nhất FORM để tính chỉ số độ tin cậy tại điểm D thì giá trị β lần lượt là β1= 8.92 tương ứng 100%, β2= 11.53 tương ứng 100%, β3= 1.77 tương ứng 96.13% và β4= 1.30 tương ứng 90.24%. Điều đó có nghĩa là mức an toàn của kết cấu chỉ đạt 90.24%. Tại điểm E, giá trị β lần lượt là β1= 8.87 tương ứng 100%, β2= 9.59 tương ứng 100%, β3= 3.05 tương ứng 99.88% và β4= 3.10 tương ứng với 99.90%. Điều đó có nghĩa kết cấu luôn đảm bảo mức độ an toàn cho trước là 99.88%. Do đó bài toán DO chỉ có thể áp dụng để giải các bài toán trong điều kiện lý tưởng (bài toán lý thuyết), khi xét đến những điều kiện thiết kế ngẫu nhiên cụ thể thì độ an toàn của kết cấu thấp, kể cả khi lựa chọn những kết cấu có khối lượng lớn. Còn bài toán tối ưu RBDO là một bài toán thực tế, kết quả bài toán luôn đảm bảo độ an toàn nhất định cho kết cấu vận hành trong những điều kiện thiết kế ngẫu nhiên cụ thể. Cũng như kết quả tối ưu hóa đa mục tiêu DO, bài toán tối ưu hóa đa mục tiêu RBDO cũng cho kết quả là một đường Pareto như Hình 3.15. Kết quả cho thấy đường Pareto P2 một phần nằm trùng với đường Pareto P1 và một phần nằm lệch sang bên phải so với đường Pareto P1. Tuy nhiên đường Pareto P2 ngắn hơn đường Pareto P1. Điều này một lần nữa khẳng định chỉ số độ tin cậy ảnh hưởng lớn đến kết quả tối ưu. Trong tập hợp các nghiệm tối ưu của Pareto P1 sẽ có nghiệm thỏa giá trị β và có nghiệm không 72 thỏa giá trị β. Còn khi xét đến độ tin cậy thì tất cả nghiệm tối ưu của P2 đều thỏa giá trị β cho trước. Qua đó một lần nữa có thể khẳng định tính logic hợp lý của bài toán tối ưu RBDO. 73 Bảng 3.12: Kết quả tối ưu hóa RBDO kết cấu dàn 72 thanh. Luận văn RBDO-NSGA-II Biến thiết kế
(diện tích cm2) A1
A2
A3
A4
A5
A6
A7
A8
A9
A10
A11
A12
A13
A14
A15
A16
Khối lượng (kg)
Chuyển vị (cm) Tần số Điểm A
17.11
7.69
1.14
0.65
12.77
7.76
0.65
0.65
7.81
7.84
0.65
0.65
3.43
8.50
0.65
0.65
325.86
0.503
f1 = 4.0
f3 = 6.0 DO-NSGA-II
Điểm B
17.28
7.70
1.17
0.65
12.80
7.76
0.65
0.65
7.85
8.27
0.68
0.65
3.52
8.50
0.72
2.28
334
0.415
f1 = 4.02
f3 = 6.04 Điểm C
20.51
9.57
2.13
0.65
16.86
10.20
0.65
0.65
10.44
10.65
0.65
0.65
4.11
10.02
0.65
0.75
412.9
0.419
f1 = 4.50
f3 = 6.76 Điểm E
20.64
9.59
2.13
0.65
20.65
10.32
0.77
0.88
11.02
10.52
0.65
0.65
2.39
10.36
3.45
5.18
441
0.26
f1 = 4.51
f3 = 6.78 βFORM Điểm D
20.65
10.01
1.15
0.65
16.12
11.90
0.67
0.66
10.03
7.09
0.65
3.11
2.75
7.95
8.16
11.28
436
0.26
f1 = 4.29
f3 = 6.33
β1= 8.92(100%)
β1= 2.25(98.77%) β1= 4.20(99.99%)
β2= 11.53(100%)
β2= 11.88(100%) β2= 11.96(100%)
β3= 0.00(42.91%) β3= 0.00(48.26%)
β3= 1.77(96.13%)
β4= 0.00(43.11%) β4= 0.004(50.16%) β4= 1.30(90.24%) β1= 4.10(99.99%) β1= 8.87(100%)
β2= 14.19(100%) β2= 9.59(100%)
β3= 3.00(99.87%) β3= 3.05(99.88%)
β4= 3.00(99.87%) β4= 3.10(99.90%) 74 Tương tự như bài toán dàn 10 thanh, trong bài toán này tác giả cũng khảo sát sự ảnh hưởng của các mức độ tin cậy (an toàn) khác nhau lên kết cấu. Kết quả khảo sát được thể hiện qua Hình 3.16. Kết quả đạt được là 4 đường Pareto tương ứng với 4 giá trị β khác nhau (β =2; β =2.5; β =3; và β =3.5) với trục đứng (trục tung) thể hiện chuyển vị và trục ngang (trục hoành) thể hiện khối lượng của hệ dàn. Kết quả tại vị trí cực tiểu khối lượng và cực đại chuyển vị của hệ dàn được trình bày trong Bảng 3.13. Bảng 3.13 cho thấy khi mức độ tin cậy yêu cầu càng tăng (hay chỉ số độ tin cậy tăng) thì khối lượng của kết cấu cũng tăng theo. Bảng 3.13: Kết quả thiết kế RBDO của bài toán kết cấu 72 thanh ở các mức độ an toàn khác nhau. βFORM β = 2
(Điểm A)
19.27
9.47
1.82
0.65
15.45
9.52
0.65
0.65
9.14
9.58
0.65
0.65
4.35
8.89
0.65
0.82
383.28
0.44
3.61
14.02
2.01
2.01 β = 2.5
(Điểm B)
20.59
10.03
1.91
0.65
15.73
9.88
0.65
0.65
9.31
9.37
0.65
0.65
4.31
9.59
0.65
0.74
396.86
0.43
3.76
13.68
2.50
2.50 β = 3
(Điểm C)
20.51
9.57
2.13
0.65
16.86
10.20
0.65
0.65
10.44
10.65
0.65
0.65
4.11
10.02
0.65
0.75
412.85
0.42
4.10
14.19
3.00
3.00 β = 3.5
(Điểm D)
20.64
10.52
2.21
0.65
17.45
10.20
0.65
0.65
11.14
10.41
0.65
0.65
4.69
10.82
0.65
0.77
428.00
0.41
4.36
14.59
3.50
3.50 Biến thiết kế
diện tích (cm2)
A1
A2
A3
A4
A5
A6
A7
A8
A9
A10
A11
A12
A13
A14
A15
A16
Khối lượng (kg)
Chuyển vị (cm)
β1
β2
β3
β4 75 Hình 3.16 cho thấy khi giá trị β càng lớn thì đường Pareto càng ngắn lại và dịch chuyển sang bên phải. Điều này sẽ làm khả năng chịu lực của kết cấu càng cao (hay độ an toàn của kết cấu sẽ càng cao), tương ứng với diện tích của các thanh dàn càng lớn, khối lượng của kết cấu càng tăng và chuyển vị của kết cấu càng nhỏ. Vì vậy, tuỳ thuộc vào tầm quan trọng của kết cấu khi thiết kế mà người thiết kế nên đưa ra những mức độ an toàn phù hợp nhằm giảm chi phí và tránh sự lãng phí không cần thiết. Hình 3.16: Kết quả thiết kế RBDO của bài toán kết cấu 72 thanh ở các mức độ an toàn khác nhau. Bài toán kết cấu dàn phẳng 200 thanh Bài toán này được tối ưu với mục đích đánh giá hiệu quả của việc áp dụng giải thuật NSGA-II và RBDO-NSGA-II cho kết cấu phức tạp. Xét kết cấu dàn phẳng 200 thanh như Hình 3.17. Kết cấu có 200 thanh và được chia làm 29 nhóm tương ứng với 29 biến thiết thế như Bảng 3.14. Bài toán thiết kế DO (Dertemined Optimization) và RBDO (Reliability Based-Design Optimization) nhằm cực tiểu khối lượng và cực tiểu chuyển vị của toàn hệ kết cấu. Kết cấu được thiết kế cho ba trường hợp tải khác nhau [19] [52]–[54]: 76 Trường hợp 1: kết cấu chịu tải trọng 4.45 kN theo phương x lên nút các 1, 6, 15, 20, 29, 43, 48, 57, 62, và 71. Trường hợp 2: kết cấu chịu tải trọng 44.5 kN theo phương y lên các nút 1, 2, 3, 4, 5, 6, 8, 10, 12, 14, 15, 16, 17, 18, 19, 20, 22, 24,. . . , 71, 72, 73, 74, và 75. Trường hợp 3: kết cấu chịu tải trọng là tổ hợp của hai trường hợp trên. Hình 3.17: Mô hình kết cấu dàn phẳng 200 thanh. 77 Bảng 3.14: Nhóm thanh diện tích trong kết cấu dàn 200 thanh. Nhóm thanh Nhóm thanh Biến diện
tích 1 1 2 3 4 82 83 85 86 88 89 91 92 103 104 Biến
diện tích
16 106 107 109 110 112 113 2 5 8 11 14 17 17 115 116 117 118 3 19 20 21 22 23 24 18 119 122 125 128 131 4 18 25 56 63 94 101 132 139 19 133 134 135 136 137 138 170 177 5 26 29 32 35 38 20 140 143 146 149 152 6 6 7 9 10 12 13 15 16 27 28 30 21 120 121 123 124 126 127 129 31 33 34 36 37 130 141 142 144 145 147 148 150 151 7 153 154 155 156 39 40 41 42 22 8 157 160 163 166 169 43 46 49 52 55 23 9 171 172 173 174 175 176 57 58 59 60 61 62 24 10 178 181 184 187 190 64 67 70 73 76 25 11 44 45 47 48 50 51 53 54 65 66 26 158 159 161 162 164 165 167 68 69 71 72 74 75 168 179 180 182 183 185 186 188 189 12 191 192 193 194 77 78 79 80 27 13 195 197 198 20 81 84 87 90 93 28 14 196 199 95 96 97 98 99 100 29 15 102 105 108 111 114 78 Bảng 3.15: Thông số bài toán kết cấu dàn phẳng 200 thanh. Thông số Giá trị Mô-đun đàn hồi E 2.06x1011 N/m2 Khối lượng riêng ρ 7933.410 kg/m3 100 kg Khối lượng bổ sung W0 = 15 Hz Giới hạn về tần số [9]–[11]; [49] = 5 Hz, = 10 Hz, Giới hạn về chuyển vị [1];[28] 1.5 cm Giới hạn ứng suất [19]; [52]–[54] 68.95x106 N/m2 Giới hạn biến thiết kế 0.6452 cm2 ≤ Ai ≤ 129.03 cm2 Thông số bài toán được cho như trong Bảng 3.15. Tương tự như bài toán dàn 10 thanh và bài toán dàn 72 thanh, bài toán tối ưu dàn 200 thanh được tính toán lần lượt qua 2 trường hợp: Không xét đến độ tin cậy (DO - Dertemined Optimization) Trong trường hợp này hai mục tiêu của bài toán tối ưu DO là cực tiểu khối lượng và cực tiểu chuyển vị. Mô hình toán học của bài toán DO được cho như Hình 3.17. Trong đó diện tích mặt cắt ngang của các thanh dàn là biến thiết kế. Giá trị cận dưới và cận trên của biến thiết, mô-đun đàn hồi E, khối lượng riêng của vật liệu ρ, khối lượng bổ sung tại các nút W0 (từ nút 1 đến nút 5), giới hạn chuyển vị, ứng suất và tần số dao động được cho như trong Bảng 3.15. Bài toán này đã được nghiên cứu trước đó bởi nhiều tác giả như Lamberti [52]; Mustafa Sonmez [53]; Lee và Geem [54]; Kaveh và Talatahari [19] v.v. Các bài toán này đều là tối ưu hóa đơn mục tiêu, với mục tiêu là cực tiểu khối lượng chịu ràng buộc ứng suất; nghiên cứu của Poldee và cộng sự (2012, 2013) [1][28] nhằm tối ưu hóa đơn mục tiêu, với mục tiêu là cực tiểu khối lượng chịu ràng buộc chuyển vị; và các nghiên cứu của Farshchin [9]; Kaveh [10]; Khatibinia [11] v.v. nhằm tối ưu hóa đơn mục tiêu, với mục tiêu là cực tiểu khối lượng chịu ràng buộc tần số dao động. Tuy nhiên, cho đến nay các bài toán tối ưu hoá đa mục tiêu kết cấu dàn 200 thanh gồm cực tiểu khối lượng, 79 cực tiểu chuyển vị chịu ràng buộc chuyển vị, ứng suất và tần số dao động vẫn chưa được thực hiện. Cũng giống như bài toán tối ưu dàn phẳng 10 thanh và dàn không gian 72 thanh, trong trường hợp này bài toán cũng được tối ưu với 3 trường hợp ràng buộc khác nhau: + Ràng buộc chuyển vị và ứng suất (tĩnh học) với giới hạn chuyển vị của hệ là 1.5 cm. + Ràng buộc chuyển vị, ứng suất (tĩnh học) và tần số dao động (động học) với giới hạn chuyển vị của hệ là 1.5 cm. + Ràng buộc chuyển vị, ứng suất (tĩnh học) và tần số dao động (động học) với giới hạn chuyển vị của hệ là 3 cm. a) Không xét đến ràng buộc tần số dao động b) Có xét đến ràng buộc tần số dao động Hình 3.17: Mô hình bài toán tối ưu DO kết cấu dàn 200 thanh. Kết quả tối ưu của bài toán trong trường hợp này được thể hiện qua Hình 3.18. Kết quả cho thấy 3 đường Pareto tương ứng với 3 trường hợp ràng buộc khác nhau với trục đứng (trục tung) thể hiện chuyển vị và trục ngang (trục hoành) thể hiện khối lượng của hệ dàn. Tai vị trí cực tiểu khối lượng, kết quả tối ưu được trình bày trong Bảng 3.16 với điểm A là kết quả tối ưu trường hợp kết cấu chịu ràng buộc tĩnh học với giới hạn chuyển vị 1.5 cm, điểm B là kết quả tối ưu trường hợp kết cấu chịu ràng buộc tĩnh học và động học với giới hạn chuyển vị 1.5 cm, điểm C là kết quả tối ưu trường hợp kết cấu chịu ràng buộc tĩnh học và động học với giới hạn chuyển vị 3 cm. Ngoài ra các kết quả 80 trên cũng được so sánh với các kết quả nghiên cứu của Lamberti [52]; Mustafa Sonmez [53]; Lee và Geem [54]; Kaveh và Talatahari [19] để đánh giá hiệu quả của bài toán. Kết quả Bảng 3.16, cho thấy tại vị trí cực tiểu khối lượng, phương pháp NSGA- II cho kết quả tốt với giá trị là 15470.20 kg (điểm A và điểm B), kết quả này không tốt bằng kết quả nghiên cứu của Lamberti [52]; Mustafa Sonmez [53]; Kaveh và Talatahari [19]. Tuy nhiên khi xét về giới hạn chuyển vị thì chuyển vị các nghiên cứu tham khảo trên không thỏa |u| ≤ umax = 1.5 cm. Ngoài ra nghiên cứu này cũng đã khắc phục được nhược điểm của nghiên cứu của Poldee và cộng sự (2012, 2013) [1][28] khi chỉ xét ràng buộc chuyển vị; và nhược điểm của nghiên cứu của Farshchin [9]; Kaveh [10]; Khatibinia [11] khi chỉ xét ràng buộc tần số dao động. Hình 3.18: Kết quả tối ưu hóa DO kết cấu dàn 200 thanh. Kết quả Bảng 3.16 cũng cho thấy khi tăng giá trị ràng buộc chuyển vị lên (từ 1.5 cm lên 3 cm) thì giá trị khối lượng và chuyển vị của hệ không thay đổi (điểm B và điểm C). Điều đó có nghĩa rằng việc tăng giá trị chuyển vị không đồng nghĩa với khối lượng của kết cấu nhỏ đi một cách tương ứng vì trong bài toán này ta cũng xét đến ba ràng buộc. Do đó kết quả này một lần nữa cho thấy rằng một bài toán khi xét nhiều ràng buộc 81 sẽ cho kết quả càng chính xác, các ràng buộc có ảnh hưởng qua lại với nhau và giải thuật NSGA-II tối ưu được cho các kết cấu phức tạp. Kết quả Hình 3.18 một lần nữa cho thấy rằng kết quả tối ưu là một tập hợp của nhiều lời giải tối ưu. Trong đó nếu lựa chọn kết cấu có khối lượng nhỏ thì chuyển vị của kết cấu sẽ lớn và ngược lại. Hình 3.18 cho thấy kết quả tối ưu là 3 đường Pareto. Trong đó cả ba đường Pareto P1, Pareto P2 và Pareto P3 trùng lên nhau. Điều này một lần nữa cho thấy sự tương quan giữa chuyển vị và khối lượng kết cấu, sự ảnh hưởng của các ràng buộc lên bài toán tối ưu. Từ đó chúng ta có thể thiết lập và giải một bài toán tối ưu hóa đa mục tiêu cho kết cấu, và tùy theo các điều kiện đầu ra mà lựa chọn một nghiệm thích hợp nằm trên đường Pareto. Tất cả các nghiệm nằm trên đường Pareto đều phải thỏa tất cả các điều kiện ràng buộc. 82 Bảng 3.16: Kết quả tối ưu hóa DO kết cấu dàn 200 thanh. Mustafa
Sonmez [53] Luận văn
DO-NSGA-II Lamberti
[52]
CMLPSA
0.95
6.06
0.65
0.65
12.52
1.91
0.65
20.03
0.65
26.48
2.60
1.23
35.02
0.65
41.48
3.70
0.86
51.43
0.65
57.89
4.55
2.71
70.09
0.65
76.52 ABC-AP 0.67
6.11
0.67
0.73
12.59
1.89
0.69
20.16
0.69
26.64
2.74
0.67
35.36
0.68
41.84
3.61
1.18
51.90
0.66
58.28
5.06
4.84
72.94
1.42
79.19 Lee và
Geem [54]
PSO
5.17
15.50
28.01
36.76
25.51
3.84
36.18
59.33
29.12
29.69
3.58
120.98
38.67
0.65
52.62
1.75
71.95
45.98
28.81
59.13
17.82
3.58
104.29
3.21
104.68 Kaveh và
Talatahari [19]
HPSACO
0.67
5.93
0.78
0.65
12.04
1.82
0.65
19.15
0.65
25.46
2.41
2.90
32.00
6.93
38.57
5.07
4.76
47.62
4.31
53.55
7.72
6.45
69.85
0.65
75.47 Điểm A(*) Điểm B(**) Điểm C(***)
1.23
7.33
0.65
0.76
17.00
0.65
2.03
25.71
0.65
28.73
5.51
4.21
44.61
2.13
51.60
7.26
2.04
65.39
0.83
75.99
7.29
2.00
88.54
13.90
94.50 1.23
7.33
0.65
0.76
17.00
0.65
2.03
25.71
0.65
28.73
5.51
4.21
44.61
2.13
51.60
7.26
2.04
65.39
0.83
75.99
7.29
2.00
88.54
13.90
94.50 1.23
7.33
0.65
0.76
17.00
0.65
2.03
25.71
0.65
28.73
5.51
4.21
44.61
2.13
51.60
7.26
2.04
65.39
0.83
75.99
7.29
2.00
88.54
13.90
94.50 Biến thiết kế diện
tích (cm2)
A1
A2
A3
A4
A5
A6
A7
A8
A9
A10
A11
A12
A13
A14
A15
A16
A17
A18
A19
A20
A21
A22
A23
A24
A25 83 Mustafa
Sonmez [53] Luận văn
DO-NSGA-II 15.10
43.89
86.22
116.55 Biến thiết kế diện
tích (cm2)
A26
A27
A28
A29
Khối lượng (kg)
Chuyển vị (cm) Tần số ABC-AP Lamberti
[52]
CMLPSA
6.67
43.11
69.75
89.32
11542.92
1.93
f1 = 5.7
f2 =15.8
f3 =20.6 9.07
33.29
64.47
94.85
11582.05
1.93
f1 = 5.7
f2 =16.0
f3 =20.8 Lee và
Geem [54]
PSO
6.48
23.29
53.99
100.41
19995.19
1.57
f1 = 3.38
f2 =11.3
f3 =17.5 Kaveh và
Talatahari [19]
HPSACO
8.85
31.95
56.77
94.61
11410.914
2.00
f1 = 5.8
f2 =15.7
f3 =21.5 Điểm A(*) Điểm B(**) Điểm C(***)
15.10
43.89
86.22
116.55
15470.20
15470.20
1.497
f1 = 6.11
f2 =13.87
f3 =22.19 15.10
43.89
86.22
116.55
15470.20
1.497
f1 = 6.11
f2 =13.87
f3 =22.19 1.497
f1 = 6.11
f2 =13.87
f3 =22.19 (**): Giới hạn chuyển vị 1.5 cm, có xét ràng buộc tần số. (***): Giới hạn chuyển vị 3.0 cm, có xét ràng buộc tần số. Chú thích: (*): Giới hạn chuyển vị 1.5 cm, không xét ràng buộc tần số. 84 Có xét đến độ tin cậy RBDO (Reliability Based-Design Optimization) Trong trường hợp này bài toán thiết kế RBDO nhằm cực tiểu khối lượng, cực tiểu chuyển vị của toàn hệ kết cấu như được phát biểu trong Hình 3.19. Mô-đun đàn hồi E, lực tác dụng P, khối lượng bổ sung tại các nút W0 (từ nút 1 đến nút 5) và khối lượng riêng của vật liệu ρ được xem là những biến ngẫu nhiên với giá trị trung bình được cho trong Bảng 3.15. Tất cả các biến ngẫu nhiên tuân theo luật phân phối Gauss (normal distribution) và được thống kê độc lập với hệ số thay đổi là 0.05. Giá trị giới hạn chuyển vị, ứng suất và tần số dao động được cho trong Bảng 3.15. Chỉ số độ tin cậy yêu cầu là βt = 2.5, tương ứng với xác suất an toàn là 99.39%. Diện tích mặt cắt ngang của các thanh dàn vừa là biến thiết kế vừa là biến ngẫu nhiên. Trong trường hợp này bài toán được tính toán với ràng buộc chuyển vị, ứng suất (tĩnh học) và tần số dao động (động học) với giới hạn chuyển vị của hệ là 1.5 cm. Hình 3.19: Mô hình bài toán tối ưu RBDO kết cấu dàn 200 thanh. Dựa trên thông tin được thu thập cho đến nay của tác giả, bài toán thiết kế RBDO cho kết cấu này cũng chưa được thực hiện Giống như kết quả tối ưu DO, kết quả tối ưu bài toán này cũng là một đường Pareto như trong Hình 3.20 với trục đứng (trục tung) thể hiện chuyển vị và trục ngang (trục hoành) thể hiện khối lượng của hệ dàn. Để đánh giá hiệu quả của bài toán trong phần này, kết quả tối ưu DO-NSGA-II trong Mục 3.2.3.1 ứng với ràng buộc chuyển vị, ứng suất (tĩnh học) và tần số dao động (động học) với giới hạn chuyển vị của hệ là 1.5 85 cm cũng được xem xét. Để đánh giá thêm sự khác biệt giữa bài toán DO và RBDO, một số nghiệm nằm trên 2 đường Pareto P1 và Pareto P2 (Hình 3.20) được sử dụng để so sánh. Kết quả được thể hiện như trong Bảng 3.17, trong đó điểm A và điểm B là hai nghiệm thuộc Pareto P1 và điểm C nghiệm thuộc Pareto P2. Kết quả Bảng 3.17 cho thấy các yếu tố ngẫu nhiên ảnh hưởng lớn đến kết quả tối ưu DO-NSGA-II. Cụ thể tại điểm A, giá trị β lần lượt là β1= 0.06 tương ứng 52.51%, β2= 0.05 tương ứng 52.05%. Điều đó có nghĩa là kết cấu chỉ đảm bảo an toàn ở mức 52.05%. Kết quả Bảng 3.17 cũng cho thấy rằng, khi xét tại cùng một vị trí khối lượng tối ưu thì độ an toàn của kết cấu khi xét đến chỉ số độ tin cậy luôn lớn hơn so với kết cấu không xét đến độ tin cậy. Cụ thể tại điểm B giá trị β nhỏ nhất là 2.26 tương ứng với 98.84 % còn tại điểm C giá trị β nhỏ nhất là 2.5 tương ứng với 99.39%. Điều đó có nghĩa khi sử dụng phương pháp RBDO-NSGA-II để tính toán tối ưu, kết cấu luôn đảm bảo mức độ an toàn cho trước là 99.39%. Bảng 3.17: Kết quả tối ưu hóa RBDO kết cấu dàn 200 thanh. Luận văn DO-NSGA-II Biến thiết kế
diện tích (cm2) Điểm A
1.23
7.33
0.65
0.76
17.00
0.65
2.03
25.71
0.65
28.73
5.51
4.21
44.61 RBDO-NSGA-II
β= 2.5
Điểm C
3.14
19.89
8.51
2.31
47.34
0.65
0.65
65.18
6.83
57.71
10.06
6.05
72.88 Điểm B
3.51
24.68
0.66
3.71
63.71
2.90
6.64
32.38
0.65
69.04
6.70
0.65
69.19 A1
A2
A3
A4
A5
A6
A7
A8
A9
A10
A11
A12
A13 86 Luận văn DO-NSGA-II Biến thiết kế
diện tích (cm2) A14
A15
A16
A17
A18
A19
A20
A21
A22
A23
A24
A25
A26
A27
A28
A29
Khối lượng (kg)
Chuyển vị (cm) Tần số Điểm A
2.13
51.60
7.26
2.04
65.39
0.83
75.99
7.29
2.00
88.54
13.90
94.50
15.10
43.89
86.22
116.55
15470.20
1.497
f1 = 6.11
f2 =13.87
f3 =22.19 RBDO-NSGA-II
β= 2.5
Điểm C
0.65
80.12
14.12
0.74
80.79
0.65
95.29
12.45
1.98
100.91
5.18
108.96
14.46
81.92
128.08
128.99
22270.1
0.98
f1 = 5.55
f2 =12.56
f3 =23.23 βFORM Điểm B
1.72
85.40
7.08
1.58
116.46
0.79
107.39
10.25
0.65
107.03
6.77
113.76
11.74
97.07
129.03
128.12
22528.9
0.96
f1 = 5.41
f2 =16.38
f3 =21.98
β1=0.06 (52.51%) β1=4.43 (99.99%) β1=4.21 (99.99%)
β2=0.05 (52.05%) β2=2.26 (98.84%) β2=2.50 (99.39%)
β3=2.38(99.14%) β3= 2.93(99.83%)
β3= 5.59 (100%)
β4= 5.75 (100%)
β4=11.85 (100%)
β4=8.32 (100%)
β5=10.69 (100%)
β5=8.62 (100%)
β5=9.56 (100%) Cũng như kết quả tối ưu hóa đa mục tiêu DO, bài toán tối ưu hóa đa mục tiêu RBDO cũng cho kết quả là một đường Pareto như Hình 3.20. Hình 3.20 cho thấy đường Pareto P2 nằm trùng với đường Pareto P1 tuy nhiên đường P2 ngắn hơn đường P1. Điều đó có nghĩa tập hợp các nghiệm của P1 là nghiệm tối ưu nhất, trong đó có cả nghiệm thỏa giá trị β và nghiệm không thỏa giá trị β. Còn khi 87 xét đến độ tin cậy thì tất cả các nghiệm của P2 là nghiệm tối ưu nhất, và đều thỏa giá trị β cho trước. Kết quả đạt được một lần nữa cho thấy tính logic hợp lý của bài toán, kết quả tối ưu dựa trên độ tin cậy RBDO có giá trị lớn hơn so với kết quả tối ưu không xét đến độ tin cậy DO. Sự gia tăng chi phí này là hiển nhiên nhằm đảm bảo độ an toàn nhất định cho kết cấu vận hành trong những điều kiện thiết kế ngẫu nhiên cụ thể. Hình 3.20: Kết quả tối ưu RBDO-NSGA-II kết cấu dàn 200 thanh. Tương tự như hai bài toán tối ưu trên, trong bài toán này tác giả cũng khảo sát sự ảnh hưởng của các mức độ tin cậy (an toàn) khác nhau lên kết cấu. Kết quả tối ưu được thể hiện qua Hình 3.21. Kết quả đạt được là 4 đường Pareto tương ứng với 4 giá trị β khác nhau (β =1.8; β =2; β =2.3; và β =2.5) với trục đứng (trục tung) thể hiện chuyển vị và trục ngang (trục hoành) thể hiện khối lượng của hệ dàn. Tai vị trí cực tiểu khối lượng của hệ dàn, kết quả được trình bày trong Bảng 3.18. Kết quả Bảng 3.18 cũng cho thấy khi mức độ tin cậy yêu cầu càng tăng (hay chỉ số độ tin cậy tăng) thì khối lượng của kết cấu cũng tăng theo. 88 Bảng 3.18: Kết quả thiết kế RBDO của bài toán kết cấu 200 thanh ở các mức độ an toàn khác nhau. β = 2.5
(Điểm D) β = 1.8
(Điểm A)
0.70 β = 2.0
(Điểm B)
1.32 β = 2.3
(Điểm C)
0.65 Biến thiết kế
diện tích (cm2)
A1 3.14 17.30 34.27 13.12 A2 19.89 0.65 0.65 13.43 A3 8.51 1.81 1.12 1.30 A4 2.31 29.00 36.90 32.80 A5 47.34 1.00 0.65 0.65 A6 0.65 5.99 0.65 6.20 A7 0.65 27.54 40.59 41.88 A8 65.18 0.65 0.65 2.88 A9 6.83 36.73 71.71 73.87 A10 57.71 7.41 5.07 6.61 A11 10.06 1.54 18.34 0.65 A12 6.05 50.92 53.80 76.78 A13 72.88 0.65 1.68 0.65 A14 0.65 59.96 70.93 66.23 A15 80.12 8.50 8.69 7.56 A16 14.12 3.35 8.09 0.65 A17 0.74 73.49 72.73 106.86 A18 80.79 0.65 0.65 0.65 A19 0.65 85.23 79.77 96.67 A20 95.29 10.11 10.79 8.61 A21 12.45 5.15 2.35 5.73 A22 1.98 100.13 98.51 100.90 A23 100.91 0.65 0.65 0.65 A24 5.18 117.61 108.48 107.55 A25 108.96 89 β = 2.5
(Điểm D) β = 1.8
(Điểm A)
9.44 β = 2.0
(Điểm B)
10.57 β = 2.3
(Điểm C)
15.68 Biến thiết kế
diện tích (cm2)
A26 14.46 70.79 101.79 89.13 A27 81.92 126.23 115.57 118.93 A28 128.08 126.01 127.74 129.03 A29 128.99 Khối lượng (kg) 18903.2 20397.2 21086.5 22270.1 Chuyển vị (cm) 1.18 1.07 1.04 2.34 3.38 3.56 0.98
4.21 β1 1.90 2.06 2.44 2.50 β2 βFORM 5.27 2.36 3.12 2.93 β3 10.08 4.93 6.57 5.75 β4 6.55 5.77 7.84 10.69 β 5 Hình 3.21 cũng cho thấy 4 đường Pareto tương ứng với 4 giá trị β khác nhau, cả 4 đường Pareto này xấp xỉ nằm chồng lên nhau, tuy nhiên giá trị β càng lớn thì đường Pareto càng ngắn lại. Điều này một lần nữa cho thấy tính logic hợp lý của bài toán RBDO, khi mức độ tin cậy yêu cầu càng tăng (hay chỉ số độ tin cậy tăng) thì khối lượng của kết cấu cũng tăng theo, và chuyển vị của kết cấu càng nhỏ lại. 90 Hình 3.21: Kết quả thiết kế RBDO của bài toán kết cấu 200 thanh ở các mức độ an toàn khác nhau. Nhận xét: Qua ba bài toán tối ưu hóa đa mục tiêu kết cấu dàn, các kết quả tối ưu cho thấy: - Một bài toán tối ưu khi xét đến nhiều ràng buộc sẽ cho kết quả càng chính xác và các ràng buộc này có ảnh hưởng qua lại với nhau. - Kết quả tối ưu RBDO luôn cho kết quả có giá trị lớn hơn kết quả tối ưu DO, nhằm đảm kết cấu vận hành theo mức độ an toàn cho trước (trong các điều kiện ngẫu nhiên cho trước). - Giải thuật NSGA-II có thể áp dụng để giải tất cả các bài toán tối ưu hóa đa mục tiêu dựa trên độ tin cậy kết cấu dàn, từ những kết cấu đơn giản đến kết cấu phức tạp. 91 4.1. Kết luận Luận văn được thực hiện nhằm hai mục tiêu chính: + Nghiên cứu áp dụng giải thuật NSGA-II trong việc giải các bài toán tối ưu đa mục tiêu kết cấu dàn nhằm khắc phục các nhược điểm của các nghiên cứu trước đây. Nghiên cứu xem xét đồng thời các ràng buộc tĩnh học và động học tác động lên kết cấu; + Áp dụng giải thuật NSGA-II để giải bài toán tối ưu hóa dựa trên độ tin cậy bằng cách kết hợp nó với giải thuật tối ưu hóa dựa trên độ tin cậy một vòng lặp đơn xác định với chỉ số độ tin cậy cho trước nhằm đảm bảo độ an toàn cho phép kết cấu dàn khi đưa vào vận hành sử dụng trong các điều kiện ngẫu nhiên cho trước. Để thực hiện hai mục tiêu trên, tác giả bước đầu tìm hiểu những ưu điểm cũng như nhược điểm của các phương pháp tối ưu hoá và tối ưu hoá dựa trên độ tin cậy, và cũng đã tìm hiểu những ưu điểm cũng như nhược điểm của các phương pháp tối ưu hoá đa mục tiêu so với tối ưu hóa đơn mục tiêu. Khảo sát, thu thập các nghiên cứu và áp dụng tối ưu hóa cho kết cấu dàn. Chi tiết nội dung của các cơ sở lý thuyết liên quan nhằm phục vụ cho việc thực hiện đề tài như: tổng quan về bài toán tối ưu hoá, phương pháp phần tử hữu hạn cho kết cấu dàn, giải thuật di truyền phân loại không trội NSGA-II, giải thuật tối ưu hoá dựa trên độ tin cậy một vòng lặp đơn xác định SLDM và giải thuật tối ưu hoá dựa trên độ tin cậy một vòng lặp đơn xác định sử dụng giải thuật di truyền phân loại không trội RBDO-NSGA-II cũng đã được trình bày. Và cuối cùng độ tin cậy của PP- PTHH, giải thuật NSGA-II và RBDO-NSGA-II được kiểm chứng thông qua 4 ví dụ số. Hiệu quả của việc áp dụng phương pháp NSGA-II và RBDO-NSGA-II cho kết cấu dàn được thực hiện qua 3 bài toán: dàn phẳng 10 thanh, dàn không gian 72 thanh và dàn phẳng 200 thanh. Dựa trên các kết quả đạt được, luận văn rút ra một số kết luận như sau: 92 + Đối với mục tiêu thứ nhất, Chương 2 đã trình bày chi tiết giải thuật tối ưu hoá đa mục tiêu NSGA-II và cách thành lập bài toán tối ưu hóa đa mục tiêu cho kết cấu dàn chịu ràng buộc tĩnh học và động học và Chương 3 đã trình bày các kết quả số liên quan. Kết quả số chỉ ra rằng giải thuật tối ưu hoá đa mục tiêu đề xuất trong luận văn là hiệu quả và cho kết quả đáng tin cậy. Cụ thể, cho cả 3 bài toán, phương pháp đề xuất cho kết quả tương đương so với các phương pháp khác đã được nghiên cứu trước đây tại vị trí so sánh. Ngoài ra, các kết quả đạt được ở Chương 3 còn khắc phục các nhược điểm sau của các nghiên cứu trước đây: Đã kết hợp được các điều kiện ràng buộc tĩnh (chuyển vị và ứng suất) và động (tần số dao động) vào chung một bài toán tối ưu hóa, điều mà các nghiên cứu trước đây ở trong nước và trên thế giới chưa thực hiện; Khắc phục được một số nhược điểm của tối ưu hóa đơn mục tiêu. Cụ thể nếu như tối ưu hóa đơn mục tiêu chỉ cho một nghiệm tối ưu duy nhất, thì tối ưu đa mục tiêu không những cho ta nhiều nghiệm tối ưu mà giữa các nghiệm này còn có mối tương quan với nhau. Hơn nữa nếu dữ liệu đầu ra có sự thay đổi như thay đổi giá trị giới hạn chuyển vị, giới hạn ứng suất hay giới hạn tần số dao động v.v. thì đối với bài toán tối ưu hóa đơn mục tiêu, ta phải thiết lập mô hình và chạy lại bài toán để tìm ra kết quả, việc làm này sẽ tốn nhiều chi phí tính toán và thời gian. Còn đối với bài toán tối ưu hóa đa mục tiêu, ta có thể dựa vào dữ liệu đầu ra để đưa ra một mô hình tổng quát của bài toán, từ đó thực hiện giải bài toán tối ưu. Kết quả đạt được là một tập hợp nghiệm tối ưu trên đường Pareto, khi đó các nhà nghiên cứu sẽ có cái nhìn tổng quan về bài toán hơn và tùy theo các điều kiện đầu ra mà chọn một phương án tối ưu thích hợp nằm trên đường Pareto. + Đối với mục tiêu thứ hai, Chương 2 đã trình bày chi tiết của giải thuật tối ưu hoá dựa trên độ tin cậy một vòng lặp đơn xác định và kết hợp giữa nó với giải thuật NSGA-II và Chương 3 đã trình bày các kết quả số liên quan. Kết quả số chỉ ra rằng, giải thuật tối ưu hoá dựa trên độ tin cậy đề xuất trong luận văn là hiệu quả và cho kết quả 93 đáng tin cậy so với các nghiên cứu trước đây. Cụ thể, trong ví dụ tối ưu kết cấu dàn phẳng 10 thanh, phương pháp đề xuất trong luận văn luôn cho kết quả tương đương với các phương pháp được so sánh. Ngoài ra ví dụ này cũng đã khắc phục được nhược điểm của các nghiên cứu trước là không xét đến ràng buộc của tần số dao động. Thêm vào đó bài toán thiết kế RBDO cho kết cấu dàn không gian 72 thanh và dàn phẳng 200 thanh chưa được thực hiện. Các kết quả số này cho thấy phương pháp đề xuất có thể áp dụng để giải tất cả các bài toán tối ưu hóa đa mục tiêu dựa trên độ tin cậy của kết cấu dàn, từ những kết cấu đơn giản (10 thanh) đến kết cấu phức tạp (200 thanh), từ dàn phẳng (10 thanh, 200 thanh) đến dàn không gian (72 thanh). Như vậy, với những kết quả đạt được, luận văn đã đóng góp thêm một công cụ tính toán mới tương đối hiệu quả và đáng tin cậy cho việc giải bài toán tối ưu hoá đa mục tiêu dựa trên độ tin cậy cho kết cấu dàn. So với các nghiên cứu trước đây, các nghiên cứu trong luận văn có những ưu điểm: + Kết hợp được các điều kiện ràng buộc tĩnh (chuyển vị và ứng suất) và ràng buộc động (tần số dao động) vào một bài toán tối ưu hóa; + Xét đến hệ số tin cậy trong quá trình tối ưu nhằm đảm bảo những giải pháp tối ưu tìm được không những thỏa mãn các mục tiêu và ràng buộc tối ưu mà còn thỏa mãn điều kiện an toàn trong những điều kiện thiết kế ngẫu nhiên cụ thể; + Kết quả tối ưu là một tập hợp các nghiệm tối ưu. Kết quả này giúp cho người thiết kế có cái nhìn tổng quát và có nhiều sự lựa chọn tùy theo mục đích sử dụng của mình; + Thuật toán đề xuất dễ dàng áp dụng cho tất cả các bài toán tối ưu hoá kết cấu dàn. Vì vậy nghiên cứu của luận văn không chỉ mang ý nghĩa về mặt lý thuyết mà nó còn mang tính ứng dụng cao. Tuy nhiên luận văn thực hiện vẫn còn một số điểm hạn chế: 94 + Giải thuật đề xuất có thời gian tính toán vẫn còn khá lâu khi giải các bài toán có độ phức tạp cao (kết cấu dàn nhiều thanh); + Giải thuật còn giới hạn cho các bài toán kết cấu có biến thiết kế liên tục, chứ chưa mở rộng cho các bài toán có biến thiết kế rời rạc. 4.2. Hướng phát triển Mặc dù đã nỗ lực rất nhiều trong suốt quá trình thực hiện luận văn, tuy nhiên luận văn khó có thể tránh khỏi những thiếu sót. Vì vậy tác giả mong nhận được những đóng góp quý giá từ quý thầy cô nhằm bổ sung những kiến thức còn thiếu sót để luận văn được hoàn chỉnh hơn. Bên cạnh đó, dựa vào những kết quả đã đạt được, tác giả cũng mong đề tài nhận được sự quan tâm và phát triển hơn nữa từ các nhà nghiên cứu. Luận văn có thể được mở rộng theo các hướng nghiên cứu sau: + Phát triển thuật toán tối ưu hoá dựa trên độ tin cậy (RBDO-NSGA-II) cho các loại kết cấu khác như: kết cấu khung, tấm, vỏ, kết cấu vật liệu composite, kết cấu dầm v.v.; + Giải thuật đề xuất cũng có thể được mở rộng để giải cho nhiều loại bài toán tối ưu trong nhiều lĩnh vực khác nhau như: các bài toán tối ưu hoá trong lĩnh vực kinh tế, các bài toán tối ưu hoá trong lĩnh vực điều khiển, v.v.; + Hoàn thiện và phát triển giải thuật để có thể giải quyết được các bài toán kết cấu có biến thiết kế rời rạc. 95 [1] D. G. Dựng, N. T. Liêm, H. H. Vịnh, V. D. Trung, and N. T. Trung (2016). Tối ưu hóa đa mục tiêu kết cấu dầm composites nhiều lớp chịu ràng buộc về tần số. Bài báo được trình bày tại Hội nghị Khoa học toàn quốc Vật liệu và Kết cấu Composite Cơ học, Công nghệ và ứng dụng, 28-29/7/2016, Đại học Nha Trang- TP Nha Trang. [2] T. Vo-Duy, D. Duong-Gia, V. Ho-Huu, H. C. Vu-Do, and T. Nguyen-Thoi (2017). “Multi-objective optimization of laminated composite beam structures using NSGA-II algorithm”. Compos. Struct., vol. 168, pp. 498–509. http://dx.doi.org/10.1016/j.compstruct.2017.02.038 (ISI) 96 [1] N. Pholdee and S. Bureerat (2012). “Performance enhancement of multiobjective evolutionary optimisers for truss design using an approximate gradient”. Comput. Struct., vol. 106, pp. 115–124. [2] K. Deb, A. Pratap, S. Agarwal, and T. Meyarivan (2002). “A fast and elitist multiobjective genetic algorithm: NSGA-II”. IEEE Trans. Evol. Comput., vol. 6, no. 2, pp. 182–197. [3] R. Storn and K. Price (1997). “Differential evolution–a simple and efficient heuristic for global optimization over continuous spaces”. J. Glob. Optim., vol. 11, no. 4, pp. 341–359. [4] M. Dorigo, M. Birattari, and T. Stutzle (2006). “Ant colony optimization”. IEEE Comput. Intell. Mag., vol. 1, no. 4, pp. 28–39. [5] J. Kennedy and R. C. Eberhart (1997). “A discrete binary version of the particle swarm algorithm”. Systems, Man, and Cybernetics, 1997. Computational Cybernetics and Simulation., 1997 IEEE International Conference on, vol. 5. IEEE, pp. 4104–4108. [6] Q. Zhang and H. Li (2007). “MOEA/D: A multiobjective evolutionary algorithm based on decomposition”. IEEE Trans. Evol. Comput., vol. 11, no. 6, pp. 712–731. [7] F. Li, T. Wu, A. Badiru, M. Hu, and S. Soni (2013). “A single-loop deterministic method for reliability-based design optimization”. Eng. Optim., vol. 45, no. 4, pp. 435–458. [8] S. Shan and G. G. Wang (2008). “Reliable design space and complete single-loop reliability-based design optimization”. Reliab. Eng. Syst. Saf., vol. 93, no. 8, pp. 1218–1230. [9] M. Farshchin, C. V Camp, and M. Maniat (2016). “Multi-class teaching–learning- based optimization for truss design with frequency constraints”. Eng. Struct., vol. 97 106, pp. 355–369. [10] A. Kaveh and M. I. Ghazaan (2015) “Hybridized optimization algorithms for design of trusses with multiple natural frequency constraints”. Adv. Eng. Softw., vol. 79, pp. 137–147. [11] M. Khatibinia and S. Sadegh Naseralavi (2014). “Truss optimization on shape and sizing with frequency constraints based on orthogonal multi-gravitational search algorithm”. J. Sound Vib., vol. 333, no. 24, pp. 6349–6369. [12] S. Pezeshk, C. V Camp, and D. Chen (2000). “Design of nonlinear framed structures using genetic optimization”. J. Struct. Eng., vol. 126, no. 3, pp. 382– 388. [13] C. V Camp, B. J. Bichon, and S. P. Stovall (2005). “Design of steel frames using ant colony optimization”. J. Struct. Eng., vol. 131, no. 3, pp. 369–379. [14] S. O. Degertekin (2008). “Optimum design of steel frames using harmony search algorithm”. Struct. Multidiscip. Optim., vol. 36, no. 4, pp. 393–401. [15] M. P. Saka (2009). “Optimum design of steel sway frames to BS5950 using harmony search algorithm”. J. Constr. Steel Res., vol. 65, no. 1, pp. 36–43. [16] A. Kaveh and S. Talatahari (2010). “An improved ant colony optimization for the design of planar steel frames”. Eng. Struct., vol. 32, no. 3, pp. 864–873. [17] A. Kaveh and H. Abbasgholiha (2011). “Optimum design of steel sway frames using Big Bang-Big Crunch algorithm”. Asian J. Civ. Eng., vol. 12, no. 3, pp. 293– 317. [18] E. Doğan and M. P. Saka (2012). “Optimum design of unbraced steel frames to LRFD–AISC using particle swarm optimization”. Adv. Eng. Softw., vol. 46, no. 1, pp. 27–34. [19] A. Kaveh and S. Talatahari (2009). “Particle swarm optimizer, ant colony strategy and harmony search scheme hybridized for optimization of truss structures”. Comput. Struct., vol. 87, no. 5–6, pp. 267–283. 98 [20] L. J. Li, Z. B. Huang, F. Liu, and Q. H. Wu (2007). “A heuristic particle swarm optimizer for optimization of pin connected structures”. Comput. Struct., vol. 85, no. 7–8, pp. 340–349. [21] K. S. Lee and Z. W. Geem (2004). “A new structural optimization method based on the harmony search algorithm”. Comput. Struct., vol. 82, no. 9–10, pp. 781– 798. [22] C. Camp, S. Pezeshk, and G. Cao (1998). “Optimized design of two-dimensional structures using a genetic algorithm”. J. Struct. Eng., vol. 124, no. 5, pp. 551–559. [23] T. Vo-Duy, D. Duong-Gia, V. Ho-Huu, H. C. Vu-Do, and T. Nguyen-Thoi (2017). “Multi-objective optimization of laminated composite beam structures using NSGA-II algorithm”. Compos. Struct., vol. 168, pp. 498–509. [24] C. A. Coello and A. D. Christiansen (2000). “Multiobjective optimization of trusses using genetic algorithms”. Comput. Struct., vol. 75, no. 6, pp. 647–660. [25] S. N. Omkar, D. Mudigere, G. N. Naik, and S. Gopalakrishnan (2008). “Vector evaluated particle swarm optimization (VEPSO) for multi-objective design optimization of composite structures”. Comput. Struct., vol. 86, no. 1, pp. 1–14. [26] S. Gurugubelli and D. Kallepalli (2014). “Weight and deflection optimization of Cantilever Beam using a modified Non-Dominated sorting Genetic Algorithm”. IOSR J. Eng., vol. 4, no. 3, pp. 19–23. [27] R. F. Coelho and P. Bouillard (2011). “Multi-objective reliability-based optimization with stochastic metamodels”. Evol. Comput., vol. 19, no. 4, pp. 525– 560. [28] N. Pholdee and S. Bureerat (2013). “Hybridisation of real-code population-based incremental learning and differential evolution for multiobjective design of trusses”. Inf. Sci. (Ny)., vol. 223, pp. 136–152. [29] J. S. Angelo, H. S. Bernardino, and H. J. C. Barbosa (2015). “Ant colony approaches for multiobjective structural optimization problems with a cardinality 99 constraint”. Adv. Eng. Softw., vol. 80, pp. 101–115. [30] C. K. Dimou and V. K. Koumousis (2009). “Reliability-based optimal design of truss structures using particle swarm optimization”. J. Comput. Civ. Eng., vol. 23, no. 2, pp. 100–109. [31] M. R. Ghasemi and M. Yousefi (2011). “Reliability-based optimization of steel frame structures using modified genetic algorithm”. Asian J. Civ. Eng., vol. 12, no. 4, pp. 449–475. [32] M. Shayanfar, R. Abbasnia, and A. Khodam (2014). “Development of a GA-based method for reliability-based optimization of structures with discrete and continuous design variables using OpenSees and Tcl”. Finite Elem. Anal. Des., vol. 90, pp. 61–73. [33] J.-O. Lee, Y.-S. Yang, and W.-S. Ruy (2002). “A comparative study on reliability- index and target-performance-based probabilistic structural design optimization”. Comput. Struct., vol. 80, no. 3, pp. 257–269. [34] Y. Su, H. Tang, S. Xue, and D. Li (2016). “Multi-objective differential evolution for truss design optimization with epistemic uncertainty”. Adv. Struct. Eng., vol. 19, no. 9, pp. 1403–1419. [35] L.-K. Song, C.-W. Fei, J. Wen, and G.-C. Bai (2017). “Multi-objective reliability- based design optimization approach of complex structure with multi-failure modes”. Aerosp. Sci. Technol., vol. 64, pp. 52–62. [36] T. V. Dần (2013). “Tối ưu hóa dựa trên độ tin cậy tấm composite laminate bằng giải thuật di truyền”. Luận văn (Thạc sĩ), Khoa Kỹ thuật xây dựng, Đại học Bách khoa TP.HCM. [37] T. B. Quốc (2015). “Tối ưu hóa kết cấu dàn với biến diện tích rời rạc sử dụng phương pháp DE cải tiến”. Luận văn (Thạc sĩ), Khoa xây dựng. Đại học Công nghệ TP.HCM. [38] N. T. Đ. Nhi (2014). “Tối ưu hóa đa mục tiêu Kết cấu tấm gấp composite Dùng 100 phương pháp NSGA – II”. Luận văn (Thạc sĩ), Khoa Kỹ thuật công trình, Đại học Tôn Đức Thắng. [39] H. H. Vịnh (2016). “Phát triển giải thuật tiến hóa DE cho bài toán tối ưu hóa dựa trên độ tin cậy với biến thiết kế rời rạc”. Luận văn (Thạc sĩ), Khoa khoa học Ứng dụng, Đại học Bách khoa TP.HCM. [40] D. G. Dựng, N. T. Liêm, H. H. Vịnh, V. D. Trung, and N. T. Trung (2016). Tối ưu hóa đa mục tiêu kết cấu dầm composites nhiều lớp chịu ràng buộc về tần số. Bài báo được trình bày tại Hội nghị Khoa học toàn quốc Vật liệu và Kết cấu Composite Cơ học, Công nghệ và ứng dụng, 28-29/7/2016, Đại học Nha Trang-TP Nha Trang. [41] B. T. V. Thái (2016). “Tối ưu hóa đa mục tiêu dựa trên độ tin cậy cho kết cấu khung thép sử dụng giải thuật một vòng lặp đơn”. Luận văn (Thạc sĩ), Khoa xây dựng dân dụng và công nghiệp, Đại học Mở TP.HCM. [42] L. T. P. Hạnh (2016). “Tối ưu hóa dựa trên độ tin cậy kết cấu khung thép sử dụng giải thuật một vòng lặp đơn xác định”. Luận văn (Thạc sĩ), Khoa xây dựng, Đại học Công nghệ TP.HCM. [43] Nguyễn Thời Trung and Nguyễn Xuân Hùng (2015). Phương pháp phần tử hữu hạn - Sử dụng Matlab. NXB Xây Dựng, Hà Nội. [44] X. Du and W. Chen (2004). “Sequential optimization and reliability assessment method for efficient probabilistic design”. J. Mech. Des., vol. 126, no. 2, pp. 225– 233. [45] N. Srinivas and K. Deb (1994). “Muiltiobjective optimization using nondominated sorting in genetic algorithms”. Evol. Comput., vol. 2, no. 3, pp. 221–248. [46] A. J. M. Ferrerai (2009). matlab_codes_for_finite_element_analysis-me. Springer, Berlin. [47] V. Ho-Huu, T. Nguyen-Thoi, T. Truong-Khac, L. Le-Anh, and T. Vo-Duy (2016). “An improved differential evolution based on roulette wheel selection for shape and size optimization of truss structures with frequency constraints”. Neural 101 Comput. Appl., pp. 1–19. [48] K. Deb, S. Gupta, D. Daum, J. Branke, A. K. Mall, and D. Padmanabhan (2009). “Reliability-based optimization using evolutionary algorithms”. IEEE Trans. Evol. Comput., vol. 13, no. 5, pp. 1054–1074. [49] A. Kaveh and A. Zolghadr (2012). “Truss optimization with natural frequency constraints using a hybridized CSS-BBBC algorithm with trap recognition capability”. Comput. Struct., vol. 102, pp. 14–27. [50] F. Erbatur, O. Hasançebi, I. Tütüncü, and H. Kılıç (2000). “Optimal design of planar and space structures with genetic algorithms”. Comput. Struct., vol. 75, no. 2, pp. 209–224. [51] R. E. Perez and K. Behdinan (2007). “Particle swarm approach for structural design optimization”. Comput. Struct., vol. 85, no. 19–20, pp. 1579–1588. [52] L. Lamberti (2008). “An efficient simulated annealing algorithm for design optimization of truss structures”. Comput. Struct., vol. 86, no. 19, pp. 1936–1953. [53] M. Sonmez (2011).“Artificial Bee Colony algorithm for optimization of truss structures”. Appl. Soft Comput., vol. 11, no. 2, pp. 2406–2418. [54] K. S. Lee and Z. W. Geem (2004). “A new structural optimization method based on the harmony search algorithm”. Comput. Struct., vol. 82, no. 9, pp. 781–798. Một số đoạn code Matlab cho bài toán tối ưu hóa kết cấu dàn 1. Code Matlab phân tích ứng xử kết cấu dàn phẳng 10 thanh
clc
clear all
close all
%-----------------------------------------------
E = 6.89e10 ;
A = [0.0203493
6.45E-05
0.017202
0.00889875
6.45E-05
6.45E-05
0.00528199
0.0132676
0.0130541
6.45E-05];
P = 444.82e3;
%%------------------------------------------------
% 1 2 3 4 5 6
gcoord = [18.288, 18.288, 9.144, 9.144, 0, 0
9.144, 0, 9.144, 0, 9.144, 0];
% 1 2 3 4 5 6 7 8 9 10
element = [5, 3, 6, 4, 4, 1, 4, 3, 2, 1
3, 1, 4, 2, 3, 2, 5, 6, 3, 4];
%------------------------------------------------
nel = length(element); % total element
nnode = length(gcoord); % total node
ndof = 2; % number of degree of freedom of one node
sdof = nnode*ndof; % total dgree of freedom of system
% plot_mesh_Spatial % plot model of problem
%**********************
% processing *
%**********************
K = zeros(sdof,sdof); % creat matrix of stiffness initial is zeros
F = zeros(sdof,1); % creat matrix of force is zeros
%---------------------------
% calculate stiffness matrix
for i=1:nel nd = element(:,i);
x = gcoord(1,nd); y = gcoord(2,nd);
% compute long of each bar
le = sqrt((x(2)-x(1))^2 + (y(2)-y(1))^2);
% compute direction cosin
l_ij = (x(2)-x(1))/le; % Eq.5.19
m_ij = (y(2)-y(1))/le; % Eq.5.19
% compute transform matrix
Te = [l_ij m_ij 0 0;
0 0 l_ij m_ij];
% compute stiffness matrix of element
ke = A(i)*E/le*[1 -1;
-1 1];
ke = Te'*ke*Te;
% compute force vector
%fe = A(i)*rho*le*9.81/2*[0 1, 0 1];
% find index assemble
index = [2*nd(1)-1 2*nd(1) 2*nd(2)-1 2*nd(2)];
% assemble ke in K
K(index,index) = K(index,index) + ke;
% assemble fe in F
%F(index,1) = F(index,1) + fe';
end
F = -F;
%*****************************
% boundary condition *
%*****************************
bcdof = [(5:6)*2-1, (5:6)*2]; % boundary condition displacement
bcval = zeros(1,length(bcdof)); % value of boundary
condition
%*********************
% Node | PX PY |
%---------------------
% 2 | 0 1e5 |
% 4 | 0 1e5 |
%---------------------
% Note: loads are in kips. 1 (kips) = 1000 (lb)
%------------------------
% consider to self weight
% Node 2
F([4 8]) = -P;
%*******************************************************
% apply boundary condition force and displacement *
%******************************************************* %--------------------------------------------------------------------------
% 1 2 3 4 5 6
gcoord = [18.288, 18.288, 9.144, 9.144, 0, 0
9.144, 0, 9.144, 0, 9.144, 0];
% 1 2 3 4 5 6 7 8 9 10
element = [3, 1, 4, 2, 3, 1, 4, 3, 2, 1
5, 3, 6, 4, 4, 2, 5, 6, 3, 4];
nel = length(element); % total element
nnode = length(gcoord); % total node
ndof = 2; % number of degree of freedom of one node
sdof = nnode*ndof; % total dgree of freedom of system
% plotModel( type,gcoord,element );
% calculate stiffness matrix
[ K,M ] = Cal_K_and_M(type,gcoord,element,A,rho,E);
% add non-structural mass
addedMass = 454; %kg
for idof = 1:sdof
M(idof,idof) = M(idof,idof) + addedMass;
end
% apply boundary
bcdof = [(5:6)*2-1, (5:6)*2]; % boundary condition displacement
% Giai phuong trinh tim tri rieng va vector rieng
[omega_2,~]=eigens(K,M,bcdof);
f=sqrt(omega_2)/2/pi;
c1 = 7/f(1) -1;
c2 = 15/f(2)-1;
c3 = 20/f(3)-1;
K = zeros(sdof,sdof); % creat matrix of stiffness initial is zeros
F = zeros(sdof,1); % creat matrix of force is zeros
%---------------------------
% calculate stiffness matrix
for i=1:nel
nd = element(:,i);
x = gcoord(1,nd); y = gcoord(2,nd);
% compute long of each bar
le = sqrt((x(2)-x(1))^2 + (y(2)-y(1))^2);
% compute direction cosin
l_ij = (x(2)-x(1))/le; % Eq.5.19
m_ij = (y(2)-y(1))/le; % Eq.5.19
% compute transform matrix
Te = [l_ij m_ij 0 0;
0 0 l_ij m_ij];
% compute stiffness matrix of element ke = A(i)*E/le*[1 -1;
-1 1];
ke = Te'*ke*Te;
% compute force vector
%fe = A(i)*rho*le*9.81/2*[0 1, 0 1];
% find index assemble
index = [2*nd(1)-1 2*nd(1) 2*nd(2)-1 2*nd(2)];
% assemble ke in K
K(index,index) = K(index,index) + ke;
% assemble fe in F
%F(index,1) = F(index,1) + fe';
end
F = -F;
%*****************************
% boundary condition *
%*****************************
bcdof = [(5:6)*2-1, (5:6)*2]; % boundary condition displacement
bcval = zeros(1,length(bcdof)); % value of boundary
condition
%*********************
% Node | PX PY |
%---------------------
% 2 | 0 1e5 |
% 4 | 0 1e5 |
%---------------------
% Note: loads are in kips. 1 (kips) = 1000 (lb)
%------------------------
% consider to self weight
% Node 2
F([4 8]) = -P;
%*******************************************************
% apply boundary condition force and displacement *
%*******************************************************
n = length(bcdof);
for i = 1:n
c = bcdof(i);
K(c,:) = 0;
K(c,c) = 1 ;
F(c) = bcval(i);
end
u = K\F;
max(abs(u));
% rang buoc chuyen vi
Cdisp = max(abs(u))/0.0508-1; % chuyen vi 5.08cm rho = 2767;
sigrho = 0.05*rho;
X = [E muyP addedMass rho xval]; % mean value of random variables
dX = 1e-6*X+1e-10;
sigma = [sigE sigP sigadd sigrho sigVal]; % standard deviation
beta = norminv(0.9987);
% % beta = 3;
[d11,d21,d31,d41,d51]=Truss10(X);
for i = 1:length(X)
xx = X;
xx(i)= xx(i)+dX(i);
[d12,d22,d32,d42,d52]=Truss10(xx);
G(1,i)= (d12-d11)/dX(i)*sigma(i);
G(2,i)= (d22-d21)/dX(i)*sigma(i);
G(3,i)= (d32-d31)/dX(i)*sigma(i);
G(4,i)= (d42-d41)/dX(i)*sigma(i);
G(5,i)= (d52-d51)/dX(i)*sigma(i);
end
x_add(1,:) = beta*sigma.*G(1,:)./norm(G(1,:));
x_new1 = X + x_add(1,:);
x_add(2,:) = beta*sigma.*G(2,:)./norm(G(2,:));
x_new2 = X + x_add(2,:);
x_add(3,:) = beta*sigma.*G(3,:)./norm(G(3,:));
x_new3 = X + x_add(3,:);
x_add(4,:) = beta*sigma.*G(4,:)./norm(G(4,:));
x_new4 = X + x_add(4,:);
x_add(5,:) = beta*sigma.*G(5,:)./norm(G(5,:));
x_new5 = X + x_add(5,:);
[c11,~,~,~,~] = Truss10(x_new1);
[~,c22,~,~,~] = Truss10(x_new2);
[~,~,c32,~,~] = Truss10(x_new3);
[~,~,~,c42,~] = Truss10(x_new4);
[~,~,~,~,c52] = Truss10(x_new5);
c = [c11 c22 c32 c42 c52];
ceq = [];
end
function [ f1, f2, f3, Cdisp, Csig] = Truss10( x ) type = '2D';
E = x(1); P = x(2); addedMass =x(3) ; rho=x(4) ; A = x(5:end);
%--------------------------------------------------------------------------
% 1 2 3 4 5 6
gcoord = [18.288, 18.288, 9.144, 9.144, 0, 0
9.144, 0, 9.144, 0, 9.144, 0];
% 1 2 3 4 5 6 7 8 9 10
element = [3, 1, 4, 2, 3, 1, 4, 3, 2, 1
5, 3, 6, 4, 4, 2, 5, 6, 3, 4];
nel = length(element); % total element
nnode = length(gcoord); % total node
ndof = 2; % number of degree of freedom of one node
sdof = nnode*ndof; % total dgree of freedom of system
% plotModel( type,gcoord,element );
% calculate stiffness matrix
[ K,M ] = Cal_K_and_M(type,gcoord,element,A,rho,E);
% add non-structural mass
for idof = 1:sdof
M(idof,idof) = M(idof,idof) + addedMass;
end
% apply boundary
bcdof = [(5:6)*2-1, (5:6)*2]; % boundary condition displacement
% Giai phuong trinh tim tri rieng va vector rieng
[omega_2,~]=eigens(K,M,bcdof);
f=sqrt(omega_2)/2/pi;
f1 = 7/f(1) -1;
f2 = 15/f(2)-1;
f3 = 20/f(3)-1;
K = zeros(sdof,sdof); % creat matrix of stiffness initial is zeros
F = zeros(sdof,1); % creat matrix of force is zeros
%---------------------------
% calculate stiffness matrix
for i=1:nel
nd = element(:,i);
x = gcoord(1,nd); y = gcoord(2,nd);
% compute long of each bar
le = sqrt((x(2)-x(1))^2 + (y(2)-y(1))^2);
% compute direction cosin
l_ij = (x(2)-x(1))/le; % Eq.5.19
m_ij = (y(2)-y(1))/le; % Eq.5.19
% compute transform matrix
Te = [l_ij m_ij 0 0;
0 0 l_ij m_ij];
% compute stiffness matrix of element ke = A(i)*E/le*[1 -1;
-1 1];
ke = Te'*ke*Te;
% compute force vector
%fe = A(i)*rho*le*9.81/2*[0 1, 0 1];
% find index assemble
index = [2*nd(1)-1 2*nd(1) 2*nd(2)-1 2*nd(2)];
% assemble ke in K
K(index,index) = K(index,index) + ke;
% assemble fe in F
%F(index,1) = F(index,1) + fe';
end
F = -F;
%*****************************
% boundary condition *
%*****************************
bcdof = [(5:6)*2-1, (5:6)*2]; % boundary condition displacement
bcval = zeros(1,length(bcdof)); % value of boundary
condition
%*********************
% Node | PX PY |
%---------------------
% 2 | 0 1e5 |
% 4 | 0 1e5 |
%---------------------
% Note: loads are in kips. 1 (kips) = 1000 (lb)
%------------------------
% consider to self weight
% Node 2
F([4 8]) = -P;
%*******************************************************
% apply boundary condition force and displacement *
%*******************************************************
n = length(bcdof);
for i = 1:n
c = bcdof(i);
K(c,:) = 0;
K(c,c) = 1 ;
F(c) = bcval(i);
end
u = K\F;
max(abs(u));
% rang buoc chuyen vi
%Cdisp = max(abs(u))/0.0508-1; %Chuyen vi 5.08cm nnode = length(gcoord); % total node
ndof = 2; % number of degree of freedom of one node
sdof = nnode*ndof; % total dgree of freedom of system
% plot_mesh_Spatial % plot model of problem
%**********************
% processing *
%**********************
K = zeros(sdof,sdof); % creat matrix of stiffness initial is zeros
F = zeros(sdof,1); % creat matrix of force is zeros
%---------------------------
% calculate stiffness matrix
for i=1:nel
nd = element(:,i);
x = gcoord(1,nd); y = gcoord(2,nd);
% compute long of each bar
le = sqrt((x(2)-x(1))^2 + (y(2)-y(1))^2);
% compute direction cosin
l_ij = (x(2)-x(1))/le; % Eq.5.19
m_ij = (y(2)-y(1))/le; % Eq.5.19
% compute transform matrix
Te = [l_ij m_ij 0 0;
0 0 l_ij m_ij];
% compute stiffness matrix of element
ke = A(i)*E/le*[1 -1;
-1 1];
ke = Te'*ke*Te;
% compute force vector
%fe = A(i)*rho*le*9.81/2*[0 1, 0 1];
% find index assemble
index = [2*nd(1)-1 2*nd(1) 2*nd(2)-1 2*nd(2)];
% assemble ke in K
K(index,index) = K(index,index) + ke;
% assemble fe in F
%F(index,1) = F(index,1) + fe';
end
F = -F;
%*****************************
% boundary condition *
%*****************************
bcdof = [(5:6)*2-1, (5:6)*2]; % boundary condition displacement
bcval = zeros(1,length(bcdof)); % value of boundary
condition
%*********************
% Node | PX PY | 6. Code Matlab giải bài toán tối ưu hóa dàn phẳng 10 thanh
clc; clear; close all
addpath NSGA-II-functions
% thong so thuat toan
options = nsgaopt(); % create default options structure
options.popsize = 50; % populaion size
options.maxGen = 5000; % tong so vong lap
options.numObj = 2; % so luong ham muc tieu
options.numVar = 10; % so luong bien thiet ke
options.numCons = 1;
nvars = 10;
options.lb = 0.645e-4*ones(1,nvars); % can duoi
options.ub = 225.806e-4*ones(1,nvars); % can tren
options.objfun = @TP_CONSTR_objfun % ham muc tieu
options.plotInterval = 5; % interval between two calls of
"plotnsga".
options.useParallel = 'no';
result = nsga2(options) % begin the optimization!
7. Code Matlab tính giá trị β
function beta = CheckReliability2(xval)
sigVal = 0.05*xval;
E = 6.89*1e10;
sigE = 0.05*E;
P = 444.82e3; V_P = 0.1;
siglogP = sqrt(log(1 + V_P^2));
muylogP = log(P)-1/2*siglogP^2;
sigP = P*siglogP;
muyP = P*(1 - log(P) + muylogP);
muy_X = [E muyP xval]; % mean value of random variables
sigma = [sigE sigP sigVal]; % standard deviation
[ beta1 ] = FORM_Reliability( muy_X,sigma,@f1 );
[ beta2 ] = FORM_Reliability( muy_X,sigma,@f2 );
beta = [beta1;beta2];
function [ Cdisp ] = f1(x)
% E = 6.89*1e10; % Young's elastic modulus (N/m^2)
% A = x; % P = 444.82e3;
E = x(1); P= x(2); A = x(3:end);
%--------------------------------------------------------------------------
% 1 2 3 4 5 6
gcoord = [18.288, 18.288, 9.144, 9.144, 0, 0
9.144, 0, 9.144, 0, 9.144, 0];
% 1 2 3 4 5 6 7 8 9 10
element = [3, 1, 4, 2, 3, 1, 4, 3, 2, 1
5, 3, 6, 4, 4, 2, 5, 6, 3, 4];
nel = length(element); % total element
nnode = length(gcoord); % total node
ndof = 2; % number of degree of freedom of one node
sdof = nnode*ndof; % total dgree of freedom of system
% plotModel( type,gcoord,element );
% calculate stiffness matrix
K = zeros(sdof,sdof); % creat matrix of stiffness initial is zeros
F = zeros(sdof,1); % creat matrix of force is zeros
%---------------------------
% calculate stiffness matrix
for i=1:nel
nd = element(:,i);
x = gcoord(1,nd); y = gcoord(2,nd);
% compute long of each bar
le = sqrt((x(2)-x(1))^2 + (y(2)-y(1))^2);
% compute direction cosin
l_ij = (x(2)-x(1))/le; % Eq.5.19
m_ij = (y(2)-y(1))/le; % Eq.5.19
% compute transform matrix
Te = [l_ij m_ij 0 0;
0 0 l_ij m_ij];
% compute stiffness matrix of element
ke = A(i)*E/le*[1 -1;
-1 1];
ke = Te'*ke*Te;
% compute force vector
%fe = A(i)*rho*le*9.81/2*[0 1, 0 1];
% find index assemble
index = [2*nd(1)-1 2*nd(1) 2*nd(2)-1 2*nd(2)];
% assemble ke in K
K(index,index) = K(index,index) + ke;
% assemble fe in F
%F(index,1) = F(index,1) + fe';
end F = -F;
%*****************************
% boundary condition *
%*****************************
bcdof = [(5:6)*2-1, (5:6)*2]; % boundary condition displacement
bcval = zeros(1,length(bcdof)); % value of boundary
condition
%*********************
% Node | PX PY |
%---------------------
% 2 | 0 1e5 |
% 4 | 0 1e5 |
%---------------------
% Note: loads are in kips. 1 (kips) = 1000 (lb)
%------------------------
% consider to self weight
% Node 2
F([4 8]) = -P;
%*******************************************************
% apply boundary condition force and displacement *
%*******************************************************
n = length(bcdof);
for i = 1:n
c = bcdof(i);
K(c,:) = 0;
K(c,c) = 1 ;
F(c) = bcval(i);
end
u = K\F;
max(abs(u));
% rang buoc chuyen vi
Cdisp = 0.0508-max(abs(u));
function [ Csig ] = f2(x)
% E = 6.89*1e10; % Young's elastic modulus (N/m^2)
% A = x;
% P = 444.82e3;
E = x(1); P= x(2); A = x(3:end);
%--------------------------------------------------------------------------
% 1 2 3 4 5 6
gcoord = [18.288, 18.288, 9.144, 9.144, 0, 0
9.144, 0, 9.144, 0, 9.144, 0];
% 1 2 3 4 5 6 7 8 9 10
element = [3, 1, 4, 2, 3, 1, 4, 3, 2, 1
5, 3, 6, 4, 4, 2, 5, 6, 3, 4]; nel = length(element); % total element
nnode = length(gcoord); % total node
ndof = 2; % number of degree of freedom of one node
sdof = nnode*ndof; % total dgree of freedom of system
% plotModel( type,gcoord,element );
% calculate stiffness matrix
K = zeros(sdof,sdof); % creat matrix of stiffness initial is zeros
F = zeros(sdof,1); % creat matrix of force is zeros
%---------------------------
% calculate stiffness matrix
for i=1:nel
nd = element(:,i);
x = gcoord(1,nd); y = gcoord(2,nd);
% compute long of each bar
le = sqrt((x(2)-x(1))^2 + (y(2)-y(1))^2);
% compute direction cosin
l_ij = (x(2)-x(1))/le; % Eq.5.19
m_ij = (y(2)-y(1))/le; % Eq.5.19
% compute transform matrix
Te = [l_ij m_ij 0 0;
0 0 l_ij m_ij];
% compute stiffness matrix of element
ke = A(i)*E/le*[1 -1;
-1 1];
ke = Te'*ke*Te;
% compute force vector
%fe = A(i)*rho*le*9.81/2*[0 1, 0 1];
% find index assemble
index = [2*nd(1)-1 2*nd(1) 2*nd(2)-1 2*nd(2)];
% assemble ke in K
K(index,index) = K(index,index) + ke;
% assemble fe in F
%F(index,1) = F(index,1) + fe';
end
F = -F;
%*****************************
% boundary condition *
%*****************************
bcdof = [(5:6)*2-1, (5:6)*2]; % boundary condition displacement
bcval = zeros(1,length(bcdof)); % value of boundary
condition
%*********************
% Node | PX PY | %---------------------
% 2 | 0 1e5 |
% 4 | 0 1e5 |
%---------------------
% Note: loads are in kips. 1 (kips) = 1000 (lb)
%------------------------
% consider to self weight
% Node 2
F([4 8]) = -P;
%*******************************************************
% apply boundary condition force and displacement *
%*******************************************************
n = length(bcdof);
for i = 1:n
c = bcdof(i);
K(c,:) = 0;
K(c,c) = 1 ;
F(c) = bcval(i);
end
u = K\F;
% TINH UNG SUAT
%
for i=1:nel
nd = element(:,i);
x = gcoord(1,nd); y = gcoord(2,nd);
% compute long and direction cosin of bar
le = sqrt((x(2)-x(1))^2 + (y(2)-y(1))^2);
% compute direction cosin
l_ij = (x(2)-x(1))/le; m_ij = (y(2)-y(1))/le;
% compute transform matrix
Te = [l_ij m_ij 0 0;
0 0 l_ij m_ij];
% compute stiffness matrix of element
Be = [-1/le 1/le];
% chi so bac tu do cua phan tu
index = [2*nd(1)-1 2*nd(1) 2*nd(2)-1 2*nd(2)];
disp_e = u(index); % chuyen vi cua phan tu
de_o = Te*disp_e;
stress_e(i) = E*Be*de_o; % tinh ung suat
end
% rang buoc ung suat
LimitSig = 173.375e6;
max(abs(stress_e));
Csig =LimitSig-max(abs(stress_e)); function beta = CheckReliability(xval)
sigVal = 0.05*xval;
E = 6.89*1e10;
sigE = 0.05*E;
addmass = 454;
sigadd = 454*0.05;
rho = 2767;
sigrho = 0.05*rho;
muy_X = [E addmass rho xval]; % mean value of random
variables
sigma = [sigE sigadd sigrho sigVal]; % standard deviation
[ beta1 ] = FORM_Reliability( muy_X,sigma,@f1 );
[ beta2 ] = FORM_Reliability( muy_X,sigma,@f2 );
[ beta3 ] = FORM_Reliability( muy_X,sigma,@f3 );
beta = [beta1;beta2;beta3];
function [ c1 ] = f1( x )
type = '2D';
E = x(1); addedMass = x(2); rho = x(3); A = x(4:end);
%--------------------------------------------------------------------------
% 1 2 3 4 5 6
gcoord = [18.288, 18.288, 9.144, 9.144, 0, 0
9.144, 0, 9.144, 0, 9.144, 0];
% 1 2 3 4 5 6 7 8 9 10
element = [3, 1, 4, 2, 3, 1, 4, 3, 2, 1
5, 3, 6, 4, 4, 2, 5, 6, 3, 4];
nnode = length(gcoord); % total node
ndof = 2; % number of degree of freedom of one node
sdof = nnode*ndof; % total dgree of freedom of system
% plotModel( type,gcoord,element );
% calculate stiffness matrix
[ K,M ] = Cal_K_and_M( type,gcoord,element,A,rho,E );
for idof = 1:sdof
M(idof,idof) = M(idof,idof) + addedMass;
end
% apply boundary
bcdof = [(5:6)*2-1, (5:6)*2]; % boundary condition displacement
% Giai phuong trinh tim tri rieng va vector rieng
[omega_2,~]=eigens(K,M,bcdof);
f=sqrt(omega_2)/2/pi;
c1 =f(1)-7; function [ c2 ] = f2( x )
type = '2D';
E = x(1); addedMass = x(2); rho = x(3); A = x(4:end);
%--------------------------------------------------------------------------
% 1 2 3 4 5 6
gcoord = [18.288, 18.288, 9.144, 9.144, 0, 0
9.144, 0, 9.144, 0, 9.144, 0];
% 1 2 3 4 5 6 7 8 9 10
element = [3, 1, 4, 2, 3, 1, 4, 3, 2, 1
5, 3, 6, 4, 4, 2, 5, 6, 3, 4];
nnode = length(gcoord); % total node
ndof = 2; % number of degree of freedom of one node
sdof = nnode*ndof; % total dgree of freedom of system
% plotModel( type,gcoord,element );
% calculate stiffness matrix
[ K,M ] = Cal_K_and_M( type,gcoord,element,A,rho,E );
for idof = 1:sdof
M(idof,idof) = M(idof,idof) + addedMass;
end
% apply boundary
bcdof = [(5:6)*2-1, (5:6)*2]; % boundary condition displacement
% Giai phuong trinh tim tri rieng va vector rieng
[omega_2,~]=eigens(K,M,bcdof);
f=sqrt(omega_2)/2/pi;
c2 =f(2)-15;
function [ c3 ] = f3( x )
type = '2D';
% E = 6.89*1e10; % Young's elastic modulus (N/m^2)
% A = x;
% rho = 2767; % density of material (kg/m^3)
% addedMass = 454; %kg
E = x(1); addedMass = x(2); rho = x(3); A = x(4:end);
%--------------------------------------------------------------------------
% 1 2 3 4 5 6
gcoord = [18.288, 18.288, 9.144, 9.144, 0, 0
9.144, 0, 9.144, 0, 9.144, 0];
% 1 2 3 4 5 6 7 8 9 10
element = [3, 1, 4, 2, 3, 1, 4, 3, 2, 1
5, 3, 6, 4, 4, 2, 5, 6, 3, 4];
nnode = length(gcoord); % total node
ndof = 2; % number of degree of freedom of one node
sdof = nnode*ndof; % total dgree of freedom of systemVÍ DỤ SỐ
KẾT LUẬN VÀ HƯỚNG PHÁT TRIỂN
MỘT SỐ KẾT QUẢ ĐÃ ĐƯỢC CÔNG BỐ
TÀI LIỆU THAM KHẢO
PHỤ LỤC
n = length(bcdof);
for i = 1:n
c = bcdof(i);
K(c,:) = 0;
K(c,c) = 1 ;
F(c) = bcval(i);
end
u = K\F;
% Tinh chuyen vi
Ux = u(1:2:end)
Uy = u(2:2:end)
% TINH UNG SUAT
%
for i=1:nel
nd = element(:,i);
x = gcoord(1,nd); y = gcoord(2,nd);
% compute long and direction cosin of bar
le = sqrt((x(2)-x(1))^2 + (y(2)-y(1))^2);
% compute direction cosin
l_ij = (x(2)-x(1))/le; m_ij = (y(2)-y(1))/le;
% compute transform matrix
Te = [l_ij m_ij 0 0;
0 0 l_ij m_ij];
% compute stiffness matrix of element
Be = [-1/le 1/le];
% chi so bac tu do cua phan tu
index = [2*nd(1)-1 2*nd(1) 2*nd(2)-1 2*nd(2)];
disp_e = u(index); % chuyen vi cua phan tu
de_o = Te*disp_e;
stress_e(i) = E*Be*de_o; % tinh ung suat
end
max(abs(stress_e))
2. Code Matlab hàm ràng buộc của bài toán tối ưu DO dàn phẳng 10 thanh
function [c,ceq] = hamrangbuoc10bars(x)
%**********************
% Input parameter *
%**********************
type = '2D';
E = 6.89*1e10; % Young's elastic modulus (N/m^2)
A = x;
rho = 2767; % density of material (kg/m^3)
P = 444.82e3;
% Cdisp = max(abs(u))/0.1016-1; % chuyen vi 10.16cm
%
% TINH UNG SUAT
%
for i=1:nel
nd = element(:,i);
x = gcoord(1,nd); y = gcoord(2,nd);
% compute long and direction cosin of bar
le = sqrt((x(2)-x(1))^2 + (y(2)-y(1))^2);
% compute direction cosin
l_ij = (x(2)-x(1))/le; m_ij = (y(2)-y(1))/le;
% compute transform matrix
Te = [l_ij m_ij 0 0;
0 0 l_ij m_ij];
% compute stiffness matrix of element
Be = [-1/le 1/le];
% chi so bac tu do cua phan tu
index = [2*nd(1)-1 2*nd(1) 2*nd(2)-1 2*nd(2)];
disp_e = u(index); % chuyen vi cua phan tu
de_o = Te*disp_e;
stress_e(i) = E*Be*de_o; % tinh ung suat
end
% rang buoc ung suat
LimitSig = 173.375e6;
max(abs(stress_e));
Csig = max(abs(stress_e))/LimitSig-1;
%% rang buoc bai toan
% c = [Cdisp Csig c4]; % Rang buoc bat phuong trinh
c = [c1 c2 c3 Cdisp Csig];
end
3. Code Matlab hàm ràng buộc của bài toán tối RBDO dàn phẳng 10 thanh
function [c,ceq] = hamrangbuoc10bars(xval)
sigVal = 0.05*xval;
E = 6.89e10;
sigE = 0.05*E;
P = 444.82e3; V_P = 0.1;
siglogP = sqrt(log(1 + V_P^2));
muylogP = log(P)-1/2*siglogP^2;
sigP = P*siglogP;
muyP = P*(1 - log(P) + muylogP);
addedMass = 454;
sigadd = addedMass*0.05;
Cdisp = max(abs(u))/0.1016-1; %Chuyen vi 10.16 cm
% TINH UNG SUAT
%
for i=1:nel
nd = element(:,i);
x = gcoord(1,nd); y = gcoord(2,nd);
% compute long and direction cosin of bar
le = sqrt((x(2)-x(1))^2 + (y(2)-y(1))^2);
% compute direction cosin
l_ij = (x(2)-x(1))/le; m_ij = (y(2)-y(1))/le;
% compute transform matrix
Te = [l_ij m_ij 0 0;
0 0 l_ij m_ij];
% compute stiffness matrix of element
Be = [-1/le 1/le];
% chi so bac tu do cua phan tu
index = [2*nd(1)-1 2*nd(1) 2*nd(2)-1 2*nd(2)];
disp_e = u(index); % chuyen vi cua phan tu
de_o = Te*disp_e;
stress_e(i) = E*Be*de_o; % tinh ung suat
end
% rang buoc ung suat
LimitSig = 173.375e6;
max(abs(stress_e));
Csig = max(abs(stress_e))/LimitSig-1;
end
4. Code Matlab hàm mục tiêu chuyển vị của bài toán dàn phẳng 10 thanh
function CV = hammuctieu_chuyenvi(x)
%**********************
% Input parameter *
%**********************
E = 6.89e10;
A = x;
P = 444.8e3; % tai tac dung tai nut
%--------------------------------------------------------------------------
% 1 2 3 4 5 6
gcoord = [18.288, 18.288, 9.144, 9.144, 0, 0
9.144, 0, 9.144, 0, 9.144, 0];
% 1 2 3 4 5 6 7 8 9 10
element = [3, 1, 4, 2, 3, 1, 4, 3, 2, 1
5, 3, 6, 4, 4, 2, 5, 6, 3, 4];
%------------------------------------------------
nel = length(element); % total element
%---------------------
% 2 | 0 1e5 |
% 4 | 0 1e5 |
%---------------------
% Note: loads are in kips. 1 (kips) = 1000 (lb)
%------------------------
% consider to self weight
% Node 2
F([4 8]) = -P;
%*******************************************************
% apply boundary condition force and displacement *
%*******************************************************
n = length(bcdof);
for i = 1:n
c = bcdof(i);
K(c,:) = 0;
K(c,c) = 1 ;
F(c) = bcval(i);
end
u = K\F;
CV=max(abs(u));
5. Code Matlab hàm mục tiêu khối lượng của bài toán dàn phẳng 10 thanh
function Weight= hammuctieu_khoiluong(section)
A = section; % area of bar (m^2)
rho = 2767; % density of material (kg/m^3)
%--------------------------------------------------------------------------
% 1 2 3 4 5 6
gcoord = [18.288, 18.288, 9.144, 9.144, 0, 0
9.144, 0, 9.144, 0, 9.144, 0];
% 1 2 3 4 5 6 7 8 9 10
element = [3, 1, 4, 2, 3, 1, 4, 3, 2, 1
5, 3, 6, 4, 4, 2, 5, 6, 3, 4];
%--------------------------------------------------------------------------
% calculate Weight matrix
Weight = 0;
for i=1:length(element)
nd = element(:,i);
x = gcoord(1,nd); y = gcoord(2,nd);
% compute long of each bar
le = sqrt((x(2)-x(1))^2 + (y(2)-y(1))^2);
Weight = Weight + rho*le*A(i);
end
% plotModel( type,gcoord,element );
% calculate stiffness matrix
[ K,M ] = Cal_K_and_M( type,gcoord,element,A,rho,E );
% add non-structural mass
% addedMass = 454; %kg
for idof = 1:sdof
M(idof,idof) = M(idof,idof) + addedMass;
end
% apply boundary
bcdof = [(5:6)*2-1, (5:6)*2]; % boundary condition displacement
% Giai phuong trinh tim tri rieng va vector rieng
[omega_2,~]=eigens(K,M,bcdof);
f=sqrt(omega_2)/2/pi;
% f(1:5)
c3 = f(3)-20;