論文 Original Paper
エンジン減衰およびクランク軸装着ダンパプーリーの 動的特性値を考慮したねじり振動波形の
モデリングと数値シミュレーション方法に関する研究
児 玉 知 明
＊1，本 田 康 裕
＊2A Study on the Numerical Simulation Method and Modeling of Torsional Vibration Waveforms Considering Engine Damping and
Dynamic Properties of Engine Crankshaft System with Rubber Damper Pulley
Tomoaki Kodama
＊1， Yasuhiro Honda
＊2Abstract: In this paper, we first describe the dynamic properties of rubber parts of rubber damper pulley that are necessary for the modeling and numerical simulation of torsional vibration. Secondly, we describe an experiment in which two crankshaft pulleys with a torsional rubber damper pulley are fitted to a 6-cylinder, high-speed diesel engine. Torsional waveforms of the rubber damper inertia ring and the pulley are measured by means of phase-shift torsiograph equipment. The measured waveforms are harmonically analyzed and the dynamic properties of the torsional stiffness and the torsional damping are investigated from an experimental viewpoint. As a results of comparisons with experimental data, certain dynamic properties of damper pulleys with a torsional rubber damper have been clarified. The model used for the numerical simulation of the torsional vibration is a multi-degree- of-freedom equivalent torsional vibration system. The tension acting on the damper pulley and the rotational resistance of the alternator, the cooling water pump, the valve train system, etc., as well as the frictional resistance of other accessories, were considered. Moreover, as a numerical simulation method of the torsional vibration waveforms, a transition matrix method is adopted. The validity of this method has been confirmed from comparison and examination of measured values of torsional vibration and simulation values results.
Key words: Rubber Damper Pulley, Engine Damping, Torsional Vibration, Forced Vibration, Free Vibration, Dynamic Stiffness, Damping Coefficient, Torsional Vibration Waveform, Experimental, Numerical Simulation, Torsional Angular Displacement, Transition Matrix Method
１．は じ め に
近年，自動車は性能向上，低燃費化に加えて，低公害 化，高品質化がこれまで以上に要求されている。これに 伴って，自動車用エンジンは，振動および騒音を低減 し，音色を改善することがエンジン設計時における重要
な課題の一つとなってきた。これへの対策として，高性 能ねじり振動ダンパ［1］−［3］，ねじりおよび曲げ振動を抑 制するゴムダンパ付クランクプーリー［4］，あるいはゴ ムねじり振動ダンパ付フライホイール［5］が自動車用エ ンジンに採用されている。これらのうち，最も広く採用 されているのはゴムねじり振動ダンパ（以下，ダンパと 呼ぶ）である［11］， ［15］， ［20］， ［21］。
しかし，このような種々の形式のダンパ付エンジンク ランク軸系のねじり振動角変位をエンジン設計時に予測 する際に問題となるのは，
［1］ エンジン部のねじり振動減衰，および
［2］ ダンパ部の動的ねじりこわさおよび減衰値の見積り
＊1 国士舘大学理工学部, 博士（工学），kodama@kokushikan.
ac.jp, Department of Science and Engineering, School of Science and Engineering, Kokushikan University. Dr. of Engineering
＊2 国士舘大学理工学部, 工学博士，honda@kokushikan.ac.jp, Department of Science and Engineering, School of Science and Engineering, Kokushikan University. Ph. D..
である［6］−［15］， ［19］−［21］。これらの特性は不明確で，正し く推定できる式がこれまで存在していないし，定量的な データも少ない。これらの値をある程度正確に見積らな いと満足のいく計算精度は得られない［16］−［18］。
そこで，本論文は，まず最初にねじり振動ダンパを装着 しないクランク軸系を対象として，著者らが提案した実験 法を採用して求めた自由振動実験結果（文献［6］，［11］，
［14］−［16］にデータのみを報告）に新たに検討を加え，エ ンジン部のねじり振動減衰特性について考察する。次に，
これまで広く採用されているゴムねじり振動ダンパ付クラ ンクプーリー（以下，ダンパプーリーと呼ぶ）を装着した エンジンクランク軸系を対象として，ねじり振動角変位 測定を行い，防振ゴムの複雑な動的特性の一端を実験的 に明らかにする。最後に，これらの実験データを採用し て，著者らが開発した「推移マトリックス法（Transition matrix method）を導入したクランク軸系のねじり振動 波形のシミュレーション法」［6］， ［11］， ［14］−［16］， ［20］， ［21］を採 用して，得られた実験結果に解析的な検討を加える。
２．エンジン減衰およびダンパプーリーの動的特 性値の測定装置および実験方法
本章は供試エンジンおよび供試ダンパプーリーの主要
諸元ならびにねじり振動波形の測定方法について述べ る。なお，3章で言及する実験データに関する実験方法 は，文献［6］，［11］，［14］−［16］，［20］，［21］で詳述 してあるので，ここでは詳しい説明を省略する。自由減 衰ねじり振動波形は，全負荷運転中に吸気管通路を閉鎖 して新たな空気の流入を止めるかあるいはデコンプ
（Decompression device）を使用して減圧して，燃焼を 急に停止させることにより発生させる［14］−［16］，［20］，［21］。 これを，自由減衰ねじり振動測定（ねじり減衰定数測 定）と定義する。
２. １ 供試エンジンおよびダンパプーリーの主要諸元 本研究に使用した供試エンジンは，直列6シリンダー 高速ディーゼルエンジンであり，その主要諸元をTable
1に示す。また，供試ダンパプーリーはFigure 1に示す
形状および寸法をもち，Table 2は諸元を示す。この表 の中には，静的試験および自由振動試験より得られたね じりばね定数および固有振動数がそれぞれ示されてい る。タイプAおよびBとは，ゴム硬度は異なるが，形状 係数の影響を除くためにゴム形状および寸法は同一であ る。また，ゴム材質はTable 2に示すように同様である。
Table 1 Main specifications of the test engine.
Figure 1 Dimensions and shapes of the test rubber damper pulley.
２. ２ ダンパプーリー装着エンジンクランク軸系のね じり振動波形の測定方法
供試エンジンはユニバーサルジョイントを介して渦電 流式電気動力計を接続し，クランク軸前端部（フライホ イールと反対側）に供試ダンパプーリーを装着してい る。ダンパ慣性リング部およびプーリー部にパルス発生 用歯車を取付けて，電磁式ピックアップにより回転数に 比例した周波数信号を取り出す。取り出された信号は，
平均角速度（中心周波数）を演算するアダプターを経て 位相差形ねじり振動計に入力され検出周波数と中心周波 数との関係からねじり角度が演算されてねじり振動波形 を得る。得られた両波形は，次数ごとに動的特性を検討 するためにF.F.T.アナライザーを用いてスペクトル解析
される。測定はエンジン回転数800 ［r/min］から3200
［r/min］の間で全負荷運転にて行う。また，後述のね じり振動解析に必要なインジケーター線図は，ピエゾ式 インジケーターを用いてプーリー側から第6番目シリン ダ ー の 内 圧 を 採 取 し て 求 め る。 実 験 装 置 の 概 略 を Figure 2に示す。
３．エンジン部のねじり振動減衰特性（ねじり減 衰係数および減衰比の算出）
３. １ 減衰定数の定義および減衰比との関係
自由減衰ねじり振動波形から減衰定数および減衰比を 求める方法について述べる。この自由減衰振動波形記録
から直接求めることができる 値を次式で定義し，
Table 2 Parameters of the test rubber damper pulley.
Figure 2 Schematic diagram for measuring of amplitude of torsional vibration angular displacement ［Measurement of relative torsional vibration］．
これを減衰定数（Decay constant）と呼ぶ。
：（Decay constant） （1）
ここで，対数減衰率： と自由減衰ねじり振動の周 期： は実測ねじり波形記録より求まるので，この 式より の値が決定される。
次に，減衰定数： と減衰比 ：および固有振
動数： との間には，次の関係が成立する。
（2）
ここで， および の値はねじり減衰実験より
既知となるので，減衰比： が求まる。
３. ２ エンジン部のねじり振動減衰特性
エンジン部のねじり振動減衰特性を検討するために，
ピストンリングの装着本数を変化させて自由ねじり振動 実験を行い，前述の減衰定数の特性を調査する。
同一のエンジンであれば，式（2）の関係から，減衰定 数の特性を調査することにより，減衰比等の減衰特性を 調べることになる。本節で対象とする直列6シリンダー エンジン（総排気量10.2 ［litter］，最高出力：143/2400
［kW/r/min］，最大トルク：667/1600 ［Nm/r/min］）は 各ピストンにコンプレッションリング3本とオイルリン グ2本が装着されている。そこで，このエンジンのすべ てのピストンについて，
［1］ 全ピストンリング付，
［2］ No.4のオイルリングを抜く，
［3］ No.2およびNo.3のコンプレッションリングを抜く，
［4］ No.2，No.3およびNo.4のリングを抜く，
の各条件について自由ねじり振動実験を行って求めた減 衰定数の結果をTable 3に示す。なお，ピストンリング はシリンダヘッド側から順に番号を付けた。この実験結 果のみについては，前述したように文献［6］，［11］，
［14］−［16］，［20］，［21］ですでに報告したが，以下に 新たな考察を加える。
この結果よりピストンリング部はエンジン減衰に大き な影響を及ぼしていることが推定できる。特に，オイル リング部の影響が大きいと推定される。そこで，オイル リング部およびコンプレッションリング部等で分担して
いる減衰仕事割合について検討する。オイルリング部
（No.4， No.5），コンプレッションリング部（No.1〜No.3）
およびリング部以外でも，おのおの一定の減衰仕事をし ていると仮定すると，Table 3に示された四つの実験結果 より各減衰仕事割合に関する代数方程式が得られる。こ の表中の上から実験［1］， ［2］， ［3］ および［4］とすると，
実験［1］と［3］，ならびに実験［2］と［4］との差は，
それぞれコンプレッションリング2本分，また実験［1］
と［2］，および実験［3］と［4］との差は，それぞれオ イルリング1本分の減衰定数の相違からでてきた結果で ある。以上を考慮してこれらの代数方程式を解くと，コ ンプレッションリングおよびオイルリングについて，そ れぞれ二つの値（ほぼ同値）が求まるので，それらの値 を平均してコンプレッションリングおよびオイルリングの 減衰定数を決定した。その結果，コンプレッションリン グ1本分の値は，4.20 ［1/s］，およびオイルリング1本分 の値は，7.20 ［1/s］となる。これによりオイルリング１ 本分の減衰仕事はコンプレッションリング1本分の減衰 仕事の約1.7倍となる。また，ピストンリング部全体で全 減衰仕事の約60.0 ［%］を分担している結果が得られる。
次に，ピストンリング部以外で分担している約40.0
［%］の減衰仕事の内容について検討する。Table 4はこ れまで提示されたヒステリシス減衰に関する実験式のう ちで比較的信頼性のあるルイス［7］の式を使用してヒス
テリシス減衰仕事： ， さらに減衰仕事割合：
， （ ：強制仕事）を計算した結 果を示す。（式（3）参照）。これよりヒステリシス減衰 仕事の割合は10.0 ［%］以下であると予測されるので，
30.0 ［%］の減衰仕事が軸受部等の機械損失などで消費 されていると推定される。
・100［%］ （3）
ところで，最近の高速ディーゼルエンジンは，ピストン リングの本数が以前と比較して少なくなりコンプレッシ ョンリング2本とオイルリング1本を装着している場合 が多い。このような場合には，前述の減衰仕事割合が異 なりピストンリング部に占める仕事割合は当然低下す る。オイルリング1本分の減衰仕事はコンプレッション Table 3 Relationship between the number of piston rings and decay constant .
［No. 1 to No. 3：compression piston rings、 No. 4 and No. 5：oil-control piston rings］．
リング1本分の減衰仕事の約1.7倍であることと，ピス トンリング部全体での減衰仕事は全減衰仕事の約60.0
［%］を分担しているという前述の結果を前提として再 計算すると，ピストンリングの本数が3本の場合には，
ピストンリング部全体での減衰仕事は全減衰仕事の約 45.0 ［%］を分担しているという結果が得られる。した がって，ピストンリング部以外で全減衰仕事の55.0 ［%］
を分担していることになる。また，ヒステリシス減衰仕 事の分担割合は15.0 ［%］以下となる。
４．ダンパプーリゴム部の動的特性
４. １ 実測ねじり振動波形2.1節で述べた2種類の供試ダンパプーリーをそれぞ れ装着して，ダンパ慣性リング部およびプーリー部にお けるねじり振動波形を測定した結果の一例として，
Figure 3は6次ねじり振動共振点近傍2635 ［r/min］に おける1周期2回転分の実測ねじり振動波形を示す。こ の図中には，波形の内容を調べるために調和解析して求 めたおもな次数の片振幅値が示されている。
6次共振点近傍で測定しているので，プーリー部での 実測波形は6次振動が大きく増幅された波形である。一 方，ゴムを介した慣性リング部の実測波形は6次成分が 一番大きいが，3次および4.5次も大きい成分のため，6 次成分だけが特に大きく増幅されたプーリー部での実測 ねじり波形とは多少傾向を異にする。
４. ２ 実測ねじり角変位振幅曲線
実測ねじり角変位振幅曲線の一例としてFigure 4お よびFigure 5は，ダンパプーリーを装着しない場合の プーリー端ならびにBタイプダンパプーリーを装着した 場合のダンパ慣性リング部とプーリー部における実測ね じり波形を調和解析して得られる各次数ごとの振幅曲線 をそれぞれ示す。これらの図中には主要次数の振幅曲線 のみを示す。Figure 4のダンパプーリーを装着しない場 合は，大きな6次振動の共振点が2524 ［r/min］にあり，
共振振幅は14.9×10−3 ［rad］である。次に，この6次振 動に着目すると，Bタイプダンパプーリーを装着した場 合は多少の振幅低減効果はあるが，動力計を含めた全エ Table 4 The ratio of hysteresis loss to exciting torque energy ［10.2 ［litter］ In-line 6 cylinder diesel engine］．
Figure 3 Measured waveform of torsional vibration angular displacement.
（a）Waveform of damper inertia ring. （b）Waveform of damper pulley.
ンジンクランク軸系では最適調整の状態でない。
４. ３ 実測波形からの動的特性値の算出
ダンパゴム部の動的ねじりこわさおよび減衰係数は，
ダンパ慣性リング部およびプーリー部の実測波形の調和 解析結果から算出される。
ダンパ慣性リング部の運動方程式は，
（4）
で示される。
そこで， ，
と表し，上式（4）に代入して ， について
整理すると次式が求めらる。
，
（5）
ここに， である。
式（5）に実測ねじり振動波形を調和解析して得られ る各次数振動ごとのダンパ慣性リングとプーリーとの振
幅比： および位相差： の値を代入すれば，
， は既知であるので，単次数振動ごとに動
的ねじりこわさ： および減衰係数： が算
出される。Table 5は求められた ， の値お よびこれらの値を使用して算出された絶対ねじりこわ
さ： と損失係数： を示す。
防振ゴムの動的特性は，ダンパゴム部の形状，材質を 同じとすると，主として，
［1］ 温度，
［2］ 振動数（周波数），
［3］ 平均ひずみとひずみ振幅，
に依存する。この3因子のうち1因子の影響を調査する ときは，ほかの2因子を一定に保って実験する必要があ る。エンジン実験では，ゴム部表面温度を313.0 ［K］に 保ったが，ほかの2因子はエンジンクランク軸系の振動 特性上，それぞれが変化し一定に保つことができない。
そこで，両因子を考慮したひずみ振幅と角振動数の積で あるひずみ速度： を次式で定義してTable 5の実験 結果を整理する。
（6）
ここに，
（7）
Figure 4 Measured amplitude curves of torsional vibration angular displacement ［without rubber damper pulley］．
（a） Measured amplitude curves of torsional vibration angular displacement ［Inertia ring of rubber damper pulley］．
Figure 5 Measured amplitude curves of torsional vibration angular displacement ［with B-type rubber damper pulley］．
（b） Measured amplitude curves of torsional vibration angular displacement ［Pulley of rubber damper］．
上式は均一円形中空断面において，相当半径： の 位置のせん断応力：（＝一定）
（＝Constant） （8）
が断面全体に均一に分布する場合と円周方向に単純せん 断変形を受ける際のせん断応力分布：
（9）
の場合とのねじりモーメントが等しいとすることにより 求まる。
Figure 6は絶対ねじりこわさ： とひずみ速
度： ，Figure 7は損失係数： とひずみ速度：
との関係を示す。非線形性をもつダンパゴム部の 動的特性を線形近似して周波数や振幅に対する特性変化 を示すことを目的として，両図はタイプAおよびBのダ
ンパゴム部の絶対ねじりこわさ： と損失係数：
をひずみ速度： で表現することを試みた図で ある。第5章の数値シミュレーション計算に必要なダン パゴム部の動的ねじりこわさ： および減衰係数：
（ でまとめる。）について
も同様に整理すると，それぞれFigure 6およびFigure 8となる。なお，絶対ねじりこわさは
（10）
の式で計算できるが， の値と比較して
の値が小さいので，
（11）
となる。したがって，図中で整理するとFigure 6に示 したように両者はほぼ一致する。
ダンパゴム部のひずみ速度： が50〜950 ［1/s］の Table 5 Strain rate ， dynamic torsional stiffness and damping coefficient ［（a） part I］．
Table 5 Absolute torsional stiffness ， loss factor and imaginary part of complex torsional stiffness ［（b） part II］．
範囲内で，絶対ねじりこわさ： はダンパプーリ ーAで1.7×105〜8.0×104， ダンパプーリーBで1.7×
105〜1.0×105 ［Nm/rad］（ひずみ速度が増加すると絶 対ねじりこわさは低下する傾向にある。（Figure 6参 照））である。また，損失係数： はおのおの2.0×
10-1〜3,5×10-1，3.3×10-1〜1.0×100（ひずみ速度が増 加すると損失係数も増加する傾向にある。（Figure 7参 照））である。さらに，複素ねじりこわさの虚部：
（ ）は3.0×104〜2.0×104，6.5×
104〜6.0×104 ［Nm/rad］（ひずみ速度が増加しても複 素ねじりこわさの虚部はほぼ一定である。（Figure 8参
照））である。
Figure 6〜Figure 8の図中のカーブフィットした曲線 に対して相関係数を求めると，ダンパAに関しては0.82
〜0.89，またダンパBは0.86〜0.94となった。
５． ダンパプーリー付クランク軸系のねじり振動 特性に関する解析的な検討
５. １ ねじり振動波形の数値シミュレーション法 著者らが提案した「推移マトリックス法を導入したク ランク軸系のねじり振動波形のシミュレーション法」は ティラー級数を間接使用したOne-step数値計算法の一 種である。この数値計算法における詳細な式の誘導は，
文献［6］，［11］，［14］−［16］，［19］−［21］に示してあ るので数値シミュレーション計算をするのに必要な最終 段階の式（ティラー級数は4階導関数まで考慮）のみを 次に示す。
（12）
ここで，マトリックスを構成している ， および
（ ）は （ ：質点数）の部分マ
トリックスである。これらの部分マトリックスは，等価 ねじり振動系の各要素の諸元とステップサイズよりなっ ているので，これらの値が与えられることにより既知と なる。また，強制トルク： およびその微分の
， はインジケータ線図および往復運動部
質量等が与えられれば既知となる。したがって，各質点
のねじり角変位： ，角速度： の初期値を与
えてやると ， が求まり以下順次連続し
て計算することができる。最初の2，3周期は立ち上が り問題があり定常状態を示さないので，計算結果はそれ Figure 6 Relationship between absolute torsional stiffness
（dynamic torsional stiffness） and strain rate of damper rubber .
Figure 7 Relationship between loss factor and strain rate of damper pulley .
Figure 8 Relationship between imaginary part of complex torsional stiffness and strain rate of damper rubber .
以降を採用する。
なお，自由減衰ねじり振動波形を計算する場合は，強 制トルクが零となるので，右辺の強制トルクに関係する 第2項および第3項を除いて計算すればよい。
５. ２ ねじり振動波形計算に必要な入力データ
Figure 9は5.1節で言及した数値解析法にしたがって
置換されたダンパプーリー付クランク軸系の等価ねじり 振動系を示す。この等価振動系は，従来の方法にしたが って置換することができる。ダンパプーリーゴム部はフ ォークトモデルを採用する。Table 6はこの等価振動系 の諸元である。ねじり振動波形計算で不明確な入力デー タとしては，前述したように，エンジン部の減衰係数お よびダンパプーリーゴム部のねじりこわさおよび減衰係 数の値である。
まず，エンジン部の減衰係数： （ と仮定）
について言及する。Figure 10は文献［6］−［21］で示した直 列小形高速ディーゼルエンジン6種類（シリンダー径：
83.0〜120.0 ［mm］，総排気量：1.9〜10.2 ［litter］）を対 象として自由ねじり振動実験から求められた減衰定数：
，固有振動数： および減衰比： の関係
である。あらかじめ神田法［10］で と の関係 を求めておき，ホルツァ法で を求めFigure 10か
ら に相当する を算出し，神田法の結果から
を逆算することにより，Figure 9のエンジン部の減 衰係数を決定する。その結果，エンジン部の減衰係数は Table 6に示すように6.86 ［Nms/rad］（＝一定）となる。
次に，ダンパプーリーゴム部のねじりこわさ：
および減衰係数： の値の見積りについて述べる。
波形計算はエンジン回転数を変えて行うので，各エンジ
ン回転数における支配的次数の角変位振幅： ，
， ：（ ）の実測結果
をもとにしてひずみ速度： を計算し，Figure 6お よびFigure 8より動的ねじりこわさ： および減
衰係数： の値を求めて入力データとして使用す
る。Figure 6およびFigure 8はすべての実験データを代
Figure 9 Equivalent torsional vibration system of an engine system with a rubber damper pulley （Numerical simulation model of test diesel engine crankshaft system with rubber damper pulley）．
表して一本の曲線でまとめたものであるので，これらの
および の値は代表値となる。数値計算シミ
ュレーション方法の詳細についてはAppendix Iを参照。
５. ３ 数値計算結果およびその検討
ダンパプーリーを装着しないクランク軸系の計算結果 として，Figure 11数値計算波形を調和解析して得られ る主要次数のねじり振動角変位振幅曲線を示す。この図 はFigure 4の実測より求められたねじり振動角変位振
幅曲線と対応する。各次数振動の共振回転数および共振 振幅は，ほぼ一致している。このことよりエンジン部の ねじり振動減衰係数の値は妥当な値であると考えられる。
次に，ダンパプーリーを装着したエンジンクランク軸 系の計算結果として，Figure 12はダンパ慣性リング部 およびプーリー部の数値計算波形を調和解析して得られ る主要次数のねじり振動角変位振幅曲線を示す。この図 はFigure 5の実測ねじり振動角変位振幅曲線と対応す る。これらの数値計算結果と対応する実測結果とを比較 Figure 10 Relation among decay constant ， natural angular frequency and damping ratio .
Table 6 Numerical calculation values of the equivalent torsional vibration system ［with a rubber damper pulley］．
すると，計算精度上多少の誤差はあるが実用上は有効な 数値計算結果であると判断できる。
６．お わ り に
ダンパプーリー付クランク軸系のねじり振動計算をす る際に問題となるエンジン部のねじり振動減衰特性およ
びダンパプーリーゴム部の動的特性について実験と数値 シミュレーションより検討を加えた。その結果は以下の とおりである。
［1］ 最近の高速小形ディーゼルエンジンはピストンリン グの本数を少なくする傾向にあるので，全減衰仕事 に占めるピストンリング部の割合は，これまでより 低下しているが，依然として大きいと推定される。
［2］ 同一ゴム材質で，温度および形状係数を一定とし た場合のダンパプーリーゴム部の絶対ねじりこわ さおよび減衰係数をひずみ速度で表現することを 試みた。
［3］ ダンパゴム部のひずみ速度： が50〜950 ［1/s］
の範囲内で，絶対ねじりこわさは低下，損失係数は 増加，複素ねじりこわさの虚部はほぼ一定である。
［4］ 推移マトリックス法を導入したねじり振動波形の 数値シミュレーション法は，プーリーダンパ付エン ジンクランク軸系のねじり振動計算に有効である。
記号の定義
本論文に使用するおもな記号を以下に示す。
：ダンパゴム部相当半径の位置におけるゴム厚さ ［m］
：ダンパゴム部の減衰係数 ［Nms/rad］
：ダンパ慣性リング部の慣性モーメント ［kgm2］
：ダンパゴム部の動的ねじりこわさ ［Nm/rad］
：ダンパゴム部の複素ねじりこわさ ［Nm/rad］，
，
：虚数）
：ダンパゴム部の絶対ねじりこわさ ［Nm/rad］，
：損失係数，
：振幅比，
：減衰ねじり振動の周期 ［s］
：ねじり振動の強制仕事 ［Nm］
：ヒステリシス減衰仕事 ［Nm］
：減衰仕事割合 ［%］，
×100
， ：ダンパゴム部の外径および内径 ［m］
：減衰定数 ［s-1］
：ひずみ速度 ［1/s］
：減衰比
： プーリー部（クランク軸前端部）のねじり振動 Figure 11 Numerical calculation amplitude curves of torsional
vibration angular displacement ［without rubber damper pulley］．
（a） Numerical calculation amplitude curves of torsional vibration angular displacement ［Inertia ring of rubber damper pulley］．
Figure 12 Numerical calculation amplitude curves of torsional vibration angular displacement ［with B-type rubber damper pulley］．
（b） Numerical calculation amplitude curves of torsional vibration angular displacement ［Pulley of rubber damper］．
角変位 ［rad］
：ダンパ慣性リングのねじり振動角変位 ［rad］
：ダンパゴム部の相対ねじり振動角変位 ［rad］，
：対数減衰率
：ダンパゴム部の相当半径 ［m］
：プーリー部とダンパ慣性リング部との位相差 ［rad］
：ねじり振動の角振動数 ［rad/s］
：固有円振動数 ［rad/s］
参考文献
［ 1 ］ 岡村・沼尻：高粘度ビスカスダンパの開発，自動車技術，
第40巻，9号 （1986） P. 1140-1146.
［ 2 ］ 引野・新井野・笠原：大型商用車用ディーゼルエンジン のトーショナルダンパについて，自動車技術，第38巻，
12号 （1984） P. 1440-1147.
［ 3 ］ 川村・山口・田坂・中窪・田中：クランク軸トーショナ ルダンパプーリの性能予測，自動車技術，第40巻，12号
（1986） P. 1581-1587.
［ 4 ］ 神谷・熱海・田坂：新型ダンパ付きクランクプーリの開 発，自動車技術，第41巻，13号 （1987） P. 1557-1561.
［ 5 ］ 加藤・角・工藤・小林：トーショナルダンパ付きフライ ホイールの開発，自動車技術，第39巻，12号 （1985） P.
1428-1432.
［ 6 ］ 若林・関：クランク軸系のねじり振動波形のシミュレー ション，日本舶用機関学会誌，第11巻，11号 （1976） P.
808-818.
［ 7 ］ Lewis, F. M., Trans. Soc. Nav. Archit, No.23 （1925） P.
109.
［ 8 ］得丸：防振ゴム，日本機械学会誌，第64巻，505号 （1961）
P. 90-100.
［ 9 ］ 関口：防振ゴムの設計法，日本機械学会論文集，第27巻，
173号 （1961） P. 85-95.
［10］ 神田：多数の粘性減衰要素を有するクランク軸系ねじり 非規準振動モードの電子計算機による研究，日本機械学 会論文集，第32巻，235号 （1966） P. 465.
［11］ Tomoaki Kodama, Yasuhiro Honda, Katsuhiko Wakabayashi, Shoichi Iwamoto：A Calculation Method for Torsional Vibration of a Crankshafting System with a Conventional Rubber Damper by Considering Rubber Form, SAE 1996 International Congress and Exposition, SAE Technical Paper Series No. 960060, （1996） P. 103- 121.
［12］ Yasuhiro Honda, Katsuhiko Wakabayashi, Tomoaki Kodama, Shoichi Iwamoto：Effect of Rubber Hardening by Secular Change on Properties of Vibration Proof Rubber for Torsional Vibration Dampers, SAE 1996 International Congress and Exposition, SAE Technical Paper Series No. 9600139, （1996） P. 1-14.
［13］ Yasuhiro Honda, Katsuhiko Wakabayashi, Tomoaki Kodama, Shoichi Iwamoto：An Experimental Study on Dynamic Properties of Rubber Test Specimens for Design of Shear-type Torsional Vibration Dampers, The Ninth International Pacific Conference on Automotive Engineering, Volume 2 No. 971400 Abstract Code 88,
（1997） P. 217-222.
［14］ Katsuhiko Wakabayashi, Yasuhiro Honda, Tomoaki
Kodama, Hiroshi Okamura：Torsional Vibration Analysis of Internal Combustion Engine Crankshaft System by Transition Matrix Method, 31st International Symposium on Automotive Technology and Automation （31st ISATA）， Simulation Virtual Reality and Supercomputing Automotive Applications, 98VR00036, （1998） P. 253-265.
［15］ Yasuhiro Honda, Katsuhiko Wakabayashi, Tomoaki Kodama, Takashige Kasai, Hiroshi Okamura：Simulation Method for Torsional Vibration of Diesel Engine Crankshaft with a Rubber Damper in Consideration of Temperature Dependency, 31st International Symposium on Automotive Technology and Automation （31st ISATA）， Simulation Virtual Reality and Supercomputing Automotive Applications, 98VR00037, （1998） P. 1-13.
［16］ Katsuhiko Wakabayashi, Yasuhiro Honda, Tomoaki Kodama, Hiroshi Okamura：Simulation Method of Torsional Vibration Waveform of Reciprocating Engine Crankshaft System by Transition Matrix Method, ASME 1999 International Design Engineering Technical Conference & The Computers and Information in Engineering Conference, 17th Biennial Conference on Mechanical Vibration and Noise Symposium on Dynamics and Vibration of Machine Systems, Session VIB-00037, Dynamics and Vibration of Machine System 2, DETC99/VIB-08124, （1999） P. 1-12.
［17］ Yasuhiro Honda, Katsuhiko Wakabayashi, Tomoaki Kodama, Hiroshi Okamura：An Experimental Study of Dynamic Properties of Rubber Specimens for a Crankshaft Torsional Vibration of Automobiles, ASME 1999 International Design Engineering Technical Conference & The Computers and Information in Engineering Conference, 17th Biennial Conference on Mechanical Vibration and Noise Symposium on Dynamics and Vibration of Machine Systems, Session VIB-00037, Dynamics and Vibration of Machine System 2, DETC99/VIB-08125, （1999） P. 1-14.
［18］ Yasuhiro Honda, Tomoaki Kodama, Katsuhiko Wakabayashi：Relationships between Rubber Shapes and Dynamic Characteristics of some Torsional Vibration Rubber Dampers for Diesel Engines, The 15th Pacific Automotive Engineering Conference – APAC15, No.0320,
（2009） P. 01-12.
［19］ Tomoaki Kodama, Yasuhiro Honda: A Study on Rubber Shapes and Dynamic Characteristics of some Torsional Vibration Shear Type Rubber Dampers for High Speed Diesel Engine Crankshaft System, Transactions of the Kokushikan University School of Science and Engineering, No.11 （2018） March P.01-10.
［20］ Tomoaki Kodama, Yasuhiro Honda: A Study on the Modeling Consideration Dynamic Properties of Vibration Damper Rubber Parts, Proceedings of 2018_The 7th International Conference on Engineering and Innovative Materials ICEMI 2018, EM028 （2018） P.01-07.
［21］ Tomoaki Kodama, Yasuhiro Honda: A Study on the Modeling and Simulation Method of Torsional Vibration Considering Dynamic Properties of Rubber Parts for Engine Crankshaft System with Rubber Damper Pulley, Selections of 2019_The 8th International Conference on Engineering and Innovative Materials ICEMI 2019, IM067 （2019） P.01-07.
Appendix I
General formula derivation A-1. Notation
The following are the main symbols used in the present paper.
：Main torque ［Nm］
：j/2-th order component of resultant excitation torque ［Nm］
：Viscous damping coefficient proportional to the absolute velocity of mass ［Nms/rad］
：Viscous damping coefficient proportional to the relative velocity between adjacent masses ［Nms/rad］
：Cylinder bore ［m］
：Resultant excitation torque acting on m-th mass ［Nm］
：Rank of Taylor series
：Order of torsional vibration in four-stroke cycle engine （j=1,2,3,4･･････）
：Mass polar moment of inertia ［kgm2］
：Torsional stiffness of shaft ［Nm/rad］
：Positive integer
：Center-to-center distance between large and small ends of connecting rod ［m］
：Mass number
：Mass of reciprocating part ［kg］
：Gas pressure in cylinder ［Pa］
：Crank radius ［m］
：Time ［s］
：Step size
：Torsional angular displacement of mass ［rad］
：Connection rod/crank radius ration （l/r）
：Firing interval between the m-th cylinder and 1-st cylinder ［rad］
：Phase angle of j/2-th order component of resultant :torque ［rad］
：Revolving angular velocity of crankshaft system ［rad /s］
：Damped free circular frequency ［rad/s］
Matrices and vectors
：Damping matrix
：Unit matrix
：Inertia matrix
：Inverse inertia matrix
：Stiffness matrix
：Partial matrix
：Resultant torque column vector
：Components of resultant torque column vector （T：Transposition）
：Angular displacement column vector
：Components of angular displacement column vector A-2. Numerical simulation method of resultant torsional vibration waveform A-2-1. Derivation of general expressions
If it is assumed that the engine crankshaft system retains a sufficient linearity in inertia and rigidity within the range of such deformations as to have an engineering significance and that the engine damping in the shaft does not have so high a degree of nonlinearity the shaft system may be linear damping by the equivalence of dissipation of energy during one period as is evident from the previously reported results of studies
From what we have just discussed, torsional vibration numerical calculations will be made by replacing the engine crankshaft system with such equivalent torsional vibration system as shown in Figure A-1. The equation of motion for the m-th mass can be expressed as follows.
（A-1）
In the case of four cycle straight type engines, the resultant excitation torque acting on the m-th mass can be expressed by the following equations.
（A-2）
（A-3）
In equation （A-3）， the resultant excitation torque has been expanded into a Fourier series to express it in the form of the sum of each order components.
Considering equations of motion for all masses （ ： number of masses） from equation （A-1） they can be expressed in matrix as follows.
（A-4）
Where, Inertia matrix
Damping matrix
Stiffness matrix
Angular displacement column vector
Resultant torque column vector
（T: Transposition）
（A-5）
Now, if the resultant torsional vibration in the engine crankshaft system at an arbitrary engine speed is to be numerical calculation from the n elements 2-nd derivative ordinary differential equation （A-4） is solved for each order and compounded, taking phase difference into considered.
［a］ Assuming that the each order components of the resultant torque, which can be obtained by equation （A-3），
act independently, equation （A-4） is solved for each order and compounded, taking phase difference into consideration.
［b］ The resultant excitation torque shown by equation （A-2） is directly applied to the torsional vibration system as external forces and the direct simulation using numerical calculation is carried out.
In the case of ［a］， it is necessary to make torsional vibration calculation for each order and moreover to sum up all of them, taking phase differences into consideration, and therefore errors are apt to occur during the process of such numerical calculations. The method ［b］ is better suited for the numerical calculation of resultant torsional vibration waveforms. In the case of ［b］， several methods of numerical calculation are conceivable. In the present paper, the transition matrix method which is a kind of step by step numerical calculation method making indirect use of Tayor series is used to solve equation （A-4）．
In this method, numerical calculation are made, taking into consideration the constant coefficient in equation （A- 4）．Therefore, it has the advantage that the time required for making numerical calculation can be so much reduced.
Equation （A-4） can be rewritten
（A-6）
If equation （A-6） is successively differentiated with respect to time, taking the constant coefficient into consideration, the n-th derived function （n>2） of can be written as follows.
（A-7）
Now if the angular displacement at the time and are replaced by and respectively and the angular velocities by and , and are expand into a Taylor series, considering up to the i-th derived function, they can be expressed by equation （A-8），
（A-8）
Where in is the step size.
If equation （A-7） is repeatedly applied to （A-8） the second and higher derivatives on the right hand side of these equations may be written in terms of the angular displacements and velocity .
（A-9）
Where,
（A-10）
A-2-2. Derivation of equations taking Taylor series up to the fourth derived function into consideration
It has been empirically shown that a considerably high degree of numerical calculation accuracy can be obtained by considering up to the fourth derived function or higher. Therefore, in this section equations will be derived, considering up to the fourth derived function. If the angular displacement and velocity at the time are expanded into the Taylor series, considering up to the fourth derived function, they can be expressed by equation （A-11） instead of equations （A-9） and （A-10）．
（A-11）
Equation （A-11） can be rearranged
（A-12）
Where and , which constitute the transition matrix, are partial matrices of . As an example, a general term of will be shown below.
（A-13）
As is analogized from equation （A-13）， the partial matrices , , and are by the step size which satisfies the condition for not diverging during the process of repeated numerical calculations, that is the stabilizing condition, and also satisfies the condition for obtaining the correct solution and by the various factors constituting the equivalent torsional vibration system in Figure A-1.
The excitation torque and the time differentials and in equation （A-12） can be determined by adopting the measured indicator diagram and the specifications of reciprocating mass parts, and the angular displacement and velocity at the time can be numerical calculated by giving initial values of and . This numerical calculation can be repeated continuously. Since no steady waveform was obtained for the first two or three periods, the subsequent periods were adopted as the final results.
Since the values of the exciting torques in equation （A-12） are equals to zero in the case of free torsional vibrations, only the first term need to be considered.
A-3. Numerical simulation flow chart
All the numerical simulation values required in making numerical simulations now have been determined as discussed above. The numerical simulation chart is shown in Figure A-2. With the initial value given, repeated numerical calculations are made. The numerical calculation capacity is approximately 100 ［kB］ in the case where the numerical calculation is possible up to ten masses.
Figure A-1 Numerical simulation model of equivalent torsional vibration of multi-degree of freedom for high-speed diesel engine crankshaft system with rubber damper pulley ［General torsional vibration numerical simulation model］
Figure A-2 Forced torsional vibration numerical simulation flow chart by the transition matrix method
