Математическая формулировка

Расположение массы $$i$$ является $$x(i)$$, с горизонтальной координатой $$x_1(i)$$ и вертикальная координата $$x_2(i)$$. Масса $$i$$ имеет потенциальную энергию от силы тяжести $$gm(i)x_2(i)$$. Потенциальная энергия весной $$i$$ является $$k(i)(d(i)-l(i){)}^2/2$$, где $$d(i)$$ - длина пружины между массами $$i$$ и масса $$i-1$$. Примите якорную точку 1 как положение массы 0 и якорную точку 2 как положение массы $$n+1$$. Предыдущий энергетический расчет показывает, что потенциальная энергия пружины $$i$$ является

$$Energy(i)=\frac{k(i)(\Vert x(i)-x(i-1)\Vert -l(i){)}^2}{2}$$.

Переформулирование этой потенциальной энергетической задачи как задачи программирования конуса второго порядка требует введения некоторых новых переменных, как описано в Лобо [1]. Создайте переменные $$t(i)$$ равным квадратному корню термина $$Energy(i)$$.

Давайте $$e$$ быть модулем вектора-столбца $$\left[\genfrac{}{}{0ex}{}{0}{1}\right]$$. Тогда $$x_2(i)=e^{T}x(i)$$. Задача становится

$$\underset{x,t}{\min }(\sum _{i}gm(i)e^{T}x(i)+\Vert t{\Vert }^2).$$               (1)

Теперь рассмотрим $$t$$ как переменная свободного вектора, не заданная предыдущим уравнением для $$t(i)$$. Включите отношения между $$x(i)$$ и $$t(i)$$ в новом наборе конусных ограничений

$$\Vert x(i)-x(i-1)\Vert -l(i)\le \sqrt{\frac{2}{k(i)}}t(i).$$   (2)

Целевая функция еще не линейна в своих переменных, как требуется для coneprog. Введите новую скалярную переменную $$y$$. Заметьте, что неравенство $$\Vert t{\Vert }^2\le y$$ эквивалентно неравенству

$$\Vert \left[\genfrac{}{}{0ex}{}{2t}{1-y}\right]\Vert \le 1+y$$.                               (3)

Теперь проблема в том, чтобы минимизировать

$$\underset{x,t,y}{\min }(\sum _{i}gm(i)e^{T}x(i)+y)$$                     (4)

удовлетворяющее ограничениям конуса на $$x(i)$$ и $$t(i)$$ приведенный в (2) и дополнительное ограничение конуса (3). Ограничение конуса (3) гарантирует, что $$\Vert t{\Vert }^2\le y$$. Поэтому задача (4) эквивалентна задаче (1).

Целевая функция и ограничения конуса в задаче (4) подходят для решения с coneprog.

Состав MATLAB ®

Задайте шесть коэффициентов упругости $$k$$, шесть констант длины $$l$$и пять масс $$m$$.

k = 40*(1:6);
l = [1 1/2 1 2 1 1/2];
m = [2 1 3 2 1];
g = 9.807;

Задайте переменные оптимизации, соответствующие математическим переменным задачам. Для простоты установите якорные точки как две виртуальные масс-точки x(1,:) и x(end,:). Эта композиция позволяет каждой пружине растягиваться между двумя массами.

nmass = length(m) + 2;
% k and l have nmass-1 elements
% m has nmass - 2 elements
x = optimvar('x',[nmass,2]);
t = optimvar('t',nmass-1,'LowerBound',0);
y = optimvar('y','LowerBound',0);

Создайте задачу оптимизации и установите целевую функцию в выражение в (4).

prob = optimproblem;
obj = dot(x(2:(end-1),2),m)*g + y;
prob.Objective = obj;

Создайте ограничения конуса, соответствующие выражению (2).

conecons = optimineq(nmass - 1);
for ii = 1:(nmass-1)
    conecons(ii) = norm(x(ii+1,:) - x(ii,:)) - l(ii) <= sqrt(2/k(ii))*t(ii);
end
prob.Constraints.conecons = conecons;

Задайте якорные точки anchor0 и anchorn. Создайте ограничения равенства, определяющие, что две виртуальные массы концов расположены в точках привязки.

anchor0 = [0 5];
anchorn = [5 4];
anchorcons = optimeq(2,2);
anchorcons(1,:) = x(1,:) == anchor0;
anchorcons(2,:) = x(end,:) == anchorn;
prob.Constraints.anchorcons = anchorcons;

Создайте ограничение конуса, соответствующее выражению (3).

ycone = norm([2*t;(1-y)]) <= 1 + y;
prob.Constraints.ycone = ycone;

Решите задачу

Формулировка задачи завершена. Решите проблему, позвонив solve.

[sol,fval,eflag,output] = solve(prob);

Solving problem using coneprog.
Optimal solution found.

Постройте график точек решения и анкеров.

plot(sol.x(2:(nmass-1),1),sol.x(2:(nmass-1),2),'ro')
hold on
plot([sol.x(1,1),sol.x(end,1)],[sol.x(1,2),sol.x(end,2)],'ks')
plot(sol.x(:,1),sol.x(:,2),'b--')
legend('Calculated points','Anchor points','Springs','Location',"best")
xlim([sol.x(1,1)-0.5,sol.x(end,1)+0.5])
ylim([min(sol.x(:,2))-0.5,max(sol.x(:,2))+0.5])
hold off

Можно изменить значения параметров m, l, и k чтобы увидеть, как они влияют на решение. Можно также изменить количество масс; код берёт количество масс из предоставляемых вами данных.

Ссылки

[1] Лобо, Мигель Соуза, Ливен Ванденберге, Стивен Бойд и Эрве Лебрет. «Приложения программирования конуса второго порядка». Линейная алгебра и ее приложения 284, № 1-3 (ноябрь 1998): 193-228. https://doi.org/10.1016/S0024-3795(98)10032-0.