Weak Limits of the Quantum Hydrodynamic Model

A numerical study of the dispersive limit of the quantum hydrodynamic equations for semiconductors is presented. The solution may develop high frequency oscillations when the scaled Planck constant is small. Numerical evidence is given of the fact that in such cases the solution does not converge to the solution of the formal limit equations.


INTRODUCTION
Semiconductor models based on classical or semi-classical mechanics (like the drift-diffusion equation, the hydrodynamic models and the semiclassical solid state physics Boltzmann equation) cannot be used to reasonably describe the performance of ultra-integrated devices, which are based on quantum effects. Typical examples for such devices are resonant tunneling diodes [5]. Recently the so-called quantum hydrodynamic model (QHD) has been introduced (see e.g. [2,5,10]). The QHD model has the advantage of dealing with macroscopic fluid-type unknowns (and consequently fluid-type boundary conditions can be incorporated) and it is able to describe quantum phenomena, such as negative differential resistance in a resonant tunneling diode.
Mathematically, the QHD system is a dispersive regularization of the hydrodynamic equations (HD). The regularization depends on the scaled Planck constant c and, formally, vanishes in the classical limit c 0. Thus, in the formal limit, the QHD model tends to the HD model. However, due to the non-linearity and the dispersive nature of the regularization, one cannot expect that the formal limit describes the correct limiting behavior in general. Due to the dispersive term, the solutions of the QHD system may develop fast oscillations which are not damped as e goes to zero and in that case the limiting system is not expected to be the HD system.
Here we present numerical evidence of the fact that the solution of the QHD system develops dispersive oscillations when the corresponding HD system exhibits a shock wave and that the weak *Tel." + 39 (0) 382 505690, e-mail: pietra@dragon.ian.pv.cnr.it *Corresponding author. Tel.: +43 (0) 732-2468 9187, e-mail: carsten.pohl@jk.uni-linz.ac.at limit of the QHD solution is not a solution of the HD equations. The understanding of the dispersive nature of the problem and of its dispersive limit is an important issue in semiconductor applications, since the scaled Planck constant e is often small (e 10-2, 10-3). In particular, it gives important hints on the choice of reliable numerical schemes.
The paper is organized as follows. Section 2 presents the derivation of the isothermal/isentropic QHD equations from the Schr6dinger equation.
Section 3 collects known results about a dispersive model problem (Korteweg-de Vries equation) and compares two low order space discretization schemes. Section 4 is devoted to rigorous (partial) results for the QHD model in the literature. Finally, the numerical study of the dispersive limit for the isothermal QHD system is presented in Section 5.
where n-I/)12 is the particle density and , is the scaled Debye length. In semiconductor applications C(x) denotes the given doping profile. By using the WKB-Ansatz Here we present the simple derivation based on a quantum mechanical pure state described by a non-linear Schr6dinger equation [13]. We consider the (scaled) Schr6dinger equation (2.7) Taking the gradient of (2.6), multiplying by n and using (2.5), we obtain the QHD flow equations for the density n and current density J coupled to the The isothermal case corresponds to p(n) n (and h(n) log (n)), and the isentropic case corresponds to p(n) n"r, < 7 < 3 (and h(n) n7-1 ).
In order to model the influence of collisions between electrons and the lattice atoms the relaxation term J with a small relaxation time -, is added to the right hand side of Eq. (2.8b). It is obtained with purely phenomenological considerations and cannot be deduced from the Schr6dinger equation.
In the Wigner-Boltzmann derivation of QHD models, the relaxation term comes from a Fokker-Planck relaxation approximation of the collision term. This is discussed in [11].

KORTEWEG-DE VRIES EQUATION
The Korteweg-de Vries (KdV) equation is the simplest model problem for dispersive phenomena.
u. Clearly, since the weak limit of (UC) 2 in general is not the square of the weak limit of uc, u is not a weak solution of (B).
This result is illustrated in Figure 1 A generalized KdV equation is considered in [3,17] +u + u -0, The Korteweg-de Vries equation is a dispersive approximation of the non-linear hyperbolic Burgers' (B) equation with the same initial fun- (3.2b) The small dispersion limit of the KdV equation has been studied in [12,19] using the complete integrability of KdV. As long as (B) has a smooth solution, the KdV solution u(x, t) tends strongly to that smooth solution. However, when exceeds the time where the smooth solution of (B) breaks down, u(x, t) develops dispersive oscillations. More Only when dispersion dominates, the weak limit is not a solution of the (B) equation.
It is clear that the choice of a discretization scheme for KdV must take into account the nature of the dispersive character of the equation. Looking at KdV as a perturbation of the Burgers' equation can be misleading. A numerical scheme of hyperbolic type for the first order term can easily give wrong results, unless the space mesh size is very small. A scheme for the Burgers' equation which guarantees convergence to the entropy solution contains a diffusion term with a diffusion coefficient going to zero with the mesh size. For instance, the classical upwind scheme introduces an artificial viscosity term of the order of the mesh size. Due to the presence of numerical viscosity, the computed solution may fail to exhibit fundamental properties of the true solution, unless the mesh size is much smaller than e.
When a standard central difference scheme is used for the first order term, no numerical dissipation is introduced. The leading term of the Taylor expansion of the truncation error is l g ( X i + l ) -U ( X t -1 ) Z h --Ux(xi) hzuxxx and smoothing of the oscillations is not expected.   For high order approximation schemes we refer to [4], where a pseudospectral discretization is proposed and analyzed also in the case of a small dispersion parameter .

WEAK LIMITS
The small dispersion limit of the QHD system is an open problem. Only partial answers are known in special cases [6,9]. These results are given in the context of the Schr6dinger equation and they apply to QHD-type problems without relaxation term. Moreover, WKB initial conditions with ni, Ji independent of are considered for the Schr6dinger equation in those references. Hence the initial flow of the QHD system is irrotational. In the following, we denote by n o and j0 the weak limits of n" and J', solution of (2.8).
A complete characterization of n o and j0 is given in [6]  Burgers' equation (4.3b), n o and j0 develop singularities and they are not anymore solutions of (4.3). In [9] a nonlinear Schr6dinger equation (without self-consistent potential) is considered. Convergence to the solution of the hyperbolic limit problem is proved as long as the solution of the isentropic HD (Euler) system is smooth. After the breakdown, the weak limits are not identified.

NUMERICAL STUDY
Here we present a numerical study of the weak  It is well known that system (5.2) can develop a shock discontinuity in the transition from a supersonic region to a subsonic region, while the transition from the subsonic to the supersonic regime is smooth (see e.g. [1]). The sound speed is c 1, thus, J < n characterizes the subsonic zone, and J > n characterizes the supersonic zone (we recall that the velocity is u J/n). With asymptotic analysis arguments, it is possible to show that the QHD system (5.1) may develop dispersive oscillations with c wave length in the supersonic zone and it does note have dispersive oscillations in the subsonic region. Both transitions (subsonic to supersonic, and supersonic to subsonic) are not fully understood yet. We refer to [18] and, for a numerical study, to [8].
In order to avoid numerical damping of the oscillations, system (5.1) is discretized with a central finite difference scheme. A damped Newton algorithm is then used to solve the non-linear discrete problem. We refer to [16] for details.
In the numerical tests presented here, a n + nn + diode is considered, with the smooth doping profile C plotted in Figure 4 (max(C)= 1, rain(C)= 0.1). The Debye length is chosen A-0.1, the relaxation time --1/8, and the applied voltage V0 6.5.
In Figure 5 we present the solution of the hydrodynamic equations (5.2). The density n has a shock discontinuity (approximately at x 0.7). The horizontal line is the (constant) current J. The shock occurs in the transition from the supersonic region (J > n) to the subsonic region (J < n). Figures 6-9 show the solutions n and JC of (5.1) for various values of c ( 0.01, 0.005, 0.0026, 0.001, respectively). Figure 10 is a zoom of Figure 9 in the neighbourhood of the       the solution of the hydrodynamic equation is plotted. The high frequency oscillations show a typical dispersive behaviour: the wavelength is of order c and the amplitude is of order 1. Moreover, as in the KdV case presented in Section 3, the oscillations of n , solution of (5.1), go beyond the shock discontinuity of n, solution of (5.2), and it is evident that the weak limit of n as c --, 0 (denoted by n) is not equal to n. Figure 11 illustrates the situation. The HD solution is plotted, together with two curves, which are obtained by taking mean values (in two different ways) over the oscillations. They can be considered as good approximations of n . A comparison between Figures 6-9 and Figure 11 suggests a conjecture on the behaviour of the QHD solution and its weak limit in the transonic transitions. The subsonicsupersonic transition, smooth in the HD case, occurs in a regular manner also in the QHD case.
The oscillations do not perturb that subsonic region and the weak limit n o seems to coincide with the HD solution n before the transition. The situation is more intricate for the supersonic-subsonic transition. For small e's, (see Fig. 8 and Fig. 9), some large amplitude oscillations become partly subsonic, however the oscillations are confined in the supersonic region of the weak limit n o (see Fig. 11). The supersonic-subsonic transition point of n and of n o are different. Moreover, the effect of the dispersive oscillations is transmitted to the subsonic zone. The weak limit n o differs from the HD solution in the supersonic region and in the subsonic region after that transition.