Nghiên cứu tương tác lưu chất và vật thể đàn hồi bằng phương pháp IBM (Immersed Boundary Method)

Nội dung text: Nghiên cứu tương tác lưu chất và vật thể đàn hồi bằng phương pháp IBM (Immersed Boundary Method)

1. NGHIÊN CỨU TƯƠNG TÁC LƯU CHẤT VÀ VẬT THỂ ĐÀN HỒI BẰNG PHƯƠNG PHÁP IBM (IMMERSED BOUNDARY METHOD) FLUID AND ELASTIC STRUCTURE INTERACTION RESEARCH BY USING IBM (IMMERSED BOUNDARY METHOD) Lê Thanh Tùng, (1) PGS. TS. Nguyễn Hoài Sơn, (2) TS. Phan Đức Huynh. (1) (2) Khoa Xây dựng và Cơ học ứng dụng, trường Đại học Sư Phạm Kỹ Thuật, tp. HCM. TÓM TẮT Để tính toán cho các bài toán về tác động giữa lưu chất và vật thể, lưu chất và lưu chất trong thực tế hiện nay đã có rất nhiều phương pháp được sử dụng như FEM, FVM, FDM . Nhìn chung các phương pháp này đều sử dụng việc chia lưới, việc chia lưới lại phải phụ thuộc vào hình dạng của vật thể. Đối với những bài toán mà hình dạng vật thể phức tạp thì việc chia lưới theo những cách thông thường sẽ càng trở nên khó khăn và chi phí tính toán lớn. Với phương pháp IBM (Immersed Boundary Method) thì việc chia lưới rất đơn giản bằng cách chia lưới lưu chất trực tiếp trên các ô vuông của lưới Đề các, các ô lưới này cố định trong suốt quá trình tính toán và không phụ thuộc hình dạng của vật thể. Từ khóa: IBM, Immersed Boundary Method, tương tác lưu chất và lưu chất, tương tác lưu chất và vật thể. Abstract Nowadays, there are many methods be applied to the problems of fluid-structure and fluid-fluid interaction such as FEM, FVM, FDM In general, most of these methods are used for meshing that still depends on object shapes. It will make our spend a large cost calculations and difficulty if object has complex shape. However, with IBM (Immersed Boundary Method), it is more simple for the mesh that be on directly on Cartesian grid system that is fixed during calculation process and to be independent on the object shapes. Key words: IBM, Immersed Boundary Method, fluid- fluid interaction, fluid-structure interaction 1. Giới thiệu Phương pháp IBM là phương pháp áp dây nhỏ mềm trong một màng xà phòng đang dụng cho việc mô phỏng số của một dòng chảy có thể biểu hiện ba trạng thái động lực chảy nhớt, không nén được qua một vật thể học; trạng thái duỗi thẳng, trạng thái vỗ, và đàn hồi không nén được trong không gian hai trạng thái ổn định kép. Khi đó mô hình thí chiều. Ở phương pháp này, việc chia lưới nghiệm này được coi như là một lá cờ một không phụ thuộc vào hình dạng của vật thể. chiều trong một dòng chảy hai chiều, phát Các phương pháp sai phân hữu hạn và một hệ hiện của họ trái ngược với quan niệm phổ thống lưới Đề các so le được dùng để giải biến trong nhiều thập kỉ qua rằng lá cờ luôn quyết sự chuyển động của chất lỏng không vỗ trong gió. Trong bài báo này, bài toán về nén được trong miền Euler. Mặt khác, một sợi dây mềm trong màng xà phòng sẽ được lưới Lagrange di động dùng để rời rạc hóa mô phỏng bằng phương pháp IBM như là một miền kết cấu. Các biến Euler và Lagrange liên phiên bản hai chiều của bài toán lá cờ bay hệ với nhau qua các phương trình tương tác, trong gió. trong đó việc áp dụng một xấp xỉ trơn của hàm delta Dirac đóng một vai trò chủ đạo. Zhang và các đồng nghiệp [1] đã thực hiện một thí nghiệm tại viện Courant đại học New York, USA, Zhang đã cho thấy một sợi 1
2. U V 0 Llc Luc L Duỗi thẳng Lưỡng ổn Vỗ I I Hình 2: Trạng thái vỗ một sợi dây F S S Trường hợp đặt hai sợi dây song song với nhau, trong thí nghiệm của Zhang [1] đã chỉ O O ra rằng khi khoảng cách d giữa hai sợi dây đủ nhỏ (d/L 0,21±0,04), thì hai sợi dây sẽ vỗ lệch pha U: thùng chứa dung dịch xà phòng, V: van và vòi phun, nhau (xem hình 3). S: hai màng nylon, F: sợi dây, D: thùng hứng dung dịch xà phòng. d/L=0,21±0,04 Zhang và các đồng nghiệp [1] đã thực Cùng pha Ngược pha nghiệm nghiên cứu động lực học của các sợi dây mềm trong màng xà phòng đang chảy Hình 3: Trạng thái vỗ hai sợi dây như một mô hình thí nghiệm thích hợp đối với một lá cờ kích thước một chiều trong một cơn Trong bài báo này sẽ sử dụng phương gió chảy tầng hai chiều. Hình 1 cho thấy các pháp IBM để mô phỏng lại trạng thái vỗ ở cả thiết lập của thí nghiệm. Một màng xà phòng hai trường hợp một sợi dây và hai sợi dây. thẳng đứng được tạo thành, dòng chảy đi xuống giữa hai lớp nylon song song và nhanh 2. Phương pháp tính toán chóng đạt vận tốc giới hạn dưới tác động của Phương pháp IBM vừa là một công thức lực hấp dẫn và lực cản không khí. Một sợi toán học vừa là một phương pháp số. Công dây nhỏ mềm được đưa vào giữa dòng chảy thức toán học của IBM sử dụng được pha trộn của kênh màng xà phòng với đầu trên của nó giữa các biến Euler và Lagrange. Những biến được cố định bằng một ống mỏng vuông góc này có liên hệ với nhau qua các phương trình với màng xà phòng. tương tác trong đó hàm delta (δ) Dirac đóng Thí nghiệm của họ đã cho thấy rằng khi một vai trò chủ đạo. Trong phương pháp số chiều dài của sợi dây ngắn hơn chiều dài dưới được hình thành từ công thức IB, các biến tới hạn Llc thì sợi dây không thể duy trì Euler được định nghĩa trên một lưới Đề các chuyển động vỗ của nó và luôn tiếp cận trạng cố định, và các biến Lagrange được định thái thẳng của nó theo hướng của dòng chảy, nghĩa trên một lưới cong di chuyển tự do qua không phụ thuộc vào độ lớn các yếu tố gây lưới Đề các cố định mà không buộc phải thích nhiễu ở bên ngoài áp đặt vào. Ngược lại, khi ứng với nó theo bất kỳ một cách nào cả. Mô chiều dài sợi dây trên chiều dài trên tới hạn tả Euler dựa trên các phương trình Navier – Luc , trạng thái duỗi thẳng biến mất và sợi dây Stokes được sử dụng cho động lực chất lỏng, luôn vỗ trong dòng chảy màng xà phòng. Họ và một mô tả Lagrange được sử dụng cho cũng phát hiện ra rằng, khi chiều dài ở giữa từng đối tượng nhúng chìm trong chất lỏng. Biên được nhúng là những đối tượng có dạng Llc và Luc , sợi dây sẽ ở một trong hai trạng như lá mỏng van tim, cánh dù, các dạng cánh thái duỗi thẳng hoặc vỗ tùy thuộc vào điều mỏng (bao gồm chim, côn trùng hoặc cánh kiện ban đầu hoặc các yếu tố gây nhiễu bên dơi), cánh diều, lá cờ và chong chóng gió ngoài (Hình 2). Ta coi thể tích bị chiếm bởi các bề mặt như 2
3. thế bằng 0 [3][4]. Phương trình chuyển động trình động lượng), phương trình (2) là phương của hệ như sau: trình liên tục. Phương trình (3) và (4) là các u phương trình tương tác giữa biên nhúng và    uug  p ufu (1) t chất lỏng, chúng chuyển từ các biến Lagrange  u 0 (2) sang Euler và ngược lại phương trình (5) là phương trình chuyển động của biên, chuyển fx,FxX tststds , δ, (3)  đổi các biến Euler sang Lagrange. Tất cả các x, t Mstdsδ,xX (4) 0 phương trình tương tác liên quan đến các toán  tử tích phân cùng với một nhân δ,xX st , X st, U st, u x, t δ, x X s t d x (5) t  tuy nhiên có sự khác biệt giữa (3), (4) và (5). Tτˆ E Các phương trình tương tác xác định ρ và f có FFF ststst,,, u (6) kus X mối liên hệ giữa các mật độ tương ứng, giá trị  X T  1 (7) bằng số của ρ và M không giống nhau tại các k s điểm tương ứng, và tương tự giá trị bằng số X τˆ s (8) của f và F không giống nhau tại các điểm X s tương ứng. Nhưng giá trị bằng số của X và 2 t 1 2X st, Edsuu  (9) ux ,t lại giống nhau tại các điểm tương ứng. 2 s2 Mặc dù x,t và fx,t không xác định Ở đây các biến Euler độc lập là toạ độ trên biên nhúng, nhưng tích phân của chúng Đề các x = ( )x , y và thời gian t, và các biến lại xác định trên bất kì một thể tích xác định Lagrange độc lập là toạ độ cong của vật thể s nào đó. Đặc biệt, tích phân của x, t trên và thời gian t; trong đó Ω là miền chất lỏng, Г một thể tích là tổng khối lượng của chất lỏng là miền biên nhúng. Các biến Euler phụ thuộc được chứa trong thể tích đó và khối lượng của là vận tốc chất lỏng u x,,,,t u x t v x t , phần biên đàn hồi được chứa trong thể tích áp suất pt x, , mật độ x, t , và mật độ lực đó. Tương tự, tích phân của fx , t trên một thể tích là lực toàn phần tác động vào chất Euler fxxx ,,,,tftft xy . Các biến lỏng bởi phần biên nhúng được chứa trong thể Lagrange phụ thuộc là vị trí của biên nhúng tích đó. Mật độ lực đàn hồi Lagrange F bằng X s,,,, t X s t Y s t , mật độ lực Lagrange 0 tại các điểm chất lỏng, và mật độ lực đàn F stFstFst,,,, , và vận tốc biên xy hồi Euler fx , t bằng 0 tại các điểm không . Và ở đây μ là độ U stUstVst,,,, gian x ngẫu nhiên nằm trong vùng chất lỏng nhớt chất lỏng, 0 là mật độ khối lượng chất tại thời điểm t. Tác động của biên nhúng lên lỏng, g là gia tốc trọng lượng và g = (0, −gy), chất lỏng xung quanh chủ yếu có được bằng  là hệ số cản không khí và  00g / v [2], cách truyền ứng suất của biên nhúng tới chất M là mật độ khối lượng Lagrange đồng nhất lỏng. Bởi vì vị trí của biên nhúng không của biên (chú ý M không phụ thuộc vào thời thường trùng với các điểm nút của lưới Đề gian, vì khối lượng được bảo toàn). các nên lực được phân bố qua một dãi các ô xung quanh mỗi điểm Lagrange (xem hình 4) Phương trình (1) là phương trình Navier– thông qua một xấp xỉ trơn của hàm delta Stokes không nén được (còn gọi là phương Dirac. 3
4. Đặt t là khoảng thời gian của mỗi bước, n là kí hiệu cho bước thời gian, chẳng hạn u un x ,nt . Như vậy nn TD ks X 1 (12) D Xn τˆ n s (13) n Ds X Ta sẽ có lực kéo – nén tại thời điểm n Hình 4: Chuyển đổi lực Fk từ điểm biên như sau Lagrange ()Xk tới các nút chất lỏng xung n n n quanh. Vùng tô đậm biểu thị mức độ phân Fks DT τˆ (14) bố lực. Và năng lượng uốn và lực uốn được rời Mật độ lực Lagrange F(s,t) gồm hai rạc theo thành phần: lực kéo – nén Fk(s,t) và lực uốn 2 Fu(s,t). T là ứng suất trong sợi dây, được tính 1 Nb 1 DDX ss (15) bằng định luật Hook (xem phương trình 7). τˆ EstuuX(,)    2 2 k 2 s là véc tơ đơn vị tiếp tuyến được xác định tại mỗi điểm của sợi dây.  là hệ số đàn hồi, llN 3212 CC k F  lbl 1 uu l 2 (16) được chọn trong tính toán sao cho sợi dây s lN bl2 C 1 không giãn nở. u là hệ số uốn được đo s2 lường trong phòng thí nghiệm. Mật độ lực uốn thu được bằng cách lấy đạo hàm Fréchet Trong đó năng lượng uốn Eu (xem phương trình 9). XXX 2 C kkk 11, k 2 (17) 3. Phương pháp số s kN 2, ,1 Phương trình Navier – Stokes nhớt không b nén được (phương trình 1) được giải để xác ở đây có các điều kiện, giống như l 3 , định vận tốc chất lỏng. Sự tương tác giữa chất điều kiện này sẽ bằng 1 khi nó đúng, và bằng lỏng và sợi dây được xử lý bằng một xấp xỉ 0 khi sai. Nb là số tổng số điểm rời rạc của trơn với hàm delta Dirac. Có nhiều phương biên nhúng. pháp số để giải phương trình Navier – Stokes không nén được, trong bài báo này phương Và ba phương trình (3), (4) và (5) được pháp sai phân hữu hạn [5] được sử dụng để rời rạc như sau giải phương trình này. Các đạo hàm bậc nhất fxFx-nnn X ()sssδ (18) theo không gian trong các phương trình (1), s  (2), (6) ,(7), (8) và (9) được tính theo một nn()xx X Mssδ (19) hàm bất kì như sau 0  s     , (10) nnn 112 D  Uuxx ss X h δ (20)  x  các đạo hàm bậc hai theo không gian thì được XXnn 1(s)( ) s Un 1(s) (21) tính theo t   2    DD  Trong đó h là bề rộng ô lưới Đề các.  2 (11) 4
5. δ là một xấp xỉ trơn hai chiều của hàm 4. Kết quả delta Dirac. Trong bài báo này hàm delta Các tham số mô phỏng trong phần kết Dirac được chọn theo [3] quả này được thể hiện trong bảng 1 và kết quả 1 xy được so sánh với kết quả trong thí nghiệm của δ x 2  (22) hhh Zhang [1] và kết quả mô phỏng bằng FEM của T. Sawada và T. Hisada [2]. ở đây x, y là thành phần Đề các của x và  được chọn  rr 0, 2 1 527124, 21rrrr 2 8 1 32144, 10rrrr 2 8 1 32144, rrrr 01 2 8 1 527124, rrrr1 2 2 8 0, 2 r Bảng 6.1 Thông số thí nghiệm của Zhang [1], mô phỏng bằng FEM [2] và mô phỏng hiện tại của một sợi dây mềm Giá trị thí Giá trị mô phỏng Giá trị mô nghiệm [1] [2] phỏng hiện tại Vận tốc dòng vào (v0) 280 cm/s 280 cm/s 280 cm/s -4 2 3 -4 2 Mật độ khối lượng xà phòng ( 0 ) 3×10 g/cm 1 g/cm 3×10 g/cm Độ nhớt động xà phòng (μ) 1,2×10-5 g/cm.s 0,12 Pa.s 1,2×10-3 g/cm.s Bề rộng kênh màng xà phòng 8,5 cm 8,5 cm 8,5 cm Chiều dài kênh màng xà phòng 120 cm 25 cm 25,5 cm Mật độ khối lượng của sợi dây (M) 2×10-4 g/cm 1,132 g/cm3 3×10-4 g/cm Đường kính sợi dây 0,015 cm - - Diện tích mặt cắt ngang sợi dây 1,767×10-4 cm2 1,767×10-4 cm2 - 4 Độ cứng uốn sợi dây(u ) 0,1 erg.cm 0,01Pa cm 0,1 erg.cm Chiều dài sợi dây(L) 3 cm 3 cm 3 cm 2 2 2 Gia tốc trọng trường (gy) 980 cm/s 980 cm/s 980 cm/s Hệ số cản không khí (λ) - 00g / v 5
6. Trong phần mô phỏng này thì chiều dài Hình 5 thể hiện các trạng thái vỗ ở các sợi dây được chọn là 3 cm và được rời rạc ra thời điểm khác nhau, các đường xoáy bắt đầu làm 20 điểm, kích thước lưới của màng xà từ đầu tự do của sợi dây và nối tiếp ra phía phòng được thực hiện là 200×600, bề rộng sau của sợi dây xuôi chiều với dòng chảy của kênh là 8,5 cm, chiều dài kênh là 25,5 cm, với chất lỏng, sau kết thúc của mỗi đường xoáy -4 thì đường xoáy tiếp theo sẽ có chiều xoáy bước thời gian ∆t=10 . Vận tốc ban đầu v0 = ngược lại so với đường xoáy trước đó. Hình 6 280 cm/s, và số Reynolds là Re = 210. Hệ số thể hiện trường áp suất trong từng thời điểm nhớt μ được chọn là 1,2x10-3 g/cm lớn hơn hệ tương ứng với các hình đường xoáy ở hình 5. số nhớt trong thí nghiệm, việc chọn hệ số Hình 7 thể hiện tọa độ dao động của đầu tự do nhớt như vậy để tránh việc số Reynolds sợi dây, ban đầu sợi dây có phương lệch một không quá cao, đảm bảo rằng phương trình góc 100 so với phương thẳng đứng về phía Navier-Stokes vẫn đúng trong quá trình tính bên phải. Mới đầu sợi dây dao động với biên toán. Mật độ khối lượng được chọn lớn hơn độ tương đối nhỏ, có xu hướng tăng dần lên, so với mật độ khối lượng sợi dây trong thí sau khoảng một thời gian ngắn (khoảng 0,13 nghiệm, vì thực tế khi sợi dây ngâm trong s) sợi dây bắt đầu đi vào trạng thái vỗ ổn màng xà phòng thì khối lượng sợi dây tăng định với biên độ tương đối đồng đều, không thêm một phần do nước xà phòng thẩm thấu có xu hướng tăng thêm nữa. Và tần số dao vào trong sợi dây cho tới khi bảo hòa. Với các động của sợi dây là 50 Hz bằng với tần số thông số đã chọn trong bảng 1 dùng để mô trong thí nghiệm của Zhang [1], và trong mô phỏng lại hiện tượng vỗ của sợi dây trong phỏng bằng FEM của T. Sawada và T. Hisada màng xà phòng (được coi như là mô hình của [2] cũng bằng là 50 Hz. Biên độ dao động đo lá cờ có kích thước 1 chiều vỗ trong gió), kết được ở đây là 1,4 cm, trong thí nghiệm là 1,5 quả mô phỏng bằng IBM đã cho thấy kết quả cm, kết quả mô phỏng nhỏ hơn 6,6% so với tốt, tương đồng với kết quả trong thí nghiệm thí nghiệm. và kết quả trong mô phỏng bằng FEM. - Hình 5 Các đường xoáy của dòng chảy 6
7. Hình 6 Trường áp suất Vị trí đầu tự do Vị trí đầu Hình 7 Tọa độ dao động đầu tự do của sợi dây 7
8. Vị trí đầu tự do Vị trí đầu Hình 8 Tọa độ dao động đầu tự do của sợi dây Ở hình 8, ban đầu sợi dây trùng với Sử dụng các thông số trong bảng 1, thực phương thẳng đứng, ta thấy rằng ở thời điểm t hiện mô phỏng với các chiều dài dây khác = 0,08s sợi dây mới bắt đầu có sự dao động nhau. Ta thu được các kết quả tần số và biên và biên độ vỗ tăng dần tới thời điểm t = 0,32s độ vỗ của sợi dây (hình 9 (c)), ở chiều dài dây thì bắt đầu đi vào trạng thái dao động ổn định. 1,5 cm và 1,8 cm sợi dây không duy trì được So với trường hợp ban đầu sợi dây nghiêng trạng thái vỗ mà luôn tiếp cận trạng thái duỗi một góc 100 ta thấy rằng khi ban đầu sợi dây thẳng. ở vị trí khác nhau thì sợi dây vẫn vỗ trong màng xà phòng. ố n s ầ T Hz cm Biên độ 1 a) Chiều dài dây cm b) Chiều dài dây cm c) Chiều dài dây cm Hình 9 Tần số và biên dộ vỗ như hàm của chiều dài sợi dây. So sánh giữa kết quả thí nghiệm của Zhang [1], kết quả mô phỏng bằng FEM [2] và kết quả mô phỏng hiện tại. (a) Kết quả thí nghiệm. (b) Kết quả bằng FEM. (c) Kết quả hiện tại. 8
9. Tần số Hz Biên độ cm Biên độ a) Chiều dài dây cm b) Chiều dài dây cm c) Chiều dài dây cm Hình 10 Ảnh hưởng của hệ số uốn sợi dây u đến tần số và biên độ vỗ như là một hàm theo chiều dài sợi dây. (a) Kết quả mô phỏng bằng FEM [15] ( ). (b) Kết quả mô phỏng hiện tại ( ). (c) Kết quả mô phỏng hiện tại ( ). Khi tăng hệ số uốn từ 0,1 lên 0,5, quan dây đối xứng nhau, trong thí nghiệm của sát hình 10 (b) và (c) ta thấy tần số vỗ của 2 Zhang [1] đã chỉ ra rằng khi khoảng cách d trường hợp gần như không khác nhau theo giữa hai sợi dây đủ nhỏ (d/L 0,21±0,04), thì hai sợi dây giữa hai trường hợp. Trường hợp hình 10 (b) sẽ vỗ lệch pha nhau. Các thông số mô phỏng sợi dây ở trạng thái duỗi thẳng khi độ dài dây trong phần này được thể hiện trong bảng 3. là 2,2 cm và 2,6 cm, trường hợp hình 10 (a) Ban đầu vị trí hai sợi dây đối xứng nhau sợi dây ở trạng thái duỗi thẳng khi độ dài dây (xem hình 11), dây bên trái lệch với phương 0 là 1,5 cm và 1,8 cm. So sánh hai trường hợp thẳng đứng một góc 1 5 và dây bên phải này cho thấy khi hệ số uốn tăng thì sợi dây sẽ 0 lệch với phương thẳng đứng một góc 2 5 . tiếp cận trạng thái duỗi thẳng ở chiều dài dây lớn hơn so với khi hệ số uốn thấp hơn (bảng d 2). Sợi dây Bảng 2: Trạng thái duỗi thẳng ở hệ số uốn và 1 chiều dài dây khác nhau 2 Hz u 0 ,5 u 0 ,1 0 2,2 cm 1,5 cm Hình 11 Vị trí ban đầu 2 sợi dây trong màng xà phòng 0 2,6 cm 1,8 cm 50 3 cm 2 cm Trong phần tiếp theo của bài báo này sẽ Hình 12 thể hiện tọa độ dao động của thực hiện mô phỏng với trường hợp hai sợi đầu tự do hai sợi dây. Kết quả này tương đồng dây. Một vài thông số mô phỏng ở đây tương với phát hiện trong thí nghiệm của Zhang [1] tự như đối với phần mô phỏng của một sợi là d/L < 0,21±0,04 thì hai sợi dây sẽ vỗ song dây ở bảng 1. Kích thước lưới thực hiện là song với nhau. 200×600, chiều dài mỗi sợi dây là 3,6 cm và được rời rạc ra làm 25 điểm. Với bước thời -4 gian ∆t=10 . Vận tốc ban đầu v0 220 cm/s, và số Reynolds là Re = 198. Ban đầu các sợi 9
10. Bảng 6.3 Thông số thí nghiệm của Zhang [1] và mô phỏng hiện tại của hai sợi dây Giá trị mô phỏng Giá trị thí nghiệm [1] hiện tại Vận tốc dòng vào (v0) 220 cm/s 220 cm/s -4 2 -4 2 Mật độ khối lượng xà phòng ( 0 ) 3×10 g/cm 3×10 g/cm Độ nhớt động xà phòng (μ) 1,2×10-5 g/cm.s 1,2×10-3 g/cm.s Bề rộng kênh màng xà phòng 8,5 cm 8,5 cm Chiều dài kênh màng xà phòng 120 cm 25.5 cm Mật độ khối lượng của sợi dây (M) 2×10-4 g/cm 3×10-4 g/cm Đường kính sợi dây 0,015 cm - Diện tích mặt cắt ngang sợi dây 1,767 cm2 - Độ cứng uốn sợi dây (b ) 0,1 erg.cm 0,1 erg.cm Chiều dài sợi dây (L) 3.6 cm 3,6 cm 2 2 Gia tốc trọng lượng (gy) 980 cm/s 980 cm/s Hệ số cản không khí (λ) - 00g / v Sợi dây bên phải Sợi dây bên trái do ự u t ầ trí đ trí ị V (cm) Hình 12 Tọa độ dao động đầu tự do hai sợi dây với khoảng cách d/L = 0,15 Mặt khác, đối với khoảng cách d/L = 0,25, hai Zhang [1] là d/L > 0,21±0,04 thì hai sợi dây sợi dây đối xứng nhau qua trục trung tâm sẽ vỗ ngược pha nhau. thẳng đứng, mỗi sợi dây lệch một góc 50 so Trên đây ta mới khảo sát đối với trường với phương thẳng đứng. Kết quả này tương hợp hai sợi dây nghiêng và đối xứng nhau. đồng với phát hiện trong thí nghiệm của Phần tiếp theo sẽ khảo sát trường hợp hai sợi dây không song song và cùng nghiêng về một 10
11. phía, sợi bên trái nghiêng và hợp với phương 5. Kết luận thẳng đứng 1 góc 100, sợi bên phải hợp 1 Trong bài báo cho thấy mô phỏng số với phương thẳng đứng 1 góc 50 (xem 2 của một sợi dây mềm vỗ trong màng xà hình 13). phòng đang chảy, thực nghiệm nghiên cứu của Zhang và các đồng nghiệp [1] như một mô hình cho một phiên bản hai chiều của một Sợi dây lá cờ trong gió. Các kỹ thuật tính toán ở đây là phương pháp IBM, mô phỏng trong luận 1 văn này có kết quả tương đồng với thí nghiệm 2 của Zhang [1] và mô phỏng bằng FEM của T. Sawada và T. Hisada [2]. Thí nghiệm của Zhang cũng phát hiện ra rằng, khi chiều dài ở giữa L và L , sợi dây sẽ ở một trong hai Hình 13 Vị trí ban đầu 2 sợi dây trong màng lc uc xà phòng trạng thái duỗi thẳng hoặc vỗ tùy thuộc vào điều kiện ban đầu hoặc các yếu tố gây nhiễu Với vị trí ban đầu của hai sợi dây như bên ngoài. Với kỹ thuật tính toán bằng trong hình 13 thì kết quả mô phỏng cho ta phương pháp IBM cũng đã thể hiện được hai thấy với d/L = 0,25 thì hai sợi dây lại vỗ cùng trạng thái trên của sợi dây trong màng xà pha với nhau, khi đó trong thí nghiệm của phòng, ngoài ra mô phỏng bằng IBM còn bổ Zhang lại cho thấy nếu d/L>0,21±0,04 thì hai sung một số kết quả sau: sợi dây sẽ vỗ ngược pha với nhau. Tiếp tục thực hiện mô phỏng mô hình Vị trí ban đầu sợi dây không ảnh hưởng trên cũng với d/L = 0,25, nhưng ta thay đổi tới trạng thái vỗ của sợi dây. 00 góc 12 0,0 , nghĩa là hai sợi dây sẽ Sự thay đổi độ lớn hệ số uốn sợi dây song song với phương thẳng đứng. Kết quả không ảnh hưởng lớn tới tần số dao cho thấy hai sợi dây vỗ ngược pha với nhau động mà thay đổi về chiều dài tới hạn đúng với phát hiện trong thí nghiệm của phân nhánh giữa hai trạng thái duỗi Zhang. thẳng và vỗ. Khi khoảng cách d/L = 0,35 ở các góc α Với hai sợi dây cạnh nhau trong màng khác nhau hai sợi dây vẫn vỗ ngược pha nhau, xà phòng đang chảy, tác giả đã tìm thấy hai ngay cả trường hợp ban đầu hai sợi dây đan trạng thái vỗ khác nhau. Hai trạng thái này chéo (hình 14). xảy ra phụ thuộc vào các tham số không thứ nguyên d/L, trong đó d là khoảng cách giữa Sợi dây hai sợi dây và L là chiều dài của sợi dây. Với d/L=0,15 cho thấy hai sợi dây vỗ cùng pha với nhau (song song với nhau). Với d/L=0,35 hai sợi dây sẽ vỗ lệch pha nhau. Các kết quả 30° 5° mô phỏng này tương đồng với kết quả thực nghiệm của Zhang [1] rằng khi khoảng cách d giữa hai sợi dây đủ nhỏ (d/L 0,21±0,04), thì hai sợi dây sẽ vỗ lệch pha nhau. Nhưng trong phần mô Hình 14 Vị trí ban đầu 2 sợi dây trong phỏng với hai sợi dây tác giả đã phát hiện ra màng xà phòng 3000 , 5 12 rằng kết luận trong thí nghiệm của Zhang vẫn chưa đúng hoàn toàn, với trường hợp d/L = 0,25 > 0,21±0,04 theo thí nghiệm của Zhang 11
12. thì hai sợi dây phải vỗ ngược pha nhau. the immersed boundary method, Journal of Nhưng kết quả mô phỏng hiện tại đã phát hiện Computational Physics 153, 1999, 509–534. rằng với d/L = 0,25 hai sợi dây có thể ở hai [5]. M. Griebel, T. Dornseifer, T. trạng thái dao động khác nhau, cùng pha hoặc Neunhoeffer, Numerical simulation in fluid ngược pha tùy thuộc vào vị trí ban đầu của hai dynamics: A practical introduction, Society sợi dây. Cụ thể là khi khoảng cách giữa hai for Industrial and Applied Mathematics. đầu tự do d đủ lớn thì hai sợi dây sẽ vỗ td Philadelphia, 1998. ngược pha, ngược lại khi dtd đủ nhỏ thì hai sợi dây sẽ vỗ cùng pha. [6]. Charles S. Peskin, The immersed boundary method, Handwritten lecture notes Như vậy, kết luận của Zhang khi (under construction) of a course taught at the khoảng cách d giữa hai sợi dây đủ nhỏ (d/L 0,21±0,04), thì hai sợi dây sẽ vỗ lệch pha nhau đã không còn cture_notes/index.html đúng. Mà hai sợi dây sẽ thể hiện ở một trong ba trạng thái dao động: cùng pha, ngược pha và cùng pha hoặc ngược pha, tùy thuộc vào khoảng cách d/L và vị trí giữa hai đầu tự do. Giả sử ta có hai khoảng cách d/L1 và d/L2, trong đó d/L1 < d/L2, tác giả có kết luận sau Khi khoảng cách d/L1 đủ nhỏ thì hai sợi dây sẽ vỗ cùng pha. Khi khoảng cách d/L2 đủ lớn thì hai sợi dây sẽ vỗ ngược pha. Khi khoảng cách d/L nằm giữa d/L1 và d/L2 thì hai sợi dây sẽ vỗ cùng pha hoặc ngược pha tùy vào vị trí ban đầu của hai đầu tự do sợi dây. TÀI LIỆU THAM KHẢO [1]. J. Zhang, S. Childress, A. Libchaber, and M. Shelley, Flexible filaments in a flowing soap film as a model for one-dimensional flags in a two-dimensional wind, Nature 408 , 2000, 835–839. [2]. T. Sawada, T. Hisada, Fluid-structure interaction analysis of the two-dimensional flag-in-wind problem by an interface-tracking ALE finite element method, Computer & Fluid 36, 2007, 136–146. [3]. Charles S. Peskin, The immersed boundary method, Acta Numer.11 (2), 2002, 479–517. [4]. Alexandre M. Roma, Charles S. Peskin, and Marsha J. Berger, An adaptive version of 12
13. BÀI BÁO KHOA HỌC THỰC HIỆN CÔNG BỐ THEO QUY CHẾ ĐÀO TẠO THẠC SỸ Bài báo khoa học của học viên có xác nhận và đề xuất cho đăng của Giảng viên hướng dẫn Bản tiếng Việt ©, TRƯỜNG ĐẠI HỌC SƯ PHẠM KỸ THUẬT TP. HỒ CHÍ MINH và TÁC GIẢ Bản quyền tác phẩm đã được bảo hộ bởi Luật xuất bản và Luật Sở hữu trí tuệ Việt Nam. Nghiêm cấm mọi hình thức xuất bản, sao chụp, phát tán nội dung khi chưa có sự đồng ý của tác giả và Trường Đại học Sư phạm Kỹ thuật TP. Hồ Chí Minh. ĐỂ CÓ BÀI BÁO KHOA HỌC TỐT, CẦN CHUNG TAY BẢO VỆ TÁC QUYỀN! Thực hiện theo MTCL & KHTHMTCL Năm học 2016-2017 của Thư viện Trường Đại học Sư phạm Kỹ thuật Tp. Hồ Chí Minh.