シミュレーションのイ​ンプットパラメータを​最小化するには？

10 次查看（过去 30 天）
Shigenori Ishihara
评论： Tohru Kikawada 2017-2-22
例えば、物体の質量・初速・射出角度を与えて、重力と空気抵抗の影響を考慮した運動計算を行い、到達高度・落下位置までの距離を出力するシミュレーションにおいて、所定の条件を満たすように質量を最小化する（設計値を目的関数に据える）ことは可能でしょうか？ ＞x=質量、初速、射出角度, 非線形制約１；到達高度　設定値以上, 非線形制約２；落下位置～射点　設定値以上,　目的関数=質量　　の設定では、初期値付近で最適化を終了してしまいます。 ＞x=質量、初速、射出角度, 非線形制約１；落下位置～射点　設定値以上,　目的関数=到達高度-設定値　として、質量の初期値の与え方を工夫するのが一般的なのでしょうか？
  1 个评论
Jiro Doke
Jiro Doke 2017-2-9
質問からは少しイメージしづらいので、試されたコードがありましたら、再現できるものを追加して下さい。

请先登录，再进行评论。

请先登录，再回答此问题。

采纳的回答

Tohru Kikawada
Tohru Kikawada 2017-2-13
编辑：Tohru Kikawada 2017-2-22
Jiro Dokeさんのコメントのとおりもう少し具体的な処理がイメージできる形で書いていただけると迅速な回答が期待できます。
参考： Tutorial: MATLAB Answersで早く的確な回答を得るためのポイント
いただいた情報を元に簡単なサンプルを作ってみました。
一般的にこのような非線形最適化問題を fmincon のような勾配ベースの探索手法を使うと局所解に陥ってしまう可能性があります。
大域的な探索を行いたい場合には Global Optimization Toolbox の最適化ソルバーをご活用ください。
function [x,f,eflag,outpt] = ans_324101
xLast = []; % Last place computeall was called
myf = []; % Use for objective at xLast
myc = []; % Use for nonlinear inequality constraint
myceq = []; % Use for nonlinear equality constraint
x_pts = [];
y_pts = [];
fun = @objfun; % the objective function, nested below
cfun = @constr; % the constraint function, nested below
options = optimoptions('fmincon','MaxFunctionEvaluations',500,...
'PlotFcn',@optimplotfval,'OutputFcn',@outfun);
% Solve Second-Order ODE with Initial Conditions
% syms y(t) m k y0 g theta v0;
% Dy = diff(y);
% y(t) = dsolve(m*diff(y,t,2) == -k*diff(y,t) - m*g,y(0)==y0, Dy(0)==v0*cos(theta))
% myf = matlabFunction(y,'file','my_y','vars',{t,g,v0,y0,k,m,theta})
% Initial values
x0 = [10 deg2rad(30) 1];
% Conditions
x_0 = 0;
y_0 = 0;
y_threshold = 3;
x_threshold = 3;
% Visualization
figure;
hopt = plot(0,0);
hold on;
plot([x_0+x_threshold; x_0+x_threshold],[y_0; y_0 + y_threshold*1.5],'r');
plot([x_0; x_threshold*1.5],[y_0 + y_threshold; y_0 + y_threshold],'r');
axis([0 x_threshold*1.5 0 y_threshold*1.5]);
% Lower and upper bounds
lb = [0 0 0];
ub = [10 deg2rad(90) 10];
% Call fmincon
[x,f,eflag,outpt] = fmincon(fun,x0,[],[],[],[],lb,ub,cfun,options);
function y = objfun(x)
if ~isequal(x,xLast) % Check if computation is necessary
[myf,myc,myceq,~,x_tmp,y_tmp] = myModel(x);
xLast = x;
x_pts = x_tmp;
y_pts = y_tmp;
end
% Now compute objective function
y = myf;
end
function [c,ceq] = constr(x)
if ~isequal(x,xLast) % Check if computation is necessary
[myf,myc,myceq,~,x_pts,y_pts] = myModel(x);
xLast = x;
end
% Now compute constraint functions
c = myc; % In this case, the computation is trivial
ceq = myceq;
end
function [f,c,ceq,t,x_pts,y_pts] = myModel(x)
% Ref: https://ja.wikipedia.org/wiki/%E6%96%9C%E6%96%B9%E6%8A%95%E5%B0%84
v_0 = x(1);
theta = x(2);
m = x(3);
t = 0:0.01:5;
k = 0.5;
g = 9.8;
x_pts = m*v_0/k*(1-exp(-k/m*t))*cos(theta)+x_0;
y_pts = -m/k*((v_0*sin(theta) + m/k*g)*(exp(-k/m*t)-1)+g*t)+y_0;
[~,idxMaxY] = max(y_pts);
c(1) = -y_pts(idxMaxY)+ y_threshold;
[~,idxYZero] = min(abs(y_pts(idxMaxY+1:end)));
c(2) = -x_pts(idxMaxY+idxYZero) + x_threshold;
f = m;
ceq = [];
end
function stop = outfun(~,~,~)
% Visualization
stop = false;
set(hopt,{'XData','YData'},{x_pts,y_pts});
end
fprintf('v0=%f[m/s], theta=%f[deg], m=%f[kg]\n',x(1),rad2deg(x(2)),x(3));
end
軌道：
質量の変化：
%
  2 个评论
Shigenori Ishihara
Shigenori Ishihara 2017-2-14
Tohru Kikawada様 有難うございます。早速提示頂いたサンプルを試してみました。 瑣末ですが、非線形不等号制約c(2)≦0を満たさずに終わっているのは　どうしてなのでしょうか？
Tohru Kikawada
Tohru Kikawada 2017-2-22
c(2) の計算方法が正しくなかったので修正しました。
制約が守られていないなど終了したときの状態を知りたい場合には fmincon の戻り値の exitflag をご覧ください。

请先登录，再进行评论。

更多回答（0 个）

类别

在 Help Center 和 File Exchange 中查找有关 Signal Processing Toolbox 的更多信息

标签

产品

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by