?  search:
summary | log | graphiclog1 | graphiclog2 | commit | commitdiff | tree | refs | edit | fork
blame | history | raw | HEAD
more fix on Ec/Ev.
[gss-tcad.git] / src / solver / mix1 / semiequ1mix.cc
blob4dc1193ae5a904947028409922f53f09cc2448c0
1 /*****************************************************************************/
2 /* 8888888 88888888 88888888 */
3 /* 8 8 8 */
4 /* 8 8 8 */
5 /* 8 88888888 88888888 */
6 /* 8 8888 8 8 */
7 /* 8 8 8 8 */
8 /* 888888 888888888 888888888 */
9 /* */
10 /* A Two-Dimensional General Purpose Semiconductor Simulator. */
11 /* */
12 /* GSS 0.4x */
13 /* Last update: March 29, 2007 */
14 /* */
15 /* Gong Ding */
16 /* gdiso@ustc.edu */
17 /* NINT, No.69 P.O.Box, Xi'an City, China */
18 /* */
19 /*****************************************************************************/
20
21 #include "zonedata.h"
22 #include "vec3.h"
23
24
25 //-----------------------------------------------------------------------------
26 // special process of electrode bcs for mix type simulation
27 //-----------------------------------------------------------------------------
28
29
30 void SMCZone::F1E_mix_ddm_ombc(int i,PetscScalar *x,PetscScalar *f, ODE_Formula &ODE_F, vector<int> & zofs, DABC &bc)
31 {
32 const VoronoiCell* pcell = pzone->davcell.GetPointer(i);
33 int z = zone_index;
34 PetscScalar e = mt->e;
35 PetscScalar kb = mt->kb;
36 PetscScalar Na = aux[i].Na;
37 PetscScalar Nd = aux[i].Nd;
38 mt->mapping(&pzone->danode[i],&aux[i],ODE_F.clock);
39 PetscScalar Vi = x[zofs[zone_index]+3*i+0]; //potential of node i
40 PetscScalar ni = x[zofs[zone_index]+3*i+1]; //electron density of node i
41 PetscScalar pi = x[zofs[zone_index]+3*i+2]; //hole density of node i
42 mt->mapping(&pzone->danode[i],&aux[i],ODE_F.clock);
43 PetscScalar nie = mt->band->nie(fs[i].T);
44 PetscScalar Nc = mt->band->Nc(fs[i].T);
45 PetscScalar Nv = mt->band->Nv(fs[i].T);
46 OhmicBC *pbc = dynamic_cast <OhmicBC * > (bc.Get_pointer(pcell->bc_index-1));
47
48 PetscScalar electron_density,hole_density;
49
50 if(Fermi) //Fermi
51 {
52 PetscScalar Ec = -(e*Vi + aux[i].affinity + mt->band->EgNarrowToEc(fs[i].T) );
53 PetscScalar Ev = -(e*Vi + aux[i].affinity - mt->band->EgNarrowToEv(fs[i].T) + mt->band->Eg(fs[i].T));
54 PetscScalar phin = pbc->Vapp;
55 PetscScalar phip = pbc->Vapp;
56 PetscScalar etan = (-e*phin-Ec)/kb/fs[i].T;
57 PetscScalar etap = (Ev+e*phip)/kb/fs[i].T;
58 f[zofs[z]+3*i+0] = Nc*fermi_half(etan) - Nv*fermi_half(etap)+ Na - Nd;
59 f[zofs[z]+3*i+1] = ni - Nc*fermi_half(etan);
60 f[zofs[z]+3*i+2] = pi - Nv*fermi_half(etap);
61 }
62 else //Boltzmann
63 {
64 f[zofs[z]+3*i+0] = Vi - kb*fs[i].T/e*asinh((Nd-Na)/(2*nie)) + kb*fs[i].T/2/e*log(Nc/Nv) + mt->band->Eg(fs[i].T)/2/e
65 + aux[i].affinity - pbc->Vapp;
66 PetscScalar electron_density,hole_density;
67 if(Na>Nd) //p-type
68 {
69 hole_density = (-(Nd-Na)+sqrt((Nd-Na)*(Nd-Na)+4*nie*nie))/2.0;
70 electron_density = nie*nie/hole_density;
71 }
72 else //n-type
73 {
74 electron_density = ((Nd-Na)+sqrt((Nd-Na)*(Nd-Na)+4*nie*nie))/2.0;
75 hole_density = nie*nie/electron_density;
76 }
77 f[zofs[z]+3*i+1] = ni - electron_density; //electron density
78 f[zofs[z]+3*i+2] = pi - hole_density; //hole density
79 }
80 }
81
82
83 void SMCZone::F1E_mix_ddm_stkbc(int i,PetscScalar *x,PetscScalar *f, ODE_Formula &ODE_F, vector<int>& zofs,DABC &bc)
84 {
85 const VoronoiCell* pcell = pzone->davcell.GetPointer(i);
86 SchottkyBC *pbc = dynamic_cast<SchottkyBC * >(bc.Get_pointer(pcell->bc_index-1));
87 PetscScalar Vi = x[zofs[zone_index]+3*i+0]; //potential of node i
88 PetscScalar ni = x[zofs[zone_index]+3*i+1]; //electron density of node i
89 PetscScalar pi = x[zofs[zone_index]+3*i+2]; //hole density of node i
90 PetscScalar L = 0.5*(pcell->ilen[0]+pcell->ilen[pcell->nb_num-1]);
91
92 //Schottky current
93 PetscScalar Fn = mt->band->SchottyJsn(ni,fs[i].T,pbc->WorkFunction-aux[i].affinity)*L;
94 PetscScalar Fp = mt->band->SchottyJsp(pi,fs[i].T,pbc->WorkFunction-aux[i].affinity)*L;
95 //Schotty Barrier Lowerring
96 PetscScalar deltaVB=mt->band->SchottyBarrierLowerring(aux[i].eps,sqrt(aux[i].Ex*aux[i].Ex+aux[i].Ey*aux[i].Ey));
97
98 f[zofs[zone_index]+3*i+0] = x[zofs[zone_index]+3*i+0] + pbc->WorkFunction - deltaVB - pbc->Vapp;
99 f[zofs[zone_index]+3*i+1] = (f[zofs[zone_index]+3*i+1]+Fn)/pcell->area;
100 f[zofs[zone_index]+3*i+2] = (f[zofs[zone_index]+3*i+2]-Fp)/pcell->area;
101
102 if(ODE_F.TimeDependent==true)
103 {
104 //second order
105 if(ODE_F.BDF_Type==BDF2&&ODE_F.BDF2_restart==false)
106 {
107 PetscScalar r = ODE_F.dt_last/(ODE_F.dt_last+ODE_F.dt);
108 PetscScalar Tn = (2-r)/(1-r)*ni-1.0/(r*(1-r))*fs[i].n+(1-r)/r*fs[i].n_last;
109 PetscScalar Tp = (2-r)/(1-r)*pi-1.0/(r*(1-r))*fs[i].p+(1-r)/r*fs[i].p_last;
110 f[zofs[zone_index]+3*i+1] += -Tn/(ODE_F.dt_last+ODE_F.dt);
111 f[zofs[zone_index]+3*i+2] += -Tp/(ODE_F.dt_last+ODE_F.dt);
112 }
113 else //first order
114 {
115 f[zofs[zone_index]+3*i+1] += -(ni-fs[i].n)/ODE_F.dt;
116 f[zofs[zone_index]+3*i+2] += -(pi-fs[i].p)/ODE_F.dt;
117 }
118 }
119 }
120
121
122 void SMCZone::F1E_mix_ddm_insulator_gate(int i,PetscScalar *x,PetscScalar *f, ODE_Formula &ODE_F, vector<int> & zofs, DABC &bc)
123 {
124 const VoronoiCell* pcell = pzone->davcell.GetPointer(i);
125 InsulatorContactBC *pbc = dynamic_cast<InsulatorContactBC * >(bc.Get_pointer(pcell->bc_index-1));
126
127 PetscScalar Vi = x[zofs[zone_index]+3*i+0]; //potential of node i
128 PetscScalar ni = x[zofs[zone_index]+3*i+1]; //electron density of node i
129 PetscScalar pi = x[zofs[zone_index]+3*i+2]; //hole density of node i
130 mt->mapping(&pzone->danode[i],&aux[i],ODE_F.clock);
131 PetscScalar L = 0.5*(pcell->ilen[0]+pcell->ilen[pcell->nb_num-1]);
132
133 PetscScalar grad_P = 0;
134
135 for(int j=0;j<pcell->nb_num;j++)
136 {
137 int nb = pcell->nb_array[j];
138 PetscScalar Vj = x[zofs[zone_index]+3*nb+0]; //potential of nb node
139 //the poisson's equation on third boundary type
140 if(j==0||j==pcell->nb_num-1)
141 {
142 PetscScalar vgate = pbc->Vapp - pbc->WorkFunction;
143 PetscScalar q = mt->e*pbc->QF; //sigma is the surface change density
144 PetscScalar Thick = pbc->Thick;
145 PetscScalar eps_ox = mt->eps0*pbc->eps;
146 PetscScalar r=q + eps_ox/Thick*vgate;
147 PetscScalar s=eps_ox/Thick;
148 grad_P += 0.5*pcell->ilen[j]*(r-0.25*s*(3*Vi+Vj));
149 }
150 }
151 f[zofs[zone_index]+3*i+0] = (f[zofs[zone_index]+3*i+0]+grad_P)/pcell->area
152 + mt->e*((pi-aux[i].Na)+(aux[i].Nd-ni));
153 f[zofs[zone_index]+3*i+1] = f[zofs[zone_index]+3*i+1]/pcell->area;
154 f[zofs[zone_index]+3*i+2] = f[zofs[zone_index]+3*i+2]/pcell->area;
155
156 if(ODE_F.TimeDependent==true)
157 {
158 //second order
159 if(ODE_F.BDF_Type==BDF2&&ODE_F.BDF2_restart==false)
160 {
161 PetscScalar r = ODE_F.dt_last/(ODE_F.dt_last+ODE_F.dt);
162 PetscScalar Tn = (2-r)/(1-r)*ni-1.0/(r*(1-r))*fs[i].n+(1-r)/r*fs[i].n_last;
163 PetscScalar Tp = (2-r)/(1-r)*pi-1.0/(r*(1-r))*fs[i].p+(1-r)/r*fs[i].p_last;
164 f[zofs[zone_index]+3*i+1] += -Tn/(ODE_F.dt_last+ODE_F.dt);
165 f[zofs[zone_index]+3*i+2] += -Tp/(ODE_F.dt_last+ODE_F.dt);
166 }
167 else //first order
168 {
169 f[zofs[zone_index]+3*i+1] += -(ni-fs[i].n)/ODE_F.dt;
170 f[zofs[zone_index]+3*i+2] += -(pi-fs[i].p)/ODE_F.dt;
171 }
172 }
173 }
174
175
176
177
178 void SMCZone::J1E_mix_ddm_ombc(int i,PetscScalar *x,Mat *jac,Mat *jtmp,ODE_Formula &ODE_F, vector<int> &zofs,DABC &bc)
179 {
180 PetscScalar A[9];
181 PetscInt index[3];
182 const VoronoiCell* pcell = pzone->davcell.GetPointer(i);
183 OhmicBC *pbc = dynamic_cast <OhmicBC * > (bc.Get_pointer(pcell->bc_index-1));
184 index[0] = zofs[zone_index]+3*i+0;
185 index[1] = zofs[zone_index]+3*i+1;
186 index[2] = zofs[zone_index]+3*i+2;
187 if(Fermi) //Fermi
188 {
189 PetscScalar e = mt->e;
190 PetscScalar kb = mt->kb;
191 PetscScalar Vt = kb*fs[i].T/e;
192 mt->mapping(&pzone->danode[i],&aux[i],ODE_F.clock);
193 PetscScalar Nc = mt->band->Nc(fs[i].T);
194 PetscScalar Nv = mt->band->Nv(fs[i].T);
195 PetscScalar Vi = x[zofs[zone_index]+3*i+0]; //potential of node i
196 PetscScalar Ec = -(e*Vi + aux[i].affinity + mt->band->EgNarrowToEc(fs[i].T) );
197 PetscScalar Ev = -(e*Vi + aux[i].affinity - mt->band->EgNarrowToEv(fs[i].T) + mt->band->Eg(fs[i].T));
198 PetscScalar phin = pbc->Vapp;
199 PetscScalar phip = pbc->Vapp;
200 PetscScalar etan = (-e*phin-Ec)/kb/fs[i].T;
201 PetscScalar etap = (Ev+e*phip)/kb/fs[i].T;
202 PetscScalar detan_dVi = 1.0/Vt;
203 PetscScalar detap_dVi = -1.0/Vt;
204 A[0] = Nc*fermi_mhalf(etan)*detan_dVi - Nv*fermi_mhalf(etap)*detap_dVi;
205 A[1] = 0;
206 A[2] = 0;
207 A[3] = -Nc*fermi_mhalf(etan)*detan_dVi;
208 A[4] = 1;
209 A[5] = 0;
210 A[6] = -Nv*fermi_mhalf(etap)*detap_dVi;
211 A[7] = 0;
212 A[8] = 1;
213 MatSetValues(*jac,3,index,3,index,A,INSERT_VALUES);
214
215 //f[zofs[z]+3*i+0] = Nc*fermi_half(etan) - Nv*fermi_half(etap)+ Na - Nd;
216 //f[zofs[z]+3*i+1] = ni - Nc*fermi_half(etan);
217 //f[zofs[z]+3*i+2] = pi - Nv*fermi_half(etap);
218 }
219 else
220 {
221 A[0] = 1; A[1] = 0; A[2] = 0;
222 A[3] = 0; A[4] = 1; A[5] = 0;
223 A[6] = 0; A[7] = 0; A[8] = 1;
224 MatSetValues(*jac,3,index,3,index,A,INSERT_VALUES);
225 }
226 }
227
228
229 void SMCZone::J1E_mix_ddm_stkbc(int i,PetscScalar *x,Mat *jac,Mat *jtmp,ODE_Formula &ODE_F, vector<int> &zofs,DABC &bc)
230 {
231 Mat3 A;
232 PetscInt index[3],col[3];
233 const VoronoiCell* pcell = pzone->davcell.GetPointer(i);
234 SchottkyBC *pbc = dynamic_cast<SchottkyBC * >(bc.Get_pointer(pcell->bc_index-1));
235 int z = zone_index;
236
237 PetscScalar VB = pbc->WorkFunction-aux[i].affinity;
238 PetscScalar deltaVB=mt->band->SchottyBarrierLowerring(aux[i].eps,sqrt(aux[i].Ex*aux[i].Ex+aux[i].Ey*aux[i].Ey));
239 PetscScalar Vi = x[zofs[z]+3*i+0]; //potential of node i
240 PetscScalar ni = x[zofs[z]+3*i+1]; //electron density of node i
241 PetscScalar pi = x[zofs[z]+3*i+2]; //hole density of node i
242
243 PetscScalar area = pcell->area;
244 PetscScalar d_Fn_dni = 0;
245 PetscScalar d_Fp_dpi = 0;
246 //--------------------------------
247 index[0] = zofs[z]+3*i+0;
248 index[1] = zofs[z]+3*i+1;
249 index[2] = zofs[z]+3*i+2;
250 //--------------------------------
251 for(int j=0;j<pcell->nb_num;j++)
252 {
253 int nb = pcell->nb_array[j];
254 col[0] = zofs[z]+3*nb+0;
255 col[1] = zofs[z]+3*nb+1;
256 col[2] = zofs[z]+3*nb+2;
257 MatGetValues(*jtmp,3,index,3,col,A.m);
258 A.m[0] = 0.0;
259 A=A/area;
260 if(j==0||j==pcell->nb_num-1)
261 {
262 d_Fn_dni += mt->band->pdSchottyJsn_pdn(ni,fs[i].T,VB-deltaVB)*0.5*pcell->ilen[j];
263 d_Fp_dpi += mt->band->pdSchottyJsp_pdp(pi,fs[i].T,VB+deltaVB)*0.5*pcell->ilen[j];
264 }
265 MatSetValues(*jac,3,index,3,col,A.m,INSERT_VALUES);
266 }
267
268 MatGetValues(*jtmp,3,index,3,index,A.m);
269 A=A/area;
270 A.m[0] = 1.0; //dfun(0)/dP(i)
271 A.m[1] = 0.0; //dfun(0)/dn(i)
272 A.m[2] = 0.0; //dfun(0)/dp(i)
273
274 A.m[4] += d_Fn_dni/area; //dfun(1)/dn(i)
275 A.m[8] += -d_Fp_dpi/area; //dfun(2)/dp(i)
276
277 if(ODE_F.TimeDependent==true)
278 {
279 //second order
280 if(ODE_F.BDF_Type==BDF2&&ODE_F.BDF2_restart==false)
281 {
282 PetscScalar r = ODE_F.dt_last/(ODE_F.dt_last+ODE_F.dt);
283 A.m[4] += -(2-r)/(1-r)/(ODE_F.dt_last+ODE_F.dt);
284 A.m[8] += -(2-r)/(1-r)/(ODE_F.dt_last+ODE_F.dt);
285 }
286 else //first order
287 {
288 A.m[4] += -1/ODE_F.dt;
289 A.m[8] += -1/ODE_F.dt;
290 }
291 }
292 MatSetValues(*jac,3,index,3,index,A.m,INSERT_VALUES);
293 }
294
295
296 void SMCZone::J1E_mix_ddm_insulator_gate(int i,PetscScalar *x, Mat *jac, Mat *jtmp, ODE_Formula &ODE_F, vector<int> & zofs, DABC &bc)
297 {
298 Mat3 A;
299 PetscInt index[3],col[3];
300 const VoronoiCell* pcell = pzone->davcell.GetPointer(i);
301 InsulatorContactBC *pbc = dynamic_cast<InsulatorContactBC * >(bc.Get_pointer(pcell->bc_index-1));
302 int z = zone_index;
303
304 PetscScalar ni = x[zofs[z]+3*i+1]; //electron density of node i
305 PetscScalar pi = x[zofs[z]+3*i+2]; //hole density of node i
306
307 PetscScalar area = pcell->area;
308 PetscScalar L = 0.5*(pcell->ilen[0]+pcell->ilen[pcell->nb_num-1]);
309 PetscScalar d_grad_P_dVi = 0;
310
311 //--------------------------------
312 index[0] = zofs[z]+3*i+0;
313 index[1] = zofs[z]+3*i+1;
314 index[2] = zofs[z]+3*i+2;
315 //--------------------------------
316 for(int j=0;j<pcell->nb_num;j++)
317 {
318 int nb = pcell->nb_array[j];
319 col[0] = zofs[z]+3*nb+0;
320 col[1] = zofs[z]+3*nb+1;
321 col[2] = zofs[z]+3*nb+2;
322 MatGetValues(*jtmp,3,index,3,col,A.m);
323 A=A/area;
324 if(j==0||j==pcell->nb_num-1)
325 {
326 PetscScalar Thick = pbc->Thick;
327 PetscScalar eps_ox = mt->eps0*pbc->eps;
328 PetscScalar s=eps_ox/Thick;
329 d_grad_P_dVi += -0.5*pcell->ilen[j]*0.25*s*3;
330 A.m[0]+= -0.5*pcell->ilen[j]*0.25*s/area;
331 }
332 MatSetValues(*jac,3,index,3,col,A.m,INSERT_VALUES);
333 }
334
335 mt->mapping(&pzone->danode[i],&aux[i],ODE_F.clock);
336
337 MatGetValues(*jtmp,3,index,3,index,A.m);
338 A=A/area;
339 A.m[0] += d_grad_P_dVi/area; //dfun(0)/dP(i)
340 A.m[1] += -mt->e; //dfun(0)/dn(i)
341 A.m[2] += mt->e; //dfun(0)/dp(i)
342
343 if(ODE_F.TimeDependent==true)
344 {
345 //second order
346 if(ODE_F.BDF_Type==BDF2&&ODE_F.BDF2_restart==false)
347 {
348 PetscScalar r = ODE_F.dt_last/(ODE_F.dt_last+ODE_F.dt);
349 A.m[4] += -(2-r)/(1-r)/(ODE_F.dt_last+ODE_F.dt);
350 A.m[8] += -(2-r)/(1-r)/(ODE_F.dt_last+ODE_F.dt);
351 }
352 else //first order
353 {
354 A.m[4] += -1/ODE_F.dt;
355 A.m[8] += -1/ODE_F.dt;
356 }
357 }
358 MatSetValues(*jac,3,index,3,index,A.m,INSERT_VALUES);
359 }
360
361
362 void SMCZone::F1E_mix_om_electrode_current(int i,PetscScalar *x,PetscScalar *f, ODE_Formula &ODE_F, vector<int> & zofs, DABC &bc, PetscScalar DeviceDepth)
363 {
364 int bc_index = electrode[i];
365 OhmicBC *pbc = dynamic_cast <OhmicBC * > (bc.Get_pointer(bc_index));
366 PetscScalar e = mt->e;
367
368 //calculate the total current of ohmic electrode
369 PetscScalar current=0;
370 for(int j=0;j<bc[bc_index].psegment->node_array.size();j++)
371 {
372 int node=bc[bc_index].psegment->node_array[j];
373 const VoronoiCell* pcell = pzone->davcell.GetPointer(node);
374 PetscScalar Vi = x[zofs[zone_index]+3*node+0]; //potential of node i
375 //conduction current
376 current += DeviceDepth*(f[zofs[zone_index]+3*node+1]-f[zofs[zone_index]+3*node+2]);
377 for(int k=0;k<pcell->nb_num;k++)
378 {
379 int nb = pcell->nb_array[k];
380 PetscScalar Vj = x[zofs[zone_index]+3*nb+0]; //potential of nb node
381 //displacement current
382 current += DeviceDepth*pcell->elen[k]*aux[node].eps*((Vi-Vj)-(fs[node].P-fs[nb].P))/pcell->ilen[k]/ODE_F.dt;
383 }
384 }
385 pbc->Set_Current_new(current);
386
387 }
388
389 void SMCZone::F1E_mix_stk_electrode_current(int i,PetscScalar *x,PetscScalar *f, ODE_Formula &ODE_F, vector<int> & zofs, DABC &bc, PetscScalar DeviceDepth)
390 {
391 int bc_index = electrode[i];
392 SchottkyBC *pbc = dynamic_cast <SchottkyBC * > (bc.Get_pointer(bc_index));
393 PetscScalar current=0;
394 for(int j=0;j<bc[bc_index].psegment->node_array.size();j++)
395 {
396 int node = bc[bc_index].psegment->node_array[j];
397 const VoronoiCell* pcell = pzone->davcell.GetPointer(node);
398 PetscScalar Vi = x[zofs[zone_index]+3*node+0]; //potential of node i
399 PetscScalar ni = x[zofs[zone_index]+3*node+1]; //electron density of node i
400 PetscScalar pi = x[zofs[zone_index]+3*node+2]; //hole density of node i
401 mt->mapping(&pzone->danode[node],&aux[node],ODE_F.clock);
402 PetscScalar deltaVB=mt->band->SchottyBarrierLowerring(aux[node].eps,sqrt(aux[node].Ex*aux[node].Ex+aux[node].Ey*aux[node].Ey));
403 PetscScalar Jsn = mt->band->SchottyJsn(ni,fs[node].T,pbc->WorkFunction-aux[node].affinity-deltaVB);
404 PetscScalar Jsp = mt->band->SchottyJsp(pi,fs[node].T,pbc->WorkFunction-aux[node].affinity+deltaVB);
405 current += -Jsn*0.5*pcell->ilen[0]*DeviceDepth;
406 current += -Jsp*0.5*pcell->ilen[0]*DeviceDepth;
407 current += -Jsn*0.5*pcell->ilen[pcell->nb_num-1]*DeviceDepth;
408 current += -Jsp*0.5*pcell->ilen[pcell->nb_num-1]*DeviceDepth;
409 for(int k=0;k<pcell->nb_num;k++)
410 {
411 int nb = pcell->nb_array[k];
412 PetscScalar Vj = x[zofs[zone_index]+3*nb+0]; //potential of nb node
413 //displacement current
414 current += DeviceDepth*pcell->elen[k]*aux[node].eps*((Vi-Vj)-(fs[node].P-fs[nb].P))/pcell->ilen[k]/ODE_F.dt;
415 }
416 }
417 pbc->Set_Current_new(current);
418
419 }
420
421
422 void SMCZone::F1E_mix_ins_electrode_current(int i,PetscScalar *x,PetscScalar *f, ODE_Formula &ODE_F, vector<int> & zofs, DABC &bc, PetscScalar DeviceDepth)
423 {
424 int bc_index = electrode[i];
425 InsulatorContactBC *pbc = dynamic_cast <InsulatorContactBC * > (bc.Get_pointer(bc_index));
426 PetscScalar current=0;
427
428 for(int j=0;j<bc[bc_index].psegment->node_array.size();j++)
429 {
430 int node = bc[bc_index].psegment->node_array[j];
431 const VoronoiCell* pcell = pzone->davcell.GetPointer(node);
432 PetscScalar Vi = x[zofs[zone_index]+3*node+0]; //potential of node i
433 for(int k=0;k<pcell->nb_num;k++)
434 {
435 int nb = pcell->nb_array[k];
436 PetscScalar Vj = x[zofs[zone_index]+3*nb+0]; //potential of nb node
437 //displacement current
438 current += DeviceDepth*pcell->elen[k]*aux[node].eps*((Vi-Vj)-(fs[node].P-fs[nb].P))/pcell->ilen[k]/ODE_F.dt;
439 }
440 }
441
442 pbc->Set_Current_new(current);
443 }
444
445
446 void SMCZone::F1E_mix_om_electrode_Load(int bc_index, double & current, PetscScalar & pdI_pdV, Vec & pdI_pdw, Vec & pdF_pdV,
447 PetscScalar *x, Mat *jtmp, ODE_Formula &ODE_F, vector<int> &zofs, DABC & bc, PetscScalar DeviceDepth)
448 {
449 OhmicBC *pbc = dynamic_cast <OhmicBC * > (bc.Get_pointer(bc_index));
450 pdI_pdV = 0;
451 PetscScalar e = mt->e;
452 PetscScalar kb = mt->kb;
453 VecZeroEntries(pdI_pdw);
454 VecZeroEntries(pdF_pdV);
455 PetscScalar * apdI_pdw;
456 PetscScalar * apdF_pdV;
457 VecGetArray(pdI_pdw,&apdI_pdw);
458 VecGetArray(pdF_pdV,&apdF_pdV);
459 current = pbc->Get_Current_new();
460 PetscScalar A1[3],A2[3];
461 PetscInt index[3],col[3];
462 for(int j=0;j<bc[bc_index].psegment->node_array.size();j++)
463 {
464 int node=bc[bc_index].psegment->node_array[j];
465 const VoronoiCell* pcell = pzone->davcell.GetPointer(node);
466 PetscScalar Vi = x[zofs[zone_index]+3*node+0]; //potential of node i
467 PetscScalar dJdisp_dVi=0,dJdisp_dVr=0;
468 index[0] = zofs[zone_index]+3*node+0;
469 index[1] = zofs[zone_index]+3*node+1;
470 index[2] = zofs[zone_index]+3*node+2;
471 for(int k=0;k<pcell->nb_num;k++)
472 {
473 int nb = pcell->nb_array[k];
474 PetscScalar Vj = x[zofs[zone_index]+3*nb+0]; //potential of nb node
475 col[0] = zofs[zone_index]+3*nb+0;
476 col[1] = zofs[zone_index]+3*nb+1;
477 col[2] = zofs[zone_index]+3*nb+2;
478
479 MatGetValues(*jtmp,1,&index[1],3,col,A1);
480 MatGetValues(*jtmp,1,&index[2],3,col,A2);
481
482 //for displacement current
483 dJdisp_dVi += aux[node].eps/pcell->ilen[k]/ODE_F.dt*pcell->elen[k];
484 dJdisp_dVr = -aux[node].eps/pcell->ilen[k]/ODE_F.dt*pcell->elen[k];
485
486 //
487 apdI_pdw[col[0]] += (A1[0]-A2[0]+dJdisp_dVr)*DeviceDepth;
488 apdI_pdw[col[1]] += (A1[1]-A2[1])*DeviceDepth;
489 apdI_pdw[col[2]] += (A1[2]-A2[2])*DeviceDepth;
490 }
491
492 MatGetValues(*jtmp,1,&index[1],3,index,A1);
493 MatGetValues(*jtmp,1,&index[2],3,index,A2);
494
495 apdI_pdw[index[0]] += (A1[0]-A2[0]+dJdisp_dVi)*DeviceDepth;
496 apdI_pdw[index[1]] += (A1[1]-A2[1])*DeviceDepth;
497 apdI_pdw[index[2]] += (A1[2]-A2[2])*DeviceDepth;
498
499 if(Fermi)
500 {
501 PetscScalar Vt = kb*fs[node].T/e;
502 mt->mapping(&pzone->danode[node],&aux[node],ODE_F.clock);
503 PetscScalar Nc = mt->band->Nc(fs[node].T);
504 PetscScalar Nv = mt->band->Nv(fs[node].T);
505 PetscScalar Vi = x[zofs[zone_index]+3*node+0]; //potential of node
506 PetscScalar Ec = -(e*Vi + aux[node].affinity + mt->band->EgNarrowToEc(fs[node].T) + kb*fs[node].T*log(aux[node].Nc));
507 PetscScalar Ev = -(e*Vi + aux[node].affinity - mt->band->EgNarrowToEv(fs[node].T) - kb*fs[node].T*log(aux[node].Nv)+ aux[node].Eg);
508 PetscScalar phin = pbc->Vapp;
509 PetscScalar phip = pbc->Vapp;
510 PetscScalar etan = (-e*phin-Ec)/kb/fs[node].T;
511 PetscScalar etap = (Ev+e*phip)/kb/fs[node].T;
512 PetscScalar detan_dVi = 1.0/Vt;
513 PetscScalar detap_dVi = -1.0/Vt;
514 apdF_pdV[index[0]] = Nc*fermi_mhalf(etan)*detan_dVi - Nv*fermi_mhalf(etap)*detap_dVi;
515 apdF_pdV[index[1]] = -Nc*fermi_mhalf(etan)*detan_dVi;
516 apdF_pdV[index[2]] = -Nv*fermi_mhalf(etap)*detap_dVi;
517 }
518 else
519 apdF_pdV[index[0]] = 1.0;
520 }
521 VecRestoreArray(pdI_pdw,&apdI_pdw);
522 VecRestoreArray(pdF_pdV,&apdF_pdV);
523 }
524
525
526 void SMCZone::F1E_mix_stk_electrode_Load(int bc_index, double & current, PetscScalar & pdI_pdV, Vec & pdI_pdw, Vec & pdF_pdV,
527 PetscScalar *x, Mat *jtmp, ODE_Formula &ODE_F, vector<int> &zofs, DABC & bc, PetscScalar DeviceDepth)
528 {
529 SchottkyBC *pbc = dynamic_cast <SchottkyBC * > (bc.Get_pointer(bc_index));
530 pdI_pdV = 0;
531 PetscScalar e = mt->e;
532 VecZeroEntries(pdI_pdw);
533 VecZeroEntries(pdF_pdV);
534 PetscScalar * apdI_pdw;
535 PetscScalar * apdF_pdV;
536 VecGetArray(pdI_pdw,&apdI_pdw);
537 VecGetArray(pdF_pdV,&apdF_pdV);
538 current = pbc->Get_Current_new();
539
540 for(int j=0;j<bc[bc_index].psegment->node_array.size();j++)
541 {
542 int node = bc[bc_index].psegment->node_array[j];
543 const VoronoiCell* pcell = pzone->davcell.GetPointer(node);
544 PetscScalar VB = pbc->WorkFunction-aux[node].affinity;
545 PetscScalar deltaVB=mt->band->SchottyBarrierLowerring(aux[node].eps,sqrt(aux[node].Ex*aux[node].Ex+aux[node].Ey*aux[node].Ey));
546 PetscScalar ni = x[zofs[zone_index]+3*node+1]; //electron density of node i
547 PetscScalar pi = x[zofs[zone_index]+3*node+2]; //hole density of node i
548 PetscScalar dJ_dni = -mt->band->pdSchottyJsn_pdn(ni,fs[node].T,VB-deltaVB)*0.5*pcell->ilen[0]
549 -mt->band->pdSchottyJsn_pdn(ni,fs[node].T,VB-deltaVB)*0.5*pcell->ilen[pcell->nb_num-1];
550 PetscScalar dJ_dpi = -mt->band->pdSchottyJsp_pdp(pi,fs[node].T,VB+deltaVB)*0.5*pcell->ilen[0]
551 -mt->band->pdSchottyJsp_pdp(pi,fs[node].T,VB+deltaVB)*0.5*pcell->ilen[pcell->nb_num-1];
552
553 //for displacement current
554 for(int k=0;k<pcell->nb_num;k++)
555 {
556 int nb = pcell->nb_array[k];
557 PetscScalar dI_dVi = DeviceDepth*pcell->elen[k]*aux[node].eps/pcell->ilen[k]/ODE_F.dt;
558 PetscScalar dI_dVr = -DeviceDepth*pcell->elen[k]*aux[node].eps/pcell->ilen[k]/ODE_F.dt;
559 apdI_pdw[zofs[zone_index]+3*node+0] += dI_dVi;
560 apdI_pdw[zofs[zone_index]+3*nb+0] += dI_dVr;
561 }
562 apdI_pdw[zofs[zone_index]+3*node+1] = DeviceDepth*dJ_dni;
563 apdI_pdw[zofs[zone_index]+3*node+2] = DeviceDepth*dJ_dpi;
564
565 apdF_pdV[zofs[zone_index]+3*node+0] = 1.0;
566 }
567
568 VecRestoreArray(pdI_pdw,&apdI_pdw);
569 VecRestoreArray(pdF_pdV,&apdF_pdV);
570 }
571
572
573 void SMCZone::F1E_mix_ins_electrode_Load(int bc_index, double & current, PetscScalar & pdI_pdV, Vec & pdI_pdw, Vec & pdF_pdV,
574 PetscScalar *x, Mat *jtmp, ODE_Formula &ODE_F, vector<int> &zofs, DABC & bc, PetscScalar DeviceDepth)
575 {
576 InsulatorContactBC *pbc = dynamic_cast <InsulatorContactBC * > (bc.Get_pointer(bc_index));
577 pdI_pdV = 0;
578 PetscScalar e = mt->e;
579 VecZeroEntries(pdI_pdw);
580 VecZeroEntries(pdF_pdV);
581 PetscScalar * apdI_pdw;
582 PetscScalar * apdF_pdV;
583 VecGetArray(pdI_pdw,&apdI_pdw);
584 VecGetArray(pdF_pdV,&apdF_pdV);
585 current = pbc->Get_Current_new();
586
587 for(int j=0;j<bc[bc_index].psegment->node_array.size();j++)
588 {
589 int node = bc[bc_index].psegment->node_array[j];
590 const VoronoiCell* pcell = pzone->davcell.GetPointer(node);
591
592 //for displacement current
593 for(int k=0;k<pcell->nb_num;k++)
594 {
595 int nb = pcell->nb_array[k];
596 PetscScalar dI_dVi = DeviceDepth*pcell->elen[k]*aux[node].eps/pcell->ilen[k]/ODE_F.dt;
597 PetscScalar dI_dVr = -DeviceDepth*pcell->elen[k]*aux[node].eps/pcell->ilen[k]/ODE_F.dt;
598 apdI_pdw[zofs[zone_index]+3*node+0] += dI_dVi;
599 apdI_pdw[zofs[zone_index]+3*nb+0] += dI_dVr;
600 if(k==0||k==pcell->nb_num-1)
601 {
602 PetscScalar Thick = pbc->Thick;
603 PetscScalar eps_ox = mt->eps0*pbc->eps;
604 PetscScalar s=eps_ox/Thick;
605 apdF_pdV[zofs[zone_index]+3*node+0] += 0.5*pcell->ilen[k]*s/pcell->area;
606 }
607 }
608 }
609
610 VecRestoreArray(pdI_pdw,&apdI_pdw);
611 VecRestoreArray(pdF_pdV,&apdF_pdV);
612 }
613
General-purpose Semiconductor Simulator