ode45を用いて二階常微分方程式を解こうとしていますが、
なぜか入力引数が不足しているという表示が出てきてしまいます。
入力した文字に対する数値の個数はあっているはずですし、色々弄ってみてはいるものの同じエラーが出てきて困っています。
下の図のようなシーソー型の振り子について、
二階常微分方程式で表される、振り子の振動を表す方程式を解こうとしています。
図左側の機械からは推力が発生し、振り子には重力方向に力が加わって変位が発生し、
図右側のカウンターウェイト（振り子の変位を復元するための重り）によって、
発生した振り子の振れ角θ（図にはありませんが…）は復元され、元の位置に戻ります。
この時推力F(t)は,時間t=0の時は0で、t>0の場合はFtなる推力を定常的に発生させます。
振動方程式は以下のように表現されます。
当初ラプラス変換で解こうとしましたが、頓挫したのでodeの使用を試みています。
以下に現状のコードを記載します。
clc
clear
close all
syms ze om I F(t) L s(t)
assume([I L s(t)]>0)
assume([L]<0.261)
assume([F(t)]>=0)
syms sita(t)
assume(abs(sin(sita(t))-sita(t))<=0.01 & abs(1-(1/2)*sita(t)^2)<=0.01)
ds2=diff(sita,t,2)
ds1=diff(sita,t,1)
eqn = ds2+2*ze*om*ds1+(om^2)*sita(t) == F(t)*L/I
[V] = odeToVectorField(I*ds2+2*ze*om*ds1+om^2*sita(t)==F(t)*L/I)
cond1=F(0)==0
cond2=sita(0)==0
cond3=ds1(0)==0
cond4=ds2(0)==0
M = matlabFunction(V,'vars', {'t','Y','I','L','om','ze','F'})
m1=1
m2=3
m3=0.6
r1=0.04
r2=0.04
L=0.26
Lcm1=0.24
Lcm2=0.200
l=0.04
J1=(1/2)*m1*r1^2+m1*Lcm1^2
J2=((m2*(3*r2^2+(l/2)^2))/12)+m2*Lcm2^2
J3=(m3*L^2)/12
I=J1+J2+J3
ze=0
g=9.8
k=m2*g*Lcm2-m1*g*Lcm1
om=sqrt(k/I)
F(t)=10^-9
sol=ode45(M,[0 50],[0 10],[I],[L],[om],[0],[10^-9]);
エラーが表示されている時、以下のような表示がでます。
入力引数が不足しています。
エラー:
symengine>@(t,Y,I,L,om,ze,F)[Y(2);-1.0./I.^2.*(-L.*F(t)+I.*om.^2.*Y(1)+I.*om.*ze.*Y(2).*2.0)]
エラー: odearguments (line 90)
f0 = feval(ode,t0,y0,args{:});
エラー: ode45 (line 115)
odearguments(FcnHandlesUsed, solver_name, ode, tspan, y0, options, varargin);
エラー: thrust_stand_ode (line 73)
sol=ode45(M,[0 50],[0 10],[I],[L],[om],[0],[10^-9]);
元々、「関数または'F'が未定義です」という表示も出ていたのですが、
色々弄っているうちに無くなっていました。
