ห้องทดลองช่องไหล: คำตอบวิเคราะห์เทียบวิธีเชิงตัวเลข

เป้าหมายและขอบเขต

จบบทนี้ผู้เรียนจะตรวจหน่วย แก้ความเร็วระหว่างแผ่น และอธิบายได้ว่าความเร็วตรงที่จุดกริดยังไม่ได้ทำให้อัตราไหลจากการอินทิเกรตตรงทั้งหมด แบบฝึกเป็น ปัญหา steady laminar 1D ที่ลดรูปจาก Navier–Stokes ไม่ใช่ solver สนามไหลทั่วไปหรือ CFD ใบพัด

สมมติของไหล Newtonian ไม่อัดตัว ρ และ μ คงที่ แผ่นนิ่งห่าง H กว้างมากเมื่อเทียบกับช่องว่าง ไม่มีผลผนังด้านข้าง สนามพัฒนาเต็มที่ u = u(y) และไม่มีแรงปริมาตรตามแนวไหล พิจารณาบริเวณห่างจากทางเข้าออก ไม่จำลองช่วงพัฒนาของโปรไฟล์ เงื่อนไข no-slip คือ u(0) = u(H) = 0 ตรงกับโจทย์ MIT plane Poiseuille flow

จากแรงดันสู่ความเร็ว

นิยาม Δp = p_in − p_out ≥ 0 และ L เป็นระยะที่ใช้วัดความดันตก จึงได้ dp/dx = −Δp/L สมการ μd²u/dy² = −Δp/L อินทิเกรตสองครั้งแล้วใส่เงื่อนไขผนังจะได้

u(y) = Δp y(H − y)/(2μL)

u_max = ΔpH²/(8μL), Ū = ΔpH²/(12μL), Q′ = ŪH = ΔpH³/(12μL)

Q′ เป็นปริมาตรไหลต่อหน่วยความกว้าง หน่วย m²/s ไม่ใช่ m³/s หากความกว้างจริง W มากพอให้ละเลยขอบข้างได้จึงประมาณ Q ≈ WQ′ ส่วน u_max/Ū = 1.5 เป็นคุณสมบัติของโปรไฟล์นี้ ความหนาแน่นไม่ปรากฏในสูตรความเร็วเมื่อกำหนด Δp แต่ยังใช้ตรวจ Reynolds ได้

ตัวอย่างอ้างอิงแม้ปิด JavaScript

ใช้ ρ = 1.225 kg/m³, μ = 0.0000181 Pa·s, H = 0.01 m, L = 1 m และ Δp = 0.01 Pa

ค่า Re ต่ำนี้สอดคล้องกับการเลือกตัวอย่างที่ความหนืดเด่น แต่การเพิ่ม Δp จนสูงมากไม่ทำให้คำตอบพาราโบลาเป็นคำตอบทางกายภาพที่เชื่อได้เสมอ เครื่องคำนวณยังแสดงสูตรเดิม จึงต้องประเมินสมมติฐานใหม่เอง

วิธีเชิงตัวเลขที่แก้อย่างอิสระ

แบ่ง H เป็น N ช่วงเท่ากัน มี N + 1 จุด ใช้ h = H/N และจุดภายใน i = 1,…,N − 1 ประมาณอนุพันธ์อันดับสองด้วย central difference:

u_i−1 − 2u_i + u_i+1 = −Δp h²/(μL)

นำค่าผนังศูนย์เข้าสู่ระบบสมการสามแนวทแยง แล้วแก้ค่าภายในโดยไม่แทนสูตรพาราโบลาเป็นคำตอบเชิงตัวเลข จากนั้นเปรียบเทียบแต่ละจุดกับสูตรวิเคราะห์ เครื่องคำนวณและสมุดงานใช้วิธีนี้เพื่อให้ตรวจทั้งการตั้งระบบและการแก้ระบบได้

มีจุดน่าสังเกต: อนุพันธ์อันดับสองแบบ central difference ให้ค่าตรงสำหรับพหุนามกำลังสองบนกริดสม่ำเสมอ กรณีนี้จึงอาจมี error ที่จุดเพียงระดับปัดเศษ แม้กริดหยาบ ไม่ใช่หลักฐานว่า solver ทั่วไปแม่นเท่านี้ การหาปริมาตรไหลด้วยกฎคางหมูยังมี error โดย Q′_trap = Q′(1 − 1/N²) สำหรับโปรไฟล์และกริดนี้

ที่ N = 20 จะได้ Q′_trap ≈ 0.0000459254 m²/s การเพิ่ม N สองเท่าลด error ของการอินทิเกรตประมาณสี่เท่า นี่เป็นผลคำนวณของแบบฝึกนี้ ไม่ใช่ผลทดลอง NASA หรือ MIT และยังไม่ใช่ validation กับของไหลจริง ตามความแตกต่างที่ NASA อธิบายใน verification

ลงมือทดลองและส่งชิ้นงาน

1. เปิดเครื่องคำนวณท้ายบท คำนวณค่าเริ่มต้นและตรวจหน่วยตามตาราง
2. เพิ่ม Δp สองเท่า แล้วตรวจว่า u และ Q′ เพิ่มสองเท่า
3. เปลี่ยน N จาก 10 เป็น 20 และ 40 เก็บ CSV และเปรียบเทียบ error ของ Q′
4. ดาวน์โหลด สมุดงาน Python และไฟล์อ้างอิง อ่าน README แล้วรัน solver และ tests ตามคำสั่งในนั้น เก็บเวอร์ชัน Python พร้อมผลที่ได้

โจทย์: หากเพิ่ม H สองเท่า แต่ตรึง μ, L และ Δp ความเร็วเฉลี่ยกับ Q′ เปลี่ยนอย่างไร?

แนวคำตอบ: Ū เพิ่มสี่เท่า และ Q′ เพิ่มแปดเท่า ส่วน Re_Dh เพิ่มแปดเท่า จึงต้องกลับไปตรวจขอบเขต laminar และเงื่อนไขช่องกว้างอีกครั้ง ก่อนตีความว่าเป็นผลของช่องจริง

ต่อไป: จากแบบฝึกหนึ่งมิติสู่งาน CFD ที่ตรวจสอบได้

บรรณานุกรม

• MIT OCW — Couette & Poiseuille Flows
• NASA NPARC — Verification Assessment

ตรวจแหล่งต้นทาง 26 กันยายน 2569