[gyfgyf]

for var1=0:0.001:15
var2=var1/R1;
[x,y]=ode23('fun0088',[0,60],[var1;var2]);
pp=(-p0.*(1+R1)+sqrt(p0^2.*(1+R1)^2+4.*R1.*y(:,1).*(p0+R2.*y(:,1)+p0.*R2)))/(2.*R1);
ppp=(-p0.*(1+R2)+sqrt(p0^2.*(1+R2)^2+4.*R2.*pp.*(p0+R1.*pp+p0.*R1)))/(2.*R2);
pppp=ppp*R2;
if abs(mean((y(:,2)-pppp)))<1;
break
end
end