22 views (last 30 days)

Research area: Chemical kinetics

problem description: So this is parallel kinetics. I want to determine rate of reaction (as dCdt) and constant rate (as k) from experimental data that I have, time (as t) and concentration (as C). There are 6 variable called as mag,dag,tag,gly,ffa, and air.

so there are 6 constant rate and rate of reaction that I need to know, using ode23 or ode45. and the equations are described below:

dCdt(mag)=k1*Ctag*Cgly+2*k2*Cdag*Cgly-k3*Ctag*Cmag+k5*Cdag*Cair+k6*Cgly*Cffa;

dCdt(dag)=k1*Ctag*Cgly-k2*Cdag*Cgly+2*k3*Ctag*Cmag+k4*Ctag*Cair-k5*Cdag*Cair;

dCdt(tag)=-k1*Ctag*Cgly-k3*Ctag*Cmag-k4*Ctag*Cair;

dCdt(ffa)= k4*Ctag*Cair+k5*Cdag*Cair-k6*Cgly*Cffa;

dCdt(gly)=-k1*Ctag*Cgly-k2*Cdag*Cgly-k6*Cgly*Cffa;

dCdt(air)=-k4*Ctag*Cair-k5*Cdag*Cair+k6*Cgly*Cffa;

all I want to do is fit my data to this function to solve for k1,k2,k3,k4,k5,k6 and find each dCdt. I've tried make a code by using fminsearch and ode23, but I dont know how to define each concentration for ode so as finding each k. Can anybody help me with the code?

Thanks in advance.

function constantrate

clc

Tspan = [0.1 0.1]

constant=fminsearch(@rate, Tspan);

function constantrate = rate (k)

Data = [

60 90.5 5.2 1.1 10.5 2 3

120 90.6 5.1 1 10.1 1.8 3.3

180 90.7 5.1 1.4 12.5 1.9 4

240 90.7 5.4 2.1 10 1.7 4.5

];

t=Data(:,1); %time

Cmag=Data(:,2); %concentration of mag

Cdag=Data(:,3); %concentration of dag

Ctag=Data(:,4); %concentration of tag

Cgly=Data(:,5); %concentration of gly

Cair=Data(:,6); %concentration of air

Cffa=Data(:,7); %concentration of ffa

[tint, Cint]=ode23(@reaction,t,Cmag,[],k);

constantrate=sum((Cmag-Cint).^2);

plot (tint, Cint, 'r','o',t,Cmag,'b')

legend('calculation','experimental')

xlabel('time (min)')

ylabel('concentration (M)')

function dCdt=reaction(Ctag,Cmag,Cdag,Cgly,Cair,Cffa,k1,k2,k3,k4,k5,k6)

dCdt(mag)=k1*Ctag*Cgly+2*k2*Cdag*Cgly-k3*Ctag*Cmag+k5*Cdag*Cair+k6*Cgly*Cffa;

dCdt(dag)=k1*Ctag*Cgly-k2*Cdag*Cgly+2*k3*Ctag*Cmag+k4*Ctag*Cair-k5*Cdag*Cair;

dCdt(tag)=-k1*Ctag*Cgly-k3*Ctag*Cmag-k4*Ctag*Cair;

dCdt(ffa)= k4*Ctag*Cair+k5*Cdag*Cair-k6*Cgly*Cffa;

dCdt(gly)=-k1*Ctag*Cgly-k2*Cdag*Cgly-k6*Cgly*Cffa;

dCdt(air)=-k4*Ctag*Cair-k5*Cdag*Cair+k6*Cgly*Cffa;

Star Strider
on 24 May 2020

Star Strider
on 25 May 2020

You need to have the Global Optimization Toolbox to use the ga function. I am running 2020a, and the code I posted (for reference) ran without error.

My objective otherwise was to present a set ot parameters (and slightly revised code) that you can use as initial estimates with fminsearch or other optimization functions to estimate the parameters. Those parameters should give a reasonable fit.

I still recommend that you be certain that the order and coding of the equations in ‘kinetics’ are correct. I generally expect a much better fit.

Star Strider
on 26 May 2020

My pleasure!

This finds the best coefficients of the differential equation that minimise the norm or least-square difference between the differential equation and the data. I do not know the system you are modeling (it appears to be lipolysis or something similar), so I do not know what the parameters mean. As I remember (from undergraduate physical chemistry back when log tables and slide rules were all the rage, and computers were behemoths that ate punchcards), the order would be given by the exponents of the reactions in kinetic equations, and I do not see any exponents here (so I assume everything is first-order in the respective reactants).

So if you want to optimise or estimate the order and rate, and the differential equations are set up to do that, optimising the fit will estimate the parameters you want.

Opportunities for recent engineering grads.

Apply Today
## 0 Comments

Sign in to comment.