clc
clear
close all
e=1.602176634*10^-19;
por=0.5;
S=0.88;
Cf=(4.3*10^-20)/0.0000000000000000033;
ep=8.85*(10^-12);
E0=0;
Tee=2.5;
dc=1.2*(10^-8);
Vg=5000;
d2=0.08;
Ef=Vg/d2;
row=3000;
Current=[1 2 5 8 10 16];
DLR1=[1.6 5 10 25 30 60]*10^4;
DLR1=[1 5 10 24 30 60]*10^4/(0.01^2);
A=[0.01;0.1;1];
Vcon=0.1*0.1*0.01;
g=0;
num=500;
count1=1;
for i=0.1:0.1:18
r=i*10^-6;
A1(count1,1)={r};
m=(4/3)*pi*r^3*row;
m_dot=A*m.*DLR1;
Vdust=(4/3)*pi*(r^3);
nd=(Vcon./Vdust)*por;
nd=nd';
md=Vdust*row;
mall=md*nd;
Cc=4*pi*ep*(r);
d1=200*r;
s=2*r;
et=(s/Tee)*sqrt(((2*Cf*(S^2))/(pi*ep*r))-2*E0^2);
Qm=-0.5*et.*Cc*Tee;
V1=Qm/Cc;
A2(count1,1:4)={r,Qm,Cc,V1};
count1=count1+1;
Fef1=Qm*Ef;
h=2/m;
a=Qm^2/(4*pi*ep*s);
b=2*pi*ep*r^2*E0^2;
c=m*g;
f=2*Cc*S^2*r*dc;
Vini=(h*(a+((b-c)*d1)-f))^0.5;
g=9.8;
y=Vini^2/(2*g);
D=0.05;
t=sqrt((Vini^2)-(2*g*D));
S2=(r^2)*pi;
e=1.602*10^-19;
V1=Qm/Cc;
En1=Qm*V1;
count2=1;
for V2=100:1:num
eV2=V2;
if abs(V1)<=abs(eV2);
A2(count2,1:4)={r,V1,V2,eV2};
Q_new=V2*Cc;
Vini=abs(((2/md)*(((Qm.^2)./(4*pi*ep.*s))+(2*pi*ep*(r.^2)*(E0^2)).*d1-(2*Cf*(S^2).*r*b)))^0.5);
dltV1=sqrt((2*abs(Qm)*Ef*d2)/md);
Vex1=abs(Vini+dltV1);
F1=Vex1*m_dot;
F11=F1(1,1:6);
F12=F1(2,1:6);
F13=F1(2,1:6);
AF11(count2,1:2)={F11,r};
AF12(count2,1:2)={F12,r};
AF13(count2,1:2)={F13,r};
dltV2=sqrt((2*(abs(Qm)+abs(Q_new))*Ef*d2)/md);
Vex2=abs(Vini+dltV2);
F21=Vex2*m_dot;
F211=F21(1,:);
F212=F21(2,:);
F213=F21(3,:);
AF211(count2,1)={F211};
AF212(count2,1:3)={F212,r,V2};
AF213(count2,1)={F213};
AFr(count2,1)={r};
AFr211(count2,1:3)={r,F211,V2};
AFV2(count2,1)={V2};
AFr212(count2,1:3)={r,F212,V2};
AFr213(count2,1:3) ={r,F213,V2};
count2=count2+1;
end
end
end
AFA2=cell2mat(A2);
AF11=cell2mat(AF11);
AF12=cell2mat(AF12);
AF13=cell2mat(AF13);
AF211=cell2mat(AF211);
AF212=cell2mat(AF212);
AF213=cell2mat(AF213);
AFr=cell2mat(AFr);
AFV2=cell2mat(AFV2);
AFr211=cell2mat(AFr211);
AFr212=cell2mat(AFr212);
AFr213=cell2mat(AFr213);