Luận văn Tính toán tấm composite cốt hạt có tính đến sự truyền nhiệt

x2 [ −D (∂ 2w ∂ x2 +ν ∂ 2w ∂ y2 ) − MT 1−ν ] −2D(1−ν) ∂ 4w ∂ x2∂ y2+ + ∂ 2 ∂ y2 [ −D (∂ 2w ∂ y2 +ν ∂ 2w ∂ x2 ) − MT 1−ν ] +Nx ∂ 2w ∂ x2 +2Nxy ∂ 2w ∂ x∂ y +Ny ∂ 2w ∂ y2 = 0 . hay D (∂ 4w ∂ x4 +2 ∂ 4w ∂ x2∂ y2 + ∂ 4w ∂ y4 ) + 1 1−ν (∂ 2MT ∂ x2 + ∂ 2MT ∂ y2 ) − −Nx ∂ 2w ∂ x2 −2Nxy ∂ 2w ∂ x∂ y −Ny ∂ 2w ∂ y2 = 0 . (1.43) phương trình (1.43) có thể viết lại như sau: D∇2∇2w+ 1 1−ν ∇ 2MT −Nx ∂ 2w ∂ x2 −2Nxy ∂ 2w ∂ x∂ y −Ny ∂ 2w ∂ y2 = 0 . (1.44) trong đó ∇2 là toán tử Laplace ∇2 = ∂ 2 ∂ x2 + ∂ 2 ∂ y2 . 1.5. Điều kiện biên Xét một số trường hợp đơn giản về điều kiện biên thường gặp trong các bài toán đối với tấm hình chữ nhật 1. Biên tựa bản lề: (Hình 1.8) chẳng hạn cạnh y = 0 tựa tự do, khi đó điều kiện biên là (w)y=0 = 0 , (My)y=0 = 0 . (1.45) hay (w)y=0 = 0 ,[ −D (∂ 2w ∂ y2 +ν ∂ 2w ∂ x2 ) − MT 1−ν ] y=0 = 0 (1.46) 17 Hệ thức (1.45) có thể viết lại như sau( −D∂ 2w ∂ y2 − MT 1−ν ) y=0 = 0 (1.47) 2. Biên bị ngàm: (Hình 1.9) ví dụ cạnh y = 0 bị ngàm, thì điều kiện biên sẽ là (w)y=0 = 0 , (∂ w ∂ y ) y=0 = 0 . (1.48) 3. Biên tự do: giả sử cạnh y = 0 tự do. Theo Poisson, điều kiện biên sẽ là (My)y=0 = 0 , (Mxy)y=0 = 0 , (Qy)y=0 = 0 . (1.49) Nhưng theo Kirchhoff, tại cạnh đó điều kiện biên sẽ là (My)y=0 = 0 , ( Qy + ∂ Mxy ∂ x ) y=0 = 0 . (1.50) Thay (1.33) vào (1.37) và (1.38), ta được biểu thức xác định lực cắt như sau Qx =−D (∂ 3w ∂ x3 + ∂ 3w ∂ x∂ y2 ) − 1 1−ν ∂ MT ∂ x , (1.51) Qy =−D (∂ 3w ∂ y3 + ∂ 3w ∂ x2∂ y ) − 1 1−ν ∂ MT ∂ y . (1.52) Khi đó,(1.50) trở thành[ −D (∂ 2w ∂ y2 +ν ∂ 2w ∂ x2 ) − MT 1−ν ] y=0 = 0 , [ −D (∂ 3w ∂ y3 +(2−ν) ∂ 3w ∂ x2∂ y ) − 1 1−ν ∂ MT ∂ y ] y=0 = 0 . (1.53) 18 Chương 2 Uốn tấm composite mỏng khi có truyền nhiệt dừng 2.1. Modun đàn hồi và hệ số dãn nở nhiệt của composite cốt hạt Giả thiết tấm composite có độn các hạt hình cầu có độ dài hai cạnh là a,b, chiều dày h, tấm tựa bản lề tại các cạnh. Ta xét trường hợp tấm có sự trao đổi nhiệt đối lưu ổn định trên bề mặt z = ±h/2, và khi có sự khác nhau lớn về nhiệt độ giữa các mặt z =±h/2, thì có gradient nhiệt độ cao thông qua bề dày của tấm gây ra uốn nhiệt. Bài toán đặt ra là hãy xác định độ uốn của tấm khi có truyền nhiệt dừng. Để giải được bài toán đó cần xác định được môdun đàn hồi cho composite. Giả thiết nền và hạt đều là vật liệu đàn hồi đồng nhất đẳng hướng, có tính đến tương tác giữa nền và hạt. Khi đó, ta có [18]: K = Km + (Kc−Km)ξ 1+(Kc−Km) /( Km + 43Gm ) , G = Gm− 15(1−νm)(Gm−Gc)ξ 7−5νm +(8−10νm) GcGm . (2.1) trong đó, ξ = N ∑ i=1 Vi V là tỷ lệ thể tích các hạt độn. Gm,Gc tương ứng là môdun trượt của pha nền và pha hạt; Km,Kc tương ứng là môdun kéo nén thể tích của pha nền, pha hạt; νm,νc tương ứng là hệ số Poisson của pha nền và pha hạt. Với αm,αc là hệ số dãn nở nhiệt của pha nền, pha hạt tương ứng, khi đó hệ số dãn nở nhiệt của composite cốt hạt được xác định như công thức sau [4]: α = αm +(αc−αm) Kc (3Km +4Gm)ξKm (3Kc +4Gm)+4(Kc−Km)Gmξ . (2.2) .2.2. Sự phân bố nhiệt độ trong tấm Do có sự trao đổi nhiệt đối lưu trên bề mặt tấm, và gradient nhiệt theo bề dày của tấm, nên phương trình truyền nhiệt dừng (1.21) có dạng [12]: d2T dz2 = 0 (2.3) Với điều kiện biên (1.25) sẽ là  k∂ T∂ z = β1 (T −T1) tại z =− h 2 , −k∂ T∂ z = β2 (T −T2) tại z = h 2 , (2.4) trong đó, T1,β1 tương ứng là nhiệt độ môi trường, hệ số truyền nhiệt đối lưu tại mặt z = −h/2; T2,β2 là nhiệt độ môi trường, hệ số truyền nhiệt đối lưu tại mặt z = h/2 tương ứng; k là hệ số dẫn nhiệt của vật liệu, và giả thiết được xác định theo công thức [3] k = (1−ξ )km+ξ kc. (2.5) Đặt : γ1 = β1 k , γ2 = β2 k (2.6) Khi đó, điều kiện biên (2.4) trở thành:  ∂ T ∂ z = γ1 (T −T1) tại z =− h 2 , −∂ T∂ z = γ2 (T −T2) tại z = h 2 , (2.7) Tích phân phương trình (2.3) ta được nghiệm T có dạng T =C0 +C1z , trong đó,C0,C1 là các hằng số được xác định từ điều kiện biên. Thay biểu thức nghiệm T này vào điều kiện biên (2.7), ta được hệ gồm hai phương trình hai ẩn C0,C1  ( 1+ γ1 h 2 ) C1− γ1C0 =−γ1T1 − ( 1+ γ2 h 2 ) C1− γ2C0 =−γ2T2 (2.8) 20 Giải hệ (2.8), ta được C0 = 1 γ1 + γ2 + γ1γ2h [ T1γ1 ( 1+ γ2 h 2 ) +T2γ2 ( 1+ γ1 h 2 )] , C1 = γ1γ2 (T2−T1) γ1 + γ2 + γ1γ2h , (2.9) Với C0,C1 xác định theo (2.9), thay vào biểu thức nghiệm của T , ta thu được T = 1 γ1 + γ2 + γ1γ2h [ T1γ1 ( 1+ γ2 h 2 ) +T2γ2 ( 1+ γ1 h 2 ) + γ1γ2(T2−T1)z ] (2.10) Khi đó, thế (2.10) vào biểu thức (1.11), dẫn tới ∆T =−T0 + 1γ1 + γ2 + γ1γ2h [ T1γ1 ( 1+ γ2 h 2 ) +T2γ2 ( 1+ γ1 h 2 )] + + γ1γ2(T2−T1) γ1 + γ2 + γ1γ2h z . (2.11) 2.3. Uốn tấm composite mỏng khi có truyền nhiệt dừng 2.3.1. Mặt giữa không biến dạng Thay (1.31) vào (1.35) và (1.36), tương ứng ta được hệ sau  Eh 1−ν2 (∂ 2u ∂ x2 +ν ∂ 2v ∂ x∂ y ) − 1 1−ν ∂ NT ∂ x + Eh 2(1+ν) (∂ 2u ∂ y2 + ∂ 2v ∂ x∂ y ) = 0 Eh 1−ν2 (∂ 2v ∂ y2 +ν ∂ 2u ∂ x∂ y ) − 1 1−ν ∂ NT ∂ y + Eh 2(1+ν) ( ∂ 2u ∂ x∂ y + ∂ 2v ∂ x2 ) = 0 hay hệ trên có thể viết lại như sau  Eh 1−ν2 ∂ 2u ∂ x2 + Eh 2(1+ν) ∂ 2u ∂ y2 + Eh 2(1−ν) ∂ 2v ∂ x∂ y − 1 1−ν ∂ NT ∂ x = 0 Eh 1−ν2 ∂ 2v ∂ y2 + Eh 2(1+ν) ∂ 2v ∂ x2 + Eh 2(1−ν) ∂ 2u ∂ x∂ y − 1 1−ν ∂ NT ∂ y = 0 (2.12) thay biểu thức ∆T ở (2.11) vào biểu thức NT ở (1.34), tích phân ta được NT = N1 = const (2.13) khi đó, hệ (2.12) trở thành  1 1−ν2 ∂ 2u ∂ x2 + 1 2(1+ν) ∂ 2u ∂ y2 + 1 2(1−ν) ∂ 2v ∂ x∂ y = 0 1 1−ν2 ∂ 2v ∂ y2 + 1 2(1+ν) ∂ 2v ∂ x2 + 1 2(1−ν) ∂ 2u ∂ x∂ y = 0 (2.14) 21 Do tấm tựa bản lề tại các cạnh nên điều kiện biên có dạng w = v = 0; Mx = 0, tại x = 0, x = a w = u = 0; My = 0, tại y = 0, y = b Nghiệm u,v thỏa mãn điều kiện biên trên có dạng u = ∞ ∑ m=1 ∞ ∑ n=1 umn cos mpix a sin npiyb v = ∞ ∑ m=1 ∞ ∑ n=1 vmn sin mpix a cos npiy b trong đó m,n là các số tự nhiên. Thay biểu thức nghiệm của u,v vào hệ (2.14), ta có:  ∞ ∑ m=1 ∞ ∑ n=1 { umn 1 1+ν ( m2 (1−ν)a2 + n2 2b2 ) + vmn mn 2ab(1−ν) } cos mpix a sin npiy b =0 ∞ ∑ m=1 ∞ ∑ n=1 { vmn 1 1+ν ( n2 (1−ν)b2 + m2 2a2 ) +umn mn 2ab(1−ν) } sin mpix a cos npiy b = 0 Do hệ đúng với mọi x,y nên ta có  umn 1 1+ν ( 1 1−ν m2 a2 + 1 2 n2 b2 ) + vmn 1 2(1−ν) mn ab = 0 vmn 1 1+ν ( 1 1−ν n2 b2 + 1 2 m2 a2 ) +umn 1 2(1−ν) mn ab = 0 (2.15) hay hệ (2.15) có thể viết lại như sau  umn 2b2m2 +a2n2 (1−ν) ab(1+ν) + vmnmn = 0 umnmn+ vmn 2a2n2+b2m2 (1−ν) ab(1+ν) = 0 (2.16) Từ phương trình đầu của hệ (2.16) suy ra: vmn =−umn 2b 2m2 +a2n2 (1−ν) abmn(1+ν) thay vào phương trình thứ hai của hệ đó, ta có umn [ mn− ( 2b2m2 +a2n2 (1−ν))(2a2n2+b2m2 (1−ν)) a2b2mn(1+ν)2 ] = 0 22 hay umnImn = 0 , (2.17) trong đó, Imn = mn− ( 2b2m2 +a2n2 (1−ν))(2a2n2 +b2m2 (1−ν)) a2b2mn(1+ν)2 Biến đổi biểu thức Imn, ta được Imn = a2b2m2n2 (1+ν)2−(2b2m2 +a2n2 (1−ν))(2a2n2 +b2m2 (1−ν)) a2b2mn(1+ν)2 = −2(1−ν)(a2n2+m2b2)2 a2b2mn(1+ν)2 Do: −1 6 ν 6 12 nên Imn < 0. Do đó, từ (2.17) suy ra: umn = 0 suy ra vmn = 0 Vậy u = 0,v = 0 Như vậy, khi nhiệt độ được truyền theo bề dày của tấm theo quy luật (2.10), (2.11) thì nó không gây ra chuyển vị tại các điểm ở mặt giữa theo phương x,y nên không có biến dạng nhiệt tại mặt giữa. 2.3.2. Biểu thức nghiệm xác định uốn tấm Mặt giữa của tấm không bị biến dạng, khi đó, (1.31) trở thành: Nx =− NT1−ν , Ny =− NT1−ν , Nxy = 0 , (2.18) Thay biểu thức NT ở (2.13) vào (2.18), ta được Nx = Ny =− N11−ν = N0 , Nxy = 0 , (2.19) tương tự, thay ∆T từ (2.11) vào biểu thức của MT ở (1.34), thực hiện tích phân, ta thu được MT = M0 = const (2.20) 23 Thay (2.20) vào biểu thức xác định mômen (1.33), ta được Mx =−D (∂ 2w ∂ x2 +ν ∂ 2w ∂ y2 ) − M0 1−ν , My =−D (∂ 2w ∂ y2 +ν ∂ 2w ∂ x2 ) − M0 1−ν , Mxy =−D(1−ν) ∂ 2w ∂ x∂ y . (2.21) trong đó, D được xác định như trong (1.34). Thế (2.19) và (2.21) vào (1.43), ta được phương trình xác định uốn tấm trong trường hợp có sự truyền nhiệt dừng là D (∂ 4w ∂ x4 +2 ∂ 4w ∂ x2∂ y2 + ∂ 4w ∂ y4 ) −N0 (∂ 2w ∂ x2 + ∂ 2w ∂ y2 ) = 0 . (2.22) hay ∇2 ( ∇2w−N0 D w ) =0 . (2.23) Do tấm tựa bản lề tại các cạnh, nên điều kiện biên có dạng w = 0, Mx =−D (∂ 2w ∂ x2 +ν ∂ 2w ∂ y2 ) − M0 1−ν = 0, tại x = 0,x = a. w = 0, My =−D (∂ 2w ∂ y2 +ν ∂ 2w ∂ x2 ) − M0 1−ν = 0, tại y = 0,y = b. (2.24) Phương trình (2.23) có thể được đưa về hệ gồm hai phương trình sau:  ∇2 f = 0 , ∇2w− N0 D w = f . (2.25) hay cụ thể là   ( ∂ 2 ∂ x2 + ∂ 2 ∂ y2 ) f = 0 ( ∂ 2 ∂ x2 + ∂ 2 ∂ y2 ) w− N0 D w = f (2.26) Mặt khác, ta có thể xem các điều kiện biên của tựa bản lề như sau: w = 0, ∂ 2w ∂ y2 = 0, tại x = 0,x = a, w = 0, ∂ 2w ∂ x2 = 0, tại y = 0,y = b. (2.27) 24 Do đó, sử dụng các điều kiện biên (2.27), ta thu được các điều kiện biên cho hàm ẩn f : f + M0 D(1−ν) = 0, tại x = 0,x = a f + M0 D(1−ν) = 0, tại y = 0,y = b Suy ra, hàm ẩn f trong phương trình đầu tiên của hệ (2.26): f =− M0 D(1−ν) , Khi đó, phương trình thứ hai của hệ (2.26) thành:( ∂ 2 ∂ x2 + ∂ 2 ∂ y2 ) w− N0 D w =− M0 D(1−ν) , (2.28) đây chính là phương trình cơ bản xác định độ võng w của tấm. Dễ dàng thấy rằng điều kiện biên w = 0 được thỏa mãn, nếu w biểu thị qua chuỗi Fourier [1]: w = ∞ ∑ m=1 ∞ ∑ n=1 wmn sin mpix a sin npiyb , (2.29) trong đó, m,n là các số tự nhiên. Ta cũng biểu diễn M0 dưới dạng chuỗi Fourier: M0 = ∞ ∑ m=1 ∞ ∑ n=1 amn sin mpix a sin npiyb , (2.30) trong đó, amn =   16M0 pi2mn , m,n = 1,3,5, ... 0 m,n = 2,4,6, ... (2.31) Thay (2.29) và (2.30) vào phương trình (2.28), và do phương trình đúng với mọi x,y, và mặt giữa khi uốn phải đối xứng, các số hạng m,n chẵn tương ứng với độ võng không đối xứng do đó chúng bằng 0, nên suy ra: wmn =   amn D(1−ν) . 1 m2pi2 a2 + n 2pi2 b2 + N0 D m,n = 1,3,5, ... 0 m,n = 2,4,6, ... 25 Vậy, ta được nghiệm giải tích của bài toán xác định độ uốn của tấm khi có truyền nhiệt dừng là: w = ∑ m=1,3,5,... ∑ n=1,3,5,... amn D(1−ν) . 1 m2pi2 a2 + n 2pi2 b2 + N0 D sin mpix a sin npiyb . (2.32) trong đó, amn = 16M0 pi2mn , D = Eh3 12(1−ν2) , N0 =− αE1−ν h/2∫ −h/2 ∆T dz =− αE 1−ν N ∗ T , M0 = αE h/2∫ −h/2 z∆T dz = αEM∗T . (2.33) Do tấm đang xét là tấm composite có độn các hạt hình cầu với các hằng số môdun kéo nén thể tích K và môdun trượt G được xác định theo công thức (2.1) và hệ số dãn nở nhiệt theo công thức (2.2), nên ta có thể biểu diễn nghiệm (2.32) theo các hằng số K,G. Ta có [1] E = 9KG 3K +G , ν = 3K−2G 6K +2G , thay vào biểu thức (2.33), ta được N0 =− 18KGα3K +4GN ∗ T , M0 = 9KGα 3K +GM ∗ T D = h3G(3K +G) 3(3K +4G) , amn = 144KGαM∗T pi2mn(3K +G) (2.34) Thế (2.34) vào (2.32), ta thu được dạng khác của nghiệm xác định độ võng của tấm w = ∑ m=1,3,5,... ∑ n=1,3,5,... wmn sin mpix a sin npiyb . (2.35) trong đó, wmn = 864KαM∗T pi2mnh3 (3K +G) . 1 w1 , w1 = m2pi2 a2 + n2pi2 b2 − 54KαN∗T h3 (3K +G) 26 2.4. Tính toán số Với biểu thức nghiệm giải tích xác định độ võng của tấm khi có truyền nhiệt dừng (2.32), ta sẽ đi tính toán độ võng của tấm trong trường hợp cụ thể. Xét tấm composite độn các hạt Titan được làm từ các vật liệu thành phần với các đặc trưng tương ứng như sau [3, 4]: Nền PVC: Em = 3.109(Pa),νm = 0.2,αm = 8.10−5/K,km = 0.16(W/m.K). Hạt Titan: Ec = 100.109(Pa),νc = 0.34,αc = 4.8.10−6/K,kc = 22.1(W/m.K). Giả sử, tấm có chiều dài a = 2.25m, chiều rộng b = 1.5m, bề dày h = 0.02m. Môi trường xung quanh tấm là không khí, với hệ số truyền nhiệt đối lưu, nhiệt độ môi trường ở mặt (z = −h/2 ) tương ứng là β1 = 60(W/m2.K), T1 = 3300K và ở mặt (z = h / 2) tương ứng là β2 = 40(W/m2.K), T2 = 3000K. T0 = 2930K. Với các số liệu trên, thay vào biểu thức (2.11), sử dụng phần mềmMatlab 7.1 (xem phụ lục 1), kết quả thu được là sự phân bố nhiệt độ trong tấm theo bề dày của tấm như hình 2.1. Tương tự, thay vào biểu thức nghiệm giải tích (2.35), tác giả đã thu được độ võng của tấm như hình 2.2. (Xem phụ lục 2). Chú ý biểu thức nghiệm xác định độ võng của tấm (2.35) là một chuỗi kép vô hạn, nhưng chỉ cần vài số hạng đầu đã cho ta kết quả tương đối chính xác, vì vậy ở đây luận văn chỉ xét với m,n = 1−3. Hình 2.1 Đồ thị phân bố nhiệt độ theo chiều dày tấm. 27 Bảng 2.1 Độ võng của tấm tại một số điểm a a a a a a a a a w (x,y) (0, 0) (0.36, 0.24) (0.9, 0.6) (1.125, 0.75) (1.8, 1.2) (2.16, 1.44) w|ξ=0.1 0 0.0028 0.0008 0.0019 0.0027 0.0003 w|ξ=0.15 0 -0.0003 -0.0003 -0.0011 -0.0001 -0.0001 w|ξ=0.2 0 0.00027 -0.00022 -0.0007 0.0004 -0.00001 w|ξ=0.3 0 -0.0022 0.0028 0.0048 -0.0025 -0.0002 Nhận xét: Nhìn vào hình 2.1, ta thấy sự truyền nhiệt theo bề dày của tấm từ mặt trên (z =−h/2 ) xuống mặt dưới (z = h/2 ) là tuyến tính và khi tỉ lệ thể tích hạt Titan tăng (ξ tăng) thì sự chênh lệch nhiệt độ giữa mặt trên và mặt dưới giảm. Từ kết quả trên hình biểu thị độ võng của tấm dọc theo chiều dài tấm với những tỉ lệ thể tích hạt Titan nhất định dưới ảnh hưởng của sự truyền nhiệt dừng theo bề dày của tấm (Hình 2.2), ta thấy rằng khi trộn hạt Titan với những tỉ lệ thể tích khác nhau, các hạt Titan đã làm ảnh hưởng đáng kể đến ứng xử uốn của tấm. Khi ta trộn hạt với tỉ lệ thể tích dưới hoặc bằng 20% (ξ ≤ 0.2) thì độ võng tại các điểm của tấm giảm và 28 giữa các điểm là đồng đều hơn khi tăng tỉ lệ thể tích hạt. Vậy, với tỉ lệ thể tích tăng dần (≤ 20% ), các hạt Titan đã làm tăng khả năng chống uốn, chống dạn nứt và kháng nhiệt của tấm. Nhưng khi tỉ lệ thể tích hạt khoảng 30% (ξ = 0.3) thì ứng xử uốn của tấm biến đổi mạnh, độ võng tại các điểm của tấm cao và giữa các điểm có sự chênh lệch lớn. Như vậy, khi được trộn với một tỉ lệ thể tích hạt hợp lí, các hạt Titan sẽ có vai trò quan trọng trong việc tăng khả năng chống uốn, chống dạn nứt và kháng nhiệt của tấm. 2.5. Kết luận Trong chương 2, luận văn đã thu được kết quả như sau: 1. Thiết lập bài toán uốn tấm khi có truyền nhiệt dừng. 2. Giải tìm được nghiệm phương trình truyền nhiệt dừng (2.10), từ đó chỉ ra rằng mặt giữa của tấm không có biến dạng nhiệt, xác định được các hằng số (2.34) và xác định được độ võng của tấm w. 3. Trên cơ sở nghiệm giải tích uốn tấm tìm được, tác giả tính toán số với các composite cốt hạt Titan nền PVC, do đó thấy được sự ảnh hưởng của các hạt tới ứng xử uốn của tấm. Khi được trộn với một tỉ lệ thể tích hạt hợp lí, các hạt Titan sẽ có vai trò quan trọng trong việc tăng khả năng chống uốn, chống dạn nứt và kháng nhiệt của tấm. 29 Chương 3 Uốn tấm composite mỏng khi có truyền nhiệt không dừng Xét tấm composite độn các hạt hình cầu với các hằng số đàn hồi được xác định theo công thức (2.1)và hệ số dãn nở nhiệt theo công thức (2.2). Giả thiết tấm hình chữ nhật có cạnh a,b, chiều dày h, tấm tựa bản lề tại các cạnh. Ta xét trường hợp có sự trao đổi nhiệt đối lưu trên bề mặt z = ±h/2, và có sự truyền nhiệt không dừng theo bề dày của tấm. Bài toán đặt ra là hãy xác định uốn của tấm khi có truyền nhiệt không dừng. 3.1. Uốn tấm composite mỏng khi có truyền nhiệt không dừng 3.1.1. Sự phân bố nhiệt độ trong tấm Do có truyền nhiệt không dừng và gradient thông qua bề dày của tấm nên phương trình truyền nhiệt (1.20) là [9]: ∂ T ∂ t = a1 ∂ 2T ∂ z2 , (3.1) Với điều kiện đầu và các điều kiện biên:  T = Ti tại t = 0 k∂ T∂ z = β1(T −T1) tại z =− h 2 −k∂ T∂ z = β2(T −T2) tại z = h 2 (3.2) trong đó, a1 = kCρ là độ khuếch tán nhiệt của vật liệu, T1,β1 tương ứng là nhiệt độ môi trường, hệ số truyền nhiệt đối lưu tại mặtz = −h/2; T2,β2 là nhiệt độ môi trường, hệ số truyền nhiệt đối lưu tại mặt z = h / 2 tương ứng; k tương ứng là hệ số dẫn nhiệt và được giả thiết xác định theo công thức (2.5), C,ρ nhiệt dung riêng, mật độ khối của vật liệu, và giả thiết được xác định bởi các công thức sau [3]: ρ = (1−ξ )ρm +ξ ρc , C = Cmρm (1−ξ )+Ccρcξρm (1−ξ )+ρcξ . Đặt : τ = a1t h2 , ζ = z h , γ1 = β1h k , γ2 = β2h k . (3.3) thì phương trình (3.1) trở thành: ∂ T ∂ τ = ∂ 2T ∂ ζ 2 , (3.4) với điều kiện đầu và các điều kiện biên (3.2) xác định như sau:  T = Ti tại τ = 0 , ∂ T ∂ ζ = γ1 (T −T1) tại ζ =− 1 2 , −∂ T∂ ζ = γ2 (T −T2) tại ζ = 1 2 , (3.5) Dùng phép biến đổi Laplace: T ∗ = ∫ ∞ 0 Te−pτdτ , suy ra ∞∫ 0 ∂ T ∂ τ e −pτdτ = pT ∗−Ti khi đó, phương trình (3.4) và các điều kiện (3.5) tương ứng có dạng như (3.6) và (3.7): d2T ∗ dζ 2 − p ( T ∗− Ti p ) = 0 , (3.6)   dT ∗ dζ − γ1 ( T ∗− T1 p ) = 0 tại ζ =−1 2 , dT ∗ dζ + γ2 ( T ∗− T2 p ) = 0 tại ζ = 1 2 , (3.7) 31 Nghiệm của phương trình (3.6) có dạng: T ∗ = Ti p +C3ch √ pζ +C4sh√pζ , (3.8) trong đó, C3,C4 được xác định từ điều kiện biên (3.7). Thay (3.8) vào điều kiện biên (3.7), ta được hệ phương trình hai ẩn C3,C4 là  −C3 (√ psh √p 2 + γ1ch √p 2 ) +C4 (√ pch √p 2 + γ1sh √p 2 ) = γ1 (Ti−T1) p C3 (√ psh √p 2 + γ2ch √p 2 ) +C4 (√ pch √p 2 + γ2sh √p 2 ) =−γ2 (Ti−T2) p Giải hệ phương trình này, ta thu được C3 =− γ1 (Ti−T1) (√pch√p2 + γ2sh√p2 )+ γ2 (Ti−T2)(√pch√p2 + γ1sh√p2 ) p√p [ (p+ γ1γ2) sh √p√p +(γ1 + γ2)ch √p ] , C4 = γ1 (Ti−T1) (√psh√p2 + γ2ch√p2 )− γ2 (Ti−T2)(√psh√p2 + γ1ch√p2 ) p√p [ (p+ γ1γ2) sh √p√p +(γ1 + γ2)ch √p ] . (3.9) Thế kết quả C3,C4 từ (3.9) vào (3.8), với chú ý (3.12), ta dẫn đến T ∗ = Ti p + A(p) B(p) (3.10) trong đó, A(p) = γ1 (T1−Ti) [ ch√p ( 1 2 −ζ ) + γ2√psh √ p ( 1 2 −ζ )] + +γ2 (T2−Ti) [ ch√p ( 1 2 +ζ ) + γ1√psh √ p ( 1 2 +ζ )] B(p) = p [ (p+ γ1γ2) sh√p√p +(γ1 + γ2)ch √ p ] (3.11) 32 Chú ý shx = e x− e−x 2 , chx = e x + e−x 2 , shxchy = 1 2 [sh(x+ y)+ sh(x− y)] , chxshy = 1 2 [sh(x+ y)− sh(x− y)] , shxshy = 1 2 [ch(x+ y)− ch(x− y)] , chxchy = 1 2 [ch(x+ y)+ ch(x− y)] . (3.12) Ở đây, muốn đưa nghiệm T ∗ ở (3.10) về nghiệm ban đầu T , ta vận dụng định lí sau [9]: Ta có phép biến đổi Laplace của hàm f (s) là: f ∗(p) = ∫ ∞ 0 f (s)e−psds , Khi đó, định lí phát biểu như sau: Nếu phép biến đổi f ∗(p) là tỉ lệ của các hàm siêu việt Φ(p) và ψ(p), tức là f ∗(p) = Φ(p)ψ(p) , thì hàm ban đầu là f (s) = m ∑ n=1 1 (kn−1)! limp→pn dkn−1 d pkn−1 [ f ∗ (p)(p− pn)kn eps ] (3.13) trong đó, pn là nghiệm của ψ(p). Nếu ψ(p) chỉ có các nghiệm đơn pn (n = 1,2, ...), thì công thức (3.13) có thể đưa về dạng f (s) = m ∑ n=1 Φ(pn) ψ ′ (pn) epns , (3.14) trong đó ψ ′ (pn) = ∂ ψ(pn) ∂ p . Giải phương trình B(p) = 0 tức là p [ (p+ γ1γ2) sh√p√p +(γ1 + γ2)ch √ p ] = 0 phương trình trên có nghiệm p1 = 0 và pn+1 = −µ2n (n = 1,2, ...) trong đó, pn+1 là nghiệm của phương trình (p+ γ1γ2) sh√p√p +(γ1 + γ2)ch √ p = 0 33 khi đó, µn = i √p n+1 được xác định bởi phương trình sau tgµ = γ1 + γ2µ2− γ1γ2 µ . (3.15) với chú ý sh(iµ) = isin µ , ch(iµ) = cos µ Áp dụng định lí trên, ta có nghiệm T ban đầu có dạng T = Ti + lim p→0 d d p [ A(p) B(p) p2epτ ] + ∞ ∑ n=1 A (−µ2n) B′ (−µ2n ) epn+1τ Thực hiện tính toán, ta được lim p→0 d d p [ A(p) B(p) p2epτ ] = A(0) B′ (0) . trong đó A(0) = γ1(T1−Ti) [ 1+ γ2 ( 1 2 −ζ )] + γ2(T2−Ti) [ 1+ γ1 ( 1 2 +ζ )] , B ′ (0) = γ1 + γ2 + γ1γ2 , (3.16) Từ biểu thức B(p) trong (3.11) thực hiện đạo hàm theo p, ta có B ′ (p) = (p+ γ1γ2) sh√p√p +(γ1 + γ2)ch √ p+ + 1 2√p {[(1+ γ1 + γ2) p− γ1γ2]sh √ p+ √ p(p+ γ1γ2)ch √ p} khi đó, với pn+1 =−µ2n và chú ý sh(iµn) =−sin µn,ch(iµn) = cos µn, suy ra B ′ (−µ2n) =− 1 2µn {[ (1+ γ1 + γ2)µ2n + γ1γ2 ] sin µn + ( µ2n − γ1γ2 ) µn cos µn } . tương tự, thế pn+1 =−µ2n vào biểu thức A(p) ở (3.11) suy ra A(−µ2n) = γ1(T1−Ti) [ cos µn ( 1 2 −ζ ) + γ2 µn sin µn ( 1 2 −ζ )] + + γ2(T2−Ti) [ cos µn ( 1 2 +ζ ) + γ1 µn sin µn ( 1 2 +ζ )] , (3.17) 34 Kết quả thu được nghiệm T như sau T = Ti+ A(0) B′(0) + ∞ ∑ n=1 A(−µ2n) B′(−µ2n) e−µ 2 n τ , (3.18) trong đó, A(0) = γ1(T1−Ti) [ 1+ γ2 ( 1 2 −ζ )] + γ2(T2−Ti) [ 1+ γ1 ( 1 2 +ζ )] , B ′ (0) = γ1 + γ2 + γ1γ2 , A(−µ2n) = γ1(T1−Ti) [ cos µn ( 1 2 −ζ ) + γ2 µn sin µn ( 1 2 −ζ )] + (3.19) + γ2(T2−Ti) [ cos µn ( 1 2 +ζ ) + γ1 µn sin µn ( 1 2 +ζ )] , B ′ (−µ2n) =− 1 2µn {[ (1+ γ1 + γ2)µ2n + γ1γ2 ] sin µn + ( µ2n − γ1γ2 ) µn cos µn } . Lấy T0 = Ti, khi đó, thay (3.18) vào biểu thức (1.11) ta được ∆T∗ = A(0) B′(0) + ∞ ∑ n=1 A(−µ2n) B′(−µ2n ) e−µ 2 n τ . (3.20) trong đó, A(0),B′(0),A(−µ2n),B ′ (−µ2n ) xác định như (3.19). Thay biểu thức của τ và ζ từ (3.3) vào (3.20), khi đó ta thu được biểu thức ∆T∗ biểu diễn thông qua z,t như sau ∆T∗ = A1(0) B′(0) + ∞ ∑ n=1 A1(−µ2n) B′(−µ2n) e−ηµ 2 n t , (3.21) trong đó, η = a1h2 A1(0) = γ1(T1−Ti) [ 1+ γ2 ( 1 2 − z 2 )] + γ2(T2−Ti) [ 1+ γ1 ( 1 2 + z 2 )] , B ′ (0) = γ1 + γ2 + γ1γ2 , A1(−µ2n) = γ1(T1−Ti) [ cos µn ( 1 2 − z 2 ) + γ2 µn sin µn ( 1 2 − z 2 )] + (3.22) + γ2(T2−Ti) [ cos µn ( 1 2 + z 2 ) + γ1 µn sin µn ( 1 2 + z 2 )] , B ′ (−µ2n) =− 1 2µn {[ (1+ γ1 + γ2)µ2n + γ1γ2 ] sin µn + ( µ2n − γ1γ2 ) µn cos µn } . 35 3.1.2. Biểu thức nghiệm xác định uốn tấm Với biểu thức xây dựng được ∆T∗ theo (3.21), tương tự như phần 2.2 chương 2, ta cũng chỉ ra được rằng đối với tấm mỏng mà nhiệt độ chỉ truyền theo bề dày của tấm thì không gây ra biến dạng mặt giữa. Khi đó, thực hiện cách làm tương tự theo phần 2.3 chương 2, ta cũng thu được biểu thức nghiệm giải tích xác định uốn của tấm khi có truyền nhiệt không dừng w = ∑ m=1,3,5,... ∑ n=1,3,5,... wmn sin mpix a sin npiyb . (3.23) trong đó, wmn = 864KαMtT pi2mnh3 (3K +G) . 1 w1 , w1 = m2pi2 a2 + n2pi2 b2 − 54KαNtT h3 (3K +G) NtT = h/2∫ −h/2 ∆T∗dz , MtT = h/2∫ −h/2 z∆T∗dz 3.2. Tính toán số Để nghiên cứu cho trường hợp cụ thể, ta xét tấm composite độn các hạt Titan được làm từ các vật liệu thành phần với các đặc trưng tương ứng như sau [3, 4]: Nền PVC: Em = 3.109(Pa),νm = 0.2,αm = 8.10−5/K,km = 0.16(W/m.K), Cm = 900(J/kg.K),ρm = 1380(kg/m3). Hạt Titan: Ec = 100.109(Pa),νc = 0.34,αc = 4.8.10−6/K,kc = 22.1(W/m.K), Cc = 523(J/kg.K),ρc = 4500(kg/m3). Giả sử, tấm có chiều dài a = 2.25m, chiều rộng b = 1.5m, bề dày h = 0.02m. Môi trường xung quanh tấm là không khí, với hệ số truyền nhiệt đối lưu, nhiệt độ môi trường ở mặt (z = −h/2) tương ứng là β1 = 60(W/m2.K), T1 = 3300K và ở mặt (z = h / 2) tương ứng là β2 = 40(W/m2.K), T2 = 3000K. T0 = 2930K. Sử dụng phần mềm Matlab 7.1, lập trình ta giải phương trình siêu việt (3.15) ứng với mỗi giá trị cụ thể của tỉ lệ thể tích hạt độn ξ (ξ = 0.1,ξ = 0.15,ξ = 0.2, ξ = 0.3 ), ta tìm được nghiệm µn (n = 1,2, ...) tương ứng (xem phụ lục 3). 36 Sau khi tìm được nghiệm µn ứng với mỗi giá trị cụ thể ξ đó, thay vào biểu thức ∆T∗ theo (3.21), tiếp tục ứng dụng phần mềm Matlab (xem phụ lục 4), ta thu được sự phân bố nhiệt độ trong tấm theo bề dày tại những thời gian cố định xem hình 3.1-3.4. Khi đã xác định được sự phân bố nhiệt độ trong tấm, nhờ biểu thức xác định độ võng của tấm (3.23), áp dụng tính toán (xem phụ lục 5, 6), tác giả thu được độ võng của tấm composite cốt hạt Titan tại thời gian cố định, và sự biến đổi độ võng của tấm theo thời gian với những tỉ lệ trộn hạt Titan nhất định xem hình 3.5, 3.6. Tương tự như phần 2.4, ở đây luận văn chỉ xét với m,n = 1−3 vì chỉ cần vài số hạng đầu của biểu thức nghiệm đã cho ta kết quả tương đối chính xác. Hình 3.1 Sự phân bố nhiệt độ theo bề dày tấm với ξ = 0.1 37 Hình 3.2 Sự phân bố nhiệt độ theo bề dày tấm với ξ = 0.15 Hình 3.3 Sự phân bố nhiệt độ theo bề dày tấm với ξ = 0.2 38 Hình 3.4 Sự phân bố nhiệt độ theo bề dày tấm với ξ = 0.3 Nhận xét: Nhìn vào đồ thị biểu diễn sự phần bố nhiệt độ trong tấm (Hình 3.1-3.4) với những tỉ lệ thể tích hạt Titan nhất định, ta thấy khi thời gian nhỏ sự truyền nhiệt trong tấm theo bề dày từ mặt trên xuống mặt dưới là phi tuyến, nhưng khi thời gian đủ lớn sự truyền nhiệt là tuyến tính. Và khi tăng tỉ lệ thể tích hạt Titan (ξ tăng) sự chênh lệch nhiệt độ giữa lớp trên và lớp dưới giảm. Bảng 3.1 Độ võng của tấm tại một số điểm với t=1200s a a a a a a a a a w (x,y) (0, 0) (0.36, 0.24) (0.9, 0.6) (1.125, 0.75) (1.8, 1.2) (2.16, 1.44) w|ξ=0.1 0 0.0028 0.0008 0.0019 0.0027 0.0003 w|ξ=0.15 0 -0.0003 -0.0003 -0.0011 -0.0001 -0.0001 w|ξ=0.2 0 0.00035 -0.00027 -0.00078 0.00048 -0.00002 w|ξ=0.3 0 -0.0008 0.0012 0.0019 -0.0009 -0.0001 39 Từ kết quả trên hình biểu thị độ võng của tấm dọc theo chiều dài tấm với những tỉ lệ thể tích hạt Titan nhất định dưới ảnh hưởng của sự truyền nhiệt không dừng theo bề dày của tấm (Hình 3.5), ta thấy rằng khi trộn hạt Titan với những tỉ lệ thể tích khác nhau, các hạt Titan đã làm ảnh hưởng đáng kể đến ứng xử uốn của