Commit: a06fe694ca95feda2bd31c7a8f52e353897fbe28 - myprojects-hg-reborn (hg)

Commit MetaInfo

• Revisão: a06fe694ca95feda2bd31c7a8f52e353897fbe28 (tree)
• Hora: 2009-04-23 01:45:58
• Autor: iselllo
• Commiter: iselllo

Mensagem de Log

I added a code to investigate numerically the condensation equation. The space-dependent part is
treated using a finite-difference scheme (central derivative) and the equation is translated into
a system of coupled ODE's.

Mudança Sumário

• add: Python-codes/growth-semi-difference.py

Diff

diff -r 3d61cd0fb648 -r a06fe694ca95 Python-codes/growth-semi-difference.py
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/Python-codes/growth-semi-difference.py	Wed Apr 22 16:45:58 2009 +0000
@@ -0,0 +1,83 @@
+#! /usr/bin/env python
+
+import scipy as s
+import numpy as n
+import pylab as p
+import scipy.integrate as si
+
+
+def simple_coupling(y,t):
+    coupling=s.zeros(len(y))
+
+    y_new=y/x_grid
+    
+    coupling[1:-1]=y_new[2:]-y_new[:-2]
+    coupling=-coupling*A/(2.*dx)
+    return coupling
+
+A=0.00016358757687
+
+#fix grid in x (careful to avoid zero!)
+
+NN=2000
+
+x_grid=s.linspace(0.01,3.,NN)
+print "x_grid is, ", x_grid
+
+dx=x_grid[2]-x_grid[1]
+
+print "dx is, ", dx
+
+#now create a constant initial state
+
+my_sel=s.where((x_grid>=0.4) & (x_grid<=0.5))
+
+# print "x_grid[my_sel] is, ", x_grid[my_sel]
+
+# y0=s.zeros(NN)
+# y0[my_sel]=1000.
+
+y0=s.exp(-(x_grid-0.5)**2./0.02)
+
+print "y0 is, ", y0
+
+t=s.linspace(0.,3500., 60)
+
+print "t is, ", t
+
+# y = integrate.odeint(coupling_optim, y0, \
+# 	t,printmessg=1,rtol=1e-10,atol=1e-10)
+
+
+y = si.odeint(simple_coupling, y0, \
+	t,printmessg=1, mxstep=600)
+
+#p.save("evolution_of_the_population.dat",y)
+
+print "y[20,:] is, ", y[20,:]
+
+fig = p.figure()
+axes = fig.gca()
+
+# p.semilogx(dp,my_initial_state,"ro")
+# p.semilogx(dp[selection],grown_state[selection],"kx")
+p.plot(x_grid,y0,"ro")
+p.plot(x_grid,y[(len(t)-1),:],"b")
+
+p.legend(("Initial state", "10 minutes"))
+
+
+p.subplots_adjust(bottom=0.2)
+
+
+
+
+p.xlabel(r'$d_p [\mu m]$', fontsize=20)
+p.ylabel(r'$n(d_p,t)$', fontsize=20)
+p.savefig("evolution_semi_difference.pdf")
+ 
+p.clf()    
+
+
+
+print "So far so good"