# PDE Concepts I need a PDF with all the answers. I’m providing the text you will follow to answer the questions. JCAMECH

Vol. 51, No. 1, June 2020, pp 21-2

PDE Concepts I need a PDF with all the answers. I’m providing the text you will follow to answer the questions. JCAMECH

Vol. 51, No. 1, June 2020, pp 21-29

DOI: 10.22059/jcamech.2019.282591.406

A theoretical model for analysis of ionic polymer metal composite

sensors in fluid environments

Mohammad Reza Salehi Kolahia, Hossein Moeinkhaha*

aDepartment of Mechanical Engineering, University of Sistan and Baluchestan, Zahedan, Iran

———

*

Corresponding author. Tel.: +54-31136447; e-mail: hmoein@eng.usb.ac.ir

1. Introduction

In recent years ionic polymer metal composites (IPMCs),

which are a unique and new class of smart materials, have been

studied by so many theorists and practitioners for their great

potential applications in sensing and actuation [1-7]. IPMC

sensors are very sensitive to low mechanical stimulus which is an

advantage over some piezoelectric sensors. Flexibility, low

power consumption, large stroke, excellent biocompatibility and

operation in wet environments are some other advantages of the

IPMCs [8, 9]. As illustrated in Figure 1, an IPMC is made of a

soft ionic polymer membrane sandwiched between thin noble

metal electrode layers [10]. Typical polymer membrane materials

are Nafion )Perfluorosulfonate) and Flemion

(Perfluorocarboxylate) [11] and for the electrodes with high

conductivity specification, metals such as platinum, silver and

gold can be used [12]. According to the requirements, these

groups of smart materials can be made in variety of shapes and

sizes [13, 14]. IPMCs are capable of converting energy between

chemical, electrical and mechanical domains. Applying a voltage

between the conductive electrodes, causes the free cations and

attracted water molecules to migrate from the anode to the

cathode side. Due to this mechanism, bending motions of the

IPMC occur. The polarity of the applied voltage determines the

bending direction of the IPMC. Noting that the anions are fixed

to the polymer backbone and are not free to move [15, 16]. On

the other hand, an applied deformation on an IPMC strip causes

redistribution of the free cations inside the polymer and produces

a short circuit current across the electrodes [17]. Thus the main

cause of both of the actuation and sensing phenomena is induced

ionic current.

To study and describe the behavior of IPMCs, up till now

several models have been proposed and applied. We can classify

these models into three categories as [18]: 1- the black box

model, 2- the gray box model and 3- the physical model. The

most complex forms among them are physical models, in which

the ionic current through the polymer membrane is related to the

mechanical deformation [19].

Comparing with the intensive works which have been done on

modeling of IPMC actuators, limited researches have been done

on IPMC sensors. Based on linear irreversible thermodynamics

De Gennes et al. [20] proposed a static model to capture both

actuation and sensing mechanisms of IPMCs. Farinholt and Leo

[21] derived the charge sensing response for a cantilevered IPMC

beam under a step change in tip displacement. Chen et al. [22]

developed a physics-based, geometrically scalable model for

IPMC sensors. Their mechanical model was simple and did not

include damping effects. They also investigated structural health

monitoring capabilities of IPMC sensors. Aureli et al. [23]

experimentally and theoretically investigated energy harvesting

capabilities of base excited IPMC sensors in fluid environments.

Their theoretical modeling approach was complicated. In 2011

A R T I C L E I N F O A B S T R A C T

Article history:

Received: 31 May 2019

Accepted: 28 June 2019

By the past two decades IPMCs have been intensively studied because of their

special capabilities for actuation and sensing. This paper presents a theoretical

physics based model for analyzing the behavior of IPMC sensors in fluid

environments. The mechanical vibration of the IPMC strip is described by the

classical Euler–Bernoulli beam theory. The model also takes in to account the

physical properties of the surrounding fluid. The resulting model is an infinite-

dimensional transfer function that relates the input tip displacement to the output

sensing current. Further the original model is reduced to a finite-dimensional one,

for pure sensing applications of IPMC sensors such as structural health

monitoring. The proposed model is verified using existing experimental data.

Then the effect of various parameters is investigated. The acoustics physics

interface in COMSOL Multiphysics software is used for coupled modal analysis

of the IPMC strip. It is shown that the effect of surrounding fluid cannot be

neglected.

Keywords:

Dynamic model

IPMC

Physical model

Smart materials

Mohammad Reza Salehi Kolahi and Hossein Moeinkhah

22

Ganly et al. [24] presented the first systematic studies on

temperature-dependent IPMC sensing dynamics. They used a

gray box model in form of a forth order transfer function. Lei and

Tan [25] developed a dynamic model for a base-excited IPMC

sensor operating in air. Using a mathematical model, an

estimation of the power harvested by the multilayered cantilever

IPMC sensors from the base excitation is reported in [26]. Most

of these papers are investigating the features and possible

applications of IPMC sensors instead of proposing a complete

model.

To the best of our prior studies and knowledge, there is no

complete mechanical model for the IPMC sensors. In this

research a dynamic, physics based model for analyzing the

behavior of the IPMC sensors is presented. This work extends

previous studies as it incorporates the effect of the surrounding

fluid on the IPMC vibration. The frequency response of an elastic

structure is very sensitive to the nature of the fluid in which it is

submerged. The deformation of the IPMC is modeled by Euler-

Bernoulli beam theory, incorporating viscoelastic damping of the

polymer membrane and the surrounding fluid damping force. The

description of ion transport dynamics is based on the governing

partial differential equation in [22]. Since the resulting model is

an infinite-dimensional transfer function, is further reduced to a

finite-dimensional one, using some approximations and Taylor

series. The resulting model is validated based on the existing

experimental data. COMSOL Multiphysics software is used for

coupled modal analysis of the IPMC sensor.

The current paper is organized as follows. In section 2, the

mathematical model for the IPMC sensor is derived. Model

reduction is discussed in section 3. The simulation results are

discussed in Section 4. Finally, the conclusions are provided in

Section 5.

Figure 1. Schematic model of an IPMC strip.

2. Model derivation

Here we consider a slim cantilevered IPMC with rectangular

cross section of thickness 2h, length L, and width b. The neutral

axis of the beam is shown by z = 0, see Figure 2. In order to

obtain the dynamic equation of motion we consider the effects

that two types of damping have on the frequency response of an

elastic structure: viscoelastic damping and damping due to the

submergence of the structure in a viscous fluid. The main

difference between polymers and other materials resides in the

viscoelastic properties of the polymers. Viscoelastic material can

be defined as the one whose stress depends on deformation rate

[27]. The mechanical response of viscoelastic materials to

mechanical excitation has traditionally been modeled in terms of

elastic and viscose components such as springs and dashpots.

Here we use kelvin-Voigt model. This model combines a spring

and a dashpot in parallel, see Figure 3 [28].

v

xx

xx xxY C

t

e

s e

¶æ ö

= +ç ÷

¶è ø

(1)

here Cv denotes the strain rate damping coefficient and Y denotes

the Young’s modulus.

Figure 2. Schematic model of a cantilever IPMC sensor.

Figure 3. Schematic representation of the Kelvin-Voigt model.

The effect of the damping can be neglected in air, but it cannot

be neglected when the cantilever is immersed in a denser

medium, like water. The viscous fluid damping force per unit

length of the beam can be expressed as [8, 29]

( ) ( )

( )

( )

( )

( ) ( )( )

( ) ( )( )

2

fluid a 2

2

a rect

2

rect

, ,

,

Re

4

Im

4

f

f

w x t w x t

f x t m C

tt

b

m

b

C

w w

r p

w w

r p

w w

¶ ¶

= +

¶¶

= G

= – G

(2)

were ω (rad/s) is the oscillation, w(x, t) is the transverse

displacement, b is the width of the IPMC beam, ρf is the density

of the fluid, and Γrect(ω) is the hydrodynamic function, ma(ω) is

the added mass and C(ω) is the viscous fluid damping

coefficient. The hydrodynamic function for a beam with

rectangular cross section can be expressed as

( ) ( )

( )( )

( )( )

1

rect corr

0

4

1

e

e e

jK j jR

jR K j jR

w y w

æ ö-

ç ÷

G = +

ç ÷

-ç ÷

è ø

(3)

Journal of Computational Applied Mechanics, Vol. 51, No. 1, June 2020

23

where K0 and K1 are modified Bessel functions of the second

type, ψcorr(ω) is a complex valued correction function that

corrects the results for a beam of circular cross section to a beam

of rectangular cross section and

2

4

f

e

f

b

R

r w

m

= (4)

where Re is the Reynolds number and μf is the viscosity of fluid.

The correction function for the rectangular cross section is given

by [30]

( ) ( )corr r ijy y w y w= + (5)

where

( )

( )

r

r

r

i

i

i

a

b

a

b

y w

y w

=

=

(6)

and

2

3 4

5 6

2

3 4

5 6

2 3

0.46842

0.12866 0.044055

0.0035117 0.00069085

1 0.56964 0.48690

0.13444 0.045155

0.0035862 0.00069085

0.024134 0.029256

0.016294 0.00010961

0.00

0.91324 0.48 74

00645

2r

r

i

a

b

a

l

l l

l l

l l

l l

l l

l

l l

l

– +

– +

– +

– +

= – +

=

=

– +

– –

+ –

+

4 5

2

3 4

5 6

77 0.00004451

1 0.59702 0.55182

0.18375 0.079156

0.014369 0.0028361 .

ib

l l

l l

l l

l l

=

–

– +

– +

– +

(7)

The quantity λ is defined as

10log eRl = (8)

Expression (3) of the hydrodynamic function can be used for a

wide range of Reynolds numbers [30]. Considering these

assumptions, the governing equation of motion is obtained as

[31]

( )

( )

( )

( )

( ) ( )

2

2

5 4

4 4

, ,

, ,

0v

w x t w x t

C

xx

w x t w x t

C I Y I

x t x

h w w

¶ ¶

+

¶¶

¶ ¶

+ + =

¶ ¶ ¶

(9)

where η(ω) = ρmA + ma(ω), ρm is the density, A is the cross

sectional area of IPMC and I is the area moment of inertia with

respect to the x axis. By converting equation (9) into the Laplace

domain we have

( )

( )

( ) ( )( ) ( )

4

2

4

,

, 0v

w x s

Y C I C s w x s

x

w h w

¶

+ + + =

¶

(10)

where s is the Laplace variable. Equation (10) can be rewritten as

( )

( ) ( )

4

4 2

4

,

4 , 0

w x s

s s w x s

x

g

¶

+ =

¶

(11)

where γ4(s) is

( )

( ) ( )4

4 ( )v

C s

s

Is Y C s

w h w

g

+

=

+

(12)

The general solution for equation (11) is as follows [32]

( ) ( ) ( ) ( )

( ) ( ) ( )

1

2

, cos cosh

cos sinh

w x s A s px px

A s px px

=

+

(13)

( ) ( ) ( ) ( ) ( ) ( )3 4sin cosh sin sinhA s px px A s px px+ +

where

p sg= (14)

For the clamped-free boundary conditions, the following

relations can be obtained.

( )

( )

( ) ( )

( )2

2

0,

0, 0 , 0

,

, , 0

w s

w s

x

w L s

w L s W s

x

¶

= =

¶

¶

= =

¶

(15)

By substituting the boundary conditions into (13), we get

( )

( )

( )

( )

,

,

A x s

w x s W s

B s

= (16)

where

( ) ( ) ( ) ( ) ( )

( ) ( ) ( )( )

, cos sinh sin cosh

sin sinh

A x s px px px px

C s px px

= –

+

( ) ( ) ( ) ( ) ( )

( ) ( )

( ) ( )

( ) ( ) ( ) ( )( )

cos sinh sin cosh

sin sinh

cos cosh

sin cosh cos sinh

B s pL pL pL pL

pL pL

pL pL

pL pL pL pL

= –

+

´ +

( )

( ) ( ) ( ) ( )

( ) ( )

sin cosh cos sinh

cos cosh

pL pL pL pL

C s

pL pL

+

=

(17)

As illustrated in Figure 2, an applied displacement w(t) leads

to the generation of the sensing current i(t) due to redistribution

of the free cations. The governing PDE which describes the

charge density distribution ρ within the polymer membrane is

given by [22]:

( ) ( )

( ) ( )

2

2

2

0

0

, , , ,

1 , , 0

x z t x z t

d

t x

F dC

C V x z t

RT

r r

r

e

¶ ¶

–

¶ ¶

+ – D =

(18)

where d is the ionic diffusion, F is Faraday’s constant, R is the

gas constant, T is the absolute temperature, ε is the effective

dielectric constant of the polymer membrane, C0 is the anion

Mohammad Reza Salehi Kolahi and Hossein Moeinkhah

24

concentration, and ΔV is the molar volume change. By taking the

Laplace transform for the of ρ(x, z, t), we can convert (18) into

the frequency domain as

( )

( ) ( )

2

2

2

, ,

, , 0

x z s

s x z s

x

r

b r

¶

+ =

¶

(19)

where β2(s) is

( ) ( )

2

2 0

0, 1

F dCs k

s k C V

d RT

b

e

+

= = – D (20)

Assuming that the charge density is symmetric relative to z =

0, the solution of (19) is obtained as

( ) ( ) ( )( )1, , 2 , sinhx z s a x s s zr b= (21)

The following equations describe the electric field E and the

electric potential φ.

( )

( ), z,s

, ,

x

E x z s

x

j¶

= –

¶

(22)

( )

( )

, ,

, ,

E x z s

x z s

x

e r

¶

=

¶

(23)

With the above equations and equation (21), we can evaluate

E(x, z, s) and φ(x, z, s) as

( ) ( )

( )( )

( )

( )

1

1

cosh

, , 2 ,

,

s z

E x z s a x s

s

b x s

b

eb

=

+

(24)

( ) ( )

( )( )

( )

( ) ( )

1 2

1 2

sinh

, , 2 ,

, ,

s z

x z s a x s

s

b x s z b x s

b

j

eb

= –

– +

(25)

The charge density ρ(±h, z, s) inside the polymer

membrane at the boundary z = ±h is proportional to the

induced stress σ(±h, z, s) [33]:

( ) ( ), , , ,x h s x h ss ar± = ± (26)

where α is the charge–stress coupling constant. Now we can

relate the induced stress σ(x, h, s) to the external stimulus W(s).

According to linear beam theory we have

( )

( )

( )

( )2

2

, ,

, z, , ,

M x s z w x s

x s M x s Y I

I x

s

¶

= =

¶

(27)

where M is the bending moment. Substituting the equation (16)

in to (27) one can obtain the stress at z = h

( )

( )

( )

( )2

2

,

, ,

Y hW s A x s

x h s

B s x

s

¶

=

¶

(28)

Using equation (26) we can get a1(x, s)

( )

( )

( ) ( )( )

( )2

1 2

,

,

2 sinh

Y hW s A x s

a x s

B s s h xa b

¶

=

¶

(29)

Now we can obtain b1(x, s) and b2(x, s) using the short circuit

boundary condition φ(x, h, s) ─ φ( x, -h, s) = 0.

( ) ( )

( )( )

( )

( )1 1 22

sinh

, 2 , , , 0

s h

b x s a x s b x s

h s

b

e b

= – = (30)

By integrating the electric displacement field D = εE over the

cross section of the beam at the z = h the total induced electric

charge can be obtained as

( ) ( )

0 0

, ,

b L

Q s E x h s dxdye= ò ò (31)

The short-circuit current i(t) in the Laplace domain is i(s) =

sQ(s). Finally we can achieve to the transfer function which

relates the mechanical input w(s) to the output sensing current

i(s) as

( )

( )

( )

( ) ( )( )( )

( )

( )

( )

2

2

2

0

Ybs coth 1

,1

L

s h s hI s

H s

W s s

A x s

B s x

b b

ab

–

= =

¶

´

¶

ò

(32)

where

( )

( ) ( )

( ) ( ) ( ) ( ) ( )( )

2

2

0

,

2 sin sinh

cos sinh sin cos

L

A x s

p pL pL

x

p C s pL pL pL pL

¶

= –

¶

+ +

ò

(33)

3. Reduced order model

The resulting model in equation (33), is infinite-dimensional

as it contains nonrational terms like sin(○), cos(○) and , etc.

For practical implementation of the model, such as structural

health monitoring, we aim to reduce the model to a finite order.

First, we take 1 – C0ΔV ≈ 1 since |C0ΔV | 1. Based on the data

listed in Table 1, for s = jω, the term

k

d is of the order of 1012,

which means

( ) 610

s k

s

d

b

+

= ³ (34)

On the other hand, due to the thickness of the sensor it can be

concluded that

( ) 10s hb (35)

Equation (34) allows us to consider the following

approximation

( )

( ) ( )( )( )

( )

( )( )

( )

( )

( )

1 2

2

coth 1

1

Y bs s h s h

H s

s

Y bs d s k dY bs s

s ks

b b

ab

b

aab

–

=

+ —

= =

+

(36)

Further we can approximate the term s K+ using its Taylor

series about s = 0. The following approximation can be achieved

by considering up to the second order.

Journal of Computational Applied Mechanics, Vol. 51, No. 1, June 2020

25

( )

( )

1

1

s

Y bs d h k d

k

H s

s ka

æ öæ ö

+ -ç ÷ç ÷

è øè ø

=

+

(37)

To simplify the second part of (32), Taylor series expansion

for terms sin(○), cos(○), sinh(○) and cosh(○) will be used.

Considering just the first three terms in each series leads to

equation (39).

( )

( ) ( ) ( )(

( ) ( ) ( )(

( ) )

( ) )

16 12 8

2 16 12 8

4

4

5 128 384

56 10560

36864 8294400

576000 27648000

pL pL pL

H s

L pL pL pL

pL

pL

– –

=

– –

+ +

+ +

(38)

The resulting finite-dimensional approximation is

( ) ( ) ( )1 2H s H s H s= × (39)

4. Results and discussion

In order to verify the obtained model, parameters in the model

need to be identified. Some parameters in the model are physical

constants (gas constant R and Faraday constant F) and some of

them such as sensor dimensions, temperature T, density ρm and

young modulus Y are directly measurable. The other parameters

must be identified through a curve fitting process in MATLAB

[8]. The MATLAB function fminsearch can be used to find the

remaining parameters. This command minimizes the squared

error between the experimental frequency response and the

proposed model prediction. Table 1 shows the proposed model

parameters. The frequency response of the proposed model is

obtained for an IPMC sensor with the same dimensions, 22 × 7 ×

0.360 mm3, as used by Chen et al. [22]. By substituting s = j2πf

in equation (32) and using MATLAB the Bode plot of the

frequency response can be obtained. Figure 4 shows the Bode

plot of the frequency response. As illustrated in Figure 4, there is

a good correlation between the proposed model and experimental

results. This indicates that the proposed model is effective in

capturing the sensing of the IPMC.

It is worth noting that Chen et al. [22] conducted their

experiments in air. As the proposed model is a physical model it

is geometrically scalable. Thus we can study the effects of

geometrical parameters. Figure 5 shows the effect of length

variation on the sensor’s behavior. The length variation affects

only the magnitude and has not any effect on the phase. It also

can be seen that as the strip gets longer the magnitude decreases.

It means that the output sensing current decreases too. Next we

investigate the effect of hydration level on the IPMC sensor

behavior. Humidity variation leads to changes in diffusion

coefficient d. Figure 6 shows and compares the predicted

frequency response for three different diffusion coefficients d.

Table 1. Model parameters

C0 1091 (mol/m3)

Cv 20000(Pa.s)

d 1.2⨯10-11 (m2/s)

F 96487 (C/mol) [22]

R 8.3143 (J/mol.K) [22]

T 300 (K) [22]

Y 571 (MPa) [22]

ε 2.2 (mF/m)

ρm 9159.3 (kg/m3) [22]

α 109 (J/C)

Figure 4. Comparison of the frequency response with experimental data.

Figure 5. The effect of length variation.

When the IPMC sensor is oprating in air, it faces with water

evaporation and its hydration level decreases and the diffusion

coefficient decreases too. Figure 7 shows the cyclic current–

displacement plot. This plot corresponds to a harmonic

mechanical tip-displacement with peak of 1 mm and frequency of

10 Hz. The overall behavior of the sensor in figure 7 is the same

as figure 6, as the diffusion coefficient decreases the output

sensing current gets smaller.

Mohammad Reza Salehi Kolahi and Hossein Moeinkhah

26

Figure 6. The effect of diffusion coefficient variation.

Figure 7. Current-displacement cyclic diagram.

The natural frequencies of a cantilever beam submerged in a

viscous fluid, can be obtained by the following equation [34]

( )

1

2

f

vac rect vacfluid

m

1

4

n n nb

h

p r

w w w

r

–

æ ö

= + Gç ÷

è ø

(40)

where ωvac is the natural frequency in vacuum. As the IPMCs

can be used in wet conditions and are biocompatible we

investigate and compare their behavior in water and human blood

as well as air. The first coupled natural frequency of the IPMC

sensor submerged in air, deionized water and human blood is

given in Table 2. Noting that the boundary conditions,

geometrical and physical parameters are the same as figure 4.

The values are captured from COMSOL Multiphysics software

and equation (40). Table 3 summarizes the characteristics of the

mentioned fluids.

Table 2. Evaluation of the first coupled natural frequency (Hz)

Fluid

Analytical

Eq. (40)

FE

Experimental

[22]

Air 29.975 32.013 ~30

Deionized

water

18.3861 21.00 ‒

Human blood 17.6907 20.282 ‒

Table 3. Fluid properties of Table 2

Fluid ρf (kg.m

-3) μf (Pa.s)

Air 1 0.0001

Deionized water 997 0.0008

Human blood 1125 0.004

The acoustics physics interface in COMSOL Multiphysics

software can be used to compute the coupled natural frequencies

of an elastic structure submerged in a fluid. We used a three

dimensional geometry for simulation. It is worth noting that the

used mesh is free triangular with 78154 elements. We selected

this type and number of mesh after several simulations.

According to the Table 2, there is a fine correlation between

analytical and FE simulation results. In particular with denser

fluid, the natural frequency is lower. The changes in the natural

frequency values enable us to realize resonant sensors. Now we

investigate the effect of the surrounding fluid on the sensors

frequency response. Figure 8 shows and compares the frequency

response of the sensor submerged in three different mediums. It

is clear that the effect of the surrounding fluid cannot be

neglected specially for the vibrating frequencies more than 1 Hz.

Figure 8. Comparison of the frequency response for three different mediums.

For some sensing applications like structural health

monitoring (SHM), our aim is to reconstruct the original

mechanical signal u(t) based on the output sensing current i(t).

The reconstructed mechanical signal can be achieved by

inverting the reduced order model.

( )

( ) ( )

inv

1W s

H

I s H s

= = (41)

Journal of Computational Applied Mechanics, Vol. 51, No. 1, June 2020

27

It shows also the wide applicability of the proposed model.

Figures 9 to 11 show the reconstructed signal for three different

sensing currents. Figure 9 shows a smooth step excitation leads

to a sensing current which rises to a peak value and then

decreases to zero. For figures 10 and 11 reconstructed signal and

sensing current are similar. It means a decay oscillatory

excitation leads to a decay oscillatory current and a multitone

oscillatory excitation leads to a multitone oscillatory current.

5. Conclusions

In current paper a theoretical model which is physics-based, is

proposed to simulate the mechano-electrical response of IPMC

sensors. Moreover, the model accommodates the effect of the

surrounding fluid on the IPMC sensing. Since the original model

is an infinite-dimensional transfer function, it is not suitable for

some applications. Further the order of the original model is

reduced and converted to a rational transfer function. It is shown

that the resulting model has a fine correlation with the existing

experimental data.

More the effect of various parameters investigated. The results

showed that as the sensor gets longer, the output sensing current

gets smaller. Results also depicted that decreasing the hydration

level leads to reduction of output sensing current. It has also been

observed that the effect of the surrounding fluid cannot be

neglected. We used acoustics physics interface in COMSOL

Multiphysics software for modal analysis. Results showed that

the IPMCs natural frequency differs from one medium to

another. This change enables one to use IPMC sensors as

resonant sensors. On the other hand the reduced order model

enables one to use IPMC sensors for structural health monitoring,

especially for underwater structures.

The …