|
| |
Linear Poisson Equation
In this chapter it is described how to set up the Poisson matrix for the
linear Poisson equation.
The following topics are addressed:
- How Poisson clusters are defined and set up.
- How the Poisson matrix is assembled.
- How boundary conditions are implemented.
- How the right hand side of the Poisson equation is assembled.
- How the Poisson matrix is stored.
Poisson clusters
Before the Poisson matrix can be created, it must
be defined
- where the Poisson equation is solved.
- which boundary conditions on domain/inner
boundaries to take.
The relevant information is written to
MODULE poisson:
[first half]
Poisson contact regions are a set of points on
material grid. All adjoining points
on the physical grid are assigned the corresponding boundary condition.
A Poisson contact region is the contact region defined by the keyword
$poisson-boundary-conditions. A Poisson contact region is NOT
the region where the Poisson equation has to be solved, in fact in a Poisson
contact region, no equations have to be solved at all.
MODULE poisson
num_poisson_clusters ...
poisson_clusterM(1..num_poisson_clusters)
has TYPE
TYPE :: poisson_cluster
INTEGER :: num_points Number of points in that cluster
INTEGER,DIMENSION(1:num_points),POINTER :: num_guenther_points
Number of points on
material grid which belong to that cluster.
INTEGER :: kind_point = 1 -->
exempt region,
set field in valueV(1:3) on boundary as
Neumann condition.
= 2 -->
set to Dirichlet value valueV(1)
= 3 -->
set to Dirichlet value (valueV(1) + built-in potential on that point)
= 4 -->
fix conduction band edge [barrier in eV]
above Fermi levels for electrons (fermi_nV)
Note: min(E_cM(i,:)) is taken.
=> used for Schottky contact
REAL(8),DIMENSION
(3) :: valueV ...
potential in V respectively barrier in [eV] are stored in ..(1)
field [V/m]
Note: The field is a vector E =
(Ex,Ey,Ez) = Ex * x + Ey *
y +
Ez * z.
The Poisson assembler takes all 8 octants for each
point.
Wherever there is a transition from an octant on a
Poisson cluster to an octant outside a Poisson cluster and a field is specified
for that cluster, the corresponding component of the field is taken.
If any of the 8 octants has the option kind_point=2,3,4
(i.e. any kind
of Dirichlet) the assembler ignores all fields and
sets the potential on that
point to the specified Dirichlet value.
For kind_point=4 --> the band edge of
the adjoining octant is taken.
Note: This option is used for
Schottky contacts, i.e. at a contact/semiconductor
interface. So one has to take the band edge in the semiconductor, not in
the contact.
If
colliding potentials are specified the cluster with the lowest value
is valid.
The same procedure for boundary points on boundary
of simulation domain:
poissonbnd_xM(1..2,j,k,4): poissonbnd_xM(1,j,k)
--> i=1
poissonbnd_xM(2,j,k) --> i=n_1
1..4 --> interfaces 9,10,11,12
poissonbnd_yM(1..2,i,k,4): poissonbnd_xM(1,i,k)
--> j=1
poissonbnd_xM(2,i,k) --> j=n_2
1..4 --> interfaces 5,6,7,8
poissonbnd_zM(1..2,i,j,4):
poissonbnd_xM(1,i,j) --> k=1
poissonbnd_xM(2,i,j) --> k=n_1
1..4 -->
Note: Numbering of interfaces (1..12)
interface 1 <--> between octant 1 and 5
interface 2 <--> between octant 2 and 6
interface 3 <--> between octant 3 and 7
interface 4 <--> between octant 4 and 8
interface 5 <--> between octant 1 and 3
interface 6 <--> between octant 2 and 4
interface 7 <--> between octant 5 and 7
interface 8 <--> between octant 6 and 8
interface 9 <--> between octant 1 and 2
interface 10 <--> between octant 3 and 4
interface 11 <--> between octant 5 and 6
interface 12 <--> between octant 7 and 8
= 1 -->
Neumann
= 2 -->
set to Dirichlet with absolute value
= 3 -->
set to Dirichlet value (value+built-in)
= 4 -->
fix conduction band edge [barrier in eV]
above Fermi levels for electrons (fermi_nV)
bndval_xM(1..2,j,k,4)
bndval_yM(1..2,i,k,4)
bndval_zM(1..2,i,j,4)
... field, if
poissonbnd_.M=1
(directed orthonormal to the surface)
potential poissonbnd_.M=2-4 ("."
stands for x, y or z)
Note: Up to now, on simulation domain boundary only field zero is implemented.
For applying
voltages one has to define contact regions.
For having
nonvanishing fields, one has to specify a semiconductor-air interface
with surface charge.
Must be allocated and defined before call of poissonmatrix3D
and is deallocated in
deallocate_poisson.
In the input file Poisson contact regions can be
specified with a contact model ($poisson-boundary-conditions):
For the Poisson problem it is important to know
- if the built-in potential (equilibrium) is
calculated, or
- if the built-in potential is already known and a
voltage is applied.
- which contact model is specified.
In
MODULE poisson_control the variable built_inL must
be defined (.TRUE. by default):
built_inL = .TRUE. --> The
routine poisson_block calculates the built-in potential.
= .FALSE. --> The routine poisson_block calculates the
potential;
the built-in potential is assumed to exist in
BuiltInPotentialV in MODULE potentials.
The Poisson clusters and contact models are
provided in the modules:
mod_num_of_ohmics
mod_num_of_schottkys
mod_num_of_dirichlets
mod_num_of_neumanns
mod_num_of_mixed
mod_apply_dirichlet_to_cluster
mod_dirichlet_grid_points
mod_potential_for_dirichlet
mod_apply_ohmic_to_cluster
mod_ohmic_grid_points
mod_voltage_for_ohmic
mod_apply_schottky_to_cluster
mod_schottky_grid_points
mod_voltage_for_schottky
mod_barrier_for_schottky
See there for more details.
With that information SUBROUTINE
define_poisson_clusters defines the Poisson problem in MODULE
poisson:
ohmic, voltage
controlled --> |
When built-in potential is
calculated, electric field = 0 is chosen as boundary condition.
When voltage is applied, Dirichlet condition
(for potential and Fermi levels) with
pot = pot_built_in + V_applied
fermi = fermi_equilibrium - V_applied. |
ohmic, current
controlled --> |
not yet implemented |
Schottky,
voltage controlled --> |
When built-in potential is
calculated, Dirichlet condition is chosen with
E_c-pot = fermi_equilibrium + barrier --> pot.
That means that the potential is chosen such
that the conduction band edge is
fixed (barrier) above the
equilibrium (constant) Fermi level.
When voltage is applied, Dirichlet condition
(for potential and Fermi levels) with
pot = pot_built_in + V_applied
fermi = fermi_equilibrium - V_applied. |
Schottky,
current controlled --> |
not yet implemented |
Dirichlet
--> |
Potential is set brute force
to value Dirichlet (pure mathematical condition). |
Neumann
--> |
Electric field is set to Neumann. |
mixed, i.e.
relation between Dirichlet
and Neumann (Cauchy) --> |
not yet
implemented |
NOTE: The grid points inside a Poisson
cluster (contact) are exempted from simulation (i.e. set to Dirichlet which is
of no importance).
NOTE: Boundary conditions on domain
boundaries cannot be specified via input yet, so they are set to
type Neumann (electric field = 0) here.
Setup of Poisson matrix
After the Poisson clusters have been defined by
SUBROUTINE
define_poisson_clusters the
Poisson matrix is set up in
SUBROUTINE poissonmatrix3D
(poissonmatrix2D,poissonmatrix1D).
The Poisson solver allows one to exempt parts of the
simulation domain from simulation. On these internal surfaces Neumann boundary
conditions can be specified (zero by
default). The corresponding information is provided by
SUBROUTINE irregular_boundary(i,j,k,ch,which_octantsVL,fieldM):
SUBROUTINE irregular_boundary(i,j,k,irr_bcC,which_octantsVL,fieldM)
(i,j,k) point on physical grid
[Input]
which_octantsVL(1..8) .TRUE. --->
octant
within simulation region [Output]
.FALSE. ---> octant outside simulation region, i.e. within a
Poisson contact region
irr_bcC = 'in' if all octants within
simulation region.
'ot' if all octants outside simulation region.
'on' if
at least one octant within simulation region, but not all.
fieldM(1..8,1..3) : fieldM(n,1) field component
in x-direction for box n
fieldM(n,2) field component in y-direction for box n
fieldM(n,3)
field component in z-direction for box n
[electric field in
V/m]
The Poisson assembler further needs the
position dependent dielectric constant which is obtained via
epsilon_poisson:
-
FUNCTION
epsilon_poisson(ch,i,j,k) result(epsi)
provides for point number (i,j,k) on
physical grid
value of epsilon on
that point.
ch = number of octant
('1','2','3', ... )
epsilon in SI units (including
vacuum_permeability)
ch = 'n' <--> normal (i.e.
single) point
Finally, the program in first SETUP always sets
Neumann boundary conditions which are set to zero if none are defined. Later on, Dirichlet boundary
conditions are set (i.e. Neumann are overwritten by Dirichlet).
The fields on the domain boundaries can be obtained
via
boundary_value:
FUNCTION
boundary_value(ch,i,j,k,whichinterface) RESULT(bv)
This may just be called for (i,j,k) on the domain boundary, otherwise stop.
Input:
(i,j,k) on physical grid
ch = 'x' or 'y' or 'z' (direction of differential)
called at any setup of Poisson
whichinterface=1..4 --> interfaces
9,10,11,12 for 'x'
5, 6, 7, 8
for 'y'
1, 2, 3, 4
for 'z'
For whichinterface the same convention as for poissonbnd_xM... is
.TRUE.
Output:
If Neumann boundary conditions are valid, the electric field is provided.
Otherwise: set to zero
Note: In the first setup always Neumann conditions are implemented.
They are overwritten with Dirichlet values (if specified) later on.
The matrix which is created at the setup is
stored in MODULE
poisson:
diagV ..
DIMENSION(1:lengvecV)
contains diagonal elements of Poisson matrix
matrixindicesM ..
DIMENSION(1:nummaelemente,1:2)
contains for nonvanishing offdiagonal matrix elements row and column
index
Because matrix is symmetric, only upper half is stored.
matrixwerteV .. DIMENSION(1:nummaelemente)
contains corresponding values
matrixposV(i) ..
DIMENSION(1..lengvecV)
matrix_pos(n) points to the position, where the the nonvanishing elements of row n
(upper triangle) are placed in the arrays.
matrixindices and matrixwerteV
If the line doesn't have any offdiagonal elements it points to the next entry
unequal to zero.
nummaelemente .. number of
offdiagonal matrix elements
bV ..
DIMENSION(1:lengvecV)
contains the right hand side of the matrix equation (but only the charge
independent parts = boundary conditions)
rhoeinfachV .. DIMENSION(1..lengvecV)
contains for each single point the factor by which
the volume density rho must be multiplied to get the charge-dependent
right hand side b of the matrix equation. At multiple points rhoeinfachV is
zero.
skarhoM ..
DIMENSION(1..size(punktartCV),1:8)
contains the vector by which the vector (rho1,...rho8) must be (scalar)
multiplied to get the right hand side for multiple points, i.e. the
volume charge density (rho1, ... ,rho8) in these eight octants
skarhoM contains the corresponding partial volumes. All
components 1..8 are always defined, even if the octant doesn't lie in the
simulation area (in this case the value is zero.)
Note: Numbering of boxes(1..8)
octant 1 <--> (i-,j-,k-)
octant 2 <--> (i+,j-,k-)
octant 3 <--> (i-,j+,k-)
octant 4 <--> (i-,j+,k-)
octant 5 <--> (i-,j-,k+)
octant 6 <--> (i+,j-,k+)
octant 7 <--> (i-,j+,k+)
octant 8 <--> (i-,j+,k+)
where 'i-' denotes the
octant between i-1,i
and 'i+'
denotes the octant between i,i+1
skasigma3D_M .. DIMENSION(1..size(punktartCV),1:12)
skarhoM but
corresponds to surface charge densities at the 12 inner interfaces. All
components 1..12 are always defined, even if the interior interfaces are not
lying in the simulation area (in this case the value is zero). The numbering
1..12 is based on the convention of numbering of interfaces.
In 2D:
Note: Numbering of interfaces (1..4)
interface 1 <--> between
octant 1 and 2
interface 2 <--> between
octant 1 and 3
3
interface 3 <--> between
octant 3 and 4
2 + 4
interface 4 <--> between
octant 2 and 4 1
At domain boundaries there are only 4 interior interfaces in each case. They are
stored in:
ska_bnd_x3DM(1..2,j,k,4) : ska_bnd_x3DM(1,j,k) --> i=1
ska_bnd_x3DM(2,j,k) --> i=n_1
1..4 --> 9,10,11,12
ska_bnd_y3DM(1..2,i,k,4) : ska_bnd_y3DM(1,i,k) --> j=1
ska_bnd_y3DM(2,i,k) --> j=n_2
1..4 --> interfaces 5, 6, 7, 8
ska_bnd_z3DM(1..2,i,j,4) : ska_bnd_z3DM(1,i,j) --> k=1
ska_bnd_z3DM(2,i,j) --> k=n_3
1..4 --> interfaces 1, 2, 3, 4
The setup is done in SUBROUTINE
poissonmatrix3D (poissonmatrix2D,poissonmatrix1D):
The program loops through all points:
1) CALL WhichKindOfBoundaryPoint3D(ch,i,j,k)
[internal subroutine, in 2D/1D this distinction is done explicitly]
-->
ch ='in', if (i,j,k) is within domain.
otherwise --> ch provides information about the position of the boundary
point (see source code).
The
variables beteiligtL_V and randL_V provide
information for box integration:
beteiligtL_V(i)=.TRUE. <--> octant i is within simulation region, so box
integration is performed over it.
randL_V(1..3): says if flux through surface in x-, y-, z-direction has to
be considered in box integration.
2) CALL irregular_boundary provides for (i,j,k) in irr_bcC the following information:
irr_bcC
= 'ot' <--> (i,j,k) lies in a
region which is exempted from simulation (i.e. a Poisson
contact region)
--> set all beteiligtL_V=0
= 'in' <--> (i,j,k) lies
within simulation region (i.e. not in a Poisson contact region)
--> do not alter beteiligtL_V
= 'on' <--> (i,j,k) lies on
internal interface separating
simulation / nonsimulation region (Poisson contact region)
--> set the corresponding octants of beteiligtL_V to 0
3) Now the setup of matrix elements for that
point:
a) if
outer point (irr_bcC=='ot') --> CALL outer_point,
writes for (i,j,k) a Dirichlet condition (=0)
b) if
inner point which is no multiple point and does not adjoin a multiple point and
not to a Poisson contact region:
the differential is defined explicitly: CALL standardfall3D
c) else: perform box integration over all boxes specified by
beteiligtL_V(:)=.TRUE.:
CALL berechne_koeff3D
4) Now allocate and define
ska_bnd_x3DM, ska_bnd_y3DM, ska_bnd_z3DM.
Technical details of subroutine berechne_koeff3D
The Laplace operator is discretized in a
computational molecule with following indices:
1 2
3 4
5 6
7 |
8
---------------------------------------------------------------------------------------------------------------------
(i,j,k-1)
(i,j-1,k) (i-1,j,k) (i,j,k)
(i+1,j,k) (i,j+1,k) (i,j,k+1) |
right hand side
This convention is valid for koeff_M(n,:), n=1,...,8.
The differential star is set up for each octant n
(set to 0 if beteiligtL_V(n)=.FALSE.).
Finally koeffsum_V is created by
simply adding all octants (integral over whole box = sum of all eight integrals
over octants).
The differentials are created via Box3D(box_number=1-8,...):
NOTE: The right hand side contains charges
which in turn depend on the potential.
--> The right-hand side contains a fixed part (due to boundary conditions, ...)
and a
variable part (i.e. charge density multiplied with volume).
In bV only the fixed part is stored.
The functions Box3D(1-8,...) carry out the volume integration over the octants 1,...,8.
For nonvanishing electric fields the flux due to the field on boundaries has to be taken into account. Up to now this is only done for points on the boundary of the simulation region but not for points on the boundaries to Poisson contact regions (the corresponding values are contained in fieldM).
-
FUNCTION Box3D(box_number,eps1,eps2,eps3,eps4,randV_L)
RESULT(erg)
INTEGER
:: box_number
REAL :: eps1,eps2,eps3,eps4
! eps1 ->
permeability in (i,j,k)
REAL ,DIMENSION(1:9) :: erg
! eps2 ->
in i-direction
LOGICAL,DIMENSION(1:3) :: randV_L
! eps3 ->
j-direction
! eps4 ->
k-direction
erg = 0.0d0
!
- Laplacian = -div(grad)
erg(1) =
-(eps4+eps1)*h_1V(i-1)*h_2V(j-1)/h_3V(k-1)/8d0
erg(2) = -(eps3+eps1)*h_1V(i-1)*h_3V(k-1)/h_2V(j-1)/8d0
erg(3) = -(eps2+eps1)*h_2V(j-1)*h_3V(k-1)/h_1V(i-1)/8d0
erg(4) = -(erg(1)+erg(2)+erg(3))
erg(9) = h_1V(i-1)*h_2V(j-1)*h_3V(k-1)/8d0
! erg(9)
contains volume of Box3D
! -->
enters ska_rhoM/rhoeinfachV: variable part of right hand side
! -->
IF (randV_L(1))
erg(8) = erg(8)+(eps1*h_1V(i-1)*h_2V(j-1)*
&
boundary_value('z',i,j,k,1)/4d0)
IF (randV_L(2))
erg(8) = erg(8)+(eps1*h_1V(i-1)*h_3V(k-1)* &
boundary_value('y',i,j,k,1)/4d0)
IF (randV_L(3))
erg(8) = erg(8)+(eps1*h_3V(k-1)*h_2V(j-1)* &
boundary_value('x',i,j,k,1)/4d0)
!
----------- Note: fieldM
not yet implemented -----------
Now the area of internal interfaces are defined
(only for multiple points):
- Consider interfaces
Note: skasigma3D_M only different from zero when the interface is within device domain.
IF (PointInfoM(i,j,k)/=0) THEN
skasigma3D_M(PointInfoM(i,j,k),:)=0.0d0
IF ((i>1) .AND.(j>1 )) skasigma3D_M(PointInfoM(i,j,k),1) =
&
h_1V(i-1)*h_2V(j-1)/4d0
IF ((i<n_1).AND.(j>1 )) skasigma3D_M(PointInfoM(i,j,k),2) =
&
h_1V(i )*h_2V(j-1)/4d0
IF ((j<n_2).AND.(i>1 )) skasigma3D_M(PointInfoM(i,j,k),3) =
&
h_1V(i-1)*h_2V(j)/4d0
IF ((i<n_1).AND.(j<n_2))
skasigma3D_M(PointInfoM(i,j,k),4) = &
h_1V(i )*h_2V(j)/4d0
IF ((i>1) .AND.(k>1 )) skasigma3D_M(PointInfoM(i,j,k),5) = &
h_1V(i-1)*h_3V(k-1)/4d0
IF ((i<n_1).AND.(k>1 )) skasigma3D_M(PointInfoM(i,j,k),6) = &
h_1V(i )*h_3V(k-1)/4d0
IF ((i>1) .AND.(k<n_3))
skasigma3D_M(PointInfoM(i,j,k),7) = &
h_1V(i-1)*h_3V(k)/4d0
IF ((i<n_1).AND.(k<n-3))
skasigma3D_M(PointInfoM(i,j,k),8) = &
h_1V(i )*h_3V(k)/4d0
IF ((j>1) .AND.(k>1 )) skasigma3D_M(PointInfoM(i,j,k),9) = &
h_2V(j-1)*h_3V(k-1)/4d0
IF ((j<n_2).AND.(k>1 ))
skasigma3D_M(PointInfoM(i,j,k),10)= &
h_2V(j )*h_3V(k-1)/4d0
IF ((j>1) .AND.(k<n_3))
skasigma3D_M(PointInfoM(i,j,k),11)= &
h_2V(j-1)*h_3V(k)/4d0
IF ((j<n_2).AND.(k<n_3))
skasigma3D_M(PointInfoM(i,j,k),12)= &
h_2V(j )*h_3V(k)/4d0
END IF
The matrix elements are written (just for
upper part of matrix (symmetric) in erzeuge_matrix3D:
It creates the elements for single points and
the first line of multiple points.
Note: For multiple points the balance
from integration is fully contained in the first line.
In Poisson equation the potential is assumed to be continuous.
The vector of Poisson equation is defined on Stefan's grid (1..lengvecV)
which contains
multiple points; so additional lines are added to account for that.
Since the matrix must be symmetric and positive definite, we use following
scheme:
Generate multiple points
Double points:
The matrix is supplemented by an identity relation for the second point. The
following structure is obtained:
2nd line: -a
+ a = 0
Matrix must be symmetric, thus add
+ a
and -a
to 1st line.
-1 0 ... 0 -1
4+a
-a -1 0
0 ... 0 -1 ...
= Q
... 0
0 -a
a 0 0 0 ... ...
= 0
-1 ...
0 0 -1 0
4 -1 0 ... 0 -1
... = Q
Eightfold point:
-1 ... -1 0 ... 0 -1
6+a -a
0 0 0 0 0 0 -1
... -1 ... -1 ... = Q
-a 7a
-a -a -a -a -a -a
= 0
0 -a a
0 0 0 0 0
= 0
0 -a 0
a 0 0 0
0
= 0
0 -a 0 0
a 0 0 0
= 0
0 -a 0 0
0 a 0 0
= 0
0 -a 0 0
0 0 a 0
= 0
0 -a 0 0
0 0 0 a
= 0
-1 ... -1 0 ... 0 -1
0 0 0 0 0 0
0 6 -1 .. -1 .. -1 .. = Q
(An alternative would be, but it didn't converge with conjugate gradient,
maybe later check again...)
-1 ... -1 0 ... 0 -1 6+7a -a -a -a -a -a -a -a -1 ... -1 ... -1
-a a
-a 0 a 0 0 0 0 0
-a 0 0 a 0 0 0 0
-a 0 0 0 a 0 0 0
-a 0 0 0 0 a 0 0
-a 0 0 0 0 0 a 0
-a 0 0 0 0 0 0 a
-1 ... -1 0 ... 0 -1 0 0 0 0 0 0 0 6 -1 .. -1 .. -1
Boundary Conditions (Dirichlet points)
After the matrix has been completely set up,
Neumann boundary conditions are overwritten with Dirichlet values:
There is an auxiliary subroutine DIRICHLET3D(i,j,k,dirval)(2D/1D)
which
defines Dirichlet conditions for a certain point.
and the
right hand side bV is defined
accordingly.
Note: Matrix must remain symmetric.
So if you overwrite one line to a Dirichlet boundary condition you also have
to change the line which shares the matrix element in order to restore symmetry.
- Subroutine
set_dirichlet_poisson
sets the boundary conditions as specified in MODULE
poisson.
Subroutine
set_dirichlet_poisson
sets Dirichlet points as specified in MODULE
poisson.
Poisson matrix must be defined.
Loops through all Poisson clusters and loops through all points of the Poisson clusters and sets the according
Dirichlet value.
The points with higher cluster-number (priority) are valid.
Right hand side of Poisson equation
If one wants to calculate the right hand side
one has to allow for the fixed part of the charge density bV.
If one wants to calculate the differential of the
right hand side with respect to the potential, one has
to neglect the fixed part of the charge denisty bV. Of
course, the differential of the charge with respect to the potential is
needed:
This distinction is done via varL.
This is done in subroutine right_hand_side_poisson3D(2D,1D).
SUBROUTINE
right_hand_side_poisson3D(rhoV,sigmaV,bnd_sigma_xM,
&
bnd_sigma_yM, &
bnd_sigma_zM,varL,bendV)
Calculates right hand side of Poisson equation on basis of data given in modules.
varL = .TRUE. : Only part that depends on charge density is considered.
.FALSE.: Whole part is taken into account.
rhoV ... DIMENSION(1:lengvecV) contains charge density at each point
(in C/m3)
Note: To calculate the Jacobian the derivative with respect to the potential is needed
--> Then varL=.TRUE. has to be used.
sigmaV ... DIMENSION(1:SIZE(punktartCV),1:12) contains surface charge density at interface points (in SI
unit: C/m2)
Contains for each multiple point the interface charge of the internal interfaces.
Note: Numbering of interfaces (1..12)
interface 1 <-->
between octant 1 and 5
interface 2 <-->
between octant 2 and 6
interface 3 <-->
between octant 3 and 7
interface 4 <-->
between octant 4 and 8
interface 5 <-->
between octant 1 and 3
interface 6 <-->
between octant 2 and 4
interface 7 <-->
between octant 5 and 7
interface 8 <-->
between octant 6 and 8
interface 9 <-->
between octant 1 and 2
interface 10 <--> between octant 3 and 4
interface 11 <--> between octant 5 and 6
interface 12 <-->
between octant 7 and 8
bnd_sigma_xM,bnd_sigma_yM,bnd_sigma_zM
--> bnd_sigma_xM(1:2,j,k,4)
1 -->
for i=1
surface charge
2 -->
i=n_1
surface charge ...
Contains surface charge density at boundaries (in SI units: C/m2).
bendV ... dimension(1:lengvecV) : Output
!---------------------------------
For varL=.TRUE. only the parts in the Poisson equation are taken into
account whose derivative with respect to phiV is nonzero.
So piezo- and pyroelectric charges do not enter.
!---------------------------------
Storage of Poisson Matrix
The Poisson matrix is stored in MODULE
poisson.
diagV ...
DIMENSION(1:lengvecV)
diagonal elements of Poisson matrix
matrixindicesM ... DIMENSION(1:nummaelemente,1:2)
contains for nonvanishing offdiagonal matrix-elements row and column
index
Because matrix is symmetric, only upper half is stored.
matrixwerteV ... DIMENSION(1:nummaelemente)
contains corresponding values
matrixposV(i) ... DIMENSION(1:lengvecV)
matrix_pos(n)
points to the position where
the the nonvanishing elements of row n
(upper triangle) are placed in the arrays matrixindicesM
and matrixwerteV.
nummaelemente
number of
off-diagonal matrix elements
NOTE: checken: - wenn Teilgebiet
ausgenommen, muss Grenzflaeche auf Interface liegen?
es muesste auch gehen wenn nicht der Fall
- fieldM in poissonmatrix nicht verwendet,d.h. el. Felder als Neumann-RB werden
automatisch als Null behandelt
- Vorzeichen bei Neumann RB (dPhi/dx) checken/ Konsistenz mit Vz. bei
Boxintegration (Randwert/Feld als RB]
|
