Ink bottle seems right

This commit is contained in:
James E McClure 2018-01-30 11:06:49 -05:00
parent 3d795facde
commit 46c22f8ac3

View File

@ -1,4 +1,6 @@
require(ggplot2)
require("ggplot2")
GGTHEME<-theme_bw()+theme(panel.grid.major = element_blank(),panel.grid.minor = element_blank())
PI=3.14159
@ -37,30 +39,52 @@ Rb=1.2
# Length of the bulb
Hb=sqrt(Rb^2-r1^2)+sqrt(Rb^2-r2^2)
time=seq(1,100)
time=seq(0,100)
time_inlet<-time
M0_inlet<-CapVolume(r1,r1)+CylVolume(r1,time*0.03)
M1_inlet<-CapArea(r1,r1)+CylArea(r1,time*0.03)
M0_inlet<-CapVolume(r1,r1)+CylVolume(r1,time*0.01*L1)
M1_inlet<-CapArea(r1,r1)+CylArea(r1,time*0.01*L1)
M2_inlet<-2*CapArea(r1,r1)/r1+CylArea(r1,time*0.03)/r1
Jwn_inlet<-2/r1
# radius of the fluid interface while contact line is pinned
rt_pore<-r1+0.01*(Rb-r1)*time
# height of the spherical cap at the inlet side
ht_pore<-(-rt_pore+sqrt(rt_pore^2+r1^2))
ht_pore<-(rt_pore-sqrt(rt_pore^2-r1^2))
# depth of penetration into the pore
xt_pore<-2*rt_pore-ht_pore
time_pore<-time+100;
M0_pore<-CylVolume(r1,L1)+CapVolume(r1,xt_pore)
M1_pore<-CylArea(r1,L1)+CapArea(rt,xt_pore)
M2_pore<-2*CapArea(r1,r1)/r1+CylArea(r1,L1)/r1+2*CapVolume(r1,xt_pore)/rt_pore
time_pore<-time+100
M0_pore<-CylVolume(r1,L1)+pi*xt_pore*(3*r1^2+xt_pore^2)/6
M1_pore<-CylArea(r1,L1)+2*pi*rt_pore*xt_pore
M2_pore<-2*CapArea(r1,r1)/r1+CylArea(r1,L1)/r1+2*pi*xt_pore
Jwn_pore<-2/rt_pore
rt_pin<-Rb+(r2-Rb)*time*0.01
ht_pin<-rt_pin-sqrt(rt_pin^2-r2^2)
M0_pin<-CylVolume(r1,L1)+BulbVolume(r1,r2,Hb)+pi*ht_pin*(3*r2^2+ht_pin^2)/6
M1_pin<-CylArea(r1,L1)+2*pi*Rb*Hb+2*pi*rt_pin*ht_pin
M2_pin<-2*CapArea(r1,r1)/r1+CylArea(r1,L1)/r1+2*pi*Hb+2*pi*ht_pin
Jwn_pin<-2/rt_pin
time_outlet<-time_pore+100;
Jwn_outlet<-2/r2
M0_outlet<-CylVolume(r1,L1)+BubVolume(r1,r2,Hb)+CapVolume(r2,r2)+CylVolume(r2,time*0.03)
M1_outlet<-CylArea(r1,L1)+CapArea(Rb,Hb)+CapArea(r2,r2)+CylArea(r2,time*0.03)
M2_outlet<-2*CapArea(r2,r2)/r2+CylArea(r2,time*0.03)/r2
M0_outlet<-CylVolume(r1,L1)+BulbVolume(r1,r2,Hb)+CapVolume(r2,r2)+CylVolume(r2,time*0.01*L2)
M1_outlet<-CylArea(r1,L1)+CapArea(Rb,Hb)+CapArea(r2,r2)+CylArea(r2,time*0.01*L2)
M2_outlet<-2*CapArea(r2,r2)/r2+CylArea(r2,time*0.01*L2)/r2
M0<-as.vector(rbind(M0_inlet,M0_pore,M0_pin,M0_outlet))
Jwn<-as.vector(rbind(Jwn_inlet,Jwn_pore,Jwn_pin,Jwn_outlet))
Vol=max(M0)
sw<-M0/Vol
p1<-ggplot()+geom_line(aes(time,M0_inlet,color="inlet"))+
geom_line(aes(time+100,M0_pore,color="pore"))+
geom_line(aes(time+200,M0_pin,color="pin"))+
geom_line(aes(time+300,M0_outlet,color="outlet"))
p2<-ggplot()+geom_line(aes(sw,Jwn))+GGTHEME+
xlab(expression(s^w))+ylab(expression(p^c))