Next Article in Journal
Convolutional Neural Networks for Challenges in Automated Nuclide Identification
Previous Article in Journal
A Deep-Learning Based Posture Detection System for Preventing Telework-Related Musculoskeletal Disorders
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

A Semi-Linear Elliptic Model for a Circular Membrane MEMS Device Considering the Effect of the Fringing Field

by
Mario Versaci
1,*,†,
Alessandra Jannelli
2,†,
Francesco Carlo Morabito
1,† and
Giovanni Angiulli
3,†
1
DICEAM Department, “Mediterranea” University, I-89122 Reggio Calabria, Italy
2
MIFT Department, University of Messina, I-98166 Messina, Italy
3
DIIES Department, “Mediterranea” University, I-89122 Reggio Calabria, Italy
*
Author to whom correspondence should be addressed.
This work is the result of a research project shared between the authors of this paper which deals with numerical methodologies for the recovering of membrane MEMS devices.
Sensors 2021, 21(15), 5237; https://doi.org/10.3390/s21155237
Submission received: 3 July 2021 / Revised: 29 July 2021 / Accepted: 29 July 2021 / Published: 2 August 2021
(This article belongs to the Section Physical Sensors)

Abstract

:
In this study, an accurate analytic semi-linear elliptic differential model for a circular membrane MEMS device, which considers the effect of the fringing field on the membrane curvature recovering, is presented. A novel algebraic condition, related to the membrane electromechanical properties, able to govern the uniqueness of the solution, is also demonstrated. Numerical results for the membrane profile, obtained by using the Shooting techniques, the Keller–Box scheme, and the III/IV Stage Lobatto IIIa formulas, have been carried out, and their performances have been compared. The convergence conditions, and the possible presence of ghost solutions, have been evaluated and discussed. Finally, a practical criterion for choosing the membrane material as a function of the MEMS specific application is presented.

1. Introduction

The development of embedded technologies in the last decade has been mainly due to microdevices, which can manage the connection between the physical layer of a particular problem and the logic of the machine language [1,2,3,4,5,6]. In this context, the interest in both static and dynamic MEMS devices is high [7,8,9]. In particular, a lot of effort have been made concerning dynamic investigations into the response of resonant microbeam [10] as much as micro-resonators [11,12]. This because they provide an excellent approximation of the human–machine interface in all the cases wherein miniaturized-integrated electromechanical systems are required [13,14,15,16,17,18,19,20,21]. Since the first batch device was produced in 1964 [22], the scientific and technological development in the field has strongly influenced the analysis and synthesis of physical-mathematical models that can describe the extremely complex underneath MEMS multi-physics [23,24,25,26,27,28,29]. However, such theoretical models hardly provide analytical solutions. Hence, it appears necessary to have conditions that can ensure both the existence and uniqueness of the solution [23,24,25,26,30] to avoid ghost solutions when the model is solved numerically [30,31,32]. Since the aforementioned analytical solutions are unobtainable, suitable numerical procedures have to be selected, and it becomes imperative to evidence the pros and cons of each method [25,26,30,31,32]. MEMS is a rampant technology employed for realizing thermo-elastic systems [1,33,34,35,36,37], and it has a wide variety of applications ranging from biomedical engineering to microfluidics [38,39,40,41,42]. Moreover, many researchers are actively engaged in the development of important experimental research works for the development and prototyping of special MEMS such as, for example, circular graphene membrane MEMS devices [43,44], SiN circular membrane MEMS devices [45,46], and CMOS MEMS-based membrane-bridge devices [47] particularly useful for industrial applications. Moreover, the scientific community is intensively working on the analysis/synthesis of multi-physical models characterized by a high degree of symmetry, because these are more easily achievable from a technological point of view [24,30,48,49,50,51]. Accordingly, we have focused our attention on a 2D circular membrane MEMS device, a kind of geometry widely used in many industrial applications [24,25,48,52]. The membrane deformation is described by its profile, u ( r ) , r [ 0 , R ] with r being the radial coordinate and R the radius of the device. As it is well known, for a circular geometry, the physical-mathematical model describing the membrane curvature assumes the following semilinear elliptic form [1,53]:
Δ u ( r ) = λ 2 ( 1 u ( r ) ) 2 u ( r ) = 0 , u ( 0 ) = 0 , 0 < u ( r ) < d
where λ 2 is a parameter related to the external electrical voltage V and d is the distance between the parallel disks. (1) is a theoretical model disregarding the fringing field phenomenon [54,55,56,57,58] that occurs when d R [17]. A more precise analytic model, from which (1) is derived, is the following [59]:
Δ u ( r ) = λ 2 ( 1 u ( r ) ) 2 + λ 2 F ( u ( r ) , u ( r ) , δ , ) u ( L ) = u ( L ) = 0 .
where F is a suitable function describing the effects close to the MEMS edges, and δ 0 is a term that weighs the effects of the fringing field. It is apparent that (1) derives directly from (2) omitting the additional term F ( u ( r ) , u ( r ) , δ , ) under the hypothesis that the effect due to its presence is negligible. However, this hypothesis turned out to be false despite using uniform approximation for the electrostatic field E . On the basis of this premise, we employed the so-called corner-corrected theory, developed in [59], for our analysis, thus exploiting the “corner-corrected model”; more precisely, following the work in [59], the term F ( u ( r ) , u ( r ) , δ , ) in (2) can be written as
λ 2 δ | u ( r ) | 2 .
Now, considering that Δ u ( r ) in (1) has only the radial component, (1) specializes as [48]
u ( r ) = 1 r u ( r ) λ 2 ( 1 + δ | u ( r ) | 2 ( 1 u ( r ) ) 2 u ( R ) = 0 , u ( 0 ) = 0 , 0 < u ( r ) < d
(It can be noticed that if δ = 0 , no fringing field occurs and (4) reduces to (1)). At this stage, we point out that the term δ | u ( r ) | 2 becomes appreciable only if the term | u ( r ) | assumes relevant values, a real possibility at the edge of the device in most cases (see Theorem 1 in [23]). Once the membrane deforms, the electrostatic capacitance of the device, C e l , varies since the distance between the membrane and the upper disk is not constant. Moreover, as physically
λ 2 ( 1 u ( r ) ) 2 | E | 2
it follows that
λ 2 ( 1 u ( r ) ) 2 = θ | E | 2 , θ R + .
As proved in [48], the field E is locally orthogonal to the straight-line tangent to the membrane profile itself, so that it is natural to consider | E | K ( r , u ( r ) ) with K ( r , u ( r ) ) mean curvature of the membrane [24]:
K ( r , u ( r ) ) = 1 2 u ( r ) + 1 r u ( r ) .
Considering both (6) and (7), (4) can be rewritten as
u ( r ) = 1 r u ( r ) 4 ( 1 u ( r ) d ) 2 θ λ 2 ( 1 + δ | u ( r ) | 2 ) u ( R ) = 0 , u ( 0 ) = 0 , 0 < u ( r ) < d ,
which is a 2D second-order differential semi-linear elliptic model without explicit nonlinearity. The term d , the so-called critical security distance, guarantees that the membrane does not touch the upper disk. Note that the factor θ , by (6), is an important parameter because it takes into account the electrical properties of the membrane. This last differential model is exceedingly difficult to solve analytically. Accordingly, we focus our attention on obtaining conditions that ensure both the existence and uniqueness of its solutions. Before proceeding with our analysis, we point out how our study is framed in a broader line of research. The first paper was that in [24] in which a 2D nonlinear second-order differential model for electrostatic circular membrane MEMS devices (wherein the singularity was not explicitly evident) was studied, starting from (1), from which the following chain of proportionality,
λ 2 / ( 1 u ( r ) ) 2 | E | 2 K 2 ( r , u ( r ) )
was proved. As the model studied in [24] did not allow to obtain an explicit analytical solution, the existence was derived through the use of two auxiliary functions satisfying particular properties. In this way, an algebraic condition on the product θ λ 2 , depending on both d and V, was provided. However, as far as uniqueness is concerned, it has not been guaranteed [24]. In [60], stable numerical approaches for recovering the membrane profile, based on the three-stage Lobatto formulas, have been exploited, thus obtaining the ranges for operative parameters and the areas of applicability of the device avoiding the ghost solutions. In contrast, in [25], the equilibrium configurations were analyzed. Finally, in [48], the model was improved and rewritten. In particular, in this study, the effect of the fringing field on membrane curvature has been considered. Furthermore, the study into the stability of the equilibrium configurations, the analysis of the range of possible values for V, and the study on the device’s optimal control was conducted.
The main results presented in this paper can be summarized as follows:
  • A new algebraic condition governing the uniqueness of the solution for (8), depending on the electromechanical properties of the membrane material, has been demonstrated. Unlike the 1D geometry, this new condition cannot guarantee both the existence and uniqueness of the solution for the model (8).
  • Shooting procedure, Keller-box scheme, and III/IV Stage Lobatto IIIa formulas have been employed, and their numerical performances, related to the membrane profile recovering task, when δ varies in the range of its possible values, have been compared. Furthermore, the values of the parameter θ λ 2 ensuring the procedures’ convergence have been determined.
  • Ghost solutions have been investigated for obtaining the values of the factor θ λ 2 that ensures the convergence of each considered numerical procedure, avoiding the ghost solutions’ computation.
  • Finally, the relationship among the numerical convergence criteria, the parameter θ λ 2 , and the intended use of the device has been highlighted.
The paper is organized as follows. In Section 2, the 2D electrostatic circular-membrane MEMS device considered in our analysis has been mentioned. Section 3 discusses the curvature-dependent | E | for modeling the 2D circular membrane MEMS device with the fringing field. In Section 4 and Section 5, we recall well-known results of the existence and uniqueness for the problem under study. Section 6 suggests a new result concerning the uniqueness of the solution, which depends on the electromechanical properties of the material constituting the membrane, whose proof is reported in the Appendix A. In Section 7, we derive a new condition ensuring both the existence and uniqueness of the solution. In Section 8, we present and discuss the numerical results conducted in our study for recovering the MEMS membrane profile. Finally, in the last section, our conclusions are given.

2. A Description of the 2D Electrostatic Circular-Membrane MEMS Device

2.1. The Point of View of the Actuator

In the usual 3D Euclidean space, R 3 , we consider a system of orthonormal Cartesian axes O x y z , where the z axis represents the vertical axis (see Figure 1). The device consists of two parallel circular metal disks with radius R, located at a distance of d from each other. The lower disk is located on the x y plane, and its center coincides with the origin of the system O x y z . The device is subjected to an external electrical voltage V, and the lower disk is fixed at the potential V = 0 . A circular membrane with radius R is anchored to the edges of the lower disk. It deforms toward the upper disk when a voltage V is applied. If we denote with r ( 0 < r R ), the radial coordinate, we have that the membrane profile, u ( r ) , results to be a function of r. To overcome its mechanical inertia, the applied voltage V must assume values such that the corresponding value of the electrostatic field E inside the device can generate an appropriate electrostatic pressure equal to p e l = 0.5 ϵ 0 | E | 2 (with ϵ 0 , the permittivity of free space). The term p e l can be translated into an equivalent electrostatic force, f e l , computable as [1,25,48,61]
f e l = 0.5 ϵ 0 π R 2 V 2 ( d u ( r ) ) 2 ,
which deflects the membrane, thus achieving a displacement in its center, u 0 , equal to R 2 p e l / 4 T (where T is the radial mechanical tension of the membrane when it is at rest) [1,48]. Electrostatically, if the membrane deforms, the field E , which depends on the distance between the membrane and the upper disk, results to be locally orthogonal to the tangent line to the membrane at the same point [2]. Moreover, the electrostatic capacitance, C e l , is variable as the distance between the membrane and the upper disk varies locally. We can also observe that the greater the | E | , the greater the curvature of the membrane. On this basis, we have that | E | can be considered locally proportional to the mean curvature K ( r , u ( r ) ) of the membrane [23,24,25,30].
The device geometry being such that L d , we have that the effect due to the fringing field cannot be neglected [48]. Accordingly [1],
λ 2 = ϵ 0 V 2 ( 2 R ) 2 2 d 3 T = ρ V 2
in which ϵ 0 is the permittivity of the free space.
Remark 1.
λ 2 , being T-dependent, is a parameter that expresses the mechanical properties of the membrane. Therefore, θ λ 2 expresses the electromechanical properties of the membrane.

2.2. The Point of View of the Sensor

In the following, it will be convenient to exploit some similarities with the circular plate MEMS transducer model subject to mechanical pressure p.
As is known, the study of the deformation of membranes in MEMS devices starts from considering metal plates subject to mechanical pressure p whose deflections u satisfy [1,62]:
ρ h u t t T Δ u + D Δ 2 u = 0
where ρ is the density of the material constituting the deformable plate, and h and T are the thickness and mechanical tension, respectively. Moreover, indicating the Young’s modulus with Y and the Poisson ratio with ν , the stiffness coefficient D assumes the following form [62]:
D = Y h 3 12 ( 1 ν 2 )
Furthermore, if the plates are circular, u only depends on r so that in the steady-state condition results [1,62],
u ( r ) = R 4 64 D 1 r R 2 2 p
However, if r = 0 , it follows that u 0 = R 4 p / 64 D , and from (14),
u ( r ) = u 0 1 r R 2 2 p
In this case, the device works as a transducer. In fact, p deforms the membrane and generates u ( r ) 0 profile (except at the device’s boundaries). Then, C e l becomes [1,62]
C e l ( u 0 ) = 0 R 2 ϵ 0 π r d 1 u ( r ) d d r , with | u 0 | d .
As both h and D are bounded values, u ( r ) becomes unobtrusive. Therefore, the distance between the two plates remains almost constant and equal to d. Moreover, exploiting the Taylor series and the value of C e l at equilibrium, C 0 = ϵ 0 π r 2 d , for p = 0 , (16) can be written as
C e l ( u 0 ) C 0 1 + u 0 3 d + u 0 2 5 d 2
by which it is possible to obtain the electrostatic charge of the membrane, the co-energy of the system, and the electrostatic force. Finally, it is also possible to obtain
| E ( r ) | V d u 0 1 r R 2 2 .
Remark 2.
The physical quantities involved here clearly depend on d, because the circular plate had a significantly high value of D; further, u ( r ) appeared extremely limited so that any dependency on d u ( r ) could be replaced by the dependence on d (which is mathematically more straightforward and, therefore, easier to manage).
If we consider a membrane instead of the deformable plate, h can be neglected. D then decreases substantially (see (13)) in the case in which a deformable plate is considered. Bearing in mind that the smaller the D value, the more flexible the membrane is, it follows that u 0 becomes higher with an increased risk of the membrane accidentally touching the upper disk. If the membrane replaces the deformable disk, u ( r ) takes the form [1]
u ( r ) = u 0 1 r R 2
and
u 0 = p R 2 4 T
with f e l as formulated in (10).
Remark 3.
To calculate f e l and p e l for evaluating the surface of the membrane, π R 2 was considered as if the membrane was at rest. This approximation can be justified because d R such that the surface of the deformed membrane is approximately equal to the surface of the membrane at rest.
The parameters p and p e l are linked to each other because applying V produces the field E in the device which, in turn, generates p e l , deforming the membrane toward the upper plates. Moreover, indicating by k 1 = R 2 4 T , from (20), we can write
u 0 = k 1 p .
We can assume that, in our case, p comes exclusively from p e l due to | E | inside the device. Thus, p depends on p e l so that the following chain of equalities have to be considered valid:
u 0 = k 1 p = k 1 k 2 p e l = k p e l .
As we will see below, k plays an important role in formulating the algebraic condition that governs the existence of the solution for (8) (see, inequality (45)).
Remark 4.
When the membrane is at rest, the distance between the membrane and the upper disk is d. Therefore, C e l along any vertical plane σ ¯ whose support is the straight line
r ¯ : x = 0 y = 0 ,
with fringing field phenomenon, is [54]
( C e l ) C = 2 ϵ 0 R d 1 + d 2 π R ln 2 π R d ,
where C represents the curve that arises from the intersection σ ¯ and the membrane. Therefore, the total C e l , ( C e l ) t o t a l can be written as follows:
( C e l ) t o t a l = 0 π B ( ϕ ) ( C e l ) C d ϕ = ( C e l ) C 0 π B ( ϕ ) d ϕ .
where B ( ϕ ) is a bounded and continuous electrostatic function based on ϕ [63]. Therefore,
0 π B ( ϕ ) d ϕ = D ¯ < + ,
so that (25) becomes
( C e l ) t o t a l = ( C e l ) C D ¯ = 2 ϵ 0 R D ¯ d 1 + d 2 π R ln 2 π R d .

3. The Mathematical Model

The model (1), with fringing field phenomenon, becomes
Δ u ( r ) = λ 2 ( 1 + δ | u ( r ) | 2 ) ( 1 u ( r ) ) 2 u ( R ) = 0 , u ( 0 ) = 0 , 0 < u ( r ) < d .
Moreover, the radial symmetry in it with respect to vertical axes r = 0 allows to write
Δ u ( r ) = 1 r u ( r ) + u ( r ) ,
and considering u ( r ) = u ( r ) , we can write (28) as
1 r u ( r ) + u ( r ) = λ 2 ( 1 + δ | u ( r ) | 2 ) ( 1 u ( r ) ) 2 u ( R ) = 0 , u ( 0 ) = 0 , 0 < u ( r ) < d .
Furthermore, putting (5) in (30), the (6) holds. Accordingly, (30) becomes
u ( r ) + 1 r u ( r ) = θ | E | 2 ( 1 + δ | u ( r ) | 2 ) u ( R ) = 0 , u ( 0 ) = 0 , θ R + , 0 < u ( r ) < d .
As observed in [24,48], E on the membrane is locally orthogonal to the straight-line tangent to the membrane’s profile of. Thus, | E | can be considered locally proportional to the mean curvature K ( r , u ( r ) ) of the membrane (see (7)). Of course, the greater the | E | , the greater the deformation of the membrane. Therefore, we can assume that | E | K ( r , u ( r ) ) ; based on this premise, we can write
| E | = μ ( r , u ( r ) , λ ) K ( r , u ( r ) )
in which μ ( r , u ( r ) , λ ) represents the function of proportionality denoted as [24,48]
μ ( r , u ( r ) , λ ) = λ ( 1 u ( r ) d )
where μ ( r , u ( r ) ) C 0 ( A ) , and A = [ R , R ] × [ 0 , 1 ) .
Remark 5.
Usually, the models closest to the physical reality of MEMS are highly complex and cannot be faced analytically. It is thus necessary to operate with simplifications in the geometry of the devices, thereby obtaining simplified models that can be studied analytically. In other words, the results obtained from studying (8) will hardly agree with the experimental data. However, they will give an excellent qualitative contribution to the device’s behavior even if it is characterized by a simplified geometry.
Remark 6.
From (6), we can observe that λ 2 is directly proportional to | E | 2 and therefore to V 2 . Then, λ 2 is a bounded parameter, because a minimum value of V is required to overcome the mechanical inertia of the membrane (see below). On the other hand, the value of V cannot increase indefinitely because the intended use of the device fixes its maximum admissible value. Then, λ 2 is also limited by the intended use of the device.
Analytically, θ has no limitations except that it must be non-zero. However, as we will see below, the θ λ 2 product will be subject to specific limitations due to problems related to the convergence of the numerical procedures used for recovering the membrane profile in the presence/absence of ghost solutions.
As far as μ ( r , u ( r ) , λ ) is concerned, no apparent limitations to the values that can be attributed have appeared.
Moreover, from (32), considering both (7) and (33), we obtain
| E | 2 = 1 4 λ 2 ( 1 u ( r ) d ) 2 u ( r ) + 1 r u ( r ) 2
so that (31) becomes
u ( r ) + 1 r u ( r ) = θ λ 2 4 ( 1 u ( r ) d ) 2 u ( r ) + 1 r u ( r ) 2 ( 1 + δ | u ( r ) | 2 ) u ( R ) = 0 , u ( 0 ) = 0 , θ R + , 0 < u ( r ) < d .
Additionally, from (35), (8) follows because u ( r ) + 1 r u ( r ) 0 (see [23]). It is a special case of the following general problem:
u ( r ) + F ( r , u ( r ) , u ( r ) ) = 0 u ( b ) = B , u ( a ) = m ,
where F C 0 ( ( a , b ] × R × R ) and B , m R . In fact, if
F ( r , u ( r ) , u ( r ) ) = 1 r u ( r ) + 4 ( 1 u ( r ) d ) 2 θ λ 2 ( 1 + δ | u ( r ) | 2 ) ,
B = m = 0
and
b = R , a = 0 ,
Equation (4) is easily achieved. The general formulation (36) is allowed to exploit two necessary lemmas (see Lemmas 1 and 2 in [24]) for achieving an algebraic inequality that, if satisfied, ensured the existence of the solution for (8). For the simplicity for reading, we report this critical result of the existence [48].

4. On the Existence of At Least One Solution

Theorem 1.
Let us consider Problem (8). In addition, let us take into account two functions continuously differentiable twice, u 1 ( r ) and u 2 ( r ) , both defined on [ 0 , R ] , in order that u 1 ( r ) < u 2 ( r ) . Furthermore, let us suppose that
u 1 ( r ) + 1 r u 1 ( r ) + 4 ( 1 u 1 ( r ) d ) 2 θ λ 2 ( 1 + δ | u 1 ( r ) | 2 ) > 0
and
u 2 ( r ) + 1 r u 2 ( r ) + 4 ( 1 u 2 ( r ) d ) 2 θ λ 2 ( 1 + δ | u 2 ( r ) | 2 ) > 0
for r ( 0 , R ) . If 1 r u ( r ) + 4 ( 1 u ( r ) d ) 2 θ λ 2 ( 1 + δ | u ( r ) | 2 ) is a continuous function (except for r = 0 ), which satisfies the Lipschitz condition in U × ( , + ) , with
U = { ( r , u ) : 0 < r < R a n d u 1 ( r ) u ( r ) u 2 ( r ) }
and if
u 1 ( 0 ) u 2 ( 0 )
and
u 1 ( R ) = u 2 ( R ) = 0
with
θ λ 2 > 2 d 2 R 2 k ϵ 0 V 2 1 + δ k ϵ 0 V 2 r d 2 R 2 2 ,
thus, Problem (8) admits at least one solution.
Proof of Theorem 1.
See [48]. □
Remark 7.
(45) represents the algebraic condition ensuring the existence of at least one solution for (8). (45) makes sense if the device is in operation ( V 0 ). If on the plane θ λ 2 r , we represent (45) in the form of an equation, and we obtain a curve of the type shown in Figure 2. Above the curve is the area of the plane that satisfies (45), while below the curve is the forbidden area. Further, if we take into account the limitation
0 δ < 2
obtained in [64], the trends of θ λ 2 , varying by δ, for 0 r R are displayed in Figure 2, which indicates that as δ increases, the curve subsides, increasing the allowed area while decreasing the forbidden area. This effect becomes more evident from the axis of symmetry toward the device’s edges where the fringing effect is most evident.
Then, following Remark 7, Figure 2 could be an indicative criterion for choosing the material for the membrane once the intended use of the device has been chosen. In other words, once V is fixed and quantified, the effect due to the fringing field (i.e., δ ) and the allowed area on the θ λ 2 r plane is immediately visible, from which it is possible to select the most suitable range of values of θ λ 2 . Therefore, by (11), T is computable. Fixed δ and all materials that allow θ λ 2 below the curve have to be excluded.
Remark 8.
Note that the criterion above mentioned is purely theoretical as the differential analytical model proposed in this work, even if starts from physical considerations found in industrial reality (curvature of the membrane locally proportional | E | ) has not yet had an experimental confirmation (not even in terms of hardware prototyping). Therefore, in the future, it would be appropriate to proceed with software simulations of materials more adherent to industrial reality such as graphene or the highly stressed silicon nitride in order to carry out a hardware prototyping of sure industry interest.
Remark 9.
From the work in [65], we can extrapolate the following important result:
H > θ λ 2 2 R 2 ( 2 δ ) ( 1 d ) 6
which makes sense if only if δ < 2 . Obviously, for δ = 1 no problem is encountered. However, for δ = 2 , it follows that H + . However, to say that H + means to say that, given the symmetry of the membrane with respect to the origin of the Cartesian axes, at the boundary of the lower plate, the membrane adhere to the lateral surface of the device (unwanted harmful effect). Therefore, the value δ = 2 , avoided mathematically, must also be physically avoided.
Remark 10.
We note that the term due to the fringing field in the model is
λ 2 δ | u ( r ) | 2 ( 1 u ( r ) ) 2 = ξ | E f r i n g i n g f i e l d | 2 , ξ R
where
λ 2 ( 1 u ( r ) ) 2 | E | 2
and the following dimensionless term,
δ | u ( r ) | 2
weighs the fringing field phenomenon. Then, the electrostatic force due solely to the fringing field holds:
[ f e l ] f r i n g i n g f i e l d = 0.5 ϵ 0 π R 2 λ 2 δ | u ( r ) | 2 ( 1 u ( r ) ) 2
considering that | u ( r ) | 2 < H 2 = 146 2 [48], we can write:
[ f e l ] f r i n g i n g f i e l d < 296174 · 10 12 R 2 λ 2 δ d 2 .
Moreover, considering condition (11) and subbstituting the usual values for each physical parameter, we easily achieve
[ f e l ] f r i n g i n g f i e l d < 5242279 · 10 24 δ d 3 T d 2
which represents the link between the electrostatic force due to the fringing field and the δ parameter which weighs the effect due to the fringing field. Finally, we observe that the greater the mechanical tension of the membrane (term T in the denominator), the lower the effect due to the fringing field will be.

5. A Well-Known Result of Uniqueness

Theorem 2.
If all the hypotheses of Theorem 1 concerning (8) are satisfied, and u 1 ( r ) and u 2 ( r ) together satisfy the assigned boundary conditions, then the uniqueness of solution u ( r ) , such that u 1 ( r ) u ( r ) u 2 ( r ) , is not guaranteed.
Proof of Theorem 2.
See [48] □
The uniqueness result proved in Theorem 2, although theoretically interesting, has a defect in that it does not depend on θ λ 2 (i.e., the electromechanical properties of the material constituting the membrane). Therefore, in this paper, we present a new algebraic condition that ensures uniqueness for (8) involving θ λ 2 .

6. A New Condition Ensuring the Uniqueness of the Solution

The following result holds:
Proposition 1.
The algebraic condition ensuring the uniqueness of the solution for (8) is the following:
θ λ 2 > 4 R ( 1 + R ) ( 1 + δ H 2 )
where H = sup r ( 0 , R ] | u ( r ) | .
Proof. 
As the proof to obtain (54) is quite extensive, please refer to its reading in the Appendix A. □
Note that sup r ( 0 , R ] | u ( r ) | = | u ( ± R ) | . This result has been proved in [23] and confirmed numerically in [30]. On the other hand, due to reasons of symmetry, the points of Ω characterized by the greatest slope of the membrane profile can only be those located on Ω .
Remark 11.
Inequality (54) is very interesting because it expresses that the uniqueness of the solution for (8) is ensured if inf { θ λ 2 } depends on both the geometry of the device (i.e., radius R) and δ by that given algebraic combination. This makes physical sense because the greater the R, the greater the V that must be applied to overcome the inertia of the membrane (see (11)). Furthermore, the higher the δ, the higher the V must be overcome this effect and move the membrane.
(45) and (54) represent the algebraic conditions that ensure, respectively, the existence and uniqueness of the solution for (8). With δ = 0 , in [64], it was proved that the condition of uniqueness weighed more than that of existence. We wonder if this result remains even in circular 2D geometry. The following section answers this question.

7. A New Algebraic Condition Ensuring Both the Existence and Uniqueness

The following results yield the following.
Proposition 2.
Inequality (45) guarantees both the existence and uniqueness of the solution for (8).
Proof. 
As both (45) and (54) are verified, the following system makes sense:
θ λ 2 > 2 d 2 R 2 k ϵ 0 V 2 1 + δ k ϵ 0 V 2 r d 2 R 2 2 ( existence ) θ λ 2 > 4 R ( 1 + R ) ( 1 + δ H 2 ) ( uniqueness ) .
We observe that
2 d 2 R 2 k ϵ 0 V 2 1 + δ k ϵ 0 V 2 r d 2 R 2 2 > 4 R ( 1 + R ) ( 1 + δ H 2 ) .
In fact, if absurdly
2 d 2 R 2 k ϵ 0 V 2 1 + δ k ϵ 0 V 2 r d 2 R 2 2 < 4 R ( 1 + R ) ( 1 + δ H 2 ) ,
we would easily obtain
H > 1 δ d R 2 ( 1 + R ) ϵ 0 V 2 1 + δ k ϵ 0 V 2 r d 2 R 2 2 1 .
Then, substituting in (58) the plausible values for each parameter, we will get H > 10 12 , which would contradict the fact that, in circular 2D geometry with fringing field, H = 146 (see in [48]). Therefore, the existence condition also ensures the uniqueness of the solution for (8) so that each numerical solution that does not satisfy (45) represents a ghost solution. □
Remark 12.
Unlike in 1D geometry with fringing field [64], in 2D circular geometry, the condition of existence also guarantees uniqueness. On the other hand, as proved in [48], even if uniqueness are not guaranteed, equilibrium configurations are obtained dangerously close to the upper disk of the device, which is fortunately stable.
The model (8), as already specified previously, does not allow for obtaining explicit solutions. Therefore, to recover the membrane profiles, we resort to suitable numerical techniques for solving nonlinear BVPs with singularities. The following section presents a quick overview of the main characteristics of the numerical techniques used in this study.

8. Numerical Results

In this section, we present and discuss the numerical results obtained solving model (8) for recovering the membrane profile and the performances of the methods used to this end. More precisely, we have selected the following numerical approaches: the shooting method, Keller–Box scheme, and II/IV Lobatto IIIa formulas as a collocation procedure. A detailed discussion of these methods can be found in [32] with the only difference that in the present work, we used the secant method to obtain the zeros in the nonlinear equation generated by the shooting procedure. All simulations have been carried out by using the MATLAB ® R2019a environment. To facilitate the discussion, the section is divided into five subsections: the first subsection discusses the criteria for avoiding the computation of ghost solutions. The second subsection addresses the convergence of the numerical procedures as a function of the parameter θ λ 2 . In the third subsection, remarks on the choice of the number of nodes are provided. A comparison of the recovered membrane profile retrieved by the selected numerical approach is provided in the four subsection. In the last subsection, we consider the values assumed by the applied voltage V to the device to overcome the inertia of the membrane, which is useful in practical applications.

8.1. Detection of Ghost Solutions

As mentioned previously, each numerical solution must satisfy (45); otherwise, it must be considered a ghost solution. Therefore, from (45), we can express
θ λ 2 > 2 d 2 R 2 k ϵ 0 V 2 1 + δ k ϵ 0 V 2 r d 2 R 2 2 = 2 d 6 R 6 k ϵ 0 V 2 ( d 4 R 4 + δ k 2 ϵ 0 2 V 4 r 2 ) ,
but, as r < R and d 4 R 4 1 , (59) becomes
θ λ 2 > 2 d 2 R 2 k ϵ 0 V 2 1 + δ k ϵ 0 V 2 r d 2 R 2 2 = 2 d 6 R 6 k ϵ 0 V 2 ( d 4 R 4 + δ k 2 ϵ 0 2 V 4 R 2 ) .
However, d 4 R 4 1 , so (60) becomes
θ λ 2 > 2 d 6 R 6 k ϵ 0 V 2 ( 1 + δ k 2 ϵ 0 2 V 4 R 2 )
Again, combining (22), (19), and (10), and taking into account that u ( r ) < d , we easily obtain
k = 2 u ( r ) ( d u ( r ) ) 2 1 r R 2 ϵ 0 V 2 < 2 d ( d u ( r ) ) 2 1 r R 2 ϵ 0 V 2
which substituted into (60), gives us
θ λ 2 d 6 R 6 V 2 ( 1 + 4 δ d 2 ( d u ( r ) ) 4 ) .
Therefore, any numerical procedure used, δ [ 0 , 2 ) and under convergence conditions, produces a corresponding membrane profile. Then, we denote by u j , δ ( r ) the profile of the membrane obtained by applying the j-th numerical procedure with a given value of δ . Then, for each numerical procedure j with a given value of δ , the value of θ λ 2 , starting from which convergence is ensured without ghost solutions, denoted by ( ( θ λ 2 ) c o n v n o g h o s t s o l u t i o n s ) j , δ , is the following:
( ( θ λ 2 ) c o n v n o g h o s t s o l u t i o n s ) j , δ = d 6 R 6 V 2 ( 1 + 4 δ d 2 ( d max r { max j { u j , δ ( r ) } } ) 4 ) .
Therefore,
[ ( ( θ λ 2 ) c o n v n o g h o s t s o l u t i o n s ) j , δ , + )
are the ranges of values ensuring convergence without ghost solutions, j and δ .

8.2. On the Convergence of the Numerical Procedures

As shown in previous works [30,31,32,60], θ λ 2 that ensured the convergence of all the utilized numerical procedures (indicated by min ( θ λ 2 ) conv ) was obtained without considering the effects of the fringing field. Noticeably, if θ λ 2 < min ( θ λ 2 ) conv , all the numerical procedures did not converge. Conversely, if θ λ 2 > min ( θ λ 2 ) conv , all numerical procedures converged, and, in some cases, ghost solutions took place. Unlike the 1D case [65], wherein a decrease in the value of θ λ 2 was highlighted, ensuring convergence as δ increases, whatever the numerical procedure used, the circular geometry 2D studied here does not show the same behavior. Rather, there is an increasing trend of these minimum values but, in any case, of extremely reduced amplitudes. Looking at the present work in more detail, starting from (8), we achieved ( θ λ 2 ) conv exploiting, as in [65], the shooting technique by ode23 and ode45 (MatLab routines), Keller–Box scheme, and III/IV Stage Lobatto IIIa formulas by bpv4c and bpv5c (MatLab routines) when δ [ 0 , 2 ) increased, according to (46). Particularly, Table 1 and Table 2 highlight, as δ increases, the minimum value of θ λ 2 ensuring convergence for each numerical procedure such that min ( θ λ 2 ) conv = 10 3 , which corresponds to the Keller–Box scheme with δ 1 . Then, δ [ 0 , 2 ) the ranges of values relative to θ λ 2 that do not allow convergence of the numerical procedures are shown in Table 3 and Table 4, which indicates that increasing the effect due to the fringing field increases the non-convergence intervals. Moreover, for each numerical procedure and δ [ 0 , 2 ) , the convergence intervals in the presence of ghost solutions (i.e., the numerical solutions that do not respect the condition (45)) have been obtained as shown in Table 5 and Table 6. Finally, Table 7 and Table 8 show the ranges of possible values for θ λ 2 , which ensure the convergence of all numerical procedures without ghost solutions (i.e., the numerical solutions obtained satisfy the condition (45)).
Remark 13.
Notably, on increasing the effect due to the fringing field, we observe a shift toward lower values of the range of possible values for θ λ 2 that ensures convergence (with/without ghost solutions). Furthermore, the numerical procedure that determines the minimum value of θ λ 2 ensuring convergence for each numerical procedure without ghost solution, in our study, turned out to be the Keller–Box scheme (see Table 7).

8.3. A Few Remarks on the Number of Nodes N

AAs far as the Keller–Box and Lobatto procedures (bpv4c and bpv5c, respectively) are concerned, a number of nodes, N, equal to 40 have been chosen not to deviate much from the number of nodes selected by ode23 (see Table 9). We observe that the Keller–Box scheme did not converge to the solution for values of N smaller. Moreover, the performance of Lobatto formula for values of N smaller did not reach the accuracy set by default, providing unreliable solutions (for example, deformed triangle-shaped membranes with a mesh consisting of an excessive number of points). Furthermore, note that the presented results were obtained by applying numerical methods of different convergence to show the efficiency and performance of the various numerical approaches in the existing literature. From a qualitative point of view, the recovered profiles are comparable even though tiny differences appear due to the different orders of accuracy of the proposed methods.
Remark 14.
(65) is more general than the range of values shown in the previous tables because it is formulated δ [ 0 , 2 ) , while in the tables, the ranges shown refer to certain values of δ.

8.4. The Recovering of the Membrane Profile: Performance of Numerical Procedures

By applying the numerical procedures to (8) to recover the profile of the membrane, interesting results were obtained, as shown in Figure 3, Figure 4, Figure 5, Figure 6 and Figure 7. Particularly, Figure 3 depicts u ( r ) achieved by each numerical procedure when δ = 0 with its corresponding value of θ λ 2 as shown in Table 1 and Table 2. Similarly, in the other figures, for δ = 0.5 , 1 , 1.5 , 1.99 , we observe that, by ode23 with δ = 0 and θ λ 2 = 1 , a displacement of the membrane is appreciated with a maximum value of about 0.4 at r = 0 , while for θ λ 2 < 1 , the solution remains almost constant (equal to 1) as if the whole membrane is moved up by an amount equal to 1. Finally, if θ λ 2 10 15 , the solution assumes values that are too high. Situations with completely similar results are obtained with ode45. However, note that a good recovering of the membrane is achieved with bpv4c when θ λ 2 = 10 5 and δ = 0 ; with bpv4c when θ λ 2 = 10 6 and δ = 1 ; and with bpv5c when θ λ 2 = 10 7 and δ = 0 .
Remark 15.
As we will see below, the simulations performed by the numerical procedure did not show the recovering of the membrane profile in a regime of small displacements as proved in 1D geometry [65]. This allows us to highlight that the studied device could not be intended for biomedical applications (such as in drug delivery systems) wherein small displacements of the membrane are required. We also observe that the numerical methods used in this work are the most suitable for solving elliptic semi-linear differential problems with ordinary derivatives with singularities 1 r . Among them, in terms of performance, the Keller–Box scheme method stands out, which pays a higher computational cost in the face of better performance. Therefore, Keller–Box method, even if it provides the best performance, could give problems in all those real-time applications where an extremely reduced membrane recovering time is required. Fortunately, such real-time applications are not numerous, so Keller–Box method is still attractive for many engineering applications.

8.5. inf { V } to Overcome the Inertia of the Membrane with Fringing Field

As specified previously, if (45) is not satisfied, the numerical solution represents a ghost solution. Therefore, from (45), considering that r < R , we obtain
θ ϵ 0 V 2 d 3 T > d 2 k ϵ 0 V 2 1 + δ k ϵ 0 V 2 r d 2 R 2 2 > > 2 d 4 R 2 k ϵ 0 V 2 ( d 2 + δ k 2 ϵ 0 2 V 2 ) > 2 d 4 R 2 k V 2 ( d 2 + δ k 2 V 2 )
using which and exploiting (11), we achieve
V > d 3 d 4 T θ ϵ 0 k ( d 2 + δ k 2 V 2 ) 4 inf { V } .
Remark 16.
(67), without fringing field ( δ = 0 ), is the same inequality obtained in [25] necessary to overcome the inertia of the membrane without fringing field. Conversely, when the fringing field phenomenon occurs, in circular 2D geometry and 1D geometry [65], the effect because of it appears in the denominator. In both kinds of geometry, the more intense the fringing field effect, the lower the V is to move the membrane because the lines of force of E near the edges of the device are noticeably curved outward, facilitating the deformation of the membrane. Thus, the movement of the membrane requires smaller values of V. This phenomenon is of great help in all those cases wherein the membrane MEMS device is a part of an electronic device for applications that, usually, are subject to reduced values of V. In these cases, the effects of the fringing field help move the membrane more easily, overcoming its inertia.

8.6. Properties of the Material Constituting the Membrane & Intended Use of the Device in Non-Convergence Conditions

We preliminarily observe that | E | 2 λ 2 ( 1 + δ | u ( x ) | 2 ( 1 u ( x ) ) 2 such that
θ | E | 2 = λ 2 ( 1 + δ | u ( x ) | 2 ) ( 1 u ( x ) ) 2 .
Moreover, considering (11), (68) becomes
θ | E | 2 = 4 ϵ 0 2 V 4 R 4 ( 1 + δ | u ( x ) | 2 ) d 6 T 2 ( 1 u ( r ) ) 2
from which, multiplying both sides by λ 2 and considering again (11), we obtain
θ λ 2 = 4 ϵ 0 2 V 4 R 4 ( 1 + δ | u ( r ) | 2 ) d 6 T 2 ( 1 u ( r ) ) 2 | E | 2 .
Furthermore, as
| E | 2 < sup { | E | 2 }
it follows that
1 | E | 2 > 1 sup { | E | 2 }
As
1 u ( r ) < 1
it follows that
1 ( 1 u ( r ) ) 2 > 1 .
Therefore, taking into account both (72) and (74), (70) becomes
θ λ 2 > 4 ϵ 0 2 V 4 R 4 ( 1 + δ | u ( r ) | 2 ) d 6 T 2 sup { | E | 2 } .
If the non-convergence conditions of each numerical procedure occur, it follows that
θ λ 2 < min j , δ { ( ( θ λ 2 ) c o n v ) j , δ } ,
where min j , δ { ( ( θ λ 2 ) c o n v ) j , δ } represents the minimum value of θ λ 2 that ensures the convergence of all numerical procedures used for any value of δ [ 0 , 2 ) . Therefore, (75) can be written as follows:
min j , δ { ( ( θ λ 2 ) c o n v ) j , δ > θ λ 2 > 4 ϵ 0 2 V 4 R 4 ( 1 + δ | u ( r ) | 2 ) d 6 T 2 sup { | E | 2 } ,
so that
T > 4 ϵ 0 2 V 4 R 4 ( 1 + δ | u ( r ) | 2 ) d 6 min j , δ { ( ( θ λ 2 ) c o n v ) j , δ sup { | E | 2 } .
Therefore, if we choose the intended use of the device, V is fixed so that, from (67), the effect due to the fringing field is not arbitrary. Furthermore, sup { | E | 2 } inside the device also has a specific value. It follows that, in the conditions of non-convergence, all materials whose mechanical tension T satisfies the (78) must be avoided. Conversely, if the device has been chosen a priori (i.e., if the material constituting the membrane has been chosen a priori), then T is fixed, and, in the conditions of non-convergence, the intended uses of the device (i.e., the pairs { V , sup { | E | 2 } ) satisfying
V 4 sup { | E | 2 } < min j , δ { ( ( θ λ 2 ) c o n v ) j , δ } d 6 T 2 4 ϵ 0 2 R 4 ( 1 + δ | u ( r ) | 2 )
are not allowed.
Remark 17.
(78), without fringing field effects ( δ = 0 ), is the condition (73) achieved in [65]. This is because the membrane in circular geometry 2D, under symmetry conditions, can be considered as generated by the rotation of a 1D c u r v e C , lying on the vertical plane x z , and rotating around the vertical axis z, making a complete rotation. Therefore, the electromechanical behavior of the membrane is the same on any vertical plane whose support is the rotation axis z. Similarly, (79), without fringing field, becomes the condition (74) in [65]. However, we note that, in 2D circular geometry, unlike 1D geometry, the effects due to the fringing field is manifested not only by the presence of min j , δ { ( ( θ λ 2 ) c o n v ) j , δ } (as manifested in [65]), but also by the presence of δ | u ( r ) | 2 .

8.7. Properties of the Material Constituting the Membrane & Intended Use of the Device in the Presence of Ghost Solutions

Indicating by min j , δ { ( ( θ λ 2 ) c o n v ) j , δ } l i m i t the value of θ λ 2 , in convergence conditions, below which ghost solutions occur, (69) with ghost solutions satisfies the following chain of inequalities:
min j , δ { ( ( θ λ 2 ) c o n v ) j , δ } < 4 ϵ 0 R 4 V 4 ( 1 + δ | u ( r ) | 2 ) d 6 T 2 sup { | E | 2 } < min j , δ { ( ( θ λ 2 ) c o n v ) j , δ } l i m i t
so that
d 6 T 2 min j , δ { ( ( θ λ 2 ) c o n v ) j , δ } 4 ϵ 0 R 4 ( 1 + δ | u ( r ) | 2 ) < V 4 sup { | E | 2 } < d 6 T 2 min j , δ { ( ( θ λ 2 ) c o n v ) j , δ } l i m i t 4 ϵ 0 R 4 ( 1 + δ | u ( r ) | 2 ) .
Therefore, once the intended use of the device is chosen (i.e., once V satisfies (67) is fixed so that sup { | E | 2 } is also fixed ), T that characterizes the material constituting the membrane must necessarily satisfy the chain of inequalities (81). Conversely, from (80), comes the following chain of inequalities:
d 6 min j , δ { ( ( θ λ 2 ) c o n v ) j , δ } sup { | E | 2 } 4 ϵ 0 R L 4 V 4 ( 1 + δ | u ( r ) | 2 ) < 1 T 2 < d 6 min j , δ { ( ( θ λ 2 ) c o n v ) j , δ } l i m i t sup { | E | 2 } 4 ϵ 0 R 4 V 4 ( 1 + δ | u ( r ) | 2 ) .
Furthermore, in this case, once T has been set (i.e., once the membrane material has been selected), the operating parameters {V , sup { | E | 2 } } must be selected so that it is the (67) that the (82) conditions are both satisfied.
Remark 18.
The range of θ λ 2 that shows ghost solutions represents an electrostatic problem of great interest. In fact, the profiles recovered numerically are not able to satisfy the analytical model, thus requiring this range to be as small as possible. Therefore, from (81), we achieve
( min j , δ { ( ( θ λ 2 ) c o n v ) j , δ } min j , δ { ( ( θ λ 2 ) c o n v ) j , δ } l i m i t ) < < V 4 sup { | E | 2 } 1 T 2 4 ϵ 0 R 4 ( 1 + δ | u ( r ) | 2 ) d 6 .
(83) has a specific physical meaning. Once the device’s geometry is fixed (i.e., both R and d are fixed), the higher the T, the lower the θ λ 2 that risk manifesting ghost will be solutions. In other words, the stiffer the membranes, the lower the risk of obtaining ghost solutions. Again, the higher the effects due to the fringing field, the greater the risk of obtaining ghost solutions. Finally, high values of the ratio V 4 sup { | E | 2 } produces an increase in the area dedicated to ghost solutions by allocating the devices only to applications wherein reduced values of V are required to reduce the risk of ghost solutions.

8.8. Properties of the Material Constituting the Membrane and Intended Use of the Device in Absence of Ghost Solutions

If T and { V , sup { | E | 2 } } (such that V satisfies (67)) satisfy
4 ϵ 0 R 4 V 4 ( 1 + δ | u ( r ) | 2 ) d 6 T 2 sup { | E | 2 } > min j , δ { ( ( θ λ 2 ) c o n v ) j , δ } l i m i t
the device works in convergence condition without ghost solutions. Therefore, fixed { V , sup { | E | 2 } } , T must satisfy
1 T 2 > d 6 min j , δ { ( ( θ λ 2 ) c o n v ) j , δ } l i m i t sup { | E | 2 } 4 ϵ 0 R 4 V 4 ( 1 + δ | u ( r ) | 2 ) .
Conversely, once T is fixed, { V , sup { | E | 2 } } must be chosen in such a way that (85) is satisfied.
Remark 19.
(84) represents an interesting limitation for the range of values of θ λ 2 . In fact, once the geometry of the device has been fixed, the membranes that are characterized by high mechanical stresses severely limit the presence of ghost solutions, allocating the device for all those applications with low electrical potential values.

9. Conclusions

In this work, the membrane profile u ( r ) of an MEMS device with 2D circular geometry subject to an external electrical V applied between the disks has been numerically recovered. In the analytical model, we have considered that the electric field E , caused by the potential V, is locally orthogonal to the membrane’s tangent plane at the same point in a way that the principle according to which | E | results locally proportional to the membrane mean curvature, K ( r , u ( r ) ) , turns out to be valid. In addition, for taking into account the effects due to the fringing field according to Pelesko and Driscoll’s theory, an addendum term has been joined in the considered analytical model. The term is weighted by a parameter depending on the square of the first derivative membrane profile amplitude. It takes into account the deformation of the lines of force of the field E near the edges of the device. Furthermore, an algebraic condition ensuring the uniqueness of the solution, which depends on the electromechanical properties of the membrane material, has also been demonstrated. Although the founded condition results to be less incisive than the condition of existence well known in the literature on these devices, it also ensures the solution’s uniqueness. The membrane profile has been recovered through several numerical procedures; we used the shooting method, the Keller–Box method, and the Lobatto IIIa stage III/IV formulas, comparing their performances. For this particular MEMS geometry, the links among (i) the electromechanical properties of the material constituting the membrane and (ii) the operational electrical parameters in conditions of convergence and non-convergence in the presence/absence of ghost solutions have been obtained. Furthermore, the ranges of the θ λ 2 parameter, which ensures the convergence in the presence/absence of ghost solutions, have been obtained. Based on these results, a possible criterion for choosing the material constituting the membrane starting from V and vice versa in all device operating conditions has been provided. To conclude, we point out that the numerical results carried out in our study showed the clear superiority of the Keller–Box method compared with the other numerical techniques considered. However, because of its higher computational effort, the Keller–Box method could be less attractive for real-time applications.
Note that the numerical recovering in our work is a means to recognize any ghost solutions and understand the link between the properties of the material constituting the membrane and the intended use of the device in the different operating conditions. Therefore, a FEM analysis of the proposed analytical model is required (and developed in the near future) because it would allow us to obtain software formulations implementable in hardware for any prototyping of industrial use. However, in this work, we relied on the aforementioned numerical techniques because they currently represent the gold standard for this type of analytical models. Finally, it is worth nothing that concerning the material of the disks, in this work, we have considered metal plates. This arises from the fact that the analytical model here proposed derives from a well-known MEMS model with parallel metal plates proposed by D. Cassani et al. and cited in [23]. However, in the next future, it is our interest to propose analytical models closer to industrial reality, so the metal making up the disks should be replaced by more suitable materials (such as polysilicon).

Author Contributions

Conceptualization, M.V. and A.J.; methodology, A.J.; software, G.A.; validation, G.A. and F.C.M.; formal analysis, M.V.; investigation, G.A.; resources, F.C.M.; writing—original draft preparation, M.V. and G.A.; writing—review and editing, A.J. and G.A.; supervision, G.A. All authors have read and agreed to the published version of the manuscript.

Funding

This research received no external funding.

Institutional Review Board Statement

Not applicable.

Informed Consent Statement

Not applicable.

Data Availability Statement

Our study does not report any data.

Conflicts of Interest

The authors declare no conflict of interest.

Abbreviations

The following abbreviations are used in this manuscript:
Ω bounded circular smooth domain
ddistance between the two parallel disks
rradial coordinate
Rradius of the device
u ( r ) profile of the membrane
Vexternal electrical voltage
Tradial mechanical tension of the membrane at rest
λ 2 parameter depending on V and T
δ parameter that weighs the fringing field effect
E electrostatic field
| E | amplitude of the electrostatic field
K ( r , u ( r ) ) mean curvature
d critical security distance
θ factor of proportionality
ϵ 0 permittivity of the free space
f e l electrostatic force
p e l electrostatic pressure
pmechanical pressure
u 0 displacement in the center of the membrane
C e l electrostatic capacitance
ρ density
hthickness
YYoung modulus
ν Poisson ratio
Dstiffness coefficient
μ ( r , u ( r ) , λ ) function of proportionality
u 1 ( r ) , u 2 ( r ) auxiliary functions
kconstant of proportionality between p and p e l
H sup r ( 0 , 1 ] | u ( r ) |
FEMFinite Element Method
t o l tolerance for Brent procedure
T O L tolerance for Keller–Box Scheme
R K Runge–Kutta Methods
G ( r , s ) Green function

Appendix A. Proof of Proposition 1

We preliminarly observe that model (8), exploiting a suitable Green’s function G ( r , s ) , can be rewritten in its equivalent integral formulation [66,67,68]:
u ( r ) = R R G ( r , s ) F ( r , u ( r ) , u ( r ) ) d s
and considering (37), (A1) becomes:
u ( r ) = R R G ( r , s ) u ( r ) r + 4 ( 1 u ( r ) d ) 2 θ λ 2 ( 1 + δ | u ( r ) | 2 ) d s .
As known, G ( r , s ) , in our problem, assumes the following form [66]:
G ( r , s ) = ( s + R ) ( R r ) 2 R
if R s r and
G ( r , s ) = ( R s ) ( r + R ) 2 R
when r < s R . In addition:
G r ( r , s ) = ( s + R ) 2 R
if R s < r and
G r ( r , s ) = ( R s ) 2 R
when r < s R . Furthermore, it is easy to prove that [66]:
0 G ( r , s ) R 2 r , s R , R ]
and
G r ( r , s ) 1 2 .
By contradiction, we assume that u 1 ( r ) and u 2 ( r ) are tho different solutions for (8), in order that u 1 ( r ) = T ( u 1 ( r ) ) and u 2 ( r ) = T ( u 2 ( r ) ) . Thus,
u 1 ( r ) = T ( u 1 ( r ) ) = R R G ( r , s ) u 1 ( r ) r + 4 ( 1 u 1 ( r ) d ) 2 θ λ 2 ( 1 + δ | u 1 ( r ) | 2 ) d s
and
u 2 ( r ) = T ( u 2 ( r ) ) = R R G ( r , s ) u 2 ( r ) r + 4 ( 1 u 2 ( r ) d ) 2 θ λ 2 ( 1 + δ | u 2 ( r ) | 2 ) d s ,
from which
u 1 ( r ) = T ( u 1 ( r ) ) = R R G r ( r , s ) u 1 ( r ) r + 4 ( 1 u 1 ( r ) d ) 2 θ λ 2 ( 1 + δ | u 1 ( r ) | 2 ) d s
and
u 2 ( r ) = T ( u 2 ( r ) ) = R R G r ( r , s ) u 2 ( r ) r + 4 ( 1 u 2 ( r ) d ) 2 θ λ 2 ( 1 + δ | u 2 ( r ) | 2 ) d s .
Therefore, taking into account also both (A7) and (A8), we can write:
| | u 1 ( r ) u 2 ( r ) | | C 1 ( [ R , R ] ) = = sup r [ R , R ] | u 1 ( r ) u 2 ( r ) | + sup r [ R , R ] | u 1 ( r ) u 2 ( r ) | = = sup r [ R , R ] | R R G ( r , s ) u 1 ( r ) r + 4 ( 1 u 1 ( r ) d ) 2 θ λ 2 ( 1 + δ | u 1 ( r ) | 2 ) d s R R G ( r , s ) u 2 ( r ) r + 4 ( 1 u 2 ( r ) d ) 2 θ λ 2 ( 1 + δ | u 2 ( r ) | 2 ) d s | + + sup r [ R , R ] | R R G r ( r , s ) u 1 ( r ) r + 4 ( 1 u 1 ( r ) d ) 2 θ λ 2 ( 1 + δ | u 1 ( r ) | 2 ) d s R R G r ( r , s ) u 2 ( r ) r + 4 ( 1 u 2 ( r ) d ) 2 θ λ 2 ( 1 + δ | u 2 ( r ) | 2 ) d s | = = sup r [ R , R ] | R R G ( r , s ) { u 1 ( r ) r + 4 ( 1 u 1 ( r ) d ) 2 θ λ 2 ( 1 + δ | u 1 ( r ) | 2 ) u 2 ( r ) r + 4 ( 1 u 2 ( r ) d ) 2 θ λ 2 ( 1 + δ | u 2 ( r ) | 2 ) } d s | + + sup r [ R , R ] | R R G r ( r , s ) { u 1 ( r ) r + 4 ( 1 u 1 ( r ) d ) 2 θ λ 2 ( 1 + δ | u 1 ( r ) | 2 ) u 2 ( r ) r + 4 ( 1 u 2 ( r ) d ) 2 θ λ 2 ( 1 + δ | u 2 ( r ) | 2 ) } d s | sup r [ R , R ] R R | G ( r , s ) | | u 1 ( r ) r + 4 ( 1 u 1 ( r ) d ) 2 θ λ 2 ( 1 + δ | u 1 ( r ) | 2 ) u 2 ( r ) r + 4 ( 1 u 2 ( r ) d ) 2 θ λ 2 ( 1 + δ | u 2 ( r ) | 2 ) | d s + + sup r [ R , R ] R R | G r ( r , s ) | | u 1 ( r ) r + 4 ( 1 u 1 ( r ) d ) 2 θ λ 2 ( 1 + δ | u 1 ( r ) | 2 ) u 2 ( r ) r + 4 ( 1 u 2 ( r ) d ) 2 θ λ 2 ( 1 + δ | u 2 ( r ) | 2 ) | d s R 2 sup r [ R , R ] R R | u 1 ( r ) r u 2 ( r ) r + + 4 θ λ 2 ( 1 u 1 ( r ) d ) 2 ( 1 + δ | u 1 ( r ) | 2 ) ( 1 u 2 ( r ) d ) 2 ( 1 + δ | u 2 ( r ) | 2 ) | d s + + 1 2 sup r [ R , R ] R R | u 1 ( r ) r u 2 ( r ) r + + 4 θ λ 2 ( 1 u 1 ( r ) d ) 2 ( 1 + δ | u 1 ( r ) | 2 ) ( 1 u 2 ( r ) d ) 2 ( 1 + δ | u 2 ( r ) | 2 ) | d s = = R 2 + 1 2 sup r [ R , R ] R R | u 1 ( r ) r u 2 ( r ) r + + 4 θ λ 2 ( 1 u 1 ( r ) d ) 2 ( 1 + δ | u 1 ( r ) | 2 ) ( 1 u 2 ( r ) d ) 2 ( 1 + δ | u 2 ( r ) | 2 ) | d s .
We note that
1 1 + δ | u ( r ) | 2 = 1 + δ | u ( r ) | 2 δ | u ( r ) | 2 1 + δ | u ( r ) | 2 1 + δ | u ( r ) | 2 δ | u ( r ) | 2 1 + δ | u ( r ) | 2 ,
therefore, (A13), taking into account (A14), becomes
| | u 1 ( r ) u 2 ( r ) | | C 1 ( [ R , R ] ) R 2 + 1 2 sup r [ R , R ] { R R 1 r | u 1 ( r ) u 2 ( r ) | d s + + 4 θ λ 2 R R ( 1 u 1 ( r ) d ) 2 ( 1 + δ H 2 ) 2 d s + 4 θ λ 2 R R ( 1 u 2 ( r ) d ) 2 ( 1 + δ H 2 ) d s } d s R 2 + 1 2 sup r [ R , R ] { R R 1 r | u 1 ( r ) u 2 ( r ) | d s + + 4 θ λ 2 ( 1 + δ H 2 ) 2 R R [ ( 1 u 1 ( r ) ) 2 + ( 1 u 2 ( r ) ) 2 ] d s } = = R 2 + 1 2 sup r [ R , R ] { 2 R r | u 1 ( r ) u 2 ( r ) | + + 8 R θ λ 2 ( 1 + δ H 2 ) [ ( 1 u 1 ( r ) ) 2 + ( 1 u 2 ( r ) ) 2 ] d s } = = R 2 + 1 2 2 R r sup r [ R , R ] | u 1 ( r ) u 2 ( r ) | + + R 2 + 1 2 8 R θ λ 2 ( 1 + δ H 2 ) sup r [ R , R ] [ ( 1 u 1 ( r ) ) 2 + ( 1 u 2 ( r ) ) 2 ] 2 ( R + 1 ) R sup r [ R , R ] | u 2 ( r ) u 1 ( r ) | + + 4 R ( R + 1 ) ( 1 + δ H 2 ) θ λ 2 sup r [ R , R ] | u 2 ( r ) u 1 ( r ) | .
From (A15), to achieve a contradiction, it is necessary to write
2 ( R + 1 ) R < 1 4 R ( R + 1 ) ( 1 + δ H 2 ) θ λ 2 < 1 .
From (A16) we observe that
4 R ( R + 1 ) ( 1 + δ H 2 ) θ λ 2 > 2 ( R + 1 ) R .
In fact, if absurdly
4 R ( R + 1 ) ( 1 + δ H 2 ) θ λ 2 < 2 ( R + 1 ) R
we easily would obtain
H < θ λ 2 2 1 1 δ
which results an absurd condition because in absence of fringing field (i.e., δ = 0 ), it follow that H + . In other words, without fringing field the membrane profiles numerically recovered are symmetric with respect the axis r = 0 (see [32]) putting in evidence that max { | u ( r ) | } corresponds to r = ± R . Moreover [48],
max { | u ( r ) | } < sup { | u ( r ) | } = H = 146 .
Therefore, on x = ± R , there is no risk of dangerous adhesion between the membrane and the vertical walls of the device that could generate unwanted electrostatic effects. This phenomenon could occur if | u ( ± R ) | + . From which, it follows that the algebraic condition of uniqueness for the solution for (8) is
θ λ 2 > 4 ( R + 1 ) R ( 1 + δ H 2 )
(i.e., the algebraic inequality (54)).

References

  1. Pelesko, J.A.; Bernstein, D.H. Modeling MEMS and NEMS; Chapman & Hall/CRC Press: Boca Raton, FL, USA; London, UK; New York, NY, USA; Washington, DC, USA, 2003. [Google Scholar]
  2. Di Barba, P.; Wiak, S. MEMS: Field Models and Optimal Design; Springer International Publishing: Berlin/Heidelberg, Germany, 2020. [Google Scholar]
  3. Pelesko, J.A. Electrostatic in MEMS and NEMS; Micromechanics and Nanoscale Effects; Springer: Dordrecht, The Netherlands, 2004. [Google Scholar]
  4. Adrian, R.J.; Agarwal, R.K.; Beheim, G.M.; Bergstrom, P.L.; Bernestein, G.H.; Beskok, A.; Bewley, T.R.; Breuer, K.S.; Chang, H.C.; Chen, L.-Y.; et al. The MEMS Handbook; Gad-el-Hak, M., Ed.; CRC Press: Boca Raton, FL, USA; London, UK; New York, NY, USA; Washington, DC, USA, 2015. [Google Scholar]
  5. Lindros, V.; Tilli, M.; Lohto, A.; Motooka, T. Handbook of Silicon Based MEMS Materials and Technologies; Tilli, M., Paulasto-Kröckel, M., Petzold, M., Theuss, H., Motooka, T., Lindroos, V., Eds.; Elsevier: Amsterdam, The Netherlands, 2020. [Google Scholar]
  6. Versaci, M.; Morabito, F.C. Fuzzy Time Series Approach for Disruption Prediction in Tokamak Reactors. IEEE Trans. Magn. 2003, 39, 1503–1506. [Google Scholar] [CrossRef]
  7. Peushuai, S.; Chaowei, S.; Zhang, M.; Zhao, Y.; He, Y.; Liu, W.; Wang, X. A Movel Piezoresistive MEMS Pressure Sensors Based on Temporary Bonding Technology. Sensors 2020, 20, 337. [Google Scholar]
  8. Setiono, A.; Bertke, M.; Nyang’au, W.O.; Xu, J.; Fahrbach, M.; Kirsch, I.; Uhde, E.; Deutschinger, A.; Fantner, E.J.; Peiner, E.; et al. In-Plane and Out-of-Plane MEMS Piezoresistive Cantilever Sensors for Nanoparticle Mass Detection. Sensors 2020, 20, 618. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  9. Wang, Y.; Fu, Q.; Zhang, Y.; Zhang, W.; Chen, D.; Yin, L.; Liu, X. A Digital Closed-Loop Sense MEMS Disk Resonator Gyroscope Circuit Design Based on Integrated Analog Front-End. Sensors 2020, 20, 687. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  10. Younis, M.I.; Nayfeh, A.H. A Study of the Nonlinear Response of a Resonant Microbeam to an Electric Actuation. Nonlinear Dyn. 2003, 31, 91–117. [Google Scholar] [CrossRef]
  11. Wang, X.; Huan, R.; Zhu, W.; Pu, D.; Wei, X. Frequency Locking in the Internal Resonance of Two Electrostatically Coupled Micro-Resonators with Frequency ratio 1:3. Mech. Syst. Signal Process. 2021, 146, 106981. [Google Scholar] [CrossRef]
  12. Abdolvand, R.; Bahreyni, B.; Lee, J.E.Y.; Nabki, F. Micromachined Resonators, A Review. Micromachines 2016, 7, 160. [Google Scholar] [CrossRef]
  13. Shao, X.; Si, H.; Zhang, W. Fuzzy Wavelet Neural Control with Improved Prescribed Performance for MEMS Byroscope Subject to Input Quantization. Fuzzy Sets Syst. 2021, 411, 136–154. [Google Scholar] [CrossRef]
  14. Ali, I.A. Modeling and Simulation of MEMS Components: Challenges and Possible Solutions, Micromachining Techniques for Fabrication of Micro and Nano Structures; Springer Nature: Singapore, 2012. [Google Scholar]
  15. Bechtold, T.; Schrag, G.; Fent, L. System-Level Modeling of MEMS; Wyley-VCH Verlag GmbH & Co. KGaA: London, UK, 2013. [Google Scholar]
  16. Zhu, J. Development Trends and Perspectives of Future Sensors and MEMS/NEMS. Micromachines 2020, 11, 7. [Google Scholar] [CrossRef] [Green Version]
  17. Quakad, H. Electrostatic Fringing-Fields Effects on the Structural Behavior of MEMS Shallow Arches. Microsyst. Technol. 2018, 24, 1394–1399. [Google Scholar]
  18. de Oliveira Hanse, R.; Mátéfi-Tempfli, M.; Safonovs, R.; Adam, J.; Chemnitz, S.; Reimer, T.; Wagner, B.; Benecke, W.; Mátéfi-Tempfli, S. Magnetic Films for Electromagnetic Actuation in MEMS Switches. Microsyst. Technol. 2018, 24, 1394–1399. [Google Scholar] [CrossRef]
  19. Mohammadi, A.; Ali, N. Effect of High Electrostatic Actuation on Thermoelastic Damping in Thin Rectangular Microplate Resonators. J. Theor. Appl. Mech. 2015, 53, 317–329. [Google Scholar] [CrossRef] [Green Version]
  20. Farhangian, F.; Landry, R., Jr. Accuracy Improvement of Attitude Determination Systems Using EKF-Based Error Prediction Filter and PI Controller. Sensors 2020, 20, 4055. [Google Scholar] [CrossRef]
  21. Zhou, N.; Jia, P.; Liu, J.; Ren, Q.; An, G.; Liang, T.; Xiong, J. MEMS-Based Reflective Intensity-Modulated Fiber-Optic Sensor for Pressure Measurements. Sensors 2020, 20, 2233. [Google Scholar] [CrossRef]
  22. Nathanson, H.; Newell, W.; Wickstrom, R.; Lewis, J. The Resonant Gate Transistor. IEEE Trans. Electron Device 1964, 14, 117–133. [Google Scholar] [CrossRef]
  23. Di Barba, P.; Fattorusso, L.; Versaci, M. Electrostatic Field in Terms of Geometric Curvature in Membrane MEMS Devices. Commun. Appl. Ind. Math. 2017, 8, 165–184. [Google Scholar] [CrossRef] [Green Version]
  24. Di Barba, P.; Fattorusso, L.; Versaci, M. A 2D Non-Linear Second-Order Differential Model for Electrostatic Circular Membrane MEMS Devices: A Result of Existence and Uniqueness. Mathematics 2019, 7, 1193. [Google Scholar] [CrossRef] [Green Version]
  25. Versaci, M.; Morabito, F.C. Electrostatic Circular Membrane MEMS: An Approach to the Optimal Control. Computation 2021, 7, 41. [Google Scholar] [CrossRef]
  26. Versaci, M.; Morabito, F.C. Membrane Micro Electro-Mechanical Systems for Industrial Applications; Handbook of Research on Advanced Mechatronic Systems and Intelligent Robotics; IGI Global: Hershey, PA, USA, 2019. [Google Scholar]
  27. Esposito, P.; Ghoussoub, N.; Guo, Y. Mathematical Analysis of Partial Differential Equations Modeling Electrostatic MEMS; American Mathematical Society: New York, NY, USA, 2010. [Google Scholar]
  28. Liu, D.; Liu, H.; Liu, J.; Hu, F.; Fan, J.; Wu, W.; Tu, L. Temperature Gradient Method for Alleviating Bonding-Induced Warpage in a High-Precision Capacitive MEMS Accelerometer. Sensors 2020, 20, 1186. [Google Scholar] [CrossRef] [Green Version]
  29. Zhang, Y.; Li, B.; Li, H.; Shen, S.; Li, F.; Ni, W.; Cao, W. Investigation of Potting-Adhesive-Induced Thermal Stress in MEMS Pressure Sensors. Sensors 2021, 21, 2011. [Google Scholar] [CrossRef]
  30. Versaci, M.; Di Barba, P.; Morabito, F.C. Curvature-Dependent Electrostatic Field as a Principle for Modelling Membrane-Based MEMS Devices. A Review. Membranes 2020, 10, 361. [Google Scholar] [CrossRef] [PubMed]
  31. Angiulli, G.; Jannelli, A.; Morabito, F.C.; Versaci, M. Reconstructing the Membrane Detection of a 1D Electrostatic-Driven MEMS Device by the Shooting Method: Convergence Analysis and Ghost Solutions Identification. Comput. Appl. Math. 2018, 37, 4484–4498. [Google Scholar] [CrossRef]
  32. Versaci, M.; Jannelli, A.; Angiulli, G. Electrostatic Micro-electro-Mechanical-Systems (MEMS) Devices: A Comparison Among Numerical Techniques for Recovering the Membrane Profile. IEEE Access 2020, 8, 125874–125886. [Google Scholar] [CrossRef]
  33. Howell, L.L.; Luon, S.M. Thermomechanical in-Plane Microactuator, (TIM). U.S. Patent N. US6734597B1, 11 May 2004. [Google Scholar]
  34. Farokhi, H.; Ghayesh, M.H. Nonlinear Thermo-Mechanical Behaviour of MEMS Resonators. Membranes 2017, 23, 5303–5315. [Google Scholar] [CrossRef]
  35. Ren, Z.; Yuan, J.; Su, X.; Mangla, S.; Nam, C.Y.; Lu, M.; Camino, F.; Shi, Y. Thermo-Mechanical Modeling and Experimental Validation for Multilayered Metallic Microstructures. Microsyst. Technol. 2020, 21, 751–783. [Google Scholar] [CrossRef]
  36. Mistry, K.K.; Mahapatra, A. Design and Simulation of a Thermo-Transfer Type MEMS Based Micro Flow Sensor for Arteria Blood Flow Measurement. Microsyst. Technol. 2012, 18, 683–692. [Google Scholar] [CrossRef]
  37. Scaccabarozzi, D.; Saggin, B.; Somaschini, R.; Magni, M.; Valnegri, P.; Esposito, F.; Molfese, C.; Cozzolino, F.; Mongelluzzo, G. “MicroMED” Optical Particle Counter: From Design to Flight Model. Sensors 2020, 20, 611. [Google Scholar] [CrossRef] [Green Version]
  38. Di Barba, P.; Mognaschi, M.E.; Sieni, E. Many Objective Optimization of a Magnetic Micro-Electro-Mechanical (MEMS) Micromirror with Bounded MP-NSGA Algorithm. Mathematics 2020, 8, 1509. [Google Scholar] [CrossRef]
  39. Akoz, G.; Gokdel, Y.D. Field-of-View Optimization of Magnetically Actuated 2D Gimballed Scanners. Turk. J. Electr. Eng. Comput. Sci. 2020, 28, 2385–2399. [Google Scholar] [CrossRef]
  40. Patterson, L.H.; Walker, J.L.; Naivar, M.A.; Rodriguez-Mesa, E.; Hoonejani, M.R.; Shields, K.; Foster, J.S.; Doyle, A.M.; Valentine, M.T.; Foster, K.L. Inertial Flow Focusing: A Case Studi in Optimixing Cellular Trajectory Through a Microfluidic MEMS Device for Timing-Critical Application. Biomed. Microdevices 2020, 22, 1–12. [Google Scholar] [CrossRef]
  41. Can Atik, A.; Ozkan, M.D.; Ozgur, E.M.; Kulah, H.; Yildirim, E. Modeling and Fabrication of Electrostatically Actuated Diaphragms for on-Chip Valving of MEMS Compatible Microfluidic Systems. J. Micromech. Microeng. 2020, 30, 115001. [Google Scholar] [CrossRef]
  42. Xu, Y.; Hu, X.; Kundu, S.; Nag, A.; Afsarimanesh, N.; Sapra, S.; Mukhopadhyay, S.C.; Han, T. Silicon-Based Sensors for Biomedical Applications: A Review. Sensors 2019, 19, 2908. [Google Scholar] [CrossRef] [Green Version]
  43. Lee, S.; Chen, C.; Deshp, V.V.; Lee, G.H.; Lee, I.; Lekas, M.; Gondarenko, A.; Yu, Y.-J.; Shepard, K.; Hone, J.; et al. Electrically Integrated SU-8 Clamped Graphene Drum Resonators for Strain Engineering. Appl. Phys. Lett. 2013, 102, 153101. [Google Scholar] [CrossRef] [Green Version]
  44. Fan, X.; Smith, A.D.; Forsberg, F.; Wagner, S.; Schröder, S.; Akbari, S.S.A.; Fischer, A.C.; Villanueva, L.G.; Östling, M.; Niklaus, F.; et al. Manufaccture and Characterization of Graphene Membranes with Suspended Silicon Proof Masses for MEMS and NEMS Applications. Microsyst. Nanoeng. 2020, 102, 1–27. [Google Scholar]
  45. Zhou, X.; Venkatachalam, S.; Zhou, R.; Xu, H.; Pokharel, A.; Fefferman, A.; Zaknoune, M.; Collin, E. High-Q Silicon Nitride Drum Resonators Strongly Coupled. Nano Lett. 2021, 21, 5738–5744. [Google Scholar] [CrossRef]
  46. Phan, A.; Truong, P.; Schade, C.; Joslin, K.; Talke, F. Analytical Modeling of an Implantable Opto-Mechanical Pressure Sensor to Study Long Term Drift. In Proceedings of the ASME 2020 29th Conference on Information Storage and Processing Systems, Virtual, Online, 2–3 June 2021. V001T06A003, ASME. [Google Scholar]
  47. Yen, Y.K.; Chiu, C.Y. A CNOS- MEMS-Based Membrane-Bridge Nanomechanical Sensors for Small Molecule Detection. Sci. Rep. 2020, 10, 2931. [Google Scholar] [CrossRef]
  48. Di Barba, P.; Fattorusso, L.; Versaci, M. A 2D Membrane MEMS Device Model with Fringing Field: Curvature-Dependent Electrostatic Field and Optimal Control. Mathematics 2021, 6, 465. [Google Scholar] [CrossRef]
  49. Lin, Z.; Li, X.; Jin, Z.; Qian, J. Fluid-Structure Interaction Analysis on Membrane Behavior of a Microfluidic Passive Valve. Membranes 2020, 10, 300. [Google Scholar] [CrossRef]
  50. Pal, M.; Lalengkima, C.; Maity, R.; Baishya, S.; Maity, N.P. Effects of Fringing Capacitances and Electrode’s Finiteness in Improved SiC Membrane Based Micromachined Ultrasonic Transducers. Microsyst. Technol. 2021. [Google Scholar] [CrossRef]
  51. Kottapalli, A.G.; Asadnia, M.; Miao, J.M.; Barbastathis, G.; Triantafyllou, M.S. A Flexible Liquid Crystal Polymer MEMS Pressure Sensor Array for Fish-Like Underwater Sensing. Smart Mater. Struct. 2012, 21, 281–299. [Google Scholar] [CrossRef] [Green Version]
  52. Versaci, M.; Mammone, N.; Ieracitano, C.; Morabito, F.C. Micropumps for Drug Delivery Systems: A New Semi-Linear Elliptic Boundary-Value Problem. Comput. Appl. Math. 2021, 40, 1–21. [Google Scholar] [CrossRef]
  53. Cassani, D.; d’O, M.; Ghoussoub, N. On a Fourth Order Elliptic Problem with a Singular Nonlinearity. Nonlinear Stud. 2009, 9, 189–209. [Google Scholar] [CrossRef]
  54. Leus, V.; Elata, D. Fringing Field Effect in Electrostatic Actuator; Technical Report ETR-2004-2; Faculty of mechanical engineering Technion—Israel Institute of Technology: Haifa, Israel, 2004. [Google Scholar]
  55. Weng, C.C.; Kong, J.A. Effects of Fringing Fields on the Capacitance of Circular Microstrip Disk. IEEE Trans. Microw. Theory Tech. 1980, 28, 98–104. [Google Scholar] [CrossRef] [Green Version]
  56. Zhang, Y.; Wang, T.; Luo, A.; Hu, Y.; Li, X.; Wang, F. Micro Electrostatic Energy Harvester with Both Broad Bandwidth and High Normalized Power Density. Appl. Energy 2018, 212, 363–371. [Google Scholar] [CrossRef]
  57. Batra, R.C.; Porfiri, M.; Spinello, D. Electromechanical Model of Electrically Actuated Narrow Microbeams. J. Microelectromech. Syst. 2006, 15, 1175–1189. [Google Scholar] [CrossRef]
  58. Batra, R.C.; Porfiri, M.; Spinello, D. Review of Modeling Electrostatically Actuated Microelectromechanical Systems. Smart Mater. Struct. 2007, 16, 1119–1131. [Google Scholar] [CrossRef]
  59. Pelesko, J.A.; Driscoll, T.A. The Effect of the Small-Aspect-Ratio Approximation on Canonical Electrostatic MEMS Models. J. Eng. Math. 2005, 53, 239–252. [Google Scholar] [CrossRef] [Green Version]
  60. Versaci, M.; Angiulli, G.; Jannelli, A. Recovering of the Membrane Profile of an Electrostatic Circular MEMS by the Three-Stage Lobatto Procedure: A Convergence Analysis in the Absence of Ghost Solutions. Mathematics 2020, 8, 487. [Google Scholar] [CrossRef] [Green Version]
  61. Jonassen, N. Electrostatics; Springer Science+Business Media: New York, NY, USA, 2002. [Google Scholar]
  62. Timoshenko, S.; Woinowsly-Krieger, S. Theory of Plates and Shells; McGraw Hill: New York, NY, USA, 1959. [Google Scholar]
  63. Fujimoto, M. Physics of Classical Electromagnetism; Springer: New York, NY, USA, 2007. [Google Scholar]
  64. Di Barba, P.; Fattorusso, L.; Versaci, M. Curvature-Dependent Electrostatic Field as a Principle for Modelling Membrane MEMS Device with Fringing Field. Comp. Appl. Math. 2021, 40, 1–30. [Google Scholar] [CrossRef]
  65. Versaci, M.; Di Barba, P.; Morabito, F.C. MEMS with Fringing Field: Curvature-Dependent Electrostatic Field and Numerical Techniques for Recovering the Membrane Profile. Comput. Appl. Math. 2021, 40, 1–28. [Google Scholar] [CrossRef]
  66. Bayley, P.B.; Shampine, L.F.; Waltman, P.E. Nonlinear Two Points Boundary Value Problems; Academic Press: Cambridge, MA, USA, 1969. [Google Scholar]
  67. Russell, R.D.; Shampine, L.F. Numerical Methods for Singular Boundary Value Problems. SIAM J. Num. Anal. 1975, 12, 13–36. [Google Scholar] [CrossRef]
  68. Quarteroni, A.; Sacco, R.; Saleri, F. Numerical Mathematics; Springer: Berlin/Heidelberg, Germany, 2007. [Google Scholar]
Figure 1. A 2D electrostatic circular membrane MEMS device whose metal plates (upper and support ones) are displayed in gray. Between them, a circular membrane, clumped to the edges of the support plate, deforms towards the upper plate without touching it to avoid unwanted electric discharges.
Figure 1. A 2D electrostatic circular membrane MEMS device whose metal plates (upper and support ones) are displayed in gray. Between them, a circular membrane, clumped to the edges of the support plate, deforms towards the upper plate without touching it to avoid unwanted electric discharges.
Sensors 21 05237 g001
Figure 2. A graphical representation of (8) when δ changes; the forbidden area is located below each curve, while the permitted area is highlighted above each curve.
Figure 2. A graphical representation of (8) when δ changes; the forbidden area is located below each curve, while the permitted area is highlighted above each curve.
Sensors 21 05237 g002
Figure 3. Recovering of u ( r ) for δ = 0 and θ λ 2 as reported in Table 1 and Table 2.
Figure 3. Recovering of u ( r ) for δ = 0 and θ λ 2 as reported in Table 1 and Table 2.
Sensors 21 05237 g003
Figure 4. Recovering of u ( r ) for δ = 0.5 and θ λ 2 as reported in Table 1 and Table 2.
Figure 4. Recovering of u ( r ) for δ = 0.5 and θ λ 2 as reported in Table 1 and Table 2.
Sensors 21 05237 g004
Figure 5. Recovering of u ( r ) for δ = 1 and θ λ 2 as reported in Table 1 and Table 2.
Figure 5. Recovering of u ( r ) for δ = 1 and θ λ 2 as reported in Table 1 and Table 2.
Sensors 21 05237 g005
Figure 6. Recovering of u ( r ) for δ = 1.5 and θ λ 2 as reported in Table 1 and Table 2.
Figure 6. Recovering of u ( r ) for δ = 1.5 and θ λ 2 as reported in Table 1 and Table 2.
Sensors 21 05237 g006
Figure 7. Recovering of u ( r ) for δ = 1.99 and θ λ 2 as reported in Table 1 and Table 2.
Figure 7. Recovering of u ( r ) for δ = 1.99 and θ λ 2 as reported in Table 1 and Table 2.
Sensors 21 05237 g007
Table 1. Ranges of θ λ 2 ensuring convergence (shooting and Keller–Box procedures).
Table 1. Ranges of θ λ 2 ensuring convergence (shooting and Keller–Box procedures).
δ Shooting (ode 23)Shooting (ode 45)Keller–Box
0 ( θ λ 2 ) conv [ 10 14 + ) ( θ λ 2 ) conv [ 10 14 + ) ( θ λ 2 ) conv [ 10 3 + )
0.5 ( θ λ 2 ) conv [ 10 10 + ) ( θ λ 2 ) conv [ 10 10 + ) ( θ λ 2 ) conv [ 10 5 + )
1 ( θ λ 2 ) conv [ 10 10 + ) ( θ λ 2 ) conv [ 10 9 + ) ( θ λ 2 ) conv [ 10 3 + )
1.5 ( θ λ 2 ) conv [ 10 9 + ) ( θ λ 2 ) conv [ 10 9 + ) ( θ λ 2 ) conv [ 10 3 + )
1.99 ( θ λ 2 ) conv [ 10 8 + ) ( θ λ 2 ) conv [ 10 7 + ) ( θ λ 2 ) conv [ 10 3 + )
Table 2. Ranges of θ λ 2 ensuring convergence (Three-/Four-Stage Lobatto IIIa).
Table 2. Ranges of θ λ 2 ensuring convergence (Three-/Four-Stage Lobatto IIIa).
δ Three-Stage Lobatto IIIa (bpv4c)Four-Stage Lobatto IIIa (bpv5c)
0 ( θ λ 2 ) conv [ 10 4 + ) ( θ λ 2 ) conv [ 10 6 + )
0.50 ( θ λ 2 ) conv [ 10 5 + ) ( θ λ 2 ) conv [ 10 5 + )
1 ( θ λ 2 ) conv [ 10 5 + ) ( θ λ 2 ) conv [ 10 5 + )
1.50 ( θ λ 2 ) conv [ 10 5 + ) ( θ λ 2 ) conv [ 10 5 + )
1.99 ( θ λ 2 ) conv [ 10 5 + ) ( θ λ 2 ) conv [ 10 5 + )
Table 3. Ranges of θ λ 2 in conditions of non-convergence (shooting and Keller–Box procedures).
Table 3. Ranges of θ λ 2 in conditions of non-convergence (shooting and Keller–Box procedures).
δ Shooting (ode 23)Shooting (ode 45)Keller–Box
0 [ 0 , 10 14 ) [ 0 , 10 14 ) [ 0 , 10 3 )
0.5 [ 0 , 10 10 ) [ 0 , 10 10 ) [ 0 , 10 5 )
1 [ 0 , 10 10 ) [ 0 , 10 9 ) [ 0 , 10 3 )
1.5 [ 0 , 10 9 ) [ 0 , 10 9 ) [ 0 , 10 3 )
1.99 [ 0 , 10 8 ) [ 0 , 10 7 ) [ 0 , 10 3 )
Table 4. Ranges of θ λ 2 in conditions of non-convergence (Three-/Four-Stage Lobatto IIIa).
Table 4. Ranges of θ λ 2 in conditions of non-convergence (Three-/Four-Stage Lobatto IIIa).
δ Three-Stage Lobatto IIIa (bpv4c)Four-Stage Lobatto IIIa (bpv5c)
0 [ 0 , 10 4 ) [ 0 , 10 6 )
0.50 [ 0 , 10 5 ) [ 0 , 10 5 )
1 [ 0 , 10 5 ) [ 0 , 10 5 )
1.50 [ 0 , 10 5 ) [ 0 , 10 5 )
1.99 [ 0 , 10 5 ) [ 0 , 10 5 )
Table 5. Ranges of θ λ 2 ensuring convergence with ghost solutions (shooting and Keller–Box procedures).
Table 5. Ranges of θ λ 2 ensuring convergence with ghost solutions (shooting and Keller–Box procedures).
δ Shooting (ode 23)Shooting (ode 45)Keller–Box
0 [ 10 14 , 0.639 ) [ 10 14 , 0.633 ) [ 10 3 , 0 . 721 )
0.5 [ 10 10 , 0.627 ) [ 10 10 , 0.625 ) [ 10 5 , 0 . 716 )
1 [ 10 10 , 0.614 ) [ 10 9 , 0.618 ) [ 10 3 , 0 . 709 )
1.5 [ 10 9 , 0.611 ) [ 10 9 , 0.612 ) [ 10 3 , 0 . 703 )
1.99 [ 10 8 , 0.599 ) [ 10 7 , 0.603 ) [ 10 3 , 0 . 694 )
Table 6. Ranges of θ λ 2 ensuring convergence with ghost solutions (Three-/Four-Stage Lobatto IIIa).
Table 6. Ranges of θ λ 2 ensuring convergence with ghost solutions (Three-/Four-Stage Lobatto IIIa).
δ Three-Stage Lobatto IIIa (bpv4c)Four-Stage Lobatto IIIa (bpv5c)
0 [ 10 4 , 0.693 ) [ 10 6 , 0.698 )
0.50 [ 10 5 , 0.686 ) [ 10 5 , 0.691 )
1 [ 10 5 , 0.679 ) [ 10 5 , 0.684 )
1.50 [ 10 5 , 0.672 ) [ 10 5 , 0.677 )
1.99 [ 10 5 , 0.668 ) [ 10 5 , 0.669 )
Table 7. Ranges of θ λ 2 ensuring convergence without ghost solutions (shooting and Keller–Box procedures).
Table 7. Ranges of θ λ 2 ensuring convergence without ghost solutions (shooting and Keller–Box procedures).
δ Shooting (ode 23)Shooting (ode 45)Keller–Box
0 [ 0.639 , + ) [ 0.633 , + ) [ 0 . 721 , + )
0.5 [ 0.627 , + ) [ 0.625 , + ) [ 0 . 716 , + )
1 [ 0.614 , + ) [ 0.618 , + ) [ 0 . 709 , + )
1.5 [ 0.611 , + ) [ 0.612 , + ) [ 0 . 703 , + )
1.99 [ 0.599 , + ) [ 0.603 , + ) [ 0 . 694 , + )
Table 8. Ranges of θ λ 2 ensuring convergence without ghost solutions (Three-/Four-Stage Lobatto IIIa).
Table 8. Ranges of θ λ 2 ensuring convergence without ghost solutions (Three-/Four-Stage Lobatto IIIa).
δ Three-Stage Lobatto IIIa (bpv4c)Four-Stage Lobatto IIIa (bpv5c)
0 [ 0.693 , + ) [ 0.698 , + )
0.50 [ 0.686 , + ) [ 0.691 , + )
1 [ 0.679 , + ) [ 0.684 , + )
1.50 [ 0.672 , + ) [ 0.677 , + )
1.99 [ 0.668 , + ) [ 0.669 , + )
Table 9. Number of bondes for each numerical technique.
Table 9. Number of bondes for each numerical technique.
δ Shooting
(ode 23)
Shooting
(ode 45)
Keller-BoxThree-Stage
Lobatto IIIa (bpv4c)
Four-Stage
Lobatto IIIa (bpv5c)
01140404040
0.51140404040
16440404040
1.563125404040
1.9958101404040
Publisher’s Note: MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Share and Cite

MDPI and ACS Style

Versaci, M.; Jannelli, A.; Morabito, F.C.; Angiulli, G. A Semi-Linear Elliptic Model for a Circular Membrane MEMS Device Considering the Effect of the Fringing Field. Sensors 2021, 21, 5237. https://doi.org/10.3390/s21155237

AMA Style

Versaci M, Jannelli A, Morabito FC, Angiulli G. A Semi-Linear Elliptic Model for a Circular Membrane MEMS Device Considering the Effect of the Fringing Field. Sensors. 2021; 21(15):5237. https://doi.org/10.3390/s21155237

Chicago/Turabian Style

Versaci, Mario, Alessandra Jannelli, Francesco Carlo Morabito, and Giovanni Angiulli. 2021. "A Semi-Linear Elliptic Model for a Circular Membrane MEMS Device Considering the Effect of the Fringing Field" Sensors 21, no. 15: 5237. https://doi.org/10.3390/s21155237

APA Style

Versaci, M., Jannelli, A., Morabito, F. C., & Angiulli, G. (2021). A Semi-Linear Elliptic Model for a Circular Membrane MEMS Device Considering the Effect of the Fringing Field. Sensors, 21(15), 5237. https://doi.org/10.3390/s21155237

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop