Can someone model a reactor for me? All data is provided- I'm desperate

Really really struggling to model my reactor for my project. I have spent days trying to do this but had no success. DESPERATE AS MY DEADLINE IS SOON. Reactions are two reversible reactions:

Ethylene+Benzene=Ethylbenzene and then Ethylene+Ethylbenzene=DiEthylbenzene

Rates are as follows: REB=((k1*p1*p2)-(k2*p3))/((1+(KB*p2)+(KE*p1)+(KEB*p3)+(KDEB*p4))^2);

RDEB=((k3*p1*p3)-(k4*p4))/((1+(KB*p2)+(KE*p1)+(KEB*p3)+(KDEB*p4))^2);

I have initial values coded as following:

%Initial values------------------------------------------------------------

ini=[673,... %T 15,... %P 0.145833333,... %Ethylene 0.729166666,... %Benzene 0,... %Ethylbenzene 0]; %DiEthylbenzene

T=x(3); p=ini(2); f10=ini(3); f20=ini(4); f30=ini(5); f40=ini(6); f0=f10+f20+f30+f40;

%Initial mole fractions---------------------------------------------------- x10=f10/f0; x20=f20/f0; x30=f30/f0; x40=f40/f0;

%Variables: Moles Ethylbenzene and DiEthylbenzene made---------------------------------- y=x(1); v=x(2);

%Mole fractions for each variable------------------------------------------ x1=(x10-v-y)/(0.875-y-v); x2=(x20-y)/(0.875-y-v); x3=(x30+y-v)/(0.875-y-v); x4=(x40+v)/(0.875-y-v);

%partial pressures--------------------------------------------------------- p1=p*x1; p2=p*x2; p3=p*x3; p4=p*x4; P=[p1;p2;p3;p4];

%Rate constants------------------------------------------------------------ k1=255.3*exp(-42496/(8.314*T)); %k1 k2=5.13*exp(-68741/(8.314*T)); %k-1 k3=52.41*exp(-53148/(8.314*T)); %k2 k4=2.667*exp(-73491/(8.314*T)); %k-2

KE=0.457*exp(7769/(8.314*T)); KB=0.00677*exp(12436/(8.314*T)); KEB=0.00319*exp(10173/(8.314*T)); KDEB=0.000854*exp(14059/(8.314*T));

Other information I have are as follows:

delH1=-105450; delH2=-57670; U=300; A=1934.215765; pb=21.025; Tc=663; Q=(((U*A)/pb)*(Tc-T));

Also these'functions that someone has done sand I have no idea what they are for or anything about them.

r(1,1)=(1/f0)*REB; r(2,1)=(1/f0)*(RDEB+REB); r(3,1)=(1/(f0*(0.875-x(1)-x(2)))*cp)*(REB*-delH1+RDEB*-delH2+Q); where cp=177850

If anyone can help me do this I would be massively appreciative. I have been trying to use ODE45 and ODE15sto solve the equations. Ultimitely I had trying been getting plots of the molar flow rates of each component (y-axis) against the length along the reactor (x axis). Thank you

[wv,xv]=ode15s('greact1',0:10:50000,x0);

M=wv(min(find(C(:,4)>90)),1); disp('Mass Catalyst:') disp(M)

[Maneet Goyal]

Hello!

This paper might help you.

Modelling and Simulation Of Benzene Alkylation Process Reactors For Production Of Ethylbenzene