Đề tài Nghiên cứu sự ổn định khoang hầm trong môi trường đá nứt nẻ bằng phương pháp Phân tích biến dạng không liên tục

1 MỞ ĐẦU 1. Tính cấp thiết của đề tài nghiên cứu Một trong những vấn đề đặt ra cho việc xây dựng công trình ngầm trong đá là nghiên cứu, đánh giá, phân tích ổn định các khoảng trống ngầm, không gian ngầm nhằm có được thiết kế hợp lý về kết cấu chống đỡ, kết cấu công trình và biện pháp thi công. Trong những năm gần đây, để khắc phục những khó khăn của các lời giải giải tích cũng như phương pháp thực nghiệm và thí nghiệm, các nhà nghiên cứu đã sử dụng nhiều phương pháp số khác nhau. Phương pháp Phân tích biến dạng không liên tục DDA (Discontinuous Defor mation Analysis) là phương pháp số được sử dụng để phân tích lực tương tác và chuyển dịch khi các khối tiếp xúc với nhau. Đối với mỗi khối, DDA cho phép xác định các chuyển dịch, biến dạng ở mỗi bước thời gian; đối với toàn bộ hệ các khối thì cho phép mô phỏng quá trình tiếp xúc, tương tác giữa các khối. Với các lí do trên, đề tài nghiên cứu của luận án được chọn là “Nghiên cứu sự ổn định khoang hầm trong môi trường đá nứt nẻ bằng phương pháp Phân tích biến dạng không liên tục”. 2. Mục đích, nội dung, phƣơng pháp, phạm vi nghiên cứu của luận án  Mục đích của luận án Xây dựng mô hình, thuật toán và chương trình để xác định các trường chuyển dịch, ứng suất và biến dạng của khối đá theo thời gian xung quanh khoang hầm trong môi trường biến dạng không liên tục. Thông qua các nghiên cứu lý thuyết và các thử nghiệm số trên máy tính, phân tích ảnh hưởng của trạng thái nứt nẻ khối đá đến tính ổn định của kết cấu công trình ngầm.  Nội dung nghiên cứu của luận án 2 1. Tìm hiểu và sử dụng phương pháp Phân tích biến dạng không liên tục DDA. 2. Xây dựng mô hình tính và thuật toán cùng việc thiết lập chương trình tính toán chuyển dịch, biến dạng và ứng suất theo DDA. 3. Tiến hành một số tính toán, thử nghiệm số phân tích chuyển dịch của khối đá nứt nẻ xung quanh khoang hầm và sự tiếp xúc, tương tác giữa công trình ngầm với môi trường đá nứt nẻ.  Phƣơng pháp nghiên cứu của luận án Nghiên cứu lý thuyết kết hợp với thử nghiệm số trên máy tính.  Phạm vi nghiên cứu của luận án Xét mô hình tính là các bài toán phẳng trong môi trường không liên tục. 3. Cấu trúc của luận án Cấu trúc của luận án bao gồm phần mở đầu, bốn chương và phần kết luận, cuối cùng là tài liệu tham khảo và phụ lục. Nội dung luận án gồm 120 trang, 19 bảng biểu, 92 hình vẽ và đồ thị, 27 tài liệu tham khảo, 05 bài báo khoa học phản ánh nội dung của luận án. Phần phụ lục trình bày mã nguồn của các chương trình đã lập trong luận án. CHƢƠNG I TỔNG QUAN Trong chương này đã tiến hành tổng quan các nghiên cứu về sự ổn định khối đá xung quanh khoang hầm và một số phương pháp số áp dụng trong môi trường không liên tục. Ứng dụng nghiên cứu này trong xây dựng công trình ngầm trong môi trường đá nứt nẻ cho phép đánh giá tương tác giữa môi trường và công trình để từ đó có những giải pháp hợp lý giúp cho việc xây dựng an toàn, hiệu quả và chất lượng. Các kết luận rút ra trong chương tổng quan là: 3  Lý thuyết về nghiên cứu ổn định công trình ngầm cũng như áp lực địa tầng tác dụng lên công trình được phát triển rất đa dạng, từ lâu. Bằng các nghiên cứu của mình các nhà khoa học đã có những đóng góp to lớn trong việc xây dựng hệ thống công trình ngầm trong các môi trường khác nhau đặc biệt là môi trường đá nứt nẻ.  Trong việc phân tích ổn định khoang hầm hiện nay có hai phương pháp chủ yếu là: phương pháp giải tích và phương pháp số. Trong đó phương pháp số là phương pháp có thể mô phỏng được điều kiện bài toán gần sát với làm việc thực tế của kết cấu và môi trường. Đối với các bài toán trong môi trường rời, nhóm theo quan điểm mô hình không liên tục có những ưu thế vượt trội so với nhóm theo quan điểm môi trường liên tục. Phương pháp phân tích biến dạng không liên tục DDA là một trong những phương pháp số nghiên cứu các bài toán cơ học biến dạng không liên tục, đặc biệt được áp dụng có hiệu quả trong các bài toán về cơ học đá. CHƢƠNG II PHƢƠNG PHÁP PHÂN TÍCH BIẾN DẠNG KHÔNG LIÊN TỤC (DDA) 2.1 Phƣơng pháp DDA và quá trình phát triển Phương pháp DDA nghiên cứu tính toán chuyển dịch, ứng suất và biến dạng các khối trong môi trường không liên tục; trong đó chú trọng nhất vào việc nghiên cứu tiếp xúc và tương tác giữa các khối với nhau trong cơ hệ. Phân tích biến dạng không liên tục do G.H Shi và R.E. Goodman [20],[21]giới thiệu vào những năm 1984, 1985. Tuy nhiên, DDA chính thức trở thành phương pháp được mọi người biết đến năm 1988 [22]. Mặc dù các tài liệu về DDA khá phổ biến trên các mạng thông 4 tin nhưng các phần mềm ứng dụng lại ít được giới thiệu. Tại Việt Nam, DDA còn ít được nghiên cứu và giới thiệu trong các chương trình giảng dạy cũng như các nghiên cứu, báo cáo khoa học. 2.2 Nội dung cơ bản của phƣơng pháp DDA 2.2.1 Chuyển dịch và biến dạng của khối đơn Xét cơ hệ trong hệ tọa độ Descartes xOy , trong trường hợp tổng quát của bài toán phẳng, trạng thái chuyển động của khối được xác định bởi 3 thành phần: hai thành phần chuyển động tịnh tiến u,v và một thành phần chuyển động quay r ; trạng thái biến dạng gồm 3 thành phần: hai thành phần biến dạng thẳng x , y và một thành phần biến dạng góc xy . Như vậy, chuyển vị (u,v) tại một điểm bất kỳ có tọa độ (x, y) của khối có thể được biểu diễn qua 6 thành phần chuyển vị và biến dạng 0 0 0 x y xy(u v r )   tại một điểm xác định (xo,yo) thuộc khối. Trong đó: 0 0(u ,v ) là chuyển vị tại một điểm cụ thể 0 0(x ,y ) của khối; 0r là góc quay của khối với tâm quay tại 0 0(x ,y ) ; x , y , xy là biến dạng thẳng và biến dạng góc của khối. Bằng việc biểu diễn chuyển dịch (u,v) tại một điểm bất kỳ (x,y) của khối bởi đa thức bậc nhất. Sau khi biến đổi ta có công thức xác định chuyển dịch (u,v) tại một điểm bất kỳ (x,y) qua 6 thành phần chuyển vị và biến dạng 0 0 0 x y xy(u v r )   tại một điểm xác định (xo,yo) thuộc khối dưới dạng ma trận như sau:          i i u T D v (2.11) trong đó: 0 0 0i 0 0 0 1 0 (y y ) (x x ) 0 (y y ) / 2 [T ] 0 1 (x x ) 0 (y y ) (x x ) / 2                 T i 0 0 0 x y xyD u v r    5 2.2.2 Hệ phƣơng trình chuyển động của cơ hệ Hệ phương trình tổng quát của DDA được xây dựng theo nguyên lý cực tiểu cơ năng toàn phần. Hệ phương trình tổng quát của DDA cho một cơ hệ bao gồm n khối được biểu diễn dưới dạng ma trận: [K][D]=[F] (2.14) Ma trận [K] được gọi là ma trận độ cứng tổng thể; ở đây, mỗi phần tử trên đường chéo chính iiK là một ma trận con ii[K ] phụ thuộc vào tính chất cơ học của khối thứ i, các ma trận con ij[K ] với i j được xác định khi khối thứ i tiếp xúc với khối thứ j;  iD là véc tơ chuyển vị của khối thứ i 1i 2i 3i 4i 5i 6id d d d d d , iF là tải trọng tác dụng lên khối thứ i (bao gồm lực quán tính, tải trọng ngoài, lực dính kết, lực khối, điều kiện tiếp xúc). Trong DDA, sau mỗi một bước tích phân, vị trí tương đối giữa các khối trong cơ hệ sẽ thay đổi, hệ lực tác dụng lên mỗi khối cũng thay đổi, vì vậy phương trình (2.14) sẽ được thiết lập lại, hay nói cách khác mỗi một hệ phương trình chuyển động chỉ được dùng cho một bước tích phân. Như vậy, hệ phương trình chuyển động cho cơ hệ sẽ được xây dựng theo hai bước: + Thiết lập phương trình chuyển động cho khối đơn. + Tiếp xúc và tương tác giữa các khối. 2.2.3 Phƣơng trình chuyển động khối đơn Phương trình chuyển động của khối đơn thứ i được biểu diễn theo công thức (2.14), lúc này ma trận ij[K ] với i j là các ma trận 0. Tổng cơ năng của hệ  được xác định theo nguyên lý cộng tác dụng . Những năng lượng này được tính riêng rẽ, sau đó được lấy đạo hàm từng phần, các ma trận con (năng lượng thành phần) thu được sẽ đưa vào thành phần của ma trận ii[K ] và véc tơ i{F} trong phương trình (2.14). Các trường hợp cụ thể được xác định như sau: 6 2.2.3.1 Ma trận con biến dạng đàn hồi Thế năng biến dạng đàn hồi của một khối thứ i là: e x x y y xy xy 1 ( )dxdy 2           ( 2.18) Đạo hàm thế năng biến dạng đàn hồi của khối theo các thành phần chuyển vị và biến dạng của khối ta nhận được:    i iiS E K sẽ được đưa vào ma trận ii[K ] trong ma trận độ cứng tổng thể [K] . E và  lần lượt là mô đun đàn hồi và hệ số Poisson của vật liệu khối. 2.2.3.2 Véc tơ tải trọng ứng với ứng suất ban đầu Thế năng tạo ra bởi ứng suất ban đầu  0 0 0x y xy   của khối thứ i : 0 0 0 x x y y xy xy( )dxdy           (2.23) Đạo hàm  theo các thành phần chuyển vị và biến dạng của khối ta được véc tơ 6 thành phần : 0 i S{ } {F}  sẽ được bổ sung vào i{F} trong phương trình (2.14). 2.2.3.3 Véc tơ tải trọng ứng với tải trọng tập trung Giả sử khối thứ i chịu tác dụng của tải trọng tập trung ( xF , yF ) tác dụng tại điểm (x,y). Thế năng được tạo ra bởi tải trọng tập trung sẽ có dạng như sau:        TTx xp x y i i y y F F F u F v u v D T (x, y) F F                          ( 2.26) Đạo hàm phương trình (2.26) cho ta véc tơ 6 thành phần mô tả lực tác dụng vào khối: T x11 12 13 14 15 16 i y21 22 23 24 25 26 Ft t t t t t {F} Ft t t t t t           (2.28) sẽ được bổ sung vào véc tơ i{F} trong phương trình tổng thể (2.14). 2.2.3.4 Véc tơ tải trọng ứng với tải trọng phân bố theo đường Giả sử khối thứ i chịu tải trọng phân bố có cường độ thay đổi 7 dọc theo đường phân bố (phương trình tham số) là: x xF F (t) , y yF F (t) 0 t 1  trên một đoạn thẳng với chiều dài l. Thế năng tạo bởi tải trọng phân bố x y(F (t),F (t)) được biểu diễn:       l l TTx x l i i y y0 0 F (t) F (t) u v ldt D T ldt F (t) F (t)                         (2.31) Đạo hàm l nhận được véc tơ 6x1:   l T x i i y 0 F (t) T ldt {F} F (t)          được bổ sung vào véc tơ i{F} trong phương trình tổng thể (2.14). 2.2.3.5 Ma trận con tạo bởi lực quán tính Lực quán tính trên đơn vị diện tích của khối thứ i được xác định qua chuyển vị theo thời gian  u(t),v(t) tại một điểm bất kỳ (x,y) và M là khối lượng trên đơn vị diện tích sẽ là: 2u(t ) 2x t 2v(t )y 2t f M f                      (2.35) Thế năng i tạo ra bởi lực quán tính được xác định:        2TT i i i i 2 D(t) M D T T dxdy t     (2.37) Bằng cách lấy đạo hàm theo thời gian, ta có được:      22 i 2 D(t) D(t)t D t 2 tt       (2.38)           2 i i 02 2 2 D(t) D(t)2 2 2 2 D D V t t tt t t             (2.39) Thay (2.39) vào (2.37) ta có được             TT i i i i i 02 2M 2M D T T dxdy( D V ) tt     (2.41) Lấy đạo hàm i theo các giá trị chuyển vị và biến dạng ta được: 8     T i i ii2 2M T T dxdy [K ] t    được đưa vào ma trận ii[K ] trong phương trình tổng quát (2.14). Lấy đạo hàm i tại giá trị 0:   T i i 0 i 2M ( [T ] [T ]dxdy) V {F } t    được đưa vào véc tơ i{F} trong phương trình tổng quát (2.14). 2.2.3.6 Véc tơ tải trọng ứng với trọng lượng bản thân của khối Giả sử x y( )f ,f là trọng lượng bản thân tác dụng lên khối thứ i, khi đó thế năng của tải trọng bản thân x y( )f ,f sẽ là:     TT x w x y i i y f (f u f v)dxdy D T dxdy f                 (2.50) Lấy đạo hàm của w sẽ cho ta rf là một véc tơ 6x1:   T x y if S f S 0 0 0 0 {F} được đưa vào véc tơ tải trọng i{F} trong phương trình (2.14). 2.2.3.7 Ma trận con tạo bởi lực cản nhớt Lực cản nhớt tỷ lệ với vận tốc cũng như diện tích của khối. Khi chuyển vị thay đổi tính theo đơn vị thời gian, lực cản nhớt sẽ là: x y f u f vt               (2.54) ở đây t là bước thời gian; u và v là chuyển dịch tính trên một đơn vị thời gian. Thế năng do lực nhớt của khối phần tử thứ i sẽ là: T T v i i i i[D ] [T ] [T ][D ]dxdy t      (2.56) Lấy đạo hàm của v sẽ cho ta là ma trận 6x6 :       T i i ii 2 T T dxdy K     đưa vào ma trận ii[K ] trong (2.14). 2.2.3.8 Ma trận con do chuyển dịch cưỡng bức tại một điểm Giả sử một khối bị ngăn cản chuyển dịch theo hai phương x và y. Khi đó, chuyển dịch (u,v) tại điểm cố định (x,y) của khối sẽ bằng 9 0. Vấn đề này được thực hiện bằng cách sử dụng hai lò xo có độ cứng p rất lớn đặt theo hai phương x và y. Thế năng biến dạng đàn hồi của lò xo là m là:           T T2 2m i i i i up p p u v ( u v ) D T T D v2 2 2           (2.61) Lấy đạo hàm theo các thông số biến dạng và chuyển vị. Kết quả nhận được là ma trận 6x6: Ti i iip[T ] [T ] [K ] được đưa vào ma trận  iiK trong phương trình tổng quát (2.14). 2.3 Tiếp xúc và tƣơng tác giữa các khối 2.3.1 Vấn đề tiếp xúc Về mặt tổng quát có 3 dạng tiếp xúc cơ bản được mô tả trên hình 2.8 bao gồm: tiếp xúc đỉnh-cạnh, đỉnh-đỉnh, cạnh-cạnh. a)Tiếp xúc đỉnh-cạnh b)Tiếp xúc đỉnh-đỉnh c)Tiếp xúc cạnh-cạnh Hình 2.8 Ba dạng khác nhau của tiếp xúc Tiếp xúc cạnh-cạnh có thể chuyển thành tiếp xúc hai góc với cạnh. Để xử lý vấn đề tiếp xúc giữa các khối với nhau, DDA sử dụng một phương pháp được gọi là phương pháp “penalty”. Nguyên tắc đặt ra khi các khối tiếp xúc với nhau là không thể xảy ra trạng thái chồng lên nhau hoặc xuyên vào nhau. Vấn đề này được gọi là “cưỡng bức không xuyên”(inter-penetration). Trong phương pháp “penalty”, khi hai khối tiếp xúc nhau, “cưỡng bức không xuyên” được thực hiện bằng cách đặt vào một tham số “penalty” giống như một lò xo có độ P3 P4 P2P1 10 cứng p tại điểm tiếp xúc, lò xo này được đặt theo phương của đỉnh xâm nhập nhằm ngăn cản việc xuyên vào nhau của các khối. 2.3.2 Liên kết tại điểm tiếp xúc Hai khối được xem là ở trong trạng thái tiếp xúc khi và chỉ khi: ij 2   ( ij là khoảng cách giữa hai khối, là chuyển dịch lớn nhất của một trong hai khối trong bước thời gian trước đó) và không xảy ra việc chồng lên nhau khi đỉnh của khối này chuyển dịch tới cạnh khối kia mà không bị xoay. Ba trạng thái tiếp xúc giữa hai khối với nhau là: trạng thái “mở”, “đóng” và “khóa”. Khi điểm tiếp xúc ở trạng thái “mở” thì không có bất kỳ một lò xo nào được đặt vào tại điểm tiếp xúc. Khi ở trạng thái “đóng”, một lò xo cứng (hay còn gọi là một khóa) được đặt theo phương vuông góc với “đường tham chiếu”,còn ở trạng thái “khóa” thì có hai lò xo có độ cứng khác nhau, lần lượt đặt theo phương pháp tuyến và phương tiếp tuyến . Quá trình thêm vào hay bỏ đi các lò xo tiếp xúc (giá trị penalty) được xem là tiêu chuẩn “mở-đóng”. 2.3.3 Quy định về khóa và sự xuyên vào nhau Trạng thái tiếp xúc được xác định dựa vào tính toán khoảng cách vuông góc d giữa đỉnh và đường tham chiếu. Giả thiết rằng độ cứng của lò xo là p và khoảng cách xuyên là d, thế năng biến dạng đàn hồi của lò xo là: 2k (1/ 2).p.d  Lấy đạo hàm của k theo các tham số dri , dsi ta nhận được:     T 1 2 3 4 5 6 1 2 3 4 5 6 iip e e e e e e e e e e e e [K ] (2.73) được đưa vào ma trận ii[K ] trong phương trình tổng quát (2.14). Lấy đạo hàm của k theo các tham số dri , dsj ta nhận được :     T 1 2 3 4 5 6 1 2 3 4 5 6 ijp e e e e e e g g g g g g [K ] (2.75) được đưa vào ma trận ij[K ] trong phương trình tổng quát (2.14). Lấy đạo hàm của k theo các tham số drj , dsi là ma trận 6x6: 11     T 1 2 3 4 5 6 1 2 3 4 5 6 jip g g g g g g e e e e e e [K ] (2.77) được đưa vào ma trận ji[K ] trong phương trình tổng quát (2.14). Lấy đạo hàm của k theo các tham số drj , dsj là ma trận 6x6:     T 1 2 3 4 5 6 1 2 3 4 5 6 jjp g g g g g g g g g g g g [K ] (2.79) được đưa vào ma trận jj[K ] trong phương trình tổng quát (2.14). Lấy đạo hàm của k theo tham số dri tại giá trị 0 là véc tơ :     T0 1 2 3 4 5 6 i pS e e e e e e F l   (2.81) được đưa vào véc tơ tải trọng i{F} trong phương trình (2.14). Lấy đạo hàm của k theo tham số drj tại giá trị 0 là véc tơ :    0 1 2 3 4 5 6 j pS g g g g g g F l   (2.83) được đưa vào véc tơ tải trọng j {F}trong phương trình (2.14). ở đây :    r 2 3 1r 1 1 3 2 2r 1 1e [ y y t (x ,y ) x x t (x ,y )] / l       r 3 1 1r 2 2 1 3 2r 2 2g [ y y t (x ,y ) x x t (x ,y )] / l       1 2 1r 3 3 2 1 2r 3 3[ y y t (x ,y ) x x t (x ,y )] / l    2.3.4 Trƣợt giữa các khối Khi thành phần pháp tuyến của lực tiếp xúc nR là lực kéo, tức là: nR pd 0   .Trường hợp này tiếp xúc ở trạng thái “mở”, lúc này sẽ không có một lò xo penalty nào được đặt vào tại điểm tiếp xúc. Khi thành phần pháp tuyến của lực tiếp xúc nR là lực nén, hai khối tiếp xúc với nhau, tức là: nR pd 0   . Lúc này sẽ sử dụng tiêu chuẩn phá hoại Mohr-Coulomb để kiểm tra việc trượt giữa các khối. Giả sử ,c là góc ma sát trong và cường độ lực liên kết trên bề mặt tiếp xúc. Khi thành phần tiếp tuyến sR của lực tiếp xúc dọc theo đường tham chiếu có giá trị đủ lớn: s nR R tan c  .Trường hợp này, tiếp xúc ở dạng trượt; khi đó một lò xo theo phương pháp tuyến với đường tham chiếu được đặt vào để không cho các khối xuyên vào nhau nhưng vẫn 12 cho phép quá trình trượt diễn ra dọc theo đường tham chiếu. Khi thành phần tiếp tuyến sR của lực tiếp xúc dọc theo đường tham chiếu có giá trị: s nR R tan c  . Lúc này, tiếp xúc ở dạng “khóa” ; khi đó điểm tiếp xúc là cố định và bị khoá bởi hai lò xo theo phương pháp tuyến và tiếp tuyến để không cho phép quá trình trượt diễn ra. 2.4 Những ứng dụng của DDA Từ khi được đề xuất cho đến nay đã qua hơn hai thập kỷ, DDA đã chứng minh tính hiệu quả của mình trong việc dự đoán các nguy cơ mất ổn định cũng như giảm thiểu các thiệt hại trong trường hợp xảy ra sự phá hoại các khối đá. Các bài toán được thực hiện như : + Ổn định của mái đá nghiêng + Chuyển động do động đất + Sự xuất hiện và lan truyền khe nứt 2.5 Những hạn chế của DDA 1-Tính chính xác của phương pháp phụ thuộc đáng kể vào các thông số đầu vào. 2-Việc nghiên cứu trạng thái trượt các khối bằng cách sử dụng tiêu chuẩn Mohr-Coulomb nhưng hệ số ma sát vẫn xem xét là hằng số. 3-Hầu hết các chương trình của DDA giới hạn cho bài toán phẳng, trong khi các vấn đề đặt ra trong thực tế thường là bài toán ba chiều. CHƢƠNG III XÂY DỰNG THUẬT TOÁN VÀ CHƢƠNG TRÌNH TÍNH 3.1 Đặt bài toán 3.1.1 Đặt vấn đề Bài toán phân tích chuyển động của các khối rời rạc được gặp tương đối nhiều trong thực tế, đặc biệt là trong lĩnh vực địa chất núi đá. Bằng việc sử dụng phương pháp “penalty” như đề cập trong 13 chương 2, lý thuyết DDA giúp chúng ta mô phỏng được quá trình tương tác, chuyển động của các khối trong hệ thông qua việc tích phân phương trình chuyển động theo thời gian để xác định giá trị chuyển dịch của khối. 3.1.2 Mô hình tính toán Giới hạn xét là bài toán phẳng, việc đưa bài toán không gian của hệ các khối thực tế về bài toán phẳng bằng cách chọn vị trí mặt cắt phẳng cần nghiên cứu đi qua (hướng của mặt phẳng tùy thuộc vào từng bài toán cụ thể); giao tuyến của mặt cắt phẳng với các khối không gian cho hình ảnh đại diện các khối của cơ hệ trong bài toán phẳng. Mô hình trong các bài toán nghiên cứu được lấy theo mô hình được trình bày trong các tài liệu của giáo sư Shi Genhua [22],[23]. 3.2 Xây dựng thuật toán và sơ đồ khối 3.2.1 Giả thiết tính toán + Giới hạn phân tích là bài toán phẳng. + Trong quá trình chuyển động các khối không được đứt gãy. + Hình dạng và kích thước các khối được xấp xỉ bằng các đa giác có số đỉnh bất kỳ và vật liệu được giả thiết là đẳng hướng trong phạm vi từng khối. 3.2.2 Xây dựng thuật toán Như đã trình bày ở chương 2, quá trình tính toán được chia thành nhiều bước thời gian; các thành phần của ma trận độ cứng [K] và véc tơ tải trọng  F đều phải được xây dựng lại tương ứng mỗi bước thời gian. Để làm được điều này, trong mỗi bước tích phân (bước thời gian tính toán) phải xác định trạng thái của các khối trong cơ hệ. Khi đã xác định được toàn bộ các thành phần của ma trận [K] và véc tơ  F của mỗi khối trong hệ, tiến hành tích phân phương trình (2.14) ta có được véc tơ chuyển vị của mỗi khối  iD . Từ đó chuyển dịch 14 của mỗi khối trong bước thời gian đó sẽ được xác định theo công thức (2.11). Giá trị xuất ra của bước thời gian này lại là giá trị đầu vào cho bước thời gian kế tiếp. Trên cơ sở thuật toán nêu trên, sơ đồ giải bài toán DDA được tóm tắt như trên hình 3.3. Hình 3.3 Sơ đồ giải bài toán DDA Khối khởi tạo bài toán bao gồm các nội dung công việc : + Xây dựng số liệu hình học các khối; + Nhập giá trị đặc trưng vật liệu: E, ,  , C,  ; + Các loại tải trọng tác dụng lên các khối ; + Số liệu thời gian: bước thời gian t , tổng số bước thời gian ; Khối xây dựng, giải phương trình chuyển động bao gồm : + Xây dựng phương trình chuyển động cho tất cả các khối đơn có trong hệ. + Xây dựng phương trình do tương tác, tiếp xúc các khối: kiểm tra điều kiện tiếp xúc các khối tại thời điểm ban đầu và trong suốt quá trình tính toán. Khối xuất kết quả bao gồm : + Chuyển dịch, biến dạng khối (hoặc ứng suất các khối); + Trạng thái tiếp xúc giữa các khối 3.3 Các thông số đầu vào theo phân tích DDA 3.3.1 Tham số vật lý Để phân tích biến dạng của khối theo phương pháp DDA, thì các tính chất vật lý và cơ học của khối cần xác định như mô đun đàn hồi E, hệ số Poatxông  và trọng lượng thể tích  . Ngoài ra, để xác định Khëi t¹o bµi to¸n X©y dùng, gi¶i p.tr×nh chuyÓn ®éng XuÊt kÕt qu¶ 15 trạng thái trượt và tách giữa các khối theo tiêu chuẩn Mohr-Coulomb thì góc ma sát  và lực dính kết c phải được xét đến. 3.3.2 Tham số điều khiển + Bước thời gian t : Trong thực tế t có thể có giá trị từ 0,0001 đến 0,01 s. + Độ cứng lò xo liên kết (theo hướng pháp tuyến và tiếp tuyến) ( , )n sk k : Nghiên cứu gần đây chỉ ra rằng nên sử dụng độ cứng lò xo liên kết có giá trị trong khoảng giới hạn nhỏ thì việc phân tích chính xác hơn. + Hệ số chuyển vị lớn nhất: Giá trị của nó được đề nghị lấy trong khoảng từ 0,001 đến 0,01 để phân tích được hội tụ. + Tiêu chuẩn “mở”,“đóng”: Các tiêu chuẩn “mở”,“đóng” thường được sử dụng là 0f với giá trị của 0f là được đề xuất theo kinh nghiệm bằng 1e-7. + Hệ số kháng nhớt 01k : Trong phân tích bài toán tĩnh, khi đó 01k = 0. Đối với các bài toán động, khuyến cáo nên lấy 01k = 0,8 [18]. 3.4 Giới thiệu chƣơng trình tính DDA.m 3.4.1 Giới thiệu chƣơng trình Trên cơ sở thuật toán nêu trên tác giả đã lập chương trình mang tên DDA.m được viết bằng ngôn ngữ lập trình MATLAB [1],[2],[10]. Số liệu đầu vào và kết quả số của chương trình DDA được lưu ở dạng file văn bản, file đồ họa. 3.4.2 Khả năng tính toán của chƣơng trình 1- Tính toán các đặc trưng hình học của đa giác có số đỉnh bất kỳ; 2- Tính toán trường chuyển vị, biến dạng, ứng suất trong môi trường biến dạng không liên tục; 3- Mô hình hóa các dạng tiếp xúc đỉnh-đỉnh, đỉnh-cạnh, cạnh-cạnh; tương tác giữa các khối khi tiếp xúc. 16 3.5 Một số thử nghiệm số 3.5.1 Bài toán chuyển động tự do của hệ khối Hệ gồm 3 khối, có vị trí ban đầu được xác định bởi tọa độ của chúng như hình 3.5. Trong đó: khối 1 được cố định, khối 2 chuyển động tự do, khối 3 chuyển động tự do nhưng chịu tác dụng của hai lực F1=1 kN và F2=3kN như hình vẽ. Hình 3.5 Sơ đồ bài toán Hình 3.6 Mô hình trong DDA Các kết quả nhận được từ chương trình DDA.m:  Đặc trưng hình học các khối Giá trị các đặc trưng hình học như: diện tích S, mô men tĩnh Sx, Sy, mô men quán tính Sxx, Syy và mô men quán tính ly tâm Sxy hoàn toàn phù hợp với kết quả tính bằng công thức giải tích.  Chuyển dịch các khối theo thời gian Hình ảnh chuyển dịch các khối theo thời gian được mô tả như trên hình 3.7 và hình 3.10. Trong đó, kết quả chuyển động khối 2 được so sánh với lời giải giải tích (bài toán vật rắn rơi tự do) cho sai số chấp nhận được. Nhận xét: Về cơ bản kết quả tính theo DDA là phù hợp với lý xo y A ED C B FG 2 1 2 3 4 5 6 7 1 2 3 4 5 6 7 IH JK 8 8 F1 F2 3 1 2 F1 F2 3 1 Fix 1(2.5,0.5) Fix 3(6.5,2.0) Fix 2(6.5,0.5) 17 thuyết, do đó chương trình DDA.m là đáng tin cậy. Giá trị tính toán chuyển dịch theo DDA nhỏ hơn so với lời giải giải tích là do hàm xấp xỉ chuyển vị chỉ là bậc 1 và bước thời gian được chọn chưa đủ nhỏ; sai số trên về sử dụng là chấp nhận được. Chuyển động của vật thể tự do chịu tác dụng của ngoại lực (khối 3) phù hợp với quy luật chuyển động. Hình 3.7 Vị trí các khối thời điểm ban đầu t = 0,000s Hình 3.10 Vị trí các khối thời điểm t = 1,050s 3.5.2 Bài toán chuyển động của mái dốc đá Giả sử có các hòn đá mồ côi được sắp xếp ổn định trên mái dốc đá, giữa các hòn đá có các chất lấp nhét. Vì lý do nào đó các chất lấp nhét giữa các hòn đá bị rửa trôi, quá trình mất ổn định diễn ra. Vấn đề đặt ra là quá trình chuyển dịch của toàn bộ các khối đá theo thời gian sẽ diễn ra như thế nào? Việc phân tích được thực hiện sau 100 bước tính toán (t=0,500 s). Kết quả tính toán chuyển dịch được mô phỏng trên hình 3.15 và hình 3.20. Chương trình cũng xác định lực tương tác khi các khối tiếp xúc, va chạm nhau. Nhận xét: Việc mô phỏng quá trình trượt của các khối đá trên mái dốc theo thời gian có ý nghĩa rất lớn trong thực tế. Kết quả của 18 mô phỏng số cho phép kiểm tra quá trình ổn định mái dốc đá khi có tác động của yếu tố tự nhiên hoặc nhân tạo. Bên cạnh đó, nó còn cho phép xác định phạm vi ảnh hưởng cũng như những tác động khi quá trình mất ổn định diễn ra. Hình 3.15 Mô hình bài toán Hình 3.20 Mái dốc tại t= 0,50s CHƢƠNG IV SỰ ỔN ĐỊNH KHOANG HẦM TRONG MÔI TRƢỜNG ĐÁ NỨT NẺ 4.1 Đặt bài toán Giả sử khoang hầm được tạo ra trong môi trường gồm các khối đá rời rạc, không có nước ngầm; giữa các khe nứt có chất lấp nhét liên kết các khối đá với nhau. Tùy thuộc vào khoảng trống được tạo ra bởi các khe nứt nẻ và kích thước khoang hầm mà chuyển dịch các khối đá theo thời gian cũng như tương tác giữa các khối trong các trường hợp khác nhau sẽ có giá trị khác nhau. 4.2 Mô hình nghiên cứu Với ý nghĩa thử nghiệm ứng dụng của DDA.m vào tính toán ổn định khoang hầm, mô hình các khối đá để làm thử nghiệm số trong chương này hoàn toàn là nhân tạo, các tính chất cơ lý của khối đá được lấy gần giống đá trầm tích trong tự nhiên. Tuy nhiên, do hạn M¸i dèc ®¸ LÊp nhÐt Khèi ®¸ 19 chế của chương trình DDA.m nên các khối đá được giả thiết có góc ma sát và lực dính C = 0. 4.3 Giới hạn miền khảo sát Trong mọi chương trình tính đều phải xác định kích thước của miền khảo sát. Trên biên của miền này cần xác định tác động của miền ngoài đối với kết quả tính toán thông qua ảnh hưởng của tải trọng hay biến dạng. Trong chương này giớí hạn miền khảo sát được lấy theo khuyến cáo của Hội Địa Kỹ thuật Đức [7]. 4.4 Bài toán khoang hầm trong môi trƣờng đá phân lớp Mô hình phân tích có giới hạn kích thước: chiều cao là 6 m và chiều rộng là 12m. Khoang hầm có dạng hình tròn với đường kính D=1,5m. Hệ khe nứt được tạo ra bao gồm một hệ các phân lớp có góc nghiêng so với phương ngang là . Các phân lớp ngang được giả định được mở rộng (có chiều dài) vô hạn, với khoảng cách trung bình là h. Độ mở rộng của hệ khe nứt được ký hiệu là . Để tính toán thuận lợi ta tiệm cận hình dáng khoang hầm về hình lục giác.Trong thử nghiệm số này, chúng ta nghiên cứu chuyển dịch tại hai vị trí biên ở bên hông A và điểm nóc B của khoang hầm trong sự phụ thuộc vào: 1. Khoảng cách giữa cáckhe nứt; 2. Chiều dày phân lớp; 3. Góc nghiêng các phân lớp so với mặt phẳng ngang. 4.4.1 Trƣờng hợp góc nghiêng của các phân lớp thay đổi Thử nghiệm nghiên cứu chuyển dịch tại hai điểm A, B trên biên khoang hầm khi thay đổi giá trị góc nghiêng các phân lớp . Chiều dày các phân lớp đều nhau và có giá trị h = 0,8m. Quá trình tính toán được tiến hành sau 500 bước tính (tổng thời gian t=2,25s). Kết quả thể hiện trên biểu đồ 4.13. Phân tích kết quả cho thấy rằng: 1. Khi góc nghiêng của các phân lớp tăng lên, chuyển dịch bên 20 hông khoang hầm có xu hướng tăng theo, trong khi đó chuyển dịch tại nóc khoang hầm có xu hướng giảm dần. 2. Giá trị các chuyển dịch đều nhỏ hơn 5cm nên theo phân loại ổn định của VNIMI (bảng 1.2, chương 1) thì các vị trí này đều được xem là ổn định. Hình 4.13 Biểu đồ quan hệ giữa chuyển dịch biên khoang hầm “U” với góc nghiêng “ ” Hình 4.25 Biểu đồ quan hệ giữa chuyển dịch biên khoang hầm “U” với khoảng cách “h” 4.4.2 Trƣờng hợp chiều dày các phân lớp thay đổi Trong thử nghiệm số này góc nghiêng các phân lớp được lấy = 45 o, độ mở rộng giữa các phân lớp =0. Quá trình tính toán được thử nghiệm sau 500 bước tính (tại thời điểm t=2,25s). Kết quả thể hiện trên biểu đồ 4.25. Phân tích chuyển dịch trên biên cho thấy rằng: 1. Khi chiều dày giữa các phân lớp tăng lên thì chuyển dịch tại điểm bên hông A và điểm nóc B đều giảm, giá trị dịch chuyển tại điểm nóc và điểm hông đều giảm xấp xỉ 5 lần khi giá trị chiều dày phân lớp tăng 2 lần( từ 0,8m đến 1,6m). 2. Khi chiều dày các phân lớp thay đổi từ 0,8m đến 1,4m thì sự biến thiên chuyển dịch rất lớn (xấp xỉ 4 lần) trong khi chiều dày phân lớp thay đổi từ 1,4m đến 1,6m thì sự thay đổi này không nhiều (xấp xỉ 1,3 lần). 4.4.3 Trƣờng hợp độ mở các phân lớp thay đổi Trong thử nghiệm số này góc nghiêng các phân lớp được lấy = 0 2 4 6 15 30 45 60 75 U, cm , độ UA UB 0 2 4 6 0.8 1 1.2 1.4 1.6 U,cm h,m UA UB 21 45 o, chiều dày của các phân lớp h=0,8m. Kết quả chuyển dịch tại A, B được thể hiện trên biểu đồ 4.37. Hình 4.29 Sơ đồ bài toán với độ mở rộng  =0,2cm Hình 4.30 Chuyển dịch các khối tại t= 2,25s(dA=5,351, dB=7,028) Hình 4.37 Biểu đồ quan hệ giữa chuyển dịch biên khoang hầm “U” với độ mở rộng khe nứt “ ” Kết quả chuyển dịch trên biên cho thấy rằng: 1. Khi khoảng cách giữa các phân lớp tăng lên thì chuyển dịch tại các điểm trên biên khoang hầm tăng lên. 2. Giá trị chuyển dịch tại điểm nóc lớn hơn chuyển dịch tại điểm hông.Tỷ số giữa chuyển dịch tại điểm nóc và chuyển dịch tại điểm hông được xem không có sự biến đổi nhiều. 3. Đại lượng độ mở rộng “” được xem có ảnh hưởng mạnh mẽ nhất tới chuyển dịch của các điểm trên biên khoang hầm. 4.5 Khoang hầm vòm tƣờng thẳng trong môi trƣờng đá nứt nẻ Khảo sát dịch chuyển tại điểm nóc khoang hầm dạng vòm tường thẳng (chuyển dịch điểm A thuộc khối 22) trong môi trường đá nứt nẻ theo thời gian như hình 4.38 . Mô hình phân tích có kích thước giới hạn với chiều rộng là 14m và chiều cao là 5,5 m. Các khối đá nứt nẻ phân bố đối xứng, phía trên nóc khoang hầm có một nêm đá. 0 10 20 30 0.1 0.2 0.3 0.4 0.5 U, cm , cm UA UB 22 Hình 4.38 Sơ đồ công trình và các khối đá trong mô hình Hình 4.39 Hình ảnh chuyển dịch tại t=1,125s (dA =4,57 cm) Chuyển dịch của khoang hầm ở đây chủ yếu là chuyển dịch của nêm đá phía trên nóc khoang hầm. Do phân bố các khối mang tính đối xứng nên kết quả chuyển dịch của nêm đá trên nóc chỉ có chuyển dịch theo phương thẳng đứng xuống dưới. 4.6 Công trình ngầm trong môi trƣờng đá nứt nẻ 4.6.1 Công trình ngầm trong môi trƣờng đá nứt nẻ đối xứng Xét công trình ngầm nằm trong môi trường đá phân lớp, nứt nẻ trong đó các khối đá phía trên được đánh số thứ tự từ 5-24 có dạng như ở hình 4.41. Do hình dáng công trình đối xứng, được đặt trong môi trường đối xứng nên chuyển dịch của các khối đá và của công trình là đối xứng. Hình 4.41 Sơ đồ công trình và các khối đá trong mô hình Hình 4.43 Các khối tại thời điểm t = 0,144 s 23 Kết quả nhận được từ chương trình DDA.m cho biết thời gian các khối va chạm vào công trình, sự biến đổi ứng suất, biến dạng trong các khối đá và trong công trình trong suốt thời gian nghiên cứu. 4.6.2 Công trình ngầm trong môi trƣờng đá nứt nẻ không đối xứng Xét công trình ngầm nằm trong môi trường đá phân lớp, nứt nẻ trong đó các khối đá phía trên có dạng như ở hình 4.45. Hình 4.45 Sơ đồ công trình và các khối đá trong mô hình Hình 4.47 Các khối tại thời điểm t = 0,144 s Để nghiên cứu tương tác giữa môi trường đá và kết cấu công trình ngầm chúng ta xác định chuyển dịch tại điểm A (xem hình 4.41 và 4.45) trong hai trường hợp tại mục 4.6.1 và 4.6.2. Bằng cách biểu thị tọa độ theo phương ngang của điểm A theo thời gian, chúng ta thấy rằng trong trường hợp địa tầng phân bố đối xứng thì chuyển dịch điểm A là rất ít; trong trường hợp địa tầng phân bố không đối xứng thì chuyển dịch điểm A thể hiện khá rõ rệt: khi bị các khối đá bên phải tác dụng thì điểm A dịch sang trái, sau khi chịu tác dụng của kháng lực do các khối đá bên trái thì điểm A lại dịch sang phải. KẾT LUẬN CHUNG I. Các kết quả chính của luận án 1. Đã tiếp cận và nghiên cứu những nội dung cơ bản của lý thuyết phân tích biến dạng không liên tục. 24 2. Trên cơ sở lý thuyết phân tích biên dạng không liên tục, đã thiết lập chương trình tính toán với các khả năng: 2.1-Tính toán trường chuyển vị, biến dạng, ứng suất trong môi trường biến dạng không liên tục đã xét đến hầu hết các lực: trọng lượng bản thân, ngoại lực tác dụng dưới dạng lực tập trung, lực phân bố, lực quán tính, lực cản, lực ma sát, lực dính (bài báo [2]); 2.2- Mô hình hóa các dạng tiếp xúc đỉnh-đỉnh, đỉnh-cạnh, cạnh-cạnh; tương tác giữa các khối khi tiếp xúc. Từ đó mô phỏng chuyển dịch, quá trình tiếp xúc, tương tác giữa các khối trong cơ hệ (bài báo [3], [4], [5]). 3. Trên cơ sở chương trình tính toán DDA.m được thiết lập, bước đầu đã nghiên cứu ổn định của khoang hầm trong khối đá nứt nẻ với các kết quả như sau: 3.1- Đã tiến hành thử nghiệm nghiên cứu chuyển dịch trên biên khoang hầm tròn trong môi trường đá phân lớp, nứt nẻ trong các trường hợp thay đổi góc nghiêng, độ dày, độ mở rộng giữa các phân lớp và so sánh kết quả đạt được với lời giải bằng giải tích; 3.2- Nghiên cứu tương tác giữa khối đá nứt nẻ với kết cấu công trình. II. Hƣớng phát triển của luận án Với mục tiêu chính của nghiên cứu này, sự ổn định của khoang hầm trong khối đá nứt nẻ bằng phương pháp Phân tích biến dạng không liên tục, hướng nghiên cứu tiếp theo cần tập trung vào các nội dung sau: + Các tính chất cơ lý khác nhau của đá. + Ảnh hưởng của góc ma sát, cường độ lực dính theo các phân lớp và khe nứt. + Biện pháp gia cố để giảm chuyển dịch trên biên khoang hầm bằng phun vữa xi măng và hệ thống neo.