Hi dear deal.ii community!
I am trying to solve 2D transient problem in semiconductor devices and I'm
stuck with application of Dirichlet boundary values via ConstrainMatrix or
AffineMatrix (in the last version of this library)
The program (which I am rewriting) solves consecutively two equation: (i)
for a given charge density it solves Poisson equation, and than (ii) the
continuity equation to find the densities in next step.
The poisson equation looks like this:
[image: eq.png]
This is a formulation of mixed finite element problem with phi as electric
potential and E = (Ex, Ey) as a vector of electric field.
Program uses a *Raviart-Thomas* finite element to describe E, and *Discontinous
Galerkin* (FE_DGQ) for potential.
Weak fe formulation:
[image: phi.png]
[image: E.png]
The geometry (as well as boundary conditions) is quite simple:
[image: bc.png]
I have choosen this kind of boundary conditions because:
- the potential is applied on the left side of device
- Interesting things (non-zero carrier densities) shound appeard only near
the interface of two domains so I assumed that Ex=0, Ey=0 on Dirichlet
boundaries
- there should be no current through top and bottom of domains, so I
decided in most simple case to have zero electric field perpendicular to
the boundary.
I hope those are reasonable bc.
My implementation of boundary conditions goes like this:
For a given FESystem fe and DoFHanfler dof_handler:
dof_handler.distribute_dofs(fe);
DoFRenumbering::component_wise(dof_handler);
// make hanging node constraints
constraints.clear();
DoFTools::make_hanging_node_constraints(dof_handler, constraints);
// add constaints to Poisson equation
const FEValuesExtractors::Vector ElectricField(0);
ComponentMask electric_field_mask = fe.component_mask(ElectricField);
DoFTools::make_zero_boundary_constraints(dof_handler,
Dirichlet, // DIRICHLET BOUNDARY INDICATOR
constraints,
electric_field_mask);
std::set<unsigned int> no_flux_boundary;
no_flux_boundary.insert(Neumann);
VectorTools::compute_no_normal_flux_constraints(dof_handler,
0,
no_flux_boundary,
constraints);
constraints.close();
during the procedure of assembling the global matrix I use in each loop:
constraints.distribute_local_to_global( data.local_matrix,
data.local_dof_indices,
Poisson_object.system_matrix);
and for right hand side:
constraints.distribute_local_to_global(data.local_rhs,
data.local_dof_indices,
Poisson_object.system_rhs);
Unfortunately when I look at the output the E value on Dirichlet boundary
only x-component Ex is 0, but not Ey. On Neumann boundary non of the
components are 0, but of course Ey should be. Additionally the value of a
potential (which comes from face integral on right hand side) on Dirichlet
boundaries is not longer as it should be (eg. I ascribed 1V and in the
output I read 1e-5V)
What did I do wrong? Is it because the Raviart-Thomas should be treated
differently? Maybe the value that I look at in Paraview is extrapolated
from the values from the middle of the cells face and that is why the
values are incorrect?
Another problem is that for the simplest case all Ey on all of the domain
should be zero. Maybe this is the problem, but I kind of don't want to
write a code especially for the simplest case when Ey=0 -> is there a
better way to do this?
I will appreciate any help!
--
The deal.II project is located at http://www.dealii.org/
For mailing list/forum options, see
https://groups.google.com/d/forum/dealii?hl=en
---
You received this message because you are subscribed to the Google Groups
"deal.II User Group" group.
To unsubscribe from this group and stop receiving emails from it, send an email
to [email protected].
To view this discussion on the web visit
https://groups.google.com/d/msgid/dealii/ff80bf4d-f269-4ad6-b131-8aab24e94f0d%40googlegroups.com.
