troubles with ode45 - Programmers Heaven

Howdy, Stranger!

It looks like you're new here. If you want to get involved, click one of these buttons!

Categories

troubles with ode45

NotYetMScNotYetMSc Posts: 1Member
Hello everybody!
I've written a program:
[code]
clear all
close all
clc

global mk Kfi g lk alfa beta Kz Jk a0 nsteps1
g=9.81;
E=7*10^10;
alfa=20*pi/180;
L=50;
N=5;
mk=300;
lk=1.5;
H=40000;
ro1=8.45;
d=0.04;
v0=0;
a0=0.5;
%**************************************************************************
roz1=ro1+(N*mk)/L;
Jk=mk*lk^2;
I=d^4*pi/64;
Kfi=(L/(N*cos(alfa)))/6*E*I;
%**************************************************************************
h=0.01; %step
%**************************************************************************
droga=L/cos(alfa);
delta=v0^2+4*(a0/2)*droga;
tc=(-v0+sqrt(delta))/a0;
Vmax=3;
V=v0+a0*tc;
if V>Vmax
t1=(Vmax-v0)/a0;
s1=v0*t1+((a0*t1^2)/2);
s2=droga-s1;
t2=s2/Vmax;
tc=t1+t2;
t01=[0:h:t1];
t02=[h:h:t2];
Lk1=cos(alfa).*((v0.*t01+(a0.*t01.^2))./2);
Lk2=Lk1(1,end)+(cos(alfa).*(Vmax.*t02));
Lk=[Lk1 Lk2];
t0=[0:h:tc];
%modyfying (x/0)
Lkk=Lk;
Lkk(Lkk<1)=1;
Lkk(Lkk>(L-1))=(L-1);
z=((roz1*g)/(2*H*(cos(alfa))^2)).*Lk.*(L-Lk);
beta=atan((roz1*g/H)*(L/2-Lk));
prosta=-tan(alfa).*Lk;
tor=prosta-z;
Kz=H*((1./Lkk)+(1./(L-Lkk)));
%beta=beta';
%beta=median(beta,1);
%Kz=Kz';
%Kz=median(Kz,1);
nsteps1=t1/h;
%nsteps2=t2/h;
tspan1 = linspace(0,t1,nsteps1);
tspan2 = linspace(t1,t2,nsteps2);
[tc1,x1]=ode45(@przysp1, tspan1, [0.02 0 0.0028 0]); [code/]
and function:
[code]
function xprime=przysp1(tc1,x1)
global mk Kfi g lk alfa beta Kz Jk a0 nsteps1
i=1;

xprime=zeros(4,1);
xprime(1)=x1(2);
xprime(2)=(mk*(Kfi*x1(1)+g*mk*lk*sin(x1(1)))*cos(2*(alfa+beta(:,i)))+(x1(2)^2)*(mk^2)*...
(lk^2)*(sin(2*(x1(1)))*cos(2*(alfa+beta(:,i))))-cos(x1(1)-alfa-beta(:,i))*...
sin(x1(1)+alfa+beta(:,i))*a0*(mk^2)*lk*sin(-beta(:,i))*sin(x1(1)+alfa+beta(:,i))-...
cos(alfa-x1(1))*cos(2*(alfa+beta(:,i)))-Kz(:,i)*x1(3)*mk*lk*sin(x1(1)+alfa+beta(:,i)))/...
(mk*cos(2*(alfa+beta(:,i)))*(mk*(lk^2)+Jk)-(mk^2)*(lk^2)*sin(x1(1)+alfa+beta(:,i))*...
sin(x1(1)-alfa-beta(:,i)));
xprime(3)=x1(4);
xprime(4)=(xprime(2)*mk*lk*sin(x1(1)-alfa-beta(:,i))-Kz(:,i)*x1(3)+mk*a0*sin(-beta(:,i))-mk*lk*...
(x1(2)^2)*cos(x1(1)-alfa-beta(:,i)))/(mk*cos(2*(alfa+beta(:,i))));
[code/]

I want that "i" will increase in each step. (in first step i=1, in secound i=2... and so on).

Now program works only with all this values. If I change a0 lk, mk... matlab show me: "Warning: Failure at t=2.802832e+000. Unable to meet integration tolerances without reducing the step size below the smallest value allowed (9.957659e-015) at time t.
(Type "warning off MATLAB:ode45:IntegrationTolNotMet" to suppress this warning.)
> In C:Matlab6p5 oolboxmatlabfunfunode45.m at line 335
In C:Documents and Settings...proba.m at line 64 "
I dont know what is wrong. Please help me... and sorry for my english
Sign In or Register to comment.