So sánh trường tốc độ của tia phun rối, khuếch tán tính theo mô hình tích phân và code CFD đa phương fluent 6.0

rối, khuếch tỏn. Tuy nhiờn để cú thể tổng quỏt húa việc ỏp dụng, mụ hỡnh cần được đỏnh giỏ bằng kết quả cho bởi cỏc phần mềm đa phương cú sẵn. Bài bỏo này so sỏnh trường tốc độ cho bởi mụ hỡnh tớch phõn và code CFD Fluent 6.0. Sai lệch giữa hai mụ hỡnh nằm trong giới hạn 10% khi số Reynolds ở miệng vũi phun nhỏ hơn 5000. ABSTRACT The integral model is simple in utilization, low CPU time calculation, suitable for a lot of pratical applications of turbulent diffusion jet. However, for a general application, the model should be assessed by the results of available multidirectional codes. This paper shows the comparison of velocity profiles given by the integral model and the CFD FLUENT 6.0 Code. The difference in results of the two models is less than 10% when the Reynolds number at the exit nozzle is lower than 5000. 1. Giới thiệu Tia phun rối và khuếch tán có rất nhiều ứng dụng trong kỹ thuật. Tr-ớc đây, việc nghiên cứu tia phun đ-ợc tiến hành chủ yếu bằng thực nghiệm và những qui luật cơ bản rút ra đ-ợc từ các nghiên cứu này đã có những ứng dụng thiết thực trong công nghiệp, đặc biệt trong lĩnh vực động cơ đốt trong. Ngày nay, với sự phát triển của các công cụ tin học, bài toán tia phun rối, khuếch tán đã đ-ợc nghiên cứu một cách t-ờng tận nhờ các phần mềm tính toán động học chất lỏng (CFD). Sự phát triển của tia phun trong những điều kiện khác nhau, kể cả những tr-ờng hợp mà tr-ớc đây thực nghiệm khó hay không thể thực hiện đ-ợc, đã đ-ợc xác định. Tuy nhiên những phần mềm nh- vậy rất phức tạp, thời gian tính toán kéo dài, đôi lúc không phù hợp với thực tiễn áp dụng. Do đó, việc xây dựng các công cụ toán học đơn giản hơn nhằm hỗ trợ cho nghiên cứu ứng dụng tia phun rối, khuếch tán là rất cần thiết. Công cụ nh- vậy cần đ-ợc thiết lập trên cơ sở hệ ph-ơng trình tích phân mô tả sự biến thiên của các đại l-ợng vật trung bình theo ph-ơng trục tia kết hợp với các qui luật thực nghiệm về diễn biến của chúng theo ph-ơng h-ớng kính [6], [7], [11]. Mô hình đơn giản mô tả tia phun rối đ-ợc thiết lập trong môi tr-ờng không khí đứng yên. Mô hình này có ý nghĩa trong kiểm chứng các điều kiện biên và tính chính xác của các hệ số thực nghiệm sử dụng. Trong thực tế, dù trong buồng cháy động cơ hay ngòai khí quyển, tia phun cũng chịu những tác động của môi tr-ờng không khí vận động. Vì vậy mô hình tia phun có tính tổng quát đ-ợc xây dựng trong điều kiện có sự t-ơng tác của môi tr-ờng [8]. Tuy kết quả cho bởi mô hình tích phân và thực nghiệm trong các tr-ờng hợp cụ thể khá phù hợp [9], [10] nh-ng để có thể tổng quát hóa cho những tr-ờng hợp áp dụng khác nhau, mô hình này cần đ-ợc đánh giá bởi những kết quả của phần mềm đa ph-ơng. Trong bài báo này, chúng tôi so sánh kết quả tr-ờng tốc độ cho bởi mô hình tích phân đã thiết lập với phần mềm đa ph-ơng FLUENT 6.0. B-ớc đầu việc đánh giá đ-ợc thực hiện trong cùng điều kiện tính tóan đối với tia phun thẳng đứng trong môi tr-ờng không khí đứng yên. Vòi phun có đ-ờng kính 2 và 3 mm. Vận tốc phun thay đổi từ 50 đến 100 m/s. Môi chất trong tia phun là khí dầu mỏ hóa lỏng LPG. 2. Hệ ph-ơng trình không chế tia phun Hệ ph-ơng trình khống chế tia phun rối, khuếch tán nghiêng một góc bất kỳ trong môi tr-ờng không khí chuyển động ngang đã đ-ợc trình bày trong [7]. Hệ ph-ơng trình bao gồm các ph-ơng trình bảo tòan viết d-ới dạng tích phân và mô hình rối k- tiêu chuẩn. Kết quả cho bởi mô hình là biến thiên của các đại l-ợng vật lý theo ph-ơng h-ớng trục. Biến thiên của chúng theo ph-ơng h-ớng kính đ-ợc xác định theo qui luật đồng dạng [10]. Mô hình đã đ-ợc đánh giá bằng số liệu thực nghiệm của tia phun rối, khuếch tán ngoài khí quyển và trong buồng cháy động cơ Diesel [7]. Hình 1a,b,c trình bày kết quả so sánh giữa mô hình và thực nghiệm một số tr-ờng hợp tiêu biểu. Số liệu thực nghiệm về tr-ờng tốc độ đ-ợc đo bằng ph-ơng pháp Laser Doppler [11]. Kết quả cho thấy mô hình tích phân cho giá trị hơi thấp hơn thực nghiệm ở các bán kính trung gian của profil. Kết quả tính toán tiêu biểu về tr-ờng tốc độ theo mô hình tích phân đ-ợc trình bày trên hình 2. Trong mô hình đa ph-ơng của code FLUENT, hệ ph-ơng trình mô tả tia phun đ-ợc viết tổng quát nh- sau: Ph-ơng trình liên tục: i i i i u u 0 t x t x             (1) Ph-ơng trình bảo tòan động l-ợng:   i ji j ji k ij i j j j j i k j u uu t x up u 2 u u u x x x x 3 x x                                 (2) Trong ph-ơng trình (8) dấu l-ợn sóng (~) để chỉ trung bình Favre. ứng suất Reynolds đ-ợc mô hình hoá nhờ giả thiết của Boussinesq [4] liên kết giữa sức căng Reynolds và các gradient vận tốc trung bình: ji i i j t t ij j i i uu 2 u u u k x x 3 x                       (3) Trong đó đại l-ợng ứng suất Reynolds đ-ợc mô hình hoá bởi mô hình k-  tiêu chuẩn gồm hai ph-ơng trình của Launder và Spalding [2]: 0 2 4 6 0 5 10 15 20 r (mm) U ( m /s ) X = 100mm 0 1 2 3 0 20 40 X = 200mm r (mm) U ( m /s ) 0 0,4 0,8 1,2 1,6 0 20 40 60 80 r (mm) U ( m /s ) X = 400mm Hình 1: So sánh profil tốc độ cho bởi mô hình tích phân và thực nghiệm Hình 2: Tr-ờng tốc độ của tia phun tính theo mô hình tích phân i k b M k i j k j k ( k) ( ku ) ( ) G G Y S t x x x                           (4)   2 i 1 k 3 b 2 i j j ( ) ( u ) ( ) C G C G C S t x x x k k                                  (5) k và  đ-ợc định nghĩa nh- sau: i i 1 k u u 2   và i i j j u u x x         (6) Trong các ph-ơng trình trên: - Gk là đại l-ợng sản sinh năng l-ợng rối do các gradient vận tốc trung bình gây ra. Theo giả thuyết Boussinesq Gk có thể đ-ợc biểu diễn gần đúng bằng biểu thức sau: 2 kG S  trong đó 2 kG S  (7) - Gb là đại l-ợng sản sinh động năng rối do tác động của nhiệt độ và lực trọng tr-ờng: t b i t i T G g Pr x      (8) Trong đó Prt là số Prandlt rối, th-ờng chọn Prt=0.85, gi là thành phần vectơ gia tốc trọng tr-ờng theo h-ớng i,  là hệ số giãn nở nhiệt cho bởi biểu thức: p 1 p T           (9) Trong khi ảnh h-ởng của lực trong tr-ờng đến k có thể xác định một chính xác thì ảnh h-ởng của nó đối với  vẫn còn nhiều ý kiến khác nhau. Trong code FLUENT ảnh h-ởng của lực trọng tr-ờng đến  đ-ợc đơn giản hoá bằng cách đặt Gb=0 trong ph-ơng trình vận chuyển của . Tuy nhiên ảnh h-ởng này lại đ-ợc tính đến trong hệ số C3e (5) đ-ợc tính theo biểu thức: 3 tanh v C u   (10) Trong đó v và u là hai thành phần vận tốc song song và vuông góc với ph-ơng của lực trọng tr-ờng, do đó C3e =1 trong tr-ờng hợp vận tốc song song với lực trọng tr-ờng và C3e =0 khi vận tốc vuông góc với ph-ơng của lực trọng tr-ờng [2]. - YM biểu thị ảnh h-ởng của sự thay đổi thể tích trong dòng chảy. Đại l-ợng này th-ờng đ-ợc bỏ qua trong các mô hình đối với dòng chảy không chịu nén. Trong dòng chảy chịu nén YM đ-ợc tính theo quan hệ Sarkar [4] 2 M tY 2 M  (11) Trong đó Mt là số Mach: t 2 k M a  , a là vận tốc âm thanh. -  là độ nhớt rối đ-ợc tính thông qua k và : Hình 3: Tính tóan tia phun thẳng đứng 2D a. Chia l-ới b. Kết quả tính tiêu biểu a. b. 2k C     (12) - C1, C2, C3 và C là các hằng số. k, , là các số Prandtl cho k và . Sk và S là các hàm số do ng-ời sử dụng định nghĩa. Trong FLUENT các giá trị hằng số sau đây đ-ợc sử dụng: C1 = 1.44, C2 = 1.92, C = 0.09, k = 1.0,  = 1.3. Hệ ph-ơng trình trên đ-ợc giải bằng ph-ơng pháp thể tích hữu hạn. Không gian của miền tính toán đ-ợc chia thành một số l-ợng hữu hạn các phần tử. Các điểm nút l-ới là trọng tâm của các khối. Các ph-ơng trình bảo toàn đ-ợc áp dụng cho mỗi phần tử và các giá trị của các biến sẽ đ-ợc tính toán tại trung tâm của khối [3]. Ph-ơng pháp nội suy sau đó sẽ đ-ợc áp dụng để tìm ra các giá trị tại bề mặt của khối. Trong nghiên cứu này, chúng tôi sử dụng ph-ơng pháp nội suy sai phân tiến [1]. 3. Chia l-ới tia phun L-ới tính 2D đ-ợc thiết lập với phần mềm Gambit 2.1. Đối với tr-ờng hợp tia phun thẳng đứng miền tính toán đ-ợc chọn đối xứng, với bán kính 500mm và chiều cao 1000mm. B-ớc l-ới dày ở phần trục tia và th-a dần về phía ngoài theo ph-ơng h-ớng kính và t-ơng tự nh- vậy theo ph-ơng h-ớng trục kể từ miệng vòi phun. L-ới tính bao gồm 2600 phần tử hình chữ nhật và 2727 nút. L-ới tính toán đ-ợc trình bày trên hình 3a. Đầu vào miền tính toán là đ-ờng kính miệng vòi phun d=0.002m và d=0.003m, vận tốc phun thay đổi từ 50m/s đến 100m/s. Điều kiện biên đ-ợc đặt là tốc độ ra khỏi miệng vòi phun với c-ờng độ rối 10%. Các mặt xung quanh là các biên đối xứng. Hình 3b trình bày kết quả tính toán tiêu biểu trong tr-ờng hợp tia phun thẳng đứng trong môi tr-ờng không khí đứng yên. 4. Kết quả và bình luậN Hình 4 trình bày kết quả so sánh profil tốc độ cho bởi mô hình tích phân và mô hình đa ph-ơng FLUENT. Các hệ số tính toán của hai mô hình đ-ợc chọn thống nhất. Profil tốc độ trong mô hình đa ph-ơng đ-ợc giải cho từng điểm theo ph-ơng h-ớng kính, trong khi đó, ở mô hình tích phân, chúng đ-ợc xác định theo qui luật đồng dạng do Ebrahini [5] đề nghị:  2cU( ) U .exp ln(2)   (13) Hình 4: So sánh profil tốc độ cho bởi mô hình tích phân và FLUENT (D=2mm, Uo=30m/s, Re=3.000) Hình 5: So sánh profil tốc độ cho bởi mô hình tích phân và FLUENT (D=2mm, Uo=100m/s, Re=10.000) Hình 6: So sánh profil tốc độ cho bởi mô hình tích phân và FLUENT (D=3mm, Uo=50m/s, Re=7.500) Trong đó Uc là tốc độ trên trục tia;  là độ dài không thứ nguyên: 0,5 r r  (14) r0,5 là bán kính tại điểm có vận tốc bằng 0,5Uc Khi tốc độ phun thấp với số Reynolds ở miệng vòi phun khỏang 3000, sai lệch của hai mô hình nhỏ, chủ yếu diễn ra ở khu vực gần trục tia (hình 4). Biên dạng của profil hầu nh- phù hợp với nhau. Điều này cho thấy qui luật đồng dạng (13) phù hợp với tia phun rối khuếch tán. Khi số Reynolds ở miệng vòi phun tăng lên, mức độ sai lệch về profil tốc độ cho bởi hai mô hình gia tăng (hình 5). Khi đ-ờng kính vòi phun 2mm và tốc độ phun 100m/s (số Reynolds 10.000), sai lệch tốc độ cực đại cho bởi hai mô hình khỏang 15%. Kết quả t-ơng tự nhận đ-ợc khi đ-ờng kính vòi phun tăng lên 3mm với số Reynolds 7.500 (hình 6). Sự khác biệt về profil tốc độ ở khu vực gần trục tia là do trong mô hình tích phân chúng ta áp dụng giả thuyết ôtop hat profileằ. Theo giả thuyết này tốc độ qua mỗi mặt cắt ngang đ-ợc giả định là hằng số sao cho tích phân l-u l-ợng và động l-ợng qua mặt cắt này đ-ợc bảo toàn. Kết quả so sánh trên cho thấy mô hình tích phân một chiều cho kết quả phù hợp với mô hình đa ph-ơng đối với tia phun rối có số Reynolds ở miệng vòi phun bé. Nếu xem sự khác biệt cực đại về tốc độ cho bởi hai mô hình là một hàm số theo số Reynolds ở miệng vòi phun thì theo kết quả nghiên cứu trên đây chúng ta thấy khi số Reynolds nhỏ hơn 5000 thì sự khác biệt này nằm trong giới hạn 10%. Giới hạn số Reynolds này phù hợp với hầu hết các hệ thống phun áp dụng trong động cơ đốt trong. 5. Kết luận Biên dạng của các profil tốc độ cho bởi mô hình tích phân một chiều phù hợp với mô hình đa ph-ơng. Giá trị vận tốc trên trục tia phun th-ờng nhỏ hơn mô hình đa ph-ơng do giả thiết top hat profile áp dụng trong tính tóan tích phân l-u l-ợng và động l-ợng. Sai lệch về giá trị vận tốc cho bởi hai mô hình nhỏ hơn 10% khi số Reynolds ở miệng vòi phun nhỏ hơn 5000. Vì vậy mô hình tích phân có thể áp dụng trong hầu hết cấu hình tia phun trong động cơ đốt trong. Tài liệu tham khảo [1] Fluent 6.0 User's Guide, Fluent, Inc., Cennterra Resource Park, Lebanon, NH 30766, 2002. [2] B. E. Launder, D. B. Spalding, Lectures in Mathematical Models of Turbulence, Academic Press, London, England, 1972. [3] J. Y. Murthy, S. R. Mathur, A Finite Volume Method For Radiative Heat Transfer Using Unstructured Meshes, AIAA-98-0860, 1998. [4] S. Sarkar, L. Balakrishnan, Application of a Reynolds-Stress Turbulence Model to the Compressible Shear Layer, ICASE Report 90-18, NASA CR 182002, 1990. [5] EBRAHINI I., KLEINE R, The nozzle fluid concentration fluctuation field in round turbulent fuel jets and jet diffusion flames, Sixteenth symposium (International) on Combustion, pp. 1711-1723, 1976. [6] Bui Van Ga, An Integral Model for Calculation of LPG Jet Development in Combustion Chamber Of Spark Ignition Engine, International Conference on HPSC, Hanoi 10-14 March 2003. [7] BUI VAN GA, PHAM XUAN MAI, Liviu GEORGESCU, A mathematical model for calculation of turbulence diffusion combustion in air and in Diesel engines, Proceedings of the VII International Conference of Motor Vehicles CAR-2000 (FISITA, SIAR), Romania, 16-17 Nov. 2000, Vol ICE, pp 8-16. [8] Bựi Văn Ga, Nhan Hồng Quang, J.M. Vignon, Calculation of turbulent diffusion jets under effects of gravity and moving surrounding air, Vietnamese Journal of Mechanics, Vol. 23, No. 2, pp. 87-94, 2001. [9] Bựi Văn Ga, Dương Việt Dũng, Trần Văn Nam, Mụ phỏng tia phun khớ dầu mỏ húa lỏng LPG trong buồng chỏy động cơ đỏnh lửa cưỡng bức, Khoa Học và Cụng Nghệ, No. 30-31, pp. 97- 103, 2001.