Các bài toán cơ bản kênh hở hình thang,tính toán theo phương pháp đối chiếu mặt cắt có lợi nhất về thủy lực và dòng chảy trong ống
Ta xét thấy: Q=f(n, i, b, h, m) Tính kênh đã biết. Bài toán 1 : khi có n, i, b, h, m ta cần tìm Q Ta tính những trị số A, C, R rồi thay vào (1-10) tìm được Q. Bài toán 2 : khi có n, Q, b, h ...
Ta xét thấy: Q=f(n, i, b, h, m)
Tính kênh đã biết.
Bài toán 1: khi có n, i, b, h, m ta cần tìm Q
Ta tính những trị số A, C, R rồi thay vào (1-10) tìm được Q.
Bài toán 2: khi có n, Q, b, h ta cần tìm i.
Ta tính những trị số A, C, R rồi thay vào (1-9) tìm được theo công thức:
i=Q2A2C2R size 12{i= { {Q rSup { size 8{2} } } over {A rSup { size 8{2} } C rSup { size 8{2} } R} } } {} (1-29)
Bài toán 3: Khi có Q, i, b, h ta cần tìm n .
Thiết kế kênh mới.
Khi thiết kế kênh, cần tính chiều rộng và độ sâu mực nước kênh (b, h), cần thu thập các số liệu sau:
- Xác định độ dốc đáy kênh i, từ tuyến kênh theo bản đồ địa hình.
- Xác định hệ số nhám n và hệ số mái dốc m, căn cứ vào vật liệu lòng dẫn.
- Xác định lưu lượng Q, căn cứ vào nhu cầu sử dụng nước hay tiêu thoát nước được xác định ở các bài toán thủy nông, thủy văn công trình, cân bằng nước, v.v...
Sau khi xác định Q, m, n, i và chọn một trong các thông số, tùy từng trường hợp, thường gặp các bài toán có cách giải khác nhau như sau :
Bài toán 1 : Chọn β.
Từ công thức (1-10), tính theo Manning ta được:
Q=AnR23i size 12{Q= { {A} over {n} } R rSup { size 8{ { {2} over {3} } } } sqrt {i} } {}, (m3/s) (1-30)
Kết hợp các công thức(1-15), (1-17) và (1-8) thay vào ta tính được:
h=nQi38β+21+m20.25β+m58 size 12{h= left ( { { ital "nQ"} over { sqrt {i} } } right ) rSup { size 8{ { {3} over {8} } } } { { left (β+2 sqrt {1+m rSup { size 8{2} } } right ) rSup { size 8{0 "." "25"} } } over { left (β+m right ) rSup { size 8{ { {5} over {8} } } } } } } {}, (m) (1-31)
b=βh, (m) (1-31a)
Bài toán 2 : Chọn R hay v.
Từ (1-14) và (1-16), ta có:
(a)
(b)
Để giải bài toán, tìm nghiệm b và h từ hệ phương trình trên, cần xác định A và P
+ Nếu biết R, từ (1-28) ta tính :
A=nQR23i size 12{A= { { ital "nQ"} over {R rSup { size 8{ { {2} over {3} } } } sqrt {i} } } } {}, (m2) (1-32)
P=AR size 12{P= { {A} over {R} } } {}, (m) (1-33)
+ Nếu biết v, từ (1-9) theo Manning ta có:
v=R32ni size 12{v= { {R rSup { size 8{ { {3} over {2} } } } } over {n} } sqrt {i} } {}, (m/s) (1-34)
Nên: R=nvi32 size 12{R= left ( { { ital "nv"} over { sqrt {i} } } right ) rSup { size 8{ { {3} over {2} } } } } {} , (m) (1-35)
A=Qv size 12{A= { {Q} over {v} } } {}, (m2) (1-36)
Từ hệ phương trình, dùng phương pháp suy ra được như (1-26), sau đó khử h, ta được phương trình bậc hai:
m0h2 - Ph + A = 0 (1-37)
ở đó: mo=21+m2−m size 12{ ital "mo"=2 sqrt {1+m rSup { size 8{2} } } - m} {}
Giải phương (1-35) ta tìm được h.
h1,2=P±P2−4m0A2m0 size 12{h rSub { size 8{1,2} } = { {P +- sqrt {P rSup { size 8{2} } - 4m rSub { size 8{0} } A} } over {2m rSub { size 8{0} } } } } {} (1-38)
Từ h1 và h2 thay vào (1-26), ta chọn nghiệm dương, chiều rộng b và độ sâu mực nước hợp lý làm nghiệm.
Chú ý : Bài toán có nghiệm khi :
- Điều kiện của (1-38) là P2 > 4m0A
- Ngoài ra ta biết rằng khi mặt cắt có lợi nhất về thủy lực, thì bán kính thủy lực và vận tốc là lớn nhất và diện tích mặt cắt là nhỏ nhất. Như vậy bài toán chỉ có lời giải khi R và v cho trước nhỏ hơn R và v lợi nhất về thủy lực.
Bài toán 3 : Chọn b (hay h). Tính h (hay b)
Từ (1-12), ta tính (1-39)
Từ (1-11) ta cũng có thể truy tìm nghiệm bằng cách lập bảng hoặc bằng đồ thị. Dùng cách lập trình trong Visual basic, Pascal hay dùng phần mềm Mathcad để tính.
Bài toán có b tìm h hay có h tìm b, thường phải giải đúng dần, cho nên việc tính toán dùng máy tính tay gặp khó khăn về thời gian và mức độ chính xác phụ thuộc người tính. Vì vậy trong phần này giới thiệu phương pháp tính của Agơrôtskin bằng cách lập bảng tra đối với mặt cắt hình thang.
Agơrôtskin đặt hệ sốđặc trưng mặt cắt hình thang, không thứ nguyên, biểu thị quan hệ giữa b, h, m, nghĩa là biểu thị hình dạng mặt cắt.
Từ đó xác định các yếu tố thuỷ lực theo đặc trưng mặt cắt, điều quan trọng mặt cắt hình thang lợi nhất về thuỷ lực, có giá trị đặc trưng mặt cắt lợi nhất bằng một.
Từ đó xác định được bán kính lợi nhất thuỷ lực, đặc biệt quan hệ mặt cắt lợi nhất về thuỷ lực và mặt cắt bất kỳ là hàm số phụ thuộc vào đặc trưng mặt cắt.
Quan hệ hình dạng mặt cắt.
Từ (1-14), đặt bề rộng trung bình hình thang:
b¯=b+mh size 12{ {overline {b}} =b+ ital "mh"} {} (1-40)
nên: A=b¯h size 12{A= {overline {b}} h} {} (1-41)
Từ (1-40) rút b thay vào (1-16) xắp xếp lại ta được :
A=b¯+m0h size 12{A= {overline {b}} +m rSub { size 8{0} } h} {} (1-42)
ở đó đặt : m0=21+m2−m size 12{m rSub { size 8{0} } =2 sqrt {1+m rSup { size 8{2} } } - m} {} (1-43)
Tính bán kính thuỷ lực theo(1-41) và (1-42), ta được
R=b¯hb¯+m0h=h1+σ size 12{R= { { {overline {b}} h} over { {overline {b}} +m rSub { size 8{0} } h} } = { {h} over {1+σ} } } {} (1-44)
ở đó đặt: σ=m0hb¯ size 12{σ= { {m rSub { size 8{0} } h} over { {overline {b}} } } } {} (1-45)
Từ các công thức trên, nếu ta biết hệ số đặc trưng mặt cắt, thì quan hệ giữa các yếu tố của mặt cắt sẽ được xác định như sau:
Từ (1-44) rút h ta được :
h=(1+ σ)R (1-46)
Từ (1-45) rút chiều rộng trung bình và thay (1-46) vào, ta được:
b¯=m0hσ=m0σ1+σR size 12{ {overline {b}} = { {m rSub { size 8{0} } h} over {σ} } = { {m rSub { size 8{0} } } over {σ} } left (1+σ right )R} {} (1-47)
Từ (1-40) rút chiều rộng và thay (1-47) vào, ta được :
b=b¯−mh=m0σ−m1+σR size 12{b= {overline {b}} - ital "mh"= left ( { {m rSub { size 8{0} } } over {σ} } - m right ) left (1+σ right )R} {} (1-48)
Từ (1-41) thay (1-46) và (1-47) tính lại diện tích theo công thức :
A=1+σ2σm0R2 size 12{A= { { left (1+σ right ) rSup { size 8{2} } } over {σ} } m rSub { size 8{0} } R rSup { size 8{2} } } {} (1-49)
Suy ra R2=A.σm01+σ2 size 12{R rSup { size 8{2} } = { {A "." σ} over {m rSub { size 8{0} } left (1+σ right ) rSup { size 8{2} } } } } {} (1-50)
Từ (1-46) và (1-48) ta tình được hệ số:
β = m 0 σ − m size 12{β= { {m rSub { size 8{0} } } over {σ} } - m} {}
hay σ=m0β+m size 12{σ= { {m rSub { size 8{0} } } over {β+m} } } {} (1-51)
Đặc trưng của mặt cắt có lợi nhất về thủy lực.
Cũng như ở 1.3, xét mặt cắt lợi nhất, theo (1-50) ta biết rằng diện tích mặt cắt và mái dốc cho trước, nên mặt cắt lợi về thủy lực khi có R lớn nhất. Để R đạt gía trị lớn nhất ta xét đạo hàm sau :
d dσ σ 1 + σ 2 = 1 + σ 2 − 2σ 1 + σ 1 + σ 4 = 0 size 12{ { {d} over {dσ} } left [ { {σ} over { left (1+σ right ) rSup { size 8{2} } } } right ]= { { left (1+σ right ) rSup { size 8{2} } - 2σ left (1+σ right )} over { left (1+σ right ) rSup { size 8{4} } } } =0} {}
Tính đạo hàm và giải ra ta được σ=1. Vậy điều kiện để có mặt cắt lợi nhất về thủy lực của hình thang là khi :
σLn=1 (1-52)
Từ (1-51) cho bằng 1, và chú ý công thức (1-43), ta sẽ tìm được công thức (1-27). Điều này cho thấy mặt cắt lợi nhất thuỷ lực hình thang có thể biểu thịquan hệ khác nhau nhưng bản chất là như nhau.
Quan hệ giữa mặt cắt có lợi nhất về thủy lực và mặt cắt bất kỳ.
Xét phương trình cơ bản, ta có:
Ta tính hệ số C theo công thức (1-5) của Pavơlôpski; còn A tính theo (1-49) thay vào công thức trên, chuý thay σLN=1 ứng với mặt cắt lợi nhất. Sau đó, tính tỉ số bán kính bất kỳ trên mặt cắt lợi nhất về thuỷ lực và rút gọn ta được:
RRln=4σ1+σ21y+2.5=fσ size 12{ { {R} over {R rSub { size 8{"ln"} } } } = left [ { {4σ} over { left (1+σ right ) rSup { size 8{2} } } } right ] rSup { size 8{ { {1} over {y+2 "." 5} } } } =f left (σ right )} {} (1-53)
Nếu xem y là hằng số, ứng với σ cho trước, ta tính được công thức (1-52). Nếu chia hai vế công thức (1-46) và (1-48) cho RLn ta được:
hRln=1+σRRln=fσ size 12{ { {h} over {R rSub { size 8{"ln"} } } } = left (1+σ right ) { {R} over {R rSub { size 8{"ln"} } } } =f left (σ right )} {} (1-54)
bRln=m0σ−mhRln=fσ,m size 12{ { {b} over {R rSub { size 8{"ln"} } } } = left ( { {m rSub { size 8{0} } } over {σ} } - m right ) { {h} over {R rSub { size 8{"ln"} } } } =f left (σ,m right )} {} (1-55)
Theo Phoocơrâyme lấy y = 0.2, ta sẽ lập bảng các quan hệ giữa các đại lượng không thứ nguyên RRLn size 12{ { {R} over {R rSub { size 8{ ital "Ln"} } } } } {}, hRLn size 12{ { {h} over {R rSub { size 8{ ital "Ln"} } } } } {}, bRLn size 12{ { {b} over {R rSub { size 8{ ital "Ln"} } } } } {} theo σ﹐ từ (1-53), (1-54), (1-55) ở (Phụ lục 1-2). Bảng này tự chúng ta cũng có thể lập bảng trên excel.
Từ phụ lục, nếu biết một trong các đại lượng, tra ra các đại lượng còn lại. Do đó, có thể tính các kích thước hình thang như b, h, R nếu biết bán kính lợi nhất vế thuỷ lực.
Xác định bán kính thủy lực.
Theo lưu lượng cho mặt cắt lợi nhất về thủy lực, ta có:
Q = ωC R Ln i = 1 + σ Ln 2 σ Ln m 0 R Ln 2 C R Ln i size 12{Q= left (ωC sqrt {R} right ) rSub { size 8{ ital "Ln"} } sqrt {i} = { { left (1+σ rSub { size 8{ ital "Ln"} } right ) rSup { size 8{2} } } over {σ rSub { size 8{ ital "Ln"} } } } m rSub { size 8{0} } R rSub { size 8{ ital "Ln"} } rSup { size 8{2} } C sqrt {R rSub { size 8{ ital "Ln"} } } sqrt {i} } {}
⇔ Q = 4m 0 R Ln 2 . 5 C Ln i size 12{ dlrarrow Q=4m rSub { size 8{0} } R rSub { size 8{ ital "Ln"} } rSup { size 8{2 "." 5} } C rSub { size 8{ ital "Ln"} } sqrt {i} } {}
⇔ 4m 0 i Q = 1 CR 2 . 5 Ln = f R Ln size 12{ dlrarrow { {4m rSub { size 8{0} } sqrt {i} } over {Q} } = left ( { {1} over { ital "CR" rSup { size 8{2 "." 5} } } } right ) rSub { size 8{ ital "Ln"} } =f left (R rSub { size 8{ ital "Ln"} } right )} {}
Agơrôtskin đã tính sẵn quan hệ:
fRln=4m0iQ size 12{f left (R rSub { size 8{"ln"} } right )= { {4m rSub { size 8{0} } sqrt {i} } over {Q} } } {} (1-56)
Trong đó hệ số Chezy được tính theo công thức của tác giả và lập thành bảng (Phụ lục 1 -1)
Nếu tính C theo công thức của Maninh hay Phoocơrâyme, thì có thể tính rút trực tiếp ra RLn:
- Theo Maninh: Rln=nQ4m0i38 size 12{R rSub { size 8{"ln"} } = left ( { { ital "nQ"} over {4m rSub { size 8{0} } sqrt {i} } } right ) rSup { size 8{ { {3} over {8} } } } } {} (1-57)
- Theo Phoocơrâyme: Rln=nQ4m0i38 size 12{R rSub { size 8{"ln"} } = left ( { { ital "nQ"} over {4m rSub { size 8{0} } sqrt {i} } } right ) rSup { size 8{ { {3} over {8} } } } } {} (1-58)
Cách vận dụng cụ thể
Bài toán 1: Tìm h khi biết: Q, m, n, i và b.
+ Trước tiên xác định bán kính lợi nhất về thuỷ lực: RLn có thể dùng các công thức (1-57), (1-58) hoặc dùng phụ lục (1-1).
+ Lập tỉ: bRLn size 12{ { {b} over {R rSub { size 8{ ital "Ln"} } } } } {} tra phụ lục (1-2) suy ra được: hRLn size 12{ { {h} over {R rSub { size 8{ ital "Ln"} } } } } {}
+ Tính h theo công thức:
h=hRLnRLn size 12{h= { {h} over {R rSub { size 8{ ital "Ln"} } } } R rSub { size 8{ ital "Ln"} } } {} (1-59)
Bài toán 2: Tìm b khi biết: Q, m, n, i và h.
+ Trước tiên xác định RLn như trên
+ Lập tỉ: hRLn size 12{ { {h} over {R rSub { size 8{ ital "Ln"} } } } } {} tra phụ lục (1-2) suy ra được: bRLn size 12{ { {b} over {R rSub { size 8{ ital "Ln"} } } } } {}
+ Tính b theo công thức:
b=bRlnRln size 12{b= { {b} over {R rSub { size 8{"ln"} } } } R rSub { size 8{"ln"} } } {} (1-60)
Bài toán 3: Tìm b và h, khi biết: Q, m, n, i và β
+ Xác định RLn như trên.
+ Tính đặc trưng mặt cắt hình thang theo công thức (1-51), tra phụ lục (1-2) suy ra được hRLn size 12{ { {h} over {R rSub { size 8{ ital "Ln"} } } } } {}, bRLn size 12{ { {b} over {R rSub { size 8{ ital "Ln"} } } } } {}
+ Tính h và b theo công thức: (1-59) và (1-60)
Bài toán 4: Tìm b và h, khi biết: Q, m, n, i và R hoặc v.
+ Xác định RLn như trên.
+ Nếu có R thì lập tỉ số, tra phụ lục (8-3) suy ra được: hRLn size 12{ { {h} over {R rSub { size 8{ ital "Ln"} } } } } {}, bRLn size 12{ { {b} over {R rSub { size 8{ ital "Ln"} } } } } {}
+ Tính h và b theo công thức: (1-59) và (1-60)
- Nếu biết v: Tính vận tốc theo Chezy, hệ số Chezy xác định theo Manning. Do đó tính bán kính thuỷ lực R theo công thức (1-35), tính ra b và h như trên.
Các yếu tố thuỷ lực
Công thức tính diện tích và chu vi mặt cắt hình tròn chảy lưng ống, tuy đơn giản nhưng ít được các tài liệu chứng minh.
Tính diện tích, xét 2 phần: diện tích cung tròn MHG và diện tích tam giác OMN, tức là:
A = A MHG + A OMG = 1 8 2θ − sin 2θ d 2 size 12{A=A rSub { size 8{ ital "MHG"} } +A rSub { size 8{ ital "OMG"} } = { {1} over {8} } left (2θ - "sin"2θ right )d rSup { size 8{2} } } {}
trong đó:
d là đường kính mặt cắt hình tròn;
θ là góc được ghi chú trên hình 3. (rad)
Diện tích cung tròn MHG: AMGH=π4d22.θ2π=θ4d2 size 12{A rSub { size 8{ ital "MGH"} } = { {π} over {4} } d rSup { size 8{2} } { {2 "." θ} over {2π} } = { {θ} over {4} } d rSup { size 8{2} } } {}
Diện tích phần tam giác OMG: AOMG=2AOMN=ON.MN=−d24sinθcosθ size 12{A rSub { size 8{ ital "OMG"} } =2A rSub { size 8{ ital "OMN"} } = ital "ON" "." ital "MN"= - { {d rSup { size 8{2} } } over {4} } "sin"θ"cos"θ} {}
Vì xét tam giác vuông OMN, ta có:
ta lại có:
ON = h − d 2 size 12{ ital "ON"=h - { {d} over {2} } } {}
Do đó:
cos θ = 1 − 2 h d size 12{"cos"θ=1 - 2 { {h} over {d} } } {}
Hay:
cosθ =1-2a (1-61)
Đặt:
a=hd size 12{a= { {h} over {d} } } {} (1-61a)
Công thức (1-65) và (1-66), giúp chúng ta thiết lập mối quan hệ giữa độ sâu mực nước chảy lưng ống với đường kính ống tròn và góc θ đã đặt, để từ tính diện tích ướt và chu vi ướt.
Diện tích: A=kAd2 size 12{A=k rSub { size 8{A} } d rSup { size 8{2} } } {}(1-62)
Đặt: kA=182θ−sin2θ size 12{k rSub { size 8{A} } = { {1} over {8} } left (2θ - "sin"2θ right )} {}(1-62a)
Chu vi ướt P=θ.d size 12{P=θ "." d} {} (1-63)
Chiều rộng mặt thoáng B=dsinθ (1-64)
Bán kính thuỷ lực R=kAθd size 12{R= { {k rSub { size 8{A} } } over {θ} } d} {} (1-65)
Công thức tính lưu lượng
Tính lưu lượng theo công thức Manning (1-30), thay (1-62) và (1-65), ta được:
Q=kA53θ23ind83 size 12{Q= { {k rSub { size 8{A} } rSup { size 8{ { {5} over {3} } } } } over {θ rSup { size 8{ { {2} over {3} } } } } } { { sqrt {i} } over {n} } d rSup { size 8{ { {8} over {3} } } } } {} (1-66)
hθ=nQi.d83=kA53θ23 size 12{h left (θ right )= { { ital "nQ"} over { sqrt {i} "." d rSup { size 8{ { {8} over {3} } } } } } = { {k rSub { size 8{A} } rSup { size 8{ { {5} over {3} } } } } over {θ rSup { size 8{ { {2} over {3} } } } } } } {} (1-67)
Mặt cắt lợi nhất về thuỷ lực
Với i, n và d cho trước, ứng độ sâu mực nước trong ống là bao nhiêu để có lưu lượng lớn nhất khi:
d dθ k A 5 3 θ 2 3 = d dθ 2θ − sin 2θ 5 3 θ 2 3 = 0 size 12{ { {d} over {dθ} } left ( { {k rSub { size 8{A} } rSup { size 8{ { {5} over {3} } } } } over {θ rSup { size 8{ { {2} over {3} } } } } } right )= { {d} over {dθ} } left [ { { left (2θ - "sin"2θ right ) rSup { size 8{ { {5} over {3} } } } } over {θ rSup { size 8{ { {2} over {3} } } } } } right ]=0} {}
Sau khi lấy đạo hàm hàm số trên, ta được phương trình:
2θ- 5θcos2θ sin2θ =0
Giải phương trình, ta được: θ=1510 hay a=0,94.
Tính vận tốc theo (1-34), thay bán kính thuỷ lực (1-64), ta được:
v=inkAθ23d23 size 12{v= { { sqrt {i} } over {n} } left ( { {k rSub { size 8{A} } } over {θ} } right ) rSup { size 8{ { {2} over {3} } } } d rSup { size 8{ { {2} over {3} } } } } {} (1-68)
Với i, n và d cho trước, ứng độ sâu mực nước trong ống là bao nhiêu để có vận tốc lớn nhất khi:
d dθ k A θ 2 3 = d dθ 2 . θ − sin 2θ θ 2 3 = 0 size 12{ { {d} over {dθ} } left [ left ( { {k rSub { size 8{A} } } over {θ} } right ) rSup { size 8{ { {2} over {3} } } } right ]= { {d} over {dθ} } left [ left ( { {2 "." θ - "sin"2θ} over {θ} } right ) rSup { size 8{ { {2} over {3} } } } right ]=0} {}
Sau khi lấy đạo hàm hàm số trên, ta được phương trình:
- 2θcos2θ sin2θ =0
Giải phương trình, ta được: θ=1290 hay a=0,81
Các bài thường gặp
Bài toán 1: Bài toán thiết kế, có Q, n và i. Xác định đường kính ống.
Giải.
Từ công thức (1-66), cho thấy Q=f(n, i, d, a), vì vậy bài toán có 2 ẩn số là d và a, nhưng chỉ có một phương trình, nên tuỳ yêu cầu thực tế ta cần lưu lượng lớn thì lấy a=0,94, còn tính theo vân tốc lớn nhất lấy a=0,81.
Khi có a ta kính được θ và kA, tính theo công thức sau:
d=n.Qi38θ14kA58 size 12{d= left ( { {n "." Q} over { sqrt {i} } } right ) rSup { size 8{ { {3} over {8} } } } { {θ rSup { size 8{ { {1} over {4} } } } } over {k rSub { size 8{A} } rSup { size 8{ { {5} over {8} } } } } } } {} (1-69)
Bài toán 2: Bài toán kiểm tra, có Q, d, n và i. Xác định độ sau mực nước.
Giải.
Từ (1-67), ta tính được:
h0θ=n.Qi.d83 size 12{h rSub { size 8{0} } left (θ right )= { {n "." Q} over { sqrt {i} "." d rSup { size 8{ { {8} over {3} } } } } } } {} (1-70)
Có 2 cách để tìm nghiệm h:
- Cách 1: Phương pháp thử dần (mò nghiệm), tự chọn a tính θ và kA, ta vào biểu thức sau:
hθ=kA53θ23 size 12{h left (θ right )= { {k rSub { size 8{A} } rSup { size 8{ { {5} over {3} } } } } over {θ rSup { size 8{ { {2} over {3} } } } } } } {} (1-71)
Tính đến khi nào h0(θ)≈ h(θ) thì gía trị a đó cần tìm.
- Cách 2: Tra bảng, từ công thức (1-61), (1-62) và (1-71) ta lập bảng tra
Từ công thức (1-70) tính được h0(θ) dựa vào bảng ta tra ra giá tri cần tím a, tính h theo công thức sau:
h=a.d (1-72)
Từ các công thức (1-61a), (1-61), (1-62a) và (1-71), tiến hành lập bảng bằng excel Phụ lục 1-3 để tra, thuân tiện trong việc tính toán bằng máy tính tay. Ta cũng thể dựa vào các công thức trên lập trình tính toán hay dùng phần mềm Mathcad.