Memory Effect and Fractional Differential Dynamics in Planar Microsupercapacitors Based on Multiwalled Carbon Nanotube Arrays

Abstract: The development of portable electronic devices has greatly stimulated the need for miniaturized power sources. Planar supercapacitors are micro-scale electrochemical energy storage devices that can be integrated with other microelectronic devices on a chip. In this paper, we study the behavior of microsupercapacitors with in-plane interdigital electrodes of carbon nanotube array under sinusoidal excitation, step voltage input and sawlike voltage input. Considering the anomalous diffusion of ions in the array and interelectrode space, we propose a fractional-order equivalent circuit model that successfully describes the measured impedance spectra. We demonstrate that the response of the investigated micro-supercapacitors is linear and the system is time-invariant. The numerical inversion of the Laplace transforms for electric current response in an equivalent circuit with a given impedance leads to results consistent with potentiostatic measurements and cyclic voltammograms. The use of electrodes based on an ordered array of nanotubes reduces the role of nonlinear effects in the behavior of a supercapacitor. The effect of the disordering of nanotubes with increasing array height on supercapacitor impedance is considered in the framework of a distributed-order subdiffusion model.


Introduction
Planar micro-supercapacitors (PMSCs) are micro-scale electrochemical energy storage devices that can be embedded in a chip and integrated with other microelectronic devices. Different materials were suggested for electrodes of PMSC including graphene, nanowires, multiwalled carbon nanotubes (MWCNT), and conducting polymers (see review [1] and references therein). With an in-plane interdigital structure (Figure 1a,b), the PMSCs have a better performance than the conventional sandwich-like analogs of similar size [1]. Due to the spatial separation of the electrodes there is no need for a separator and organic binders. This simplifies the design and simulation of the device [2]. Interdigitized PMSCs possess excellent mechanical and electrical characteristics [3]. Due to design, the decreased ionic diffusion path between electrodes leads to small ion transport resistance, that results in higher power capability [3]. The high electrical conductivity of electrodes and their porosity can mitigate the ion diffusion resistance caused by the electrolyte viscosity [4].
An important problem is the development of reliable models for engineering calculations and diagnostics of PMSC. The simplest way is to extend the effective models of standard supercapacitors taking the geometry and specifics of PMSCs into account. One of the successful directions in phenomenological modeling of conventional supercapacitors is based on the use of fractional ( jω) −α . Bertrand, Sabatier et al., and Allagui et al. [8,[16][17][18] have shown that fractional order models can effectively describe the dynamical behavior of electrochemical power sources and can be used to build SOC (state of charge) and/or SOH (state of health) estimators. Fractional differential models of supercapacitors were used to describe their operation in pulse power systems [19], to describe the impedance in a frequency range [5,7,20], to study temporal response of supercapacitors [11,[21][22][23], to identify parameters from transient characteristics [24,25], to control supercapacitor systems using fractional state models [7], or fractional nonlinear systems [8,26].
Energies 2020, xx, 5 2 of 14 calculus [5][6][7][8][9][10][11][12]. Usually, these models are formulated in terms of equivalent circuits containing a constant phase element (CPE) sometimes called a capacitor of fractional order or a fractal element [13][14][15]. This element has an impedance with a frequency dependence proportional to (jω) −α . Bertrand, Sabatier et al., and Allagui et al. [8,[16][17][18] have shown that fractional order models can effectively describe the dynamical behavior of electrochemical power sources and can be used to build SOC (state of charge) and/or SOH (state of health) estimators. Fractional differential models of supercapacitors were used to describe their operation in pulse power systems [19], to describe the impedance in a frequency range [5,7,20], to study temporal response of supercapacitors [11,[21][22][23], to identify parameters from transient characteristics [24,25], to control supercapacitor systems using fractional state models [7], or fractional nonlinear systems [8,26]. In this paper, for the first time, a fractional differential model for PMSC is proposed. We show that the response of the studied supercapacitors with carbon nanotube-based electrodes is linear and the system is time-invariant. The proposed model well describes the behavior of PMSC that is confirmed by various experimental methods. We assume that the use of electrodes with an array of rod-like nanostructures should reduce the role of nonlinear effects in PMSC behavior. The effect of weak disordering of CNT array for increased height is considered in the framework of a distributed-order subdiffusion model.

Carbon Nanotubes in Planar Microsupercapacitors
The studied PMSC samples were fabricated on silicon or borosilicate glass substrates by means of standard microelectronics technologies. In the case of silicon substrates, they were pre-oxidized to form a layer of high-quality dielectric. The oxide thickness was~0.3 µm. Next, a layer of metal In this paper, for the first time, a fractional differential model for PMSC is proposed. We show that the response of the studied supercapacitors with carbon nanotube-based electrodes is linear and the system is time-invariant. The proposed model well describes the behavior of PMSC that is confirmed by various experimental methods. We assume that the use of electrodes with an array of rod-like nanostructures should reduce the role of nonlinear effects in PMSC behavior. The effect of weak disordering of CNT array for increased height is considered in the framework of a distributed-order subdiffusion model.

Carbon Nanotubes in Planar Microsupercapacitors
The studied PMSC samples were fabricated on silicon or borosilicate glass substrates by means of standard microelectronics technologies. In the case of silicon substrates, they were pre-oxidized to form a layer of high-quality dielectric. The oxide thickness was~0.3 µm. Next, a layer of metal Energies 2020, 13, 213 3 of 14 (titanium) was applied to the plates. The measurements were carried out under normal atmospheric conditions. During measurement, the samples were placed into the shielded volume to protect against electromagnetic interference. Contacting the electrodes on the chip was carried out using silver-plated clamps. Multiwalled carbon nanotubes were synthesized by catalytic PECVD process (see details in [27]) ( Figure 1). The capacities of a supercapacitor with a solid electrolyte based on polyvinyl alcohol (PVA) were investigated. Phosphoric acid H 3 PO 4 or lithium chloride LiCl were used as the ion-conducting agent.
The model developed here is applied to description of PMSC with electrodes composed of MWCNT forest. Two samples of PMSC are analyzed. They are labeled as D33 and A53, and differ in the width of stripes with MWCNTs: 30 µm for D33 and 50 µm for A53. Otherwise, they are the same, including the area. The measured voltammograms ( Figure 2) for fabricated PMSCs demonstrate shapes indicating the absence of electrochemical reactions on electrode-electrolyte interface. These shapes are typical for supercapacitors, where electric double layer plays the main role in charging kinetics. As shown below, these samples are characterized by linear response uniform in time that make these devices very convenient for description and prediction of their behavior under various conditions. Energies 2020, xx, 5 3 of 14 (titanium) was applied to the plates. The measurements were carried out under normal atmospheric conditions. During measurement, the samples were placed into the shielded volume to protect against electromagnetic interference. Contacting the electrodes on the chip was carried out using silver-plated clamps. Multiwalled carbon nanotubes were synthesized by catalytic PECVD process (see details in [27]) ( Figure 1). The capacities of a supercapacitor with a solid electrolyte based on polyvinyl alcohol (PVA) were investigated. Phosphoric acid H 3 PO 4 or lithium chloride LiCl were used as the ion-conducting agent. The model developed here is applied to description of PMSC with electrodes composed of MWCNT forest. Two samples of PMSC are analyzed. They are labeled as D33 and A53, and differ in the width of stripes with MWCNTs: 30 µm for D33 and 50 µm for A53. Otherwise, they are the same, including the area. The measured voltammograms ( Figure 2) for fabricated PMSCs demonstrate shapes indicating the absence of electrochemical reactions on electrode-electrolyte interface. These shapes are typical for supercapacitors, where electric double layer plays the main role in charging kinetics. As shown below, these samples are characterized by linear response uniform in time that make these devices very convenient for description and prediction of their behavior under various conditions.

Linear Response and Anomalous Diffusion Model
In this section, we consider the linear response in terms of the integral equation describing a renewal generally non-Markovian process. By using the Mellin tranform of its kernel, it can be shown that a time-invariant dynamical process obeys a fractional differential equation of a complex order without additional assumptions [28]. The process is characterized by the order (spectral) distribution, so this representation incorporates both the standard Debye process and the non-Debye process family depending on choice of the spectral distribution.

Linear Response and Anomalous Diffusion Model
In this section, we consider the linear response in terms of the integral equation describing a renewal generally non-Markovian process. By using the Mellin tranform of its kernel, it can be shown that a time-invariant dynamical process obeys a fractional differential equation of a complex order without additional assumptions [28]. The process is characterized by the order (spectral) distribution, so this representation incorporates both the standard Debye process and the non-Debye process family depending on choice of the spectral distribution. Generally speaking, the object under investigation in the linear response theory (LRT) is an abstract system with two points: input and output. The input signal v(t) and output signal u(t) are linked via some functional relation Here, F is a response functional and a denotes an initial moment of input signal appearance: v(t) = 0 for t < a. The lower limit a of the functional can be taken as −∞, which allows to include in consideration signals appeared in remote past or composed of several segments separated by silence intervals. For short, we will write F t instead of −∞ F t .
In the frame of LRT, the response functional F should satisfy the following two requirements.
• It is linear: for any pair of functions v 1 (t) and v 2 (t) of some family of functions v(t) and any pair of constants c 1 and c 2 , the following equality takes place: This functional can be represented in the following integral form, • If the kernel K is invariant with respect to shift, K(t, t ) = G(t − t ), and G(t) = 0 for t < 0, we deal with the Standard Linear Response Theory (SRLT).
The response formula (1) takes the form of a convolution integral An interesting behavior of discharging current in Panasonic supercapacitors (0.22, 0.47 and 1 F) was revealed in [22,23]. This behavior conforms neither of known models based on equivalent circuit because the response is nonuniform in time. Here, we will show that the response of PMSC-MWCNT is linear and the system is time-invariant.
The invariance of response kernel with respect to shift in time opens a wide field for using such power analytical tools as integral transforms. The Mellin transform casts additional light on the sense of fractional operators. Representing function G(t) via its Mellin transform and inserting the result into integral term of Equation (2): Here, is an integral operator of complex order µ = 1 − σ. Passage from the integral to fractional derivatives can be realized say by means of regularization procedure (derivation of a finite part under µ < 0). After multiplying on the weight factor and integrating over σ, this term yields the fractional derivative of distributed order according to the spectral function w(s) [29,30]. As a result, Equation (2) takes the form The spectral function w(σ) defines the LRT kinetics and links it with dynamics of an open system [28].
The application of fractional equations to the analysis of supercapacitors is usually justified by heterogeneity and complexity of porous electrodes [16,23,[31][32][33]. These equations are often interpreted in terms of anomalous diffusion processes. Anomalous diffusion is usually associated with power law expansion of the diffusion packet, ∆(t) ∝ t α/2 . The case 0 < α < 1 corresponds to subdiffusion, and the case α > 1 to superdiffusion.
Different specific mechanisms can lead to anomalous diffusion processes [33,34]. The features of anomalous diffusion are often described within the continuous time random walk (CTRW) model. The decoupled version of CTRW model implies that different jump lengths and waiting (trapping) times between two successive jumps are independent random variables. Fractional-order diffusion equation can be obtained as the asymptotic equation of CTRW with power law distributions of random trapping times τ: P{τ > t} ∝ t −ν (0 < ν < 1). Considering CTRW with space dependent jump probabilities, Barkai et al. [35] derived a time-fractional Fokker-Planck equation (FFPE) for the case, when the mean waiting time diverges. FFPE describes anomalous diffusion in an external force field, it is derived in the following form [35,36], where c(x, t) is particle concentration, V(x) is an external potential, K is a generalized diffusion coefficient, kT is the thermal energy. Here, is the Riemann-Liouville derivative of fractional order 1 − ν [37]. Fundamental solutions to Equation (3) can be found in [34,35,38].
The fractional diffusion equation for the simplest case V(x) = 0 has the following form [39], The Fourier transform of Equation (4) leads to the expression The flux density for diffusion described by (4) takes the form Energies 2020, 13, 213 6 of 14 The subdiffusive generalization of Warburg's impedance for a semi-infinite medium can be written in the form [40] Z Here, B is a frequency-independent parameter. Based on fractional Equations (3) and (4), different generalizations of diffusion impedances can be obtained for different geometries and boundary conditions (see [32,[40][41][42] and references therein). In our recent paper [43], effects produced by ion subdiffusion in lithium-ion cell are described within the fractional-order generalization of the single-particle model. Here, we use relation of fractional-order impedance to anomalous ion diffusion to interpret the behavior of PMSC.

Equivalent Circuit Model
Considering the anomalous diffusion of ions in planar electrodes of PMSC and electrolyte in the interelectrode space, we obtained the following simple equivalent circuit of PMSC ( Figure 3). It contains CPEs corresponding to ion subdiffusion in the MWCNT array and, in the interelectrode space, a capacitor C DL and a resistor R DL related to the electric double layer capacity and electrolyte resistance. The corresponding impedance is Fitting measured impedance spectra, we determined parameters of equivalent circuit for samples D33 and A53 at room temperature (Table 1). To fit the data, we used the program "EIS Spectrum Analyzer" [44]. An example of such fitting for D33 sample is shown in Figure 4. The value of parameter β can be explained by the weak influence of ion diffusion in the interelectrode space on the impedance spectrum; therefore, the interelectrode region can be described using a series-connected resistor and an almost "ideal" capacitance. Parameter α is related to ion subdiffusion in MWCNT array. According to Equation (7), value α = 0.7 corresponds to subdiffusion parameter ν = 2(1 − α) = 0.6.
Energies 2020, xx, 5 6 of 14 Based on fractional Equations (3) and (4), different generalizations of diffusion impedances can be obtained for different geometries and boundary conditions (see [32,[40][41][42] and references therein). In our recent paper [43], effects produced by ion subdiffusion in lithium-ion cell are described within the fractional-order generalization of the single-particle model. Here, we use relation of fractional-order impedance to anomalous ion diffusion to interpret the behavior of PMSC.

Equivalent Circuit Model
Considering the anomalous diffusion of ions in planar electrodes of PMSC and electrolyte in the interelectrode space, we obtained the following simple equivalent circuit of PMSC ( Figure 3). It contains CPEs corresponding to ion subdiffusion in the MWCNT array and, in the interelectrode space, a capacitor C DL and a resistor R DL related to the electric double layer capacity and electrolyte resistance. The corresponding impedance is Fitting measured impedance spectra, we determined parameters of equivalent circuit for samples D33 and A53 at room temperature (Table 1). To fit the data, we used the program "EIS Spectrum Analyzer" [44]. An example of such fitting for D33 sample is shown in Figure 4. The value of parameter β can be explained by the weak influence of ion diffusion in the interelectrode space on the impedance spectrum; therefore, the interelectrode region can be described using a series-connected resistor and an almost "ideal" capacitance. Parameter α is related to ion subdiffusion in MWCNT array. According to (7), value α = 0.7 corresponds to subdiffusion parameter ν = 2(1 − α) = 0.6.      To find electric current as a response of PMSC, we use the numerical inverse Laplace transformation of the following relation:Ĩ(s) =Ṽ(s)[Z(s)] −1 . Note that usually this relation leads to disagreement of numerical results with measured data, because Z(s) depends on operating voltage and it can not be considered as a constant characteristic of the device for all voltages. Moreover, nonlinear effects can take place, whereas Z(s) is obtained for the linearized process.
We use the following input voltage signals and their Laplace transforms. For potentiostatic charging V(t) and its Laplace transform are where H is the Heaviside function. For the discharging process, we have where θ is a charging time. To plot discharging curves in log-log scale, we shift obtained I(t) in time to charging duration θ. For cyclic voltammetry, the following voltage signal is used Its Laplace transform has the form V(s) = 2V 0 1 We  To find electric current as a response of PMSC, we use the numerical inverse Laplace transformation of the following relation:Ĩ(s) =Ṽ(s)[Z(s)] −1 . Note that usually this relation leads to disagreement of numerical results with measured data, because Z(s) depends on operating voltage and it can not be considered as a constant characteristic of the device for all voltages. Moreover, nonlinear effects can take place, whereas Z(s) is obtained for the linearized process.
We use the following input voltage signals and their Laplace transforms. For potentiostatic charging V(t) and its Laplace transform are where H is the Heaviside function. For the discharging process, we have where θ is a charging time. To plot discharging curves in log-log scale, we shift obtained I(t) in time to charging duration θ. For cyclic voltammetry, the following voltage signal is used Its Laplace transform has the form We   Table 1 (points).  Table 1 (points).
Effective capacities of the studied PMSC are estimated to be 102 µF for D33 and 54 µF for A53. D33 is characterized by a smaller electrode stripe width (30 µm) in comparison with A53 (50 µm). The difference in capacities of these samples can be explained by the difference in surface area of contact between the MWCNT array and electrolyte. Despite the fact that the MWCNT array is porous and saturated with electrolyte, it is quite dense and the surface layers of the MWCNT array-electrolyte contact play the main role.
The deviation of the tails of theoretical discharge curves from the experimental data can be explained by variations in the parameter of the constant phase element considered in the next section. These variations can be related to disordering of the array with the height of MWCNT nanotubes.

Distributed-Order Fractional Model
To determine the dependence of PMSC capacitance on nanotube height, the MWCNT arrays with similar morphology, but of different heights, were obtained by varying the synthesis time. The dependence of array height on the synthesis time is shown in Figure 7 (left panel). Next, the capacitances of PMSC with PVA+LiCl electrolyte were measured. Figure 7 (right panel) shows the dependence of specific capacity on the height of the array. The capacitive efficiency of using the MWCNT array height h, or ratio C/h, decreases with h. The use of hydrophilic treatments can increase the PMSC capacity by activating the unused area of MWCNT array, however, the power of the device does not increase. Such behavior can be explained by the small variation of the morphology of MWCNT forest for increased heights. These changes can modify the type of ion diffusion inside electrode at different heights (e.g., decrease the subdiffusion exponent). Based on this assumption,   Table 1 (points).  Table 1 (points).
Effective capacities of the studied PMSC are estimated to be 102 µF for D33 and 54 µF for A53. D33 is characterized by a smaller electrode stripe width (30 µm) in comparison with A53 (50 µm). The difference in capacities of these samples can be explained by the difference in surface area of contact between the MWCNT array and electrolyte. Despite the fact that the MWCNT array is porous and saturated with electrolyte, it is quite dense and the surface layers of the MWCNT array-electrolyte contact play the main role.
The deviation of the tails of theoretical discharge curves from the experimental data can be explained by variations in the parameter of the constant phase element considered in the next section. These variations can be related to disordering of the array with the height of MWCNT nanotubes.

Distributed-Order Fractional Model
To determine the dependence of PMSC capacitance on nanotube height, the MWCNT arrays with similar morphology, but of different heights, were obtained by varying the synthesis time. The dependence of array height on the synthesis time is shown in Figure 7 (left panel). Next, the capacitances of PMSC with PVA+LiCl electrolyte were measured. Figure 7 (right panel) shows the dependence of specific capacity on the height of the array. The capacitive efficiency of using the MWCNT array height h, or ratio C/h, decreases with h. The use of hydrophilic treatments can increase the PMSC capacity by activating the unused area of MWCNT array, however, the power of the device does not increase. Such behavior can be explained by the small variation of the morphology of MWCNT forest for increased heights. These changes can modify the type of ion diffusion inside electrode at different heights (e.g., decrease the subdiffusion exponent). Based on this assumption, Effective capacities of the studied PMSC are estimated to be 102 µF for D33 and 54 µF for A53. D33 is characterized by a smaller electrode stripe width (30 µm) in comparison with A53 (50 µm). The difference in capacities of these samples can be explained by the difference in surface area of contact between the MWCNT array and electrolyte. Despite the fact that the MWCNT array is porous and saturated with electrolyte, it is quite dense and the surface layers of the MWCNT array-electrolyte contact play the main role.
The deviation of the tails of theoretical discharge curves from the experimental data can be explained by variations in the parameter of the constant phase element considered in the next section. These variations can be related to disordering of the array with the height of MWCNT nanotubes.

Distributed-Order Fractional Model
To determine the dependence of PMSC capacitance on nanotube height, the MWCNT arrays with similar morphology, but of different heights, were obtained by varying the synthesis time. The dependence of array height on the synthesis time is shown in Figure 7 (left panel). Next, the capacitances of PMSC with PVA+LiCl electrolyte were measured. Figure 7 (right panel) shows the dependence of specific capacity on the height of the array. The capacitive efficiency of using the MWCNT array height h, or ratio C/h, decreases with h. The use of hydrophilic treatments can increase the PMSC capacity by activating the unused area of MWCNT array, however, the power of the device does not increase. Such behavior can be explained by the small variation of the morphology of MWCNT forest for increased heights. These changes can modify the type of ion diffusion inside electrode at different heights (e.g., decrease the subdiffusion exponent). Based on this assumption, we propose distributed-order fractional circuit model and derive corresponding impedances for different distribution of subdiffusion parameter on height. Consider two cases of distributed-order fractional models for PMSC. The first model implies subdiffusive behavior of fractional order homogeneously distributed over entire volume of electrodes. The second case assumes changing of subdiffusion exponent (fractional order) with height of the array. We show that these models lead to different PMSC impedances and the second case is more preferable to describe our measurements.

Subdiffusion Exponent Distributed Over Electrodes
For distributed-order fractional model, the subdiffusion equation can be written in the following form [33], where ρ(ν) is a probability density function for subdiffusion exponent 0 < ν ≤ 1, K is a generalized diffusion coefficient. Equation (9) can be derived from the continuity equation with the following relation between current density and concentration: Using the Fourier transformation of this distributed-order fractional Fick's law, and solution of the diffusion equation in the Fourier domain, Figure 7. Dependence of MWCNT height on duration of array synthesis (a) and the efficiency of using the MWCNT array height for the device capacitance (b). Circles correspond to experimental data, lines are approximation by a cubic polynomial function (a) and a linear function (b).
Consider two cases of distributed-order fractional models for PMSC. The first model implies subdiffusive behavior of fractional order homogeneously distributed over entire volume of electrodes. The second case assumes changing of subdiffusion exponent (fractional order) with height of the array. We show that these models lead to different PMSC impedances and the second case is more preferable to describe our measurements.

Subdiffusion Exponent Distributed Over Electrodes
For distributed-order fractional model, the subdiffusion equation can be written in the following form [33], where ρ(ν) is a probability density function for subdiffusion exponent 0 < ν ≤ 1, K is a generalized diffusion coefficient. Equation (9) can be derived from the continuity equation with the following relation between current density and concentration: Using the Fourier transformation of this distributed-order fractional Fick's law, and solution of the diffusion equation in the Fourier domain, we obtain the following impedance, Energies 2020, 13, 213 10 of 14 For the case of weak dependence of K ν on ν and uniform distribution of order between ν min and ν max , the impedance takes the form where B is a constant parameter.

Variation of Subdiffusion Order with MWCNT Height
Let us assume the presence of variation in morphology of CNT forest with height of array. Subdiffusion perpendicular to MWCNT inside the electrodes is described by Equation (4), but the order ν depends on height. A separate layer of thickness dh can be characterized by impedance Equation (7). These layers are assumed to operate in parallel, and the equivalent circuit presented in Figure 8 should be used instead of the previous model shown in Figure 3.
we obtain the following impedance, For the case of weak dependence of K ν on ν and uniform distribution of order between ν min and ν max , the impedance takes the form where B is a constant parameter.

Variation of Subdiffusion Order with MWCNT Height
Let us assume the presence of variation in morphology of CNT forest with height of array. Subdiffusion perpendicular to MWCNT inside the electrodes is described by Equation (4), but the order ν depends on height. A separate layer of thickness dh can be characterized by impedance (7). These layers are assumed to operate in parallel, and the equivalent circuit presented in Figure 8 should be used instead of the previous model shown in Figure 3. As a result, we obtain the following impedance of anomalous diffusion in electrodes, If we assume uniform distribution of ν between ν min and ν max and constant b ν , we have In Figure 8b, Nyquist plots are presented for four cases of subdiffusion index values. Points denoted as ν max and ν min corresponds to impedance of the circuit presented in Figure 8 with single α-CPE, characterized by ν min = 0.5 or ν max = 0.9. Points denoted by V (volume) and H (height) correspond to the similar circuit including impedances described by Equation (10) or Equation (12), respectively. The shapes of these Nyquist plots look similar, but slope strongly depend on the type of distribution of ν. Figure 9 demonstrates charging and discharging curves and voltammograms for the circuit presented in Figure 8 with ν min = 0.5, ν max = 0.9.
Other parameters are As a result, we obtain the following impedance of anomalous diffusion in electrodes, If we assume uniform distribution of ν between ν min and ν max and constant b ν , we have In Figure 8b, Nyquist plots are presented for four cases of subdiffusion index values. Points denoted as ν max and ν min corresponds to impedance of the circuit presented in Figure 8 with single α-CPE, characterized by ν min = 0.5 or ν max = 0.9. Points denoted by V (volume) and H (height) correspond to the similar circuit including impedances described by Equation (10) or Equation (12), respectively. The shapes of these Nyquist plots look similar, but slope strongly depend on the type of distribution of ν. Figure 9 demonstrates charging and discharging curves and voltammograms for the circuit presented in Figure 8 with ν min = 0.5, ν max = 0.9. Other parameters are R = 1275 Ohm, C β = 0.101 × 10 −3 s β /Ohm, β = 0.91, R DL = 0, C DL = 12.5 µF. Discharging curves indicate the presence of the memory effect, i.e., dependence on the charging prehistory. Cyclic voltammograms of the distributed-order model have the form similar to the single order model. The difference can be seen by variation of ν min . Figure 10 shows two IV-curves for ν min = 0.5 and ν min = 0.7. Parameters Energies 2020, 13, 213 11 of 14 are the same as for Figure 9 except ν min . Due to effect of disordering, smaller values of ν min can be related to higher MWCNT array. For increased height of array, the capacity of PMSC grows, but the dependence of C on h is not direct proportionality. The shape of cyclic voltammogram for the same scanning rate is changed. This shape is closer to rectangle for higher values of ν min . This is qualitative explanation of the behavior observed in Figure 7. R = 1275 Ohm, C β = 0.101 × 10 −3 s β /Ohm, β = 0.91, R DL = 0, C DL = 12.5 µF. Discharging curves indicate the presence of the memory effect, i.e., dependence on the charging prehistory. Cyclic voltammograms of the distributed-order model have the form similar to the single order model. The difference can be seen by variation of ν min . Figure 10 shows two IV-curves for ν min = 0.5 and ν min = 0.7. Parameters are the same as for Figure 9 except ν min . Due to effect of disordering, smaller values of ν min can be related to higher MWCNT array. For increased height of array, the capacity of PMSC grows, but the dependence of C on h is not direct proportionality. The shape of cyclic voltammogram for the same scanning rate is changed. This shape is closer to rectangle for higher values of ν min . This is qualitative explanation of the behavior observed in Figure 7.

Conclusions
We studied the behavior of PMSC with electrodes of MWCNT array under sinusoidal excitation, step voltage input, and sawlike voltage input. Studied PMSC samples were fabricated on glass substrates by means of standard microelectronics technologies. MWCNTs were synthesized by a catalytic PECVD process. For the first time, a fractional differential model for planar supercapacitor is proposed. The corresponding equivalent circuit model, taking the anomalous diffusion of ions in the array and interelectrode space into account, successfully describes the measured impedance spectra. We demonstrate that the response of the investigated micro-supercapacitors is linear and uniform in time. The numerical inversion of the Laplace transforms for electric current response in an equivalent circuit with a given impedance leads to results consistent with potentiostatic measurements and cyclic voltammograms. Thus, the use of electrodes based on an ordered array of nanotubes reduces the role of nonlinear effects in the behavior of a supercapacitor. In addition, the effect of the disordering of nanotubes with increasing array length on the impedance of a supercapacitor is discussed in the framework of a distributed-order subdiffusion model. The dependence of PMSC capacitance on Figure 9. Theoretical charging-discharging curves (a) and cyclic voltammograms (b) predicted by the distributed-order fractional model.
Energies 2020, xx, 5 11 of 14 R = 1275 Ohm, C β = 0.101 × 10 −3 s β /Ohm, β = 0.91, R DL = 0, C DL = 12.5 µF. Discharging curves indicate the presence of the memory effect, i.e., dependence on the charging prehistory. Cyclic voltammograms of the distributed-order model have the form similar to the single order model. The difference can be seen by variation of ν min . Figure 10 shows two IV-curves for ν min = 0.5 and ν min = 0.7. Parameters are the same as for Figure 9 except ν min . Due to effect of disordering, smaller values of ν min can be related to higher MWCNT array. For increased height of array, the capacity of PMSC grows, but the dependence of C on h is not direct proportionality. The shape of cyclic voltammogram for the same scanning rate is changed. This shape is closer to rectangle for higher values of ν min . This is qualitative explanation of the behavior observed in Figure 7.

Conclusions
We studied the behavior of PMSC with electrodes of MWCNT array under sinusoidal excitation, step voltage input, and sawlike voltage input. Studied PMSC samples were fabricated on glass substrates by means of standard microelectronics technologies. MWCNTs were synthesized by a catalytic PECVD process. For the first time, a fractional differential model for planar supercapacitor is proposed. The corresponding equivalent circuit model, taking the anomalous diffusion of ions in the array and interelectrode space into account, successfully describes the measured impedance spectra. We demonstrate that the response of the investigated micro-supercapacitors is linear and uniform in time. The numerical inversion of the Laplace transforms for electric current response in an equivalent circuit with a given impedance leads to results consistent with potentiostatic measurements and cyclic voltammograms. Thus, the use of electrodes based on an ordered array of nanotubes reduces the role of nonlinear effects in the behavior of a supercapacitor. In addition, the effect of the disordering of nanotubes with increasing array length on the impedance of a supercapacitor is discussed in the framework of a distributed-order subdiffusion model. The dependence of PMSC capacitance on

Conclusions
We studied the behavior of PMSC with electrodes of MWCNT array under sinusoidal excitation, step voltage input, and sawlike voltage input. Studied PMSC samples were fabricated on glass substrates by means of standard microelectronics technologies. MWCNTs were synthesized by a catalytic PECVD process. For the first time, a fractional differential model for planar supercapacitor is proposed. The corresponding equivalent circuit model, taking the anomalous diffusion of ions in the array and interelectrode space into account, successfully describes the measured impedance spectra. We demonstrate that the response of the investigated micro-supercapacitors is linear and uniform in time. The numerical inversion of the Laplace transforms for electric current response in an equivalent circuit with a given impedance leads to results consistent with potentiostatic measurements and cyclic voltammograms. Thus, the use of electrodes based on an ordered array of nanotubes reduces the role of nonlinear effects in the behavior of a supercapacitor. In addition, the effect of the disordering of nanotubes with increasing array length on the impedance of a supercapacitor is discussed in the framework of a distributed-order subdiffusion model. The dependence of PMSC capacitance on the height of MWCNT arrays in electrodes is measured. The capacity grows with h, while the ratio C/h decreases. The explanation of such behavior is given within the fractional differential model assuming modification of the array morphology for the increased height that leads to a variation of the subdiffusion exponent.
Funding: This work is supported by the Ministry of Science and Higher Education of the Russian Federation through the assignment of 2020 year for Scientific-Manufacturing Complex "Technological Centre". Renat Sibatov acknowledges the financial support from the Russian Science Foundation (project 19-71-10063).

Conflicts of Interest:
The authors declare no conflict of interest.