基于Runge-Kutta积分法的系统：微分方程嵌入代码方法咨询

Got it, let's walk through how to embed these coupled differential equations into your code step by step. I'll use Python with SciPy's solve_ivp (a modern, flexible tool for initial value problems) since it’s straightforward for systems like this.

Step 1: Frame the System as a Single Function #

First, we need to package your two equations into a function that describes the rate of change of the entire state vector. Let’s define our state as z = [x, y]—this function will return [dx/dt, dy/dt] given the current time and state.

Step 2: Full Working Code Example #

Here’s a complete implementation that solves the equations and visualizes the results:

import numpy as np
from scipy.integrate import solve_ivp
import matplotlib.pyplot as plt

# Define the differential equation system
def ode_system(t, z):
    x, y = z  # Unpack the state vector into x and y
    dx_dt = -x * (2 - y)  # Your first differential equation
    dy_dt = y * (1 - 2*x)  # Your second differential equation
    return [dx_dt, dy_dt]

# Set up initial conditions and time range
initial_state = [1, 2]  # x(0)=1, y(0)=2
t_start = 0
t_end = 10  # Adjust this to solve for longer/shorter time periods
time_span = [t_start, t_end]

# Solve the ODE system
# t_eval ensures we get evenly spaced time points for clean plotting
solution = solve_ivp(ode_system, time_span, initial_state,
                     t_eval=np.linspace(t_start, t_end, 100))

# Extract results
t_values = solution.t
x_values = solution.y[0]
y_values = solution.y[1]

# Visualize the results
plt.figure(figsize=(10, 4))

# Plot x(t) and y(t) over time
plt.subplot(1, 2, 1)
plt.plot(t_values, x_values, label='x(t)')
plt.plot(t_values, y_values, label='y(t)')
plt.xlabel('Time t')
plt.ylabel('Value')
plt.title('x(t) and y(t) Over Time')
plt.legend()

# Plot phase portrait (x vs y)
plt.subplot(1, 2, 2)
plt.plot(x_values, y_values)
plt.xlabel('x')
plt.ylabel('y')
plt.title('Phase Portrait: x vs y')

plt.tight_layout()
plt.show()

Fixing Your Missing F_xy Section #

If your existing code references an F_xy function, that’s exactly the ode_system function we defined above. It takes the current time t and state z, computes the derivatives using your equations, and returns them as a list. You can drop this function directly into your code wherever F_xy was intended to be.

Alternative: Using odeint (Legacy Method) #

If you prefer the older, widely used odeint function, here’s the adapted code:

from scipy.integrate import odeint

# odeint uses a slightly different function signature: state first, then time
def ode_system_odeint(z, t):
    x, y = z
    dx_dt = -x * (2 - y)
    dy_dt = y * (1 - 2*x)
    return [dx_dt, dy_dt]

# Solve
t_values = np.linspace(t_start, t_end, 100)
solution = odeint(ode_system_odeint, initial_state, t_values)
x_values = solution[:, 0]
y_values = solution[:, 1]

# Plotting code remains identical to the example above

Run either version, and you’ll get the time evolution of x and y, plus a phase portrait showing their relationship. Tweak t_end if you need to solve for a longer or shorter duration.

内容的提问来源于stack exchange，提问作者Matt Simpson