提交: cef3a753d86068239da13ef704b67722e66848e7 - myprojects-hg-reborn (hg)

Commit MetaInfo

• 修订版: cef3a753d86068239da13ef704b67722e66848e7 (tree)
• 时间: 2008-08-15 19:40:52
• 作者: iselllo
• Commiter: iselllo

Log Message

I modified the code adding the possibility of entering by hand the
radius of gyration of certain small clusters. I also am first defining
the number of monomers and THEN the volume grid (otherwise, the
statement if n_mon=2. then... can fail since n_mon could be slightly
different from two).

更改概述

• modified: Python-codes/smoluchowski_adjusted_kernel.py

差异

diff -r 0576b292dec2 -r cef3a753d860 Python-codes/smoluchowski_adjusted_kernel.py
--- a/Python-codes/smoluchowski_adjusted_kernel.py	Tue Aug 12 19:50:47 2008 +0000
+++ b/Python-codes/smoluchowski_adjusted_kernel.py	Fri Aug 15 10:40:52 2008 +0000
@@ -13,9 +13,13 @@
 
 print "the monodisperse volume is, ", v_mono
 
-x=s.logspace(s.log10(v_mono), s.log10(v_mono*2500.),200)   # volume grid
 
-n_mon=x/v_mono #volume of solids expressed in terms of number of monomers
+n_mon=s.logspace(s.log10(1.), s.log10(2500.),200)   # volume grid
+
+x=n_mon*v_mono #volume of solids expressed in terms of number of monomers
+
+
+
 
 print "n_mon is, ", n_mon
 
@@ -53,7 +57,7 @@
 
 rho=0.01  #monomer density in the box
 
-choice_kernel= 5  #0: standard continuum kernel (and standard diffusion_coefficient)
+choice_kernel= 10  #0: standard continuum kernel (and standard diffusion_coefficient)
                   #1: constant kernel
                   #2: continuum kernel, 2 df's and special treatment of monomers (special diffusion coeff)
                   #3: continuum kernel, 1 df and special treatment of monomers (special diffusion coeff)
@@ -69,10 +73,11 @@
                   #(standard diffusion coeff)
                   #9: continuum kernel, 2 df's, Fuchs beta without special treatment of monomers, but
                   # Fuchs correction is applied to monomers only.
+                  #10: as in 5, but this time we introduce by hand the radius of gyration of the smallest
+                  #clusters
                   
 
 
-
 choice_slope = 2 #2: 2 different df's
                  #1: 1 different df
 
@@ -624,6 +629,102 @@
 	return kern_mat
 
 
+
+def Brow_ker_cont_optim_diffusion_adjusted_and_monomer_beta_fuchs_special_rules(Vlist):
+	#monomer volume
+	#r_mon=0.5 #monomer radius
+	#v_mon=4./3.*s.pi*r_mon**3.
+	#v_mono already defined as a global variable 
+	#n_mon=Vlist/v_mono #number of monomers in each aggregate
+
+	#print "n_mon is, ", n_mon
+	
+	Diff=k_B*T_0/(n_mon*m_mon*beta) #vector with the cluster diffusion coefficients
+	#which just depends on the cluster size.
+
+	R_list=s.zeros(len(Vlist)) #I initialize the array which will contain the
+	#radia of gyration
+
+	#threshold=15.
+	
+# 	small=s.where(n_mon<=threshold)
+# 	large=s.where(n_mon>threshold)
+
+# 	a_small=0.36
+# 	df_small=2.19
+
+# 	a_large=0.241
+# 	df_large=1.62
+
+
+	R_list[small]=a_small*n_mon[small]**(1./df_small)
+
+	R_list[large]=a_large*n_mon[large]**(1./df_large)
+
+        #Now I actually introduce some special rules for the radius of gyration of the smallest clusters
+
+        sel=s.where(n_mon == 1.)
+
+        R_list[sel]=3.87e-1
+
+
+        sel=s.where(n_mon == 2.)
+
+        R_list[sel]=6.39e-1
+
+        sel=s.where(n_mon == 3.)
+
+        R_list[sel]=7.18e-1
+
+        sel=s.where(n_mon == 4.)
+
+        R_list[sel]=7.48e-1
+
+
+        print "len(R_list is, )", len(R_list)
+        print "R_list is, ", R_list
+
+
+        
+
+#	R_list[0]=0.5   #special case for the monomer radius
+
+#	m=rho_p*Vlist # this holds in general (i.e. for Vlist !=3.)
+	## as long as Vlist is the volume of solid
+	##and not the space occupied by the agglomerate structure
+
+        #m=Vlist #since rho = 1.
+
+	c=(8.*k_B*T_0/(s.pi*m))**0.5
+	#print 'c is', c
+	l=8.*Diff/(s.pi*c)
+	#print 'l is', l
+	diam_seq=2.*R_list
+	
+	g=1./(3.*diam_seq*l)*((diam_seq+l)**3.-(diam_seq**2.+l**2.)**1.5)-diam_seq
+	
+	beta_fuchs=(\
+	(diam_seq[:,s.newaxis]+diam_seq[s.newaxis,:])   \
+	/(diam_seq[:,s.newaxis]+diam_seq[s.newaxis,:]\
+	+2.*(g[:,s.newaxis]**2.+g[s.newaxis,:]**2.)**0.5)\
+	+ 8.*(Diff[:,s.newaxis]+Diff[s.newaxis,:])/\
+	((c[:,s.newaxis]**2.+c[s.newaxis,:]**2.)**0.5*\
+	(diam_seq[:,s.newaxis]+diam_seq[s.newaxis,:]))\
+	)**(-1.)
+	
+	
+	## now I have all the bits for the kernel matrix
+#	kern_mat=Brow_ker_cont_optim(Vlist)*beta
+
+
+	kern_mat=4.*s.pi*(Diff[:,s.newaxis]+Diff[s.newaxis,:])* \
+		  (R_list[:,s.newaxis]+R_list[s.newaxis,:])*beta_fuchs
+
+
+	return kern_mat
+
+
+
     
 def beta_fuchs_standalone(Vlist):
 
@@ -751,6 +852,8 @@
 elif (choice_kernel == 9):
     kernel=Brow_ker_cont_optim_diffusion_adjusted_and_monomer_beta_fuchs_single_value(x)
 
+elif (choice_kernel ==10):
+    kernel=Brow_ker_cont_optim_diffusion_adjusted_and_monomer_beta_fuchs_special_rules(x)
 
 
 f_garrick=mycount_garrick(x)
@@ -996,11 +1099,12 @@
 #I have to re-define a number of important quantities
 
 
-x=s.linspace(v_mono, (v_mono*2500.), 2500)  #new volume grid
+n_mon=s.linspace(1., 2500., 2500)  #new volume grid
 
-#I have to redefine this
 
-n_mon=x/v_mono #volume of solids expressed in terms of number of monomers
+x=n_mon*v_mono  #new volume grid
+
+
 
 print "n_mon is, ", n_mon
 
@@ -1046,6 +1150,8 @@
 elif (choice_kernel == 9):
     kernel=Brow_ker_cont_optim_diffusion_adjusted_and_monomer_beta_fuchs_single_value(x)
 
+elif (choice_kernel ==10):
+    kernel=Brow_ker_cont_optim_diffusion_adjusted_and_monomer_beta_fuchs_special_rules(x)
 
 
 y_new = si.odeint(coupling_optim, y0, \