+++ Моделирование кинетических уравнений

[Няшка]

Всем доброго времени суток. А также, пользуясь случаем, поздравляю всех с прошедшими и наступающими праздниками!
Вопрос у меня касается моделирования кинетических уравнений, поэтому решила написать в раздел "Физическая химия". Так уж получилось, что мне нужно создать простейшую математическую модель для двух химических реакций:
А+В=С (k1=1) и пусть для простоты стадия будет необратимой, r1=k1*CA*CB;
В+С=Р (k2=1) и тоже пусть для простоты стадия будет необратимой, r2=k2*CB*CC.
У меня получилось, что:
dСА/dt=-k1*A*B=-r1;
dСВ/dt=-k1*A*B-k2*B*C=-r1-r2;
dСС/dt=k1*A*B-k2*B*C=r1-r2;
dСР/dt=k2*B*C=r2;

Соответственно в *.m-файл для расчета в MATLAB я записала:

function dydt=react(t,y)
A=y(1);
B=y(2);
C=y(3);
P=y(4);
global k1 k2
clc
clf
clear

dydt=zeros(4,1);
k1=1;
k2=1;

dydt(1)=-k1*A*B;
dydt(2)=-k1*A*B-k2*B*C;
dydt(3)=k1*A*B-k2*B*C;
dydt(4)=k2*B*C;

[t,y] = ode45(@react,[0 20],[0.0025 0.125 0 0]);
plot(t,y(:,1),'r',t,y(:,2),'-o',t,y(:,3),'b',t,y(:,4),'k')
end

А в ответ все время получается:

??? function dydt=react(t,y)
|
Error: Function definitions are not permitted in this context.

/

Хотя должны получиться четыре кинетические кривые.
В общем-то занимаюсь этим в первый раз, поэтому вообще не уверена в правильности.
Сам файл *.m-файл сохранен c:\Users\User\Documents\MATLAB\.
Может кто с таким сталкивался или хорошо в MATLAB разбирается, помогите, пожалуйста!

P.S. Знаю, что подобные модели можно делать в Excel, может тоже кто знает как? А то я Excel пользуюсь раз в год....

Всем заранее спасибо!

[amge]

Из самых общих соображений: обращение к функции не должно быть внутри определения функции (если нет рекурсии, а ее, очевидно, нет). Для начала попробуйте перенести end парой строк выше. И все глобальные действия типа clc вынести из определения функции.

Кстати, если цель - не самописная программы на матлаб, то есть изрядное количество специализированных кинетических программ. Одна из самых простых и понятных - KINET. Для вашей задачи - более чем...

[Няшка]

amge, спасибо за ответ. Я пробовала и перенести end парой строк выше. И все глобальные действия типа clc вынести из определения функции тоже, но фраза:

"??? function dydt=react(t,y)
|
Error: Function definitions are not permitted in this context."

все равно возникает, правда еще стало появляться окно графика, но самого графика в нем нет!

Про программу Kinet я знаю, но, к сожалению, у меня и в ней в данном случае ничего не получается (точнее не получается даже того, что можно было бы ожидать исходя из общих представлений)! Хотя ранее она мне помогала несколько раз. С программой MATLAB решила попробовать в рамках повышения квалификации и самообразования. Скачала пару книг из сети:
1. L. F. Shampine, I. Gladwell, S. Thompson-Solving ODEs with MATLAB-Cambridge University Press (2003).
2. Ronald Gieschke, Daniel Serafin (auth.)-Development of Innovative Drugs via Modeling with MATLAB_ A Practical Guide-Springer-Verlag Berlin Heidelberg (2014).
Но даже приведенные там примеры у меня заканчиваются фразой;

"??? function ....
|
Error: Function definitions are not permitted in this context."

Я уже, если честно, в отчаянии. В моей практике первый раз, когда я за неделю не могу разобраться с такой "простенькой" проблемой!

[Cherep]

по-моему тема должна называться "помогите с MATLAB-ом"

UPD врятли ответ по теме, но когда я был студентом вначале 00-х, мы такие задачки решали в MathCAD, который ЕМНИП был гораздо проще в освоении. Один мой знакомый, который в ту пору учился на ВМК МГУ после моих слов сказал "а чё ты хотел? маткад попсовее чем матлаб"

И ещё, есть замечательный ресурс, на который лично я выходил почти всегда, когда встречался гемор с latex. Там есть раздел и про матлаб. http://stackoverflow.com/questions/tagged/matlab

[alien308]

Да Mathcad полегче для простых задач. Попробуйте одно уравнение, систему из двух ... Как правило помогает сначала решить простую задачу.

[Няшка]

Всем доброго времени суток.
Задачу свою я решила. На всякий случай, вдруг еще кому понадобится, пишу ход решения.

1. в *.m-файл-функцию я записала:
function dydt=react(t,y)

dydt=zeros(4,1);
k1=1; k2=1;

A=y(1);
B=y(2);
C=y(3);
P=y(4);

dydt(1) = -k1*A*B;
dydt(2) = -k1*A*B-k2*B*C;
dydt(3) = k1*A*B-k2*B*C;
dydt(4) = k2*B*C;

end

2. Далее в командном окне MATLAB ввела:
>>[t,y] = ode45('react',[0 100],[0.0025 0.125 0 0]); % это вычисление самой функции.
>>length(t) % число точек, откладываемое по оси ОХ (в принципе можно и без него.
>>size(y) % ---||---
>>plot(t,y(:,1),'r',t,y(:,2),'-o',t,y(:,3),'b',t,y(:,4),'k') % вывод итогового графика.
Но поскольку меня, в основном, интересует изменение концентрации вещества С, то я использовала
>>plot(t,y(:,3),'b')

Всем добра!

[Cherep]

я извиняюсь за занудство, однако неплохо было бы поменять название темы и вкратце упомянуть, в чем был косяк(и)

[Няшка]

Косяк, как я понимаю, был в том, что не надо путать понятия Script M-Files и Function M-Files в MATLAB.
На счет смены названия я не уверена, поскольку в тексте моего первого сообщения упоминались и MATLAB и Excel, а коллеги посоветовали еще и Kinet. Честно говоря, я бы еще задала пару вопросов по MathCad и Maple по той же теме, но поскольку задача уже решена, думаю, что не стоит этого делать.

[Cherep]

спасибо!