BRNO UNIVERSITY OF TECHNOLOGY VYSOKÉ UČENÍ TECHNICKÉ V BRNĚ

BRNO UNIVERSITY OF TECHNOLOGY

VYSOKÉ UČENÍ TECHNICKÉ V BRNĚ

FACULTY OF CHEMISTRY

FAKULTA CHEMICKÁ

INSTITUTE OF MATERIALS SCIENCE

ÚSTAV CHEMIE MATERIÁLŮ

MODELLING OF BIOELECTRONIC DEVICES

MODELOVÁNÍ PRVKŮ PRO BIOELEKTRONIKU

MASTER'S THESIS

DIPLOMOVÁ PRÁCE

AUTHOR

AUTOR PRÁCE

Bc. Jan Truksa

SUPERVISOR

VEDOUCÍ PRÁCE

doc. Ing. Ota Salyk, CSc.

BRNO 2018

Fakulta chemická, Vysoké učení technické v Brně / Purkyňova 464/118 / 612 00 / Brno

Zadání diplomové práce

Číslo práce: FCH-DIP1196/2017

Ústav: Ústav chemie materiálů

Student: Bc. Jan Truksa

Studijní program: Chemie, technologie a vlastnosti materiálů Studijní obor: Chemie, technologie a vlastnosti materiálů Vedoucí práce: doc. Ing. Ota Salyk, CSc.

Akademický rok: 2017/18

Název diplomové práce:

Modelování prvků pro bioelektroniku

Zadání diplomové práce:

1. Rešerše literatury na téma organický elektrochemický tranzistor 2. Základy programu COMSOL pro modelování fyzikálních polí.

Termín odevzdání diplomové práce: 7.5.2018

Diplomová práce se odevzdává v děkanem stanoveném počtu exemplářů na sekretariát ústavu. Toto zadání je součástí diplomové práce.

Bc. Jan Truksa student(ka)

doc. Ing. Ota Salyk, CSc.

vedoucí práce

prof. RNDr. Josef Jančář, CSc.

vedoucí ústavu

V Brně dne 31.1.2018 prof. Ing. Martin Weiter, Ph.D.

děkan

3

ABSTRACT

The topic of this thesis is computer modelling of an organic electrochemical transistor (OECT). To create such a model, the electric field and ion concentration were numerically computed, using the finite element method. The electric potential on top of the OECT channel, the changes in conductivity and the output current of the OECT were computed. To carry out the computation, a standard personal computer and the commercial software COMSOL Multiphysics were utilized. Due to a lack of computational power, the model had to be split into parts and drastically simplified. The presented results differ from those in literature, as the saturation of the transistor is not modelled correctly. This deviation from real OECT behaviour is likely caused by the simplification of the model.

ROZŠÍŘENÝ ABSTRAKT

Organické elektrochemické tranzistory (OECT) jsou v současnosti předmětem zájmu v oboru bioelektroniky pro sledování a ovlivňování biologických procesů v živých organismech a tkáních jako senzorů a sond. Toto zařízení je založeno na třech elektrodách. Dvě z nich, spojené vrstvou vodivého organického materiálu (molekulárního nebo častěji polymerního), jsou nazvány Source (S) – zdrojová elektroda, a Drain (D) – odtoková elektroda. Zmíněná vrstva vodivého materiálu je nazvána kanál. Třetí elektroda, Gate (G) – řídící elektroda, je oddělená od kanálu vrstvou elektrolytu. Funkce OECT je založena na propustnosti vrstvy organického vodivého materiálu pro ionty elektrolytu.

Nosičem náboje v organických materiálech jsou nejčastěji díry (jedná se tedy o vodiče typu p). Díry jsou ve vodivém polymeru generované anionty. Tyto ionty mohou být přidány do materiálu při výrobě. Látka, obsahující tyto ionty, přidaná do polymeru, se nazývá dopant.

Druhou možností je stabilizace děr anionty elektrolytu, pokud je vodivá látka bez dopantu. Na základě použitého materiálu lze rozdělit OECT na dopované a dedopované.

Na elektrodu D je obvykle přikládáno záporné napětí a elektroda S je uzemněna.

V dopovaném OECT (je v něm použito polymeru s dopantem) v tomto stavu teče kanálem elektrický proud. Pokud je na elektrodu G přivedeno kladné napětí, dochází k pohybu kationtů a jejich průniku do kanálu. Uvnitř kanálu poté kationty elektrolytu stabilizují ionty dopantu, čímž dochází k rekombinaci děr a jejich úbytku, zatímco nábojová neutralita je zachována. Zároveň s poklesem koncentrace děr dochází k poklesu vodivosti materiálu, protože kationty, které je nahradily, mají řádově menší pohyblivost, a proud, protékající kanálem klesá. V kanálu dedopovaného OECT je zanedbatelná koncentrace nosičů náboje, a pokud je pod napětím pouze elektroda D, elektrický proud neteče. Pokud je na elektrodu G přivedeno záporné napětí, dochází k průniku aniontů dovnitř kanálu, kde tyto ionty generují další díry. Proud, protékající kanálem, poté roste.

Pomocí procesů dopování a dedopování tedy lze transformovat změny napětí na elektrodě G, nebo změny koncentrace iontů v elektrolytu na změnu elektrického proudu. Tuto vlastnost lze využít například při výrobě senzorů.

Nejčastěji používaným vodivým materiálem pro výrobu OECT je v současnosti poly(ethylen dioxythiofen) dopovaný polystyren sulfonovou kyselinou (PEDOT:PSS). Výhodou tohoto polymeru je jeho dostupnost, nízká cena, snadné nanášení a biokompatibilita. Díky možnosti

4 kultivace buněk na povrchu PEDOTu:PSS je OECT použitelný v oblasti bioelektroniky.

Nevýhodou je nutnost pracovat při velmi nízkém napětí, jinak hrozí nevratná redoxní změna polymeru. Bohužel je tranzistor limitován rychlostí pohybu iontů, což ztěžuje jeho aplikaci ve výpočetní technice nebo zpracování signálu.

Cílem této práce je počítačové modelování OECT pro optimalizaci jeho geometrie.

K úspěšnému vytvoření modelu je třeba vypočítat rozložení elektrického pole a koncentrace iontů. Elektrické pole je definováno Poissonovou rovnicí, koncentrace iontů v elektrolytu je definována Nernst-Planckovou rovnicí. Koncentrace iontů ovlivňuje elektrický náboj uvnitř elektrolytu, elektrické pole naopak řídí iontový drift, čímž jsou obě rovnice spojeny. Celkový elektrický náboj iontů na povrchu kanálu, které prochází do polymerní vrstvy, poté ovlivňuje koncentraci nosičů náboje v kanálu. Dle Ohmova zákona poté dojde ke změně vodivosti a tedy i výstupního proudu.

Uvedené parciální diferenciální rovnice je velmi obtížné řešit analyticky. Z tohoto důvodu je v práci použito počítačové modelování metodou konečných prvků, pomocí komerčního softwaru COMSOL Multiphysics. Kvůli nedostatečnému výkonu použitého počítače musel být model rozdělen na tři části. V první části je modelováno elektrické pole v elektrolytu nad všemi elektrodami zároveň a v polymerní vrstvě mezi elektrodami D a S bez vlivu elektrody G. Výsledkem tohoto modelu jsou tedy dvě rozložení elektrického potenciálu – uvnitř polymerní vrstvy a nad polymerní vrstvou.

Druhá část modelu je zaměřená na modelování elektrické dvojvrstvy na rozhraní polymer- elektrolyt. Na základě rozložení potenciálu, získaného v předchozí části modelu, jsou modelovány změny elektrického pole a koncentrace iontů v elektrické dvojvrstvě. Podle změn koncentrace ve dvojvrstvě lze určit změny koncentrace nábojů a vodivosti v kanálu. Ve třetí části modelu je poté na základě změn vodivosti a potenciálu na povrchu kanálu vypočítáno rozložení elektrického pole uvnitř kanálu. Elektrický proud je poté dán Ohmovým zákonem.

Získané hodnoty výstupního proudu jsou v řádu mA, což je v souladu s literaturou. Výstupní charakteristiky tranzistoru se od reálných OECT odlišují tím, že u většiny z nich není možné pozorovat ustálení proudu na konstantní hodnotě pro taková napětí na elektrodách D a G, kdy je polymerní vrstva plně saturována. Některé křivky ovšem saturaci vykazují, pravděpodobně by tedy bylo třeba provést výpočet pro vyšší hodnoty napětí na D a G, což znamená problém pro praktickou aplikaci modelu. Odchylka od reálného chování OECT je pravděpodobně způsobena zjednodušením modelu oproti skutečnosti, ke které muselo dojít kvůli nedostatečnému technickému vybavení. Pokud by tento problém byl odstraněn, bylo by možné vytvořený model použít například pro návrhy geometrie tranzistorů.

5

KEYWORDS

Organic electrochemical transistor, organic electronics, electric field, ion concentration, electric current, partial differential equations, modelling, COMSOL Multiphysics

KLÍČOVÁ SLOVA

Organický elektrochemický tranzistor, organická elektronika, elektrické pole, koncentrace iontů, elektrický proud, parciální diferenciální rovnice, modelování, COMSOL Multiphysics

6 TRUKSA, J. Modelling of Bioelectronic Devices. Brno: Vysoké učení technické v Brně, Fakulta chemická, 2018. 55 s. Vedoucí diplomové práce doc. Ing. Ota Salyk, CSc..

DECLARATION

I declare that the diploma thesis has been worked out by myself and that all the quotations from the used literary sources are accurate and complete. The content of the diploma thesis is the property of the Faculty of Chemistry of Brno University of Technology and all commercial uses are allowed only if approved by both the supervisor and the dean of the Faculty of Chemistry, BUT.

…………..………….

student's signature

I would like to thank the head of my thesis, doc. Ing. Ota Salyk, CSc. for advice and patience during the making of my master's thesis.

7

CONTENTS

1. INTRODUCTION 8

1.1. Aim of This Work 8

2. THEORETICAL PART 9

2.1. Conductive Polymers 12

2.1.1. Structure of Conductive Polymers 12

2.1.2. Charge Transfer Description 16

2.2. Recent OECT Applications 18

2.3. Fundamental Physics of OECT 19

2.3.1. Bulk Electrolyte 19

2.3.2. Electric Double Layer 19

2.3.3. Charge Transfer in the Channel 21

2.4. Finite Element Method 23

2.4.1. COMSOL Multiphysics 23

3. EXPERIMENTAL PART 25

3.1. Bulk Electrolyte 25

3.2. Electric Double Layer 29

3.3. Polymer Channel 32

4. RESULTS AND DISCUSSION 35

4.1. Bulk Electrolyte Model 35

4.2. Electric Double Layer Model 38

4.3. Polymer Channel Model 42

5. CONCLUSIONS 46

6. LIST OF ABBREVIATIONS AND SYMBOLS 48

7. LITERATURE 50

8

1. INTRODUCTION

The organic electrochemical transistor (OECT) has gained considerable attention since it was first reported in 1984. This device is based on a conductive organic material (monomolecular or polymeric), which connects two electrodes, while a third one is separated by an electrolyte.

The polymer layer is permeable for ions, which interact with charge carriers in the polymer, thus changing its conductivity. Since the drift of ions can be controlled by an electric field, it is possible to change the electric current flowing through the conductive layer of the device by superimposing an electric field from the third electrode – this electrode is called the gate (G). Several types of sensors based on the OECT have been created. These sensors typically utilize a reaction of an analyte with the gate electrode, which results in a change of the electric field [1, 2, 3, 4].

Besides the ability to modulate output current, conducting polymers are a good substrate for cultivating cells, which can also modulate the output current by their natural ion exchange.

This allows the OECT to be applied in the growing field of organic electronics and in medicine as an analytic device [5, 6].

Despite the great interest in OECTs, only a small amount of mathematical modelling has been done. This is likely because such a model demands the utilization of partial differential equations on two- and three- dimensional domains, which are often impossible to solve analytically. Numerical modelling is therefore ideal for such a task. In this work, a numerical model of an OECT is described, utilizing the finite element method of computation, which has been made considerably simpler to use by the software COMSOL Multiphysics. The electric field and electric current is modelled by the Poisson equation and Ohm's law respectively. The ionic drift is modelled by the Nernst-Planck equation. The conductivity of the polymer layer is controlled by the charge of ions drifting in and out of the layer.

1.1. Aim of This Work

The aim of this work is to create a virtual model of an OECT, which will compute the electric field present in the transistor, the ionic drift to and from the polymer layer and subsequently the changes of conductivity and output current. The values of ion concentration and by extension conductivity and current are expected to be dependent on changes of the electric field in the model. The results should be in agreement with literature, which would demonstrate validity of the model.

9

2. THEORETICAL PART

The OECT is composed of three electrodes and a thin layer of an electronically conductive polymer (ECP), called the channel. The polymer layer connects two of the electrodes, which are named Source (S) and Drain (D). The third electrode, not in contact with the ECP, is the Gate. The G electrode and the CP layer are in contact with an electrolyte. An example OECT configuration is shown on figs. 1 and 2. Typically, the S electrode is grounded, and the D electrode is biased with a negative voltage VD (V). There are two modes of OECT operation – dedoped and doped mode [1, 2, 3]. An illustration of doped mode is shown on fig. 1, and dedoped mode is illustrated on fig. 2.

In doped mode, bias on D causes charge carriers, most often holes, in the polymer layer to move, creating an electric current measured at the S electrode and denoted ID (A); the transistor is said to be in the ON state. Holes in the channel are generated by the presence of negative ions (such as sulfonate groups) in the ECP material – this structure is referred to as the dopant. When the G electrode is biased with a positive electric potential VG (V), the potential difference in the electrolyte causes cations to drift into the ECP layer and compensate some of the dopant ions, effectively decreasing the number of holes, which in turn lowers the ID current. This process is referred to as dedoping the channel. When a sufficient amount of cations is injected into the channel, current ID approaches zero – the transistor is said to be in the OFF state. Furthermore, there may be a negligible concentration of ions present in the polymer when VG = 0 V due to diffusion. A negative value of VG may lead to extracting these ions, which increases conductivity [2, 3].

In dedoped mode, the ECP is used without a dopant, so there are only a few holes in the channel. The transistor is normally in the OFF state. When the G electrode is biased with a negative voltage, anions drift into the channel, where they serve as a dopant, generating holes. This leads to an increase in the number of holes, and the current ID. This process is referred to as doping the channel [3, 4].

By utilizing the doping and dedoping processes, OECTs are able to convert (transduce) changes of VG, or chemical composition of the electrolyte (i.e. changing the concentration) into changes in the ID current. This ability is described by a transfer curve – the dependence of ID on VG. The efficiency of the OECT at a given potential VG can be calculated as the derivative of the transfer curve (eq. 1):

G D

m d

d V

g  I , (1)

where gm (S) is the transconductance. The steeper the transfer curve, the more efficient the transduction. This implies that high values of transconductance are desirable [1, 3].

10 Fig. 1 – An illustration of doped mode OECT operation. The red area signifies dedoping.

S D

G

Dopant Electrolyte ions

Holes

VG = 0 V

VG > 0 V

VG < 0 V VD < 0 V

VD < 0 V Dedoped polymer boundary

VD < 0 V

11 Fig. 2 – An illustration of dedoped mode OECT operation. The green area signifies doping.

Electrolyte ions Holes

VG = 0 V

VG > 0 V

VG < 0 V VD < 0 V

VD < 0 V

VD < 0 V

Dedoped polymer boundary

12 2.1. Conductive Polymers

Development of electrochemically active polymers has been going on since the 1970s.

Electrochemically active polymers can be divided into two groups. The first are the Redox Polymers – structures with localized oxidized or reduced sites, between which electrons can transfer through hopping. An example of this type is poly(styrene sulfonate), often used as an ion exchanger. The second group are the Electronically (Intrinsically) Conducting Polymers, which contain a conjugated system of electrons. It is also possible to form copolymers and composites, using metals, metal oxides or carbon nanotubes [5]. However, ECPs are the most relevant to the topic of OECTs, and this work will focus on them.

Their advantages over inorganic conductive materials and semiconductors lie in their biocompatibility, flexibility and ionic conductivity [1, 3, 6]. Furthermore, it is possible to fine-tune their molecular structure and properties [5,7]. They can be processed by relatively cheap technologies, such as spin coating, inkjet printing or screen printing [2, 8]. On the other hand, polymers form polycrystalline or amorphous domains, which inhibit charge transport.

Due to this fact, ECPs are most often used in the form of a thin film [9].

2.1.1. Structure of Conductive Polymers

As mentioned above, ECPs contain a polymer backbone of sp² hybridized carbon atoms, with overlapping pz orbitals. On a linear chain, the Highest Occupied Molecular Orbital (HOMO) forms the valence band, and the Lowest Unoccupied Molecular Orbital (LUMO) forms the conduction band. The structure of ECP very often contains aromatic or heterocyclic compounds, such as aniline, thiophene, or pyrrole. Commonly used ECPs include poly(3,4- ethylenedioxythiophene), stabilized by poly(styrene sulfonate) (PSS) – the material is denoted PEDOT:PSS. Other polymers are poly(3-hexylthiophene) (P3HT), polyaniline and polypyrrole. Examples of ECP structure are shown on fig. 3 [5]. The aromatic and heterocyclic structures stabilize conjugated bonds and facilitate inter- and intra-molecular interactions of polymer chains, such as π-π stacking, Van der Waals interactions and hydrogen bonds [5, 7, 9].

Conductive polymer structures are most often prepared by oxidative or reductive polymerization. Oxidation removes electrons from the HOMO orbital, producing polarons and bipolarons, which can be viewed as holes. This process effectively creates a p-doped semiconductor. On the other hand, reduction adds electrons to the LUMO orbital, creating an n-doped semiconductor. It should be noted, that the redox system used in this reaction is a dopant in the traditional sense of semiconductor physics (i.e. changing the number of electrons in the valence or conduction band). Meanwhile, a dopant in the sense of organic electronics is simply a charge stabilizer [5, 6, 9]. An example of the synthesis of PEDOT:PSS is shown in fig 4. The dopant must be added during synthesis. The use of a polymer anion (such as PSS) results in an aqueous dispersion of pH 1-2. The use of a monomolecular anion (such as p-toluenesulphonic acid) results in a blue powder [10].

13 Fig. 3 – structures of various conductive polymers [5]. a) PEDOT:PSS, b) P3HT, c) polypyrrole,

d) polyaniline.

Fig. 4 – Example of PEDOT:PSS synthesis [10].

S +

O O

n

SO ₃ ^-

n

N H

n

N

N

n

S

C ₆ H ₁₃

n

a) b)

c)

d)

S

O O

+

Fe₂(SO₄)₃ Na2S2O8 PSS

n 2n

n/3 +

S O O

S O O

S O O

SO₃^-

n

14 The morphology of ECPs is also an object of research, as it influences charge carrier mobility μ (on the order of 1-10 cm² V⁻¹ s⁻¹ in a crystalline layer, versus 10⁻²-10⁻⁷ cm² V⁻¹ s⁻¹ in an amorphous layer [9]), water intake and swelling, and ionic mobility, which has a direct effect on the doping speed [6, 11, 12]. A problem in designing ECP structures for optimal morphology is that while electronic transport is facilitated by a highly crystalline structure, ion intake is impeded by it [5, 12]. The speed of ion intake is a critical parameter, dictating the speed of switching between the ON and OFF states [12, 14, 15].

However, increased swelling of the ECP layer may lead to dopant diffusion out of the layer, especially in the case of non-polymeric dopants, such as tosylate ions [6]. Furthermore, the polymer layer may delaminate during use if a crosslinking agent, such as 3-glycidoxypropyl- trimethoxysilane, is not used. On the other hand, use of such an agent further impedes device speed [15]. Unsubstituted π-conjugated molecules typically crystallize into a layered herringbone structure [13]. By adding side chains, it is possible to enhance crystallization, solubility, ion intake and light absorption, and enlarge the conjugated system, thus facilitating charge transport [14].

For dedoped mode OECTs, either pure conjugated polymers, or conjugated polyelectrolytes are used. These polymers contain a conjugated backbone, and a hydrophilic side chain (such as a sulfonate group). Furthermore, the side chain can function as a dopant, generating charge carriers in the polymer backbone. Such a polyelectrolyte is said to be self-doping [4, 11].

Giovannitti et. al. [14] have demonstrated, that by adding glycolated side chains to a polythiophene backbone, it is possible to reach higher levels of ion intake and gain a dedoped mode OECT with a high transconductance and fast switching. Zeglio et. al. [11]

have proposed to solve the problem of water solubility of self-doped polyelectrolytes by adding counterions in the form of ammonium salt complexes. By using these materials, it is possible to achieve an OECT with excellent stability and high transconductance near zero gate bias. While most ECPs are p-doped, Giovannitti et. al. [16] have manufactured a dedoped mode OECT using the n-doped polymer poly(2,6-dibromonaphthalene-1,4,5,8-tetracarboxylic diimide). Polymers, prepared by the respective groups are shown in fig. 5. In a doped mode OECT, ion intake is helped by the presence of a polar dopant. Inal et. al. [17] have shown, that the backbone of the dopant plays a major role in OECT performance. Molecular weight of the dopant and type of the counter ion are comparatively less important.

15 Fig. 5 – Polymer structures prepared by a) Giovannitti and b) Zeglio, where R1 and R2 are

ethyleneglycole oligomers [11, 16].

The most used material for a doped mode OECT is PEDOT:PSS, due to its accessibility, ease of processing and conductivity in the order of 100 S cm⁻¹, which can be further improved by heat treatment or the use of processing co-solvents, such as dimethyl sulfoxide [3, 18, 19].

A layer of PEDOT:PSS is generally known to consist of flat grains (thickness in the order of 1 nm and diameter in the order of 10 nm) of PEDOT-rich phase, which are connected in the direction parallel with the substrate. In the perpendicular direction, these grains are separated by lamellae of PSS-rich phase. An illustration of this structure is shown in fig. 6. Nardes et.

al. [20] have compared the cross section of PEDOT:PSS structure to lasagne. Unfortunately, this phase separation leads to formation of electric double layers (EDL) inside the film, which function as capacitors (i.e. store electrical charge). As a result, the polymer layer has a volumetric capacitance C^* of around 40 F cm⁻³ [21], which can be changed by altering the concentration of the cross-linker [22].

The enhancement of PEDOT:PSS conductivity through heat treatment and co-solvent usage has been attributed to growth of the PEDOT-rich phase by thermal relaxation at temperatures above 175 °C, or interaction with polar groups of the solvent and subsequent straightening of PEDOT chains. Both mechanisms lead to an increase in hole mobility, while hole concentration remains relatively unchanged [18, 19, 23]. Alternative dopants for PEDOT include tosylate, (trifluoromethylsulfonyl)sulfonylimide (TFSI) on a polystyrene (PS) or polymethylmethacrylate backbone [17, 24].

Besides doped mode OECTs, PEDOT:PSS may be used as an antistatic in photographic films, manufacture of polyethylene or polyester film, or a protective layer in LCD production due to its transparency and effectivity [10]. In addition, Wei et. al. [24] have demonstrated the possibility to use PEDOT:PSS as a thermoelectric material.

N

N O

R

O

O O

R

S R

S R

n

S

O O

O SO

n N⁺

C₉H₂₈

H H

H

a) b)

16 Fig. 6 – PEDOT:PSS structure. The PEDOT-rich parts of the structure are light green, the PSS-rich

parts are brown. The electrolyte is light blue and the substrate is dark blue.

2.1.2. Charge Transfer Description

In ECPs, charge transfer occurs by the motion of delocalized electrons or holes through the LUMO or HOMO of a conjugated system, or in polaron energy bands. Additionally, charge carriers may transfer between chains, defects, and localized states by hopping [5, 9].

Tessler et. al. [25] have published a detailed review of models describing hopping. In the case of PEDOT:PSS, the electric charge is transferred by bipolarons, created by oxidation during the synthesis (fig. 4) [5, 6]. A structure of such a bipolaron is shown in fig. 7.

The redox state of an ECP can be further changed by the applied voltage. In the case of PEDOT:PSS, Marzocchi et. al. [26] have found, that at an applied potential of − 0.9 V, PEDOT is in a neutral state, while increasing the potential increases the number of polarons and bipolarons. The polymer is completely oxidized at a potential of 0.8 V. A schematic of this redox process is shown in fig. 7. However, PEDOT can be oxidized by an aqueous medium, if it is present. The colour of PEDOT:PSS also changes with its redox state – the polymer is dark blue in its reduced form, and becomes lighter in colour with oxidation.

Holes PSS^– Electrolyte ions

17 Fig. 7 –redox changes of PEDOT with potential [26]. The polymer becomes more oxidized with

positive electric potential.

Besides electrons and holes, ionic transport is also important, as it greatly influences the speed of OECT operation [15]. If an ionic species is involved in the redox process, it can be described by equation 2:

PEDOT⁺:A⁻ + e⁻ + M⁺ ↔ PEDOT⁰ + M⁺:A⁻, (2) where A⁻ is a dopant ion (e.g. PSS or tosylate), e⁻ is an electron, and M⁺ is a positive electrolyte ion [1]. To measure the transport efficiency of both electrons/holes and ions in an ECP, Inal et. al. [26] have introduced the product of charge carrier mobility, μ, and volumetric capacitance, C^*, of the channel. In their model, μ describes electronic charge transport across the device, and C^* describes ionic transport and charge storage in the polymer backbone.

Finally, the function of an OECT may be influenced by contact resistance on the ECP- electrode interface. This phenomenon is caused by a potential barrier at the contact, which leads to poor charge transfer across the interface. The height of this barrier is given by the difference of electrode work function wel (the amount of energy needed to free an electron from the solid state) and the HOMO or LUMO level of the polymer [27].

S

S

S

S S

O O

O O

O O

O O

O O

S

S

S C⁺

S S

O O

O O

O O

O O

O O

S

S

S C⁺

S S

O O

O O

O O

O O

O O

+

Positive potential

Positive potential

Polaron

Bipolaron Reduced PEDOT

18 2.2. Recent OECT Applications

At the present time, OECTs see most use in bioelectronics. For example, it is possible to record brain activity by placing an OECT on top of the brain of a rat. Conversely, OECTs can be used to stimulate neurons by locally supplying electrical current. An electrocardiogram can be recorded by simply placing an OECT on top of human skin. Monitoring cell culture growth and health is also a possibility, which has been demonstrated on cardiomyocytes [3, 28]. The reduced form of PEDOT is more compatible with living cells, which was demonstrated by Löffler et. al. [6] for canine kidney cells. This was explained by redox-state dependent conformation changes of the surface protein fibronectin.

Furthermore, OECTs can be used as biosensors for metabolites. This is done by reaction of a selective redox enzyme with the metabolite, and subsequent electron transfer to the gate electrode, which then modulates ID. In this way, highly selective sensors can be created and used for analysis of sweat, breath, saliva or cell cultures, giving the OECT a clinical application [3]. Bihar et. al. [29] have created a sensor of ethanol in human breath. This was done by immobilizing the enzyme alcohol dehydrogenase with the cofactor nicotineamide andenin dinucleotide on the surface of an inkjet-printed OECT with an electrolyte gel.

The presence of ethanol causes an enzymatic reaction, which causes a spike in the current ID. Zhang et. al. [30] have demonstrated OECTs sensing chirality in the amino acids tyrosine and tryptophan. This was achieved by using a material composed of the respective amino acids and o-phenylene diamine, deposited by electro-polymerization on a golden electrode. These materials were used as the gate electrode. The OECT then reacted by a significant current change when exposed to the correct amino acid (L-tyrosine electrode reacted to the presence of L-tyrosine etc.).

Wang et. al. [31] have created an OECT based on nylon fibres, polypyrrole nanowires and nanofibers of a copolymer of poly(vinyl acetate) (PVA) and polyethylene (PE). First, PVA and PE were co-polymerized and dispersed in water. Then the copolymer was deposited on the nylon fibre and polypyrrole was polymerized on its surface. The OECT was then made by crossing two strands of this fibre, one of them serving as the channel and the other as the gate electrode. A drop of electrolyte, composed of PVA and HCl, was placed on the junction. This OECT was then used to detect the presence of Pb²⁺ ions in a solution, demonstrating a useful application in environmental analysis. Furthermore, the ability to create an OECT composed of fibres shows potential in the development of wearable electronics.

It is possible to use OECTs as parts of electrical circuits, due to their high transconductance and low-voltage operation. Logic circuits and memory devices have been realized by cost- efficient techniques such as screen or inkjet printing. Capacitors and energy storage may be other possible applications for electrolyte-based devices. Systems using multiple channels or gate electrodes have been prepared for simulating brain and visual system functions.

A limiting factor for OECTs is the relatively low switching speed due to the dependence on the ion mobility. Because of the slow switching, OECTs only have a limited application in computing or signal processing [3].

19 2.3. Fundamental Physics of OECT

From a physical point of view, the OECT can be split into three parts. First is the bulk electrolyte between the G electrode and the channel, together with the S and D electrodes.

This part of the model describes the electric field between the electrodes. Next is the electric double layer on the electrolyte-channel or electrolyte-gate interface. This part of the model describes the electric field at the electrolyte interfaces, as well as changes in ion concentration, based on the electric field computed in the first part. This part of the model provides information on the doping of the polymer or interactions at the G electrode, if there is an analyte present. Finally, charge transport inside the polymer layer is modelled, based on the results of the EDL model. From this part, a current-voltage characteristic of the modelled device is obtained.

2.3.1. Bulk Electrolyte

Electric field is governed by the Poisson equation (3) [32]:

2 2 2 2 2 2

r 0

v ,

ε x y z

V 

 



 



 





 

 , (3)

where ε0 is the permittivity of vacuum (8.854 ∙ 10⁻¹² F m⁻¹), εr is the relative permittivity of the used material, ρv is space charge density (C m⁻³), x, y, z are Cartesian coordinates (m).

Space charge density in an electrolyte made up of N species is given by equation 4:





i N

1

i i

A

v ^qN



 , (4)

where q is the elementary charge (1.602 ∙ 10⁻¹⁹ C), NA is the Avogadro constant (6.022∙ 10²³ mol⁻¹), zi is the charge number of an ion, ci is its concentration (mol dm⁻³) [33].

In the bulk electrolyte, the condition of electroneutrality must apply. This means that the concentrations of cations and anions are balanced. Therefore z+c+ = z−c− and the space charge density becomes zero. In this case, equation (3) reduces to the Laplace equation (5):

2 0

2 2 2 2

2 



 



 





z y

x . (5)

2.3.2. Electric Double Layer

Electroneutrality applies in most of the electrolyte, except the EDL, where ion concentrations change significantly. Therefore, the potential is almost constant in most of the electrolyte and most changes occur close to the EDL, where the electric field is governed by the Poisson equation (3). Concentration changes in the EDL occur, because ions are attracted to an electrode with opposite charge. Electric charge on the electrode interacts with opposite ions and allows the EDL to form. An illustration of the EDL is shown in fig. 8 [33, 34].

20 Fig. 8 – the electric double layer, where the grey area is the electrode, the red line signifies electric potential changes in the respective layers. The x axis is plotted in arbitrary units (a. u.), however, the

size of the EDL is typically in the order of 1 nm.

The electric double layer consists of a compact or Stern Layer (SL) closest to the electrode.

In this layer, species (mostly solvent molecules) in this layer are generally localized on the surface (specifically adsorbed). This layer is expected to have a thickness of the ion radius. Another layer is the Diffuse Layer (DL), expected to have a thickness consistent with the Debye screening lenght λD (eq. 6):

 

2 R

NaCl 2 sol

D c

 T

  , (6)

where R is the universal gas constant (8.314 J mol⁻¹ K⁻¹), T is the temperature (K) and F is the Faraday constant (96,485.34 C mol⁻¹). The SL/DL boundary is called the Inner Helmholtz Plane. In the DL, solvated ions are distributed (said to be nonspecifically adsorbed). Even specifically adsorbed ions generally cannot pass the Inner Helmholtz Plane, and therefore there is no space charge distribution in the SL, so the Laplace equation applies. Furthermore, the permittivity of solvent in the SL is significantly smaller than in the bulk electrolyte and DL, due to surface localization [34, 35, 36].

Electric potential V (V)

x (a. u.) Diffuse layer

Stern layer

Bulk electrolyte

i

21 Electric field in the DL and bulk electrolyte is governed by the Poisson equation (3). In the SL it is governed by the Laplace equation (5). In addition, ion migration in the electric field takes place. This is governed by the Nernst-Planck (N-P) equation for each ion (7):

i i N



 t

c , i j k

z y

x 

 



 



 

 (7)

Where t is time (s), i, j, k are base unit vectors, and N is the ion flux given by equation 8 [36]:





i i i

i molm s

R

D  FD  ^{ }



 V

T c c z

N , (8)

where D is the diffusion coefficient (m² s⁻¹). The first term in equation 8 corresponds to diffusion; the second term corresponds to migration in an electric field. There is no need to account for convection in this case.

2.3.3. Charge Transfer in the Channel

The electric field inside the polymer layer is governed by Poisson's equation (3). In the case of PEDOT:PSS, the ECP structure is considered to be divided into the ionic and electronic phase. In the ionic phase, the space charge density is given by equation 9; while in the electronic phase, it is described by equation 10:



 



  





 N

1

i i i fix

v q p zc c

 , (9)

 





q p V V C

v   

 , (10)

Where p is hole concentration (mol m⁻³), Vp and Vc are the electric potentials in the electronic and ionic phase respectively (V), cfix is the concentration of fixed charges. It can be said, that the ionic phase corresponds to the PSS rich phase, while the electronic phase corresponds to the PEDOT rich phase, as seen on fig. 6.

Additionally, ion transfer in the layer is governed by the Nernst-Planck equation (7 and 8), while hole transfer is governed by the modified N−P equation 11:

d , , d

q R

D_p F _p ^p _p

p t

V p T p

p   









 



 



 









 j

j 

(11) where jp is the ionic flux (mol m⁻¹ s⁻¹), and μp is the chemical potential of holes (J mol⁻¹).

This model provides results, which match experimentally gained characteristics very well.

It is also possible to model the moving front of dedoped polymer, with time-dependent characteristics [22].

i j k

22 An older model developed by Bernards and Malliaras [37] considers changes to hole concentration as a function of total charge Q (C) of cations injected from the electrolyte into the polymer layer (equation 12):



 



 

 p v

p Q p

0

0 1 q , (12)

where v a volume of polymer material (m³). Negative ions are assumed to have no effect.

The cations are expected to be injected into the polymer due to the influence of a potential difference between the channel and G electrode. The charge Q is then a function of this difference (equation 13):

 



 



Q  _d d _G  , (13)

where cd is capacitance per unit area (F m⁻²) and W is the width of the channel (m).

Equation 13 defines charge Q at steady state. The time dependence of Q is given by equation 14:

   











 

 ^{ i}

t

x Q t

Q 1 e ^ , (14)

where τi is the ionic transit time (s), e is the Euler’s number (2.718). The value of p from eq. 12 can then be used to determine electric current density J (A m⁻²) in the CP from the Ohm's law (15):





 , V, qp

 E E

J , (15)

where E is the electric field intensity (V m⁻¹) and σ is the conductivity (S m⁻¹).

Advanced models of charge transport in the ECP layer are based on the density of states (DOS) distribution in the polymer and electrolyte. It has been shown, that charge carriers in higher energy states have a larger mobility than those at lower states, due to the exponential distribution of DOS with respect to energy – there are more states available at higher energies, between which charge carriers can move. When there are enough charge carriers to fill all lower energy states, the remaining carriers are forced to occupy higher energy states, where they are more mobile. Based on this, Friedlein et. al. [38] have proposed a relation between hole concentration and mobility (equation 16):

 

, k⁰

0

0   

 



 

T A E p

p p

A



 , (16)

where E0 is the energetic width of the DOS, k is the Boltzmann constant (1.380,648,52 ∙ 10⁻²³ m² kg s⁻² K⁻¹). They have further shown, that the parameter A changes minimally with the length/width ratio of the channel and for most geometries A = 2.

23 2.4. Finite Element Method

The finite element method is a method of numerical analysis, which provides approximate solutions of partial differential equations on some geometric domain. First, the domain is divided into a collection of relatively simple subdomains (such as triangles in 2D or tetrahedrons in 3D). These subdomains are called finite elements, and their collection is referred to as the finite element mesh. Points in the geometry, which are common for multiple elements, are called nodes.

Second, approximation functions are derived on the entire mesh. This part is based on the idea that any continuous function can be represented by a linear combination of algebraic functions. These functions are referred to as shape or basis functions. Shape functions are derived based on a special integro-differential form of the equation, often called the weak form. For example, the weak form of Poisson's equation (3) in a one-dimensional case is given by equation 17:





1 2

1

d 0 d d

d d

d

v x

x x

x x

w V dx w x dx

V x

w    , ε₀_r (17)

where w is an arbitrary test function. After obtaining the required set of shape functions, a matrix system is formed, consisting of a matrix of shape functions for each node, a vector of solutions and a vector, which corresponds to the imposed boundary conditions.

Third, relations between the nodes are computed, satisfying the governing differential equations, boundary conditions and initial conditions in the case of a time-dependent problem.

Simply put, the governing equation is solved across the nodes, and the solution is approximated everywhere in between.

The finite element method is a powerful tool in numerical analysis, as it allows for discretization of complex geometries into several types of elements, and inclusion of multiple materials. Errors in this method are created by the approximation of the domain and solution, and general computational errors such as numerical integration and round-off errors [32, 39].

2.4.1. COMSOL Multiphysics

COMSOL Multiphysics is a computational program, designed to solve problems in physics and physical chemistry. These problems are most often described by partial differential equations, but it is possible to solve ordinary differential equations and algebraic equations.

After the desired equations, boundary conditions and initial conditions are set, the system is numerically solved. The finite element method is the most frequently used, but it is possible to use others, such as the finite difference method or boundary element method [40].

Besides analysing a single physical field or function at a time, it is also possible to study the behaviour of multiple fields in the same system, such as an electric field and the distribution of ions in an electrolyte. This is done by coupling the respective differential equations. In the above example, this would mean coupling the Poisson equation (3) and Nernst-Planck equation (7, 8) by taking values of electric potential to compute changes in ionic concentration due to ionic drift, while in turn taking values of positive and negative ion

24 concentration to compute changes in the space charge density (4). This type of approach is referred to as multiphysics [41, 42].

The required equations can be put into COMSOL manually, or the user can use one of the many modules, which are available to buy from the creators of COMSOL. These modules offer pre-set systems of equations, used in the respective field of physics (such as the Maxwell equations in the case of the Electrostatics and Electric Currents modules), which allow more accurate analysis. It can be said, that COMSOL itself is a program, designed to solve systems of differential equations, with a user interface. Databases of physical equation must be purchased together with the base program, and each cover only one field of physics, such as electromagnetics, chemical engineering or structural mechanics [40].

Using the modules offers some simplification for the user, especially if they are not an expert in mathematical physics. However, it is not always clear which equations are being used in the model. Therefore, it is always necessary to understand the underlying physics and closely watch the results of the model.

COMSOL is further capable of postproduction, such as applying functions and transformations to the results [40]. While it is possible to create one- to three-dimensional plots in COMSOL, its primary function is not that of a graphics processor, and it is not easy to produce plots ready for publication.

25

3. EXPERIMENTAL PART

The present model consists of three major parts – the macroscopic bulk electrolyte model, the EDL model and the PEDOT:PSS channel model. It is necessary to operate these as largely independent components, due to the differences in scale (bulk electrolyte model size in the order of 1 mm, EDL model size in the order of 10-100 nm). While the finite element method is well prepared for modelling domains of different sizes, a scale difference of around 10⁻⁶ m is too large even for such a powerful tool. In fact, due to the low computational power of the hardware used (a standard personal computer), even the individual model components have to be further simplified, in order to be able to finish the modelling. Unfortunately, these simplifications are almost certain to result in lower quality of the results.

3.1. Bulk Electrolyte

The first step of modelling the OECT is to compute the electric field in the bulk electrolyte and inside the PEDOT:PSS layer between the S and D electrodes. A circular OECT geometry is considered for this model, based on a real OECT configuration, shown in fig. 9.

Dimensions of the OECT are based on the chimney-well microplates, used as housing for the OECTs. Using such microplates allows for arranging the OECTs into arrays, as shown in fig. 10. For simplification, a fully circular G electrode with the channel in the centre is considered, with a PEDOT:PSS channel of length L, width W and thickness H. When Cartesian coordinates are used, L corresponds to coordinate x, W to coordinate y and H to coordinate z. Such a system has axial symmetry, and therefore it is possible to analyse only one half of the geometry, as shown in fig. 11. The complete geometry can then be analysed by mirroring the result. Parametres used in the geometry of all models are summarized in table 1.

It is important, that the radius of the transistor (2.5 mm) and height of the electrolyte cell (10.9 mm) are many orders of magnitude larger than the thickness of the polymer layer (200 nm). Due to this, it is not possible to model the electrolyte cell and the polymer layer at the same time, since the finite element mesh would require element sizes around 100 nm, leading to billions of elements in the electrolyte cell area – far too many for data storage and computation.

Essentially, the model has to be split into two parts. First, the electric field in the electrolyte above the transistor geometry is modelled. The electrodes in this part are approximated as infinitely thin layers. The G electrode is expected to have uniform potential on its surface, and the S and D electrodes are connected by an infinitely thin layer of PEDOT:PSS; they are expected to create a potential gradient. In the real OECT, an insulator (Resist SU-8 in this case) is used to mask the PEDOT layer, to make sure the electrolyte does not come into contact with the electrodes. However, it does not play a significant role in describing the transistor function, apart from obstructing a small part of the electrode surface, and therefore it is not modelled here.

Second, the thin layer of PEDOT:PSS with the S and D electrodes is modelled. This part of the model has a very simple geometry, and the potential field gradient is expected to be the same in the entire channel. Therefore, this model is effectively a 1D case, which is trivial to compute in COMSOL and even has an analytical solution. This part is shown in fig. 12.

26 Materials considered in the whole model are PEDOT:PSS, silver paste and the electrolyte PBS (phosphate buffered saline). The main components of PBS are sodium and chlorine ions, which will be considered in this model. The other ions (K⁺ and phosphates) have a much lower concentration and are omitted for simplicity. Material constants used in all models are summarized in table 2. For the macroscopic model, only relative permittivity values need to be considered. The material properties of silver paste differ based on manufacturer and mixture, but since silver is only ever present in the form of an infinitely thin layer with a homogeneous potential distribution and current density in the electrode is not studied, there is no need to consider the values of relative permittivity and conductivity. The temperature was set to T = 298 K for the computation.

Fig. 9 – Top view of the real OECT configuration.

27 Fig. 10 – Schematic of an OECT array in chimney-well microplates.

Table 1 – parametres used in the geometry

Table 2 – Material constants used in the model.

Constant Value Description

εPEDOT 2.2 [19] PEDOT:PSS relative permittivity

εsol 78.3 [43] NaCl solution relative permittivity

p0

μ0

1 ∙ 10²¹ m⁻³ [19]

2 m² V⁻¹ s⁻¹ [26]

PEDOT:PSS hole concentration PEDOT:PSS hole mobility DNa 1.33 ∙ 10⁻⁹ m² s⁻¹ [44] Na⁺ diffusion coefficient DCl 2.03 ∙ 10⁻⁹ m² s⁻¹ [44] Cl⁻ diffusion coefficient aNa

cNa

cCl

C^*

0.72 ∙ 10⁻⁹ m [45]

137 mol m⁻³ [46]

137 mol m⁻³ [46]

40 F cm⁻³ [21]

Na⁺ ion diameter Na⁺ ion concentration Cl⁻ ion concentration

PEDOT:PSS volumetric capacitance

Parameter Value (mm) Description

h 10.9 Chimney-well microplate height

d 2.5 Chimney-well microplate radius

WG 0.3 G electrode width

L W

0.1 1

Channel length Channel width Lel

H λD

λS

0.1 200 ∙ 10⁻⁶ 8.2 ∙ 10⁻⁶ 3.6 ∙ 10⁻⁶

S and D el. Length

PEDOT:PSS layer thickness Debye length

Stern layer thickness

28 Fig. 11 – the macroscopic model geometry from the side (top) and from the top view (bottom).

The channel and electrodes are shown with a non-zero height for clarity, but they are actually modelled as infinitely thin layers. Sizes are in millimetres. The height of the model is not to scale.

2.5

10.9

0.3

0.1 0.1

0.3

G D S G

Electrolyte PEDOT:PSS

0.5

29 Fig. 12 – An illustration (not to scale) of the electric field modelled in the PEDOT:PSS layer as a part

of the macroscopic model. Due to the very low thickness of the layer in relation to its length, this model reduces to a one-dimensional case. The red line signifies expected potential changes.

Due to the condition of electro neutrality in this case, that is cNa = cCl and therefore ρv = 0 C m^–3, the electric field is governed by the Laplace equation (5). The potential VD was set to vary between − 1 V and 0 V with a step of 0.2 V. The potential VG was set to vary between – 1 V and 1 V, with a step of 0.2 V. The S electrode was set to potential VS = 0 V.

The results of the macroscopic models can later be used in the EDL model. More specifically, the potential on top of the channel, computed from the bulk electrolyte model, denoted Ve, and the potential inside the channel, denoted Vm will be used.

3.2. Electric Double Layer

This part of the model aims to compute the electric potential on top of the ECP layer, and the change of electric charge inside the PEDOT:PSS layer. Afterwards, the new hole concentration and the mobility can be computed. In this part of the model, it is necessary to work with the Nernst-Planck equation. This is problematic, because computing the diffusion and the ionic drift requires element sizes of 1 nm or less [22, 47]. Therefore, the computer memory requirements for analysing a 3D electrolyte domain of dimensions in the order of 1 mm are extremely high. Far higher, than what is reasonable to expect from the personal computer (PC) used to run this model. However, the EDL may be reduced to a 1D case, at a significant cost to accuracy. Taking potential values at multiple points along the channel geometry and computing the 1D EDL for each of them allows for creating an approximation of the electric potential and ion concentration distribution on the top of the channel.

The area on the top of the channel is effectively split into a certain amount of subdomains, each of them having its own EDL. An illustration of this division is shown in fig. 13. It is possible to increase accuracy of this approximation, at the cost of increasing computation time, by adding multiple points, at which to analyse the EDL.

V (V)

Channel length x (mm) VD

0

0 0.1

D S

30 After each of the EDL models is computed, functions Vi = f (z), which describe the changes of electric potential and functions cNa = f (z), describing ion concentration in the EDL, are gained for given pairs of electric potential values Ve and Vm and the initial electrolyte concentration, as shown in fig. 14. These functions correspond to potential drops at the respective points xi as seen on fig. 13. The potential distribution on the top of the channel, V = f (x), and the change of ion concentration close to the channel, can be constructed from these functions. An illustration of computing the electric potential across the channel is shown in fig. 14. The ion concentration change ΔcNa from the base value of 137 mol m⁻³ is evaluated from the function cNa = f (z) for an element of volume λs ∙ H ∙ Lx (see fig. 14) on the DL-SL interface, effectively approximating concentration changes in the element of volume by a single number.

The electric field is governed by the Poisson equation (3) in the DL, where significant changes in concentration can be expected. In the channel area, the SL, and the bulk electrolyte, the electric field is governed by the Laplace equation (5). The concentration of each ion in the diffuse layer and bulk electrolyte is governed by the Nernst-Planck equation (7, 8). The equations are coupled to each other by the relation for space charge density (4) and the ionic drift term as described in paragraph 2.4.1.

Fig. 13 – division of the polymer channel into subdomains for the purpose of EDL modelling. From each subdomain, a value of V_e and V_m is evaluated at a chosen point x_i. This evaluation then serves to approximate the EDL as a 1D case for the subdomain of length L_x. The blue line signifies potential on

top of the channel.

VD

Channel length x (mm)

Electric potential V (V)

0

0 0.1

Vm1

Vm2

Vm3

Vm4

Vm5

Vm6

Ve6

Ve5

Ve4

Ve3

Ve2

Ve1

x1 Lx

EDL evaluation

31 Boundary conditions for electric potential are taken from the previous part of the model.

The potential Ve, corresponds to the influence of the G electrode on the EDL. The potential Vm corresponds to the potential inside the channel. As initial concentration, c0 = cNa = cCl was set. The solvent in the SL is expected to have a relative permittivity of 8 due to heavy polarization of molecules [34, 35, 36].

Fig. 14 – An illustration of computing the potential changes on the top of the channel. The individual potential function, gained from modelling the EDL are shown in red, yellow and dark green along the

z axis. The light grey planes signify the boundaries of the Stern and Diffuse layers. The green plane signifies the channel boundary. As shown on the green curve, each EDL is assigned a point on the x axis. The computed potential value on the top of the channel is then extracted from this curve. From the potential values and assigned x coordinates, a function V = f (x) can be plotted in the x-V plane, either by interpolation or approximation by an analytical function. This function is represented by the

dark blue line. The x axis is plotted in arbitrary units (a.u.), however, the channel length x is on the order of 1 mm in a typical OECT .The z axis plotted in nm. The y axis on the right is plotted in a.u. and corresponds to the width of the channel. The light blue line signifies the EDL area; the light green line

signifies the channel area. The element of volume λs ∙ W ∙ Lx is shown by the light red lines.

It should be noted, that the 1D EDL model does not exactly represent reality, because it fixes the potentials Ve and Vm at a chosen point in the channel and electrolyte, effectively treating them like fixed metal electrodes. In the real OECT, the potential values are only fixed on the electrodes, not at any chosen point xi. Unfortunately, there is no other way of computing this 1D model than simply fixing the potential values by setting them as boundary conditions.

Electric potential V (V)

x (a.u.) z (nm)

Diffuse layer Stern layer

Bulk electrolyte V1 = f (z)

V2 = f (z)

V3 = f (z)

V = f (x)

x3

y (a. u.)

λs

W Lx

0.00 3.60 8.21