A STEADY-STATE MATHEMATICAL MODEL FOR AN EOS CAPACITOR: THE EFFECT OF THE SIZE EXCLUSION

. In this paper we present a suitable mathematical model to describe the behaviour of a hybrid electrolyte-oxide-semiconductor (EOS) device, that could be considered for application to neuro-prothesis and bio-devices. In particular, we discuss the existence and uniqueness of solutions also including the eﬀects of the size exclusion in narrow structures such as ionic channels or nanopores. The result is proved using a ﬁxed point argument on the whole domain.Ourresults provide information about the charge distribution and the potential behaviour on the device domain, and can represent a suitable framework for the development of stable numerical tools for innovative nanodevice mod-elling.


1.
Introduction. In recent decades we have witnessed an increase in our knowledge in the fields of biology and applied physics. We are now able to produce devices with nano-metrical dimensions and to interface them with complex biological structures, such as cell membranes. The applications of these new technologies range from computer science, to medicine and robotics. Research on the interface between brain and electrical devices started in 1924, when Berger first succeeded in recording human brain activity. Nowadays we are able to interface both peripheral and central nervous systems with biomedical devices, and this is now widely considered as a valid therapeutic opportunity in many pathological situations, such as spinal cord injury, epilepsy, Parkinson's disease along with others (see [20]).
Many of the hybrid devices taken into consideration in the previous discussion require the study of the interface between a semiconductor, either organic or inorganic, and an electrolyte solution. We call electrolyte each compound which reacts forming ions in water or in any other suitable solvent, and electrolyte solution the corresponding solution. The electrochemical properties of the electrolyte solution depend on the dissolved material, which is usually salt, as in our case. The simplest case is the so called 1 : 1 solution, where a salt (NaCl for example) is decomposed into two ions with unitary charge, one positive and one negative.
The characteristics and the performance of a two-way communication between semiconductor chips and biological structures have been studied extensively in recent years. See, for example, [2], [32], [36] and [37].

FEDERICA DI MICHELE, BRUNO RUBINO AND ROSELLA SAMPALMIERI
The properties of a hybrid chip composed of a neuron from a rat brain and transistors in silicon have been discussed by Fromherz in [15]. In particular he describes the interface between individual cultured neurons and silicon microstructures, from an electronic point of view. In [16] and in [17] the author extends his study to include the electrical coupling of ion channels and silicon chips. The feasibility of totally implantable devices in the brain is discussed in [12], [28], and [31] (see also [20] and [29] for a complete review). In order to reduce the impact of these devices on the patients' quality of life, it is crucial to reduce their dimensions until the range of few micrometers, maintaining optimal performance.
We remark that the electrolyte-semiconductor interface technology is also applicable to ultra-integrate chemical sensors. In [6] the authors discuss, from a theoretical point of view, an electrolyte gate AlGaN/GaN field effect transistors, whereas in [9] a more complex case of a silicon-on-insulator sensor for charged proteins detection is examined. Other similar applications are investigated in [7], [30], [34] and [39]. A mathematical and computational description of bio-sensors is discussed, for example, in [4], [5], [22], [23], [24], and [33]. Finally in [21] an hybrid semiconductor-electrolyte system is used to model solar cell. Recently many models have been proposed to describe the behaviour of a system of charged particles. For example in [25] the authors study the behaviour of a mixture of gases containing an arbitrary number of particles. Whereas in [14] the authors present a new model for an electrolyte solution accounting for a thermodynamically consistent coupling between mechanic and diffusion process.
On the other hand, to our knowledge, a mathematical model describing the interaction of nano-sized semiconductors and electrolyte solutions is not currently available, although it could have many interesting applications.
As a first approach to this set of complex problems, we consider a very simple device, the so-called EOS capacitor (Fig. 1). It consists of an electrolyte solution (Ω E ) and a doped semiconductor layer (Ω S ), separated by a thin layer of oxide (Ω O ). The interface electrolyte-oxide and oxide-semiconductor are labelled with ∂Ω EO and ∂Ω OS , respectively, whereas, for each of the three sub-domains (electrolyte, oxide and semiconductor), the boundary of the domain is labelled with ∂Ω E , ∂Ω O and ∂Ω S , respectively. See [35] or [39] for a general review, and [38] for the details about the insulator-electrolyte interface. If we compare the EOS capacitor to the Figure 1. Electrolyte Oxide Semiconductor device EOS: general case more well known MOS (metal-oxide-semiconductor), the main difference is that, in our case, it is the electrolyte solution which acts as a conductor. The external voltage is applied to the device, by means of a standard metallic contact, on the right-hand side (semiconductor), and by an electrode immersed in the solution, on the left-hand side (electrolyte solution).
The transport of ions in a neutrally charged solute is described by the so called Poisson-Nernst-Planck equations (PNP). This model was used to describe the behaviour of ions under an external applied potential (see for example [10], [11], and [19]).
Although the ions present in the electrolyte solution have, in general, different chemical and physical properties, here we are only considering an electrolyte solution 1:1. Similar results can also be obtained assuming different classes of ions. However, the mathematical computations involved in those cases become very complex.Under the previous assumptions, the Poisson-Nernst-Planck system, in steady state conditions, reads To be precise, the last equation is usually called Poisson's Equation and it describes the behaviour of the electrical potential V .
Here ρ ± are the densities of positive and negative ions in the solution. Moreover, µ ± , R, 0 , E , C E are the mobilities, the generation and recombination rates, the electrical permeability of free space, the relative electrical permeability of the electrolyte, and the permanent charge density (where C E = 0 if there are no fixed charges in Ω E ), respectively. Finally J ± are the current densities, and q, K B , T are the electron charges, the Boltzmann constant and the (constant) temperature, respectively.
The charge distribution in a semiconductor can be described using the Drift Diffusion system (DD) coupled to the Poisson's equation, which is formally equivalent to the Poisson-Nernst-Planck system described above (P-DD). For a complete review about the mathematical problem and industrial application of these models, see, for example, [26] and [27].
It is well known that the DD system is not able to describe the behaviour of ultraintegrated devices, where quantum effects become non-negligible. The characteristic dimensions of the bio-device are often in the order of a few micrometers, and then the quantum effects can influence the dynamic of the system. To include these phenomena, it is necessary to add a correction term to the electrostatic potential, called Bohm potential where and m * are the reduced Planck's constant and the equivalent electron mass respectively. The so modified DD system is called quantum drift-diffusion (QDD) system. In this paper we consider the steady state QDD system in the following form where µ n , µ p , J n , J p , n, p, S , C S , B are the electron mobility, the holes mobility, the electron current density, the holes current density, the electron density, the holes density, the semiconductor permeability, the doping profile in Ω S and the Bohm potential, respectively.
From the mathematical point of view we study two systems of ordinary differential equations, one describing the electrolyte behaviour and the second one describing the quantum semiconductor, coupled by using a suitable set of interface conditions. The existence and uniqueness of the solutions follow by standard fixed point arguments.
Here the outline of the paper follows. In Section 2 we introduce a mathematical model to describe the behaviour of an EOS device. We account for the effect of the size exclusion in the domain of the electrolyte and discuss a suitable set of transmission conditions linking the three domains. Section 3 is devoted to the mathematical analysis of the model, in particular we prove the existence and uniqueness of solutions for our model. Finally, in Section 4 we test the model numerically, in order to quantify the effects of the size exclusion on the charge density and on the electrical potential behaviour.
2. EOS capacitor: Mathematical model. In this section we examine the system of equations applied in this paper in further detail. Systems (P-DD) and (P-QDD) (modelling the electrolyte solution and the semiconductor, respectively) are simplified in order to allow a mathematical treatment of the EOS capacitor. First of all, we limit our analysis to the steady state regime and we assume that the oxide is a perfect insulator and there is no charge flow through it, which means J i = J n = J p = 0, i = ±. Moreover, we consider a unipolar n-type semiconductor, leaving out the contribution of the holes. Qbviously all the results presented in this paper remain valid for a p-type semiconductor.
Under these approximations P-DD and P-QDD become respectively Here and in what follows ρ± is defined in Ω E , V on the whole domain Ω and n in Ω S .
The oxide can be simply described using the Poisson's equation, that is where the left-hand side is equal to zero due to the absence of both free and fixed charges.
The equations above can be written in dimensionless form, introducing a suitable set of reference variables where L is the total device length,ρ is the concentration of ions in the electrolyte solution at equilibrium and the label s indicates the dimensionless quantities.
In view of (2), we obtain where Here η is the scaled Planck's constant and (P-DDs) is written in conservative form. From now on: (A1) We divide the domain Ω ⊂ R into three parts, Ω E , Ω O , Ω S , which are the electrolyte, the oxide and the semiconductor domain, respectively.
2.1. Boundary and interface conditions. From now on all the variables taken under consideration will be assumed as already scaled, although we will omit the index s for simplicity.
In this subsection we are going to briefly summarize the boundary and interface conditions for our hybrid system. In order to understand the way in which the geometrical properties of Ω E influence the dynamic of the system, we consider two different configurations: I. The electrolyte domain is supposed to be quite large in comparison with the total volume occupied by the ions. This condition is verified, for example, in the EOS capacitor, which is the main component of the EOSFET transistor. II. The electrolyte domain is a narrow channel (nano-pore), and the size exclusion effects must be included. This case models the theoretical behaviour of a single (open) ionic channel. In both cases we can assume that the bulk equilibrium condition holds at least in a point x B ∈ Ω E . This condition is achieved for a large enough distance from both the charged surfaces (reference electrode and oxide-electrolyte interface) [39]. On the other hand it seems reasonable to assume that the system reaches its equilibrium condition in the external bath in the case of nano-pore device. Thus we assume that there exists x B such that Without loss of generality we set V = 0 as reference potential. If the electrode has zero potential, the bulk equilibrium condition is archived on ∂Ω E . Otherwise, if V |∂Ω E = V E = 0 the behaviour of the electrical potential around the electrode (Gouy-Chapman layer) can be computed exactly by solving the system P-DDs with the following boundary condition: We remark that the original theory due to Gouy-Chapman, despite being extremely simple and fairly realistic, overestimates the value of the interface charge densities. Therefore the Stern layer must be taken into consideration if we want to obtain more realistic results. Indeed, it models the fact that the distance between ions and the electrode surface cannot be smaller than the diameter of the ions. See [39] for a complete review about the EOS capacitor in the classical case. Without loss of generality, we assume that the bulk condition holds in correspondence of the electrode (V E = 0), and we remark that all the other cases can be reduced to this one as described above.
Thus, we can prescribe the values of the electrical charge density and the value of the electrical potential on the boundary of the domain. We refer to these quantities as physical boundary conditions: Moreover, the quantum effects are considered negligible on the metallic contacts: B |∂Ω S ≈ 0. Finally, we have to prescribe the interface conditions. We assume that the electric potential and the normal component of the electric displacement vector (x)∇V are continuous on both the interfaces: where [X] |x=x0 indicates the jump of the function X at x = x 0 . The last interface condition required for this purpose is From the theoretical viewpoint, it should be n| ∂Ω SO = 0, since the electrons cannot penetrate into the oxide, due to the presence of the quantum barrier. However, to avoid numerical instabilities, we consider 0 < n |∂Ω SO = n O << 1 as suggested in [13]. This parameter can be read as a measure of free carrier transport through the oxide barrier (tunnel effect).

2.2.
From the physical problem to the mathematical model. In this section we rewrite (P-DDs) and (P-QDDs) in a more compact way to simplify the mathematical analysis, which will be the subject of the next section.
In view of the boundary conditions previously introduced, the first two equations in (P-DDs) and the first one in (P-QDDs) can be integrated over Ω E and Ω S , respectively. In this way we obtain: The first two equations in (9) are the well known Boltzmann distributions.
By substituting them into the Poisson's equation in (P-DDs), and making some simple calculations, we obtain which is usually called Boltzmann-Poisson equation.
For small applied potentials, the Debye-Huchel approximation holds, and (10) assumes the linearized form (see [39]) A more general case of (10) includes the size exclusion effect, for ionic flow through a narrow structure such as a nano-channel: See [11] and its list of references for details about its derivation.
However, for small applied potentials, we have which is again in the form (11). From now on, just to simplify the notation, we assume E = α E and refer to equations (10), (11), (12) and (13) as where with size exclusion, In the last line we assume κ 2 = 2, without the size exclusion, and κ 2 = 2/3 otherwise.
Solving (14), we obtain the behaviour of the electrical potential V . The charge densities can be also determinated by using the first two equations in (9).
Considering the Poisson's equation in (P-QDDs), (3), and (14), we get We can solve the system above associating the following boundary and interface conditions The electrons behaviour in the semiconductor domain is described by the first equation in (P-QDDs). Until now, we have used the isothermal form of the enthalpy, namely h(s) = ln s (for any positive defined function s), in order to get a more readable physical derivation and to underline the relationship between (P-QDDs) and (P-DDs). In fact, in that case, it is possible, at least formally, to write a Boltzmann-like distribution also for electrons in a quantum system, as shown in the last equation of (9). However, from the mathematical point of view, a more general enthalpy function seems more appropriate. To this end, we introduce a more general form for the enthalpy h(s), that is Finally, in order to simplify the mathematical approach in the next section, we introduce the new variable w = √ n, and from now on we prefer switching to the new unknowns (w, V ), instead of (n, V ). We remark that, by definition, the charge density n(x) is a non negative function, therefore w is well defined and non-negative.
The first equation of (P-QDDs) and the correspondent boundary condition can be rewritten as where w OS = √ n OS .
Limiting our analysis to the one-dimensional case, from (18), we get the following description of the charge concentration (see Fig. 2 for the notation): The last equation in (19) models the zero-charge distribution in the oxide layer. Analogously, by using (16), for the electrical potential, we have: in Ω E , We would like to remark that the 1-D case can be derived directly from the general 3-D system, by integration with respect to x and assuming constant device cross sections. Despite their simplicity, these results can be directly correlated with the measured quantities, such as the surface charge distributions at the interfaces. In our opinion this model may represent a good starting point for the development of more complex multidimensional models.
3. Existence of solution to the hybrid electrolyte-semiconductor system. In this section we discuss the existence and the uniqueness of solutions for the electrolyte-oxide-semiconductor problem (19) and (20). We consider a bounded 1-D domain, as in Fig. 2. Basically, we are going to use the same approach proposed in [18]. We start with the following lemma, where we discuss some a priori estimates, which will be used to establish the existence of solutions (Theorem 3.3).

Charge concentration
Electrical Potentials in Ω E , where w K = max(0, min(w, K)), and w OS , w S > 0 and 0 < λ, η << 1. Then, there exists K 0 > 0 such that for all K ≥ K 0 we have where the constant c depends only on the parameters of the problem.
Proof. Let us assume w 0 a regular function which verifes the boundary and the interface conditions for w, for example the linear interpolation By multiplying the first equation in (21) by (w − w 0 ) and integrating on Ω S , after some manipulations, we get

FEDERICA DI MICHELE, BRUNO RUBINO AND ROSELLA SAMPALMIERI
Now we take V 0 as a given regular function, with compact support in Ω (verifying the interface conditions (22) 5 ) and multiply the sum of the first, the second and the third equations in (22) by (V − V 0 ). In view of the boundary condition and the transmission conditions at the interface, after integrating on Ω, we get: Using the Cauchy-Schwarz inequality (with a suitable parameter) to estimate the first term on the right hand side of the equation (24), we obtain From this point onwards c and c i are constants independent from η and K.
We are now going to discuss in details the last integral in the inequality above. In fact, the behaviour of the function f (V ) depends strictly on the geometrical properties of Ω E and on the applied potential, as explained in the previous section and summarized in (15).
Let's denote 1. When f (V ) = 2 sinh(V ), we have: where we observe that (sinh(y) − sinh(y 0 ))(y − y 0 ) ≥ 0, for all y ∈ R. By using the Cauchy-Schwarz inequality (with parameter) and the Poincaré inequality, we get: Our second case assumes f (V ) = sinh(V ) 1 2 +cosh(V ) . Observing that |f (y)| < 1 for all y ∈ R, we get where c depends on the maximum of the function f (V ) and c P is an ad hoc Poincaré constant. Lastly, working as in (26), we get 3. Our third case considers a small applied potential. Then f (V ) = κ 2 V , corresponding to the last case in (15). In this case we have By adding (23) and (25), in view of the computations above and after some manipulations, we get which holds for each function f (V ) in (15). The last two terms in (27) can be estimated, using the Poincaré inequality, as follows: where c P is, as usual, the Poincaré constant.
Subsequently, taking into account (28), the estimate (27) reduces to Following the same arguments as in [18], there exist K 0 > 0 and M 0 > 0 such that, for all K ≥ K 0 , we have From (29) there exist c Ω , c V , c > 0 depending on the boundary values, on (x) and on the doping profile C(x) such that and In order to construct our fixed point procedure, we need the following elementary result.
Proof. Clearly u is given by By using the a priori estimate discussed above, we prove the existence of a unique solution to our problem. In the following theorem we discuss the existence of solutions applying the Leray-Schauder fixed-point theorem.  (20) admit at least one solution (w, V ) such that w ∈ H 2 (Ω S ) ∩ L ∞ (Ω S ) and V ∈ H 1 (Ω) ∩ L ∞ (Ω). Moreover, there exist three constants V ± and ω + such that Proof. We start defining the map T : M −→ M , where that will be shown to verify the assumptions of Leray-Schauder fixed point theorem.
For each (u, U ) ∈ L 2 (Ω S ) × L 2 (Ω), let V ∈ H 1 (Ω) be the unique solution to the following linearized system in Ω S , where, from now on, we take σ ∈ We introduce the function D(x), called electrical displacement, which can not jump at the interfaces, where its value is From (19) and (20) it is not difficult to show that V (x EO ), V (x OS ) and D verify the following conditions where The unique solution of the linear system (34) can be expressed in the form where A i , B i , C i , i = 1, 2, 3 are constants. Moreover, by the previous estimates, we can find two constants K 1 , K 2 such that This allows us to define three closed and bounded intervals [  Once we have V , we obtain w ∈ H 1 (Ω S ) as the unique solution of Here u K is defined as in Lemma 3.1. Then T maps L 2 (Ω) × L 2 (Ω S ) into itself and it is well defined. Now we have to check that T is continuous and compact. The continuity derives from the continuity of the right hand side of the problems (33) and (36). Concerning the compactness of T , take (u, U ) ∈ B, where B is a bounded set of M . For (ω, V ) = T (u, U ) we can obtain the estimate likewise for (31). This tells us that T (B) is a bounded subset of H 1 (Ω S ) × H 1 (Ω) and then it is a pre-compact set in L 2 (Ω S ) × L 2 (Ω). Moreover, for σ = 0 we have T ((u, U )) = (0, 0) and for (w, V ) ∈ M such that (w, V ) = T (w, V ) we have obtained similarly to (31), where c 1 does not depend on σ ∈ [0, 1] and K. Finally, we can apply the Leray-Schauder fixed-point theorem to get a weak solution to problem (22)- (21). The uniqueness of this solution will be proved in Theorem 3.4. Once (w, V ) = T (u, U ) has been obtained, we can update the values of I E and I S and then the value of the electric potential at the interfaces by means of (35). The same procedure is applied to the electrical displacement D, taking the initial guess in [D − , D + ].
This allows us to define a continuous application T 1 : ). Then one can apply to the map T 1 the Schauder's fixed point theorem to obtain the final value of V and V x at the interfaces.
In order to find the solution of the original problems (19)- (20) we have to show that w ∈ L ∞ (Ω). This simply follows from the immersion H 1 (a, b) → L ∞ (a, b). In the same way, (32) 1 follows directly from (38). It remains to show that w ≥ 0. Proceeding as in [18], we multiply (21) by w − = min(0, w) ∈ H 1 (Ω S ), and integrate by parts getting which according to the definition of w K , implies w ≥ 0. Finally w ∈ H 2 (Ω S ) from (19), in view of the estimates above.
In the following theorem we present the uniqueness result. Proof. Concerning the bound of the charge density, it is possible to show that for all x ∈ Ω S the following inequalities hold: where w and w are constants not depending on η. See [18] (page 52 and the following) for the complete proof.
To get the uniqueness of the solution, we work by contradiction. We discuss the details of isothermal case, where h(s) = ln(s). We assume there exist (w 1 , V 1 ) and (w 2 , V 2 ) solving the problems (19)- (20). Subtracting the system solved by (w 1 , V 1 ) and (w 2 , V 2 ) we get: We multiply the first equation in (39) by (w 1 − w 2 ) and obtain It is not difficult to show that therefore, one has As a consequence of L ∞ bounds for V and (H1), one has Moreover, the sum of the Poisson's equations after multiplication by (V 1 − V 2 ) gives us By summing (42) and (43), in view of the monotonicity of the functions ln(w) and f (V ), we get The uniqueness follows for η 2 large enough, by using the Poincaré inequality. Finally, we observe that in the isentropic case, that is h(s) = γ γ−1 s γ−1 and γ > 1, the uniqueness can be proved in the same way, since 4. Numerical simulations. Recently some new models describing devices based on interfacing nano-sized semiconductors and electrolyte solutions have been developed (see for example [8], [7] and reference therein). Here the PNP is used in its standard form without size exclusion effects.
The main purpose of this section is to understand how (and if) the size exclusion in the electrolyte domain modifies the electrons distribution and the behaviour of the electrical potential. Then, before testing our model on more complex devices, for which experimental results are available, we start with a simple toy model.
The numerical simulations presented in this section are produced by COLNEW, a SCILAB function for boundary value problems which use collocation at Gaussian points [3] (see also [1], where the old package called COLSYS is described). In [18] COLSYS has been shown to provide good results in approximating quantum hydrodynamical equations.
We consider the same simplified domain introduced in the theoretical section, but the approach used to get solutions is a little bit different from the one presented in Theorem 3.3 in order to simplify the numerical approach. Indeed, in this section we solve the three problems each in its own domain, while, in the Theorem 3.3, first we solve globally the equations for the potential and then the electron density is evaluated in Ω S .
The potential V and the function w = √ n are the solutions of the following boundary value problems Boundary value problem in Ω S Boundary value problem in Ω E whereas the potential in the oxide domain is given by Systems (44) and (45) are solved by using COLNEW starting with two initial guesses V 0 (x EO ) and V 0 (x OS ). The value of the potential on Ω O is given by (46). This provides the values of V on whole domain and of w in Ω S . Then we can compute the new values V 1 (x EO ) and V 1 (x OS ) by using (34) and iterate the process until convergence of the interface potentials V (x EO ) and V (x OS ). In particular, we assume that the convergence is achieved if: where k indicates the iteration step.
The convergence of this iteration can be proved using a similar approach as in Theorem 3.3.
For the simulation, we assume C S = 10 17 m −3 , L = 20 nm, 0 = 8.8510 −12 F/m, T = 300 K, O = 3.9, S = 11.9, E = 80, m * = 10 −31 Kg, = 6.6 · 10 −34 J/s and q = 1.6 · 10 −19 C. The scaled electron concentration at the interface between semiconductor and oxide is taken w 2 OS = 0.001. It seems reasonable that the ions confinement effects are important when studying high ions concentrationρ in the electrolyte domain. Thus, in order to quantify the effect of the size exclusion, we consider different values of the parameter α = C S /ρ. Obviously the value ofρ is uniquely determined by α and C S .
Assuming α = 1 and V S = −0.6 V (corresponding to the scaled value V S s = −23.076923), the presence of the size exclusion does not modify the electrical potential and the charge density behaviour, as summarized in Fig. 3 and Fig. 4.
To observe the effects of the size exclusion, we reduce the value of the parameter α until 10 −7 . As a consequence, the electrical potential decreases at the interface electrolyte-oxide (Fig. 5). At the interface between oxide and semiconductor, the potential does not change remarkably (V (x OS ) = −0.1355491 V and V (x OS ) = −0.1325997 V with and without size exclusion, respectively), but the effects of the size exclusion on the electron concentration at the interface are evident as shown in Fig. 6. On the contrary, the effects on the potential on whole domain are not remarkable for this experimental parameters (Fig. 7). Finally we reduce the value of the applied potential until V S = −0.06 V (corresponding to the scaled value V S s = −2.3076923). In this case we expect that both the ions and the electrons are very close to the equilibrium conditions and the size exclusion effects does not play any role as showed in Fig. 8 and Fig. 9.

5.
Conclusion and future works. In this paper we have proposed a mathematical model to describe an EOS capacitor working in a quantum regime. The quantum effects are introduced, as usual, by means of the Bohm potential. We have proved the existence of a unique solution in the 1D steady state case and we have discussed how the presence of a narrow channel in the domain Ω E influences the behaviour of the device. The model has been validated on a test device, by using COLNEW and a simple algorithm based on a fixed point argument, discussed in the theoretical part. Our simulations show that, as expected, the effects of the size exclusion depend on the valueρ, which is the ions concentration at equilibrium.ρ is the main parameter to consider in order to control the behaviour of the device, at least in the 1D approximation.
We point out that we do not take into consideration the electrons penetration in the oxide and electrolyte domain, although this effect can strongly modify the behaviour of the system. To include such phenomena, a more general time dependent model has to be developed.
In conclusion, we remark that, only the single channel case has been examined here. In the more complex (and probably more realistic) case the electrolyte moves in a porous or nanoporous media. To describe this case, many other effects must be considered, such as the properties of nanoporous structure, the interaction between ions and the internal side of the nanoporous surfaces and the reduction of the porosity due to the presence of ions. However, the model presented in this paper seems to the authors a good starting point for developing a more complex model, which could be able to describe exhaustively most of the modern bio-devices.