1 /*****************************************************************************/
2 /* 8888888 88888888 88888888 */
5 /* 8 88888888 88888888 */
8 /* 888888 888888888 888888888 */
10 /* A Two-Dimensional General Purpose Semiconductor Simulator. */
13 /* Last update: May 11, 2005 */
15 /* Gong Ding gdiso@ustc.edu */
16 /* Xuan Chun xiaomoyu505@163.com */
17 /* NINT, No.69 P.O.Box, Xi'an City, China */
19 /*****************************************************************************/
23 //-----------------------------------------------------------------------------
24 // special process of electrode bcs for mix type simulation
25 //-----------------------------------------------------------------------------
27 void ISZone::F2_mix_ddm_gatebc_segment(int i
,PetscScalar
*x
,PetscScalar
*f
, ODE_Formula
&ODE_F
, vector
<int> & zofs
,DABC
&bc
)
29 const VoronoiCell
* pcell
= &pzone
->davcell
[i
];
30 GateBC
*pbc
= dynamic_cast<GateBC
* >(bc
.Get_pointer(pcell
->bc_index
-1));
31 PetscScalar Ti
= x
[zofs
[zone_index
]+2*i
+1];
33 mt
->mapping(&pzone
->danode
[i
],&aux
[i
],ODE_F
.clock
);
34 for(int j
=0;j
<pcell
->nb_num
;j
++)
36 int nb
= pcell
->nb_array
[j
];
37 PetscScalar Tj
= x
[zofs
[zone_index
]+2*nb
+1]; //potential of nb node
38 PetscScalar T_mid
= 0.5*(Tj
+Ti
);
39 PetscScalar dTdx_mid
= (Tj
-Ti
)/pcell
->ilen
[j
];
40 PetscScalar kapa
= mt
->thermal
->HeatConduction(T_mid
);
41 grad_T
+= kapa
*pcell
->elen
[j
]*dTdx_mid
/pcell
->area
;
42 if(j
==0||j
==pcell
->nb_num
-1)
44 PetscScalar h
= pbc
->heat_transfer
;
45 PetscScalar r
= h
*pbc
->T_external
;
46 grad_T
+= 0.5*pcell
->ilen
[j
]*(r
-0.25*h
*(3*Ti
+Tj
))/pcell
->area
;
50 f
[zofs
[zone_index
]+2*i
+0] = x
[zofs
[zone_index
]+2*i
+0] - pbc
->Vapp
+ pbc
->WorkFunction
;
51 f
[zofs
[zone_index
]+2*i
+1] = grad_T
;
53 if(ODE_F
.TimeDependent
==true)
56 if(ODE_F
.BDF_Type
==BDF2
&&ODE_F
.BDF2_restart
==false)
58 PetscScalar r
= ODE_F
.dt_last
/(ODE_F
.dt_last
+ODE_F
.dt
);
59 PetscScalar Tt
= (2-r
)/(1-r
)*Ti
-1.0/(r
*(1-r
))*fs
[i
].T
+(1-r
)/r
*fs
[i
].T_last
;
60 PetscScalar HeatCapacity
= mt
->thermal
->HeatCapacity(Ti
);
61 f
[zofs
[zone_index
]+2*i
+1] += -aux
[i
].density
*HeatCapacity
*Tt
/(ODE_F
.dt_last
+ODE_F
.dt
);
65 PetscScalar HeatCapacity
= mt
->thermal
->HeatCapacity(Ti
);
66 f
[zofs
[zone_index
]+2*i
+1] += -aux
[i
].density
*HeatCapacity
*(Ti
-fs
[i
].T
)/ODE_F
.dt
;
70 //just for fill the extra equ entry
72 int size
= pzone
->davcell
.size();
74 for(int j
=0;j
<electrode
.size();j
++)
76 if(electrode
[j
]==pcell
->bc_index
-1) {gate_equ
=j
;break;}
78 f
[zofs
[zone_index
]+equ_num
*size
+gate_equ
]=0.0;
82 void ISZone::F2_mix_ddm_gatebc_interface(int i
,PetscScalar
*x
,PetscScalar
*f
, ODE_Formula
&ODE_F
, vector
<int> & zofs
,DABC
&bc
,
85 const VoronoiCell
* pcell
= &pzone
->davcell
[i
];
86 GateBC
*pbc
= dynamic_cast<GateBC
* >(bc
.Get_pointer(pcell
->bc_index
-1));
87 PetscScalar Ti
= x
[zofs
[zone_index
]+2*i
+1];
89 mt
->mapping(&pzone
->danode
[i
],&aux
[i
],ODE_F
.clock
);
90 for(int j
=0;j
<pcell
->nb_num
;j
++)
92 int nb
= pcell
->nb_array
[j
];
93 PetscScalar Tj
= x
[zofs
[zone_index
]+2*nb
+1]; //potential of nb node
94 PetscScalar T_mid
= 0.5*(Tj
+Ti
);
95 PetscScalar dTdx_mid
= (Tj
-Ti
)/pcell
->ilen
[j
];
96 PetscScalar kapa
= mt
->thermal
->HeatConduction(T_mid
);
97 grad_T
+= kapa
*pcell
->elen
[j
]*dTdx_mid
;
99 const VoronoiCell
* ncell
= pz
->pzone
->davcell
.GetPointer(n
);
100 pz
->mt
->mapping(&pz
->pzone
->danode
[n
],&pz
->aux
[n
],ODE_F
.clock
);
101 for(int j
=0;j
<ncell
->nb_num
;j
++)
103 int nb
= ncell
->nb_array
[j
];
104 PetscScalar Tj_n
= x
[zofs
[pz
->pzone
->zone_index
]+2*nb
+1]; //potential of nb node
105 PetscScalar kapa
= pz
->mt
->thermal
->HeatConduction(0.5*(Ti
+Tj_n
));
106 grad_T
+= kapa
*ncell
->elen
[j
]/ncell
->ilen
[j
]*(Tj_n
-Ti
);
110 int size
= pzone
->davcell
.size();
112 for(int j
=0;j
<electrode
.size();j
++)
114 if(electrode
[j
]==pcell
->bc_index
-1) {gate_equ
=j
;break;}
116 f
[zofs
[zone_index
]+2*i
+0] = x
[zofs
[zone_index
]+2*i
+0] - pbc
->Vapp
+ pbc
->WorkFunction
;
117 f
[zofs
[zone_index
]+2*i
+1] = grad_T
;
119 if(ODE_F
.TimeDependent
==true)
122 if(ODE_F
.BDF_Type
==BDF2
&&ODE_F
.BDF2_restart
==false)
124 PetscScalar r
= ODE_F
.dt_last
/(ODE_F
.dt_last
+ODE_F
.dt
);
125 PetscScalar Tt
= (2-r
)/(1-r
)*Ti
-1.0/(r
*(1-r
))*fs
[i
].T
+(1-r
)/r
*fs
[i
].T_last
;
126 PetscScalar HeatCapacity1
= mt
->thermal
->HeatCapacity(Ti
);
127 PetscScalar HeatCapacity2
= pz
->mt
->thermal
->HeatCapacity(Ti
);
128 f
[zofs
[zone_index
]+2*i
+1] += -aux
[i
].density
*HeatCapacity1
*Tt
/(ODE_F
.dt_last
+ODE_F
.dt
)*pcell
->area
129 -pz
->aux
[n
].density
*HeatCapacity2
*Tt
/(ODE_F
.dt_last
+ODE_F
.dt
)*ncell
->area
;
133 PetscScalar HeatCapacity1
= mt
->thermal
->HeatCapacity(Ti
);
134 PetscScalar HeatCapacity2
= pz
->mt
->thermal
->HeatCapacity(Ti
);
135 f
[zofs
[zone_index
]+2*i
+1] += -aux
[i
].density
*HeatCapacity1
*(Ti
-fs
[i
].T
)/ODE_F
.dt
*pcell
->area
136 -pz
->aux
[n
].density
*HeatCapacity2
*(Ti
-fs
[i
].T
)/ODE_F
.dt
*ncell
->area
;
140 //just for fill the extra equ entry
141 f
[zofs
[zone_index
]+equ_num
*size
+gate_equ
]=0.0;
145 void ISZone::F2_mix_gate_electrode_current(int i
,PetscScalar
*x
,PetscScalar
*f
, ODE_Formula
&ODE_F
, vector
<int> & zofs
, DABC
&bc
, PetscScalar DeviceDepth
)
147 int bc_index
= electrode
[i
];
148 GateBC
*pbc
= dynamic_cast <GateBC
* > (bc
.Get_pointer(bc_index
));
149 PetscScalar current
=0;
150 for(int j
=0;j
<bc
[bc_index
].psegment
->node_array
.size();j
++)
152 int node
=bc
[bc_index
].psegment
->node_array
[j
];
153 const VoronoiCell
* pcell
= pzone
->davcell
.GetPointer(node
);
154 PetscScalar Vi
= x
[zofs
[zone_index
]+2*node
+0]; //potential of node i
156 for(int k
=0;k
<pcell
->nb_num
;k
++)
158 int nb
= pcell
->nb_array
[k
];
159 PetscScalar Vj
= x
[zofs
[zone_index
]+2*nb
+0]; //potential of nb node
160 //displacement current
161 current
+= DeviceDepth
*pcell
->elen
[k
]*aux
[node
].eps
*((Vi
-Vj
)-(fs
[node
].P
-fs
[nb
].P
))/pcell
->ilen
[k
]/ODE_F
.dt
;
164 pbc
->Set_Current_new(current
);
167 //-------------------------------------------------------------------------------------------------
169 void ISZone::J2_mix_ddm_gatebc_segment(int i
,PetscScalar
*x
,Mat
*jac
, ODE_Formula
&ODE_F
, vector
<int> & zofs
,DABC
&bc
)
171 const VoronoiCell
* pcell
= &pzone
->davcell
[i
];
172 GateBC
*pbc
= dynamic_cast<GateBC
* >(bc
.Get_pointer(pcell
->bc_index
-1));
173 PetscScalar Ti
= x
[zofs
[zone_index
]+2*i
+1];
174 PetscScalar d_grad_T_dTi
=0;
175 mt
->mapping(&pzone
->danode
[i
],&aux
[i
],ODE_F
.clock
);
176 for(int j
=0;j
<pcell
->nb_num
;j
++)
178 int nb
= pcell
->nb_array
[j
];
180 PetscScalar Tj
= x
[zofs
[zone_index
]+2*nb
+1]; //potential of nb node
181 PetscScalar T_mid
= 0.5*(Tj
+Ti
);
182 PetscScalar kapa
= mt
->thermal
->HeatConduction(T_mid
);
183 PetscScalar d_grad_T_dTj
= kapa
*pcell
->elen
[j
]/pcell
->ilen
[j
]/pcell
->area
;
184 d_grad_T_dTi
+= -kapa
*pcell
->elen
[j
]/pcell
->ilen
[j
]/pcell
->area
;
185 if(j
==0||j
==pcell
->nb_num
-1)
187 PetscScalar h
= pbc
->heat_transfer
;
188 d_grad_T_dTi
+= -0.5*pcell
->ilen
[j
]*0.25*h
*3/pcell
->area
;
189 d_grad_T_dTj
+= -0.5*pcell
->ilen
[j
]*0.25*h
/pcell
->area
;
191 MatSetValue(*jac
,zofs
[zone_index
]+2*i
+1,zofs
[zone_index
]+2*nb
+1,d_grad_T_dTj
,INSERT_VALUES
);
193 MatSetValue(*jac
,zofs
[zone_index
]+2*i
+0,zofs
[zone_index
]+2*i
+0,1.0,INSERT_VALUES
);
195 if(ODE_F
.TimeDependent
==true)
198 if(ODE_F
.BDF_Type
==BDF2
&&ODE_F
.BDF2_restart
==false)
200 PetscScalar r
= ODE_F
.dt_last
/(ODE_F
.dt_last
+ODE_F
.dt
);
201 PetscScalar HeatCapacity
= mt
->thermal
->HeatCapacity(Ti
);
202 d_grad_T_dTi
+= -aux
[i
].density
*HeatCapacity
*(2-r
)/(1-r
)/(ODE_F
.dt_last
+ODE_F
.dt
);
206 PetscScalar HeatCapacity
= mt
->thermal
->HeatCapacity(Ti
);
207 d_grad_T_dTi
+= -aux
[i
].density
*HeatCapacity
/ODE_F
.dt
;
210 MatSetValue(*jac
,zofs
[zone_index
]+2*i
+1,zofs
[zone_index
]+2*i
+1,d_grad_T_dTi
,INSERT_VALUES
);
213 //just for fill the matrix entry
215 int size
= pzone
->davcell
.size();
217 for(int j
=0;j
<electrode
.size();j
++)
219 if(electrode
[j
]==pcell
->bc_index
-1) {gate_equ
=j
;break;}
221 MatSetValue(*jac
,zofs
[zone_index
]+equ_num
*size
+gate_equ
,zofs
[zone_index
]+equ_num
*size
+gate_equ
,1.0,INSERT_VALUES
);
225 void ISZone::J2_mix_ddm_gatebc_interface(int i
,PetscScalar
*x
,Mat
*jac
, ODE_Formula
&ODE_F
, vector
<int> & zofs
,DABC
&bc
,
228 const VoronoiCell
* pcell
= &pzone
->davcell
[i
];
229 GateBC
*pbc
= dynamic_cast<GateBC
* >(bc
.Get_pointer(pcell
->bc_index
-1));
230 PetscScalar Ti
= x
[zofs
[zone_index
]+2*i
+1];
231 PetscScalar d_grad_T_dTi
=0;
232 mt
->mapping(&pzone
->danode
[i
],&aux
[i
],ODE_F
.clock
);
233 for(int j
=0;j
<pcell
->nb_num
;j
++)
235 int nb
= pcell
->nb_array
[j
];
236 PetscScalar Tj
= x
[zofs
[zone_index
]+2*nb
+1]; //potential of nb node
237 PetscScalar T_mid
= 0.5*(Tj
+Ti
);
238 PetscScalar kapa
= mt
->thermal
->HeatConduction(T_mid
);
239 PetscScalar d_grad_T_dTj
= kapa
*pcell
->elen
[j
]/pcell
->ilen
[j
];
240 d_grad_T_dTi
+= -kapa
*pcell
->elen
[j
]/pcell
->ilen
[j
];
241 MatSetValue(*jac
,zofs
[zone_index
]+2*i
+1,zofs
[zone_index
]+2*nb
+1,d_grad_T_dTj
,INSERT_VALUES
);
243 const VoronoiCell
* ncell
= pz
->pzone
->davcell
.GetPointer(n
);
244 pz
->mt
->mapping(&pz
->pzone
->danode
[n
],&pz
->aux
[n
],ODE_F
.clock
);
245 for(int j
=0;j
<ncell
->nb_num
;j
++)
247 int nb
= ncell
->nb_array
[j
];
248 PetscScalar Tj_n
= x
[zofs
[pz
->pzone
->zone_index
]+2*nb
+1];
249 PetscScalar kapa
= pz
->mt
->thermal
->HeatConduction(0.5*(Ti
+Tj_n
));
250 d_grad_T_dTi
+= -kapa
*ncell
->elen
[j
]/ncell
->ilen
[j
];
251 PetscScalar d_grad_T_dTj
= kapa
*ncell
->elen
[j
]/ncell
->ilen
[j
];
252 MatSetValue(*jac
,zofs
[zone_index
]+2*i
+1,zofs
[pz
->pzone
->zone_index
]+2*nb
+1,d_grad_T_dTj
,INSERT_VALUES
);
255 MatSetValue(*jac
,zofs
[zone_index
]+2*i
+0,zofs
[zone_index
]+2*i
+0,1.0,INSERT_VALUES
);
257 if(ODE_F
.TimeDependent
==true)
260 if(ODE_F
.BDF_Type
==BDF2
&&ODE_F
.BDF2_restart
==false)
262 PetscScalar r
= ODE_F
.dt_last
/(ODE_F
.dt_last
+ODE_F
.dt
);
263 PetscScalar HeatCapacity1
= mt
->thermal
->HeatCapacity(Ti
);
264 PetscScalar HeatCapacity2
= pz
->mt
->thermal
->HeatCapacity(Ti
);
265 d_grad_T_dTi
+= -aux
[i
].density
*HeatCapacity1
*(2-r
)/(1-r
)/(ODE_F
.dt_last
+ODE_F
.dt
)*pcell
->area
266 -pz
->aux
[n
].density
*HeatCapacity2
*(2-r
)/(1-r
)/(ODE_F
.dt_last
+ODE_F
.dt
)*ncell
->area
;
270 PetscScalar HeatCapacity1
= mt
->thermal
->HeatCapacity(Ti
);
271 PetscScalar HeatCapacity2
= pz
->mt
->thermal
->HeatCapacity(Ti
);
272 d_grad_T_dTi
+= -aux
[i
].density
*HeatCapacity1
/ODE_F
.dt
*pcell
->area
273 -pz
->aux
[n
].density
*HeatCapacity2
/ODE_F
.dt
*ncell
->area
;
276 MatSetValue(*jac
,zofs
[zone_index
]+2*i
+1,zofs
[zone_index
]+2*i
+1,d_grad_T_dTi
,INSERT_VALUES
);
278 //just for fill the matrix entry
280 int size
= pzone
->davcell
.size();
282 for(int j
=0;j
<electrode
.size();j
++)
283 if(electrode
[j
]==pcell
->bc_index
-1)
285 MatSetValue(*jac
,zofs
[zone_index
]+equ_num
*size
+gate_equ
,zofs
[zone_index
]+equ_num
*size
+gate_equ
,1.0,INSERT_VALUES
);
289 void ISZone::F2_mix_gate_electrode_Load(int bc_index
, double & current
, PetscScalar
& pdI_pdV
, Vec
& pdI_pdw
, Vec
& pdF_pdV
,PetscScalar
*x
, Mat
*jac
, ODE_Formula
&ODE_F
, vector
<int> &zofs
, DABC
& bc
, PetscScalar DeviceDepth
)
291 GateBC
*pbc
= dynamic_cast <GateBC
* > (bc
.Get_pointer(bc_index
));
293 PetscScalar e
= mt
->e
;
294 VecZeroEntries(pdI_pdw
);
295 VecZeroEntries(pdF_pdV
);
296 PetscScalar
* apdI_pdw
;
297 PetscScalar
* apdF_pdV
;
298 VecGetArray(pdI_pdw
,&apdI_pdw
);
299 VecGetArray(pdF_pdV
,&apdF_pdV
);
300 current
= pbc
->Get_Current_new();
302 for(int j
=0;j
<bc
[bc_index
].psegment
->node_array
.size();j
++)
304 int node
=bc
[bc_index
].psegment
->node_array
[j
];
305 const VoronoiCell
* pcell
= pzone
->davcell
.GetPointer(node
);
306 for(int k
=0;k
<pcell
->nb_num
;k
++)
308 int nb
= pcell
->nb_array
[k
];
309 //for displacement current
310 PetscScalar dJ_dVi
= aux
[node
].eps
/pcell
->ilen
[k
]/ODE_F
.dt
*pcell
->elen
[k
];
311 PetscScalar dJ_dVr
= -aux
[node
].eps
/pcell
->ilen
[k
]/ODE_F
.dt
*pcell
->elen
[k
];
312 apdI_pdw
[zofs
[zone_index
]+2*node
] += dJ_dVi
*DeviceDepth
;
313 apdI_pdw
[zofs
[zone_index
]+2*nb
] += dJ_dVr
*DeviceDepth
;
315 apdF_pdV
[zofs
[zone_index
]+2*node
] = 1.0;
318 VecRestoreArray(pdI_pdw
,&apdI_pdw
);
319 VecRestoreArray(pdF_pdV
,&apdF_pdV
);