2013-05-27 2 views
1

У меня есть небольшая проблема. У меня есть 2 уравнения движения «ph» и «ph2» Я не знаю, как установить ODE, чтобы остановить вычисление «ph», когда x (1)> 0.111, а затем снова начинает вычислять «ph2» только до 0.111, после что график «ph» + «ph2» на один график зависит от времени «w», я думаю, что мне нужно установить некоторые ограничения по времени, но не знаю, как это сделать. Я использую HELP, но не имею никакой пользы для меня.Условия запуска/остановки MATLAB ODE

[t,y] = ode45(@ph,[0,w_max],[0,0]); 

function dx = ph(tt,x) 
global F1 c m_c Ff p w s ln f_t sig dstr Ren pn Fex Fzmax xz xn l Fz mn 
Fpp = F1 + c*x(1); 

if pn<0 
pn=abs(pn); 
end 

if x(1)<ln 

    pn=spline(w,p,tt)-((2*sig)/dstr*Ren);  
    Fex=3.1416.*f_t.*pn.*(ln-x(1)); 
end 

if x(1)<42e-5 
    Fz = Fzmax*(1-(1/xz)*(x(1)+l)); 
end 

if x(1)>44e-3 
    m_c=m_c-mn; 
end 
dx=[x(2);((spline(w,p,tt)*s)-Fpp-Ff-Fex-Fz)./m_c]; 


[t2,y2] = ode45(@ph2,[0,w_max],[0,0]); 

function dx=ph2(tt,x) 

    global Fv m_c g f alfa Fzp c m_nbp 

    Ft=m_c*g*f; 
    Fv = 2*f*(Fzp/cos(alfa)); 

    if x(1)>0.44 

    m_c=m_c+m_nbp 

    end 

    dx = [x(2);((x(1)*c)-Ft-Fv)/m_c]; 
+0

Вы хотите сказать, что вы действительно хотите построить кусочную функцию как ваш ODE? – wakjah

+0

Я не знаю лучшего решения: o/ – user2401142

ответ

4

Для остановки расчета тел при x(1) > 0.111 вы можете использовать Event МЕСТОПОЛОЖЕНИЕ (manual page и example о том, как использовать его). На практике это функция, оцениваемая на каждом временном шаге, и если возвращаемое значение равно 0, то ode45 останавливает интеграцию.

Добавить функцию

function [value,isterminal,direction] = events(t,y) 
% Locate the time when y passes through 0.111 in all 
% directions and stop integration. 
value = y(1) - 0.111; % Detect y=0.111 
isterminal = 1;  % Stop the integration 
direction = 0;   % All direction 

и поставил его options = odeset('Events',@events) перед вызовом

[t,y] = ode45(@ph,[0,w_max],[0,0],options); 

Но учитывая, что рН выходы dx=[x(2); ...], для того, чтобы сделать проверку на х (1), необходимо вывести также эта переменная - что-то вроде dx=[x(1); x(2);...]

Надеюсь, это поможет.

+0

Это может помочь посмотреть демоверсию 'ballode', включенную в Matlab. Введите «edit ballode» в командном окне, чтобы увидеть код. Обратите внимание, что в этом примере они высевают последующие итерации с информацией из предыдущей, что может значительно повысить эффективность: 'options = odeset (опции, 'InitialStep', t (nt) -t (nt-уточнение), 'MaxStep', t (нт) -t (1)); '. – horchler

+0

Кроме того, использование глобальных переменных не является хорошей идеей. Вы должны передать свои параметры с помощью [анонимных функций] (http://www.mathworks.com/help/matlab/matlab_prog/anonymous-functions.html). См. Мой ответ на [этот вопрос] (http://stackoverflow.com/questions/16598558/finding-the-maxima-of-a-function-using-ode45/16639006#16639006) для примеров того, как передавать параметры на ваш функции интеграции и функции событий, используя 'ode45'. – horchler

+0

Thx много парней. Красиво сделано! Мне нужно работать, после чего я могу отлаживать и оптимизировать это. – user2401142

Смежные вопросы