more fix on Ec/Ev.
[gss-tcad.git] / src / solver / ddm1ac / insulequ1ac.cc
blob1ea68848ddf006cd25c4bd94fed002f7c7b1fac7
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: Nov 19, 2006 */
14 /* */
15 /* Gong Ding */
16 /* gdiso@ustc.edu */
17 /* NINT, No.69 P.O.Box, Xi'an City, China */
18 /* */
19 /*****************************************************************************/
20 #include "zonedata.h"
22 void ISZone::AC1_ddm_inner(int i, PetscScalar omega, PetscScalar *x, Vec *r, Mat *jac, vector<int> &zofs)
24 const VoronoiCell* pcell = pzone->davcell.GetPointer(i);
25 int N = zofs[zofs.size()-1];
26 PetscScalar d_grad_P_dVi = 0;
27 int IndexRe = zofs[zone_index]+i;
28 int IndexIm = N+zofs[zone_index]+i;
29 for(int j=0;j<pcell->nb_num;j++)
31 int nb = pcell->nb_array[j];
32 int ColRe = zofs[zone_index]+nb;
33 int ColIm = N+zofs[zone_index]+nb;
34 PetscScalar d_grad_P_dVj = pcell->elen[j]/pcell->ilen[j]/pcell->area;
35 MatSetValue(*jac,IndexRe,ColRe,d_grad_P_dVj,INSERT_VALUES);
36 MatSetValue(*jac,IndexIm,ColIm,d_grad_P_dVj,INSERT_VALUES);
37 d_grad_P_dVi += -pcell->elen[j]/pcell->ilen[j]/pcell->area;
39 MatSetValue(*jac,IndexRe,IndexRe,d_grad_P_dVi,INSERT_VALUES);
40 MatSetValue(*jac,IndexIm,IndexIm,d_grad_P_dVi,INSERT_VALUES);
41 VecSetValue(*r,IndexRe,0,INSERT_VALUES);
42 VecSetValue(*r,IndexIm,0,INSERT_VALUES);
46 void ISZone::AC1_ddm_electrode_insulator_interface(int i, PetscScalar omega, PetscScalar *x, Vec *r, Mat *jac, vector<int> &zofs,
47 DABC & bc,ElZone *pz, int n)
49 int N = zofs[zofs.size()-1];
50 MatSetValue(*jac,zofs[zone_index]+i,zofs[zone_index]+i,1.0,INSERT_VALUES);
51 MatSetValue(*jac,N+zofs[zone_index]+i,N+zofs[zone_index]+i,1.0,INSERT_VALUES);
52 MatSetValue(*jac,zofs[zone_index]+i,zofs[pz->zone_index]+n,-1.0,INSERT_VALUES);
53 MatSetValue(*jac,N+zofs[zone_index]+i,N+zofs[pz->zone_index]+n,-1.0,INSERT_VALUES);
54 VecSetValue(*r,zofs[zone_index]+i,0.0,INSERT_VALUES);
55 VecSetValue(*r,N+zofs[zone_index]+i,0.0,INSERT_VALUES);
59 void ISZone::AC1_ddm_semiconductor_insulator_interface(int i, PetscScalar omega, PetscScalar *x, Vec *r, Mat *jac, vector<int> &zofs,
60 DABC & bc,SMCZone *pz, int n)
62 int N = zofs[zofs.size()-1];
63 MatSetValue(*jac,zofs[zone_index]+i,zofs[zone_index]+i,1.0,INSERT_VALUES);
64 MatSetValue(*jac,N+zofs[zone_index]+i,N+zofs[zone_index]+i,1.0,INSERT_VALUES);
65 MatSetValue(*jac,zofs[zone_index]+i,zofs[pz->zone_index]+3*n+0,-1.0,INSERT_VALUES);
66 MatSetValue(*jac,N+zofs[zone_index]+i,N+zofs[pz->zone_index]+3*n+0,-1.0,INSERT_VALUES);
67 VecSetValue(*r,zofs[zone_index]+i,0.0,INSERT_VALUES);
68 VecSetValue(*r,N+zofs[zone_index]+i,0.0,INSERT_VALUES);
72 void ISZone::AC1_ddm_insulator_insulator_interface(int i, PetscScalar omega, PetscScalar *x, Vec *r, Mat *jac, vector<int> &zofs,
73 DABC & bc,ISZone *pz, int n)
75 int N = zofs[zofs.size()-1];
76 if(zone_index > pz->pzone->zone_index)
78 MatSetValue(*jac,zofs[zone_index]+i,zofs[zone_index]+i,1.0,INSERT_VALUES);
79 MatSetValue(*jac,N+zofs[zone_index]+i,N+zofs[zone_index]+i,1.0,INSERT_VALUES);
80 MatSetValue(*jac,zofs[zone_index]+i,zofs[pz->zone_index]+3*n+0,-1.0,INSERT_VALUES);
81 MatSetValue(*jac,N+zofs[zone_index]+i,N+zofs[pz->zone_index]+3*n+0,-1.0,INSERT_VALUES);
82 VecSetValue(*r,zofs[zone_index]+i,0.0,INSERT_VALUES);
83 VecSetValue(*r,N+zofs[zone_index]+i,0.0,INSERT_VALUES);
84 return;
86 PetscScalar d_grad_P_dVi = 0;
87 int IndexRe = zofs[zone_index]+i;
88 int IndexIm = N+zofs[zone_index]+i;
89 const VoronoiCell* pcell = pzone->davcell.GetPointer(i);
90 for(int j=0;j<pcell->nb_num;j++)
92 int nb = pcell->nb_array[j];
93 int ColRe = zofs[zone_index]+nb;
94 int ColIm = N+zofs[zone_index]+nb;
95 PetscScalar d_grad_P_dVj = aux[i].eps*pcell->elen[j]/pcell->ilen[j];
96 MatSetValue(*jac,IndexRe,ColRe,d_grad_P_dVj,INSERT_VALUES);
97 MatSetValue(*jac,IndexIm,ColIm,d_grad_P_dVj,INSERT_VALUES);
98 d_grad_P_dVi += -aux[i].eps*pcell->elen[j]/pcell->ilen[j];
101 const VoronoiCell* ncell = pz->pzone->davcell.GetPointer(n);
102 for(int j=0;j<ncell->nb_num;j++)
104 int nb = ncell->nb_array[j];
105 int ColRe = zofs[pz->zone_index]+nb;
106 int ColIm = N+zofs[pz->zone_index]+nb;
107 PetscScalar d_grad_P_dVj = pz->aux[n].eps*ncell->elen[j]/ncell->ilen[j];
108 MatSetValue(*jac,IndexRe,ColRe,d_grad_P_dVj,INSERT_VALUES);
109 MatSetValue(*jac,IndexIm,ColIm,d_grad_P_dVj,INSERT_VALUES);
110 d_grad_P_dVi += -pz->aux[n].eps*ncell->elen[j]/ncell->ilen[j];
112 MatSetValue(*jac,IndexRe,IndexRe,d_grad_P_dVi,INSERT_VALUES);
113 MatSetValue(*jac,IndexIm,IndexIm,d_grad_P_dVi,INSERT_VALUES);
115 VecSetValue(*r,IndexRe,0,INSERT_VALUES);
116 VecSetValue(*r,IndexIm,0,INSERT_VALUES);
120 void ISZone::AC1_ddm_gatebc(int i, PetscScalar omega, PetscScalar *x, Vec *r, Mat *jac, vector<int> &zofs,DABC & bc)
122 const VoronoiCell* pcell = pzone->davcell.GetPointer(i);
123 int N = zofs[zofs.size()-1];
124 int IndexRe = zofs[zone_index]+i;
125 int IndexIm = N+zofs[zone_index]+i;
126 MatSetValue(*jac,IndexRe,IndexRe,1.0,INSERT_VALUES);
127 MatSetValue(*jac,IndexIm,IndexIm,1.0,INSERT_VALUES);
129 int equ_num = 1;
130 int size = pzone->davcell.size();
131 int gate_equ;
132 for(int j=0;j<electrode.size();j++)
133 if(electrode[j]==pcell->bc_index-1)
134 {gate_equ=j;break;}
135 MatSetValue(*jac,IndexRe,zofs[zone_index]+equ_num*size+gate_equ,-1.0,INSERT_VALUES);
136 MatSetValue(*jac,IndexIm,N+zofs[zone_index]+equ_num*size+gate_equ,-1.0,INSERT_VALUES);
137 VecSetValue(*r,IndexRe,0.0,INSERT_VALUES);
138 VecSetValue(*r,IndexIm,0.0,INSERT_VALUES);
142 void ISZone::AC1_ddm_chargebc(int i, PetscScalar omega, PetscScalar *x, Vec *r, Mat *jac, vector<int> &zofs,DABC & bc)
144 const VoronoiCell* pcell = pzone->davcell.GetPointer(i);
145 int N = zofs[zofs.size()-1];
146 int IndexRe = zofs[zone_index]+i;
147 int IndexIm = N+zofs[zone_index]+i;
148 MatSetValue(*jac,IndexRe,IndexRe,1.0,INSERT_VALUES);
149 MatSetValue(*jac,IndexIm,IndexIm,1.0,INSERT_VALUES);
151 int equ_num = 1;
152 int size = pzone->davcell.size();
153 int charge_equ;
154 for(int j=0;j<electrode.size();j++)
155 if(electrode[j]==pcell->bc_index-1)
156 {charge_equ=j;break;}
157 MatSetValue(*jac,IndexRe,zofs[zone_index]+equ_num*size+charge_equ,-1.0,INSERT_VALUES);
158 MatSetValue(*jac,IndexIm,N+zofs[zone_index]+equ_num*size+charge_equ,-1.0,INSERT_VALUES);
159 VecSetValue(*r,IndexRe,0.0,INSERT_VALUES);
160 VecSetValue(*r,IndexIm,0.0,INSERT_VALUES);
163 void ISZone::AC1_gate_electrode(int i, PetscScalar omega, PetscScalar *x, Vec *r, Mat *jac, vector<int> &zofs,
164 DABC & bc, PetscScalar DeviceDepth )
166 int equ_num = 1;
167 int size = pzone->davcell.size();
168 int z = zone_index;
169 int N = zofs[zofs.size()-1];
170 int bc_index = electrode[i];
171 GateBC *pbc = dynamic_cast <GateBC * > (bc.Get_pointer(bc_index));
172 PetscScalar R = pbc->R;
173 PetscScalar C = pbc->C;
174 PetscScalar L = pbc->L;
176 complex <PetscScalar> Z1(R,omega*L);
177 complex <PetscScalar> Y2(0,omega*C);
178 complex <PetscScalar> A,I(0,1);
179 PetscInt IndexRe,IndexIm,ColRe,ColIm;
180 PetscInt EquRe = zofs[z]+equ_num*size+i;
181 PetscInt EquIm = N + zofs[z]+equ_num*size+i;
182 MatAssemblyBegin(*jac,MAT_FLUSH_ASSEMBLY);
183 MatAssemblyEnd(*jac,MAT_FLUSH_ASSEMBLY);
184 for(int j=0;j<bc[bc_index].psegment->node_array.size();j++)
186 int node=bc[bc_index].psegment->node_array[j];
187 const VoronoiCell* pcell = pzone->davcell.GetPointer(node);
188 complex <PetscScalar> dJdisp_dVi=0,dJdisp_dVr=0;
189 IndexRe = zofs[z]+node;
190 IndexIm = N + zofs[z]+node;
191 for(int k=0;k<pcell->nb_num;k++)
193 int nb = pcell->nb_array[k];
194 ColRe = zofs[z]+nb;
195 ColIm = N + zofs[z]+nb;
197 //for displacement current
198 dJdisp_dVi += omega*aux[node].eps*pcell->elen[k]/pcell->ilen[k]*I;
199 dJdisp_dVr = -omega*aux[node].eps*pcell->elen[k]/pcell->ilen[k]*I;
201 A = dJdisp_dVr*DeviceDepth*Z1;
202 MatSetValue(*jac,EquRe,ColRe,A.real(),ADD_VALUES);
203 MatSetValue(*jac,EquRe,ColIm,-A.imag(),ADD_VALUES);
204 MatSetValue(*jac,EquIm,ColRe,A.imag(),ADD_VALUES);
205 MatSetValue(*jac,EquIm,ColIm,A.real(),ADD_VALUES);
207 A = dJdisp_dVi*DeviceDepth*Z1;
208 MatSetValue(*jac,EquRe,IndexRe,A.real(),ADD_VALUES);
209 MatSetValue(*jac,EquRe,IndexIm,-A.imag(),ADD_VALUES);
210 MatSetValue(*jac,EquIm,IndexRe,A.imag(),ADD_VALUES);
211 MatSetValue(*jac,EquIm,IndexIm,A.real(),ADD_VALUES);
213 MatAssemblyBegin(*jac,MAT_FLUSH_ASSEMBLY);
214 MatAssemblyEnd(*jac,MAT_FLUSH_ASSEMBLY);
216 complex <PetscScalar> K = Z1*Y2 + PetscScalar(1.0);
217 MatSetValue(*jac,EquRe,EquRe,K.real(),INSERT_VALUES);
218 MatSetValue(*jac,EquRe,EquIm,-K.imag(),INSERT_VALUES);
219 MatSetValue(*jac,EquIm,EquRe,K.imag(),INSERT_VALUES);
220 MatSetValue(*jac,EquIm,EquIm,K.real(),INSERT_VALUES);
222 VecSetValue(*r,EquRe,pbc->Vac, INSERT_VALUES);
223 VecSetValue(*r,EquIm,0.0, INSERT_VALUES);
227 void ISZone::AC1_charge_electrode(int i, PetscScalar omega, PetscScalar *x, Vec *r, Mat *jac, vector<int> &zofs, DABC & bc)
229 int equ_num = 1;
230 int size = pzone->davcell.size();
231 int z = zone_index;
232 int N = zofs[zofs.size()-1];
233 int bc_index = electrode[i];
234 ChargedContactBC *pbc = dynamic_cast <ChargedContactBC * > (bc.Get_pointer(bc_index));
236 PetscInt IndexRe,IndexIm,ColRe,ColIm;
237 PetscInt EquRe = zofs[z]+equ_num*size+i;
238 PetscInt EquIm = N + zofs[z]+equ_num*size+i;
239 MatAssemblyBegin(*jac,MAT_FLUSH_ASSEMBLY);
240 MatAssemblyEnd(*jac,MAT_FLUSH_ASSEMBLY);
241 for(int j=0;j<bc[bc_index].psegment->node_array.size();j++)
243 int node=bc[bc_index].psegment->node_array[j];
244 const VoronoiCell* pcell = pzone->davcell.GetPointer(node);
245 PetscScalar dP_dVi=0;
246 IndexRe = zofs[z]+node;
247 IndexIm = N + zofs[z]+node;
248 for(int k=0;k<pcell->nb_num;k++)
250 int nb = pcell->nb_array[k];
251 ColRe = zofs[z]+nb;
252 ColIm = N + zofs[z]+nb;
253 PetscScalar dP_dVj=aux[node].eps/pcell->ilen[k]*pcell->elen[k];
254 MatSetValue(*jac,IndexRe,ColRe,dP_dVj,ADD_VALUES);
255 MatSetValue(*jac,IndexIm,ColIm,dP_dVj,ADD_VALUES);
256 dP_dVi+=-aux[node].eps/pcell->ilen[k]*pcell->elen[k];
258 MatSetValue(*jac,IndexRe,IndexRe,dP_dVi,ADD_VALUES);
259 MatSetValue(*jac,IndexIm,IndexIm,dP_dVi,ADD_VALUES);
261 MatAssemblyBegin(*jac,MAT_FLUSH_ASSEMBLY);
262 MatAssemblyEnd(*jac,MAT_FLUSH_ASSEMBLY);
264 VecSetValue(*r,EquRe,0.0, INSERT_VALUES);
265 VecSetValue(*r,EquIm,0.0, INSERT_VALUES);
269 void ISZone::AC1_gate_electrode_current(int i, PetscScalar omega, PetscScalar *x, Vec *v, vector<int> &zofs,
270 DABC & bc, PetscScalar DeviceDepth, PetscScalar &IacRe, PetscScalar &IacIm)
272 int equ_num = 1;
273 int size = pzone->davcell.size();
274 int z = zone_index;
275 int N = zofs[zofs.size()-1];
276 int bc_index = electrode[i];
277 GateBC *pbc = dynamic_cast <GateBC * > (bc.Get_pointer(bc_index));
278 PetscInt IndexRe,IndexIm,ColRe,ColIm;
279 PetscScalar *vv;
280 VecGetArray(*v,&vv);
281 IacRe = 0;
282 IacIm = 0;
283 //for displacement current
284 for(int j=0;j<bc[bc_index].psegment->node_array.size();j++)
286 int node=bc[bc_index].psegment->node_array[j];
287 const VoronoiCell* pcell = pzone->davcell.GetPointer(node);
288 IndexRe = zofs[z]+node;
289 IndexIm = N + zofs[z]+node;
290 PetscScalar ViRe = vv[IndexRe];
291 PetscScalar ViIm = vv[IndexIm];
293 for(int k=0;k<pcell->nb_num;k++)
295 int nb = pcell->nb_array[k];
296 ColRe = zofs[z]+nb;
297 ColIm = N + zofs[z]+nb;
299 PetscScalar VjRe = vv[ColRe];
300 PetscScalar VjIm = vv[ColIm];
302 IacRe += -omega*aux[node].eps*(ViIm-VjIm)/pcell->ilen[k]*pcell->elen[k]*DeviceDepth;
303 IacIm += omega*aux[node].eps*(ViRe-VjRe)/pcell->ilen[k]*pcell->elen[k]*DeviceDepth;
306 VecRestoreArray(*v,&vv);