动力学 matlab,用Matlab求反应动力学参数

中间先拟合为一段一级和一段二级，并求出各点的导数，最后拟合出需要的参数。拟合过程中出现导数为负数，求各位大神指导！

代码如下：

clc;

clear all;

format short g

%%%分段求解  一级

expE1=[0        2054.794521

300        903.6835616

600        449.3047945];

t1=expE1(:,1);

C0=expE1(1,2);

C1=expE1(:,2);

LCO=log(C0);

LC=log(C1);

d=LCO-LC;

K0=[0.02];

K=lsqcurvefit(@FunS,K0,t1,d);

mdcdt1=K.*C1;

%%%%分段求解 二级

expE2=[900        200.3150685

1200        107.0773973

1500        78.67123288

2100        57.00136986

2700        47.67328767

3300        38.26164384

3900        28.4

4500        21.78219178];

C02=expE2(1,2);

C2=expE2(:,2);

t2=expE2(:,1);

RC2=1./C2;

d2=RC2-1/C02;

K20=[0.02];

K2=lsqcurvefit(@FunS2,K20,t2,d2);

mdcdt2=K2.*(C2.^2);

%%%%合并求参数

mdcdt=[mdcdt1;mdcdt2];

C=[C1;C2];

%%%% T=250C下，mnts初值P0

P0=[89 0.0328 110 2284.62585];

P=lsqcurvefit(@FunP,P0,C,mdcdt)

mdcdt3=(P(1)+P(2).*C)./(P(3)+1/P(4)./C./C)

figure

plot(C,mdcdt,'r-');

hold on;

plot(C,mdcdt3,'b*')

function f=FunS(K,t1)

f=K(1).*t1;

function f=FunS2(K2,t2)

f=K2(1).*t2;

function f=FunP(P,C)

f=(P(1)+P(2).*C)./(P(3)+1/P(4)./C./C);