Để lập mô hình con lắc đơn giản, chúng tôi cần ba gói. Numpy như np đến. matplotlib. pyplot dưới dạng plt để cho phép chúng tôi vẽ kết quả của mình và từ scipy. tích hợp nhập odeint để cho phép chúng tôi giải quyết các ODE
Một trong những điều thực sự khiến vật lý “nhấp chuột” đối với tôi là học cách mô phỏng số các hệ thống mà tôi đã học trên lớp. Viết ra các phương trình chuyển động cho một bài toán nhất định là một chuyện – thực sự hình dung chúng trông như thế nào trên máy tính lại là một chuyện khác. Đó là mối liên kết còn thiếu giữa trừu tượng và thực tế
Trong bài đăng này, tôi sẽ giải bài toán con lắc đơn – một bài toán đồ chơi phổ biến trong vật lý – và chỉ ra cách mô phỏng nó trên máy tính. Để hiểu các khái niệm trong bài viết này, bạn cần có kiến thức cơ bản về giải tích và làm quen với logic lập trình
Bước đầu tiên để mô phỏng bất kỳ loại bài toán vật lý nào là viết ra các phương trình giải tích của chuyển động. Chúng ta sẽ sử dụng phương pháp Lagrangian, một cách đặc biệt hiệu quả để phân tích các hệ thống phức tạp. Sử dụng nó cho trường hợp của một con lắc đơn giản có lẽ là quá mức cần thiết, nhưng tôi thích sự sang trọng trong cách tiếp cận của nó. Bước đầu tiên của chúng tôi là xác định vấn đề của chúng tôi. Một sơ đồ của một con lắc đơn giản được trình bày dưới đây
Tuyên bố vấn đề là sau đó. Tìm phương trình chuyển động của con lắc có chiều dài l và khối lượng m. Cách tiếp cận Lagrangian để làm điều này như sau
- Xác định tọa độ tổng quát của hệ thống, .
- Viết Lagrange, trong đólà động năng vàlà thế năng, tính theo . .
- Các phương trình chuyển động của hệ khi đó sẽ được cho bởi phương trình Euler-Lagrange, .
Sẽ luôn có nhiều tọa độ tổng quát tương ứng với số bậc tự do trong hệ thống của bạn. Cách tôi muốn nghĩ về tọa độ tổng quát là tưởng tượng tôi đang nói chuyện điện thoại với ai đó và tôi cần mô tả trạng thái hệ thống của mình cho họ. Tập hợp thông tin tối thiểu mà tôi có thể nói với họ để họ có thể sao chép hoàn toàn hệ thống của tôi là gì? . Một nhỏ hơn hai, và do đó tọa độ tổng quát cho con lắc đơn của chúng ta đơn giản là
Bước tiếp theo là viết năng lượng của chúng ta theo góc . Động năng được viết dễ dàng bằng cách ghi nhớ hệ thức
Đối với thế năng, chúng ta định nghĩa “không” của mình là khi khối lượng hoàn toàn thẳng đứng, i. e. khi
Bây giờ chúng ta có thể viết Lagrangian
và áp dụng phương trình Euler-Lagrange để có được phương trình chuyển động của chúng ta
Bây giờ chúng ta có phương trình chuyển động giải tích cuối cùng. Lưu ý rằng các số hạng khối lượng triệt tiêu, cho thấy chuyển động của con lắc không phụ thuộc vào khối lượng của nó. Ở giai đoạn này, nhiều khóa học vật lý nhập môn sẽ lấy xấp xỉ góc nhỏ
Để giải số một hệ, bước đầu tiên là rời rạc hóa nó. Cách thông thường để làm điều này là viết ra chuỗi Taylor cho một hàm liên tục và cắt ngắn nó ở một số hạng. Trong trường hợp của chúng tôi, chúng tôi muốn tìm một giá trị gần đúng của
Viết lại, ta có
và cắm cái này vào hệ phương trình vi phân của chúng tôi, chúng tôi nhận được
Tuy nhiên, có một câu hỏi – đối với vế phải, chúng ta nên tính giá trị
Bây giờ chúng ta có một phương pháp để xác định và
omega = omega_old - [g/l]*sin[theta_old]*tau theta = theta_old + omega_old*tau
trong đó chúng tôi đã xác định
import numpy as np from math import * def simplePendulumSimulation[theta0,omega0,tau,m,g,l, numSteps,plotFlag]: # initialize vectors time_vec = [0]*numSteps theta_vec = [0]*numSteps omega_vec = [0]*numSteps # set initial conditions theta = theta0 omega = omega0 time = 0 # begin time-stepping for i in range[0,numSteps]: omega_old = omega theta_old = theta # update the values omega = omega_old - [g/l]*sin[theta_old]*tau theta = theta_old + omega_old*tau # record the values time_vec[i] = tau*i theta_vec[i] = theta omega_vec[i] = omega
Chạy mã này cho các điều kiện ban đầu
Tổng năng lượng đang tăng lên không giới hạn. Nói cách khác, năng lượng của hệ không được bảo toàn. Bạn có thể thấy điều này bằng sự dịch chuyển ngày càng tăng của con lắc của chúng ta, nó dao động ngày càng cao một cách nghịch lý với mỗi vòng quay. Vấn đề nằm ở phép tính xấp xỉ bậc nhất thô sơ mà chúng tôi đã thực hiện trước đó khi rời rạc hóa giải pháp của mình. Rốt cuộc, phương trình chuyển động của chúng ta là phi tuyến tính, làm cho các xấp xỉ tuyến tính đặc biệt không hiệu quả trong việc giải nó. Thay vì sử dụng các phương pháp bậc cao, chúng ta chỉ cần áp dụng một biến đổi đơn giản gọi là phương pháp Euler-Cromer, đảm bảo tiết kiệm năng lượng. Chúng tôi chỉ cần sử dụng giá trị cập nhật của vận tốc góc khi nó có sẵn
omega = omega_old - [g/l]*sin[theta_old]*tau theta = theta_old + omega*tau
Bạn có thể thấy rằng ở dòng thứ hai đã được thay đổi đơn giản thành . Chạy mô phỏng lại ta được kết quả.
Lúc này con lắc trở về đúng độ dời ban đầu sau mỗi chu kỳ. Chúng ta thấy rằng mặc dù tổng năng lượng vẫn dao động một chút [do phép tính gần đúng bậc nhất của chúng ta], nhưng hiện tại nó đã được giới hạn đúng. Để giảm dao động năng lượng, chúng ta chỉ cần giảm bước thời gian. Tập lệnh Python để thực hiện mô phỏng bằng phương pháp Euler-Cromer và tạo ảnh được đưa ra tại đây. quả lắc. py