Ode where 'x' is a column vector
조회 수: 1 (최근 30일)
이전 댓글 표시
I’m trying to solve the following ODE:
[m]x’’+[c]x’+[k]x = F
where [m] is a diagonal matrix, [c] and [k] are 7x7 matrices and both F and x are 7x1 matrices which vary with time.
I initially tried to solve it using an anonymous function but I get a ‘Matrix dimensions must agree error’ at least in part because it makes 'z' a 14x1 instead of a 2x7.
m=840; mf=53; mr=76;
Ix=820; Iy=1100;
a1=1.4; a2=1.47;
b1=0.7; b2=0.75;
w=b1+b2;
kf=10000; kr=13000;
ktf=200000; ktr=ktf;
kR=25000;
cf=10000; cr=12000;
v=20; d1=20; d2=0.1; wex=(2*pi*v)/d1;
M=zeros(7,7);
M(1,1)=m; M(2,2)=Ix; M(3,3)=Iy; M(4,4)=mf;
M(5,5)=mf; M(6,6)=mr; M(7,7)=mr;
K=zeros(7,7);
K(1,1)=2*kf+2*kr;
K(2,1)=b1*kf-b2*kf-b1*kr+b2*kr; K(1,2)=K(2,1);
K(3,1)=2*a2*kr-2*a1*kf; K(1,3)=K(3,1);
K(2,2)=kR+(b1^2)*kf+(b2^2)*kf+(b1^2)*kr+(b2^2)*kr;
K(3,2)=a1*b2*kf-a1*b1*kf-a2*b1*kr+a2*b2*kr; K(2,3)=K(3,2);
K(4,2)=-b1*kf-(1/w)*kR; K(2,4)=K(4,2);
K(5,2)=b2*kf+(1/w)*kR; K(2,5)=K(5,2);
K(3,3)=2*kf*(a1^2)+2*kr*(a2^2);
K(4,4)=kf+ktf+(1/(w^2))*kR; K(5,5)=K(4,4);
K(1,4)=-kf; K(1,5)=K(1,4); K(1,6)=-kr; K(1,7)=K(1,6);
K(2,6)=b1*kr; K(2,7)=-b2*kr;
K(3,4)=a1*kf; K(3,5)=K(3,4); K(3,6)=-a2*kr; K(3,7)=K(3,6);
K(4,1)=K(1,4); K(4,3)=K(3,4); K(4,5)=-kR/(w^2);
K(5,1)=K(1,5); K(5,3)=K(3,5); K(5,4)=K(4,5);
K(6,1)=K(1,6); K(6,2)=K(2,6); K(6,3)=K(3,6); K(6,6)=kr+ktr;
K(7,1)=K(1,7); K(7,2)=K(2,7); K(7,3)=K(3,7); K(7,7)=kr+ktr;
C=zeros(7,7);
C(1,1)=2*cf+2*cr;
C(2,1)=b1*cf-b2*cf-b1*cr+b2*cr; C(1,2)=C(2,1);
C(3,1)=2*a2*cr-2*a1*cf; C(1,3)=C(3,1);
C(2,2)=(b1^2)*cf+(b2^2)*cf+(b1^2)*cr+(b2^2)*cr;
C(3,2)=a1*b2*cf-a1*b1*cf-a2*b1*cr+a2*b2*cr; C(2,3)=C(3,2);
C(3,3)=2*cf*(a1^2)+2*cr*(a2^2);
C(1,4)=-cf; C(1,5)=C(1,4); C(1,6)=-cr; C(1,7)=C(1,6);
C(2,4)=-b1*cf; C(2,5)=b2*cf; C(2,6)=b1*cr; C(2,7)=-b2*cr;
C(3,4)=a1*cf; C(3,5)=C(3,4); C(3,6)=-a2*cr; C(3,7)=C(3,6);
C(4,1)=C(1,4); C(4,2)=C(2,4); C(4,3)=C(3,4); C(4,4)=cf;
C(5,1)=C(1,5); C(5,2)=C(2,5); C(5,3)=C(3,5); C(5,5)=cf;
C(6,1)=C(1,6); C(6,2)=C(2,6); C(6,3)=C(3,6); C(6,6)=cr;
C(7,1)=C(1,7); C(7,2)=C(2,7); C(7,3)=C(3,7); C(7,7)=cr;
iM=inv(M);
odefun=@(t,z) [z(2,:); iM*([0;0;0;(d2/2)*sin(wex*t)*ktf;(d2/2)*sin(wex*t)*ktf;(d2/2)*sin(wex*t)*ktf;(d2/2)*sin(wex*t)*ktf]-C*z(2,:)-K*z(1,:))];
tspan=0:0.01:10; ic=zeros(2,7);
[t,z]=ode45(odefun,tspan,ic);
I then tried calling it from a separate script and faced the same problem:
[t,z]=ode45(@Txt_func,[0,10],zeros(2,7));
Note that once I get past this stage, I’m aiming to optimise (probably genetic) the ‘k_’ and ‘c_’ values so, as I understand it, an anonymous solution would be more straightforward.
Also I'd like to vary the components in 'F' from
[0;0;0;(d2/2)*sin(wex*t)*ktf;(d2/2)*sin(wex*t)*ktf;(d2/2)*sin(wex*t)*ktf;(d2/2)*sin(wex*t)*ktf]
to the following
v=20; d1=20; d2=0.1; wex=(2*pi*v)/d1;
tau=(a1+a2)/v;
if t<=2
y1=(d2/2)*sin(wex*t);
else
y1=0;
end
if (t>=0.2 && t<=2.2)
y2=(d2/2)*sin(wex*(t-0.2));
else
y2=0;
end
if (t>=tau && t<=2+tau)
y3=(d2/2)*sin(wex*(t-tau));
else
y3=0;
end
if (t>=0.2+tau && t<=2.2+tau)
y4=(d2/2)*sin(wex*(t-tau-0.2));
else
y4=0;
end
F=zeros(7,1);
F(4,:)=y1*ktf;
F(5,:)=y2*ktf;
F(6,:)=y3*ktr;
F(7,:)=y4*ktr;
But 't' is only defined within ODE.
Thanks in advance.
댓글 수: 0
채택된 답변
Jan
2016년 1월 3일
편집: Jan
2016년 1월 3일
The problem gets much simpler, if you write the function to be integrated as an M-file. Then determine the parameters using an anonymous function as described here: http://www.mathworks.com/matlabcentral/answers/1971
For the problem of the discontinuities in y see: http://www.mathworks.com/matlabcentral/answers/59582#answer_72047. Matlab's ODE integrators are not specified to handle discontinous functions, such the the result is numerically instable. The reliable solution is to break the integration into time intervals with differentiable parameters.
Do not multiply the inverse of a matrix, but use the \ operator.
추가 답변 (0개)
참고 항목
카테고리
Help Center 및 File Exchange에서 Numerical Integration and Differential Equations에 대해 자세히 알아보기
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!