Full Analytical 3D Ion Transport Model for Large Deformable IPMC Soft Actuator
1 Abstract -- Ionic Polymer Metal Composite is a well-known soft electroactive polymer composite that it's promising features tell us that it has adequate potential to be a utilizable and applicative soft actuator in the practical applications, especially in the small size applications. But this smart material is still immature, and one of the reasons that lead to its immaturity is lack of a valid and proper physics-based model for large deformation situations. In practical online and large deformation applications, the inverse non-autoregressive identification based models are the proper choices but if we want to know how IPMC works and what is the physics of its behavior in the large deformation situation the inverse identification based models are deeply blind, and we have to use physical and multi-physical approaches. It is our main aim in this paper, and for the first time, we want to present a fully analytical and physics-based ion transport 3D and non-Linear model for large deformable IPMC. In this direction, based on three dimensional Nernst-Plank PDE we will find a well-defined and valid relationship between input voltage and output tip displacement of IPMC for large deformation situation and with four provable pieces of evidence it will be proven that proposed model has chosen a proper way and it is more complete than previous benchmark and well-known physics-based models for IPMC, and also it is valid and accurate enough for large deformation modeling of IPMC.
Index Terms —IPMC, Physics-based model, Nernst-Plank PDE, Large deformation, Ion Transport. N OMENCLATURES h Half the thickness of the Nafion membrane W Width of the Nafion membrane L Length of the Nafion membrane t h The thickness of the Pt electrode D (x,y,z,t) Electrical displacement E (x,y,z,t) Electric field ( ) x, y, z,t Electric potential ( ) x, y, z,t Electric charge density e The effective dielectric constant of Nafion F Faraday constant C + (x,y,z,t) C − J (x,y,z,t) Concentration of cations Concentration of anions Spatial vector of the cations d Ionic diffusivity R Gas constant T Absolute temperature p Fluid pressure v The free solvent velocity field V Volumetric change of the Nafion k Hydraulic permeability coefficient V I (t) Input applied voltage to IPMC ( ) σ x,y,z,t
Mechanical stress Constant coupling coefficient M (x,y,z,t) Global bending moment r Position Vector M zy (z,t) Effective bending moment for IPMC Y Young's modulus of Nafion I Moment of inertia ( ) , z t Curvature shape function of IPMC ( ) , Math z t The curvature of ( ) , z t (Mathematical definition) ( ) , Mech z t The curvature of ( ) , z t (Mechanical definition) ( ) Tip t Tip displacement function of IPMC I. I NTRODUCTION oft-robotics is a new promising subfield of robotics that its main target is constructing robots using soft materials in order to make more flexible, more bio-inspired, and safer robots. One of the serious choices as a soft actuator is Ionic polymer-metal composites (IPMCs). Due to particular properties of IPMC such as very large stimulus strain, low density, high toughness and lightness [1-7], this soft actuator is more suitable for micro and small size bioapplications. IPMC has a thin polymer membrane (usually Nafion) that is coated by two noble metal electrodes (usually Pt) on both sides (Anod a cathode). When we apply a low voltage (less than 5 V) to IPMC, it will be bent toward the anode side. Inversely, a low voltage will be generated between the two electrodes in response to bending of IPMC [5]. Hence we can say that IPMC is an actuator as well as sensors [1-5, 8-11]. We can find that in the vast variety of application we need largely bendable IPMCs but IPMS’s behavior in large deformation situations is nonlinear [5, 12] and besides that as it has been mentioned in the [5, 12, 13], acquiring output feedback is not feasible in dominant practical applications of IPMCs. Hence it is a sensible expectation that we should choose a non-autoregressive (non-feedback) and large bendable approach in the modeling of IPMC, and we will choose and propose it in this paper. Generally, we can divide IPMC’s modeling approaches into two groups, the first group is analytical modeling, and the second one is predictive identification methods. Analytical models like using distributed Resistor-Capacitor (RC) equivalent circuits [13-19]. Or physical and multi-physical approaches using PDEs, Finite element methods, COMSOL modeling, etc [20-23]. And also there are some other miscellaneous uncategorized analytical models like [24-26]. In the predictive identification paradigms also we have two approaches, classical and intelligent methods. In classical methods, all papers except [5] have used autoregressive structures [27-30]. In intelligent methods also, all are autoregressive [2, 3, 31, 32] and except [2] and [3] are not valid for large deformation situation. Hence we can find that all mentioned methods except [5] in two above main categories are invalid for large deformation situations or they
Full Analytical 3D Ion Transport Model for Large Deformable IPMC Soft Actuator
Mohsen Annabestani, Nadia Naghavi, and Mohammad Maymandi Nejad S presented an improved version of Tadokoro model. In these three groups, NP models propose the most straightforward way to explain ion transport [33] that this approach is the most compatible approach with the nature of operation of ionic content into IPMC. The first time NP PDE was presented for charge distribution in an IPMC by Nemat-Nasser [9] and then improved his model in 2002[8], after that Farinholt [41, 42] used this method for representation of actuation and sensing response of IPMC and also for Modeling the electrical impedance response of IPMC [43]. This approach has been frequently used by other researchers, for example, Chen [44, 45] and Choonghee [46] have used NP PDEs for modeling of the electroactive deformation and sensing behavior of IPMC. NP PDEs have also been widely used by many others [24, 47-54] in it's modified, simplified, or incomplete version to show the process of cation migration in IPMC sensors and actuators. As mentioned before NP PDEs based models can propose the most straightforward solution to explain ion transport behavior of IPMC, we want to use and solve this PDE to model ion transport process of IPMC and then modeling of its curvature and mechanical deformation in response to the input voltage. We believe the ionic content of IPMC migrate through all directions and we need to solve this problem in the three dimensions. All above mentioned NP-based method except [38] are one dimensional, [38] also is a two dimensional, non-analytical and semi-numerical method. By the way we can say that our proposed method is the first dynamic full analytical 3D ion transport based model that can cover all linearity and non-linearity of IPMC’s behavior even in large deformation situation and it is a promising result that for the first time we can find a mathematical well-defined relationship between input voltage and IPMC curvature in large deformation situation. The rest of the paper consists of three parts, part II as the main part of this paper is about proposed 3D physics-based Ion transport model and its details. After that in part III, we will estimate the unknown parameters of our model and then assess its accuracy. Finally, in part IV, we will briefly talk about the achievements of our proposed model. II. P HYSICS -B ASED I ON T RANSPORT MODEL FOR
IPMC A. Preliminary equations and obtaining the governing PDE of electric potential
Let assume we have an IPMC like Fig.1:
Fig. 1.
The assumed parametric dimension of IPMC. ( ) , , , x y z t D , ( ) , , , x y z t E , ( ) , , , x y z t and ( ) , , , x y z t represent electrical displacement, electric field, electric potential, and electric charge density, respectively. According to accepted principles in physics of electricity and magnetism, these components can be linked by the following relationships: ( ) ( ) , , ,, , , e x y z tx y z t = DE (1) ( ) ( ) , , , , , , x y z t x y z t = − E (2) ( ) ( ) . , , , , , , x y z t x y z t = D (3) Which ( ) , , , 0, 0, x y z h h W L − + and e is dielectric constant and because we want to use these equations in describing the behavior of IPMC, this coefficient is the effective dielectric constant of Nafion membrane containing sodium ions, and ( ) , , , x y z t is also defined as follows[8, 9, 41, 42, 44, 45]: ( ) ( ) ( ) , , , , , , x y z t F C x y z t C + − = − (4) Where F is the Faraday constant and ( ) +C x,y,z,t and -C are respectively the concentration of cations and anions in the membrane. Based on the continuity equation for cations, we also can find the following equation: 3 ( ) ( ) , , ,, , , C x y z tx y z t t + = − J (5) ( ) , , , x y z t J is the spatial vector of the cations flux, and it is obtained from the solution of the following PDE, known as the Nernst-Planck PDE. C F C Vd C p CRT RT + ++ + = − + + + J v (6) In the above equation, the phenomena of diffusion, ion migration, and convection have been considered that d is the ionic diffusivity, R is the gas constant, T is the absolute temperature, p is the fluid pressure, v is the free solvent velocity field and finally V represents the volumetric change of the membrane[45]. Note:
In this model, we solve the equations in a three dimensions time-variant space; that is, the equations governing the model are four variables, and they have ( ) , , , x y z t mathematic argument. But from Eq (6), in order to avoid increasing the size of equations and save space, we avoid writing this argument. Unless in those equations that lack argument makes reader mislead . Also, with the help of the modified version of Darcy’s law [55], we can link the free solvent velocity field with the gradient of fluid pressure and the electric field as follows[41, 42]: ( ) k C F p − = − v E (7) Where k is the hydraulic permeability coefficient. The second term of the right side of equation (6) is the term of convection, but since we know in the practical application of IPMC the gamut of temperature changes is not very varied and also we know the IPMC's performance is not very sensitive to temperature an especially in its optimum operation frequency (Less than 1 Hz) we can find that IPMC's performance is almost independent of temperature [56], hence convection term is ignorable, and we choose it equal to zero. To achieve zero convection term, since +C is the opposite of zero, we need to approach v to zero. If we approach v to zero, the fluid pressure gradient will be obtained according to the equation (8): if p C F − → = v E (8) Now by combination Eq (1) to (4), (6) and (8), we will have the following formula for J : ( ) ( )( ) e e d C VF RTFC C VRT −− − = − − − − − − J E E EE (9) In order to reach the J , the first we apply the divergence operator ( ) in both sides of Eq (9), then based on Eqs (1) and (3) and fundamental principles of vector calculus we substitute e − , and ( ) e − + E instead of E , ( ) and ( ) E respectively. Hence we will have: ( ) e etermterm FCd K KF − = − − + − J E (10) Where K is: ( ) K C VRT − = − (11) In the equation (10) we can find that the ratio between the term 2 and term 1 is on the order of 5×10 -10, and it means that term 2 is too smaller than term 1 and so we can ignore term 2 [41], and as a result, Eq (10) becomes the following: e FCd KF − = − − J (12) Then using a combination of equations (12), (5) and (4) we can find a PDE for calculation of as follows: d Kt = − (13) Where K is e F C d K − − and defined as (14): ( ) e F C d C VK RT − − − = (14) Now if we combine Eqs (1) to (3) the following relations will be obtained between electric potential and electric charge density: e = − (15) And finally, if we substitute (15) in (13) the following time variant 3D PDE will be resulted for ( ) , , , x y z t . In the next part, we want to solve this PDE. ( ) ( ) ( ) d Kt = − (16)
Evidence 1:
The above equation is the more complete version of PDE that Nemat-Nasser first time found for ( ) , x t in 2000 [9]. If we extend Eq (16) just for x variable, easily it will be changed to (17) that this equation is exactly the same PDE that Nemat-Nasser found for ( ) , x t (Please see Eq 37 of [9]). d ax t x x = − (17) Where a is: ( ) e F C C Va RT − − − = (18) As the first evidence, it can assure us that we have chosen a proper and more common way for our model. B. Calculation of ( ) , , , x y z t The answer of Eq (16) is the four-variable scalar function ( ) , , , x y z t .To solve this PDE, we need proper boundary conditions. The ( ) , , , x y z t represents the electric potential of the IPMC, and it is obvious that the electric potential of the IPMC in the clamping area on the surface of 4 the electrodes is the same IPMC applied input voltage, which we call it ( ) I V t . Thus, the governing boundary conditions of this boundary problem will be expressed as (19). ( ) ( )( ) ( )( ) ( ) ( ) ( )( ) ( ) ( ) ( ) ( ) ( )( ) ( )( ) ( ) ( ) ( )( ) ( ) ( ) ( ) ,0,0, t , ,0, t2 2,0,0, t , ,0, t2 2,0, , t , , , t2 2,0, , t , , , t2 2
I II II IL I L II IL I L I
V t V th h WV t V th h WV t V th L a V t h W L a V tV t V th L a V t h W L a V t + = + = − = − − = − + = + = − = − − = − (19) In these conditions, the ( ) ( ) a V tL I function is called longitudinal surface voltage attenuation function. Ideally, if IPMC electrodes are considered as perfect conductors and the fractal penetration of the electrodes into the membrane is ignored, the voltage of the clamp region is constant across the whole surface of electrodes, and there will be no any attenuations in its longitudinal direction [26]. However, since that the electrodes do not have the conductivity of a perfect conductor and also due to electrode deposition process and fractal penetration of electrodes inside the membrane, the applied voltage will be attenuated along the z -axis from clamping area of IPMC proportional to amplitude, frequency and other components of the input voltage. We describe this attenuation with a coefficient derived from ( ) V tI as ( ) ( ) a V tL I function. This function is a continuous function that produces an integer in the interval ( , which 1 means that the electrodes have perfect conductivity and means that we haven't considered the surface resistance. As a result, the measured voltage at the IPMC's tip has not been attenuated and is equal to the voltage of clamp area. Naturally, if this number is less than 1, means that we have seen the longitudinal attenuation of voltage in the problem and consider the surface resistance, which this state is more like reality. The method that we have chosen to solve Eq (16) is the separation of variables method that the calculated ( ) , , , x y z t using this method is: ( ) ( ) ( ) , , , , , , , , t x y z t k e U x y z t U x y z − = + (20) Where ( ) , , , U x y z t and ( ) , , U x y z are defined as (21) and (22): ( ) ( ) ( ) ( ) ( ) ( ) ( )( ) ( ) ( )( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( )
11 12 2 22 333 3 ˆ, , , ˆ2ˆ1 cosh ˆ ˆcoshˆ ˆˆˆ ˆ
IL I
Sinh t xV tU x y z t Sinh t hW Sinh y ySinh Wa V t cosh t Lsinh t Lsinh t z cosh t z = − + − + (21) ( ) ( )( ) ( ) ( ) ( ) , , , ,1 1 1, , sin sin sin2 mnpm n p
U x y z a x y zm n px y z x y zh W L = = = = + + + = (22) k and also are constant coefficients and ( ) ˆ t , ˆ and ( ) ˆ t define as follows: ( ) ( ) ( )
221 2 2
L I t K Ln LnW L a V t L = − − + (23)
1ˆ 1
LnW = − (24) ( ) ( ) ( ) L I t LnL a V t L = (25) Where is W and is a constant coefficient that it should be bigger than 1. And also ( ) ˆ t , ˆ and ( ) ˆ t have to always satisfy the following condition. ( ) ( ) ˆˆ ˆ ˆ t t K + + = (26) mnp a also, define by (27): ( )( ) ( ) ( ) ˆ8 , ,21 1 1 ˆ2 h W Lhmnp K x y z dx dy dzhWLa m n p Kh W L +− = + + + + + + (27) Where ( ) ( )( ) ( ) ( ) , , , ,1 1 1sin sin sin2 x y z x y zm n px y zh W L = + + + (28) And ( ) , , 2 2 2 x y z x y z xy xzyz xyz x y xz = + + + + ++ + + + + + (29) Where to and to are constant coefficients. Evidence 2:
If we extend Eq 20 just for small amount of x variable (Means small deformation situation) and ignore y and z variables, it will be changed to the following equation that it is the exact equation that Farinholt [41, 42] found for electric potential (Please see Eq 2.42 of [41]). 5 ( ) ( ) ( ) ( ) ˆ ˆ, sinh sinhˆ ˆsinh te xx t x h eh h − = − (30) Where is the boundary value of ( ) , x t in x=h . As the second evidence, it can assure us again that we have chosen a proper way for our model and proves that the Farinholt approach is just a special state of our model that only valid for small deformation situations. C. Vector proper physics matched coupling relation
When the electrical voltage is applied to IPMC, the hydrated sodium cations, ions that a couple of water molecules bind to them, move toward the IPMC cathode side, which causes the concentration of water molecules to increase gradually there then leads to Nafion swelling at the cathode side. The result of this transport process is induced stress that causes mechanical bending of IPMC [26, 37]. Thus, in fact, migration and concentration of cations cause mechanical stress and bending of the IPMC. Hence we must first calculate the concentration of cations (Electrochemical component), and then, by defining a coupling relation corresponding to the physics of the problem, we obtain the stress (Mechanical component). The process that has been used so far in many articles [14, 41-45, 57, 58] is using the idea that has been proposed by Nemat Nasser [9]. The Nemat Naser's idea process is that they first obtain the electric charge density ( ), then calculate the mechanical stress ( ) using = linear relationship. Which is a constant coefficient and plays the role of coupling between charge density and mechanical stress. According to the Eq (4), the charge density and the concentration of cations have linear relationships with each other. Therefore, it can be said that according to the coupling relationship of Nemat Naser, the relationship between the cationic concentration (The electrochemical factor that makes swelling in IPMC's membrane) and mechanical stress (The mechanical factor that makes swelling in IPMC's membrane) is followed by the following relation, which it is a linear relationship. ( ) ( ) ( ) ( )( ) ( ) ( ) ( ) , ,, ,, , x t x tx t FC x t FCx t F C x t FC + −+ − == − ⎯⎯⎯⎯⎯⎯⎯⎯→→ = − (31) Mathematically, this relation in the Nemat Nasser approach means that the behavior of concentration of cations ( C + ) is similar to the behavior of mechanical stress ( ), just with a difference in their magnitudes, while the physics of the problem tells us the opposite of this similarity. Because that base on governing mechanics in a cantilever beam (IPMC here) for making a bending in the direction of the x -axis (Main direction of IPMC’s bending) we need to have one component of induced stress in the direction of z -axis (Fig.2). On the other hand, we know the concentration of cations in the Nemat Nasser model is a function of x and t variables, and it makes an ambiguity because the stress function should have a z variable and it is in the opposite of Eq (31). Hence we can say that this relation is not a proper coupling relation and is not match with the physics of the problem at least for large deformation situations of IPMC, and we should find a new proper physics-based coupling relation. Fig. 2.
For making a bending in the direction of x -axis (Main direction of IPMC’s bending) we need to have induced stress in direction of z -axis As we mentioned before the most of papers that have presented a physical model for IPMC have assumed that cations move in the direction of x -axis between anode and cathode and it means that they have assumed the concentration of cations are one dimensional and defined as a function of x and t variables [8, 9, 14, 17, 20, 21, 37, 41-45, 57-64]. But nature is different, and ions move to all directions (Of course the dominant of them move in the direction of x -axis). This ionic movement is depended on ionic conductivity, and ionic conductivity indirectly is depended on electrical conductivity, and electrical conductivity is inversely depended on electrical resistance. Hence instead of the assumption that says ions move just in the direction of x -axis, it would be better that we assume that ions move toward a spatial vector like A that the components of this vector are depended on electric resistance of IPMC in all axis that defined as follows : ( ) ( ) ( ) , , , x M y ew z el r t r t r t = + + A i j k (32) Where ( ) , x M r t , ( ) , y ew r t and ( ) , z el r t are the functions of M r , ew r and el r , and M r , ew r and el r are membrane resistance and transverse and longitudinal resistances of electrodes respectively. As we were saying, migration and changing the concentration indirectly make mechanical stress and the bending of IPMC, now we know that this changing of the concentration is done in direction of vector A and naturally ( ) , x M r t , ( ) , y ew r t and ( ) , z el r t determine the optimum direction of A. Mathematically when we want to calculate the changes of multivariable function in direction of a vector, we should calculate the directional derivative of the function in direction of the desired vector. For example here since we want to find the changes of ion’s concentration in direction of optimum movement’s vector, we should calculate the directional derivative of ( ) , , ,
C x y z t + in direction of vector A that it is shown by ( ) , y, z, t C x + A and defines as follows: . CC C ++ + = = A aA (33) 6 Where a is the unit vector of A and defined by = Aa A . And C + is also the gradient of ( ) , , , C x y z t + . We know that in the pure and basic IPMC (An IPMC without using any controller and any modification on it), the bending displacement is attenuated over time, and it is related to the effects of resistance. It means that the effects of resistance will be increased over time and make obstacles for current flow on electrodes and ion movement on the membrane. Hence it is sensible that we choose time-ascending functions for modeling the effects of resistance, the functions like the ramp, sigmoid, tangent hyperbolic, etc, because we want to increase the effects of resistance over time. Here we have chosen the ramp function, now if we use these ramp functions for ( ) , x M r t , ( ) , y ew r t and ( ) , z el r t like (34) to (36) and define mentioned resistances as the functions of dimensions (h, W, and L). We can calculate C + A as (37): ( ) ( ) ( ) , (t) x M M r h t r h t u = (34) ( ) ( ) ( ) , (t) y ew ew r W t r W t u = (35) ( ) ( ) ( ) , (t) z el el r L t r L t u = (36) ( ) ( ) ( ) ( )
1ˆ , ,
M ew el
C C CC r h r W r Lr h W L x y z + + ++ = + + A (37) Where u(t) is the Heaviside step function and ( ) ˆ , , r h W L is: ( ) ( ) ( ) ( ) ˆ , , M ew el r h W L r h r W r L = + + (38) By having C + A we can find the changes of ion’s concentration in the direction of optimum movement’s vector that these changes induce mechanical stress, semi-directly and proper physics matched. Hence we can say that there is a direct linear relationship like (39) between C + A and mechanical stress: ( ) ( ) , , , , , , x y z t C x y z t + = A (39) Where is a constant coupling coefficient that we should estimate it in the model estimation stage. D. Mechanical Stress
To find induced stress using Eq (39). The first we should calculate C + A and for calculation of C + A we need to know C + . It is so simple to find C + by the combination of (1) to (4) as (40): e C C F + − = − (40) Now we can calculate C + A using (37) and finally will be obtained by (41): ( ) ( ) ( ) ( )
40 1 2 , , , , , , , ,ˆ , W, L e ii x y z t A x y z t A x y zF r h = = − + (41) Where ( ) ( ) ( ) ˆ M ew el
U U UA K r h r W r Lx y z = + + (42) ( ) M U U UA r h x x y x z = + + (43) ( ) ew U U UA r W y x y y z = + + (44) ( ) el U U UA r L z x z y z = + + (45) Where U and U already were obtained by (21) and (22), and we can calculate A to A easily using (42) to (45). E. Bending moment
Based on beam theory, IPMC strip is a cantilever beam, and we can consider it as a cantilever beam to calculate it’s bending moment. Basically, in a cantilever beam, the global bending moment is described by a vector like M that defined by (46): ( ) y x zx zy zzc c dx dy = + + Μ r i j k (46) Where r is the position vector and is defined by x y = + r i j . If we expand (46) we will have (47): zx zy zz M M M = + +
Μ i j k (47) Where y x zx zzc c
M y dx dy = (48) y x zy zzc c M x dx dy = (49) ( ) y x zz zy zxc c M x y dx dy = − (50) As has been depicted in Fig 3, M zx is absolute zero, and M zz is almost zero for IPMC, and in fact, only the M zy induces to the IPMC and makes IPMC bend in the direction of the x -axis. As described above we can say that the only effective bending moment for IPMC is M zy that if we extend it for the dimension of IPMC it could be defined as (51): ( ) ( ) , , , , W hzy zzh
M z t x x y z t dxdy +− = (51) Where ( ) , , , x y z tzz is the renamed version of the same ( ) , , , x y z t that we defined it in (41). Now by combination of (51) and (41) to (45) and solving integration terms and some mathematical efforts we can find M zy as (52): ( ) ( ) ( , ) ( , ) ( )ˆ , W, L ezy zy zy M z t M z t M zF r h = − + (52) Where ( ) ( )( )( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( )
31 21 2 21 13 33 ˆ ˆ2( , ) ˆ ˆ ˆˆ ˆ1 coth tanh( )2ˆ ˆcsch cothˆcosh elzyL I I
K t r LM z t t Wt t ha V t t L t Lt z V t = − − − (53) 7
32 2 5 6 mnpzy h W b WM z W z = − + + (54)
Fig. 3.
Three possible bending moments for the IPMC strip. M zx is absolute zero (a), and M zz is almost zero for IPMC (b). But M zy bending moment that makes IPMC bend in x -direction (c). That mnp b is a constant coefficient that related to numbers of harmonics that we use in calculation of U and ˆ , ˆ and ˆ are constant coefficients that they have created by the combination of to and to and we should predict them in the model prediction stage. By the way, if we want to show M zy by a unit equation we can use the following: ( ) ( ) ( ) ( ) ( ) ˆ( , ) cosh z zy M z t t t z t t = + + (55) Where ( ) t to ( ) t have been defined by (56) to (58) respectively: ( ) ( )
301 2 5 e mnp h W Wt t bF r h = + (56) ( ) ( ) ( )
302 6 e mnp h Wt t W bF r h = (57) ( ) ( ) ( ) ( )( )( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( )
303 21 221 13 3 ˆ ˆ2ˆ , W, L ˆ ˆˆˆ ˆ1 coth tanh( )2ˆ ˆcsch coth eleL I I
K t r Lt tF r h t Wt t ha V t t L t L V t = − − (58) F. Calculation of IPMC tip displacement in large deformation situation
To calculate the IPMC bending, the first we need to have curvature shape function of IPMC. Then using this function we can find its tip displacement. But in the most of paper that model IPMC behavior using physics-based approaches, there is a common mistake in this direction and it is applying linear Euler-Bernoulli (LEB) beam theory to calculate tip displacement of IPMC (For example [14, 41, 42, 44, 45, 57, 58, 65]). LEB beam theory only can consider the infinitesimal strains and it means that by this theory we can only model small deformation of a beam (Here IPMC) [66]. Based on LEB beam theory the relation of curvature shape function of IPMC, ( ) , z t , and bending moment, ( ) zy M z,t , is defined as follows: ( ) ( ) ,, zy M z tz tz Y I = (59) Where Y is Young's modulus of Nafion membrane and I is moment of inertia and calculated by the following relation: ( ) L W hh
I x z dx dy dz I h WL +− = + = (60) As mentioned before, Eq (59) is valid just for infinitesimal strains, and we cannot use it for large deformation situations of IPMC. Hence we should find a valid relation for large deformation that this relation is main relation of curvature. Mathematically the curvature of a function like ( ) , z t will be calculated by ( ) , Math z t where defined as following: ( ) ( )( ) ,, ,1 Math z tzz t z tz = + (61) Mechanically also the curvature of a beam that ( ) zy M z,t bending moment has been induced on it is defined by (62): ( ) ( ) ,, zyMech M z tz t Y I = (62) Now our aim is the finding a dynamic function like ( ) , z t that its shape is match with the shape of a cantilever beam (Here IPMC) that a bending moment like ( ) zy M z,t has been induced on it. To find this function it is only enough to set ( ) , Math z t and ( ) , Mech z t equal and solve resulted PDE. If we do this, we will obtain following nonlinear PDE: 8 ( ) ( )( ) ( ) ( )
32 22 2 , ,,, ,1 0
Math Mechzy if z t z tM z tz t z tz Y I z = − + = (63)
Evidence 3:
In the above PDE, if we ignore term ( ) z tz , it will be changed to governing PDE of LEB beam theory (Eq 59). Hence it will be the third evidence and assure us again that our method has chosen a proper way and has presented a complete approach to physics-based modeling of IPMC. Eq (63) is a nonlinear PDE, but it is solvable analytically. If we solve this PDE, ( ) , z t will be obtained as follows: ( ) ( ) ( ) ( ) ( )( )( ) ( ) ( ) ( ) , ( ) ,, 1 ( ) 0,0, f t Ln P z t S t q z tz t g tf t Ln S t q tg t P t += + − (64) Where: ( ) ( ) ( )( ) ( ) ,tan ( )2, ,tan ( )2 q z tf t S t Y IP z t q z tf t S t Y I + − = − − (65) ( ) ( ) ( ) ( ) ( ) S t YI f t YI f t = − + (66) ( ) ( ) ( ) ,, sin m z t f tq z t Arc Y I + = (67) ( ) f tf f t k e − = (68) ( ) g tg g t k e − = (69) That f k , g k , f and g are constant coefficients that they will be estimated in the model estimation stage. ( ) , z t gives us the shape of IPMC, but we want to find the position of its tip overtime or the same dynamic tip displacement. If we name the function of IPMC tip displacement ( ) Tip t , it is obvious that the following equation, that has been used in the several papers like[45], is absolutely wrong, and it is just an approximation of IPMC tip displacement in the very small deformation situation. ( ) ( ) , Tip t L t = (70) In the Eq (70), L is the length of IPMC, and this relation is valid only when IPMC is open circuit and is completely off and don’t have any movement. To solve this big mistake we should use the following formula: ( ) ( ) ( ) , Tip L t z t t = (71) Where ( ) L z t makes the dynamic coordination of the z component of IPMC’s tip using Eq (72): ( ) ( ) ,1 L z t z t dz Lz + = (72) The source of the above relation is that we know the arc length of IPMC is always constant and equal to L, and the above formulation is this equilibrium relationship, and the integration term is the arc length of IPMC that we have modeled it. Evidence 4:
As you can follow in equations (70) to (72), it was shown that the common method to find IPMC tip displacement from ( ) , z t (Eq (70)) is absolutely wrong and it is just an approximation in the very small deformation situation. But using our proposed relation (Eq (71)) we can find accurate tip displacement of IPMC over time in all deformation situations (small to large deformation situations). This evidence also tells us again that our method is more complete and more accurate approach than previous models. III. P ARAMETERS E STIMATION AND M ODEL V ALIDATION
In this part, we want to validate our claim to find an accurate, well-defined relationship between input applied voltage ( ( ) I V t ) and output tip displacement of IPMC ( ( )
Tip t ). In this direction, we have used experimental data that we have collected from IPMC. The method of measuring tip displacement is based on camera and using image processing algorithms where the detail of this system and algorithm and also details of required hardware apparatus have been described in our previous papers [2-5], and we don’t talk about them here. A. Specification
About our sample, we should say that its dimension is almost 28 ; also the electrolyte is Na +, and the type of membrane is DuPont™ Nafion 117 PFS that coated with two thin layers of Pt. In this model and basically in all physics-based models we have some specified physical parameters and some unknown parameters. The specified physical parameters have been reported in the references and have been collected in Table.1, but unknown parameters should be estimated in parameter estimation stage. Table. 1.
Specified physical parameters
Cmol
F = 96458 Fe m = − Pa Y = molm C = − ms d = − . Jmol K
R = K T =
28 6 0.0889 mm L W h = B. Parameter estimation and model validation
The method of validation is that we test the IPMC for three signals as input voltages, i.e., sine wave, chirp, and PRBS. That means we estimate the unknown parameters using half of each of these signals and test it with the second half. For example, if we measure the tip displacement of IPMC for 80 sec in response an arbitrary input voltage, the first we estimate 9 the unknown parameters using the first 40 sec of the dataset and then we use the estimated parameters for the model and test it using the second 40 second of the dataset. The method of estimation is not the goal of this research, and we can find a variety of classic and intelligent methods to estimate unknown parameters that here we have used a genetic algorithm-based optimization method to find unknown parameters. The estimated parameters using genetic algorithms have been collected in Table.2.
Table. 2.
Estimated parameters using a genetic algorithm kf ˆ 6 -3 f rM k g rew g rel ˆ 2 ˆ 5 Bmnp
As you can see in Figs 4 to 6, this model is accurate enough, and it can follow the actual output precisely. But we need to assess the accuracy of model numerically that as a proper well known numerical measure we can propose the Normalized Mean Square Error (NMSE). The NMSE of our model in response to Sine wave, PRBS, and Chirp input signals are 0.07, 0.025, and 0.0047, respectively, and their average is 0.0332 that these numbers tell us proposed model is accurate enough and it can work properly even in large deformation situations and practical applications.
Fig. 4.
Model validation. (Up) applied Chirp voltage to the IPMC (with a peak voltage of 2.3 V). (Down) the actual and estimated tip displacement of IPMC in x- direction axis in response to the Chirp input voltage. Fig. 5.
Model validation. (Up) applied PRBS voltage to the IPMC (with a peak voltage of 2.3 V). (Down) the actual and estimated tip displacement of IPMC in the x- direction axis in response to the PRBS input voltage. Fig. 6.
Model validation. (Up) applied Sine wave voltage to the IPMC (with a peak voltage of 2.3 V and a frequency of 0.16 Hz). (Down) the actual and estimated tip displacement of IPMC in x- direction axis in response to the Sinewave input voltage. IV. C ONCLUSION
It this paper as the first time we proposed a fully analytical and physics-based ion transport 3D and non-Linear model for large deformable IPMC. In this direction, based on three dimensional Nernst-Plank PDE, we found a well-defined and valid relationship between the input voltage and output tip displacement of IPMC that it is valid for large deformation situation. Also, using four provable pieces of evidence we explained that proposed model had chosen a proper way and it is more complete than the previous benchmark and well- 10 known physics-based models. And finally, with some experimental examples, we showed that our model is accurate enough and is valid for large deformation situations. R
EFERENCES [1] M.Shahinpoor, Y.BarCohen, J. O. Simpson, and J. Smith, "Ionic Polymer-MetalComposites (IPMC) As Biomimetic Sensors, Actuators & Artificial Muscles- A Review,"
J. of Smart Materials, pp. 1-27, 1998. [2] M. Annabestani and N. Naghavi. Nonuniform deformation and curvature identification of ionic polymer metal composite actuators [3] M. Annabestani and N. Naghavi, "Nonlinear identification of IPMC actuators based on ANFIS–NARX paradigm,"
Sensors and Actuators A: Physical, vol. 209, pp. 140–148, 2014. [4] M. Annabestani, M. Maymandi-Nejad, and N. Naghavi, "Restraining IPMC Back Relaxation in Large Bending Displacements: Applying Non-Feedback Local Gaussian Disturbance by Patterned Electrodes,"
IEEE Transactions on Electron Devices, vol. 63, no. 4, pp. 1689-1695, 2016, doi: 10.1109/TED.2016.2530144. [5] M. Annabestani, N. Naghavi, and M. Maymandi-Nejad, "Non-autoregressive Nonlinear Identification of IPMC in Large Deformation Situations Using Generalized Volterra-based Approach,"
IEEE Transaction on Instrumentation and Measurement, vol. PP, no. 99, pp. 1-7, 2016. [6] M. Annabestani, N. Naghavi, and M. Maymandi-Nejad, "From modeling to implementation of a method for restraining back relaxation in ionic polymer–metal composite soft actuators,"
Journal of Intelligent Material Systems and Structures, vol. 29, no. 15, pp. 3124-3135, 2018, doi: 10.1177/1045389x18783082. [7] M. Annabestani, H. Mohammadzadeh, A. Aghassizadeh, S. Azizmohseni, and M. Fardmanesh, "Active microfluidic micromixer design using Ionic polymer-metal composites," in
Electrical Engineering (ICEE), Iranian Conference on , 30 April-2 May 2019 2019. [8] S.Nemat-Nasser, "Micromechanics of actuation of ionic polymer-metal composites,"
Journal of Applied physics, vol. 92, No. 5, pp. 2899-2915, 2002. [9] S.Nemat-Nasser, and J.Y.Li, "Electromechanical response of ionic polymer metal composites,"
Journal of Applied physics, vol. 87, No. 7, pp. 3321-3331, 2000. [10] M. Shahinpoor, and K.J Kim "Ionic polymer–metal composites: I. Fundamentals,"
Institute of Physics Publishing, Journal of Smart Materials and Structures, pp. 819–833, 2001. [11] Y.Bar-Cohen, "Electroactive polymer (EAP) actuators as artificial muscles: reality, potential, and challenges,"
SPIE Press, pp. 139-191, 2001. [12] R. Dong, and X.Tan, "Modeling and open-loop control of IPMC actuators under changing ambient temperature,"
Proceedings of the SPIE,
Smart Mater. Struct, vol. 25, p. 115037 2016. [14] H. Moeinkhah, J. Rezaeepazhand, and A. Akbarzadeh, "Analytical dynamic modeling of a cantilever IPMC actuator based on a distributed electrical circuit,"
Smart Materials and Structures, vol. 22, pp. 1-10, 2013. [15] P. Brunetto, L. Fortuna, P. Giannone, S. Graziani, and S. Strazzeri, "Static and Dynamic Characterization of the Temperature and Humidity Influence on IPMC Actuators,"
IEEE Transactions On Instrumentation And Measurement, vol. 59, no. 4, pp. 893-908, 2010. [16] A. J. McDaid, K. C. Aw, E. Hämmerle, and S. Q. Xie, "A nonlinear scalable model for designing ionic polymer-metal composite actuator systems," in
Second International Conference on Smart Materials and Nanotechnology in Engineering , J. Leng, A. K. Asundi, and W. Ecke, Eds., 2009, vol. 7493. [17] A.Punning, U.Johanson, M.Anton, A.Aabloo, and M. Kruusmaa, "A Distributed Model of Ionomeric Polymer Metal Composite,"
Journal of Intelligent Material Systems and Structures, vol. 0, pp. 1-14, 2009. [18] Y.-c. Chang and W.-j. Kim, "An Electrical Model with Equivalent Elements in a Time-variant Environment for an Ionic-polymer- metal-composite System,"
International Journal of Control, Automation and Systems, vol. 15, no. 1, pp. 45-53, 2017. [19] H. Liu, K. Xiong, K. Bian, and K. Zhu1, "Experimental study and electromechanical model analysis of the nonlinear deformation behavior of IPMC actuators,"
Acta Mech. Sin, vol. 33, no. 2, pp. 382–393, 2017. [20] D. Pugal, K. J. Kim, and A. Aabloo, "An explicit physics-based model of ionic polymer-metal composite actuators,"
Journal of Applied Physics, vol. 110, 2011. [21] R. Caponetto, V. D. Luca, G. D. Pasquale, Salvatore Graziani, F. Sapuppo, and E. Umana, "A Multiphysics Frequency-Dependent Model of an IP C Actuator,"
IEEE Transactions On Instrumentation And Measurement, vol. 63, no. 5, pp. 1347-1355, 2014. [22] D. Vokoun, Q. He, L. Heller, M. Yu, and Z. Dai, "Modeling of IPMC Cantilever’s Displacements and Blocking Forces,"
Journal of Bionic Engineering, vol. 12, pp. 142-151, 2015. [23] Z. Zhu et al. , "Application-oriented simplification of actuation mechanism and physical model for ionic polymer-metal composites,"
Journal of Applied Physics, vol. 120, p. 034901 2016. [24] D.Pugal, K.J.Kim, A.Punning, H.Kasemägi, M.Kruusmaa, and A.Aabloo, "A self-oscillating ionic polymer-metal composite bending actuator,"
Journal of Applied Physics, vol. 103, p. 084908, 2008. [25] D.Bandopadhya, "Derivation of Transfer Function of an IPMC Actuator Based on Pseudo-Rigid Body Model,"
Journal of Reinforced Plastics and Composites, pp. 372-390, 2010. [26] Y.Gong, J.Fan, C.Y. Tang, and C.P.Tsui, "Numerical Simulation of Dynamic Electro-Mechanical Response of Ionic Polymer-Metal Composites,"
Journal of Bionic Engineering, pp. 263–272, 2011. [27] K.Yun, and W.J. Kim, "System identification and microposition control of ionic polymer metal composite for three-finger gripper manipulation,"
Journal of Systems and Control Engineering, pp. 539-551, 2006. [28] K.Yun, and W.J.Kim, "Microscale position control of an electroactive polymer using an electroactive polymer using an anti-windup scheme,"
Journal of Smart Materials and Structures, pp. 924–930, 2006. [29] M.Vahabi, E.Mehdizadeh, M.Kabganian and F.Barazandeh, "Experimental identification of IPMC actuator parameters through incorporation of linear and nonlinear least squares methods,"
Sensors and Actuators A: Physical, vol. 168, no. 1, pp. 140-148, 2011. [30] R. Dong and Y. Tan, "A model based predictive compensation for ionic polymer metalcomposite sensors for displacement measurement,"
Sensors and Actuators A: Physical, vol. 224, pp. 43-49, 2015. [31] D.Q.Truong, and K.K. Ahn, "Design and Verification of a Non-linear Black-Box Model for Ionic Polymer Metal Composite Actuators,"
Journal of Intelligent Material Systems and Structures,
Sensors and Actuators A: Physical, vol. 183, pp. 105-114, 2012. [33] Z. Zhu, H. Chena, L. Chang, and B. Li, "Dynamic model of ion and water transport in ionic polymer-metal composites,"
AIP ADVANCES vol. 1, 2011. [34] P.G.d.Gennes, K.Okumura, M.Shahinpoor, and K.J.Kim, "Mechanoelectric Effects in Ionic Gels,"
Europhysics Letters, vol. 50, no. 4, pp. 513-518, 2000. [35] J. W. Paquette, K. J. Kim, and D. Kim, "Low temperature characteristics of ionic polymer–metal composite actuators,"
Sensors and Actuators A: Physical, vol. 118, no. 1, pp. 135-143, 2005/01/31/ 2005, doi: https://doi.org/10.1016/j.sna.2004.08.006. [36] K. M. Newbury and D. J. Leo, "Electromechanical Modeling and Characterization of Ionic Polymer Benders,"
Journal of Intelligent Material Systems and Structures, vol. 13, no. 1, pp. 51-60, 2002. [37] S.Tadokoro , S.Yamagami and T.Takamori, "An actuator model of ICPF for robotic applications on the basis of physico-chemical hypotheses,"
Proceeding of IEEE International Conference on Robotics and Automation-San Francisco-USA, pp. 1340–1346, 2000. [38] Y. Toi and S. S. Kang, "Finite element analysis of two-dimensional electrochemical–mechanical response of ionic conducting polymer–metal composite beams," Computers & Structures, vol. 83, no. 31–32, pp. 2573–2583, 2005. [39] P. J. C. Branco and J.A.Dente, "Derivation of a Continuum Model and Its Electric Equivalent-Circuit Representation for The Ionic Polymer-Metal Composite (IPMC) Electromechanics,"
Smart Materials and Structures vol. 15, no. 2, 2006. [40] T. Yamaue, H. Mukai, K. Asaka, and M. Doi, "Electrostress Diffusion Coupling Model for Polyelectrolyte Gels,"
Macromolecules, vol. 38, no. 4, pp. 1349–1356, 2005. [41] K. M. Farinholt and D. J. Leo, "Modeling and characterization of ionic polymer transducers for sensing and actuation," Doctor of Philosophy, Mechanical Engineering, Virginia Polytechnic Institute State University, Blacksburg, Virginia, 2005. [42] K. Farinholt and D. J. Leo, "Modeling of electromechanical charge sensing in ionic polymer transducers,"
Mechanics of Materials, vol. 36, no. 5, pp. 421-433, 2004/05/01/ 2004, doi: https://doi.org/10.1016/S0167-6636(03)00069-3. [43] K. M. Farinholt and D. J. Leo, "Modeling the electrical impedance response of ionic polymer transducers,"
Journal of Applied Physics, vol. 104, p. 014512 2008. [44] Z. Chen, D. R. Hedgepeth, and X. Tan, "A nonlinear, control-oriented model for ionic polymer–metal composite actuators,"
Smart Materials and Structures, vol. 18, no. 5, pp. 1-9, 2009. [45] Z. Chen and X. Tan, "A Control-Oriented and Physics-Based Model for Ionic Polymer–Metal Composite Actuators,"
IEEE/ASME Transactions On Mechatronics, vol. 13, no. 5, pp. 519-529, 2008. [46] J. Choonghee, H. E. Naguib, and R. H. Kwon, "Modeling and optimization of the electromechanical behavior of an ionic polymer–metal composite,"
Smart Materials and Structures, vol. 17, p. 065022, 2008. [47] T. Johnson and F. Amirouche, "Multiphysics modeling of an IPMC microfluidic control device,"
Microsystem Technologies, vol. 14, no. 6, pp. 871–879, 2008. [48] S. Nemat-Nassera and S. Zamani, "Modeling of electrochemomechanical response of ionic polymer-metal composites with various solvents,"
Journal of Applied Physics, vol. 100, p. 064310 2006, doi: 10.1063/1.2221505. [49] M. Porfiri, "Charge dynamics in ionic polymer metal composites,"
Journal of Applied Physics vol. 104, p. 104915 2008, doi: 10.1063/1.3017467. [50] M. Aureli, W. Lin, and M. Porfiri, "On the capacitance-boost of ionic polymer metal composites due to electroless plating: Theory and experiments,"
Journal of Applied Physics, vol. 105, p. 104911 2009, doi: 10.1063/1.3129503. [51] T. Wallmersperger, D. J. Leob, and C. S. Kothera, "Transport modeling in ionomeric polymer transducers and its relationship to electromechanical coupling,"
Journal of Applied Physics, vol. 101, p. 024912 2007, doi: 10.1063/1.2409362. [52] L. Zhang and Y. Yang, "Modeling of an ionic polymer–metal composite beam on human tissue,"
Smart Materials and Structures, vol. 16, pp. S197–S206, 2007. [53] E. T. Enikov and G. S. Seo, "Analysis of water and proton fluxes in ion-exchange polymer–metal composite (IPMC) actuators subjected to large external potentials,"
Sensors and Actuators A: Physical, vol. 122, no. 2, pp. 264-272, 2005/08/26/ 2005, doi: https://doi.org/10.1016/j.sna.2005.02.042. [54] V. Volpini, L. Bardella, A. Rodella, Y. Cha, and M. Porfiri, "Modelling compression sensing in ionic polymer metal composites,"
Smart Materials and Structures, vol. 26, p. 035030, 2017. [55] P. E. Grimshaw, J. H. Nussbaum, and A. J. Grodzinsky, "Kinetics of electrically and chemically induced swelling in polyelectrolyte gels,"
The Journal of Chemical Physics vol. 93, no. 6, 1998. [56] R. Dong and X. Tan, "Modeling and open-loop control of IPMC actuators under changing ambient temperature,"
Smart Materials And Structures, vol. 21, 2012. [57] H.Moeinkhah, J-Y. Jung, J-H Jeon,A.Akbarzadeh, J.Rezaeepazhand, K.C.Park,and I-K.Oh, "How does clamping pressure influence actuation performance of soft ionic polymer–metal composites?,"
Smart Materials and Structures,
IEEE/ASME Transactions On Mechatronics, vol. PP, no. 99, pp. 1-9, 2014. [59] S.Tadokoro , M.Fukuhara , Y.Maeba , M.Konyo , T.Takamori and K.Oguro "A dynamic model of ICPF actuator considering ion-induced lateral strain for Molluskan robotics,"
Proceedings of the 2002 IEEE/RSJ International Conference on Intelligent Robots and Systems-Lausanne-Switzeland, pp. 2010–2017, 2002. [60] Z. Chen, S. Shatara, and X. Tan, "Modeling of Biomimetic Robotic Fish Propelled by An Ionic Polymer–Metal Composite Caudal Fin,"
IEEE/ASME Transactions On Mechatronics, vol. 15, no. 3, pp. 448-459, 2010. [61] Y. Toi and W. S. Jung, "Finite element modeling of electrochemical-poroelastic behaviors of conducting polymers,"
Computers & Structures, vol. 85, no. 19–20, pp. 1453–1460, 2007. [62] R. Caponetto, V. D. Luca, S. Graziani, and F. Sapuppo, "An optimized frequency-dependent multiphysics model for an ionic polymer–metal composite actuator with ethylene glycol as the solvent,"
Smart Materials And Structures, vol. 22, pp. 1-14, 2013. [63] P. Brunetto, L. Fortuna, S. Graziani, and S. Strazzeri, "A model of ionic polymer–metal composite actuators in underwater operations,"
Smart Materials and Structures, Volume vol. 17, no. 2, p. 025029, 2008. [64] Y.Cha, M.Aureli, and M.Porfiri, "A physics-based model of the electrical impedance of ionic polymer metal composites,"
Journal of Applied Physics,
Journal of Applied Physics, vol. 114, 2013. [66] J. N. Reddy,