MATLAB: корректность записи системы уравнений для ode45

[electric_retard]

Пытаюсь найти ошибку в программе, которая решает систему уравнений:
$\frac{dz}{dt}=\beta P$
$\frac{d\psi}{dt}=\frac{dz}{dt}\frac{SURF}{U}$
$\frac{dW}{dt}=\frac{m}{\rho}\frac{d\psi}{dt}-cm\frac{d\psi}{dt}$
$\frac{dP}{dt}=\frac{1}{W}\left\lbrace f\omega\frac{d\psi}{dt}-kp\frac{dW}{dt}\right\rbrace$
вот кусок функции, которая вычисляет производные для ode45

Используется синтаксис Matlab M

    function u_dot = dudtK(~,u)
        P = u(1);                %Давление
        W = u(2);                %Объём
        psi = u(3);              %Доля прогоревшей части заряда
        Vel = u(4);              %Скорость
        x = u(5);                %Координата
        z = u(6);                %Толщина прогоревшего свода
        
        z_dot = .110519e0*P/1e6;
        [SURF,~,U]=PRF017(.508e-1,.508e-1,1.07209,.280722,3.175,NaN,NaN,z,7,NaN,NaN); %По толщине прогоревшего свода возвращает текущую площадь горения, задано семиканальное зерно
        psi_dot = z_dot*SURF/U;
        
        Pb = P/(1+m/(3*q));
        W_dot = -psi_dot*c*m + psi_dot*m/r + S*Vel;
        P_dot = 1/W*(f*m*psi_dot - k*P*W_dot);
        Vel_dot = Pb*S/q*(P>Inf);
        x_dot = Vel;
        u_dot = [P_dot W_dot psi_dot Vel_dot x_dot z_dot]';
        u_dot = u_dot*(x<=L);
    end

Проблема в том, что при интегрировании системы уравнений конечное давление не совпадает с давлением по формуле для манометрической бомбы.
Вопрос: корректно ли в коде вычисляется производная по давлению $\frac{dP}{dt}$ относительно системы уравнений?

[Tane4ka]

У вас в коде для вычисления W_dot имеется слагаемое $S*Vel$, которого нет в вышеуказанной системе уравнений. Возможно в этом проблема?