はりの計算式を解くやり方をいろいろと紹介しています。
本稿の弾性曲線式は
V''''(x) = p(x), x | 0 => L (technical theory of beam, Euler–Bernoulli beam theory)
V''''(x) = p(x) - κ * EI / GA * p''(x), x | 0 => L (Timoshenko beam theory)
と改めます。上は普通の梁で、下はせん断変形を考慮した梁の式です。
はりの計算式は
V''''(x) = p(x) = -w (1)
V''' (x) = Q(x) = C1 + ∫ p(x) dx (2)
V'' (x) = M(x) = C2 + ∫ Q(x) dx (3)
r(x) = -M(x) (4)
G(x) = κ * EI / GA * Q(x) (4-1)
V' (x) = T(x) = C3 + ∫ r(x) dx (5)
V (x) = V(x) = C4 + ∫ T(x) + G(x) dx (6)
ただし
κ せん断変形の形状係数
EI 曲げ剛性(一定)
GA せん断剛性(一定)
x 位置
L スパン長
p(x) 分布荷重
w 等分布荷重
Q(x) せん断力
M(x) 曲げモーメント
r(x) 曲率(δ''=φ=r(x)/EI)
G(x) せん断角(γ=G(x)/EI、鉛直方向のせん断歪は全断面で一定とし、材軸方向は０とする)
T(x) たわみ角(θ=T(x)/EI)
V(x) たわみ(δ=V(x)/EI)
C1 積分定数(Q(x=0))
C2 積分定数(M(x=0))
C3 積分定数(T(x=0))
C4 積分定数(V(x=0))
とします。G(x) = 0 とすれば普通の梁です。
この計算式は 荷重 p(x) を 与えると 機械的に解けます。パソコンを利用すれば 境界条件 C1, C2, C3, C4 を計算することは それほど難しくはありません。数式処理ソフトを使えば、手計算のように 静定 と 不静定を区別したり 仮想仕事や相反作用の原理を 持ち出す必要がないからです。また、普通の梁でもせん断変形を考慮した梁でも同じ計算式が使えるので便利です。
本稿のアイデアは たわみ角 T(x) と たわみ V(x) の与え方にあります。曲げ剛性とせん断剛性が一定のとき、T(x) と V(x) を たわみ角とたわみとみなすことで 弾性曲線式は積分を繰り返すとその都度 Q(x)、M(x)、T(x)、V(x) となります。たわみ角とたわみは T(x)/EI と V(x)/EI とすればもとに戻せます。
今回は 強制変位について考えます。
Mery で以下のスクリプトを書いて ファイル名 "000-8.bat" で保存して実行すれば強制変位が処理されます。
単純支持として、はりの中間に強制変位を与えます。
: 強制変位
@echo off
path C:\maxima-5.47.0\bin;%path%
ruby -x %~f0 5 2 -1.rad -50*0 0.15 0.3 6.5 6.5/15 1.2
pause
goto:eof
#!ruby
def pt_qm_tp_3(l, a, st, sv, t, h, e, g, κ)
=begin
sv SV は(x=a)の強制たわみ mm 下向きが負、上向きが正(単位を m に補正して EI で除して計算)
st ST は(x=a)の強制回転角 rad 右回りが負、左回りが正(EI で除して計算)
=end
f = open "| maxima --very-quiet", "w"
f.print <<~Maxima
[p1, p2]: -[0, 0]$ /*(M1)*/
[Q1, Q2]: [C1 + integrate(p1, x), C5 + integrate(p2, x)]$ /*(M2)*/
[M1, M2]: [C2 + integrate(Q1, x), C6 + integrate(Q2, x)]$ /*(M3)*/
[r1, r2]: -[M1, M2]$ /*(M4)*/
[G1, G2]: κ * EI / GA * [Q1, Q2]$ /*(M4-1)*/
[T1, T2]: [C3 + integrate(r1, x), C7 + integrate(r2, x)]$ /*(M5)*/
[V1, V2]: [C4 + integrate(T1 + G1, x), C8 + integrate(T2 + G2, x)]$ /*(M6)*/
PT: float(#{st})$
PV: float(#{sv})$
co: if PV = 0.0 and PT # 0.0 then
[ev(Q1, x = a) = ev(Q2, x = 0), ev(V1, x = a) = ev(V2, x = 0),
ev(T1, x = a) = ST * EI, ev(T2, x = 0) = ST * EI]
elseif PV # 0.0 and PT = 0.0 then
[ev(M1, x = a) = ev(M2, x = 0), ev(V1, x = a) = ev(V2, x = 0),
ev(V1, x = a) = SV * EI/1000, ev(V2, x = 0) = SV * EI/1000]
elseif PV # 0.0 and PT # 0.0 then
[ev(T1, x = a) = ST * EI, ev(T2, x = 0) = ST * EI,
ev(V1, x = a) = SV * EI/1000, ev(V2, x = 0) = SV * EI/1000]$ /*(M7-co)*/
s: solve(append(
[ev(M1, x = 0) = 0, ev(V1, x = 0) = 0, ev(M2, x = L - a) = 0, ev(V2, x = L - a) = 0], co)
, [C1, C2, C3, C4, C5, C6, C7, C8])$ /*(M7)*/
[Q1, M1, T1, V1, G1]: ev([Q1, M1, T1*1000/EI, V1*1000/EI, G1*1000/EI], s)$ /*(M8)*/
[Q2, M2, T2, V2, G1]: ev([Q2, M2, T2*1000/EI, V2*1000/EI, G2*1000/EI], s)$
Q: if x < a then Q1 else ev(Q2, x = x - a)$ /*(M8-1)*/
M: if x < a then M1 else ev(M2, x = x - a)$
T: if x < a then T1 else ev(T2, x = x - a)$
V: if x < a then V1 else ev(V2, x = x - a)$
G: if x < a then G1 else ev(G2, x = x - a)$
[L, a, ST, SV, t, D, e, g, κ]: [#{l}, #{a}, #{st}, #{sv}, #{t}, #{h}, #{e}, #{g}, #{κ}]$ /*(M9)*/
[EI, GA]: [e*t*D^3/12, g*t*D] * 10^6$ /*(M10)*/
plot2d([Q, M, T, V], [x, 0, L], [legend, "Q", "M", "T", "V"])$ /*(M11)*/
?sleep(5)$ /*(M12) ５秒表示*/
quit()$ /*(M13) 終了*/
Maxima
f.close
end
class Numeric
def rad
self * Math::PI / 180.0
end
end
l, a, st, sv, t, h, e, g, κ = ARGV.map { |x| eval(x) } # (R14) m, m, rad, mm, m, m, kN/mm2, kN/mm2
pt_qm_tp_3(l, a, st, sv, t, h, e, g, κ) # (R15)
__END__
(M11)で
と表示されます。-1.rad は 1000倍に拡大しているため -17.453 となっています。
座標系は右手右ねじ系としています。
(x) は材軸方向
(y) は上下方向で上向きが正、下向きが負
(z) は回転方向で左回りが正、右回りが負
y
|
|___________x
／
z
慣用の座標系も右手右ねじ系ですが向きが違います。
(x) は材軸方向
(y) は上下方向で下向きが正、上向きが負
(z) は回転方向で右回りが正、左回りが負
z
／
／__________x
|
|
y
次回は 両端変位について取り上げます。
