>> % 定义常数和初值 a = 10;k = 10;tspan = [0, 10];x0 = [0, 0];% 定义外力 F = @(t) 5*sin(t/2) + 2*sin(t);% 采用欧拉法求解微分方程 c = 3;theta = 0.5;[t1, x1] = euler_wilson(theta, @(t, x) [x(2); ... (F(t) - c*x(2) - k*x(1))/a], tspan, x0); c = 20;theta = 0.5;[t2, x2] = euler_wilson(theta, @(t, x) [x(2); (F(t) - c*x(2) - k*x(1))/a], tspan, x0); % 绘制时间-位移曲线图 figure;plot(t1, x1(:, 1), 'r', t2, x2(:, 1), 'b'); xlabel('Time (s)'); ylabel('Displacement (m)'); legend('c = 3', 'c = 20');% 欧拉法求解微分方程的函数 function[t, x] = euler_wilson(theta, f, tspan, x0) h = 0.01; t = tspan(1):h:tspan(2); N = length(t); x = zeros(N, length(x0)); x(1, :) = x0; for n = 1:N-1 x_half = x(n, :) + h/2 * f(t(n), x(n, :)); x(n+1, :) = x(n, :) + h * (1-theta) * f(t(n), x(n, :)) ... + h * theta * f(t(n) + h/2, x_half); end end 无法执行赋值，因为左侧的大小为 1×2，右侧的大小为 2×2。 出错 euler_wilson (第 9 行) x(n+1, :) = x(n, :) + h * (1-theta) * f(t(n), x(n, :)) ...

这个错误发生在 euler_wilson 函数的第9行，因为左侧的 x(n+1,:) 是一个1x2的向量，而右侧的 h*(1-theta)*f(t(n),x(n,:)) 和 h*theta*f(t(n)+h/2,x_half) 都是2x1的向量，无法直接相加。可能的解决方法是将左侧的 x(n+1,:) 改为 x(n,:)，因为欧拉法中每一步的结果都是基于前一步的。