• 検索結果がありません。

不均一性を反映した構成式による 45 ◦ 数珠繋ぎ構造モデルの解析

第 5 章 分子鎖の不均一性を反映させた構成式の構築 63

5.2 不均一性を反映した構成式による 45 ◦ 数珠繋ぎ構造モデルの解析

(a) Nominal stress-stretch relations (b) Hysteresis loss

Hysteresis loss [J/m3 ]

0 0.05 0.1 0.15

Nominal Stress Σn22[MPa]

Stretch λ

2

Experiment

Homogeneous model

Pseudo-heterogeneous model

1 1.2 1.4

0 0.5 1 1.5 2

Fig.5.7 Comparison of(a) Nominal stress-stretch relations and (b) Hysteresis loss of silica-filled rubber by simulation and experiment (45bunching structure model).

前節で作成したゴムマトリクス部およびゲル相の構成式を用いて45数珠繋ぎ構造 モデルの変形シミュレーションを行った時の(a)公称応力―ストレッチ関係,(b) ヒス テリシスロスを図5.7に示す.図中黒線は3章の結果で,赤線は実験結果である.不均 一性を導入するとゴム相・ゲル相どちらも最大応力・ヒステリシスロスともに増加す るのでいずれも3章のそれより増大している.このため実験結果からさらに外れた結 果となっているが,45数珠繋ぎ構造と実際のシリカ粒子のランダム分布の差もあるた めここでは問題としない.

図5.8〜図5.10に最大引張(λ2 = 1.5)における分子鎖ストレッチλc分布,引張方向 応力σ22分布,セグメント数の変化量∆N分布をそれぞれ示す.図では3章の結果と比 較して示している.図中上部に記載している数字はユニットセル内のゴムマトリクス 部(R),ゲル相(G)それぞれにおけるλc,σ22,∆N の平均である.図5.8の分子鎖ス トレッチλc分布は不均一性を考慮した場合もほとんど変化がないが,図5.9の引張方

向応力σ22分布では,粒子を繋ぎ止めているゲル相の変形集中領域とその周囲のゴム マトリクス部で分子鎖の配向硬化が顕著になり応力が上昇している.一方,図5.10の セグメント数の変化量∆N分布では,ゴム部の∆ ¯NRは変わらないがゲル相の平均値

∆ ¯NGは僅かに増加している.これは,4.2.2節で述べたゲル相中の軟らかい部分(セ グメント数Nの大きい部分)での局所的な非アフィン変形によってゲル相全体の平均 セグメント数変化量∆ ¯Nが増加しているためである.このような変形集中部のゲル相 における顕著な応力σ22の上昇とセグメント数の変化量∆N の増大が,図5.7(b)に示 すヒステリシスロスの増大をもたらしたと考えられる.

1.66 1.33 1.00 1.99 2.32 2.65 λc

(a)Homogeneous model (b)Pseudo-heterogeneous model λRc=1.23

λGc=1.49

λRc=1.23 λGc=1.48

Fig.5.8 Comparison of distributions of morecular chain stretch λc (45 bun-ching structure model).

3.28 1.66 0.04 4.91 6.53 8.15 σ22[MPa]

(a)Homogeneous model (b)Pseudo-heterogeneous model σR22=1.38

σG22=3.54

σR22=1.58 σG22=3.89

Fig.5.9 Comparison of distributions of tensile stressσ22 (45 bunching struc-ture model).

5.24 2.62 0.00 7.86 10.49 13.11

N

(a)Homogeneous model (b)Pseudo-heterogeneous model

NR=0.12

NG=2.59

NR=0.11

NG=2.78

Fig.5.10 Comparison of distribution of change in the segment ∆N (45 bun-ching structure model).

6 結論

シリカ充填ゴムは,CB充填ゴムとは異なる粒子分布形態やゴム部と界面との相互作 用により幅広い力学特性を発現する可能性を秘めている.本研究では,実験により確 認されたシリカ充填ゴムの数珠繋ぎ構造についての基礎的な検討として,種々の実験 事実を反映させたシミュレーションモデルの構築,ならびに充填構造やゲル相の物性 が巨視的な力学的特性に与える影響を検討した.以下に得られた結果をまとめて示す.

第2章では分子鎖網目理論に基づく粘弾性8鎖モデルを概説し,ゴム粘弾性体の構 成式について説明した.さらに,粘弾性構成式を絡み点数が変形とともに変化するこ とを許容する非アフィン分子鎖網目理論により一般化し,更新ラグランジュ法に基づ く有限要素均質化法を定式化した.

第3章では,数珠繋ぎ構造を呈するシリカ充填ゴムの局所的な構造を模擬したユニッ トセルモデルを作成し,粒子の配置形態,ならびにゲル相の特性変化がシリカ充填ゴ ムの応答に及ぼす影響を検討した.ゴム相,ゲル相の実験に基づきモデル化した構成 式を用いて,ユニットセル中の2個のシリカ粒子の配置を0,30,45に変化させた 系で変形シミュレーションを行った.ただし,ゲル相にはセグメント数変化を生じな いアフィンモデルを適用している.ユニットセル中に直列(0)にシリカ粒子を配置 した場合,ユニットセルの中心軸付近の粒子を繋ぎ止めているゲル相に変形が集中し て計算不能となったが,シリカ粒子配置が30,45の系では,系全体の変形λ2に応じ てシリカ粒子は位置を左右にずらすことが可能であるため応力上昇は緩やかになり,1 サイクルの応答を得ることが可能であった.変形後期に強い配向硬化を示すアフィン ゲル相モデルの導入により,実験結果にみられる変形後期の顕著な応力上昇(S字状 の応力ーストレッチ曲線)も再現可能であることを示した.ただし,アフィンモデル を仮定した解析ではヒステリシスロスに変化がない.そこで,非アフィンモデルによ

るゲル相のモデル化を行い,数珠繋ぎ構造の応答に与える影響を評価した.初期セグ メント数Nsの同じアフィンゲル相での解析結果と比較すると,変形後期においてゲ ル相内のセグメント数Nsが増加し,分子鎖の配向硬化が弱くなるため系の最大応力 は低下するが,除荷時の応力低下も大きくなりヒステリシスロスが増大することが示 された.さらに,ゲル相の初期セグメント数Nsの効果について検討し,初期セグメ ント数Nsが小さいほど負荷時の応力上昇挙動ならびに負荷から除荷に移行した際の 応力低下がより顕著になり,ヒステリシスロスの増大をもたらすことを示した.これ により,非アフィンゲル相モデルを導入した数珠繋ぎ構造モデルにおいて,Nsを小さ くすることにより実験結果に見られるカップリング剤の増加に伴う変形抵抗とヒステ リシスロスの増大の傾向を再現可能であることを示した.

第4章では,実験により報告されているゴムを硫黄架橋した際に生じる絡み点の不 均一性について,ゴム相・ゲル相の各有限要素メッシュにセグメント数N をランダム に変化させたモデルにより検討した.いずれも,不均一分布を導入することで,変形 しにくい硬い部分(セグメント数Nの小さい部分)が,系の中で変形を担うゴム・ゲ ルの見かけの体積減少(=ゴムにCBやシリカ等の弾性体粒子を加えた時と同様の効 果)をもたらすため,最大応力が上昇する.ゲル相ではさらに,変形しやすい軟らか い部分(セグメント数Nsの大きい部分)への変形集中・非アフィン変形の発生によ り,局所的な∆N の値が増大し,除荷時の応力を低下させヒステリシスロスが増大す る.さらに,ゲル単相の平均セグメント数N0sを小さくする(硬いゲル相を想定する)

と,3章で述べた応力上昇・ヒステリシスロスの増大が,局所的な不均一性によって強 められることを示した.

第5章では,第4章の解析から得られた絡み点の不均一性をもつゴム相・ゲル相の 応答を,粘弾性8鎖モデルによる構成式で再度パラメータフィッティングしてスケー ルアップを行い,第3章の45数珠繋ぎ構造の解析を行って,絡み点の不均一性がシリ カ充填ゴムの応答に及ぼす影響を評価した.その結果,不均一性の導入によりゴム相・

ゲル相共に分子鎖の配向硬化が顕著になり,第3章で解析した不均一性をもたないモ デルの結果よりも最大応力・ヒステリシスロスともに大きくなることを示した.

本研究で用いたシリカ充填ゴムの直鎖状数珠繋ぎ構造モデルは,実際のシリカ充填 ゴムの局所的な一部分を取り出し単純化を行ったモデルであるため,実際のシリカ充 填ゴムの自由度を大きく制限したものとなっており,実験結果に見られる以上に顕著 な応力上昇および大きなヒステリシスロスを示している.実験結果における変形後期

の顕著な応力上昇挙動は,シリカ粒子の凝集構造の変化に起因する可能性も否定でき ない.今後,実際のシリカ充填ゴムに見られる巨視的初期等方性応答を呈する,より 大きなスケールのモデルによる検討が待たれる.

参 考 文 献

(1) 深堀美英, 設計のための高分子の力学, (2000), 技報堂出版.

(2) Mullins, L., Thixotropic Behavior of Carbon Black in Rubber,Rubber Chemistry and Technology, Vol.23, (1948), pp.281-300.

(3) Bergst¨om, J.S. and Boyce, M.C., Constitutive Modeling of the Large Strain Time-Dependent Beahvior of Elastomers, Journal of Mechanical Physics and Solids, Vol.46, No.5(1998), pp.931-954.

(4) 住友ゴム工業株式会社, TEM画像.

(5) 柴山充弘,日本ゴム協会誌, 84-1,(2011),14-20.

(6) 中島健, 岩蕗仁, 浦部匡史, 伊藤万喜子, 藤波想, 西敏夫, 日本ゴム協会誌, 85-7,(2012),46.

(7) D. Wang, S. Fujinami, K. Nakajima, S. Inukai, H. Ueki, A. Magario, T. Noguchi, M. Endo, T. Nishi, Visualization of nanomechanical mapping on polymer nano-composites by AFM force measurement, Polymer, Vol.51, (2010), pp.2455-2459.

(8) 池田裕子, シンクロトロンX線による架橋天然ゴムの伸長結晶化解析, 日本ゴム 協会誌, 84-1,(2011),pp.29-36

(9) Gao, J. and Weiner, J. H., Computer Simulation of Viscoelasticity in Polymer Melts, Macromolecules,25, (1992), 1348-1356.

(10) Ananyo, Bandyopadhyay, Pavan, K. Valavala, Thomas, C. Clancy, Kristopher, E. Wise. and Gregory M. Odegard., Molecular modeling of crosslinked epoxy po-lymers: The effect of crosslink density on thermomechanical properties,Polymer, 52, (2011), 2445-2452.

(11) Ying, Li, Martin, Kroger. and Wing, Kam, Liu., Primitive chain network study on uncrosslinked and crosslinked cis-polyisoprene polymers,Polymer,52, (2011), 5867-5878.

(12) Yashiro, K., Naito, M., Minagawa, Y. and Tomita, Y., On the Hysteresis of Polyethylene and Polybutadiene: A Molecular Dynamics Study, Proceeding of ICCES2005, Advances in Computational and Experimental Engineering and Science, CDROM

(13) 屋代如月, 内藤正登, 皆川康久, 冨田佳宏, 非晶性高分子材料のヒステリシスに関 する分子動力学的研究,日本機械学会論文集,A72, (2006), 277-284.

(14) O’Brien, J. et al., An NMR Investigation of the Interaction between Carbon Black and cis-Polybutadiene, Macromolecules, Vol.9, No.4(1976), pp.653-660.

(15) . Gert Heinrich, Manfred Kl¨uppel, Thomas A. Vilgis., Reinforcement of ela-stomers, Current Opinion in Solid State and Materials Science, Vol.6, (2002), pp.195-203.

(16) A. R. Payne., The Dynamic Properties of Carbon Black-Loaded Natural Rubber Vulcanizates. Part I , Journal of Applied Polymer Science, Vol.6, (1962), pp.57-63.

(17) J. A. C. Harwood, L. Mullins, and A. R. Payne., Stress Softening in Natural Rubber Vulcanizates.Part II . Stress Softening Effects in Pure Gum and Filler Loaded Rubbers, Journal of Applied Polymer Science, Vol.9, (1965), pp.3011-3021.

(18) G. Kraus, C. W. Childers, K. W. Rollmann., Stress Softening in Carbon Black Reinforced Vulcanizates. Strain Rate and Temperature Effects,Journal of Applied Polymer Science, Vol.39, (1966), pp.1530-1543.

(19) Dannenberg, E.M., The Effects of Surface Chemical Interactions on the Properties of Filler-Reinforced Rubbers, Rubber Chemistry and Technology, Vol.48, (1975), pp.410-444.

(20) F. Bueche., Molecular Basis for the Mullins Effect, Journal of Applied Polymer Science, Vol.4, (1960), pp.107-114.

(21) G. Marckmann, E. Verron, L. Gornet, G. Chagnon, P. Charrier, P. Fort., A theory of network alteration for the Mullins effect,Journal of the Mechanics and Physics of Solids, Vol.50, (2002), pp.2011-2028.

(22) 内藤陽子,伊藤眞義, シリカ充てん加硫ゴム材料の応力軟化に関する研究,日本ゴ ム協会誌, 82-9,(2009),pp.394-399

(23) 古谷泰大, 陸偉, 内藤正登, 冨田佳宏, 分子鎖網目モデルによるカーボンブラック 充填ゴムの粘弾性変形挙動シミュレーション, 日本機械学会第17回計算力学講演 会講演論文集, No.04-40(2004-11), pp.17-18.

(24) 古谷泰大, 内藤正登, 陸偉, 冨田佳宏, カーボンブラック充填ゴムの繰り返し変形 挙動の評価, 日本機械学会論文集A編, Vol.71, No.708(2005), pp.1109-1115.

(25) Tomita,Y.,Lu,W.and Furutani,Y., Micro- to Macroscopic Deformation Beha-vior of Carbon Black-Filled Rubber Under Monotonic and Cyclic Straining, Proc.CIMTEC2004,PartB,(2004),pp.121-132.

(26) Tomita,Y.,Azuma,K.,and Naito,M.,Strain-Rate-Dependent Deformation Beha-vior of Carbon-Black-Filled Rubber under Monotonic and Cyclic Straining, Key Engineering Materials,Vol.340,No.341(2007),pp.1017-1024.

(27) Wolff,S., Silanes in Tire Coupling After Ten Years - A Review, Tire Society and Technology,Vol.15,No.4(1987),pp.276-294.

(28) 松沢憲治,シリカ/シランフィラーシステムの化学とゴム補強性, 日本ゴム協会 誌,Vol.78,No.6(2005),pp.211-217.

(29) 北村真瑠久, 神戸大学修士論文, シリカ充填ゴムの粘弾性挙動のモデル化と変形 挙動の評価, (2010)

(30) 望月利紀, 神戸大学修士論文, シリカ充填ゴムの内部微視構造のモデル化と変形 応答のシミュレーションによる評価, (2011)

(31) 田中文彦, 高分子の物理学, (1994), 裳華房.

(32) 冨田佳宏, ガラス状ポリマーの分子鎖網目理論による構成式と変形挙動のシミュ レーション, 塑性と加工, Vol.37, No.424(1996), pp.485-491.

(33) Kuhn, W. and Grun, F., Beeziehuugen Zwischen Elastischen Konstanten und Dehuungsdoppelbrechung Hochelastischer Stoffe, Kollooidzeitschrift, Vol.101, (1942), pp.248-271.

(34) 土井正男, 小貫明,高分子物理・相転移ダイナミクス, (1992), 岩波書店.

(35) James, H.M. and Guth, E., Theory of the Elastic Properties of Rubber, J.Chem.Phys., 11 (1934), 455-481.

(36) James, H.M. and Guth, E., Theory of the elastic properties of rubber,Journal of Chemical Physics, Vol.11, (1943), pp.455-481.

(37) Wang, M.C. and Guth, E., Statistical Theory of Networks of Non-Gaussian Frex-ible Chains, Journal of Chemical Physics, Vol.20, No.7(1952), pp.1144-1157.

(38) Treloar, L.R.G., The Elasticity of a Network of Long-Chain Molecules-III, Tran-sactions of the Faraday Society, Vol.42, (1946), pp.83-94.

(39) Flory, P.J. and Rehner, J., Statistical Mechanics of Cross-Linked Polymer Net-works, Journal of Chemical Physics, Vol.11, No.11(1943), pp.512-520.

(40) Arruda, E.M. and Boyce, M.C., A Three-Dimensional Constitutive Model for the Large Stretch Behavior of Rubber Elastic Materials, Journal of Mechanical Physics and Solids, Vol.41, No.2(1993), pp.389-412.

(41) Green, A.E. and Zerna, W., Theoretical Elasticity, (1968), Oxford University Press.

(42) Tomita,Y.,Azuma,K.and Naito,M.,Int.J.Mech.Sci., Vol50, No.5,(2008),pp856-868.

(43) Truesdell, C. and Noll, W., The Non-Linear Field Theories of Mechanics, Springer-Verlag, (1965).

(44) Boyce, M.C., Weber, G.G. and Parks, D.M., On the Kinematics of Finite Strain Plasticity, International Journal of Plasticity, Vol.37, No.5(1989), pp.647-665.

(45) 北川浩, 弾・塑性力学−非線形解析のための基礎理論−, (1987), 裳華房.

(46) 久田俊明, 非線形有限要素法のためのテンソル解析の基礎, (1992), 丸善.

(47) Boyce, M.C., Parks, D.M. and Argon, A.S., Inelastic Deformation of Glassy Poly-mers. Part I :Rate Dependent Constitutive Model,Mechanics of Materials, Vol.7, (1988), pp.15-33.

(48) Doi, M. and Edwards, S.F., The Theory of Polymer Dynamics, (1986), Oxford University Press.

(49) Gennes, P.G. de., Reptation of a Polymer Chain in the Presence of Fixed Obsta-cles, Journal of Chemical Physics, Vol.55, (1971), pp.572.

(50) Raha, S. and Bowden, P.B., Birefringence of Plastically Deformed Polymethyl Methacrylate, Polymer, Vol.13, (1972),pp.174-183.

(51) Botto,P.A., Duckket,R.A. and Ward,I.M., The Yield and Thermoelastic Proper-ties of Oriented Polymethyl Methacrylate,Polymer, Vol.28, (1987),pp.257-262.

(52) Tomita, Y. and Tanaka, S., Prediction of Deformation Behavior of Glassy Poly-mers Based on Molecular Chain Network Model,International Journal of Solids and Structures, Vol.32, No.23(1995), pp.3423-3434.

(53) Tomita, Y., Adachi, T. and Tanaka, S., Modelling and Application of Constitu-tive Equation for Glassy Polymer Based on Nonaffine Network Theory,European Journal of Mechanics, A/Solids, Vol.16, No.5(1997), pp.745-755.

(54) Higa, Y. and Tomita, Y., Computational Prediction of Mechanical Properties of Nickel-based Superalloy with Gamma Prime Phase Precipitates, Proceedings of ICM8,Advance Materials and Modeling of Mechanical Behavior, Vol.III, (1999), pp.1061-1066.

(55) 比嘉吉一,冨田佳宏,粒子強化型複合材の均質化法による変形挙動のモデル化とシ ミュレーション, 日本機械学会論文集A編, Vol.66, No.648(2000), pp.1441-1446.

(56) 冨田佳宏, 数値弾塑性力学, (1990), 養賢堂.

(57) Guedes, J.M. and Kikuchi, N., Preprocessing and Postprocessing for Materials Based on the Homogenization Method with Adaptive Finite Element Methods,

Computer Methods in Applied Mechanics and Engineering., Vol.83, No.2(1990), pp.143-198.

(58) 内田真, 冨田佳宏, 非晶性高分子材料の微視組織の 不均一性が変形挙動に及ぼす 影響,日本機械学会論文集A編, 70-698, (2004), 1484-1490.

(59) 内藤正登,東圭佑, 冨田佳宏, フィラー充填ゴムの粘弾性挙動のモデル化と数値シ ミュレーション,日本機械学会第19回計算力学講演会講演論文集, No.06-9(2006-11), pp.397-398.

付録 A

非圧縮性ゴム粘弾性体の構成式の速 度形式表示

式(2.10)の速度形を導出する.

σi = 1 3

{ CαR

NαL−1 ( λc

√Nα )

+CβRNβ 1

λγL−1 (

λβ

Nβ )}

λ2i λc +1

3 {

CαBR

NαBL1

( λcB

√NαA )}λ2i

λc

−p (A.1)

式(A.1)を左Cauchy-Green変形テンソルAijを用いて表すと,

σij = 1 3

{ CαR

NαL1 ( λc

√Nα

)

+CβRNβ 1

λγL1 (

λβ

Nβ )}

λ2i λc

+1 3

{ CαBR

NαBL1

( λcB

√NαA )}λ2i

λc

−pδij (A.2) となる.続いて,下記の関係式

σ= ˙ σW σ+σW, I˙1 = (trA)·= 2A·D, A˙ = (D+W)A+A(DW) を用いることにより,式(A.2)の速度形式表示は次式になる.

σij = 1 3

[{

CαRNα

( ζ

√Nα L λc

)

+ CβRNβ λγ

( ζ λγ

Nβ L λc

)}

AijAkl/Amm

+ {

LCαR Nα

λc +LCβRNβ λc

}

ikAjl+Aikδjl} ]

˙ εkl

CβRNβλ˙γ λ2γ

3Amm (

L+ λβζ

Nβ )

Aij +1 3

[{

CαBRNαB

( ζ′′

√NαB L′′

λcB )}

AijAkl/Amm

+ L′′CαBR NαB λcB

{δikAjl+Aikδjl}]

( ˙εkl−ε˙pkl)−pδ˙ ij (A.3) ただし,σij はCauchy応力のJaumann速度で,I1 は左Cauchy-Green変形テンソル Aij の第1主不変量で,I1 =λ21+λ22+λ23となる.Dは変形速度テンソルで,W はス ピンテンソルである.trはトレースを表す.

最後に,体積一定の条件下で,式(A.3)に示すCauchy応力のJaumann速度σijを Kirchhoff応力のJaumann速度Sijで置き換えても本質的な差はないことと,変形速度 テンソルDをひずみ速度テンソルε˙で置き換えることにより,非圧縮性ゴム弾性体の 構成式の速度形式表示は次式になる.

Sij = 1 3

[{

CαRNα

( ζ

√Nα L λc

)

+CβRNβ λγ

( ζ λγ

Nβ L λc

)}

AijAkl/Amm

+ {

LCαR Nα

λc + LCβRNβ λc

}

ikAjl+Aikδjl} ]

˙ εkl

CβRNβλ˙γ

λ2γ 3Amm

(

L+ λβζ

Nβ

)

Aij + 1 3

[{

CαBRNαB

( ζ′′

√NαB L′′

λcB )}

AijAkl/Amm

+ L′′CαBR NαB λcB

{δikAjl+Aikδjl}]

( ˙εkl−ε˙pkl)−pδ˙ ij (A.4)

ここで,

N は分子鎖の限界伸び比を表す.ζ = dxdL1(x) x=

I1 3N

= 1β2βcsch2 2β である.

付録 B

[ϕ],[B],[E], { ψ } の具体形

P(x,y) A e1 A e2

A e3

3(x ,y ) 3 3

2(x ,y ) 2 2

1(x ,y ) 1 1

Fig.B.1 Triangle Element

要素は図B.1に示すような三角形1次要素を用いるので,形状マトリクス[ϕ]及び {ψ}は次のようになる.

{T}={ψ}T{θ} , {ψ}T =1 ψ2 ψ3}T (B.1) {v}=

{vx

vy }

= [ϕ]˙} (B.2)

[ϕ] = [ [ϕ1][ϕ2][ϕ3] ] (B.3)

i] =

[ψi 0 0 ψi

]

(B.4)

ここで,ψiは形状関数である.図B.1のように節点1,2,3の座標をそれぞれ(x1 , y1), (x2 , y2),(x3 , y3)とし,要素内の任意の点Pの座標を(x, y)とすると,全体の面積Ae, Ae1,Ae2,Ae3 は次のように表せる.

2Ae = det

1 x1 y1 1 x2 y2 1 x3 y3

, 2Ae1 = det

1 x y 1 x2 y2 1 x3 y3

 (B.5)

2Ae2 = det

1 x y 1 x3 y3 1 x1 y1

, 2Ae3 = det

1 x y 1 x1 y1 1 x2 y2

 (B.6)

上式を用いて形状関数ψ123 は次式のようになる.

ψ1 = Ae1

Ae , ψ2 = Ae2

Ae , ψ3 = Ae3

Ae (B.7)

次にマトリクス[B],[E]は平面問題では次のようになる.

˙}=



˙ εxx

˙ εyy 2 ˙εxy



= [B]˙} (B.8)

[B] = [ [B1][B2][B3] ] (B.9) [Bi] =

ψi,x 0 0 ψi,y ψi,y ψi,x

 (B.10)

{q}=







vx,x

vy,y vx,y vy,x







= [E]˙} (B.11)

[E] = [ [E1][E2][E3] ] (B.12)

[Ei] =



ψi,x 0 0 ψi,y ψi,y 0

0 ψi,x



 (B.13)