仮想ばね法により基本境界条件を緩和させたB-spline Ritz法に設定する最適な仮想ばね定数の値について（その２）自重の作用による均質な等断面長柱の線形座屈問題

(1)

仮想ばね法により基本境界条件を緩和させたB-spline Ritz法

に設定する最適な仮想ばね定数の値について（その２）

－自重の作用による均質な等断面長柱の線形座屈問題－

山本寧音

1

_{・名木野晴暢}

2†

_{・足立忠晴}

3

1_{長岡技術科学大学大学院，}2_{都市・環境工学科，}3_{豊橋技術科学大学} †_{連絡著者 / Corresponding author (E-mail: [email protected])}

本稿では，線形化された梁－柱理論に従う自重の作用による均質かつ真っすぐな長柱の線形座屈問題を仮想ばね法により基本境界条件を緩和させたB-spline Ritz法により解析し，座屈荷重に与えるspline次数，解析領域内に配置するノットの数，仮想ばね定数の値および支持条件の影響を詳細に調べた．数値解に影響を与えない範囲の無次元仮想ばね定数の最大値は，前報と同様の鉛直ばねおよび回転ばね共に1010_であった．この値を用いて，spline次数を5次に設定したときの座屈荷重は，支持条件にかかわらずノットの数の増大にともなって有効桁数八桁から九桁で収束した．また，片持長柱の座屈荷重は有効桁数八桁から九桁まで厳密解に一致し，これ以外の支持条件については先行研究よりも精度の高い上界と思われる値を示した．

キーワード : Ritz法，正規化されたB-spline基底関数，仮想ばね法，

線形座屈問題，軸方向分布力，線形化された梁－柱理論

１．緒言

第一報1)_{では，線形化された梁－柱理論}2)_{に従う中心軸圧} 縮力を受ける均質かつ真っすぐな長柱（以下，長柱と略す）の線形座屈問題を仮想ばね法により基本境界条件を緩和させたB-spline Ritz法（以下，本手法と略す）により解析し，座屈荷重の精度，離散化条件，仮想ばね定数の値および支持条件の関係を詳細に調べた．調査の結果，支持条件と spline次数にかかわらず，ノットの増大にともなって数値的に安定した近似解を与える仮想ばね定数の最大値は，鉛直ばねと回転ばね共に1010_{であった．また，関数空間に5次} の正規化されたB-spline基底関数（以下，B-spline基底関数と略す）を張ったときの本手法の座屈荷重はノットの数の増大にともなって上から厳密解に近づき，収束値は有効桁数八桁から九桁まで厳密解に一致するという高い精度を示した．この問題の強形式は，定数係数の四階線型の斉次常微分方程式の境界値問題で表されるので，閉じた形の厳密解を求めることができ，その数量化は比較的容易である．これに対して，構造物の合理的な設計や材料の不均質性の活用の見地から用いられる変断面長柱3)_{や軸方向傾斜機}

能材料(Axially Functionally Graded Materials: AFGMs)4)_の線

形座屈問題は，線形化された梁－柱理論2)_{に基づくとすれ} ばその強形式は変数係数の四階線型の斉次常微分方程式の境界値問題で表される．これは，柱の断面諸量や材料特性が部材軸方向に変化することによる．文献5)には，幾つかの変断面長柱と軸方向に不均質な長柱の線形座屈問題の厳密解が纏められている．これは第一種および第二種 Bessel関数により表現されるため，その数量化は容易でなく，ある程度の専門知識を必要とする．また，構造設計や材料設計の観点に立てば，柱の軸方向に任意に変化する断面諸量や材料特性の取り扱いが求められる．こうなると，厳密解を求めることはほぼ不可能になり，何らかの手段により強形式を近似的に解かざるを得ない．近似解には，対象としている工学的諸問題を解決できるだけの精度が要求されることは言うまでもない．他方で，渡辺6)_が述べているように，何らかの手法により数値的に安定して高精度な近似解が得られるようになれば，それは厳密解や級数解に代わって有限要素の性能評価や各種近似解析手法が強形式を正しく解いているかの検証へ活用できるであろう．このような背景の下，本稿では変数係数の四階線型の斉次常微分方程式の斉次境界値問題で記述される“自重による均質かつ真っすぐな長柱の線形座屈問題”を対象として，第一報1)_{と同様の調査を実施し，ある程度の有益な情報を} 得ることができたので報告する．

(2)

―39―

２．数学モデル（強形式）

(1) 座標系と仮定図-1 直線直交座標系O-xyzと長柱図-1に示すような右手直交直線座標系O-xyzにある長柱の単位長さ当たりの軸方向分布力p (x)の作用による分岐座 屈問題を考える．ここで，x軸，y軸，およびz軸の設定は前 報1)_{と同様であり，}_u_x_{(x, y, z), u}_y_{(x, y, z), u}_z_{(x, y, z)は，それ} ぞれ，長柱のx, y, z方向の変位である．柱の長さLはyz断面 の寸法に比べて卓越しているとし，yz断面の断面積をA，中 立軸に関する断面二次モーメントをIで表す（正しくはIyである）．長柱はx軸方向に一様断面とし，かつ同方向に均質 な材料からなるとする．本稿では，軸方向分布力p (x)が作用することにより長柱 が僅かにy軸まわりに曲げ変形を生じてつり合った状態を 考える．この不安定現象の数学モデルには，線形化された梁－柱理論2)_{を用いる．この理論に設けられる仮定，およ} びその仮定に基づく変位場とひずみ場についても前報1)_を参照されたい． (2) 一般化応力の平衡方程式と境界条件式単位長さ当たりの軸方向分布力p (x)を受ける長柱が僅か に曲げ変形を生じてつり合ったとき，長柱のある点の近傍における局所的な力のつり合いを考える．領域 = (0, L)にある長柱の任意の点の近傍の微小な線 分要素におけるx軸方向とz軸方向の力のつり合い条件式， およびxz平面内の任意の点に関する力のモーメントのつ り合い条件式は，それぞれ，次のように表される． d ( ) _{( )} ₀ d N x _{p x} x + = , (1a) d d ( ) ( ) ( ) 0 d d w x S x N x x x ì ü ï ï ï ₊ ï₌ í ý ï ï ï ï î þ , _(1b) d ( ) ( ) 0 d M x S x x - = in W . (1c) ここで，N (x)は僅かに曲がった長柱の断面の法線方向に生 じる軸力，S (x)は同断面の接線方向に生じるせん断力であ り，M (x)はy軸まわりの曲げモーメントである（正しくは Myである）．これらは，yz断面に分布する垂直応力とせん断 応力を軸線に関する一般化応力（断面力）として表したものである．式(1b)の下線部は，軸力N (x)が影響する非線形項であることを意味している．この項を無視すると，式(1b) は初等梁理論のz軸方向の力のつり合い条件式に一致する． 長柱の左端点(x = 0)と右端点(x = L)での境界条件式は，次のような組み合わせにより表される． u= on u u or n N_x =N on , (2a) w= on w w or n Q Q_x = on Q, (2b) q=q on  or n M_x =M on M. (2c) ここで， = {0, L}は一次元領域の境界を表しており，下付 き文字のu, w, は変位の境界条件を，N, Q, Mは一般化応力 の境界条件を意味している．nxは線部材の境界に立てた外向き法線の方向余弦のx方向の成分であり，左端点(x = 0)の 境界でnx = – 1, 右端点(x = L)の境界でnx = + 1を取る．は曲げたわみ曲線wの接線と軸線とがなす角を表すたわみ角 であり，これは微小角(| | << 1)であるとすれば，次のように表すことができる． d ( ) ( ) d w x x x q = . (3) また，Qは z軸に平行なせん断力であり，次のように表さ れる． d ( ) d ( ) ( ) ( ) d d M x w x Q x N x x x = + . (4) これは，式(1b)の{ }の中の式に相当する．僅かに曲がった長柱の断面の接線方向に生じるせん断力S (x)とは異なるこ とに注意する．なお，各条件式の右辺のバーを付した変位成分u, w，たわみ角および一般化応力N, Q, M は，与える 既知の値であることを意味している．長柱の支持条件の数学モデルは，これらの境界条件式を組み合わせて，次のように理想化される．ピン支点（ヒンジ支点）：

0 w =

on w, M = on 0 M. (5) 固定端：

0 w =

on w, q = on 0 . (6) 自由端： 0 M = on M, Q = on 0 Q. (7) これらは，斉次の境界条件式である． (3) 一般化応力と変位で表された構成方程式長柱に対して，本稿のように右手直交直線座標系O-xyz を設定すると，yz断面に設定した図心軸まわりの断面一次 モーメントは零になる．このことから，一般化応力と変位 p (x) O L y, uy z, uz x, ux

(3)

で表された構成方程式に含まれる軸方向変位u (x)と面外変 位w (x)に関する項は非連成化され，微小ひずみ・微小回転 の有限変位に基づく梁－柱理論の構成方程式は，次のように表される． 2 d ( ) 1 d ( ) ( ) d 2 d u x w x N x EA x x ì ü ï _é _ù ï ï ï ï ê ú ï = í_ï + _ê _ú ý_ï ï ë û ï ï ï î þ , _(8a) 3 3 d ( ) d ( ) ( ) ( ) d d w x w x Q x EI N x x x = - + , (8b) 2 2 d ( ) ( ) d w x M x EI x = - . (8c) ここで，EAは伸び剛性，EIは曲げ剛性であり，Eは縦弾性 係数である．なお，式(8a)の下線を付した右辺の第二項を無視したものが，線形化された梁－柱理論の構成則になる． (4) 支配方程式と境界条件による強形式式(1c)を式(1b)に代入し，更に一般化応力と変位で表された構成方程式の式(8c)を代入すれば，変位で表された平衡方程式である長柱の線形座屈問題の支配方程式（基礎方程式）が次のように得られる． 2 2 2 2 d d ( ) d d ( ) ( ) 0 d d d d w x w x EI N x x x x x é ù é ù ê_- ú₊ ê ú₌ ê ú ê ú ê ú ë û ë û in W . (9) さて，x軸方向のつり合いの微分方程式(1a)は，線部材の単 純引張理論の軸力N (x)に関する平衡方程式を表しており， 曲げとは非連成である．よって，境界条件式(2a)が与えられれば，N (x)の厳密解は容易に求められる．長柱の自重を 単位長さ当たりに作用する軸方向の分布力p (x)で表すと， p (x) = p0である．ただし，p0 = gAは分布力の荷重強度， は密度であり，gは重力加速度である．これを式(1a)に代入 し，長柱の右端点(x = L)で自由端の条件(N = 0 on N)の下で解けば，自重の作用により生じる長柱の軸力N (x)が次のよ うに求められる． 0 ( ) ( ) N x = -p L x- ( 0£ £ ). x L (10) 式(10)を式(9)に代入すれば，自重の作用を受ける均質かつ等断面長柱の線形座屈問題の支配方程式は，次のように表される． 4 0 4 d ( ) d d ( ) ( ) 0 d d d w x w x EI p L x x x x ì ü ï ï ï ï - - í_ï - ý_ï= ï ï î þ in W . (11) 境界条件式を面外変位wで表すと， ピン支点（ヒンジ支点）：

0 w =

on w, 2 2 d 0 d w EI x - = on M. (12) 固定端：

0 w =

on w, d 0 d w x = on . (13) 自由端： 2 2 d ₀ d w EI x - = on M , 3 0 3 d d ( ) 0 d d w w EI p L x x x - - - = on Q. (14) 本問題は，変数係数の四階線型の斉次常微分方程式の斉次境界値問題の自明でない解を求める固有値問題である．求められる固有値が座屈荷重qcrであり，これに対応する固有関数wが座屈したときの曲げ変形を表す．

３．強形式に等価な汎関数の停留値問題

仮想ばね法により基本境界条件を緩和させた変分原理の汎関数（全ポテンシャルエネルギ）_{は，次のように表} される． 2 2 0 0 ₀ 2 2 0 0 ₀ ( ) ( ) 2 2 ( ) ( ) 2 2 w wL L x x L L L x x L U V w w w w q q a a b b q q q q = = = = P = + + - + -+ - + -. (15) ここで，Uは長柱に蓄積されるひずみエネルギ，Vは軸方向 の等分布力（自重）p (x) = p0による外力ポテンシャルであり，それぞれ，式(16)と式(17)によって与えられる． 2 2 2 0 d _d 2 d L EI w U x x æ _ö÷ ç _÷ ç = ç_ç ÷_÷_÷ è ø ó ô ô õ . (16) 2 0 0 d ( ) d 2 d L p w V L x x x æ _ö÷ ç _÷ ç = - - _ç _÷_÷÷ çè ø ó ô ô õ . (17) また，w0, wLは，それぞれ，長柱の左端点(x = 0)と右端点 (x = L)の面外変位wに抵抗する仮想の鉛直ばね定数，0, L は，それぞれ，長柱の左端点(x = 0)と右端点(x = L)のたわみ角に抵抗する仮想の回転ばね定数，式中の右辺のバーを付した面外変位wとたわみ角は境界での既知の値であり，その下付き文字の0とLは，それぞれ，x = 0, x = Lでの条件であることを意味する．式(12)と式(13)の基本境界条件は斉次であるので，本稿で考えるべき汎関数は，次のように表される． 2 2 2 0 2 0 0 2 2 0 0 2 2 0 0 d d d ( ) d 2 d 2 d 2 2 2 2 L L w wL x x L L x x L p EI w w x L x x x x w w q q a a b b q q = = = = æ ö_÷ æ ö_÷ ç ÷ ç ÷ ç ç P = ççè ÷÷÷ø - - ççè ÷_÷÷ø + + + + ó ó ô ô ô ô_õ õ . (18) この汎関数が停留するような微分可能，かつ連続な関数 wを見出せば，それは強形式の弱解を得ることに等しい．

(4)

―41―

４．仮想ばね法により基本境界条件を緩和させた

B-spline Ritz法による離散化された座屈方程

式と座屈条件式の定式化

式(18)の汎関数の要請をすべて満足する微分可能，かつ連続な変関数（長柱の面外変位）w (x)を見出すことは困 難を極める．よって，変関数wを次のような有限自由度の 試行関数で近似し，汎関数が停留するような近似解を探求する． CP , 1 , ( ) ( ) n i i p i i p w x AN x = @ =

å

N A , (19) ここで， CP CP , [ 1,( ), 2,( ), , 1, ( ), ,( )] i p = N px N px Nn - p x Nn p x N  , ₍₂₀₎ CP CP T 1 2 1 { , , ,A A A_n _-,A_n } = A  , ₍₂₁₎ であり，Tは転置を意味する．Ni, p (x)は次数pのB-spline基底 関数，iは基底関数の番号，Aiはi番目の基底の重み係数であ り，nCPは基底関数の数である．領域に1の分割(partition of unity)を満たすB-spline基底を張るために必要なノットの数は，(m + 2p)個である．ただし，mは領域内（境界を含む）に設定するノットの数であり，2pは領域の外側（境 界を含まない）に設ける付加ノットの数である．基底関数の数nCPと境界を含む領域内に設定するノットの数m との間には，次の関係が成立する． CP 1 n =m+ - . p (22) 式(19)を式(18)の汎関数に代入すると，は有限の自由度で近似される．この近似された汎関数の停留条件： CP CP 0 ( 1,2, , 1, ) i i n n A ¶P = = -¶  , (23) より，次の離散化された座屈方程式を得る． 0 P 0 G [( + )-p ] = K K K K A 0  . (24) ここで， 2 T 2 , , 0 2 2 0 d d d d d L i p j p EI x x x = óô õ N N K , (25) は自由境界を有する長柱の剛性行列， T T P 0 , , ₀ , , T T , , , , 0 0 d d d d d d d d w i p j p_x wL i p j p_{x L} i p j p i p j p L x x L x x x x q q a a b b = = = = = + + + K N N N N N N N N , ₍₂₆₎ は基本境界条件を緩和するための仮想ばねの剛性行列であり，K = K0 + KPは長柱の剛性行列である．また， T , , G 0 0 d d ( ) d d d L i p j p p L x x x x = óô -õ N N K , (27) (a) P-P (b) C-F (c) C-C (d) C-P 図-2 数値実験に用いる長柱の支持条件 O z, w x, u L z, w x, u L O z, w x, u L O z, w x, u L O

(5)

は長柱の幾何剛性行列であり，0は零の列ベクトルである． 剛性行列K0と幾何剛性行列KGは，(p + 1) 点の Gauss-Legendreの数値積分により求めた．式(24)の自明でない解は，次式の座屈条件式： 0 P 0 G det[(K +K )-pK ]=0 , (28) から固有値p0を求め，それを式(24)に代入して固有値p0に対応する固有ベクトルAを決定することで得られる．

５．数値実験および考察

ここでは，軸方向等分布力p (x) (= p0)が作用する均質かつ等断面長柱の線形座屈問題を本手法により解析し，座屈荷重pcrに与えるspline次数p，境界を含む領域内に設定するノットの数（以下，ノットの数と略す）m，仮想ばね定数 の値および長柱の支持条件の影響を詳細に調べる．本稿では長柱の座屈荷重pcr，仮想の鉛直ばね定数w0, wL および仮想の回転ばね定数0, Lを，次のように無次元化して表している． 3 * p Lcr p EI = . (29) 3 3 * w0L wLL EI EI a a a = = , (30) * 0L LL EI EI q q b b b = = . (31) また，前報の知見1)_{を参考にして，数値実験に用いる無次} 元仮想ばね定数*_{を次のように定義した．} * * * k =a =b . (32) 数値実験に設定する支持条件は，図-2に示す四ケースとし，例えばC-Fのように表す．これは柱の下端(x = 0)で固定 Cであり，上端(x = L)で自由Fであることを意味する．また， ピン支点（ヒンジ支点）はPで表す．なお，(a)と(c)は対称の支持条件，(b)と(d)は非対称の境界条件であり，(a)と(b) は静定，(c)と(d)は不静定である．通常，工学分野の諸問題の解決を目的とするならば，強形式の近似解は有効数字三桁から四桁程度の精度があれば十分であるが，本稿は前報と同様の観点1)_{の下，無次元} 座屈荷重p*_{を有効桁数十桁程度で整理した．なお，数値計} 算には64 bit OSのpersonal computerを用い，Fortranの倍精度で計算した． (1) Timoshenkoのエネルギ法に基づく座屈荷重の評価 Timoshenokoのエネルギ法によれば，軸方向等分布力p (x) (= p0)の作用による均質かつ等断面長柱の座屈荷重pcrは，次式によって求めることができる7)_． 2 2 2 0 cr ₂ 0 0 1 d d 2 d d ( ) d 2 d L L w EI x x p p _{L x} w _x x æ _ö÷ ç _÷ ç _÷ ç _÷÷ çè ø = æ _ö÷ ç _÷ - ç_ç _÷ ÷ çè ø ó ô ô õ ó ô ô õ . (33) ここで，分子は座屈による曲げ変形によって長柱に貯えられるひずみエネルギ（内力仕事）であり，分母は長柱の曲げ座屈変形により生じる軸圧縮変位と軸方向等分布力とがなす外力仕事を意味する．式(33)の右辺の表示は，線形化された梁－柱理論に基づく長柱の線形座屈問題におけるRayleigh商である7)_{．これより，座屈荷重}_p crは，外力仕事に対するひずみエネルギの大きさとして評価されると解釈することができよう．なお，許容である試行関数を仮定したRitz法による近似解は，厳密解の上界を与える8), 9)_ので，次の関係が成り立つ． Exact cr cr p ³p . (34) ここで，不等式の左辺はRitz法による近似解であり，右辺は厳密解を意味する．ただし，前報1)_{によれば，関数空間に} 高次のB-spline基底を張ったとき，本手法の座屈荷重はノットの増大にともなって有効桁数九桁の収束値を得るが，数値計算上の誤差と思われる原因により，その値は厳密解の下界を与えることがある．さて，前報1)_{では，四つの支持条件について全て閉じた} 形式の厳密解を求めることが可能なため，本手法の近似解の上界と下界を明確にすることができた．本稿が対象としている問題についても厳密解は示されているが，前報のように閉じた形式で与えられるのはC-F（片持）の支持条件に限定されているようであり5), 10)_{，その他の支持条件につ} いては境界条件を満足するべき級数解の形式（以下，べき級数解と略す）で与えられている11) – 13)_{．そのため，上界と} 下界の明確な区別が難しい．そこで，本稿では，次の手順により，本手法の近似解に影響を与えない範囲の最大の無次元仮想ばね定数*_{の値と，これを設定したときの本手法} の近似解の精度を調べることにした． [1] 厳密解が得られるC-F（片持）について前報と同様の調査を行い，結果を整理する． [2] 残り三ケースの支持条件については，前報と[1]の本手法の近似解の収束性と計算精度を判断材料とし， Huang・Dareing11)_{のべき級数解と杉山ら}12)_{のべき級数} 解の収束性と有効桁数六桁の収束値と比較することで，本手法の近似解の収束性と精度を考察する．なお，Huang・Dareing11)_と杉山ら12)_{は，自重を受ける長柱の} 支配方程式の面外変位を式(35)のようなべき級数で近似し，支配方程式と境界条件から未定係数bnに関する連立方程式を導出している．これを数値的に解くことで，座屈荷重と，これに対応する座屈モードの数値を得ている．

(6)

―43―

表-1 自重の作用による長柱の無次元座屈荷重p_{に与えるspline次数p, ノットの数mと仮想ばね定数}__{の影響: C-F}

p _ _m _{Mode 1} _{Mode 2} _{Mode 3} _{Mode 4}

4 108_{11 7.83734735} _55.9771972 _148.517011 _285.588375 21 7.83734721 55.9770304 148.508401 285.451764 31 7.83734720 55.9770282 148.508302 285.450336 41 7.83734720 55.9770280 148.508295 285.450236 51 7.83734720 55.9770280 148.508294 285.450221 61 7.83734720 55.9770280 148.508294 285.450217 1010_{11 7.83734756} _55.9771987 _148.517015 _285.588383 21 7.83734744 55.9770320 148.508406 285.451772 31 7.83734744 55.9770299 148.508307 285.450344 41 7.83734744 55.9770297 148.508299 285.450244 51 7.83734744 55.9770297 148.508298 285.450229 61 7.83734743 55.9770297 148.508298 285.450226 1012_{11 7.83734834} _55.9772044 _148.517031 _285.588418 21 7.83734706 55.9770293 148.508398 285.451758 31 7.83734747 55.9770301 148.508307 285.450346 41 7.83734723 55.9770282 148.508296 285.450237 51 7.83734749 55.9770301 148.508299 285.450231 61 7.83734741 55.9770295 148.508298 285.450225 5 108_{11 7.83734720} _55.9770297 _148.508592 _285.463595 21 7.83734720 55.9770280 148.508294 285.450238 31 7.83734720 55.9770280 148.508294 285.450216 41 7.83734720 55.9770280 148.508294 285.450216 51 7.83734720 55.9770280 148.508294 285.450215 61 7.83734721 55.9770280 148.508294 285.450215 1010_{11 7.83734747} _55.9770317 _148.508597 _285.463604 21 7.83734744 55.9770297 148.508299 285.450247 31 7.83734743 55.9770296 148.508298 285.450225 41 7.83734744 55.9770297 148.508298 285.450224 51 7.83734743 55.9770297 148.508298 285.450224 61 7.83734744 55.9770297 148.508298 285.450224 1012_{11 7.83734943} _55.9770458 _148.508634 _285.463663 21 7.83734782 55.9770325 148.508307 285.450264 31 7.83734737 55.9770292 148.508297 285.450222 41 7.83734742 55.9770296 148.508298 285.450224 51 7.83734731 55.9770288 148.508296 285.450219 61 7.83734745 55.9770297 148.508298 285.450224 Exact5) _{7.83734744 55.9770297 148.508298 285.450224}

Huang and Dareing (1969)11) _{7.8373 55.977 148.51}

Sugiyama et al. (1977)12) _{7.83735 55.9770}

The structural mechanics handbook (1986)15) _7.8373

Wang and Ang (1988)14) _7.84

Duan and Wang (2008)13) _7.8373

Note: 杉山ら12)_の_{Mode 1のp}*_{は下から収束し，}_{N = 24としたときに有効桁数六桁の収束値を得ている．Mode 2は上か} ら収束し，N = 36としたときに有効桁数六桁の収束値を得ている．

(7)

0 ( ) N n n n wx bx = =

å

. (35) ここで，は無次元座標の変数であるが，その定義は文献 11) と文献 12) で異なることに注意されたい． Huang ・ Dareing11)_{が示した数値は，近似項数}_{Nに対する座屈荷重の} 収束性が不明であるが，一次から三次までの座屈荷重を有効桁数五桁で与えてくれている．杉山ら12)_{は，近似項数}_N に対する一次と二次の座屈荷重の変化の情報を有効桁数六桁で提供してくれている．この他にもDuan・Wang13)_の変数rn/3_{(n = 0, 1, 2, 3)と一般化された超幾何関数F (a}₁_{, a}₂_{; b}₁_, b2, b3; r)との積を基底とした境界条件を満足する厳密解（厳密解は四つの項で表され，一般化された超幾何関数の定義には無限級数が含まれる）も提案されているが，文献11)と文献12)よりも数値解に関する情報が少ない． (2) 片持(C-F)長柱表-1は，支持条件がC-Fである自重の作用による長柱の一次から四次までの無次元座屈荷重p_に与える_{spline次数p,} ノット数mと無次元仮想ばね定数*_{の影響を表形式で整理} したものである．分岐座屈荷重に相当するのは最低次の座屈荷重であるが，二次以上の座屈荷重とこれに対応する座屈モードは，座屈荷重を高めるための重要な情報を提供してくれる1)_{ことから，本稿でも一次から四次までの座屈荷} 重について検討した．無次元仮想ばね定数*_は，108_{, 10}9_, 1010_{, 10}11_{, 10}12_{の五ケースについて検討したが，本稿では}₁₀8_, 1010_{, 10}12_{の結果を提供する．B-spline基底関数の次数pは，} 前報と同様の理由により，p = 4, 5に設定した1)_．_Spline次数 p = 4のとき，面外変位wは三階導関数まで連続になるので， せん断力Q (x)は一次の区分的多項式で近似される．当然の ことではあるが，spline次数p = 5に設定すると，面外変位w は四階導関数まで連続であるので，せん断力Q (x)の近似式 は二次の区分的多項式になる．ノットの数mは，11から61 まで変化させ，有効桁数九桁で表した無次元座屈荷重p_の値の変化を調べた．なお，表のboldの数値は近似解が収束し，かつ厳密解と一致したことを意味する．また，近似解は収束したが，その収束値が厳密解と一致しないものには， boldの数値に下線を付した．先にも述べたとおり，ここでは厳密解と比較することで，近似解に影響を与えない無次元仮想ばね定数*_{の最大の値の候補と，その精度を検証す} る．表には厳密解の数値が記載されていれば十分であるが，参考までにHuang・Dareing11)_と杉山ら12)_{のべき級数解，}

Duan・Wang13)_{の厳密解，}_Wang・Ang14)_の_{Ritz法による近似}

解，および構造力学公式集15)_{の数値も併記した．} 厳密解の数値は，文献5)と文献10)を参考にして数値的に求めた．ここで，自重を受ける片持長柱の座屈方程式は，次式で与えられる5), 10)_． 1/3( ) 0 J_- z = . (36) ただし， 3 cr 2 3 q L z EI = , (37) である．ここで，J– 1/3 (z)は(– 1/3)次の第一種Bessel関数である．この方程式を解析的に解くことは困難であるため，これを有効桁数十五桁の精度で数値的に解くこととし，次の四つの零根を得た． 1 2 3 4 1.86635085887389, 4.98785323143515, 8.12426538193969, 11.2635148254276. z z z z = = = = (38) ここで，zの下付添字は，座屈モードの次数を意味する．式 (38)を式(37)に代入すれば，一次から四次までの無次元座屈荷重p*_{の厳密解の数値を求めることができる．} 表-1より，spline次数pにかかわらず，近似解に影響を与えない最も大きな無次元仮想ばね定数*_{の値は，前報}1)_と同様に1010_である．_*_{= 10}8_{の近似解は柔らかく，}_*_{= 10}12 のそれは数値的に安定していない．Spline次数p = 5に設定 すると，一次の無次元座屈荷重p_{は有効桁数八桁で，二次} から四次までのそれは有効桁数九桁での収束値を得ている．ただし，二次の無次元座屈荷重p_{は，ノットの数}_mを21 から41に変化させる過程で数値が振動しており，かつm = 31のときに下界である．上界定理から判断すると，九桁目の数値の信頼性は低いと考えられる．なお，一次と二次の無次元座屈荷重p_{は厳密解と有効桁数八桁まで一致して} おり，三次と四次のそれは有効桁数九桁まで一致する．前報1)_{と比較すると，一次と二次の座屈荷重の精度が若干低} 下している．これは，長柱が自重によって座屈するときの曲げ変形の形が前報よりも複雑になっているためであると考えられる． (3) 支持条件がP-P, C-CおよびC-Pである長柱表-2から表-4は，それぞれ，支持条件をP-P, C-C, C-Pとして表-1と同様の検討を行ったものである．ただし，本節では厳密解を有効桁数九桁程度で数値化したものがないため，上界と下界を明確に区別することができず，また，近似解の収束値が何に収束したのかを正確に判断することが難しい．よって，ここでは近似解が有効桁数九桁で収束した数値をboldで表すこととし，Huang・Dareing11)_の一次から三次までの座屈荷重のべき級数解，および杉山ら12)_のべき級数解の収束状態とその収束値と比較することで，本手法の近似解の精度を判断する．表-2より，spline次数pにかかわらず，近似解に影響を与えない範囲の無次元仮想ばね定数*_{の最大値は}₁₀10_である．また，spline次数p = 5かつ無次元仮想ばね定数*_{= 10}10_に設定したときの本手法の近似解は，ノットの数mの増大にと もない数値的に安定した収束状態を示しており，有効桁数

(8)

―45―

表-2 自重の作用による長柱の無次元座屈荷重p_{に与えるspline次数p, ノットの数mと仮想ばね定数}__{の影響: P-P}

p _ _m _{Mode 1} _{Mode 2} _{Mode 3} _{Mode 4}

4 108_{11 18.5687256} _86.4315979 _196.318759 _353.186380 21 18.5687239 86.4308408 196.291038 352.751647 31 18.5687239 86.4308316 196.290779 352.748588 41 18.5687239 86.4308309 196.290760 352.748378 51 18.5687239 86.4308308 196.290757 352.748346 61 18.5687239 86.4308308 196.290757 352.748339 1010_{11 18.5687266} _86.4316031 _196.318775 _353.186414 21 18.5687249 86.4308460 196.291055 352.751680 31 18.5687248 86.4308368 196.290796 352.748622 41 18.5687248 86.4308361 196.290777 352.748412 51 18.5687248 86.4308360 196.290774 352.748380 61 18.5687248 86.4308359 196.290773 352.748372 1012_{11 18.5687259} _86.4316024 _196.318774 _353.186413 21 18.5687247 86.4308460 196.291055 352.751680 31 18.5687248 86.4308368 196.290796 352.748622 41 18.5687248 86.4308361 196.290777 352.748412 51 18.5687248 86.4308360 196.290774 352.748380 61 18.5687248 86.4308360 196.290773 352.748373 5 108_{11 18.5687239} _86.4308564 _196.292612 _352.791084 21 18.5687239 86.4308308 196.290761 352.748427 31 18.5687239 86.4308308 196.290756 352.748338 41 18.5687239 86.4308308 196.290756 352.748335 51 18.5687239 86.4308307 196.290756 352.748335 61 18.5687239 86.4308308 196.290756 352.748335 1010_{11 18.5687248} _86.4308615 _196.292629 _352.791118 21 18.5687248 86.4308360 196.290777 352.748460 31 18.5687248 86.4308359 196.290773 352.748371 41 18.5687248 86.4308359 196.290773 352.748369 51 18.5687248 86.4308359 196.290773 352.748369 61 18.5687248 86.4308359 196.290773 352.748369 1012_{11 18.5687245} _86.4308614 _196.292629 _352.791118 21 18.5687250 86.4308362 196.290777 352.748461 31 18.5687248 86.4308360 196.290773 352.748372 41 18.5687248 86.4308360 196.290773 352.748369 51 18.5687249 86.4308360 196.290773 352.748369 61 18.5687248 86.4308360 196.290773 352.748369 Huang and Dareing (1969)11) _{18.569 86.431 196.29}

Sugiyama et al. (1977)12) _{18.5687 86.4309}

Wang and Ang (1988)14) _18.58

Duan and Wang (2008)13) _18.5687

Note: 杉山ら12)_{のMode 1のp}*_{は上から収束し，}_{N = 30としたときに有効桁数六桁の収束値を得ている．Mode 2は下か} ら収束し，N = 42としたときに有効桁数六桁の収束値を得ている．

(9)

表-3 自重の作用による長柱の無次元座屈荷重p_{に与えるspline次数p, ノットの数mと仮想ばね定数}__{の影響: C-C}

4 108_{11 74.6290737 157.044090 325.758330 495.467677} 21 74.6285668 157.032875 325.516057 494.004251 31 74.6285605 157.032742 325.513599 493.989459 41 74.6285601 157.032732 325.513431 493.988483 51 74.6285600 157.032731 325.513405 493.988339 61 74.6285600 157.032731 325.513400 493.988306 1010_{11 74.6290822 157.044139 325.758384 495.467860} 21 74.6285755 157.032924 325.516112 494.004432 31 74.6285692 157.032792 325.513654 493.989639 41 74.6285687 157.032782 325.513486 493.988663 51 74.6285687 157.032780 325.513460 493.988519 61 74.6285686 157.032780 325.513454 493.988487 1012_{11 74.6290932 157.044167 325.758436 495.468037} 21 74.6285710 157.032915 325.516094 494.004406 31 74.6285705 157.032794 325.513659 493.989647 41 74.6285668 157.032778 325.513478 493.988652 51 74.6285694 157.032782 325.513463 493.988525 61 74.6285684 157.032780 325.513454 493.988486 5 108_{11 74.6285675 157.033185 325.542366 494.427873} 21 74.6285600 157.032731 325.513442 493.988741 31 74.6285600 157.032730 325.513398 493.988303 41 74.6285600 157.032730 325.513397 493.988291 51 74.6285600 157.032730 325.513397 493.988290 61 74.6285600 157.032730 325.513397 493.988290 1010_{11 74.6285766 157.033235 325.542422 494.428056} 21 74.6285687 157.032781 325.513498 493.988923 31 74.6285686 157.032780 325.513453 493.988483 41 74.6285686 157.032780 325.513452 493.988472 51 74.6285686 157.032780 325.513452 493.988471 61 74.6285686 157.032780 325.513452 493.988471 1012_{11 74.6286030 157.033287 325.542508 494.428094} 21 74.6285738 157.032792 325.513522 493.988964 31 74.6285671 157.032777 325.513447 493.988475 41 74.6285683 157.032779 325.513451 493.988471 51 74.6285673 157.032777 325.513446 493.988463 61 74.6285688 157.032780 325.513453 493.988473 Huang and Dareing (1969)11) _{74.629 157.03 325.51}

Sugiyama et al. (1977)12) _{74.6286 157.031}

Wang and Ang (1988)14) _78.96

Duan and Wang (2008)13) _74.6286

(10)

―47―

表-4 自重の作用による長柱の無次元座屈荷重p_{に与えるspline次数p, ノットの数mと仮想ばね定数}__{の影響: C-P}

4 108_{11 52.5007890} _129.059252 _276.364762 _439.725890 21 52.5006621 129.054333 276.245006 438.864630 31 52.5006604 129.054274 276.243769 438.855769 41 52.5006603 129.054270 276.243682 438.855176 51 52.5006603 129.054269 276.243669 438.855088 61 52.5006603 129.054269 276.243666 438.855068 1010_{11 52.5007917} _129.059270 _276.364784 _439.725970 21 52.5006648 129.054352 276.245029 438.864709 31 52.5006632 129.054293 276.243792 438.855848 41 52.5006631 129.054288 276.243705 438.855255 51 52.5006631 129.054288 276.243692 438.855167 61 52.5006630 129.054287 276.243689 438.855147 1012_{11 52.5007961} _129.059284 _276.364814 _439.726070 21 52.5006626 129.054346 276.245016 438.864690 31 52.5006633 129.054293 276.243793 438.855851 41 52.5006618 129.054285 276.243698 438.855245 51 52.5006634 129.054289 276.243694 438.855171 61 52.5006630 129.054287 276.243689 438.855147 5 108_{11 52.5006615} _129.054416 _276.254325 _439.061301 21 52.5006603 129.054269 276.243684 438.855290 31 52.5006603 129.054269 276.243665 438.855065 41 52.5006603 129.054269 276.243665 438.855059 51 52.5006603 129.054269 276.243665 438.855058 61 52.5006603 129.054269 276.243665 438.855058 1010_{11 52.5006645} _129.054435 _276.254349 _439.061382 21 52.5006631 129.054288 276.243707 438.855370 31 52.5006630 129.054287 276.243688 438.855144 41 52.5006630 129.054287 276.243688 438.855138 51 52.5006630 129.054287 276.243688 438.855138 61 52.5006630 129.054287 276.243688 438.855137 1012_{11 52.5006776} _129.054467 _276.254408 _439.061413 21 52.5006656 129.054295 276.243723 438.855399 31 52.5006627 129.054286 276.243686 438.855141 41 52.5006630 129.054287 276.243687 438.855138 51 52.5006622 129.054286 276.243683 438.855131 61 52.5006631 129.054288 276.243688 438.855139 Huang and Dareing (1969)11) _{52.501 129.05 276.24}

Sugiyama et al. (1977)12) _{52.5009 129.054}

Wang and Ang (1988)14) _53.91

Duan and Wang (2008)13) _52.5007

(11)

九桁での収束値を得る．また，その値は先行研究とよく一致している．前報1)_{と本稿の片持長柱の数値実験結果から} 判断すると，本手法の近似解は，先行研究よりも精度の高い上界の値を与えていると考えられるが，二次の座屈荷重についてはもう少し詳細な検討が必要である．表-3からは，次のことがわかる．  Spline次数pにかかわらず，数値的に安定し，近似解に影響を与えない最大の無次元仮想ばね定数*_の値は，やはり1010_である．  C-FおよびP-Pの支持条件に比べると，近似解の収束性が良くない．これは，前報1)_{と同様である．}  Spline次数p = 5かつ無次元仮想ばね定数*_{= 10}10_に設定したときの本手法の近似解は，ノットの数mの増大にと もない数値的に安定した収束状態を示しており，有効桁数九桁での収束値を得る．その値は先行研究とよく一致している．  P-Pの支持条件と同様，本手法の近似解は，先行研究よりも精度の高い上界の値を与えていると考えられるが，二次の座屈荷重については検討が必要である．表-4についてもこれまでとほぼ同様の結果が得られている．異なるのは，次の二点であった．  Spline次数p = 5かつ無次元仮想ばね定数*_{= 10}10_に設定したときの本手法の一次から三次の無次元座屈荷重p はノットの数mの増大にともなって有効桁数九桁で収 束するが，四次のそれは有効桁数八桁である．  本手法の近似解は，先行研究よりも精度の高い上界の値を与えていると考えられるが，一次の座屈荷重についてはもう少し詳細な検討が必要である．

６．結言

本稿では，線形化された梁－柱理論に従う自重の作用による均質かつ真っすぐな長柱の線形座屈問題を仮想ばね法により基本境界条件を緩和させたB-spline Ritz法により解析し，座屈荷重の精度，spline次数，解析領域内に配置するノットの数，仮想ばね定数の値および支持条件の関係を詳細に調べた．本調査から，次の知見が得られた． [1] 支持条件とspline次数pにかかわらず，近似解に影響を与えない範囲の仮想ばね定数の最大値は，無次元仮想鉛直ばね定数*_{，無次元仮想回転ばね定数}_*_共に 1010_である． [2] Spline次数p = 5かつ無次元仮想ばね定数*_{= 10}10_に設定したときの長柱の一次から四次までの無次元座屈荷重p*_{は，支持条件にかかわらず，ノットの数}_mの増 大にともなって上から正解に近づくような収束状態を示し，有効桁数八桁から九桁での収束値を得る． [3] [2]で得られた収束値は，有効桁数八桁から九桁の精度で，先行研究よりも精度の高い上界の値を与えていると考えられる．なお，固有関数とその導関数の精度についても検討しているが，これについては別の機会に報告する予定である．前報1)_{と本稿の調査から汎関数の変関数をB-spline基底} 関数により近似するとき，仮想ばね法により基本境界条件を緩和させることによっても有効桁数八桁までは厳密解に一致する，または一致していると思われる高い解析精度が得られることを示した．今後は，基本境界条件を正確に満足する許容な試行関数を仮定したB-spline Ritz法により同様の調査を実施し，弱形式の近似解にB-spline基底関数を用いたときの基礎的な情報を提供したいと考えている．なお，本稿は同研究分野の専門家による査読を経ずに公開されるが，本稿で示した表の数値は，読者の判断で種々の離散化手法の検証などに使っていただければ幸いである．謝辞：本研究の一部は，LIXIL住生活財団調査研究助成（2013年度）および前田記念工学振興財団研究助成（平成 27年度）を受けて行われました．ここに記して謝意を示します．参考文献 1) 名木野晴暢，山本寧音，足立忠晴：仮想ばね法により基本境界条件を緩和させたB-spline Ritz法に設定する最適な仮想ばね定数の値について（その１）－中心軸圧縮力を受ける均質な等断面長柱の線形座屈問題－, 大分工業高等専門学校紀要, 第56号, pp.24-37, 2019. 2) 土木学会構造工学委員会計算力学とその応用に関する研究小委員会：構造工学における計算力学の基礎と応用 (構造工学シリーズ 7), 土木学会, pp.118-154, 1997. 3) 崎山毅，栗林和夫：変断面梁柱の挫屈解法について，長崎大学工学部研究報告, 第8号, pp.25-31, 1977. 4) 山本寧音，名木野晴暢，足立忠晴，樋口理宏：軸方向傾斜機能材料からなる柱部材の線形座屈特性に与える縦弾性係数の変化の影響，第5回九州橋梁・構造工学研究会シンポジウム論文集 (USB), R4-1 (6 pages), 2017.

5) C.M. Wang, C.Y. Wang, J.N. Reddy: Exact solutions for buckling of structural members, CRC Press, pp.9-75 (chapter 2), 2004. 6) 渡辺力：p-Ritz法における試行関数の改良に関する一 考察，構造工学論文集A, 61A巻, pp.1-10, 2015. 7) 福本秀士：構造物の座屈・安定解析 (新体系土木工学 9), 技報堂出版, pp.37-98, 1982. 8) T.R.トーカート (著)．岩本卓也 (翻訳)：構造力学とエネルギ原理 (理工学海外名著シリーズ32), ブレイン図書出版, pp.76-77, 1979. 9) C.L.ディム，I.H.シャームス (共著)，砂川恵 (監訳)：材料力学と変分法 (理工学海外名著シリーズ26), ブ

(12)

―49―

レイン図書出版, pp.94-151, 1977.

10) S.P. Timoshenko and J.M. Gere: Theory of Elastic Stability (Dover Civil and Mechanical Engineering), Dover Publications (2nd ed. edition), pp.100-113, 2009.

11) T. Huang and D.W. Dareing: Buckling and Frequencies of Long Vertical Pipes, Journal of the Engineering Mechanics Division (ASCE), Vol. 95(1), pp. 167-182, 1969.

12) 杉山吉彦，芦田幸逸，川越治郎：自重による長柱の座屈, 日本機械学會論文集, 43 巻, 376号, pp.4435-4443, 1977.

13) W.H. Duan and C.M. Wang: Exact solution for buckling of columns including self-weight, Journal of Engineering Mechanics (ASCE), Vol. 134 (1), pp. 116-119, 2008. 14) C.M. Wang and K.K. Ang: Buckling capacities of braced

heavy columns under an axial load, Computers & Structures, Vol. 28 (5), pp.563-571, 1988.

15) 土木学会構造工学委員会構造力学公式集改訂委員会編：構造力学公式集 (昭和61年版), 土木学会, p.118, 1986.

仮想ばね法により基本境界条件を緩和させたB-spline Ritz法 に設定する最適な仮想ばね定数の値について（その２） 自重の作用による均質な等断面長柱の線形座屈問題