[已解决]求解使用python与matlab解二阶非线性微分方程组

臧雁德

给出初始设置
m1 = 1.0
m = 2.0
l = 1.0
g = 10
F = 1.0
Dx = 0.1
Dtheta = 0.1

c = np.cos(theta)
s = np.sin(theta)

求解以下二阶非线性微分方程组
x_ddot * (m+m1) + Dx*x_dot - m*l*theta_dot*s + m*l*theta_ddot*c = F
theta_ddot * (4/3*m*(l**2)) - m*g*l*s + Dtheta*theta_dot + m*l*x_dot*c

最佳答案

月排行榜 / 总排行榜

isdkz

2023-9-23 22:53:49

对于这个问题，我们可以使用Python的SciPy库中的odeint函数来解决。首先，需要将二阶微分方程组转化为一阶微分方程组。我们可以定义状态向量x和时间t，其中x包括x，x_dot，theta，和theta_dot。然后，我们可以定义一个函数来计算状态向量的导数。最后，使用odeint函数来解决微分方程组。

以下是使用Python解决二阶非线性微分方程组的代码：

import numpy as np
from scipy.integrate import odeint

m1 = 1.0
m = 2.0
l = 1.0
g = 10
F = 1.0
Dx = 0.1
Dtheta = 0.1

def derivs(x, t):
    x_dot = x[1]
    theta = x[2]
    theta_dot = x[3]
    
    c = np.cos(theta)
    s = np.sin(theta)
    
    x_ddot = (F - Dx*x_dot + m*l*theta_dot*s - m*l*theta_ddot*c) / (m + m1)
    theta_ddot = (m*g*l*s - Dtheta*theta_dot - m*l*x_dot*c) / (4/3*m*(l**2))
    
    return [x_dot, x_ddot, theta_dot, theta_ddot]

# 初始设置
x0 = [0.0, 0.0, 0.1, 0.0]  # x, x_dot, theta, theta_dot
t = np.linspace(0, 10, 1000)

# 解决微分方程组
x = odeint(derivs, x0, t)

# 绘制结果
import matplotlib.pyplot as plt
plt.plot(t, x[:, 0], label='x')
plt.plot(t, x[:, 2], label='theta')
plt.legend()
plt.show()

对于Matlab，我们可以使用ode45函数来解决微分方程组。以下是使用Matlab解决二阶非线性微分方程组的代码：

m1 = 1.0;
m = 2.0;
l = 1.0;
g = 10;
F = 1.0;
Dx = 0.1;
Dtheta = 0.1;

function x_dot = derivs(t, x)
    x_dot = zeros(4, 1);
    x_dot(1) = x(2);
    theta = x(3);
    theta_dot = x(4);
    
    c = cos(theta);
    s = sin(theta);
    
    x_ddot = (F - Dx*x(2) + m*l*theta_dot*s - m*l*theta_ddot*c) / (m + m1);
    theta_ddot = (m*g*l*s - Dtheta*theta_dot - m*l*x(2)*c) / (4/3*m*(l^2));
    
    x_dot(2) = x_ddot;
    x_dot(3) = theta_dot;
    x_dot(4) = theta_ddot;
end

% 初始设置
x0 = [0.0, 0.0, 0.1, 0.0];  % x, x_dot, theta, theta_dot
tspan = [0, 10];

% 解决微分方程组
[t, x] = ode45(@derivs, tspan, x0);

% 绘制结果
plot(t, x(:, 1), 'r', t, x(:, 3), 'b');
legend('x', 'theta');

这两个代码都可以解决这个问题，但是Matlab的代码可能会更快一些，因为它使用了更高效的ode45函数。

跳转到最佳答案楼层

isdkz

对于这个问题，我们可以使用Python的SciPy库中的odeint函数来解决。首先，需要将二阶微分方程组转化为一阶微分方程组。我们可以定义状态向量x和时间t，其中x包括x，x_dot，theta，和theta_dot。然后，我们可以定义一个函数来计算状态向量的导数。最后，使用odeint函数来解决微分方程组。

以下是使用Python解决二阶非线性微分方程组的代码：

import numpy as np
from scipy.integrate import odeint

m1 = 1.0
m = 2.0
l = 1.0
g = 10
F = 1.0
Dx = 0.1
Dtheta = 0.1

def derivs(x, t):
    x_dot = x[1]
    theta = x[2]
    theta_dot = x[3]
    
    c = np.cos(theta)
    s = np.sin(theta)
    
    x_ddot = (F - Dx*x_dot + m*l*theta_dot*s - m*l*theta_ddot*c) / (m + m1)
    theta_ddot = (m*g*l*s - Dtheta*theta_dot - m*l*x_dot*c) / (4/3*m*(l**2))
    
    return [x_dot, x_ddot, theta_dot, theta_ddot]

# 初始设置
x0 = [0.0, 0.0, 0.1, 0.0]  # x, x_dot, theta, theta_dot
t = np.linspace(0, 10, 1000)

# 解决微分方程组
x = odeint(derivs, x0, t)

# 绘制结果
import matplotlib.pyplot as plt
plt.plot(t, x[:, 0], label='x')
plt.plot(t, x[:, 2], label='theta')
plt.legend()
plt.show()

对于Matlab，我们可以使用ode45函数来解决微分方程组。以下是使用Matlab解决二阶非线性微分方程组的代码：

m1 = 1.0;
m = 2.0;
l = 1.0;
g = 10;
F = 1.0;
Dx = 0.1;
Dtheta = 0.1;

function x_dot = derivs(t, x)
    x_dot = zeros(4, 1);
    x_dot(1) = x(2);
    theta = x(3);
    theta_dot = x(4);
    
    c = cos(theta);
    s = sin(theta);
    
    x_ddot = (F - Dx*x(2) + m*l*theta_dot*s - m*l*theta_ddot*c) / (m + m1);
    theta_ddot = (m*g*l*s - Dtheta*theta_dot - m*l*x(2)*c) / (4/3*m*(l^2));
    
    x_dot(2) = x_ddot;
    x_dot(3) = theta_dot;
    x_dot(4) = theta_ddot;
end

% 初始设置
x0 = [0.0, 0.0, 0.1, 0.0];  % x, x_dot, theta, theta_dot
tspan = [0, 10];

% 解决微分方程组
[t, x] = ode45(@derivs, tspan, x0);

% 绘制结果
plot(t, x(:, 1), 'r', t, x(:, 3), 'b');
legend('x', 'theta');

这两个代码都可以解决这个问题，但是Matlab的代码可能会更快一些，因为它使用了更高效的ode45函数。

Mike_python小

Python解二阶非线性微分方程组

import numpy as npfrom scipy.integrate import odeint# 定义微分方程组def dydt(y, t):    x, x_dot, theta, theta_dot = y    c = np.cos(theta)    s = np.sin(theta)        x_ddot = (F - Dx * x_dot - m*l*theta_dot*s + m*l*theta_ddot*c) / (m + m1)    theta_ddot = (m*g*l*s - Dtheta*theta_dot - m*l*x_dot*c) / (4/3*m*(l**2))        return [x_dot, x_ddot, theta_dot, theta_ddot]# 定义初值和时间范围y0 = [0.0, 0.0, 0.0, 0.0]  # 初值t = np.linspace(0, 10, 101)  # 时间范围# 求解微分方程组sol = odeint(dydt, y0, t)# 获取解的结果x = sol[:, 0]theta = sol[:, 2]# 打印解的结果for i in range(len(t)):    print(f"t = {t[i]:.2f}: x = {x[i]:.4f}, theta = {theta[i]:.4f}")

在上述代码中，我们首先定义了微分方程组，并将其表示为函数dydt(y, t)的形式。然后，利用odeint函数求解微分方程组的解。在微分方程组中，我们使用给定的参数和变量计算出 $x_{\ddot}$ 和 $\theta_{\ddot}$。最后，我们打印出每个时间点的 $x$ 和 $\theta$ 的值。

Matlab解二阶非线性微分方程组

[/code]matlab
% 定义微分方程组
function dydt = myodefun(t, y)
    x = y(1);
    x_dot = y(2);
    theta = y(3);
    theta_dot = y(4);
    c = cos(theta);
    s = sin(theta);
   
    x_ddot = (F - Dx*x_dot - m*l*theta_dot*s + m*l*theta_ddot*c) / (m + m1);
    theta_ddot = (m*g*l*s - Dtheta*theta_dot - m*l*x_dot*c) / (4/3*m*(l^2));
   
    dydt = [x_dot; x_ddot; theta_dot; theta_ddot];
end

% 定义初值和时间范围
tspan = [0 10];
y0 = [0.0; 0.0; 0.0; 0.0];

% 求解微分方程组
[t, y] = ode45(@myodefun, tspan, y0);

% 获取解的结果
x = y(:, 1);
theta = y(:, 3);

% 打印解的结果
for i = 1:length(t)
    fprintf('t = %.2f: x = %.4f, theta = %.4f\n', t(i), x(i), theta(i));
end
[/code]

在上述代码中，我们首先定义了微分方程组，并将其表示为函数myodefun(t, y)的形式。然后，利用ode45函数求解微分方程组的解。在微分方程组中，我们使用给定的参数和变量计算出 $x_{\ddot}$ 和 $\theta_{\ddot}$。最后，我们打印出每个时间点的 $x$ 和 $\theta$ 的值。

以上是Python和Matlab求解二阶非线性微分方程组的示例代码，你可以根据实际情况调整参数和初值，并按照需要进行结果的处理和输出。

球一个最佳答案谢谢啦！这对我非常重要！
回答可能有误，多多谅解！