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?

  • 2
    $\begingroup$ FiPy (ctcms.nist.gov/fipy) can probably handle this. $\endgroup$ Commented Nov 6, 2018 at 14:58
  • $\begingroup$ @BillGreene thank you! I will try to use this first $\endgroup$
    – wintercap
    Commented Nov 7, 2018 at 4:45
  • $\begingroup$ I know it is not your question, but since you mentioned this is for pedagogical purposes: coding a program that does just that was a big part of my master thesis and I considered it an extremely instructive experience. It still provides me with insights today (i'm a phd student) when interpreting results from the commercial software we use at my lab. $\endgroup$
    – Alexis
    Commented Dec 6, 2018 at 17:11

2 Answers 2


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++.


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


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.


Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.