On a nonlinear electromechanical model of nerve
aa r X i v : . [ q - b i o . N C ] F e b EPJ manuscript No. (will be inserted by the editor)
On a nonlinear electromechanical model of nerve
Alain M. Dikand´e a Laboratory of Research on Advanced Materials and Nonlinear Sciences (LaRAMaNS), Department of Physics, Faculty of Science,University of Buea, P.O. Box 63 Buea, CameroonReceived: date / Revised version: date
Abstract.
The generation of action potential brings into play specific mechanosensory stimuli manifestin the variation of membrane capacitance, resulting from the selective membrane permeability to ionsexchanges and testifying to the central role of electromechanical processes in the buildup mechanism ofnerve impulse. As well established [See e.g. D. Gross et al, Cellular and Molecular Neurobiology vol. 3,p. 89 (1983)], in these electromechanical processes the net instantaneous charge stored in the membraneis regulated by the rate of change of the net fluid density through the membrane, orresponding to thedifference in densities of extacellular and intracellular fluids. An electromechanical model is proposed forwhich mechanical forces are assumed to result from the flow of ionic liquids through the nerve membrane,generating pressure waves stimulating the membrane and hence controlling the net charge stored in themembrane capacitor. The model features coupled nonlinear partial differential equations: the familiarHodgkin-Huxley’s cable equation for the transmembrane voltage in which the membrane capacitor is nowa capacitive diode, and the Heimburg-Jackson’s nonlinear hydrodynamic equation for the pressure wavecontrolling the total charge in the membrane capacitor. In the stationary regime, the Hodgkin-Huxley cableequation with variable capacitance reduces to a linear operator problem with zero eigenvalue, the boundstates of which can be obtained exactly for specific values of characteristic parameters of the model. Inthe dynamical regime, numerical simulations of the modified Hodgkin-Huxley equation lead to a variety oftypical figures for the transmembrane voltage, reminiscent of action potentials observed in real physiologicalcontexts.
PACS.
XX.XX.XX No PACS code given
The generation of nerve impulse is one of most actively in-vestigated problems in the history of Neuroscience [1,2,3,4,5,6,7,8,9,10,11,12,13,14,15,16]. Study of the problemis motived by the crucial need for a good understand-ing of characteristic properties of the action potential, as-sumed to be a propagating form of transmembrane volt-age along the axon. Measurements of action potentialsin several physiological contexts have generated a wealthof data that triggered a great deal of theoretical interestin the phenomenon. As pioneer in this theoretical inter-est, the Hodgkin-Huxley model [1,2] rests on a picture bywhich the nerve impulse is an electric voltage propagat-ing in form of an asymmetric pulse along the nerve fiber.Originally the Hodgkin-Huxley model was introduced toexplain data obtained from measurements of conductiveparameters of a nerve fiber, and particularly to show howthese data could be used to directly calculate both theshape and velocity of an action potential on the squid gi-ant axon [17]. a Email address: [email protected]
According to the Hodgkin-Huxley model [1], the nerveimpulse is a self-regenerative wave associated with theelectrochemical activity of the nerve cell, and due to theflow of ion currents (Na + and K + ) through specific ionchannels. This wave propagates with a constant shape,through a mechanism that can be summarized as follow:During the generation and transmission of the nerve im-pulse, the leading edge of the depolarization region of thenerve triggers adjacent membranes to depolarize, causinga self-propagation of the excitation related to the trans-membrane voltage down the nerve fiber [1,18,19]. Hodgkinand Huxley suggested that a convenient way to describethe propagation of this transmembrane voltage is to re-gard the nerve fiber as an electric cable. Thus, in its mostconventional formulation, the Hodgkin-Huxley model as-sumes currents in intracellular and extracellular fluids tobe ohmic such that the net transmembrane current is thesum of ionic and capacitive currents. In this picture theconservation law for currents passing through the mem-brane can be written [1]: C m ∂V∂t = D ∂ V∂x − F ( V ) , (1) Alain M. Dikand´e: On a nonlinear electromechanical model of nerve where V is the transmembrane voltage, C m is the mem-brane capacitance, D is the diffusion coefficient and F accounts for contributions from some ion currents.Still, besides the indisputable electrical activity of theaxonal membrane, experiments have also pointed out [20,21,22,23,24] the existence of mechanical constraints re-lated to pressures due to flows of fluids through the mem-brane. Concretely these mechanical constraints are elec-tromechanical forces that are responsible for mechanotrans-duction processes [25], physiological processes in whichmechanical forces such as pressures exerted by ionic fluidson cell membranes and tissues, can trigger excitations ofelectrical natures playing important role in the control ofvarious stimuli-responsive organs, in homeostasis of livingorganisms and so on [25].Taking advantage of experiments suggesting sizablethermodynamic phenomena preceeding and following theaction potential, and specifically the liquid-gel transitionobserved at some critical temperature [26,27,28,29], He-imburg and Jackson [7,30,31] suggested that mechanicalforces related to pressure waves could play a major rolein nerve membrane excitation and subsequently in thebuildup of nerve impulse. In this respect they postulatedthat pressure waves associated with propagation of thedensity difference between fluids flowing through the nervemembrane, could actually be a mechanical manifestationof the action potential. Most recently there have been fewother attempts to revisit the mathematical description ofthe nerve impulse, with a main aim to combine the contri-butions of electrical and mechanical processes [25,26,27,29,32].In this work we propose a model describing the elec-tromechanical process of generation of the action poten-tial. The model combines the Hodgkin-Huxley cable modeland the pressure-wave model proposed by Heimburg andJackson [7]. Our model assumes that the membrane capac-itance changes instantanously with the difference in densi-ties of fluids through the membrane, leading to a modifiedHodgkin-Huxley equation where the membrane capacitornow behaves like a ”feedback” component (i.e. like a ca-pacitive diode).In sec. 2 we present the model which consists of twononlinear partial differential equations, namely the modi-fied Hodgkin-Huxley equation for the action potential andthe Boussinesq equation for the density-difference wave [7,30]. In sec. 3 we first consider the stationary regime ofthe action-potential (or modified Hodgkin-Huxley) equa-tion. In this purpose we use the exact soliton solutionto the Korteweg-de Vries (KdV) equation derived fromBoussinesq’s equation, to recast the modified Hodgkin-Huxley equation into a linear operator problem with zeroeigenvalue. Three exact bound-state solutions to this lin-ear operator problem are obtained analytically, for specificvalues of characteristic parameters of the model. In sec.4, numerical simulations of the modified Hodgkin-Huxleyequation are carried out assuming the three stationary so-lutions as initial profiles of the action potential. In sec. 5we conclude the study. The axon can be regarded as a long cylinder with wallsmade of cell membrane surrounded by intracellular and ex-tracellular fluids [2,32]. The intracellular fluid stands for aconductive liquid with a high concentration of potassiumions but a low concentration of sodium and chlorine ions,while the axonal cell membrane acts like a barrier pre-venting ions in the intracellular liquid from mixing withexternal solutions. Due to the difference in ion concen-trations in intracellular and extracellular fluids, a restingpotential is expected to set up through the membrane.If the nerve is depolarized, e.g. due to the presence of astimulus of any kind, the axon membrane will become se-lectively permeable to ionic currents which flow rapidlyinto the cell, reversing the polarity of the action potential[1,2].In general, for a fixed number of charged lipids aroundthe cell membrane, the charge density will be different be-cause the respective lipid areas are different [25]. Thereforewe can expect changes in the electrostatic potential of themembrane during a propagating pressure wave, indicatinga possibility of an electromechanical coupling between thenet fluid density and the electrostatic potential on the cellmembrane. This electromechanical coupling, first reportedby Petrov [33] and widely observed in recent experimentsin neurophysiology [25,26,27,34,35,36], can also be linkedwith changes in membrane capacitance as a result of vari-ation of the fluid density through the membrane.The model proposed in this study retains the key in-gredients [1] of the Hodgkin-Huxley cable model, exceptfor the self-regulatory function of the membrane capaci-tance now assumed to vary instantaneously with the netfluid density on the membrane. With this consideration,the system dynamics can be described by the following setof two nonlinear space-time partial differential equations:
D ∂ V∂x = ∂∂t (cid:18) C m ( x, t ) V (cid:19) , (2) ∂ U∂t = c ∂∂x (1 − U ) ∂U∂x ! − h ∂ U∂x . (3)In eq. (3) we introduced a dimensionless variable U = ∆ρ A /ρ to represent the density difference ∆ρ A = ρ A − ρ A , note that physical meanings of parameters ∆ρ A , ρ A and ρ are discussed in detail in refs. [7,21,31,30].To describe the instantaneous change of the membranecapacitance C m due to variation of the ion-carrying fluiddensity [22,23], we postulate that when the nerve is activethe rate of change of the membrane capacitance is propor-tional to the net density of ion-carrying fluid ∆ρ A on themembrane i.e.: ∂C m ( x, t ) ∂t = κ∆ρ A , (4)where κ is assumed positive. Using eq. (4) the modifiedHodgkin-Huxley equation (2) becomes: C m ( x, t ) ∂V∂t = D ∂ V∂x − κ∆ρ A ( x, t ) V, (5) lain M. Dikand´e: On a nonlinear electromechanical model of nerve 3 where the membrane capacitance C m ( x, t ) is given by: C m ( x, t ) = C + κ Z ∆ρ A ( x, t ) dt. (6)Instructively the value κ = 0 reproduces the standardHodgkin-Huxley model [1,2], however for nonzero valuesof κ eq. (5) turns to a modified Hodgkin-Huxley equa-tion whose solution depends on the spatio-temporal pro-file of the density-difference wave ∆ρ A ( x, t ). In the nextsection, using the exact one-soliton solution to eq. (3),we seek for possible analytical solutions to the modifiedHodgkin-Huxley equation (2). In this respect we shall seethat the modified Hodgkin-Huxley equation is analyticallytractable only in the steady-state regime. Indeed in thisregime the modified Hodgkin-Huxley equation reduces toa linear-operator problem with zero eigenvalue, the boundstates of which are Legendre polynomials [37]. By introducing new coordinates; U ( x, t ) = ψ ( ξ, T ) , ξ = cc ( x − c t ) , T = hc t, (7)and integrating once with respect to the new variable ξ ,eq. (3) reduces to the KdV equation [38]: ∂ψ∂T = αψ ∂ψ∂ξ − β ∂ ψ∂ξ , (8)where α and β are constants depending on c , h and c . Theparameters α and β can be set to any values through judi-cious coordinate transformations, however we shall retainthe most widely used values of these parameters namely α = 6 and β = − U ( x, t ) = − sech ( x − t ) , (9)which is a localized wave of depression.With the help of the one-pulse solution (9) we can re-express the modified Hodgkin-Huxley equation (5) as: C m ( x, t ) ∂V∂t = D ∂ V∂x − εU ( x, t ) V, (10)with ε = κρ . Equation (10) needs to be fully solved inorder to gain a consistent picture of the spatio-temporalevolution of the action potential V ( x, t ). Unfortunatelythis equation is complex as it stands, and no exact so-lution can be obtained except via numerical simulations.Neverthless we remark that at steady state, this equationreduces to an eigenvalue problem for which exact analyti-cal solutions can be found for specific values of µ . To this last point, in steady-state regime eq. (10) reduces to thezero-eigenvalue linear operator problem: ϑ ( x ) − ∂ ∂x ! V ( x ) = 0 , ϑ ( x ) = − µ sech x, (11)in which µ = εD , where we used the simplified notation V ( x ) = V ( x, τ = tanh x , eq. (11) can betransformed into a Legendre equation of order n i.e. [37]: ddτ (cid:26) (1 − τ ) dVdτ (cid:27) + n ( n +) V = 0 , (12)where n ( n + 1) = µ , n being a positive integer. Boundedsolutions to the linear operator equation (11), for an ar-bitrary n , are the Legendre polynomials: V n ( τ ) = 12 n n ! d n dτ n ( τ − n , n = 1 , , , · · · . (13)For illustration, below we list the three lowest boundedmodes: V ( x ) = tanh x, D = 2 ǫ, (14) V ( x ) = 12 (3 tanh x − , D = 6 ǫ, (15) V ( x ) = 12 (5 tanh x −
3) tanh x, D = 12 ǫ, (16)which are sketched in fig. (1). In the next section we shallproceed to numerical simulations of the modified Hodgkin-Huxley equation (10), using the above three bounded modesas input profiles V ( x,
0) of the action potential (i.e. as ini-tial conditions).
The modified Hodgkin-Huxley equation (10) is an initial-value problem, as such it can be solved numerically using afinite-difference algorithm. In our case we adopt a finite-difference scheme that combines a central-difference ap-proximation for the time derivative and a forward-differenceapproximation for the second-order derivative in space[39]. As we are interested more in a qualitative analy-sis than a quantitative description of the problem, valuesof the diffusion coefficient D , the bare membrane capaci-tance C , the electromechanical coupling coefficient κ andthe quantity ǫ will be arbitrary.Graphs in figs. 2, 3 and 4 represent profiles of the trans-membrane voltage V ( x, t ) at six different times t , gener-ated numerically from the modified Hodgkin-Huxley equa-tion (10) for the three distinct initial conditions (14) (fig.2), (15) (fig. 3) and (16) (fig. 4). Values of parameters are D = 6, C = 2, κ = 0 . ǫ = 3 . Alain M. Dikand´e: On a nonlinear electromechanical model of nerve -1-0.5 0 0.5 1-10 -5 0 5 10 V x -0.6-0.4-0.2 0 0.2 0.4 0.6 0.8 1 1.2-10 -5 0 5 10 V x -1-0.5 0 0.5 1 -4 -2 0 2 4 V x Fig. 1.
Sketches of the first three bounded modes of the zero-eigenvalue equation (11). From left to right: n = 1, 2, 3. - - -
10 0 10 20 30 - V t = x - - -
10 0 10 20 30 - V t = x - - - -
10 0 10 20 30 - V t = x - -
20 0 20 40051015202530 x V t = - -
20 0 20 40 60020406080100120140 x V t = V t = Fig. 2. (Color online) Profiles of the transmembrane voltageat different times t , obtained from numerical simulations ofeq. (10) with the bound state eq. (14) used as initial solution: D = 6 . C = 2 . κ = 0 . ǫ = 3 . some propagation time, the kink-shaped input is modu-lated and stabilizes permanently in a typical pulse shapecharacteristic of the action potential [41,42,43]. Graphs offigs. 3 and 4 exhibit the same feature as fig. 2, meaningthat the two other initial solutions will also stabilize ina pulse shape similar to fig. 2 after a transient propaga-tion time. To gain a global view of the spatio-temporalevolution of the three different solutions shown in the pre-vious figures, they were represented in three dimensions asdepicted in figs. 5, 6 and 7. The three-dimensional repre-sentation clearly indicate that as they propagate along theaxon, the three different initial profiles always modulateinto the same pulse pattern. - - -
10 0 10 20 30 - V t = x - - -
10 0 10 20 300.00.20.40.60.81.01.21.4 x V t = - - - -
10 0 10 20 300.00.51.01.52.02.53.03.5 x V t = - -
20 0 20 40 6005101520 x V t = - -
20 0 20 40 60 80020406080100120140 x V t = V t = Fig. 3. (Color online) Profiles of the transmembrane voltageat different times t , obtained from numerical simulations of eq.(10) with the bound state eq. (15) used as initial condition: D = 6 . C = 2 . κ = 0 . ǫ = 3 . The mechanism by which the nerve impulse is generatedand transmitted along the axon has been a fundamentalproblem not only in neurophysiology, but also in mathe-matical physics[11,12,13,14]. In their pioneer model, Hodgkinand Huxley [1] suggested a picture based on an electriccable according to which the nerve impulse would be apropagating form of a depolarization, due to ion exchangesthrough the axon membrane. However in view of the unde-niable role of membrane in selectively passing ion from in-tracellular to extracellular fluids and vice-versa, Hodgkin-Huxley’s electric-cable picture was subsequently improved lain M. Dikand´e: On a nonlinear electromechanical model of nerve 5 - -
10 0 10 20 - - - V t = - -
10 0 10 20 - - - V t = - - - -
10 0 10 20 30 - - - V t = - -
20 0 20 40 600246810 x V t = - -
20 0 20 40 60 80051015202530 x V t = V t = Fig. 4. (Color online) Profiles of the transmembrane voltageat different times t , obtained from numerical simulations of eq.(10) with the bound state eq. (16) used as initial condition: D = 6 . C = 2 . κ = 0 . ǫ = 3 . Fig. 5. (Color online) Spatio-temporal shape of the transmem-brane voltage V ( x, t ), obtained numerically with the initialprofile eq. (14): ǫ = 2, κ = 0 . C = 2 . Fig. 6. (Color online) Spatio-temporal shape of the transmem-brane voltage V ( x, t ), obtained numerically with the initialprofile eq. (15): ǫ = 2, κ = 0 . C = 2 . Fig. 7. (Color online) Spatio-temporal shape of the transmem-brane voltage V ( x, t ), obtained numerically with the initialprofile eq. (16): ǫ = 2, κ = 0 . C = 2 . by regarding the membrane capacitance as a feedback or-gan [8]. Observations of thermodynamic phenomena in thenerve activity, such as the heat release during liquid-geltransition with subsequent generation of acoustic wavesalong the axon [7,30], motivated a distinct picture involv-ing mechanical processes related to a variable density dif-ference of liquids flowing through the nerve membrane.This new picture [7] led to the idea that pressure wavescould be a manifestation of the action potential. Alain M. Dikand´e: On a nonlinear electromechanical model of nerve
For the electrical [1] and mechanical [7] pictures, takenseparately, enable only partial descriptions of the processof nerve impulse generation, the need for a descriptiontaking simultaneously into consideration the electrical andmechanical activities of the cell membrane, was an imper-ative. In this study we exploited experimental evidencesof electromechanical phenomena, and their specific mani-festations in some neurophysioical contexts [25], to intro-duce a model which rests on the Hodgkin-Huxley electri-cal model, but assumes the membrane capacitance to bedetermined by the difference in densities of ion-carryingfluids flowing accross the membrane. The proposed modelcombines the KdV equation already present in the solitonmodel of Heimburg and Jackson [7], and the Hodgkin-Huxley’s electric cable equation with a feedback capacitor(i.e. a capacitive diode). By postulating a mathematicalexpression describing the relationship between the mem-brane capacitance and the density difference, we foundthat in steady-state regime the action potential equationreduces to a zero-eigenvalue linear operator problem. Thislinear operator problem can be transformed into the Leg-endre equation, the solutions of which is the family of Leg-endre polynomials. The three lowest bound states of thislinear operator equation were obtained, and used as initialprofiles in numerical simulations of the full partial differ-ential equation describing the spatio-temporal evolution ofthe transmembrane voltage. From numerical simulationsit turned out that the three distinct initial profiles alwaysdecay into the a common pulse profi;e after a transientpropagation time.To end, let us underline that the model proposed in thepresent study, can be improved to account several relevantaspects of the process that we neglected. For instance inthe modified Hodgkin-Huxley equation (2), we ingored thecontribution of ion currents F ( V ) yet this term plays animportant role the original Hodgkin-Huxley model [1,2].Also it is known that the Boussinesq equation describ-ing the dynamics of pressure waves, is actually obtainedfrom a Tyalor expansion of the fluid velocity field with re-spect to the density difference ∆ρ A (or U ). It is well estab-lished that carrying out the expansion beyond the linearterm, leads to higher-order or modified KdV equations [44]which admit soliton solutions distinct from eq. (9). There-fore including these quantities in the present model willundoubtedly enrich qualitatively the physics of the processunder study. In particular having distinct soltion solutionsfor the density wave implies distinct bound-state spectra,and hence new profiles for the action potential that mightpossibly be more close to the reality. Conflict of interest
The author declares that he has no conflict of interest.
References
1. A. L. Hodgkin and A. F. Huxley, J. Physiol. , 500 (1952)2. A. C. Scott, Rev. Mod. Phys. , 487 (1975) 3. R. Fitzhugh, Biophys. J. (1961)4. J. Nagumo, S. Arimoto and S. Yoshizawa, Proc. IRE ,2061 (1962)5. J. Engelbretcht, Proc. R. Soc. London , 195 (1981)6. I. Tasaki and G. Matsumoto, Bull. Math. Biol. , 1069(2002)7. T. Heimburg and A. D. Jackson, PNAS , 9790 (2005)8. A. M. Dikand´e and B. Ga-Akeku, Phys. Rev. E , 041904(2009)9. R. D. Keynes and D. J. Aidley, Nerve and Muscle , 3rd edn.(Cambridge University Press, Cambridge, Massachusetts,1982)10. R. R. Poznanski, L. A. Cacha, Y. M. S. Al-Wesabi, J. Ali,M. Bahadoran, P. P. Yupapin and J. Yunus, Scientific Re-ports , 2746 (2017)11. R. R. Poznanski and L. A. Cacha, J. Integ. Neurosc. ,417 (2012)12. R. R. Poznanski and J. Integ. Neurosc. , 267 (2004)13. R. R. Poznanski, Modeling in the Neurosciences: FromIonic Channels to Neural Networks (Harwood AcademicPublishers, Amsterdam, Netherlands, 1999)14. R. R. Poznanski, K. A. Lindsay, J. R. Rosenberg andO. Sporns,
Modeling in the Neurosciences: From BiologicalSystems to Neuromimetic Robotics , 2nd edn. (Taylor andFrancis, NW, USA, 2005)15. G. F. Achu, F. M. Moukam-Kakmeni and A. M. Dikand´e,Phys. Rev. E , 012211 (2018)16. J. Engelbrecht, T. Peets, K. Tamm, M. Laasmaa andM. Vendelin, Proc. Eston. Acad. Sc. , 28 (2018)17. K. S. Cole and H. J. Curtis, J. Gen. Physiol. , 649 (1939).18. E. N. Warman, W. M. Grill and D. Durand, IEEE Trans.Bio-Med. Electron. , 1244 (1992)19. W. M. Grill, IEEE Trans. Bio-Med. Electron. , 918(1999)20. J. K. Mueller and W. J. Tyler, Phys. Biol. , 01 (2014)21. T. Heimburg, Biochim. Biophys. Acta-Biomembr. ,147 (1998)22. I. Tasaki, K. Kusano and M. Byrne, Biophys. J. , 1033(1989)23. I. Tasaki and P.M. Byrne, Biophys J. , 633 (1990)24. R. Blunck and Z. Batulan, Frontiers in Pharmacology ,166 (2012)25. D. Gross, W. S. Williams, J. A. Connor, Cell. Mol. Neu-robiol. , 89 (1983)26. A. El Hady and B. B. Machta, Nature Commun. , 6697(2015).27. M. Mussel and M. F. Schneider, Scientific Report , 2467(2019).28. I. Sasaki, Ferroelectrics , 305 (1998).29. H. Chen, D. Garcia-Gonzalez and A. J´erusalem, Phys.Rev. E , 032406 (2019).30. T. Heimburg, Biochim. Biophys. Acta Biomem. , 113(2004)31. T. Heimburg, Thermal Biophysics of Membranes (Wiley-VCH Verlag Gmbh, 2007)32. J. Engelbrecht, T. Peets and K. Tamm, Biomech. Model.Mechanobio. , 1771 (2018).33. A. G. Petrov, Biochim. Biophys. Acta , 1 (2001)34. T. Heimburg, A. Blicher, L. D. Mosgaard and K. Zecchi,J. Phys.: Conf. Ser. , 012018 (2014)35. M. Plaksin, S. Shoham and E. Kimmel, Phys. Rev. X ,011004 (2014)lain M. Dikand´e: On a nonlinear electromechanical model of nerve 736. A. Kamkin, I.K. (eds.), Mechanosensitivity of the NervousSystem (Springer, Berlin, 2009)37. M. Abramowitz, I. A. Stegun,
Handbook of MathematicalFunctions , 10 th edn. (National Bureau of Standards, 1972)38. C. S. Gardner, J. M. Greene, M. D. Kruskal and R. M.Miura, Phys. Rev. Lett. , 1095 (1967)39. R. H. Landau, M. J. P´aez and C. C. Bordeianu, Compu-tational Physics (2 nd Revised and Enlarged Edition, Wiley,Wienheim, 2007).40. See e.g. M. Renganathan, H. Wei and Y. Zhao,
CardiacAction Potential Measurement in Human Embryonic StemCell Cardiomyocytes, for Cardiac Safety Studies Using Man-ual Patch-Clamp Electrophysiology . In: M. Clements and L.Roquemore (eds)
Stem Cell-Derived Models in Toxicology:Methods in Pharmacology and Toxicology (Humana Press,New York, NY, 2017).41. B. C. Carter and B. P. Bean, Neuron (Cell Press) , 898(2009).42. Y Yao, C. Su and J. Xiong, Physica A , 121734 (2019).43. A. E. Casale and D. A. McCormick, J. Neurosc. , 18289(2011).44. See e.g. J. Engelbrecht, K. Tamm and T. Peets, Biomech.Model. Mechanobio.14