接触・摩擦

2014年2月20日 (木)

ゴムバンドの摩擦力による動力伝達(テンショナーの追加)

これまでゴム部品の解析事例として、主にゴムバンドを用いた摩擦力による動力伝達をテーマとして進め、2次元モデルによる動力伝達、および3次元モデルによるクラウン効果の検証について、Slome-Mecaによる解析を行う事ができました。

それらの解析は、駆動軸と従動軸という2つの軸しか考慮していませんでしたが、今回は、2次元モデルに、更に張力を保持するためのスプリング式のテンショナーを追加しました。

まずは、結果の動画からご覧下さい。左側が駆動ローラー(Driving roller)、右側が従動ローラー(Driven roller)です。

テンショナーの動きを確認するため、負荷トルクを変動させています。

駆動ローラー(Driving roller)はTime=1で1回転しますが、回転スタート時にはトルクがゼロで、そこから直線的に増加して、半回転後にトルク最大になり、再び直線的に減少して1回転後にゼロに戻ります。この事例では駆動ローラーを2回転しているので、2サイクルの負荷変動が有ります。

Rband1201

スプリング式のテンショナーは、ゴムバンドが負荷によって緩む側のスパンに付けられています。先端のアイドラー(抵抗なしで空回りするローラー)はバネ力によってゴムバンドに押し付けられていますが、負荷の変動によって上下方向の位置が変化します。また、それに伴い、バネ力も変化しています。

負荷がゼロの状態です。
(ゴムバンドの色分けはMises応力、アイドラーの色分けは上下変位です)。
Rband1202

負荷最大の状態です。
Rband1203

負荷最大時には、上側スパンと下側スパンの応力の差が大きくなる事がわかります。また、テンショナー線先端のアイドラーの変位は見ただけではよく分かりませんが、色を見ると、上側へ移動しているようです。つまり、負荷が大きくなると下側スパンが緩もうとするため、その分テンショナーのバネが伸びてアイドラーの位置が上に移動するということです。

解析結果から、負荷トルクと、アイドラーの上下変位、テンショナーのバネ力をグラフにしました。

Rband1204

負荷トルクが大きくなると、ゴムバンド下側スパンが緩み、バネが伸びてバネ力(スプリング力)は弱くなり、アイドラーの位置は上に移動していることがわかります。

なお、ゴムバンドとアイドラーの間に摩擦力はありません。

##########################################################

さて、以降はこの解析モデルの説明です。

まず、駆動および従動ローラー、ゴムバンドの形状・寸法や物性値、張力や負荷トルクの与え方などは従来と同じですので省略します。従来の解析は下記URL「ゴムバンドの摩擦力による動力伝達(前編)」から始まる一連の記事をご参照下さい。

http://salome-meca.cocolog-nifty.com/blog/2012/12/post-e977.html

今回の追加部分であるスプリング式のテンショナーですが、下図の通り、「ベース」、「バネ(スプリング)」、「アイドラー」の3つの部品から構成されています。

Rband1205

これらはSalomeで別のパーツとして作成し、メッシュも別々に作成しています。そのため、各パーツ間での共有節点はありません。

Rband1206

そのため、下図のように各部品の部分にグループ名を設定し、Code_Asterで「LIAISON_UNIF(自由度の連結)」を使って、上下方向(Y方向)を連結しています。LIAISON系列は大変形・大回転(GROT_GDEP)には対応していないようなのですが、並進移動だけなら問題無いようなので使っています。これについては「LIAISON_SOLIDEによる剛体は大回転出来ない(http://salome-meca.cocolog-nifty.com/blog/2012/10/liaison_solide-.html)」をご参照下さい。

Rband1207

Rband1208

また、これらの部品は全て横方向(X方向)には動かないので、上のように各点のX方向変位は拘束しています。

要素は、ベースとアイドラーはゴムバンドやローラーと同じD_PLAN(平面ひずみ要素)で、バネは2D_DIS_Tとしました。バネ定数はK_T_D_Lで設定しています(各方向、0.1[N/mm])。

なお、アイドラー部分は空回りするローラーなので、本来はゴムバンドによって回転させられます。しかし今回はアイドラーの負荷が無いため、まじめに回転させる必要は無く、ゴムバンドとの摩擦をゼロにした上で固定してしまっています。つまり、ゴムバンドはアイドラーの表面で滑ります。そのほうが設定が簡単で、おそらく計算時間も短縮できます。

解析手順ですが、従来のモデルと同様に、「(1)ローラーに軸荷重を与えて。ゴムバンドを張る」、「(2)従動ローラーにトルクを与えながら、駆動ローラーを回転させる」という2つの解析(2つのcommファイル)に分けています。具体的には、従来の事例をご参照下さい。

そのうち、(1)の段階(TIme=0~1)にて、ベースの下辺グループ「basel」に上向きに移動する変位を与えてテンショナーを移動させ、ゴムバンドに押し付けています。

Rband1209

Rband1210

また、接触設定については、アイドラー外面とゴムバンド外面のZONE設定を追加しています。ゴムバンド外面には「bandc2」というグループを追加しています。

Rband1211

「LISSAGE="OUI"」は接触面の法線を滑らかにする役目があるようですが、これを入れておいたほうが収束性が良さそうな感じだったので、とりあえず入れています。

(1)に続いて、ローラーを回転させる(2)の段階ですが、これは従来どおり別のcommファイルで「POURSUITE」を用いて継続計算させています。

ただし前述したようにトルクを変動させるように変更しています(従来は、1/4回転で所定値までトルクを上げた後は一定)。またトルクの最大値も従来の2倍に増やしました(従来は1.16[N・m])

Rband1212

後は従来どおり駆動ローラーを回転させるだけなのですが、どうも(2)のJOBの最初のステップで計算が収束しませんでした。どうやら接触問題を2つの「STAT NON LINE」に分けて継続計算する場合、接触状態の引継ぎが上手くいきにくい傾向にあるようです。これは、2つのJOB(commファイル)に分けてJOB自体を分けた場合でも、1つのJOB内で「STAT NON LINE」を分けた場合でも同様のようです。

実は、(1)のテンショナーを押し付ける動作も、当初はゴムバンドを張った後の後続として「STAT NON LINE」を分けて実施しようとしたのですが収束せず、接触設定のペナルティー係数を下げるなど色々試行錯誤しても収束しなかったので、ゴムバンドを張る段階と同時進行としました。

一方、(2)の段階の最初のステップで収束しない問題については、接触設定のペナルティー係数「COEF_PENA_CONT」および「COEF_PENA_FLOT」を下げる事で解決出来ました。これらは(1)の段階では両方とも500だったのですが、(2)の段階では、CONTの方を200、FLOTの方を400まで下げています。

Rband1213

上記の(1)でテンショナーを押し付ける方についてはペナルティー係数の調整でも上手くいかなかったのですが、なぜ(2)の方では上手くいったのかわかりません。このあたりはなかなか奥が深くて難しそうです。

あと、同じく初回ステップの収束性向上に効果が見られる「SEUIL_INIT」も0.1に設定しています。

######################################################

今回の解析は、仮想PC上のDebian7.2.0に、最近リリースされたばかりのSalome-Meca2014.1を入れて実施しています。

Salomeは7.3.0で言語を日本語に設定しています。

Code_Asterは12.1を用いました。
従来の11からの変更点として、STAT_NON_LINEからCOMP_INCRとCOMP_ELASの設定が消え、「COMPORTEMENT」になっています。記述する内容は同様のようなのですが、COMP_INCRとCOMP_ELASの区別をしなくて良くなった分、分かり易くなったのではないかと思います。

Rband1214

######################################################

今回のファイルは下記にzipで圧縮しています。

rband12.zip

rband12.mmed:メッシュファイル
rband12a.comm:前半(ゴムバンドを引っ張るまで)
rband12b.comm:後半(ローラーを回転させる)

2013年9月 7日 (土)

ゴムバンドの摩擦力による動力伝達(クラウン効果について)

前回、ゴムバンドの摩擦力による動力伝達(3次元)を実施し、ゴムバンドが回転につれてローラーの太い側へ寄っていく「クラウン効果」が再現できる事を確認しました。このクラウン効果についてもう少し考えてみます。

駆動ローラー1回転毎のバンドの動きをGIFアニメーションにしたものを、もう一度下に示します。Time=1で駆動ローラーが1回転しますが、その間にゴムバンドが奥側へ寄っていくのが分かります。

ゴムバンドがどのようにしてローラーの太い側へ寄っていくのか確認するため、下図のように、回転始めにおいてゴムバンド上の1点(点A)を指定しました。

Rband08201


この点AのZ座標(ローラー軸方向)が、回転につれてどう変化するかプロットしました。
また、軸荷重は前回実施の30[N]と、今回新たに50[N]でも実施して付け加えました。

Rband08202


ローラーの外径が約100mm、ゴムバンドの内径が約220mmなので、駆動ローラーが約2.2回転でゴムバンドが1回転する事になります。

そのゴムバンド1回転の間を、上図中の点線でa~eの領域に区切っています。点Aは回転始めに駆動ローラーの中央(領域a)に位置しており、1回転してほぼ同じ位置(領域e)に戻ります。

このプロットから、ゴムバンド上の点はローラー内(領域a,c,e)においてはローラーの細い側へ若干移動し、ローラー間(領域b,d)においてローラーの太い側へ大きく移動していることがわかります。そしてゴムバンド1回転トータルではローラーの太い側へ移動しています。

これを上側から見た図に描くと、下図のようになります(フリーハンドで描いたので、見難くてすみません)。点線は下側(駆動ローラーから従動ローラーへ向かう)、実線は上側(従動ローラーから駆動ローラーへ向かう)です。ゴムバンド上の点は8の字を描くようにしてローラーの太い側へ寄っていきます。

Rband08203_2


何故このような動きをするかというと、ゴムバンドがローラーに進入する際、ローラー軸に対して直角ではなく、直角からある程度の角度傾いて進入するためです。

上図の点線の円で囲んだ部分に注目してください。

Rband08204


A点は、ローラー間ではローラーの太い側へ向かって進んでいますが、ローラーに進入する直前に細い側へ進む向きを変えます。これはローラー面に角度が付いているためです。このため、上図の通りローラー軸直交方向に対してαの角度だけ傾いて進入します。このαを進入角と言います。この挙動はプロットのグラフでも確認できます。

A点はローラー上においてローラーの細い側へ若干滑りますが、その量はそれほど大きくはありません。そして進入角があるため、A点のすぐ後ろの後続点は、A点よりもローラーの太い側において進入します。更に次の後続点はもっと太い側で進入します。このようにして、ローラーへ進入する位置がだんだん太い側へ寄っていくので、ゴムバンド全体が回転につれて太い側へ寄ります。これがゴムバンド移動のメカニズムになります。

また計算結果では、ゴムバンドを張る際の軸荷重が大きいほうが、ゴムバンドの移動速度が遅い事が分かります。これは、ゴムバンドが強く引っ張られると進入角が小さくなるためです。

このように、クラウン効果のメカニズムをSalome-Mecaにて検証することが出来ました。

【参考文献】

(1)大倉・澤田・塚本,ベルト片寄りの基本メカニズムに関する一考察, MPTシンポジウム(伝動装置)講演論文集 巻:2004、 ページ:251-254 (2004)

(2)徳田・大倉,心線撚りによるベルト片寄り現象の解析, 日本機械学会機素潤滑設計部門講演会講演論文集 巻:6th、 ページ:31-34 (2006)

2013年8月31日 (土)

ゴムバンドの摩擦力による動力伝達(3次元)

■概要

これまでずっと取り組んできた「ゴムバンドの摩擦力による動力伝達」ですが、Salome-Meca2013.2を用いて2次元モデルでの問題が解決したので、次は3次元モデルでの試行に取り組みました。

ただし、従来の2次元モデルを厚み方向に伸ばして3次元モデルにするだけでは、3次元化の意味があまり無くて面白くありませんので、ローラーの外径を回転軸方向に変化させ、「クラウン効果」を検証してみることにしました。

「クラウン効果」は、平ベルトやバンドを太さが不均一なプーリーに巻きかけて回転させると、徐々に「太い方」へ寄っていくという現象です。この効果を利用して平ベルトの走行位置を中央に安定させるため、中央部が太くなったプーリー(クラウンプーリ)を用いる事があります。

今回の利用ソフトウェアは、Salome-Meca2013.2(Salome6.6.0 + Code_Aster11.4)、OSはDebian7.1です。インテルPC上にVirtualBoxで仮想PCを構築し、その上にDebian7.1を入れています。

■解析モデル(ジオメトリ)

下図の通りです。基本的には従来の2次元モデルを(サイズは変更していますが)厚み方向に伸ばした形状にしています。ただし、駆動ローラー、従動ローラー共に外径を軸方向で変化させています。
ローラーの幅は100mmで、奥側の外径が102mm、手前側が100mmと、奥側が太くなっています。

Rband0801


Rband0802


Rband0803


Rband0804

図では裏になって見えませんが、「dr_z1」と「dn_z1」は両ローラーの奥側の端面です。これはz=0の境界条件を与える面です

ローラーの内径は88mmです。
ゴムバンドの内径は220mm、外径は230mm(つまりゴムバンドの厚みは5mm)です。

従動ローラーの回転軸付近に小さなブロックを追加していますが、これについては後述します。

■メッシュ

Rband0805


ゴムバンドは六面体2次要素、ローラーは四面体1次要素としています。
ローラーについては応力や変形を評価しないため、精度が必要ないので1次要素としました。
ゴムバンドは2次要素ですがかなり粗いメッシュになっています。本来はもっと細かくするべきですが、今回は試解析ということで、計算時間の短縮を優先させています。

■解析の設定

従来の2次元の解析にほぼ準じているため、今回はcommファイルの詳細説明は省略します。2次元の解析説明の記事をご参照ください。

(2012/12/3)ゴムバンドの摩擦力による動力伝達(前編)
http://salome-meca.cocolog-nifty.com/blog/2012/12/post-e977.html

(2012/12/9)ゴムバンドの摩擦力による動力伝達(後編)
http://salome-meca.cocolog-nifty.com/blog/2012/12/post-02d3.html

(2013/1/3)ゴムバンドの摩擦力による動力伝達(問題点)
http://salome-meca.cocolog-nifty.com/blog/2013/01/post-5b53.html

(2013/7/24)ゴムバンドの摩擦力による動力伝達(問題点が解決)
http://salome-meca.cocolog-nifty.com/blog/2013/07/post-ae04.html

■解析初期における、変位制御と荷重制御の切り替え

今回、解析初期段階における従動側ローラーの制御に、小さなブロックを用いました。この部分が従来の2次元モデルと異なります。

下図の(A): 解析の初期段階においては、従動側ローラーを右に動かして、ゴムバンドを両ローラー間に張っていきます。しかしこの段階ではゴムバンドの剛性がほとんど無いため、すこしの力の増加で大きな変位が生じてしまい、(弱いバネを入れていても)荷重制御では収束しにくいです。
この段階では変位制御が望ましいです。

下図の(B): 一方、ゴムバンドを張り終えるときには、ゴムの伸びによって大きな力がローラーにかかっており、小さな変位で反力が大きく変化するため、変位制御では収束しにくいです。また、所定の軸荷重でゴムバンドを張るためにも、荷重制御が望ましいです。

Rband0806

この変位制御->荷重制御の切り替えは、STAT_NON_LINEを2段階に分けることで行えるのかもしれません。しかし、切り替え時に、変位制御で既に生じている反力を後続の荷重制御にどうやって移行したらよいか分かりません。(逆はOKです。荷重制御の解析の後、後続のSTAT_NON_LINEの変位設定で「TYPE_CHARGE='DIDI'」を用います)。

そこで今回、小さなブロックを従動ローラーのアーム端(ローラーの回転中心)に接して配置し、アーム端と接触を設定しました。そして解析の初期ではブロックに変位条件を与え、それで従動ローラーのアームを押して、変位制御としてゴムバンドを張っていきました。

そして、ある程度反力が生じたところでブロックの移動を停止し、次に従動ローラーのアームに荷重を与えて、所定の荷重値なるまで張りました(荷重制御)。
今回の解析では、ブロックの移動停止と荷重制御への切り替えを、time=0.5の時点で行っています。

この方法により、ゴムバンドを張る段階での収束性が向上しました。

この様子の動画です。

静止画です。

Rband0807

■条件等

今回、Time=2~5の間に、駆動ローラーを3回転させています。つまり、time=1あたり駆動ローラー1回転です。

軸荷重は30[N]、トルクは0.315[N・m]で、これはほとんど無負荷に近い弱いトルクです。
従来の2次元モデルと同様、Time=2~2.25の間にトルクを従動荷重でdn_trq面に与えています。

■結果の動画です。

静止画です。

Rband0808

クラウン効果によるゴムバンドの移動を分かりやすくするため、time=2~5のいて、駆動ローラー1回転ごとの静止画をアニメーションGIFにしました。

ゴムバンドが手前から奥へ移動しているのがわかります。
今回は駆動ローラーを3回転しかしていないので移動量はわずかですが、数十回転もすればローラーの奥の端へ到達してしまうと思います。

また、従動ローラーがゴムバンドの滑りによって駆動ローラーより若干遅れつつあるのが分かります。

クラウン効果については、後日もう少し詳しく評価したいと思います。

今回のデータは下記にzipでまとめています。

rband08.mmed: メッシュ
rband08a.comm: 前半commファイル(ゴムバンドを張るまで)
rband08b.comm: 後半commファイル(ローラー回転)

rband08.zip

2013年7月24日 (水)

ゴムバンドの摩擦力による動力伝達(問題点が解決)

The probrem of "Structural analysis of rubber products with contact, friction, and large rotation" was solved by Salome-Meca 2013.2(STA11.4).

In fact, the problem had been solved in STA11.3.22. However, I did not have the skill to update STA11.3 to STA11.3.22. Therefore, I was waiting for the release of salome-meca2013.2. It is easy to install.

半年以上前に実施して、問題点がなかなか解決できなかった事例、

「ゴムバンドの摩擦力による動力伝達」
http://salome-meca.cocolog-nifty.com/blog/2013/01/post-5b53.html

ですが、今回リリースされたSalome-Meca2013.2を用いたところ、問題が解決しました。

まず、【問題点1】(応力分布にムラがあり、ゴムバンドの動きに合わせて動いてしまう。
)ですが、これについてはCode_Asterの問題のようで、実はSTA11.3.22の段階で修正されていたようです。解説は下記URLのテキスト中にある、「Bad computation of the modulus of compressibility in ELAS_HYPER」です。

http://www.code-aster.org/documents/histor/en/11.3.22.txt

 しかし既存のSalome-Meca2013.1のCode_AsterをUpdateする事が私には難しくて出来なかったため、インストールが簡単なSalome-Mecaの最新バージョン公開を待っていました。

 そして先日公開されたSalome-Meca2013.2(Code_Aster11.4)で試したところ、問題が修正されている事が確認できました。

動画です。

負荷トルクが一定になるTime=2.25以降、ゴムバンドの応力や変形状態が定常状態になって回転している事がわかります。以前のような応力ムラや変形のガタツキはありません。

静止画です。
Rband0401

定常状態で回転し、安定しています。

次に【問題点2】(従動ローラーの回転が、駆動ローラから遅れている。)
ですが、これは商用ソフト(Marc)で試したところ、同じような結果になりました。これは問題点ではなく、正しい結果のようです。

以上のように、問題となっていた応力ムラ問題が解決したので、これからゴムと金属の接触、回転物の問題などを継続して進めて行きたいと思います。

2013年7月13日 (土)

金属シャフトとリング間へのゴム圧入

今回は、下図のように金属(鉄)シャフトとリングの間に、ゴムのリングを圧入する(押し込む)解析を検討しました。

Ring0101


ここで、上の3つのパーツの寸法は下記の通りです。

金属シャフト: 内径 - 、外径 20mm、高さ 50mm
ゴムリング : 内径 18mm、外径 62mm、高さ 10mm
金属リング : 内径 60mm、外径 70mm、高さ 20mm

金属シャフトとゴム間で直径2mm分、ゴムと金属リング間も直径2mm分オーバーラップしている事になります。このゴムを2つの金属パーツ間へ押し込みます。

このような問題のFEM解析を行う手法は複数あります。商用ソフトでは、要素を空間内で重ねて作成しておけば、この重なり(オーバーラップ分)を解析の初回インクリメントで自動的にゼロになるように調整してくれる機能や、金属パーツを剛体として扱い、その剛体の大きさを変化させる機能など、様々な工夫があります。

Salome-Meca(Code_Aster)では剛体面が使えないので、全てのパーツを有限要素で作成する必要が有ります。また、要素同士のオーバーラップを調整して圧入応力を生じさせる方法については、有るのかも知れませんが、私はまだ良く分かっていません。

そこで今回は、「熱膨張」を利用して金属シャフトと金属リングのサイズを変更する方法を用いる事にしました。これは、材料に線膨張係数を設定し、温度境界条件を与える事で、パーツの大きさを変える方法です。

ゴムリングの大きさはそのままとして、金属シャフトと金属リングがゴムリングに干渉しない大きさになるよう、初期形状を定めました。作業としては、SalomeのGeometryで本来の寸法で作成した後、Scale Transformで3方向均等に拡大(または縮小)しました。

金属シャフト: 倍率:0.869565(15%膨張させると本来のサイズになる)。
金属リング : 倍率:1.052632( 5%収縮させると本来のサイズになる)。

ここで、本来であれば高さ方向のサイズ変更は必要無い(むしろ避けるべき)だと思うのですが、異方性の線膨張係数を設定する方法が分からなかったので、今回は3方向とも均等にサイズ変更しました。

使用ソフトウェアは、debian6 + Salome-Meca2013.1(SALOME6.6 + Code_Aster STA11.3)の利用です。

また、今回の熱と構造の連成解析方法については、下記サイトの情報を参考にさせて頂きました。

Code_Aster & Salome-Meca 日本語解説

熱構造解析の部分

解析条件ですが、圧入だけでは面白くないので、圧入後、リングは固定したまま、シャフトの上部をX方向へ7mm移動させ、シャフトを傾けました。これにより、ゴムの変形が軸対称ではなくなります。

■結果

まずは結果からご覧下さい。

Time=0~1で圧入
Time=1~2でシャフト傾き

となっています。

静止画です。

初期形状(Time=0)
Ring0102

圧入(Time=1)
Ring0103

シャフト傾き(Time=2)
Ring0104

ゴムと金属間は接着ではないため、シャフトを傾けた結果、シャフトとゴムの間に隙間が生じています。

今回、圧入した時点でゴムと金属を接着する方法も検討したのですが、結局出来ませんでした。シャフトの傾きにより大回転が生じるので、「LIAISON」系は使えません(平行移動だけで大回転さえ生じなければ、LIAISON系の結合が使えるかと思います)。ContactでGLISSIERE='OUI'を指定することにより、接触面が離れないようには出来たのですが、面内でスライドしてしまいます。そこで大きな摩擦係数を設定しましたが、計算が収束しませんでした。商用ソフトのように、大回転に対応した「ボンド結合(貼り合わせ)」が有ると良いと思うのですが、Code_Asterでは良く分かりませんでした。

さて、ここから後は解析設定の詳細説明です。

この解析は「熱解析」と「構造解析」の2段階で行いますが、1ステップ(1インクリメント)毎に熱と構造を交互に行いつつ進めるのではありません。
はじめに熱解析をTime=0から最終時間まで行って、そのデータを一旦baseフォルダに保存します。次に構造解析でそれを読み込んで、温度データを参照しながら再びTime=0から最終時間まで計算します。

ジオメトリとメッシュは熱と構造で共通になります。

■ジオメトリ

○金属シャフト
下図の様になります。サイズは本来の0.869565倍となります。

Ring0105

グループ設定は下記の通りです。
高さ方向がz、リング・シャフトの半径方向がx,yとなります。

sft : シャフト本体
sft_c : 接触する面
sft_lx : x方向への移動を拘束する線
sft_ly : y方向への移動を拘束する線
sft_pt : シャフトの傾きを制御する点(x方向変位)
sft_pb: z方向への移動を拘束する点

○金属リング
下図の様になります。サイズは本来の1.052632倍となります。

Ring0106


ring : 金属リング本体
ring_c : 接触する面
ring_fz : z方向への移動を拘束する面
ring_lx : x方向への移動を拘束する線
ring_ly : y方向への移動を拘束する線

○ゴムリング
下図の様になります。サイズは本来の値そのままです。

Ring0107


rubr : ゴムリング本体
rubr_c1 : 接触する面(内側)
rubr_c1 : 接触する面(外側)
rubr_fz : z方向への移動を拘束する面
rubr_lx : x方向への移動を拘束する線
rubr_ly : y方向への移動を拘束する線

■メッシュ

全てのパーツを四面体2次要素にてメッシュ分割しました。

Ring0108


各パーツそれぞれメッシュ作成した後、「Mesh」-「Build Compound」で一つにまとめて、mmedファイルにエクスポートしました。

Ring0109


これは上から見たところ(z軸に平行な視点)ですが、3つのパーツが初期状態で離れている事が分かります。

■commファイル内容(熱解析)

commファイルの内容の説明です。まずは、第一段階の熱解析です。
ファイル名は「ring01a.comm」です。

Ring0110

AFFE_MODELEにて、要素の種類を3次元ソリッド要素「3D」、解析の種類を熱解析「THERMIQUE」にしました。

物性値の設定(DEFI_MATERIAU)は熱伝導率(LAMBDA)と比熱(RHO_CP)になるのですが、今回の解析は温度を強制的に指定して変化させるだけですので、これらの値はダミーで構わない事になります。これらの値は「DEFI_FONCTION」で温度依存値に指定できますが、今回は「1」一定にしました。

AFFE_MATERIAUにて、この物性値を全ての要素に適用しました。

###

Ring0111

温度境界条件です。この後の「THER_NON_LINE」にて、時間(INST)の関数「fonc」と、各条件の値を掛けたものが温度になります。

「ther0」はゴムの温度ですが、初期値0℃で変化させません。
そもそも、ゴム材料(ELAS_HYPER)には線膨張係数を設定できないので、寸法の変化を計算できません。

「ther1」は金属シャフトの温度で、Time=0のときの初期値を0℃として、Time=1で15℃、Time=2でも同じく15℃にしています。これは第2段階の構造解析で、線膨張係数を0.01(1℃上昇すると1%膨張する)に設定する予定であり、シャフトは本来の形状になるのに15%膨張する必要があるためです。

「ther2」は金属リングの温度で、こちらは収縮させる必要が有るので、Time=0で0℃から、Time=1で-5℃、Time=2でも-5℃としています。

###

Ring0112

pas1aは計算の時間刻み幅です。
Time=0~2の2秒間を、均等に40分割しています。
つまり、金属の寸法変化による圧入過程が20ステップ、シャフトの傾き過程が20ステップとなります。
温度を強制的に変化させるだけなので、収束の失敗は無いと思い、自動時間分割は使っていません。
出力の刻み幅(output)も同じにしています。

###

Ring0113

非線形の熱伝導解析の設定です。
温度境界条件では、ther0が固定なので、FONC_MULT無し、他の二つは時間によって温度変化するのでFONC_MULT=foncにしています。

■astkによるJOB実行

Ring0114

astkによるJOB実行時に、解析結果を保存するbaseフォルダを指定しました。
(「ring01.base」という名前の空フォルダを、現在のbath pathの下に作成しておきます。)

■熱解析結果

結果の動画です。

Time=0から1の間で温度変化をして、金属シャフトは15℃に、金属リングは-5℃になっています。ゴムは0℃のままです。
その後、Time=2までは変化ありません。

静止画です(Time=2の最終状態のみ)。

Ring0115

###

■構造解析

続いて、第2段階である構造解析のcommファイル説明です。
ファイル名は「ring01b.comm」です。

Ring0116

継続実行を示す「POURSUITE」を記述しています。

「MODI_MAILLAGE」は接触面の向きを揃えるものであり、接触を扱う場合には基本的に入れておきます。

「AFFE_MODELE」で、要素を3次元ソリッド要素(3D)、解析の種類を構造解析(MECANIQUE)にしました。

###

Ring0117

材料物性は、金属(鉄)は線形弾性「ELAS」として、線膨張係数(ALPHA)に「0.01」を指定しました。これは、温度1℃の変化で、長さが1%変化します。
ゴムはムーニーリブリン材料です。

これを「AFFE_MATERIAU」の「AFFE」で各パーツに設定しました。
また、「AFFE_VARC」にて、熱解析の結果「RESU」を参照しています。この値と、 VALE_REFの値(今回は0)の差分に線膨張係数を掛けた値が熱による膨張率となります。
(この部分は、こちらを参考にしています)。

###

Ring0118

接触の設定です。
今回は簡単に済ませるため、摩擦は無しにしました。
また、ALGO_REZO_GEOMはPOINT_FIXEにしました。Code_Aster STA11.3のデフォルトはNEWTONだったと思いますが、私の環境ではPOINT_FIXEのほうが高速です。

また、材料の硬い金属の方をマスター面(MAIT)にしました。

###

Ring0119

計算時間幅(時間刻み)の設定です。
接触を使用するので、自動時間増分機能を使っています。

pas2a:Time=0~1の時間刻み(均等20分割)
pas2b:pas2aに自動増分を適用して再定義したもの
pas3a:Time=1~2の時間刻み(均等20分割)
pas3b:pas3aに自動増分を適用して再定義したもの

###

Ring0120

境界条件の設定です。
ジオメトリでグループ設定した面、辺、点にそれぞれ拘束を設定しました。

時間(INST)の関数として変化するのは「char1」のsft_ptのみです、この点はTime=1~2の間にx方向に7mm移動し、シャフトが傾きます。

###

Ring0121

非線形構造解析の設定ですが、今回はTime=1の時点で二つに分割しています(STAT_NON_LINEを2つ重ねている)。
これには下記の事情があったためです。

始めはSTAT_NON_LINEをTime=0~2の間で一つにしていました。そのとき、金属パーツ(RELATION=ELAS)のDEFORMATIONは、シャフトの傾きによる大回転を考慮し「GROT_GDEP」にしていました。ところが、この設定ではTime=1時点の熱変形後の金属パーツのサイズが、本来の寸法と異なった値になっていました。シャフトは正確に15%膨張、リングも正確に5%縮小してほしかったのですが、それぞれ若干異なった値になってしまいました。

この原因ですが、DEFORMATION=PETITにしなければ、予定の寸法にならないようです。大回転への対応は譲れないところなので、大回転可能なGDEF_LOGやSIMO_MIEHEも試しましたが、どれも寸法が異なってしまいました。PETITしかないようです。

そのため、回転の生じない熱変形過程のTime=0~1においてはPETIT、シャフトが大回転するTime=1~2においてはGROT_GDEPを用いる事にしました。
このような理由で、STAT_NON_LINEを二つに分けることになりました。
(ゴムは両方ともGROT_GDEPです)。

なお、STAT_NON_LINEの結果は、「RESU2」として、熱解析の結果「RESU」と分けています。

###

Ring0122


結果出力の部分は、特に変わったことはしていません。

■astkによるJOB実行

Ring0123

astkによるJOB実行時に、熱解析の結果データである「ring01a.base」フォルダを読み込み専用(Dのみチェック)で指定しました。Rにチェックすると、JOB完了時に上書きされるので注意が必要です。

■構造解析の結果

 結果は一番初めに示したとおりです。

■今回のデータ

下記にzipでひとまとめにしています。

ring01.zip

ring01a.comm(第一段階・熱)
ring01b.comm(第二段階・構造)
ring01.mmed (メッシュ)

今回は以上です。

2013年6月22日 (土)

分離せずにスライドする面同士の接触(GLISSIERE)

Code_asterの接触機能の一つに、「分離せずにスライドする」というオプションがあります。これは「DEFI_CONTACT」の「ZONE」指定の部分にて、「GLISSIERE='OUI'」を指定すると実現できます。
ちなみに「GLISSIERE」は、英語で言うと「Slide」の意味だそうです。

私は普段Marcを使っているのですが、Marcにも同様の機能があり、良く使っています。
例えば、長い物の一部分を切り出して解析する際、対称面が屈曲などで移動するような条件を与えなければならないケースが出てきます。そのような場合、移動する対称面に剛体壁を接触させ、分離せずにスライドする(接触面内での変形も可能)と設定します。Marcの場合は剛体面が使えますが、Code_asterは有限要素を用いなければなりません。また、剛体結合(LIAISON_SOLIDE)は大回転できないため、回転させる場合は変形体として扱う必要が有ります。

今回、ゴムチューブの一部を切り出して、屈曲させる解析の例で試してみました。

利用ソフトウェアは、DEXCS2012-Salome-D1-B1-64日本語版(SALOME 6.6.0 + Code_Astr 11.3)です。Salome-Meca 2013.1に相当します。

まずは結果の動画からご覧下さい。

ゴムチューブと奥の壁は分離しませんが、固着もしておらず、面内でスライドおよび変形をしている事が分かります。

変形前の静止画です。メッシュは2次要素を用いています。

Slide301

変形後の静止画です。

Slide302

ゴムチューブの手前側の切断面は、単なる対称面として面内(X方向)に拘束しています。
奥側の切断面には、鉄の壁を当てて配置しています。この鉄の壁とゴムチューブの奥側の切断面に「GLISSIERE='OUI'」での接触を定義し、鉄の壁を傾けて、ゴムチューブを屈曲させていきます。

この方法により、両方の切断面でゴムチューブは面内変形する事ができ、切断面の拘束による変な応力も発生しません。

下図は、鉄の壁に接触しているゴムチューブの切断面の変形状態です。手前側の切断面同様に変形している事が分かります。

Slide303

ここからは解析内容の詳細をご説明します。

まずSalomeによる全体のモデルです。

Slide304

ゴムチューブの一部分を、鉄の壁に沿わせて配置します。

ゴムチューブのグループ設定です。

Slide305


手前の切断面(tube_fx)は対称面とします。奥の切断面(tube_c)は、鉄の壁に接触する面とします。
tube_p1は、初期の剛体移動を防ぐための弱いバネを設定するポイントです。
tube_p2は、チューブの回転を抑えるため、Z方向の変位を同期させる点のグループです。

鉄の壁のグループ設定です。

Slide306_2



手前の面(wall_c)は、ゴムチューブと接触する面です。
下側のライン(wall_l1)と上側のライン(wall_l2)は、この壁を傾けていくための変位境界条件を与えるラインです。
また、wall_p1の点二つは、壁がY方向へ移動したり、回転するのを防ぐための拘束を与えます。

ここから、commファイルの内容をご説明します。

Slide307


まず、複数の面で構成された接触面の方向を整えるための「ORIE_PEAU_3D」です。これは接触を扱うときには基本的に入れておきます。

Slide308


「CREA_MAILLAGE」にて、剛体移動を防ぐための弱いバネを、節点グループ「tube_p1」に付けます。そのバネのグループ名を「tube_sp」という名前にしています。

「AFFE_MODELE」で要素のタイプを指定しますが、ゴムチューブ(tube)と鉄の壁(steel)は3次元ソリッド要素(3D)です。バネは「DIS_T」としており、「AFFE_CARA_ELEM」でタイプを「K_T_D_N」として、バネ定数を(0,0,1e-4,)に設定しました。つまり、x,y方向はバネ定数ゼロで、z方向に1e-4 N/mmとなります。

Slide309


材料物性ですが、鉄の壁は「steel」として鉄の物性値、ゴムチューブは「rubber」としてムーニーリブリン材料として設定しました。

Slide310


接触の設定は上記の通りです。
摩擦はありません(対称面としてスライドさせるため)。
今回使用したCode_asterのバージョンは11.3なので、接触については「Generalized Newton method」を用いる事が出来ます。しかし、「ALGO_RESO_GEOM」に関してはあえて「POINT_FIXE」を用いています。ここを「NEWTON」にすると、計算時間がかなり長くなってしまうためです。
理由は分かりませんが、GEOMに関しては「POINT_FIXE」の方が速く計算できるという結果になっています。

そして今回の主題ですが、「ZONE」の設定において、「GLISSIERE='OUI'」を設定しています。これによって「分離せずにスライドする」を実現しています。

また、両面は初期状態で接触しているので、「CONTACT_INIT='OUI'」を入れています。

Slide311


pas1aは計算の時間ステップ定義ですが、time=0~1の間を、50分割しています。
outputは出力の時間ステップですが、これもtime=0~1の間を、50分割しています。

Slide312


pas1bはpas1aに自動時間ステップ調整を適用したものです。

Slide313


「AFFE_CHAR_MECA」は境界条件の設定です。

char0は、時間によって変化しない条件です。
tube_p1にはz方向に弱いバネを設定していますが、y方向には動かないようにDY=0を設定しています。
同じく鉄の壁のwall_p1も、y方向に動かないようにDY=0としています。
鉄の壁の下の辺であるwall_l1は、上下に動かないよう、DZ=0としています。

tube_fxの面は対称面なので、「FACE_IMPO」で「DX=0」として面内に拘束しています。

tube_p2の2点について、「LIAISON_UNIF」にてZ方向の変位を連結することで、ゴムチューブが軸中心に回転することを防いでいます。

char1は時間によって変化する条件です。
鉄の壁の上下の辺(wall_l1とwall_l2)をX方向に互いに逆に動かし、壁を傾けています。
壁が傾く事で、ゴムチューブの断面も傾き、屈曲します。
ゴムチューブの切断面は両方とも面内を自由にスライドできますので、X平面と鉄の壁の間で、安定する位置でゴムチューブが静止します。

Slide314


「STAT_NON_LINE」は非線形解析の設定です。
特に変わったことはしていませんが、収束判定(CONVERGENCE)については、絶対値による判定「RESI_GLOB_MAXI」にしています。相対値による判定では、解析の初期段階において反力が小さいことで、収束しにくくなるためです。

また、収束しやすくするために、ニュートン法のイテレーション毎に接線剛性を計算する設定「REAC_ITER=1」を入れています。

Slide315


出力に関しても、特別な設定はありません。
Code_asterのバージョン10で使っていた「CALC_ELEM」は使えないので、「CALC_CHAMP」を使っています。

説明は以上となります。

今回のメッシュファイル(.mmed)とcommファイルは、下記のリンクです(zipでまとめています)。

slide3.zip

2013年1月 3日 (木)

ゴムバンドの摩擦力による動力伝達(問題点)

「ゴムバンドの摩擦力による動力伝達」について、前々回(前編)前回(後編)の2回に分けてご紹介しましたが、また問題点が残っていました。

これについて改良を試みましたが、まだ解決出来ていません。現在の状況をご報告します。

利用ソフトウェアは、Salome-Meca2012.2(SALOME 6.5 + Code_Astr 11.2.0)で、OSはdebian6です。

【問題点1】

ゴムバンドの応力分布にムラがあり、ゴムバンドの動きに合わせて動いてしまう。

ゴムバンドが二つのローラーに掛けられて定常状態で回転しているとき、応力や歪はメッシュ(つまりゴムバンド)について動いていくのではなく、空間上の「位置」によって定まるようになるはずです。これは、パイプを流れる液体の流体解析において、液体の状態がその「位置」によって決まるのと同じです。

ただし、流体解析は実際に空間に対してメッシュを作成するのですが、今回は構造解析ですので、ゴムバンドについてメッシュを作成しています。それでも、定常状態における応力は位置によって決まるはずなのです。つまり、ゴムバンドを基準として見た場合、ある決まった応力分布が、ゴムバンド上を一定速度で移動していくことになります。

しかし、前回の解析では、応力の大きい部分がゴムバンドの移動につれて動いています。

Rband0301


【問題点2】

従動ローラーの回転が、駆動ローラから遅れている。

駆動ローラーはtime=2~4の間に2回転するので、従動ローラーも2回転するはずですが、下図の通り、やや遅れています。すべりが大きいのかもしれません。

Rband0302


------------------------------------------------

これらの問題を解決するため、色々試行錯誤してみました。

まず【問題点1】ですが、一番疑わしいのは、メッシュが不均質なために応力にムラが生じているのではないかという事です。そこで、ゴムバンドは出来るだけ形の揃った四辺形(2次要素)でメッシュを作成しなおしました。
また、ローラーについても、表面の要素辺サイズが揃うようにしました。

Rband0303

それから怪しいのが、剛体移動防止のために設定している「弱いバネ」です。これにゴムバンドが引っ張られているのかもしれません。
前回のモデルでは、バネ定数は全て10のマイナス3乗でしたが、これを10のマイナス6乗まで小さくしました(それ以上小さくすると、不安定になり収束しませんでした)。

Rband0304


次に【問題点2】です。
前回は、摩擦力の計算アルゴリズム(ALGO_CONT)を「STANDARD」にしていました。これでは、ローラーを回転させる段階において、「COEF_FROT」(デフォルト値:100)を、「10」まで下げないと計算が収束しませんでしたので、その値にしています。これが小さすぎて、摩擦力が小さめになっている可能性があります。

そこで今回は、ALGO_FROTを「PENALISATION」に変えました。これまでの経験から、すべりの大きな摩擦計算では、「STANDARD」よりも「PENALISATION」のほうが収束しやすかったためです。その上で、「COEF_PENA_FROT」を「500」にしました。ゴムの物性値がヤング率換算で5[MPa]程度ですので、二桁大きな値になります。

Rband0305

また、初回インクリメントでなかなか収束しにくい事が多かったのですが、それに関しては「SEUIL_INIT」に「0.1」を設定する事でかなり改善しました。これは初期の接触面圧を与えるもので、剛体移動の抑制に効くようです。この値は大きいほど効果があるようですが、反面どのような悪影響(副作用)があるか分かりませんので、出来るだけ小さな値で済むようにしました。今回は0.1程度より大きくすると収束性改善効果が感じられましたので、0.1にしました。

これについては、公式ドキュメント「u2.04.04」の4.2節に解説があります。
http://www.code-aster.org/V2/doc/default/en/man_u/u2/u2.04.04.pdf

これによると、継続して実行するSTAT_NON_LINEにおいては「CONTACT_INIT='INTERPENETRE'」も指定したほうが良さそうなので、これも入れていますが、これの効果はいまひとつ分かりません。

なお、計算のステップは下記の通りで、

(1)time 0 -> 1
 駆動ローラーを固定(drp1, drp2を固定)し、従動ローラー中心(center)に右向きの荷重を与え、ゴムバンドを所定の軸荷重で張っていく。

(2)time 1 -> 2
 ゴムバンドを張った状態で、従動ローラーをcenter位置でピン固定し、軸荷重は除去。

(3)time 2 -> 4
 駆動ローラーを反時計回りに2回転させる(time=1あたり1回転)。
 このとき、回転始め~1/4回転迄の間で、従動ローラーのトルクを増加させて所定の値にする。

前回は(2)と(3)の間でJOBを前編と後編に区切っていたのですが、(2)の最初で収束しにくい状況だったため、試行の効率を考え、今回は(1)と(2)の間で区切りました。

---------------------------------------------------

以上のような変更を行って計算したのですが、最初に書いた通り、残念ながら解決には至りませんでした。

まず【問題点1】ですが、やはりゴムバンドの応力にはムラが生じており、ゴムバンドと一緒に移動してしまいます。

Rband0306


このゴムバンドの「応力ムラ」についてよく見ると、どうやら初期状態からバンドを引っ張った状態(Time=1~2)のときの応力状態を「引きずって」いるようです。

Rband0307


このときに応力の大きい部分が、回転しても応力が大きいままです。
しかしこれはおかしい事です。
塑性変形しているなら分かりますが、これは超弾性材料なので、変形途中の応力状態に影響されるはずがありません。あくまでも初期状態(真円)のときが応力ゼロとして計算されるべきです。

なにか、塑性変形的な設定をしてしまっているのかと思い、STAT_NON_LINEの中の「COMP_INCR」を「COMP_ELAS」に変えても結果は同じでした(塑性解析では違いが出ます)。

また、摩擦の影響かと思って、Time=0~2までを「摩擦なし」で計算させてみましたが、それでも同じでした。
低減積分要素を用いても変わらず。

ということで行き詰りました。

どうも、STAT_NON_LINEを連続して実行するときの、前段階からの状態の引継ぎ設定【ETAT_INIT=_F(EVOL_NOLI=RESU,),】あたりが怪しそうなのですが・・・。

それから【問題点2】についても、前回の結果と同じで、従動ローラーの回転が遅れています。

Rband0308

しかしこれは、もしかしたら正しいのかもしれません。せん断変形などによって、実際にずれていくということも考えられます。もう少し薄いバンドで、屈曲変形の影響を排除して確かめる必要が有りそうです。

以上のように、まだまだ検討の余地が有ります。

今回のファイルは下記です。

rband03.zip

(ZIPファイルの内容は下記の通りです)

メッシュファイル: rband03.mmed
commファイル(前編): rband03a.comm
commファイル(後編): rband03b.comm

2012年12月 9日 (日)

ゴムバンドの摩擦力による動力伝達(後編)

前回(前編)の続きで、2本のローラー間にゴムのバンドを掛けて、その摩擦力でローラーの回転と動力を伝達させるという事例です。

(前編の記事掲載以後、様々な改良を行っていますが、今回はとりあえず改良前のモデルの後編の説明です)

利用ソフトウェアは、Salome-Meca2012.2(SALOME 6.5 + Code_Astr 11.2.0)で、OSはdebian6です。

前編では、ローラーに所定の軸荷重を与えてゴムバンドを張り(time=0~1)、その位置で従動ローラの中心を固定し、軸荷重および安定化の為の拘束を除去したところまで(time=1~2)でした。

後編では、time=2~4で、駆動ローラーを反時計回りに2回転させます。つまり、time=1あたり1回転になります。また、回転はじめから1/4回転までの間で、従動ローラーの内側突起部の側面(trq)に負荷トルクとしての荷重を与えていき、1/4回転からは一定のトルク値を保持します。

以下、後編のcommファイル(rband01b.comm)の説明です。

■計算の継続

Rband0117


前回のデータを読み込み、継続して計算をする場合には「POURSUITE」を用います。
Eficasで新規にcommファイルを作成する際に「POURSUITE」を追加すると、前回のcommファイルを指定するようにメッセージが表示され、ファイルを選択するBOXが開きますので、そこで前回の「rband01a.comm」を選択します。
また、「POURSUITE」から始まるcommファイルを開いたときも同様で、前回のcommファイルを指定します。
この作業によって、前回のcommファイルで定義済みの変数や関数、パラメータなどがそのまま使えるようになります。

■駆動ローラーの回転変位設定(drrot)

Rband0118



time=2~4で、駆動ローラー上の2ヶ所の点(drp1, drp2)に変位を与え、駆動ローラーを反時計回りに回転させる設定です。time=1あたり1回転とします。
2ヶ所の点のx,yの変位を、INST(time)の関数として、それぞれ関数(FORMULE)として定義し、「AFFE_CHAR_MECA」で境界条件「drrot」として設定しました。

このように変形体を回転させる方法についての詳細は「回転体とゴム板の接触」をご参照ください(参照先の制御節点はX軸上に並んでいますが、今回はY軸上に並んでいますので、式の内容は異なります)。

2ヶ所の点は駆動プーリの内径ライン上に位置するので、振幅は内径の半径である「44mm」になります。

■負荷トルクの設定(trq)

Rband0119_2


従動ローラーの内径側に出ている突起の側面(trq面)に面圧荷重を与え、トルクを発生させます(突起の詳細など、モデルについては前編をご覧下さい)。

従動ローラー中心からtrq面中央までの距離が43[mm]、圧力値が1.16[N/mm^2]、trq面の面積が2[mm^2]なので、トルクは、
    43 x 1.16 x 2 = 100 [N・mm] となります。

厳密に言えば、trq面への荷重は半径方向に分布しているので、上記のような計算式ではなく積分計算をする必要がありますが、ここはtrq面が微小であるとして、単純な計算で済ませています。

■接触の設定

Rband0120



接触の設定「DEFI_CONTACT」は既に前編のcommファイル(rband01a.comm)で設定済みなので、それと同じ設定で良ければ、このcommファイルで再設定する必要はありません。しかし今回の試行では、前編の設定のままでは後編の回転部分においてどうしても計算が収束しませんでした。そのため、前編(cont)と異なる名前(cont2)にて再設定しました。

前編と違うところは、「COEF_FROT」の値を「10」にした所です。前編はデフォルト値の「100」でしたので、その1/10になっています。「COEF_FROT」は公式ドキュメント「U4.44.11」でもあまり詳しい記述が無いのですが、おそらく、接触する両面間に相対変位が生じたときに、摩擦力が変位と共に立ち上がるときの傾きのようなものではないでしょうか(これについては自信がありません)。もしそうだとすると、この値が小さくなると、相対変位に対する摩擦力の増加が緩やかになり、収束しやすくなる反面、現実と比べて十分な摩擦力が得られない可能性もあります。

■タイムステップの設定(time=2~4)

Rband0121



pas3aは、回転開始のtime=2から、1/4回転してトルクが所定値に達するtime=2.25までを180分割と細かくし、それからトルクが一定値を保持するtime=4までは315分割と粗くしています。

output3は結果の出力間隔ですが、これは2回転の間を均等に360分割しています。つまり、駆動ローラーが2°回転する毎に結果が出力されます。

pas3bは時間ステップ幅の自動調整です。これは前編とあまり変わりませんが、ステップ幅を増やす判断に用いるNewton法の収束回数基準値を、9回と多めにしました。

■非線形解析の設定(time=2~4)

Rband0122



境界条件は、今回定義した「drrot(駆動ローラー回転変位)」と「trq(従動ローラートルク荷重)」を追加しています。従動ローラーセンターのピン固定は継続です。

ここで、「trq」については追従荷重(TYPE_CHARGE='SUIV')を指定し、ローラーが回転しても、常にtrqは突起側面に垂直になるようにしています。この設定で、常に一定のトルクが作用するようにしています。下の動画をご参照ください。

動画を見られない方は、下の静止画をご覧ください。

Rband0124

また、追従荷重の詳細については、「追従荷重」をご参照ください。

それから、収束条件を残差の相対値「REGI_GLOB_RELA」で、値を「1e-4」にしました。
デフォルトは「1e-6」なのですが、それでは収束しなかったため、緩めの判定値にしました。

■結果の計算と出力

Rband0123



これは特に変わった事はしていません。

以上の設定にて、astkから計算を実行しました。

Rband0125



baseについてですが、前編でデータを「rband01a.base」フォルダに書き出していました。
したがって後編ではそこから読み出しても良いのですが、それでは後編の計算が終わったときに上書きされてしまいます。
後編の計算を何度も試行錯誤したかったので、「rband01b.base」フォルダを別に作り、試行のたびに「rband01a.base」の内容を上書きコピーしてから計算をスタートさせるという方法にしました。

以下は結果の図です(Mises応力)。ローラーの応力は不要なので非表示にしています。

●time=2

Rband0126


引っ張っただけなので、上下対称の変形・応力状態です。

●time=2.25

Rband0127


駆動ローラーが反時計回りに1/4回転した状態で、トルクが所定値に達しています。
両ローラー間でゴムバンドが張られている所(スパン)では、トルクによって、上側が強く張られていて応力も大きい事が分かります。

●time=3

Rband0128


駆動ローラーが反時計回りに1回転した状態です。
しかし従動ローラーは1回転に若干足りていない事が分かります。
これは、ローラーとゴムバンドの「すべり」が原因と思われます。

●time=4

Rband0129


駆動ローラーが反時計回りに1回転した状態です。
従動ローラーは更に遅れています。トルクが一定状態でも「すべり」があるようです。

それに、応力の分布状態に「ムラ」があります。これについては改善が必要です。

改めて結果の動画を載せます。

このように、一応この種の解析が可能である事が分かりました。
しかし、まだ下記のような問題が残っており、実用にはそれらの解決が必要です。

  • 応力や変形が若干いびつになっている。これはバネの影響が残っているか、メッシュが粗すぎるなどの原因が考えられる。
  • 計算が不安定で、タイプステップが制限値まで細かく分割されても収束せず、JOBが止まる事が多い。
  • 計算時間が長く掛かる。これは、ローラーを剛体面扱いできず、変形体として解いている事が大きく影響していると思われる。剛体サーフェスの導入か、せめて剛体結合時にも大回転が出来るようにしてほしい。
  • 接触・摩擦計算については更なる高速化と、安定性が必要。
  • これらの問題については、開発元の対応がない限りどうしようもない事もありますが、自分で改善できそうなところは、出来るだけ詰めていきたいと思います。

今回のデータは下記です(前編と同じです)。

rband01.zip

メッシュファイル: rband01.mmed
commファイル(前編): rband01a.comm
commファイル(後編): rband01b.comm

2012年12月 3日 (月)

ゴムバンドの摩擦力による動力伝達(前編)

今回は、2本のローラー間にゴムのバンドを掛けて、その摩擦力でローラーの回転と動力を伝達させるという事例です。設定内容が複雑なため、前編と後編の2回に分けることにします。

利用ソフトウェアは、Salome-Meca2012.2(SALOME 6.5 + Code_Astr 11.2.0)で、OSはdebian6です。接触・摩擦計算が高速化されたCode_aster11.2を使いたかったため、いつものDEXCS-Salomeではなく、この環境になりました。

また、モデル化は2次元(平面ひずみ)としました。はじめから3次元で行うと、試行錯誤に大変な時間がかかるためです。

まず、ローラーとゴムバンドのモデルです。

Rband0101


ゴムバンド(band)は真円の状態が初期形状です。その中に、中空のローラー2個を配置します。向かって左側が駆動ローラー(dr)、右側が従動ローラー(dn)です。drはdriving、dnはdrivenの略です。

寸法は、
・ゴムバンド(外径:260mm、内径:240mm)
・ローラー (外径:100mm、内径:88mm)

従動ローラーの拡大図です。

Rband0102


センター位置をピン固定するためのアームが1本出ているのと、トルクを与えるための突起(高さ2mm)が、内側へ出ています。

上図のように各グループを設定しましたが、それぞれの役割は下記の通りです。

【面グループ】
・band:ゴムバンドの領域設定
・dr, dn:各ローラーの領域設定

【線グループ】
・bandc:ゴムバンドの内面(接触面の定義)
・drc, dnc:各ローラーの外面(接触面の定義)
・trq:従動ローラーに負荷トルクを与えるための圧力荷重面

【点グループ】
・bandp1, bandp2, dnp1:不安定性対策の弱いバネ、および拘束の設定箇所
・drp1, drp2:駆動ローラーの回転変位制御点
・center:従動ローラーの軸荷重点、および回転中心

このモデルについて、下記の手順で計算します。

(1)time 0 -> 1
 駆動ローラーを固定(drp1, drp2を固定)し、従動ローラー中心(center)に右向きの荷重を与え、ゴムバンドを所定の軸荷重で張っていく。

(2)time 1 -> 2
 ゴムバンドを張った状態で、従動ローラーをcenter位置でピン固定し、軸荷重は除去。

(3)time 2 -> 4
 駆動ローラーを反時計回りに2回転させる(time=1あたり1回転)。
 このとき、回転始め~1/4回転迄の間で、従動ローラーのトルクを増加させて所定の値にする。

文章で説明しても分かりにくいので、まずは最終計算結果の動画をご覧下さい。

実は、このような解析が、私がSALOME-Meca導入当初からやりたかった事です(実際はもっと複雑ですが)。のんびりやっていたため、導入・試行を始めてから1年以上かかってしまいました。
今年の春頃から、SALOMEや Code_Asterの様々な機能を検証してきたのですが、その多くがこの解析のための前準備です。

次に、メッシュ分割です。

Rband01mesh_2

全て三角形2次要素で分割しました。
band, dr, dn はそれぞれ別にMesh作成を行い、Compoundでまとめています。
ただし、座標が重なる節点についてのマージは行っていません(回転によってズレていくため)。
もちろん、各グループ名をジオメトリから設定しています。

さて、ここからが今回(前編)の設定内容なのですが、これまでのようにEficasの画面スナップショットを用いると画像が膨大な量になってしまうため、commファイルの内容を用いる事にします(それでもかなり多いですが)。
前編のcommファイル名は「rband01a.comm」です。

前編では、上記の解析手順の(2)まで、つまり、バンドを張った状態で従動ローラーを固定し、それまで与えていた軸荷重を抜くところまでです。残りの回転部分については、後編にて継続計算(POURSUITE)で行います。

■メッシュの修正(エッジの方向)

Rband0103

接触や面圧荷重を受けるエッジ(線)の方向を整える設定です。

■弱いバネ要素の作成

Rband0104


初期の不安定な状態において計算の収束性を向上させるため、ゴムバンド2ヶ所、従動ローラー1ヶ所に弱いバネを作成します。駆動ローラー近くで変位が小さい所(bandp1)をspring1、従動ローラー近くで大きく変位する所(bandp1,dnp1)をspring2と二つに分け、バネ定数などの設定を分けられるようにしています。

■要素の種類の設定

Rband0105

ゴムバンド(band)、ローラー(dr, dn,trq)は平面ひずみ「D_PLAN」、バネ(spring1, spring2)は「2D_DIS_T」としました。ここで、負荷トルクのための圧力荷重面(trq)も「D_PLAN」に含める必要がありました(その理由は分かりません)。

■バネ定数の設定

Rband0106

弱いバネのバネ定数の設定です。種類は「K_T_D_N」です。
ここの値の設定は、何度も値を変えての試行錯誤が必要で、かなり苦労しました。
バネが弱すぎると初期の不安定性を抑えきれないし、強すぎるとゴムバンドが大きく伸びたときにバネの影響が出てしまいます。最終的に、spring1, spring2共に“1e-3”で落ち着きました。

私が使っている商用ソフトでは、バネ定数をTimeの関数にして、バネを徐々に消滅させていく事が出来るのですが、Code_Asterにはそのような機能が無いようです(“VALE_F”はあったのですが、Timeの関数には出来ませんでした)。この機能はぜひ入れて欲しいところです。

■材料物性の設定

Rband0107

ゴムバンドはMooney-Rivlin材料としました(「超弾性材料の検証」を参照)。

また、ローラーは剛体にしたかったのですが、「LIAISON_SOLIDEによる剛体は大回転出来ない」の結果から、やむなく変形体として、大きなヤング率(1000[GPa])を設定しました。

■接触の設定

Rband0108

Code_aster 11.0.28以降では、接触・摩擦解析が高速化」の結果から、Code_Aster11.2.0を用い、
 ・ALGO_RESO_FROT=’NEWTON’
 ・ALGO_RESO_CONT=’NEWTON’
を設定しました。

また、ゴムバンドを「マスター」、2つのローラーを「スレーブ」としています。
硬い材料であるローラーをマスターにしたかったのですが、複数のZONEに共通している面がある場合、それをマスターにしないといけないようです(今回の場合はbandです)。

摩擦係数は0.4で、他にはあまり調整したところはありません。

■タイムステップの設定(time=0~1)

Rband0109


pas1aは、従動ローラーを右(X方向)に引っ張り、ゴムバンドを張っていく段階(time=0~1)のタイムステップ定義です。解析初期(time=0.5まで)の不安定なところの刻みを多くしています。

output1は、出力の時間刻み設定です。
これに関しては「結果出力の間隔設定」をご参照ください。
今回は、time=0~1を均等に50分割しています。

pas1bは、pas1aを基本として、タイムステップの自動調整を追加したものです。これについては「時間ステップ幅の自動調整機能」をご参照ください。
今回のような解析は非常に不安定で収束性が悪く、この機能が無ければ、現実的な時間で計算する事が難しかったと思います。

最小時間幅(PAS_MINI)の値は小さめにしています。

また、時間幅を増やす場合の増加割合(PCENT_AUGM)のデフォルトは100(%)なのですが、時間幅が元の2倍になると再び収束しなくなって、もう一度分割される、という非効率な現象が良く起こっていたので、これを50(%)にしました。

■タイムステップの設定(time=1~2)

Rband0110

time=1~2では境界条件の切り替えを行うだけなので、時間刻み幅を「1」としてます。出力の刻み幅も「1」です。
一応、収束しなかったときのために自動調整機能も入れていますが、実際は動作せず、時間刻み幅は「1」のままでした。

■境界条件の設定

Rband0111


・fix_dr:駆動ローラーの変位固定
・fix_y_c:従動ローラーのセンター位置を、上下(y方向)に固定
・force_dn:従動ローラーを右に引っ張る軸荷重

       10[N]としています。
       平面ひずみ要素を用いているので、厚み(Z方向)が1mmあたりの値になります。

・fix_x_c:従動ローラーのセンター位置を、左右(x方向)に固定
・fix_y_b:弱いバネを設定した節点を、y方向に固定
      バネがあるので本来この条件は不要なはずですが、これを入れないと収束しませんでした。

・force:force_dnを増減させるための関数です。time=0~1では増加します。特に解析初期であるtime=0.5までの増加量を少なくしています。
time=1~2でゼロになります。

■非線形解析の設定(time=0~1)

Rband0112


time=0~1では、駆動ローラーを固定し、従動ローラーに軸荷重(force_dn)を与えてゴムバンドを引っ張っていきます。そのため、「fix_x_c」は使いません。

COMP_INCRでは、ローラー(ELAS)とゴムバンド(ELAS_HYPER)を異なるRELATIONに設定しました。また「異なる「DEFORMATION」設定の混在」も可能なのですが、今回は両方とも「GROT_GDEP」にしています。

NEWTON=_F(REAC_ITER=1,), は、ニュートン法の求解過程において、イテレーション毎に接線剛性を更新するものです。デフォルトではインクリメント(タイムステップ)毎になっているので、非線形性の強い計算では収束性が悪くなります。詳しくは「追従荷重(その2)」を参照ください。

■非線形解析の設定(time=1~2)

Rband0113

time=1~2では、ゴムバンドを引っ張った状態で従動ローラーのセンター位置を固定(fix_x_cの条件を追加)し、軸荷重(force_dn)を除去します。また、後編でのローラー回転に備え、安定化のためのy方向固定(fix_y_b)も除去します。

ここで、従動ローラーのセンター一固定(fix_x_cの条件)については、初期位置から変位した場所で固定するため、「TYPE_CHARGE='DIDI'」の設定が必要になります。これに関しては「荷重を与えて変形させた状態で固定する方法(境界条件の切り替え)」を参照ください。

また、time=1からの継続計算なので、「ETAT_INIT=_F(EVOL_NORI=RESU,)」を設定しています。

■結果の計算と出力設定です。

Rband0114


これは特に問題ないと思います。

--------------------------------------------------------

以上の設定で、ASTKを用いて計算を実行しました。

Rband0115


この後に後編で続きの計算(ローラー回転)を行うので、データを「rband01a.base」フォルダに保存する設定にしておきます。

time=2の時点での結果です(Mises応力)。

Rband0116


ゴムバンドが伸ばされて応力が生じている事が分かります。
現段階では単に引っ張っただけですので、変形や応力状態は上下対称になっています(メッシュによる誤差を除いて)。

なお、ローラーも変形体なので応力が出力されているのですが、評価対象ではないため非表示にしています。そのほうが、ゴムバンドの応力が見やすくなります。

前編はここまでです。

今回のデータは下記です
rband01.zip

メッシュファイル: rband01.mmed
commファイル(前編): rband01a.comm
commファイル(後編): rband01b.comm

2012年11月11日 (日)

Code_aster 11.0.28以降では、接触・摩擦解析が高速化

Code_asterの公式サイトの記事(下記URL)にて、version 11.0.28以降において接触・摩擦の計算が高速化できるアルゴリズムが導入されている事が書かれています。

http://www.code-aster.org/V2/spip.php?article612

version10までは、接触・摩擦を含む非線形解析において、その計算の中で下記のような5重のループがありました。

ループの外側から順に

1) Loop over time steps
2) Loop on geometry (matching and local base)
3) Loop on the threshold of friction (transformation of the Coulomb problem into a series of Tresca problems)
4) Loop on the status of contact
5) Newton loop (linearization of the problem)

です(技術的に正確に訳す自信が無いので、英語のまま引用しました)。

そのうち、2,3,4が接触・摩擦の計算に関するものですが、version 11.0.28 からは「FORMULATION = ’CONTINUE’」の場合において、3と4のループを省略することが可能になりました。これは、「partial generalized Newton’s method」という新しいアルゴリズムによるものだそうです。

更に、version 11.2.4からは「total Newton generalized method」によって2のループも省略できるらしいのですが、最新のSalome-Mecaである「2012.2」に含まれているのは「version 11.2.0」なので、Code_Asterを単独でインストールする必要があります。

なお、この「generalized Newton」を適用するには、「DEFI_CONTACT」にて下記の設定を行います。

ALGO_RESO_GEOM=’NEWTON’
ALGO_RESO_FROT=’NEWTON’
ALGO_RESO_CONT=’NEWTON’

現在のversion11ではこれがデフォルトになっているそうですが、一応、確認しておくほうが良いと思います。

さて、前回Code_Aster10.5で実施した「回転体とゴムの接触」について、今回の新アルゴリズムにて再計算してみました。

前回使用したDEXCS2011-Salome-A3(64)(SALOME6.3 + Code_Aster 10.5)と同じマシン(仮想PC)に、debian6とSalome-Meca2012.2(SALOME 6.5 + Code_Astr 11.2.0)を入れて、それを用いました。

計算に用いたメッシュファイル(mmed)は前回と同じですが、commファイルはCode_aster11に合わせて多少修正しています。

まず、今回の本題である接触ですが、下図のように「generalized Newton」の使用を明示的に記述しています。

Rigid03201

ただし、この問題は摩擦が無いので、ALGO_RESO_FROT=’NEWTON’はありません(摩擦を含む問題の場合は、これも必要になります)。
また、Code_Asterのバージョンが11.2.0なので、ALGO_RESO_GEOM=’NEWTON’の指定もありません。結局有効なのは、ALGO_RESO_CONT=’NEWTON’のみになります。

また、「CALC_NO」と「CALC_ELEM」をまとめて「CALC_CHAMP」で置き換えました。

Rigid03202


さて、計算時間ですが、下記の通りになりました。

前回(version10.5)   : 13141 (秒)
今回(version11.2.0):  4897 (秒)

同じマシンで計算し、計算時間が前回の37%程度まで短縮されました。
これはかなり大きな効果と言えます。

前回(version10.5) のループのログの一部です(messファイル)。

Rigid03203

表を左から見て、GEOM, CONT, NEWTONのループが回っている事が分かります。

今回(version11.2.0)の方です。

Rigid03204

GEOMとNEWTONのループだけで、CONTのループがありません。
その分、速く計算出来ているようです。

Code_Asterのバージョンは、10が安定版で11がテスト版との位置づけですが、これだけの速さの違いがあるのならば、接触・摩擦問題については安定版にこだわらず、特に問題が無い限りバージョン11を使ったほうが良さそうです。

今回のメッシュファイル(mmed)とcommファイルは下記です。
rigid03-2.zip

より以前の記事一覧