Here are a few examples on Mathematica syntax to help with your work: You 
need to learn how to set up the ODE solver. Good luck.

Parameter assignment:

    alphal = 89.47;
    alphar = 28.46;
    cas = 0.01202;
    cap = 0.02707;
    cl = 0.01289;
    cr = 0.06077;
    cvp = 0.1204;
    cvs = 0.6630;
    Rp = 1.965;
    gammal = 37.33;
    gammar = 11.88;


Expressions:

    f[s_, p_] := 0.5(s + p) - 0.5((p - s)^2 + 0.001)^(0.5); 
    Ql[t_] := H[t]*(al[t]*f[Sl[t], Pas[t]]*cl*Pvp[t])/
       (al[t]*Pas[t] + kl[t]*f[Sl[t], Pas[t]]); 
    Qr[t_] := H[t]*(ar[t]*f[Sr[t], Pap[t]]*cr*Pvs[t])/
       (ar[t]*Pap[t] + kr[t]*f[Sr[t], Pap[t]]); 
    Pap[t_] := 
     1/cap*(V - cas*Pas[t] - cvs*Pvs[t] - cvp*Pvp[t]); 

    td[t_] := Sqrt[60/H[t]]*(Sqrt[60/H[t]] - 0.4); 
    kl[t_] := Exp[-td[t]/(cl*Rl)];
    al[t_] := 1 - kl[t]; 
    kr[t_] := Exp[-td[t]/(cr*Rr)];
    ar[t_] := 1 - kr[t]; 

    Vstrl[t_] := (al[t]*f[Sl[t], Pas[t]]*cl*Pvp[t])/
       (al[t]*Pas[t] + kl[t]*f[Sl[t], Pas[t]]); 
    Vstrr[t_] := (ar[t]*f[Sr[t], Pap[t]]*cr*Pvs[t])/
       (ar[t]*Pap[t] + kr[t]*f[Sr[t], Pap[t]]); 
    
    Rs[t_] := Apesk[t]*CvO2[t];
    Fs[t_] := (Pas[t] - Pvs[t])/Rs[t]; 
    Fp[t_] := (Pap[t] - Pvp[t])/Rp;