##^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^Lance Simms, Stanford University 2009 #FDGrid_2d.py # #A class that creates the grid to be used for the finite differencing scheme # #The diodes lay along the y axis and the detector surface lies along the x axis # #INPUTS # TPX: int # The number of points in the p region along the x direction. The total amount # of points in the x direction for the whole detector will be 4*TPX-1 # TPY: int # The number of points along the y direction for one individual diode. # The total number of points for the y direction will be TPY*NumDiodes # Simple: int # 0 - Approximate diode entirely with lightly doped n material # 1 - Use p+ n n+ regions with guassian type grading for all columns # 2 - Approximate p+ pixels on n type bulk # Mobility: int # 0 - Use constant mobilities throughout diode # 1 - Use doping concentration dependent mobilities # BC_X: int # 0 - Periodic boundary conditions: add an extra dx between MM and 0 # N_a_max: float # Doping concentration in p region # N_d_max: float # Doping concentration in n region # N_bulk: float # Doping concentration in bulk region # #EXAMPLES # from FDGrid_2d import * # Grid2d = FDGrid_2d(30,100,1.e24,1.e18,1.e24,1, Simple=1) # ############################################# from numpy import * from pylab import * import matplotlib.pylab as mplot import string import pdb class FDGrid_2d: #INSTANTIATION###################### def __init__(self, TPX, TPY, N_a_max, N_bulk, N_d_max, \ TvFlag = 0, Simple = 0, BC_X=0 ): #Define some constants to be used by objects self.e_0 = 8.85419e-12 # F * m^-1 - 8.85419 e-14 F*cm^-1 self.e_si = 11.8 # Relative Permittivity for Si self.q = double(1.6022e-19) # Coulombs per electron self.k = double(1.38066e-23) # Joules per Kelvin self.T = double(300) # Kelvin self.Theta = double(self.q/(self.k*self.T)) # Factor relating energies self.n_i = 1e16 # intrinsic concentration m^-3 #Set up basic parameters and indices self.TPY = TPY # Total number of points in the P region and N along y self.AllPY = 4*TPY-1 # Total number of points in entire region along y self.TPX = TPX # Total number of points in entire region along x self.AllPX = TPX # Total number of points in entire region along x AllPY = self.AllPY; NN = AllPY-1 # Total number of points in entire region along y N_i = 1.e16 # Intrinsic concentration Num_P = double(TPY) # Total number of points in P+ Num_I = double(TPY) # Total number of points in I Num_N = double(TPY) # Total number of points in N+ Num_Pf = double(TPY) # Total number of points in P region Num_PIf = double(TPY/8) Num_Ic = double(3*TPY/4) Num_INf = double(TPY/8) Num_Nf = double(TPY) #Define the grid self.Dy = 49e-6 # Half the distance along y self.dy_min = 1e-9 # The smallest dy in meters self.Dx = 18e-6 # An 18 micron pixel self.dx_min = self.Dx/self.AllPX # The smallest dx in meters self.BC_X = BC_X # Periodic boundary conditions in X? self.DefineGrid2d(TvFlag=TvFlag) self.N_d = zeros([self.AllPY,self.AllPX], dtype=double) self.N_a = zeros([self.AllPY,self.AllPX], dtype=double) self.n = zeros([self.AllPY,self.AllPX], dtype=double) self.p = zeros([self.AllPY,self.AllPX], dtype=double) #DC is a 2d array that holds the doping profile self.DC = zeros([self.AllPY,self.AllPX], dtype=double) #SIMPLE = 0 =======================Approximate junctions as being abrupt ##SIMPLE =0 NOT YET IMPLEMENTED if Simple == 0: #Donor concentration self.N_d[0:TPY] = double(0) self.N_d[TPY:3*TPY-1] = double(N_bulk) self.N_d[3*TPY-1:4*TPY-1] = double(N_d_max) #Acceptor concentration self.N_a[0:TPY] = double(N_a_max) self.N_a[TPY:3*TPY-1] = double(0) self.N_a[3*TPY-1:4*TPY-1] = double(0) #Now define the hole and electron densities self.n[0:TPY] = N_i**2/N_a_max self.n[TPY:3*TPY-1] = double(N_bulk) self.n[3*TPY-1:4*TPY-1] = double(N_d_max) self.p[0:TPY] = N_a_max self.p[TPY:3*TPY-1] = N_i**2/N_bulk self.p[3*TPY-1:4*TPY-1] = N_i**2/N_d_max #Difference self.DC = self.N_d-self.N_a self.n[NN] = self.DC[NN]/2*(1+(1+(2*N_i/self.DC[NN])**2)**.5) self.p[0] = -self.DC[0]/2*(1+(1+(2*N_i/self.DC[0])**2)**.5) self.n[0] = N_i**2/self.p[0] self.p[NN] = N_i**2/self.n[NN] #SIMPLE = 1 =======================Use Gaussian-type doping function elif Simple == 1: #Grade the junction according to p. 55 of Kurata x_a = 1.e-6; #Thickness of p+ region x_b = 1.e-6; #Thickness of n+ region w = 100.e-6; #Total width of diode m1 = (1./x_a**2)*log(N_a_max/N_bulk) m2 = (1./x_b**2)*log(N_d_max/N_bulk) #Now define the hole and electron densities one column at a time for i in arange(self.TPX): #DC is a 1d column array self.DC[:,i] = -N_a_max*exp(-m1*self.Y**2) +\ N_bulk +\ N_d_max*exp(-m2*(w-self.Y)**2) self.n[0:TPY,i] = N_i**2/N_a_max self.n[TPY:3*TPY-1,i] = double(N_bulk) self.n[3*TPY-1:4*TPY-1,i] = double(N_d_max) self.p[0:TPY,i] = N_a_max self.p[TPY:3*TPY-1,i] = N_i**2/N_bulk self.p[3*TPY-1:4*TPY-1,i] = N_i**2/N_d_max self.n[NN,i] = self.DC[NN,i]/2*(1+(1+(2*N_i/self.DC[NN,i])**2)**.5) self.p[0,i] = -self.DC[0,i]/2*(1+(1+(2*N_i/self.DC[0,i])**2)**.5) self.n[0,i] = N_i**2/self.p[0,i] self.p[NN,i] = N_i**2/self.n[NN,i] #SIMPLE = 2 =======================Generate Doping profile for p+ implant on n bulk elif Simple == 2: #Grade the junction according to p. 55 of Kurata x_pp_bou = 6e-6 ; #p+ implant is 12 microns/ half is 6 y_pp_bou = 1e-6 ; #p+ implant is 1 micron/ boundary along y x_a = 1.e-6; #Thickness of p+ region x_b = 1.e-6; #Thickness of n+ region w = 100.e-6; #Total width of diode m1 = (1./x_a**2)*log(N_a_max/N_bulk) m2 = (1./x_b**2)*log(N_d_max/N_bulk) #DC is a 2d array that holds the doping profile self.DC = zeros([self.AllPY,self.AllPX], dtype=double) #Now define the hole and electron densities one column at a time for i, x in zip(arange(self.AllPX/2+1), self.X): for j, y in zip(arange(self.AllPY), self.Y): #Right Half if i+self.AllPX/2 < self.AllPX: self.DC[j,i+self.AllPX/2] = \ -N_a_max*(exp(-m1*y**2)*(1-long(x/x_pp_bou)))\ -N_a_max*(exp(-m1*y**2)*exp(-m1*(x-x_pp_bou)**2)*long(x/x_pp_bou))\ +N_bulk \ +N_d_max*exp(-m2*(w-y)**2) #Left Half if self.AllPX/2-i >= 0: self.DC[j,self.AllPX/2-i] = \ -N_a_max*(exp(-m1*y**2)*(1-long(x/x_pp_bou)))\ -N_a_max*(exp(-m1*y**2)*exp(-m1*(x-x_pp_bou)**2)*long(x/x_pp_bou))\ +N_bulk \ +N_d_max*exp(-m2*(w-y)**2) if i+self.AllPX/2 < self.AllPX: #Right Half self.n[0:TPY,i+self.AllPX/2] = N_i**2/N_a_max self.n[TPY:3*TPY-1,i+self.AllPX/2] = double(N_bulk) self.n[3*TPY-1:4*TPY-1,i+self.AllPX/2] = double(N_d_max) self.p[0:TPY,i+self.AllPX/2] = N_a_max self.p[TPY:3*TPY-1,i+self.AllPX/2] = N_i**2/N_bulk self.p[3*TPY-1:4*TPY-1,i+self.AllPX/2] = N_i**2/N_d_max self.n[NN,i+self.AllPX/2] = self.DC[NN,i+self.AllPX/2]/2*(1+(1+(2*N_i/self.DC[NN,i+self.AllPX/2])**2)**.5) self.p[0,i+self.AllPX/2] = -self.DC[0,i+self.AllPX/2]/2*(1+(1+(2*N_i/self.DC[0,i+self.AllPX/2])**2)**.5) self.n[0,i+self.AllPX/2] = N_i**2/self.p[0,i+self.AllPX/2] self.p[NN,i+self.AllPX/2] = N_i**2/self.n[NN,i+self.AllPX/2] if self.AllPX/2-i >= 0: #Left Half self.n[0:TPY,self.AllPX/2-i] = N_i**2/N_a_max self.n[TPY:3*TPY-1,self.AllPX/2-i] = double(N_bulk) self.n[3*TPY-1:4*TPY-1,self.AllPX/2-i] = double(N_d_max) self.p[0:TPY,self.AllPX/2-i] = N_a_max self.p[TPY:3*TPY-1,self.AllPX/2-i] = N_i**2/N_bulk self.p[3*TPY-1:4*TPY-1,self.AllPX/2-i] = N_i**2/N_d_max self.n[NN,self.AllPX/2-i] = self.DC[NN,self.AllPX/2-i]/2*(1+(1+(2*N_i/self.DC[NN,self.AllPX/2-i])**2)**.5) self.p[0,self.AllPX/2-i] = -self.DC[0,self.AllPX/2-i]/2*(1+(1+(2*N_i/self.DC[0,self.AllPX/2-i])**2)**.5) self.n[0,self.AllPX/2-i] = N_i**2/self.p[0,self.AllPX/2-i] self.p[NN,self.AllPX/2-i] = N_i**2/self.n[NN,self.AllPX/2-i] ######################################################################## #StretchFunction2d_y #PURPOSE: Generate a 2d grid with fine spacing near the P+ or N+ regions # and smooth transitions to the coarse spacing over the I region # # The stretching function is given by # x = Dy*(exp(ky*(I-1)/(IL-1))/(exp(ky)-1)) #INPUTS: # I - The starting index for x in the stretch # IL - The maximum index for x in the stretch # Dy - The total distance (49 microns) # dy_min - The finest spacing in the y mesh # def StretchFunction2d_y(self,I,IL,Dy,dy_min, TvFlag = 0): TPY = self.TPY #Total amount of points in P, N regions dy_min = 1e-6/TPY #Minimum spacing in mesh IFac = double((I-1))/(IL-1) #The value for the x-exponent #Initial guess for the stretching coefficient ky=1 for i in arange(10): f = dy_min-Dy*(exp(ky*IFac)-1)/(exp(ky)-1) f_prime = -Dy*((exp(ky)-1)*IFac*exp(ky*IFac)-(exp(ky*IFac)-1)*exp(ky))/\ (exp(ky)-1)**2 ky = ky - f/f_prime #Plot the solution to the grid if TvFlag == 1: i = frange(IL-I+1)+1 IFac = double((i-1)/(IL-1)) f = Dy*(exp(ky*IFac)-1)/(exp(ky)-1) print 'Maximum x :'+str(f.max()) print 'Minimum x :'+str(f.min()) mplot.figure(0) mplot.clf() mplot.vlines(f, 0, 1) #Create the coefficient as a member of the class self.ky = ky ######################################################################## #StretchFunction2d_y #PURPOSE: Generate a 2d grid with fine spacing near the P+ or N+ regions # and smooth transitions to the coarse spacing over the I region # # The stretching function is given by # x = Dy*(exp(ky*(I-1)/(IL-1))/(exp(ky)-1)) #INPUTS: # I - The starting index for x in the stretch # IJ - The maximum index for x in the stretch # Dx - The total distance (49 microns) # dx_min - The finest spacing in the x mesh # def StretchFunction2d_x(self,I,IJ,Dx,dx_min, TvFlag = 0): TPX = self.TPX #Total amount of points in P, N regions dx_min = 1e-6/TPX #Minimum spacing in mesh IFac = double((I-1))/(IJ-1) #The value for the x-exponent #Initial guess for the stretching coefficient kx=1 for i in arange(10): f = dx_min-Dx*(exp(kx*IFac)-1)/(exp(kx)-1) f_prime = -Dx*((exp(kx)-1)*IFac*exp(kx*IFac)-(exp(kx*IFac)-1)*exp(kx))/\ (exp(kx)-1)**2 kx = kx - f/f_prime #Plot the solution to the grid if TvFlag == 1: i = frange(IJ-I+1)+1 IFac = double((i-1)/(IJ-1)) f = Dx*(exp(kx*IFac)-1)/(exp(kx)-1) print 'Maximum x :'+str(f.max()) print 'Minimum x :'+str(f.min()) mplot.figure(0) mplot.clf() mplot.vlines(f, 0, 1) #Create the coefficient as a member of the class self.kx = kx def DefineGrid2d(self, TvFlag=0): TPY = self.TPY #Total amount of points in P, N regions along a column AllPY = self.AllPY #Total amount of points in entire grid along a column TPX = self.TPX #Total amount of points in P, N regions along a row AllPX = self.AllPX #Total amount of points in P, N regions along a column dy_min = 1e-6/TPY #Minimum spacing in mesh I = 2 #The second index in the mesh i = frange(TPY-I+1)+1 #indices for the stretched portion of the grid IFac = double((i-1)/(TPY-1)) #The exponent in the stretching function IFacR = IFac[::-1] #The reversed set of exponents #Generate the stretching coefficient that solves for the grid self.StretchFunction2d_y(2,TPY, self.Dy, self.dy_min, TvFlag=TvFlag) ky = self.ky #Create the actual grid values in microns Yp = frange(0,1e-6-dy_min, dy_min) Yn = frange(99e-6+dy_min,100e-6, dy_min) Ypi = self.Dy*(exp(ky*IFac)-1)/(exp(ky)-1)+1e-6 Yin = 100e-6-Ypi[::-1] Y = concatenate([Yp, Ypi, Yin[1:], Yn]) self.Y = Y self.Ypi = Ypi #KURATA self.dy = Y[1:AllPY]-Y[0:AllPY-1] #h(M) self.dy_p = (self.dy[1:AllPY-1]+self.dy[0:AllPY-2])/2 #h'(N) self.m_a_y = self.dy[1:AllPY-1]/(2*self.dy_p) #h(M)/2h'(N) self.m_b_y = self.dy[0:AllPY-2]/(2*self.dy_p) #h(M-1)/2h'(N) #Create the actual grid values in microns X = frange(0,self.Dx, self.dx_min) self.X = X self.dx = X[1:AllPX]-X[0:AllPX-1] #h(M) self.dx_p = (self.dx[1:AllPX-1]+self.dx[0:AllPX-2])/2 #h'(N) self.m_a_x = self.dx[1:AllPX-1]/(2*self.dx_p) #h(M)/2h'(N) self.m_b_x = self.dx[0:AllPX-2]/(2*self.dx_p) #h(M-1)/2h'(N) #Add an extra element for the periodic boundary condition case if self.BC_X == 1: self.dx = append(self.dx , self.dx[0]) self.dx_p = append(self.dx_p, self.dx_p[0]) self.m_a_x = append(self.m_a_x, self.m_a_x[0]) self.m_b_x = append(self.m_b_x, self.m_b_x[0])