用MATLAB，已知两个固定点的坐标，两点之间拉一条悬链线，考虑重力和弹性形变，求悬链线的形状

求解悬链线的形状需要先建立悬链线的方程，并且考虑重力和弹性形变的影响。可以利用拉格朗日方程求解得到悬链线的方程。

具体步骤如下：

1. 建立拉格朗日方程：

定义悬链线的弧长为`s`，悬链线长度为$L$，重力加速度为$g$，弹性模量为$E$，悬链线的形状可用弧长$s$表示，则悬链线的势能为：

$$
V(s) = mg\int_0^s\sqrt{1 + \left(\frac{dy}{dx}\right)^2}dx
$$

式中，$m$为悬链线的线密度，$y(x)$为悬链线的形状。悬链线的势能式可以看作是悬链线的重力势能，因此可以利用拉格朗日方程建立如下的式子：

$$
\frac{d}{ds}\left(\frac{\partial L}{\partial \frac{dy}{dx}}\right) - \frac{\partial L}{\partial y} = 0
$$

式中，$L$为拉格朗日量。

2. 建立拉格朗日量：

悬链线的拉格朗日量可以表示为：

$$
L = T - V
$$

式中，$T$为悬链线的动能，$V$为悬链线的势能。悬链线的动能可以表示为：

$$
T = \frac{1}{2}m\int_0^s\left(\frac{dy}{dt}\right)^2dt
$$

将上式代入拉格朗日量中，可得到：

$$
L = \frac{1}{2}m\int_0^s\left(\frac{dy}{dt}\right)^2dt - mg\int_0^s\sqrt{1 + \left(\frac{dy}{dx}\right)^2}dx
$$

3. 求解拉格朗日方程：

将拉格朗日量带入拉格朗日方程中，可得到：

$$
\frac{d}{ds}\left(\frac{\partial L}{\partial \frac{dy}{dx}}\right) - \frac{\partial L}{\partial y} = 0
$$

代入拉格朗日量的表达式，可得到如下的微分方程：

$$
\frac{d^2y}{dx^2}\left(1 + \left(\frac{dy}{dx}\right)^2\right)^{\frac{1}{2}} + \frac{mg}{E}y = 0
$$

将上式变形，可得到如下的常微分方程：

$$
\frac{d^2y}{dx^2} = -\frac{mg}{E}\frac{y}{\sqrt{1 + \left(\frac{dy}{dx}\right)^2}}
$$

上式可以通过MATLAB中的ode45函数求解得到悬链线的形状。

完整代码如下：

% 悬链线上的两个固定点
p1 = [0, 0]
p2 = [10, 0]

% 悬链线的长度
L = 5

% 重力加速度
g = 9.81

% 弹性模量
E = 1

% 线密度
m = 1

% 定义常微分方程
odefun = @(x, y) [y(2); -g/m*y(1)/sqrt(1+y(2)^2)/E]

% 定义边界条件
bcfun = @(ya, yb) [ya(1); yb(1)-L]

% 初始值
xinit = [0; L/2]
yinit = [p1(2); 0]

% 求解常微分方程
[x, y] = ode45(odefun, [0, L], yinit, odeset('Events', bcfun), xinit)

% 将悬链线的形状绘制出来
plot(x, y(:, 1))
axis equal

运行以上代码即可求解得到悬链线的形状。