Xác định các thông số của xe chạy trên cầu từ tập số liệu đo dao động

là thời điểm khối lượng di động mi bắt đầu đi vào dầm. Gọi ui(t) là dịch chuyển theo phương thẳng đứng của khối lượng mi so với vị trí cân bằng tĩnh của nó trong hệ toạ độ tổng thể. Phương trình vi phân dao động của khối lượng di động mi có kể đến chuyển động của phần tử dầm tại vị trí ξi có dạng: (97) Trong đó: g – gia tốc trọng trường, mig – trọng lượng của khối mi, Qi – tải trọng động tác dụng lên khối lượng di động, δi – biến dạng tĩnh của lò xo do trọng lượng mig, tức là: (98) Phương trình khi không kể đến biến dạng tĩnh của lò xo có dạng: (99) Tải trọng từ khối lượng di động mi tác dụng lên dầm tại toạ độ ξi gồm lực đàn hồi và lực cản của hệ liên kết lò xo – cản: (100) Từ ta có: (101) Ri(t) là tải trọng tập trung tác dụng vuông góc với dầm. Theo , tải trọng Ri(t) quy về thành tải trọng nút của phần tử theo công thức: (102) Từ , phương trình dao động của phần tử dầm chịu tác dụng của lực các lực Pi(t) có dạng: (103) Trong đó: k – ma trận độ cứng của phần tử dầm, m – ma trận khối lượng của phần tử dầm, c – ma trận cản nhớt của phần tử dầm, – véc tơ gia tốc chuyển động tại nút của phần tử dầm, – véc tơ vận tốc chuyển động tại nút của phần tử dầm q – véc tơ chuyển vị nút của phần tử dầm, F – véc tơ tổng tải trọng quy nút của các thành phần tải trọng tác dụng lên phần tử dầm: (104) p – véc tơ tải trọng quy nút của các tải trọng tác dụng trên phần tử dầm. Lưu ý tới , chuyển vị và vận tốc theo phương thẳng đứng của dầm tại toạ độ ξ được nội suy từ các véc tơ chuyển vị và vận tốc tại các nút dầm: (105) Trong đó, Nv – ma trận các hàm nội suy độ võng của phần tử dầm tại tọa độ ξi từ các chuyển vị tại nút của phần tử: (106) Thay vào nhận được: (107) Ta thấy rằng, phương trình vi phân dao động của phần tử dầm cho bởi với véctơ tải trọng quy nút chứa thành phần gia tốc chuyển động đứng của các khối lượng di động, mặt khác, phương trình vi phân dao động của khối lượng di động lại chứa các véc tơ chuyển vị q và véc tơ vận tốc của dầm, điều đó thể hiện tính chất tương tác của hệ dầm – khối lượng di động. Từ và , mở rộng các ma trận phần tử hữu hạn cho dầm chịu uốn với số bậc tự do tăng lên bằng số khối lượng di động trên dầm, hệ phương trình dao động mở rộng cho dầm và các khối lượng di động trên dầm có dạng: (108) Trong đó: – véc tơ chuyển vị mở rộng của dầm có kể đến các bậc tự do của các khối lượng di động trên dầm: (109) – véc tơ vận tốc mở rộng của dầm có kể đến các bậc tự do của các khối lượng di động trên dầm: (110) – véc tơ gia tốc mở rộng của dầm có kể đến các bậc tự do của các khối lượng di động trên dầm: (111) – ma trận khối lượng mở rộng của dầm có kể đến các bậc tự do của các khối lượng di động trên dầm: (112) – ma trận cản nhớt mở rộng của dầm có kể đến các bậc tự do của các khối lượng di động trên dầm: (113) – ma trận độ cứng mở rộng của dầm có kể đến các bậc tự do của các khối lượng di động trên dầm: (114) px – véc tơ tải trọng quy nút mở rộng của dầm có kể đến các bậc tự do của các khối lượng di động trên dầm: (115) Biểu thức là hệ phương trình vi phân tuyến tính có các hệ số phụ thuộc thời gian vì chứa đại lượng ξ(t). Từ công thức của các ma trận phần tử hữu hạn trên, có thể nhận thấy rằng các ma trận này có dạng không đối xứng. Điều này sẽ gây khó khăn đến quá trình giải bài toán dao động sau này. Các ma trận phần tử hữu hạn trên được xây dựng trong hệ tọa độ cục bộ của phần tử. Để tập hợp được các ma trận tổng thể của toàn kết cấu, các ma trận này phải được chuyển về hệ tọa độ tổng thể theo các công thức đã dẫn trong chương trước. 1.3.2. Phương trình dao động của hệ cầu dầm có xe chạy trên Phương trình dao động của toàn hệ kết cấu cầu dầm có xe chạy trên trong hệ tọa độ tổng thể khi chưa xét đến các điều kiện biên có dạng: (116) Trong đó: – ma trận khối lượng toàn hệ, tập hợp từ các ma trận khối lượng mx của phần tử, – ma trận cản nhớt toàn hệ, tập hợp từ các ma trận cản cx của phần tử, – ma trận độ cứng toàn hệ, tập hợp từ các ma trận độ cứng kx của phần tử, – véc tơ lực quy nút tác dụng lên cầu do tải trọng xe chạy và các tải trọng khác đặt trên cầu tập hợp từ các véc tơ tải của phần tử và các tải trọng đặt trực tiếp trên các nút của hệ, – véc tơ gia tốc nút của toàn kết cấu có kể đến bậc tự do của các khối lượng di động trên cầu. – véc tơ vận tốc nút của toàn kết cấu có kể đến bậc tự do của các khối lượng di động trên cầu. – véc tơ chuyển vị nút của toàn kết cấu có kể đến bậc tự do của các khối lượng di động trên cầu. Các ma trận này được viết trong hệ toạ độ tổng thể của kết cấu. Trong hệ phương trình trên, số ẩn số của phương trình bao gồm số bậc tự do của kết cấu và tổng số bậc tự do của các xe (khối lượng di động) trên cầu. Cụ thể, đối với bài toán xe chạy trên cầu dầm liên tục, số ẩn số của phương trình được tính theo công thức: (117) Phương trình chuyển động của hệ hỗn hợp cầu và xe chạy trên cầu sau khi đã đưa các điều kiện biên vào có dạng tổng quát: (118) 1.3.3. Phương pháp Newmark giải hệ phương trình chuyển động Có nhiều phương pháp để giải hệ phương trình vi phân chuyển động của hệ kết cấu rất hiệu quả như phương pháp chồng mode, phương pháp tích phân trực tiếp, ... Như đã trình bày trong phần trước, hệ phương trình có các ma trận hệ số K, M, C không đối xứng. Do vậy, không thể sử dụng được phương pháp chồng mode. Trong luận án chọn phương pháp tích phân trực tiếp để giải hệ phương trình , cụ thể là chọn phương pháp tích phân Newmark. Theo phương pháp Newmark nghiệm của phương trình được tìm như sau. Giả sử đã biết chuyển động của hệ tại thời điểm t là , ,. Phương trình chuyển động viết tại thời điểm có dạng: (119) Khai triển theo chuỗi Taylor của và ta có: (120) Loại bỏ các thành phần bậc cao trong các phương trình và viết lại theo dạng sau: (121) Trong đó α và δ là các tham số được xác định trước để đảm bảo độ chính xác và ổn định của lời giải. Theo Newmark gia tốc của hệ được giả thiết là tuyến tính trong phạm vi mỗi bước thời gian, theo đó: (122) Thay biểu thức vào sẽ nhận được: (123) Các tham số α và δ được chọn như sau: (124) Để đảm bảo sự ổn định và độ chính xác của lời giải theo phương pháp tích phân trực tiếp Newmark, chọn α=0,25 và δ=0,5. Đưa vào một số ký hiệu mới: (125) Từ biểu thức ta có thể biểu diễn gia tốc và vận tốc tại thời điểm qua chuyển vị tại thời điểm và các đại lượng đã biết tại thời điểm t như sau: (126) Thay các giá trị và từ biểu thức vào phương trình ta có hệ phương trình chuyển động tại thời điểm chỉ chứa các ẩn số là chuyển vị tại thời điểm đó: (127) hoặc viết gọn lại dưới dạng: (128) trong đó: – ma trận độ cứng hiệu quả: (129) – véc tơ tải trọng hiệu quả: (130) Giải ta tìm được . Thay giá trị này vào sẽ tìm được và . 1.3.4. Thuật toán giải bài toán dao động của cầu dầm liên tục chịu tác dụng của đoàn tải trọng di động Từ các công thức tính toán đã trình bày ở trên, trong mục này xây dựng thuật toán giải bài toán dao động của cầu dầm liên tục chịu tác dụng của đoàn tải trọng di động. Các bước tính toán được thực hiện theo trình tự sau: Bước 1: Nhập số liệu Nhập các thông số chung cho chương trình Các hằng số vật lý Các tham số cản tỷ lệ của kết cấu α và β Thời gian tính toán: Bước thời gian Δt Tổng thời gian tính toán Nhập các số liệu của cầu bao gồm: Số liệu nút: Số lượng nút Tọa độ nút trong hệ tọa độ tổng thể Số liệu phần tử Số lượng phần tử Số liệu cho mỗi phần tử: Liên kết phần tử (nút đầu, nút cuối) Đặc trưng vật liệu của phần tử: E, ρ Đặc trưng tiết diện của phần tử: J, F Số liệu liên kết: Số lượng liên kết Bậc tự do của các liên kết Độ cứng của các liên kết Nhập các số liệu của xe bao gồm: Số lượng xe Số liệu cho mỗi xe: Khối lượng xe Độ cứng liên kết lò xo của xe Hệ số cản nhớt liên kết cản của xe Vận tốc xe chạy Bước 2: Xây dựng các ma trận phần tử hữu hạn cho các phần tử dầm Xây dựng các ma trận m, c, k của các phần tử Xây dựng các véc tơ tải trọng quy nút do các lực tập trung và phân bố tác dụng lên phần tử Đưa các ma trận và véc tơ của phần tử về hệ tọa độ tổng thể Bước 3: Tập hợp các ma trận phần tử hữu hạn cho toàn kết cấu Tập hợp ma trận M, C, K của kết cấu dầm liên tục trên cơ sở các ma trận m, c, k của phẩn tử Xây dựng véc tơ tải trọng quy nút của kết cấu trên cơ sở các véc tơ tải trọng quy nút của phần tử và các lực tập trung tác dụng trực tiếp lên nút của kết cấu. Bước 4: Đưa các điều kiện biên vào phương trình Xử lý các điều kiện biên: Gối cứng: loại bỏ các phương trình tương ứng với bậc tự do bị cản trở Gối đàn hồi: Thêm độ cứng của các bậc tự do tương ứng với liên kết vào ma trận độ cứng kết cấu K Thiết lập phương trình dao động của kết cấu với các ma trận M, C, K đã khử điều kiện biên Bước 5: Giải bài toán dao động của cầu dầm liên tục khi có đoàn tải trọng di động trên cầu Xuất phát với t=0 Xác định ảnh hưởng của từng xe lên kết cấu: Xác định vị trí xe trên cầu, tìm phần tử dầm mà xe đang chạy qua Mở rộng số bậc tự do cho các ma trận m, c, k của phần tử dầm Kể thêm ảnh hưởng của các xe chạy trên phần tử vào các ma trận m, c, k của phần tử dầm Tập hợp lại các ma trận m, c, k mới của phần tử dầm vào các ma trận tổng thể tương ứng Xây dựng véc tơ tải trọng quy nút của phần tử dầm do tải trọng của xe truyền xuống dầm thông qua liên kết lò xo và cản Tập hợp véc tơ tải trọng quy nút mới của phần tử dầm vào véc tơ tải trọng tổng thể Lặp lại các bước tính cho tất cả các xe Tính các tham số chuyển động của kết cấu và các xe theo công thức tích phân Newmark với các ma trận M, C, K tổng thể của kết cấu vừa cập nhật. Tính nội lực phần tử sau khi biết các tham số chuyển động của kết cấu Lặp lại ở bước thời gian tiếp theo cho đến khi kết thúc Sơ đồ khối của thuật toán được thể hiện trong hình 14. 1.3.5. Chương trình máy tính Trên cơ sở các thuật toán trình bày ở trên đã xây dựng chương trình giải bài toán dao động của cầu chịu tải trọng của đoàn tải trọng di động bằng ngôn ngữ MATLAB. Chương trình được đặt tên là VBIS (Vehicle-Bridge Interactive Simulation). Chương trình này sẽ được sử dụng để mô phỏng các số liệu đo dao động của kết cấu cầu khi có xe chạy trên cầu. Từ kết quả đó, kết hợp với số liệu đo đạc và thuật toán nhận dạng trình bày trong chương sau để xác định các tham số của xe chạy trên cầu. Hình 14 Sơ đồ khối chương trình 1.3.6. Các ví dụ tính toán Ví dụ 1. Trong bài báo của Gou và Xu (2000) [8] đã khảo sát một mô hình khối lượng di động trên dầm đơn giản. Liên kết giữa khối lượng di động và dầm là liên kết đàn hồi (hình 15). Để kiểm tra chương trình đã lập, khảo sát lại mô hình của Gou và Xu. Hình 15 Mô hình khảo sát của Gou & Xu Các số liệu xuất phát Số liệu dầm: Chiều dài L = 25 m; Mô đun đàn hồi E = 2,87 Gpa = 2,87.109 N/m2; Mô men quán tính tiết diện J = 2,9 m4; Khối lượng trên một đơn vị dài của dầm ρ = 2303 kg/m; Số liệu xe: Khối lượng xe mv = 5750 kg; Độ cứng lò xo kv = 1595 kN/m = 1,595.106 N/m; Vận tốc di chuyển v=100 km/h=27,778 (m/s). Kết quả tính toán Kết quả tính toán của chương trình được thể hiện trong hình 16 và hình 17. So sánh 2 kết quả chúng ta thấy sai khác là rất nhỏ. Trong kết quả khảo sát của Gou & Xu (scan từ bài báo), các tác giả cũng có so sánh với kết quả khảo sát trước đó của Yang & Yau (1997) và kết quả tính toán bằng giải tích (Analytic) Kết quả tính độ võng tại điểm giữa cầu bằng chương trình MLB Kết quả tính độ võng tại điểm giữa theo Gou & Xu Hình 16 So sánh độ võng tại điểm giữa dầm với kết quả của Gou & Xu Kết quả tính dịch chuyển thẳng đứng của khối lượng Mv bằng chương trình MLB Kết quả tính dịch chuyển thẳng đứng của khối lượng Mv theo Gou & Xu Hình 17 So sánh chuyển vị của khối lượng với kết quả của Gou & Xu Ví dụ 2. Khảo sát kết cấu cầu có mô hình là dầm đơn giản liên kết khớp 2 đầu, tải trọng có liên kết đần hồi và cản nhớt với dầm (hình 18). Đối với dầm đơn giản, thiết diện không đổi như trong ví dụ này đã có lời giải chính xác [5]. Hình 18 Xe chạy trên dầm đơn giản Các số liệu xuất phát Số liệu dầm: Chiều dài dầm L = 20 m. Độ cứng chống uốn EJ=4157.106 Nm2, Khối lượng trên đơn vị chiều dài dầm r=2277kg/m. Số liệu xe: Khối lượng xe m=4333kg, Độ cứng lò xo k=902000N/m, Hệ số cản c=11016Ns/m, Vận tốc di chuyển v=20m/s. Kết quả tính toán Chia kết cấu thành 20 phần tử bởi 21 nút. Chia khoảng thời gian xe chay qua cầu là 100 bước tính. Mỗi bước Δt=0,01s. Biểu đồ độ võng tại giữa nhịp (đường cong u) và dịch chuyển thẳng đứng của khối lượng di động (đường cong w) phụ thuộc vào toạ độ di chuyển của khối lượng m từ 0 đến 20m được cho trên hình 19. Sai số của kết quả tính toán được trình bày trong ví dụ này so với nghiệm chính xác là không đáng kể, biên độ dao động sai khác dưới 4%. Hình 19: Độ võng giữa nhịp và chuyển vị của tải trọng 1.4. Kết luận chương Trong chương này đã đạt được các kết quả chính sau: Đã trình bày cơ sở của phương pháp phần tử hữu hạn cho bài toán dao động của hệ dầm liên tục. Đã xây dựng các ma trận phần tử hữu hạn cho phần tử dầm chịu uốn thuần túy và các ma trận phần tử hữu hạn cho kết cấu dầm liên tục trên cơ sở tập hợp các ma trận tổng thể từ các ma trận phần tử. Đã xây dựng được phương trình vi phân dao động của hệ dầm liên tục bằng phương pháp phần tử hữu hạn. Đã xây dựng được phương trình vi phân dao động của hệ tương tác: phần tử dầm chịu uốn – tải trọng di động, trong đó xét đến khối lượng phân bố của dầm, khối lượng của tải trọng, các liên kết đàn hồi và cản nhớt giữa dầm và tải trọng di động. Từ đó, thiết lập được các ma trận khối lượng, ma trận cản, ma trận độ cứng mở rộng cho dầm có kể đến các khối lượng di động trên dầm và véctơ tải trọng nút. Các ma trận này chứa các thông số liên quan đến bản thân dầm và các tải trọng di động trên dầm (khối lượng m, hệ số cứng k của lò xo, hệ số cản c, vận tốc v), chúng đều phụ thuộc thời gian t và không đối xứng. Đã xây dựng phương trình dao động của kết cấu dầm liên tục chịu tác dụng của tải trọng di động. Thiết lập các ma trận tổng thể cho toàn hệ trên cơ sở tập hợp từ các ma trận của phần tử. Đã triển khai thuật toán và lập chương trình bằng ngôn ngữ MATLAB để tính dao động kết cấu dầm liên tục chịu tải trọng động của đoàn xe di động bằng phương pháp phần tử hữu hạn. Các ví dụ tính toán thể hiện độ tin cậy của thuật toán và chương trình đã lập XÁC ĐỊNH CÁC THÔNG SỐ CỦA XE CHẠY TRÊN CẦU TỪ TỆP SỐ LIỆU ĐO DAO ĐỘNG Thông số của xe chạy trên cầu được xác định bằng phương pháp nhận dạng trạng thái làm việc của kết cấu từ các tập số liệu đo đạc, giám sát dao động của cầu kết hợp với giải bài toán mô phỏng tương tác động lực học của hệ cầu – xe trên máy tính. 2.1. Số liệu đo đạc, giám sát dao động của cầu khi có xe chạy qua. Hệ thống giám sát tình trạng làm việc của cầu là hệ thống đo đạc hiện đại với các cảm biến được gắn trên toàn bộ cầu, giám sát liên tục nhiều thông số trạng thái kỹ thuật của cầu như ứng suất - biến dạng trong kết cấu cầu, độ võng, gia tốc dao động theo các phương của kết cấu nhịp, v.v.... Trong các số liệu đo đạc được của hệ thống giám sát, số liệu đo gia tốc dao động theo phương thẳng đứng của kết cấu nhịp là số liệu phản ánh rõ ràng nhất ảnh hưởng của các phương tiện chạy trên cầu. Các kết quả đo đạc giám sát được sử dụng để xác định tham số của xe chạy trên cầu là bản ghi số liệu gia tốc theo phương thẳng đứng tại các trạm đo trong suốt thời gian phương tiện di chuyển qua cầu. 2.2. Xác định thông số của xe chạy trên cầu Bài toán xác định các thông số của xe từ các số liệu đo dao động thực chất là bài toán nhận dạng các thông số của xe di chuyển trên cầu từ các số liệu đo đạc giám sát. Với kết cấu được giám sát, các thông số về tình trạng kết cấu, phản ứng của kết cấu khi có xe chạy qua luôn được xác định, chỉ còn các thông số của xe là chưa biết. Ứng dụng lý thuyết nhận dạng hệ thống trong trường hợp này để xác định các thông số của mô hình xe chạy qua. Từ đó, xác định được trọng lượng xe. Các thông số của xe cần xác định (nhận dạng) là khối lượng m, hệ số cản c độ cứng lò xo k của xe và vận tốc xe chạy v. Trong luận văn này, giới hạn với giả thiết xe qua cầu với vận tốc v là không đổi. Như đã nói ở trên, gia tốc dao động theo phương thẳng đứng là đại lượng phản ánh rõ ràng nhất ảnh hưởng của xe tới kết cấu cầu khi xe di chuyển qua cầu. Vì vậy, các thông số của xe chạy trên cầu được xem là nhận dạng đúng khi kết quả đo gia tốc theo phương thẳng đứng tại các trạm đo và kết quả tính toán mô phỏng bằng máy tính với các thông số nhận dạng được của xe là bằng nhau hặc có sự sai lệch đủ nhỏ. Như vậy, khi xác định các thông số của xe chạy trên cầu, sự phù hợp giữa kết quả đo đạc và mô phỏng được đánh giá bằng sai lệch giữa bản ghi gia tốc đo được tại vị trí quan sát và các giá trị gia tốc tính toán mô phỏng từ các tham số của xe cần xác định. Khi sự sai lệch này là tối thiểu, các tham số của xe tính toán mô phỏng sẽ là các tham số của xe thực. Gọi ei là véc tơ sai lệch của các giá trị gia tốc đo được và gia tốc tính toán mô phỏng theo thời gian: (131) trong đó: aij,mp – giá trị tại thời điểm thứ j của vector gia tốc tại trạm đo thứ i tính toán mô phỏng trên máy tính aij,do – giá trị tại thời điểm thứ j của vector gia tốc tại trạm đo thứ i lấy từ bản ghi kết quả giám sát dao động của cầu. Với một cây cầu có Q điểm đo và tại mỗi điểm đo bản ghi có R giá trị của gia tốc theo thời gian thì hàm mục tiêu có thể được viết dưới dạng chuẩn của vector sai lệch như sau: (132) Các tham số của xe nhận được khi giải bài toán tối ưu với hàm mục tiêu đạt giá trị cực tiểu: (133) Trong hàm mục tiêu , các giá trị của aij,do nhận được từ hệ thống giám sát, aij,mp nhận được từ phương trình tương tác của hệ cầu xe và phụ thuộc và các thông số của xe di chuyển trên cầu (m,c,k,v). Do vậy, các ràng buộc của bài toán tối ưu là mối quan hệ tương tác của xe di chuyển trên cầu và dao động của cầu : (134) Trong đó, các ma trận M, C, K và véc tơ tải trọng P của kết cấu thay đổi phụ thuộc vào các thông số và vị trí của xe chạy trên cầu. Ràng buộc phản ánh mối liên hệ giữa các thông số của xe cần nhận dạng và giá trị của gia tốc theo thời gian tại các trạm đo. Bài toán tối ưu với các ràng buộc đẳng thức dạng (135) có nhiều phương pháp giải tuy nhiên 2.3. Phương pháp giải bài toán quy hoạch phi tuyến bằng thuật giải di truyền Những khái niệm chung Thuật toán di truyền là một thuật toán được sử dụng trong lập trình tiến hoá. Thuật ngữ chương trình tiến hoá là khái niệm dùng để chỉ các chương trình máy tính có sử dụng các thuật toán tìm kiếm và tối ưu dựa trên quá trình tiến hoá tự nhiên. Các chương trình tiến hoá chia làm ba dạng: Quy hoạch tiến hoá – EP (Evolutionary Programming – D. B. Pogel), EP cho một lớp các phương án khả dĩ giải quyết được từng phần của vấn đề. Dựa vào quy luật tiến hoá, tìm một phương án liên hợp đủ khả năng giải quyết trọn vẹn vấn đề đó. Chiến lược tiến hoá – ES (Evolutionary Strategy – Baeck; F. H. Hofmeister; H. P. Schwefel), Trên cơ sở một số chiến lược ban đầu, áp dụng luật tiến hoá để tạo ra những chiến lược mới phù hợp tốt nhất với môi trường thực. Thuật toán di truyền – GA (Genetic Algorithm – D.E. Goldberg và sau này được phát triển bởi L. Davis và Z. Michalevicz) Thuật toán di truyền cũng như các thuật toán tiến hoá nói chung được hình thành dựa trên quan niệm: Quá trình tiến hoá tự nhiên là quá trình hoàn hảo nhất, hợp lý nhất và tự nó đã mang tính tối ưu. Tính tối ưu của quá trình tiến hoá thể hiện ở chỗ thế hệ sau bao giờ cũng tốt hơn (phát triển hơn, hoàn thiện hơn và phù hợp với môi trường hơn) thế hệ trước. Xuyên suốt quá trình tiến hoá, các thế hệ mới được sinh ra để bổ xung, thay thế thế hệ cũ, trong quá trình này cá thể nào phát triển hơn 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. Với quan điểm của quá trình tiến hoá tự nhiên thì việc giải bài toán tối ưu thực chất chỉ là: “Tìm một cá thể tốt nhất trong các thế hệ đã xem xét thoả mãn mọi điều kiện của môi trường”. Trong di tryền học mỗi sinh vật (cá thể) được đặc trưng bằng chuỗi nhiễm sắc thể nào đó. Trải qua quá trình chọn lọc tự nhiên có những cá thể tồn tại, phát sinh và phát triển, bên cạnh đó có những cá thể bị đào thải. Quá trình lai ghép của thế hệ trước (cha-mẹ) tạo ra thế hệ sau (con) có thể phát huy những tính chất tốt phù hợp hơn với môi trường – cũng có thể không. Với xác xuất nhỏ hơn, một vài cá thể nhận sự đột biến tạo ra một cá thể hoàn toàn mới. Thuật toán di truyền mô phỏng lại toàn bộ quá trình trên và là quá trình lặp của ba toán tử cơ bản như trình bầy trong hình dưới. Bắt đầu t = 0; Khởi tạo quần thể ban đầu P(t) Tính độ thích nghi cho mỗi cá thể thuộc quần thể P(t) Khi (điều kiện dừng chưa đạt) lặp: Tái sinh (Reproduction) P’(t) từ P(t) Lai ghép (Crossover) Q(t) từ P(t-1); Đột biến (Mutation) R(t) từ P(t-1); Chọn lọc (Selection) P(t) từ P(t-1) È Q(t) È R(t) È P(t) Hết lặp Kết thúc Hình 21 Thuật toán di truyền đơn giản Để kết nối bài toán cần giải với các toán tử trên cần mã hoá các biến của bài toán. Đơn giản nhất là mã hoá bằng các chuỗi nhị phân. Như vậy mỗi chuỗi nhị phân với chiều dài xác định sẽ tương ứng với một lời giải của bài toán đang khảo sát. Độ thích nghi của mỗi cá thể trong quần thể được đánh giá thông qua hàm lượng giá là giá trị hàm mục tiêu ứng với các biến tìm kiếm mà chuỗi gen của cá thể đó biểu diễn. Nếu ký hiệu hàm mục tiêu là f(x) thì: (136) Trong quần thể có n cá thể thì khả năng tồn tại của mỗi cá thể sang thế hệ sau đựơc xác định bằng tỷ số: ; (137) Với bài toán cực đại cá thể có độ thích nghi cao ứng với pi lớn và ngược lại. Nền tảng lý thuyết của thuật giải di truyền dựa trên biểu diễn chuỗi nhị phân và lý thuyết về sơ đồ. Một sơ đồ được định nghĩa là một chuỗi có chiều dài bằng chuỗi nhiễm sắc thể, các thành phần của nó có thể nhận một trong các giá trị trong tập ký tự biểu diễn gen hoặc một ký tự đại diện (*). Mỗi sơ đồ biểu diễn một không gian con trong không gian tìm kiếm. Không gian con này là tập tất cả các chuỗi trong không gian lời giải mà với mọi vị trí trong chuỗi, giá trị của gen trùng với giá trị của sơ đồ. Ký tự đại diện ‘*’ có thể trùng khớp với bất kỳ ký tự biểu diễn gen nào. Ví dụ chuỗi có chiều dài là 10 thì sơ đồ: (*111100100) sẽ khớp với hai chuỗi (0111100100) và (1111100100). Mỗi sơ đồ có hai thuộc tính quan trọng là: bậc và chiều dài xác định của sơ đồ. Bậc của sơ đồ S - ký hiệu là O(S) – là số các vị trí 0 và 1, nói cách khác là bậc chính là chiều dài chuỗi trừ đi số ký tự đại diện. Bậc của sơ đồ giúp tính xác xuất sống còn của nó do ảnh hưởng của quá trình đột biến. Chiều dài xác định của sơ đồ S – ký hiệu d(S) – là khoảng cách giữa hai vị trí cố định ở đầu và cuối của chuỗi. Nó định nghĩa độ nén thông tin của chuỗi chứa trong sơ đồ. Khái niệm về chiều dài xác định của sơ đồ giúp tính xác xuất sống còn của sơ đồ do ảnh hưởng của phép lai. Về mặt lý thuyết, người ta đã chứng minh được rằng: Các sơ đồ ngắn, bậc thấp, trên trung bình nhận được số chuỗi tăng theo luỹ thừa trong các thế hệ tiếp theo của thuật toán di truyền. Trong kỹ thuật thực hiện thuật toán di truyền sau khi xác định được chiều dài chuỗi nhiễm sắc thể, các quá trình còn lại: tạo quần thể ban đầu, lai ghép và đột biến là quá trình hoàn toàn ngẫu nhiên (thường chọn theo luật phân bố đều). Như Z. Michalewicz nhận xét: “ Thuật toán di truyền thuộc lớp các thuật toán xác xuất nhưng lại rất khác các thuật toán ngẫu nhiên vì chúng kết hợp được các phần tử tìm kiếm trực tiếp và ngẫu nhiên. Điều khác biệt quan trọng của thuật toán di truyền với các thuật toán khác là thuật toán di truyền luôn duy trì và xử lý một tập các lời giải (quần thể) trong khi các phương pháp khác chỉ xử lý một điểm trong không gian tìm kiếm. Chính vì vậy mà thuật toán di truyền mạnh hơn và hiệu quả hơn các phương pháp khác”. Các thuật toán cơ bản của thuật toán di truyền theo mã nhị phân Mã hóa các biến và xây dựng quần thể ban đầu a) Mã hóa biến thành chuỗi nhị phân Một biến thực với bL – giới hạn trên và bU- là giới hạn dưới và ta muốn xi được biểu diễn với độ chính xác là nx chữ số sau dấu thập phân. Rõ ràng để đạt được điều đó thì miền D sẽ được chia thành miền con. Gọi mi là số nguyên nhỏ nhất sao cho: (138) Như vậy mỗi biến thực xi có thể biểu diễn được bằng một chuỗi nhị phân có độ chính xác bất kỳ. Đoạn chương trình dưới cho phương pháp xác định chiều dài (số bit) cần thiết cho biến xi. b) Quần thể ban đầu Trong thuật toán di truyền, cùng một lúc xét nhiều phương án, tổng số phương án kí hiệu là npop_size. Nếu bài toán có n biến độc lập và mỗi biến được biểu diễn bởi mi(j) bít thì mỗi phương án cần bít để biểu diễn tương ứng toàn bộ bài toán cần (npop_size×nx) bít. Để tiện xử lý trên máy quần thể được lưu trữ dưới dạng mảng. Chọn lọc Như trên đã tình bầy, quá trình chọn lọc các cá thể trong GA cần thỏa mãn xác suất dựa trên độ thích nghi của chúng. Để tính xác suất lựa chọn thực hiện các bước sau: Tính độ thích nghi eval(vi) cho mỗi cá thể; vi (i = 1, 2, , npop_size) Tính tổng giá trị thích nghi cho toàn bộ quần thể (139) Tính xác suất lựa chọn pi cho mỗi cá thể vi: Tính vị trí xác suất qi cho mỗi cá thể; vi (i = 1, 2, , npop_size) (140) Tiến trình chọn lọc được thực hiện bằng cách quay bánh xe Rulet npop_size lần; mỗi một lần sẽ chọn được một cá thể đễ đưa vào quần thể mới. Các bước thực hiện: Phát một số ngẫu nhiên phân bố đều trong khoản [0,1] ; Nếu (r < qi) thì cá thể đầu tiên được chọn (v1) , ngược lại sẽ chọn cá thể thứ i (vi) thỏa mãn điều kiện: (qi-1 < r £ qi). Như vậy sẽ có một số cá thể được chọn hơn một lần và có cá thể đã bị loại bỏ. Lai ghép Thông số có ý nghĩa đối với việc lai ghép là xác suất lai ghép pc, tham số này cho biết số cá thể sẽ tham gia lai ghép trong quần thể là (pc×npop_size). Thông thường pc=0,25. Chọn ngẫu nhiên (pc×npop_size) cá thể trong quần thể, tiếp tục lại chọn ngẫu một cặp cá thể trong số các cá thể sẽ lai ghép. Chọn ngẫu nhiên điểm lai ghép và tiến hành lai ghép theo sơ đồ như hình vẽ Hình 22 Sơ đồ lai ghép trong thuậtt toán di truyền Đột biến Một thông số khác điều khiển quá trình di truyền là xác xuất đột biến pm, tham số này cho biết số bit sẽ đột biến là (pm×nx×npop_size). Cách thực hiện như sau: Phát (pm×nx×npop_size) lần số ngẫu nhiên r phân bố đều trong khoảng [1, (nx×npop_size)]. Nếu r trùng với vị trí nào sẽ tiến tiến hành đột biến bít đó, có nghĩa là nếu giá trị của bít đó bằng 0 thì chuyển thành 1 hoặc ngược lại. Ứng dụng phần mềm MATLAB giải bài toán tối ưu phi tuyến Thư viện tối ưu hoá của MATLAB bao gồm các hàm và các công cụ mạnh để giải quyết các bài toán tối ưu theo các thuật toán hiện đại. Các thuật toán này có thể giải quyết nhiều lớp bài toán tối ưu trong kỹ thuật như: quy hoạch tuyến tính, quy hoạch phi tuyến không ràng buộc, quy hoạch phi tuyến có ràng buộc,Các thủ tục của tối ưu hoá đưa ra một cách chọn lựa các chiến lược tìm kiếm dựa trên các nhóm thuật toán: Medium – Scale Algorithms (MSA): Khi tìm cực tiểu hoá các hàm không ràng buộc chủ yếu sử dụng phương pháp tìm kiếm Nelder – Mead và phương pháp Newton – BFGS. Khi tìm cực tỉểu hoá các hàm có ràng buộc, cực tiểu tuyệt đối, quy hoạch đa mục tiêu và cực tiểu hoá trên bán không gian vô hạn, chủ yếu sử dụng các biến thể khác nhau của quy hoạch toàn phương liên tiếp. Khi giải quyết các vấn đề liên quan đến bài toán phương pháp bình phương nhỏ nhất phi tuyến thường sử dụng phương pháp Gauss – Newton và phương pháp Levenberg – Marquardt . Large – Scale Algorithms (LSA): Tất cả các phương pháp trong LSA (trừ quy hoạch tuyến tính) đều là các phương pháp miền tin cậy (Trust-region). Các bài toán có ràng buộc biên thường được giải theo phương pháp Newton. Các bài toán có ràng buộc dạng phương trình thường dựa trên các phương pháp lặp gradient. Thuật giải di truyền (Genetic Algorithm – GA): Giải pháp thích hợp cho các bài toán tối ưu tổ hợp (combinatorial optimization). Giải thuật di truyền là một phân ngành của giải thuật tiến hóa vận dụng các nguyên lý của tiến hóa như di truyền, đột biến, chọn lọc tự nhiên, và trao đổi chéo. Sử dụng các thuật toán có sẵn trong MATLAB giúp làm giảm được nhiều khối lượng công việc lập trình và tăng tốc độ tính toán, đặc biệt là khi các giải bài toán tối ưu phức tạp. Trong luận văn sẽ sử dụng thư viện tối ưu hóa của MATLAB để lập trình giải bài toán. Thuật toán giải bài toán dao động của cầu dầm liên tục chịu tác dụng của đoàn tải trọng di động Từ các công thức tính toán đã trình bày ở trên, trong mục này xây dựng thuật toán giải bài toán dao động của cầu dầm liên tục chịu tác dụng của đoàn tải trọng di động. Các bước tính toán được thực hiện theo trình tự sau: Bước 1: Nhập số liệu Nhập các số liệu mô tả kết cấu cầu, bao gồm: Số liệu nút: Số lượng nút Tọa độ nút trong hệ tọa độ tổng thể Số liệu phần tử Số lượng phần tử Số liệu cho mỗi phần tử: Liên kết phần tử (nút đầu, nút cuối) Đặc trưng vật liệu của phần tử: E, ρ Đặc trưng tiết diện của phần tử: J, F Số liệu liên kết: Số lượng liên kết Bậc tự do của các liên kết Độ cứng của các liên kết Nhập các số liệu đo giám sát dao động của cầu Bước 2: Xây dựng các ma trận phần tử hữu hạn cho các phần tử cầu Xây dựng các ma trận m, c, k của các phần tử Xây dựng các véc tơ tải trọng quy nút do các lực tập trung và phân bố tác dụng lên phần tử Đưa các ma trận và véc tơ của phần tử về hệ tọa độ tổng thể Bước 3: Tập hợp các ma trận phần tử hữu hạn cho toàn kết cấu Tập hợp ma trận M, C, K của kết cấu dầm liên tục trên cơ sở các ma trận m, c, k của phẩn tử Xây dựng véc tơ tải trọng quy nút của kết cấu trên cơ sở các véc tơ tải trọng quy nút của phần tử và các lực tập trung tác dụng trực tiếp lên nút của kết cấu. Bước 4: Đưa các điều kiện biên vào phương trình Xử lý các điều kiện biên: Gối cứng: loại bỏ các phương trình tương ứng với bậc tự do bị cản trở Gối đàn hồi: Thêm độ cứng của các bậc tự do tương ứng với liên kết vào ma trận độ cứng kết cấu K Thiết lập phương trình dao động của kết cấu với các ma trận M, C, K đã khử điều kiện biên Bước 5: Giải bài toán tối ưu Mô tả các tham số của xe cần xác định (m,c,k,v) Xây dựng hàm mục tiêu cho bài toán tối ưu. Xây dựng phương trình ràng buộc của bài toán tối ưu: Điều kiện ràng buộc mô tả tương tác của cầu và xe được xác định bằng chương trình VBIS đã lập trong chương 1 Giải bài toán tối ưu . Chương trình VBIS được sử dụng trong quá trình giải lặp của bài toán tối ưu. 2.4. Chương trình máy tính Trên cơ sở các thuật toán trình bày ở trên đã xây dựng chương trình nhận dạng thông số của xe di chuyển trên cầu từ các số liệu giám sát dao động cầu bằng ngôn ngữ MATLAB. Chương trình được đặt tên là IVPMB (Identification of Vehicle Parameters Moving on Bridge) Do điều kiện thời gian và tính phức tạp của bài toán, chương trình giới hạn cho bài toán nhận dạng thông số của xe chạy trên cầu cho trường hợp chỉ có một xe qua cầu. Sơ đồ thuật toán được trình bày trong hình 23. Hình 23 Sơ đồ khối chương trình IVPMB 2.5. Kết luận chương Trong chương này đã đạt được các kết quả chính sau: Đã trình bày cơ sở của phương pháp nhận dạng thông số xe chạy trên cầu từ các số liệu đo đạc giám sát dao động của cầu. Các thông số cần nhận dạng bao gồm khối lượng xe m, hệ số cản c và độ cứng k của hệ thống treo, vận tốc xe chạy trên cầu v. Thông số của xe xác định bằng cách giải bài toán tối với hàm mục tiêu là sai lệch giữa gia tốc chuyển động của cầu đo được khi xe chạy qua và gia tốc mô phỏng trên máy tính với các tam số xe giả định là nhỏ nhất. Đã trình bày phương pháp giải bài toán tối ưu bằng thuật giải di truyền. Đây là phương pháp hiện đại và có tốc độ tính toán nhanh. Đã triển khai thuật toán và lập chương trình bằng ngôn ngữ MATLAB để tính dao động kết cấu dầm liên tục chịu tải trọng động của đoàn xe di động bằng phương pháp phần tử hữu hạn. Các thử nghiệm số Trong khuôn khổ luận văn cao học, không có điều kiện để tiến hành các thí nghiệm tại hiện trường. Để có được các số liệu nghiên cứu thuật toán, đã tiến hành mô phỏng kết quả đo bằng cách sử dụng chương trình đã lập trong chương 1 để tính toán dao động của cầu tại các trạm đo giả định khi xe đi qua. Sau đó, sử dụng số liệu này như kết quả đo để nhận dạng các thông số của xe chạy qua cầu. Do hai quá trình này là độc lập với nhau nên kết quả nhận được là đáng tin cậy. Do hạn chế về thời gian, trong chương này chỉ khảo sát trường hợp chỉ có một xe di chuyển trên cầu với vận tốc không đổi. Chọn đối tượng khảo sát là cầu treo Đakrông. Đây là cầu treo dây xiên có thể mô hình hoá như cầu dầm liên tục đặt trên các gối tựa cứng và gối tựa đàn hồi. Hình 31 Cầu treo dây xiên Đakrông Giới thiệu cầu treo Đakrông Cầu treo dây xiên Đakrông có sơ đồ như trên hình 31. Các thông số của cầu như sau: Kết cấu nhịp: Dầm chủ cấu tạo từ hai dầm thép I 910, mặt cầu là bản BTCT M300, Độ cứng chống uốn của tiết diện EJ= 2,018.1011 N/m2, Chiều dài toàn bộ cầu cầu L=128,9 m, Khối lượng phân bố trên một đơn vị chiều dài cầu ρ=10472 kg/m. Dây cáp treo: Dây cáp treo được làm từ cáp mềm theo ГОСТ3079-55 Cáp 1: tương đương với gối đàn hồi k1= 5,57.108 N/m, Cáp 2: tương đương với gối đàn hồi k2= 3,90.108 N/m, Cáp 3: tương đương với gối đàn hồi k3= 2,76.108 N/m, Cáp 4: tương đương với gối đàn hồi k4= 2,76.108 N/m, Cáp 5: tương đương với gối đàn hồi k5= 3,90.108 N/m, Cáp 6: tương đương với gối đàn hồi k6= 5,57.108 N/m. Các hệ số cản tỷ lệ của kết cấu: Hệ số cản tỷ lệ với khối lượng α=0,5%. Hệ số cản tỷ lệ với độ cứng β=0,5%. Sự tương tác của cầu Đakrông với tải trọng xe chạy được mô hình hoá cầu bằng mô hình dầm liên tục đặt trên các gối tựa cứng và gối tựa đàn hồi và có các khối lượng di động trên cầu liên kết với cầu bằng các lò xo đàn hồi và cản nhớt tuyến tính (hình 32). Hình 32 Mô hình tính toán cầu treo dây xiên Đakrông Trường hợp 1: Trên cầu chỉ có 1 trạm đo dao động Đo dao động trên cầu do xe chạy qua Trong trường hợp nghiên cứu này, chỉ bố trí một trạm đo tại vị trí giữa cầu như hình 33. Hình 33 Sơ đồ bố trí trạm đo trên cầu Giả thiết, xe di chuyển trên cầu có các thông số như sau: Khối lượng m=10000 kg, Độ cứng liên kết lò xo k=9,02.105 N/m, Hệ số cản nhớt c=1,1.105 Ns/m, Vận tốc di chuyển v=50 km/h, Sử dụng chương trình VBIS đã lập để tính toán mô phỏng bản ghi giá trị đo gia tốc theo phương thẳng đứng tại vị trí trạm đo khi xe đi qua cầu. Kết quả tính toán mô phỏng gia tốc dao động của cầu được biểu diễn trên hình 34. Hình 34 Gia tốc dao động đo được tại vị trí giữa cầu Xác định các thông số của xe chạy trên cầu từ số liệu đo Tiếp tục, sử dụng kết quả tính toán mô phỏng gia tốc dao động của cầu tại vị trí trạm đo làm đầu vào cho chương trình IVPMB. Chạy chương trình nhận được kết quả tính như bảng sau Bảng 31 Kết quả tính toán TT Thông số xe Đơn vị Số liệu thực Số liệu tính được Sai số (%) Khối lượng m () Kg 10000 9986 1,40 Độ cứng liên kết lò xo k N/m 9,02.105 8,97.105 0,54 Hệ số cản nhớt c Ns/m 1,1.105 1,09.105 0,82 Vận tốc v km/h 50 49,38 1,24 Nhận xét Kết quả tính toán cho thấy với trường hợp nghiên cứu này kết quả xác định các thông số của xe có thể xem là chính xác. Từ đó, nhận thấy rằng phương pháp xác định các thông số của xe được lựa chọn là hiệu quả và chương trình đã lập là đáng tin cậy. Trường hợp 2: Trên cầu có 3 trạm đo dao động Đo dao động trên cầu do xe chạy qua Trong trường hợp nghiên cứu này, bố trí ba trạm đo tại các vị trí trên cầu như hình 33. Hình 35 Sơ đồ bố trí trạm đo trên cầu Giả thiết, xe di chuyển trên cầu có các thông số như sau: Khối lượng m=10000 kg, Độ cứng liên kết lò xo k=9,02.105 N/m, Hệ số cản nhớt c=1,1.105 Ns/m, Vận tốc di chuyển 50 km/h, Sử dụng chương trình VBIS đã lập để tính toán mô phỏng bản ghi giá trị đo gia tốc theo phương thẳng đứng tại vị trí trạm đo khi xe đi qua cầu. Kết quả tính toán mô phỏng gia tốc dao động của cầu được biểu diễn trên các sau. Hình 36 Gia tốc dao động đo được tại vị trí trạm đo 1 Hình 37 Gia tốc dao động đo được tại vị trí trạm đo 2 Hình 38 Gia tốc dao động đo được tại vị trí trạm đo 3 Xác định các thông số của xe chạy trên cầu từ số liệu đo Tiếp tục, sử dụng kết quả tính toán mô phỏng gia tốc dao động của cầu tại vị trí trạm đo làm đầu vào cho chương trình IVPMB. Chạy chương trình nhận được kết quả tính như bảng sau Bảng 32 Kết quả tính toán TT Thông số xe Đơn vị Số liệu thực Số liệu tính được Sai số (%) Khối lượng m () Kg 10000 9938 0,62 Độ cứng liên kết lò xo k N/m 9,02.105 9,00.105 0,22 Hệ số cản nhớt c Ns/m 1,1.105 1,1.105 0,20 Vận tốc v km/h 50 49,63 0,74 Nhận xét Kết quả tính toán cho thấy với trường hợp nghiên cứu này kết quả xác định các thông số của xe sát hơn với thông số thực của xe. Từ đó, nhận thấy rằng khi số lượng trạm đo tăng lên thì kết quả nhận dạng cũng chính xác hơn. Kết luận chương Trong chương này đã đạt được các kết quả chính sau: Sử dụng các chương trình đã lập để mô phỏng dao động của cầu khi có xe chạy qua Lấy các số liệu gia tốc dao động của cầu theo phương thẳng đứng tại vị trí các trạm đo làm đầu vào, nhận dạng trở lại các thông số của xe, kết quả cho thấy phương pháp nhận dạng và chương trình đã lập là hiệu quả và đáng tin cậy. Qua các thí dụ tính toán cho thấy, khi số trạm đo dao động tăng lên thì độ chính xác khi nhận dạng các thông số của xe qua cầu cũng tăng lên. Kết luận Các kết quả chính đã đạt được trong luận án gồm: Đã thiết lập được phương trình vi phân mô tả dao động của dầm chịu tác dụng của nhiều tải trọng di động đồng thời. Đã xây dựng các ma trận phần tử hữu hạn cho phần tử dầm mang nhiều tải trọng di động. Đã thiết lập phương trình dao động của kết cấu cầu dầm liên tục chịu tác dụng của đoàn tải trọng di động. Đã lập chương trình tính toán mô phỏng dao động của kết cấu cầu dầm liên tục chịu tác dụng của đoàn tải trọng di động bằng ngôn ngữ MATLAB. Chương trình đảm bảo độ tin cậy Đã lập chương trình nhận dạng thông số xe chạy trên cầu bằng ngôn ngữ MATLAB. Sử dụng chương trình đã lập để nhận dạng thông số xe chạy trên cầu. Kết quả đảm bảo tin cậy. Hướng nghiên cứu tiếp theo: Nghiên cứu ảnh hưởng của nhiễu trong kết quả đo dao động đến kết quả nhận dạng thông số xe chạy trên cầu. Nghiên cứu ảnh hưởng của sự thay đổi vận tốc xe chạy, lực phanh xe đến các các kết quả đo đạc. Nghiên cứu nhận dạng thông số xe chạy trên cầu với trường hợp đoàn xe di chuyển trên cầu. Tài liệu tham khảo Tiếng Việt 1. Hoàng Hà, Trần Quang Vinh (1999), Nghiên cứu dao động uốn của kết cấu nhịp cầu dây văng dưới tác dụng của đoàn tải trọng ô tô di động, Đề tài nghiên cứu khoa học cấp bộ, Đại học GTVT, Hà Nội. 2. Hoàng Hà (2002), Nghiên cứu dao động uốn phi tuyến của cầu dây văng chịu tác dụng của đoàn tải trọng ô tô, Đề tài nghiên cứu khoa học cấp bộ, Đại học GTVT, Hà Nội. 3. Đỗ Xuân Thọ (1996), Tính toán dao động uốn của dầm liên tục chịu tác dụng của vật thể di động, Luận án PTS KHKT, Hà Nội. 4. Hồ Anh Tuấn, Trần Bình (1978), Phương pháp phần tử hữu hạn, Nhà xuất bản KHKT, Hà Nội. Người dịch: Phạm Hồng Giang. 5. Tạ Hữu Vinh (2005), Nghiên cứu dao động của hệ thanh chịu tải trọng di động bằng phương pháp số, Luận án tiến sĩ kỹ thuật, Học viện KTQS, Hà Nội. 6. Nguyễn Mạnh Yên (2000), Phương pháp số trong cơ học kết cấu, Nhà xuất bản KHKT, Hà Nội. Tiếng Anh 7. Bathe K. J., Wilson E. L. (1976), Numerical Methods in Finite Element Analysis, Prentice - Hall, INC. 8. Gou W.H. and Xu Y.L. (2000), Direct Assembling Matrix Method for Dynamic Analysis of Coupled Vehicle-Bridge System, Advances in structural dynamics, Vol. 1, Hongkong. 9. W.F Chen and L. Duan (2000), Bridge Engineering Handbook, CRC Press. 10. S.H Ju, H.T Lin, C.C Hsueh and S.L Wang (2006), A simple finite element model for vibration analyses induced by moving vehicles, International Journal for Numerical Methods in Engineering. 11. William T. Thomson (1988), Theory of Vibration With Applications, Prentice Hall, New Jersey. Phụ lục Chương trình tính % khai_bao_chung.m %% Dinh nghia cac so lieu global t x phan_tu xe global lket alpha beta % Gia toc trong truong g = 9.81; %m/s2; %% Cac hang so can ty le M va K mac dinh alpha = 0; beta = 0; %% Dinh nghia cau truc du lieu s_ptu = struct('nut_1',0,'nut_2',0,'EJ',0,'rho',0,'p',0); s_xe = struct('m',0,'k',0,'c',0,'v',0,'x0',0,'t0',0); s_lket = struct('nut',0,'chieu',0,'k',0); % DaoDongCau.m %% Mo ta so lieu % 1. So lieu nut % x=[x1 x2 ... xN] % 2. So lieu phan tu % phantu = [ptu1, pt2,...] % ptu = struct(I, J, EJ, rho, p); % 3. So lieu xe % xe = struc (m, k, c, v, x0, t0); % khoi luong - do cung - can - van toc - vi tri bat dau - thoi diem xe bat dau chay % 5. Danh sach lket dan hoi % lket = [g1 g2 ... gN] % g = struct(nut, k); % neu k=0 thi lket la lket cung %% Khai_bao_chung; %% Phan tich so lieu % Tim so luong nut va phan tu so_nut = length(x); so_ptu = length(phan_tu); so_xe = length(xe); % Danh so bac tu do cua ket cau va xe so_btd = so_nut*2+so_xe; so_btd_kc = so_nut*2; btd_kc = (1:so_btd_kc); btd_ptu = zeros(so_ptu,4); for i=1:so_ptu nut_1 = phan_tu(i).nut_1; nut_2 = phan_tu(i).nut_2; btd_ptu(i,:)=[ (nut_1-1)*2+[1 2] (nut_2-1)*2+[1 2] ]; end; btd_xe=so_nut*2+(1:so_xe);% Bac tu do cua khoi luong m la N+1 tro di % Tinh toan chieu dai phan tu for i=1:so_ptu nut_1 = phan_tu(i).nut_1; nut_2 = phan_tu(i).nut_2; phan_tu(i).L=abs(x(nut_2)-x(nut_1)); end; %% Dinh nghia cac bien cua bai toan % Cac ma tran tong the cua ket cau ko ke den xe M0 = zeros(so_btd); C0 = zeros(so_btd); K0 = zeros(so_btd); F0 = zeros(so_btd,1); if exist('Ptt','var') for i=1:length(Ptt) F0((Ptt(i).nut-1)*2+Ptt(i).chieu) = Ptt(i).P; end; end; %% Xay dung cac ma tran phan tu huu han cho phan tu dam for i=1:so_ptu EJ = phan_tu(i).EJ; L = phan_tu(i).L; rho= phan_tu(i).rho; Ke = EJ/L^3* [ 12 6*L -12 6*L; 6*L 4*L^2 -6*L 2*L^2; -12 -6*L 12 -6*L; 6*L 2*L^2 -6*L 4*L^2]; Me = rho*L/420*[ 156 22*L 54 -13*L ; 22*L 4*L^2 13*L -3*L^2 ; 54 13*L 156 -22*L; -13*L -3*L^2 -22*L 4*L^2]; Ce = alpha*Me + beta*Ke; M0(btd_ptu(i,:),btd_ptu(i,:)) = M0(btd_ptu(i,:),btd_ptu(i,:)) + Me; C0(btd_ptu(i,:),btd_ptu(i,:)) = C0(btd_ptu(i,:),btd_ptu(i,:)) + Ce; K0(btd_ptu(i,:),btd_ptu(i,:)) = K0(btd_ptu(i,:),btd_ptu(i,:)) + Ke; end; %% Bo sung cac ma tran cua xe vao ma tran tong the for ix=1:so_xe i = btd_xe(ix); K0(i,i)=xe(ix).k; C0(i,i)=xe(ix).c; M0(i,i)=xe(ix).m; end; %% Xu ly dieu kien bien lket cung va lket dan hoi for i=1:so_lket nut = lket(i).nut; chieu = lket(i).chieu; k = lket(i).k; btd_lket = (nut-1)*2+chieu; if k>0 K0(btd_lket,btd_lket) = K0(btd_lket,btd_lket) + k; else btd_kc(btd_kc==btd_lket)=[]; end; end; %% Phan tich dao dong cua cau khi co xe chay Mx = zeros(5); Cx = zeros(5); Kx = zeros(5); nt = length(t); btd = [btd_kc btd_xe]; for i=1:nt-1 % Lay cac gia tri ban dau cua ma tran M = M0; C = C0; K = K0; if exist('Ft','var') F = F0*Ft(i); else F = F0; end; % Thiet lap cac ma tran cua phan tu xe va tap hop vao MTTT for ix=1:so_xe % xac dinh vi tri xe if t(i)<xe(ix).t0 xi = xe(ix).x0; else xi = xe(ix).x0+xe(ix).v*(t(i)-xe(ix).t0); end; % tim phan tu co xe chay tren xe_tren_cau = 0; for ipt=1:so_ptu; nut_1 = phan_tu(ipt).nut_1; nut_2 = phan_tu(ipt).nut_2; if x(nut_1) <= xi && xi <= x(nut_2) xe_tren_cau = 1; break; end; end; if xe_tren_cau % Vi tri xe trong he toa do phan tu xi = xi-x(nut_1); % ham dang phan tu tai vi tri xi L = phan_tu(ipt).L; N = [ 1-3*xi^2/L^2+2*xi^3/L^3 ... xi-2*xi^2/L+xi^3/L^2 ... 3*xi^2/L^2-2*xi^3/L^3 ... -xi^2/L^2+xi^3/L^3 ]; % Chi so bac tu do cua he xe+dam btd_dx = [btd_ptu(ipt,:) btd_xe(ix)]; % Cac ma tran phan tu xe tai vi tri xi Mx(:,5) = [N'*xe(ix).m; 0]; Cx(5,:) = xe(ix).c*[-N 0]; Kx(5,:) = xe(ix).k*[-N 0]; % Tap hop ma tran phan tu xe vao ma tran tong the M(btd_dx,btd_dx) = M(btd_dx,btd_dx) + Mx; C(btd_dx,btd_dx) = C(btd_dx,btd_dx) + Cx; K(btd_dx,btd_dx) = K(btd_dx,btd_dx) + Kx; F(btd_dx) = F(btd_dx)+[-N'*xe(ix).m*g; 0]; end; end; % tap hop xong so lieu xe % Giai he phuong trinh dao dong bang phuong phap New-Mark if i==1 U = zeros(so_btd, nt); V = zeros(so_btd, nt); A = zeros(so_btd, nt); u0 = zeros(so_btd, 1); v0 = zeros(so_btd, 1); a0 = zeros(so_btd, 1); a0(btd_kc) = M(btd_kc,btd_kc)\(F(btd_kc)-C(btd_kc,btd_kc)*v0(btd_kc)-K(btd_kc,btd_kc)*u0(btd_kc)); U(:,1) = u0; V(:,1) = v0; A(:,1) = a0; dt = t(2)-t(1); end; Keff = M + 0.5*dt*C + 0.25*(dt^2)*K; A(btd,i+1) = Keff(btd,btd)\(F(btd) - C(btd,btd)*(V(btd,i) + 0.5*dt*A(btd,i)) ... - K(btd,btd)*(U(btd,i) + dt*V(btd,i) + 0.25*(dt^2)*A(btd,i))); V(btd,i+1) = V(btd,i) + 0.5*dt*A(btd,i) + 0.5*dt*A(btd,i+1); U(btd,i+1) = U(btd,i) + dt*V(btd,i) + dt^2*0.25*A(btd,i) + dt^2*0.25*A(btd,i+1); end; %Vi_du_Gou&Xu clear variables; close all; clc; %% Khai_bao_chung; %% 1. So lieu nut L_cau=25; so_nut = 55; x=(0:L_cau/(so_nut-1):25); %% 2. So lieu phan tu so_ptu = so_nut-1; phan_tu = repmat(s_ptu,so_ptu,1); for i=1:so_ptu phan_tu(i).nut_1 = i; phan_tu(i).nut_2 = i+1; phan_tu(i).EJ = 2.87e9*2.9; % N/m2 x m4 phan_tu(i).rho = 2303; % (kg/m); end; %% 3. So lieu xe so_xe=30; xe = repmat(s_xe,so_xe,1); % So lieu chung cho cac xe for i=1:so_xe xe(i).m = 5750; % kg xe(i).k = 1.595e6; % N/m xe(i).c = 0; xe(i).v = 10*1000/60/60; % m/s = 100 km/h xe(i).x0 = 0; %(i-1)*-25; % m ; Vi tri xuat phat cua xe so voi diem dau cau xe(i).t0 = (i-1)*L_cau/xe(i).v; end; %% 4. Danh sach lket tua so_lket=2; lket = repmat(s_lket,so_lket,1); lket(1).nut=1; lket(1).chieu=1; lket(1).k=0; lket(2).nut=so_nut; lket(2).chieu=1; lket(2).k=0; %% 5. Thoi gian tmax = (x(end)-x(1)-xe(end).x0)/xe(1).v; t = (0:0.05:20); %% Tinh dao dong cua dam DaoDongCau; %% Ve cac do thi % Ve_do_thi_chuyen_vi_nut; nut = round(so_nut/2+.25); chieu = 1; dof = (nut-1)*2+chieu; figure(1); plot(t,U(dof,:)*1000); grid on; title(['Do vong tai diem giua cau, nut ' num2str(nut)]); xlabel('t (s)'); ylabel('u (mm)'); % % Ve_do_thi_chuyen_vi_xe; figure(2); plot(t*xe(1).v,U(btd_xe(:),:)*1000); grid on; title('Chuyen vi cua khoi luong \itm'); xlabel('\xi (m)'); ylabel('w (mm)'); % % Ve Hinh dang cau; figure(3); plot(x,U(1:2:so_nut*2,10:10:100)*1000); grid on; title('Duong cong bien dang cua cau'); xlabel('x (m)'); ylabel('u_m (mm)'); %Vi_du_1 clear variables; close all; clc; %% Khai_bao_chung; %% 1. So lieu nut L_cau = 20; so_nut = 21; x=(0:L_cau/(so_nut-1):20); %% 2. So lieu phan tu so_ptu = so_nut-1; phan_tu = repmat(s_ptu,so_ptu,1); for i=1:so_ptu phan_tu(i).nut_1 = i; phan_tu(i).nut_2 = i+1; phan_tu(i).EJ = 4157e6; phan_tu(i).rho = 2277; end; % 3. So lieu xe so_xe=1; xe = repmat(s_xe,so_xe,1); xe.m = 4333; xe.k = 902000; xe.c = 11016; xe.v = 20; xe.x0 = 0; xe.t0 = 0; % 4. Danh sach lket tua so_lket=2; lket = repmat(s_lket,so_lket,1); lket(1).nut=1; lket(1).chieu=1; lket(1).k=0; lket(2).nut=so_nut; lket(2).chieu=1; lket(2).k=0; % 5. Thoi gian t = (0:0.001:1); ft = sin(50*t); % 6.Cac hang so can ty le M va K alpha = 0; beta = 0; %% Tinh dao dong cua dam DaoDongCau; %% Ve cac do thi % Ve_do_thi_chuyen_vi_nut; nut = round(so_nut/2+.25); chieu = 1; dof = (nut-1)*2+chieu; figure(1); plot(t,[U(dof,:); U(btd_xe(ix),:)]*1000); grid on; title(['Do vong tai diem giua cau (nut ' num2str(nut) ') va khoi luong \itm']); xlabel('t (s)'); ylabel('u, w (mm)'); %Vi du cau Dakrong clear variables; close all; %clc; %% Khai_bao_chung; %% 1. So lieu nut x=[0 3.68 7.36 11.04 14.72 18.4 ... 22.06 25.72 29.38 33.04 36.7 ... 40.36 44.02 47.68 51.34 55 ... 58.78 62.56 66.34 70.12 73.9 ... 77.56 81.22 84.88 88.54 92.2 ... 95.86 99.52 103.18 106.84 110.5 ... 114.18 117.86 121.54 125.22 128.9]; so_nut = length(x); %% 2. So lieu phan tu so_ptu = so_nut-1; phan_tu = repmat(s_ptu,so_ptu,1); for i=1:so_ptu phan_tu(i).nut_1 = i; phan_tu(i).nut_2 = i+1; phan_tu(i).EJ = 2.018e11; % N/m2 x m4 phan_tu(i).rho = 10472; % (kg/m); end; %% 3. So lieu xe so_xe=5; xe = repmat(s_xe,so_xe,1); % So lieu chung cho cac xe for i=1:so_xe xe(i).m = 20000; % kg xe(i).k = 9.0e10; % N/m xe(i).c = 1.10e5; % Ns/m xe(i).v = 10*1000/60/60; % m/s = 10 km/h xe(i).x0 = (i-1)*-10; % m ; Vi tri xuat phat cua xe so voi diem dau cau xe(i).t0 = 0; end; %% 4. Danh sach lket tua so_lket=8; lket = repmat(s_lket,so_lket,1); for i=1:so_lket lket(i).nut=(i-1)*5+1; lket(i).chieu=1; lket(i).k=0; end; lket(2).k=5.57e8; % N/m lket(3).k=3.90e8; % N/m lket(4).k=2.76e8; % N/m lket(5).k=2.76e8; % N/m lket(6).k=3.90e8; % N/m lket(7).k=5.57e8; % N/m %% 5. Thoi gian nt = 800; tmax = 1.2*(x(end)-xe(end).x0)/xe(1).v+xe(end).t0; t=(1:tmax/nt:tmax); nt = length(t); %% 6.Cac hang so can ty le M va K alpha = 0; beta = 0.001; %% Tinh dao dong cua cau DaoDongCau; %% Ve cac do thi % Ve_do_thi_chuyen_vi_nut; nut = round(so_nut/2+.25); chieu = 1; dof = (nut-1)*2+chieu; figure(1); %plot(t*xe(1).v,U(dof,:)*1000); plot(t,U(dof,:)*1000); grid on; title('Do vong tai diem giua cau'); xlabel('t (s)'); ylabel('u (mm)'); % % Ve_do_thi_chuyen_vi_xe; ix = 1; figure(2); plot(t,U(btd_xe(:),:)*1000); grid on; title('Chuyen vi cua khoi luong \itm'); xlabel('t (s)'); ylabel('u_m (mm)'); % Ve Hinh dang cau; figure(3); i = 200; plot(x,U(1:2:so_nut*2,i)*1000); grid on; title(['Duong cong bien dang cua cau tai lan tinh thu ' num2str(i)]); xlabel('x (m)'); ylabel('u_m (mm)'); Mục lục Danh mục các hình vẽ Hình 11 Phần tử dầm chịu uốn thuần tuý 14 Hình 21 Mô hình tương tác của khối lượng di động (xe) với phần tử dầm 32 Hình 23 Mô hình khảo sát của Gou & Xu 45 Hình 24 So sánh độ võng tại điểm giữa dầm với kết quả của Gou & Xu 46 Hình 25 So sánh chuyển vị của khối lượng với kết quả của Gou & Xu 47 Hình 26 Xe chạy trên dầm đơn giản 48 Hình 27: Độ võng giữa nhịp và chuyển vị của tải trọng 49 Danh mục các bảng Bảng 31 Kết quả tính toán 67 Bảng 31 Kết quả tính toán 70