How to simulate basic semiconductor models using the Drift-diffusion model on Python?

Question

I'm trying to simulate basic semiconductor models for pedagogical purposes--starting from the Drift-diffusion model. Although I don't want to use an off-the-shelf semiconductor simulator--I'll be learning other (common, recent or obscure) models, I do want to use an off-the-shelf PDE solver.

But even for the simple 1D case, the drift-diffusion model consists of a number of coupled nonlinear PDEs.

Current density equations:

$$J_n = q n(x) \mu_n E(x) + q D_n \nabla n$$ $$J_p = q p(x) \mu_p E(x) - q D_p \nabla p$$

Continuity equations:

$$\frac{\partial{n}}{\partial{t}} = \frac{1}{q} \nabla \cdot J_n + U_n $$ $$\frac{\partial{p}}{\partial{t}} = \frac{1}{q} \nabla \cdot J_p + U_p $$

Poisson equation:

$$\nabla \cdot (\epsilon \nabla V) = -(p - n + N_D^+ - N_A^-) $$

and a number of boundary conditions.

Is there a package (preferably open source) that would take these equations, in that form, and solve them? Or perhaps is the variational form required by the tools is not as hard? In any case, what are my options?

Answer 1

I can not provide you with a solver that will solve the equation for you.

I will, however, provide a warning:

If you want to employ finite element methods, especially so called continuous Galerkin methods, you will need a weak form of the equations. For the equations you have provided that form is not trivial to obtain, but it isn't especially difficult to obtain either. If you opt for finite element methods an do derive the weak form of the equations you should be able to use FEniCS( https://fenicsproject.org/ ) for everything else.

Solving the resulting coupled systems of linear equations is not difficult (use Newton or Newton-Kantorovich methods), but will require some thought if you want to have acceptable execution times. So in summary: A bespoke implementation of a solver will be quite time consuming (you MUST test it) so try to avoid it unless you really wan to understand what can go wrong with your solver.

My background: I have advised a friend (in numerical mathematics) during his master's thesis on how to interpret the numerical results he was obtaining from solving an optimal control problem for the Nernst-Planck-Poisson equation (which is essentially the equation you have provided). He was solving it in 2D though (or 2+1 D if you count time, too) and we were using my C++.

Answer 2

I develop an open source that uses a Python front end to specify the PDEs in a scripted form.

https://devsim.org

A key point in all of this is that the Scharfetter-Gummel approach has historically been the most successful means of solving the drift diffusion equations. This approach is most compatible with using Finite Volume methods. There are hybrid approaches which use a combined form of Finite Element and Finite Volume.

Please see accepted answer here:

PDE solvers for Drift-diffusion and related models

Diode simulation results are presented here:

https://devsim.com/1d-diode-junction/ https://devsim.com/1d-diode-junction-part-ii/

Whether you want to write your solution, or try my project, please feel free to reach out on the project mailing list.