・長文です。申し訳ありません。
・ChatGPT先生に、
(Beamモジュールで)SingularityFunction がでないようにする方法を教えて下さい。
私は、l > a の設定がわかりませんでした。
・ChatGPT先生に、
カスティリアノの定理のdummy load? dummy moment? について聞きました。
・ページ分割の予定はありません。
・aを右端からださなくても。左端からデモ? 右端からだと計算しやすいカモ。
問題:公式ホームページ
・(検索中)
sympyで(文献を参考に)
・(勉強中) 微分方程式で。
sympyで(片持ちばりの公式で)
参考
(本日の)ChatGPT先生に教えてもらいました。
・私の質問
SingularityFunction がでないようにする方法を教えて下さい。
# ver1.1
from sympy import *
from sympy.physics.continuum_mechanics.beam import Beam
var('E,I,W')
var('l,a',positive=True)
var("RA,MA,x")
def my_Beam(l,a,E,I,W):
b=Beam(l,E,I)
b.apply_load(MA,0 ,-2)
b.apply_load(RA,0 ,-1)
b.apply_load(W ,l-a,-1)
b.bc_deflection.append((0,0))
b.bc_slope .append((0,0))
b.solve_for_reaction_loads(RA,MA)
return b
b=my_Beam(l,a ,E,I,W)
print("#",b.slope() .subs({x:l}).simplify()) # W*(l-a)**2 /(2*E*I)
print("#",b.deflection().subs({x:l}).simplify()) # W*(l-a)**2*(2*l+a)/(6*E*I)
# W*(a**2 - 2*a*l + l**2 - SingularityFunction(0, -a + l, 2))/(2*E*I)
# W*(a**3 - 3*a*l**2 + 2*l**3 - 3*l*SingularityFunction(0, -a + l, 2) - SingularityFunction(0, -a + l, 3))/(6*E*I)
↑↑↑↑↑↑↑↑↑↑↑↑↑↑↑↑↑↑↑↑↑↑↑↑↑↑↑↑↑↑↑↑↑↑
・私の質問だけを、markdown形式で出力してもらいました。
1. SingularityFunction がでないようにする方法を教えて下さい。
2. 方法①でできません。
3. 全体で
4. Piecewise もださないで
5. Heavisideがでています。
6. 全体のソースコードをお願いします。
7. x をださないで
8. 分子も簡略化して
9. x が
10. 分子を簡略化して
11. 本チャット、この会話で出てきた全ての質問を、コピーできるいつもの Markdown形式 のテキストで一覧にしてください。追番も。空白行不要です。
「SingularityFunction も出さない、Piecewise も出さない、Heaviside すら出さずに、最終的に 普通の式（多項式だけ） を出したい」ですね。
(途中省略) def my_Beamにコメントをつけてもらって、ウレシイ。
# ver1.2
from sympy import *
from sympy.physics.continuum_mechanics.beam import Beam
# 変数定義
E,I,W,l,a,RA,MA = symbols('E I W l a RA MA', positive=True)
# ビーム定義関数
def my_Beam(l,a,E,I,W):
b = Beam(l,E,I)
b.apply_load(MA,0,-2) # 固定端モーメント
b.apply_load(RA,0,-1) # 支点反力
b.apply_load(W ,l-a,-1) # 集中荷重
b.bc_deflection.append((0,0))
b.bc_slope.append((0,0))
b.solve_for_reaction_loads(RA,MA)
return b
# SingularityFunction を多項式に変換（Heaviside を残さない）
def SF_to_poly(expr):
SF = SingularityFunction
def _repl(arg):
xx, aa, n = arg.args
n = int(n)
return ((xx-aa)**n/factorial(n)) # Heaviside 無し
return expr.replace(lambda e: e.is_Function and e.func==SF, _repl)
# ビーム解析
b = my_Beam(l,a,E,I,W)
# slope, deflection at x=l のみ
slope_at_l = SF_to_poly(b.slope().subs({'x':l}))
defl_at_l = SF_to_poly(b.deflection().subs({'x':l}))
# 分子も整理（factorで共通因子をまとめる）
slope_final = simplify(factor(slope_at_l))
defl_final = simplify(factor(defl_at_l))
# 出力
print("slope(l) =", slope_final)
print("deflection(l) =", defl_final)
# slope(l) = W*(a - l)**2/(4*E*I)
# deflection(l) = W*(a - l)**2*(5*a + 4*l)/(36*E*I)
sympyで。ならば、ｌ_ACとl_CB(=a)で。作図も。
・別に,SingularityFunction も Heaviside を使わなくてもです。
# ver1.3
from sympy import *
from sympy.physics.continuum_mechanics.beam import Beam
var('E,I,W')
var('l,l_AC,a',positive=True)
var('RA,MA,x')
def my_Str_al2la(mySiki):
return str(mySiki).replace("-a + l","l - a").replace("a + 2*l","2*l + a")
def my_Beam(l_AC,a,E,I,W):
b=Beam(l_AC+a,E,I)
b.apply_load(MA,0 ,-2)
b.apply_load(RA,0 ,-1)
b.apply_load(W ,l_AC,-1)
b.bc_deflection.append((0,0))
b.bc_slope .append((0,0))
b.solve_for_reaction_loads(RA,MA)
return b
def my_Sakuzu(l_AC,a,E,I,W):
b=Beam(l_AC+a,E,I)
b.apply_load(MA,0 ,-2)
b.apply_load(RA,0 ,-1)
b.apply_load(W ,l_AC,-1)
b.bc_deflection.append((0,0))
b.bc_slope .append((0,0))
b.solve_for_reaction_loads(RA,MA)
b.plot_shear_force()
b.plot_bending_moment()
b.remove_load(MA,0,-2)
b.remove_load(RA,0,-1)
b.apply_support(0,'fixed')
b.draw().show()
return
b =my_Beam(l_AC,a ,E,I,W)
rep1={x :l_AC+a}
rep2={l_AC:l-a}
print("#",my_Str_al2la((b.slope() .subs(rep1).simplify().subs(rep2))))
print("#",my_Str_al2la((b.deflection().subs(rep1).simplify().subs(rep2))))
# W*(l - a)**2/(2*E*I)
# W*(l - a)**2*(2*l + a)/(6*E*I)
#
my_Sakuzu(0.7 ,0.3 ,1,1,1)
・アメリカ式 Moment を微分したら、Force。いいですね。
私はドイツ式を理解できていません。
・ver1.4 (ver1.3 と同じ。最後で、lとaにする。ｌ>a を考えなくていいです。)
# ver1.4
from sympy import *
from sympy.physics.continuum_mechanics.beam import Beam
var('E,I,W')
var('l,a,l_AC,l_CB',positive=True)
var('RA,MA,x')
def my_Str_al2la(mySiki):
return str(mySiki).replace("-a + l","l - a").replace("a + 2*l","2*l + a")
def my_Beam(l_AC,l_CB,E,I,W):
b=Beam(l_AC+l_CB,E,I)
b.apply_load(MA,0 ,-2)
b.apply_load(RA,0 ,-1)
b.apply_load(W ,l_AC,-1)
b.bc_deflection.append((0,0))
b.bc_slope .append((0,0))
b.solve_for_reaction_loads(RA,MA)
return b
b =my_Beam(l_AC,l_CB ,E,I,W)
rep1={x :l_AC+l_CB}
rep2={l_AC+l_CB:l,l_AC:l-a,l_CB:a}
print("#",my_Str_al2la((b.slope() .subs(rep1).simplify().subs(rep2))))
print("#",my_Str_al2la((b.deflection().subs(rep1).simplify().subs(rep2))))
# W*(l - a)**2/(2*E*I)
# W*(l - a)**2*(2*l + a)/(6*E*I)
(本日の)ChatGPT先生に教えてもらいました。
・私の質問
以下のように、カスティリアの定理で、仮想荷重、仮想モーメントを同時に追加したら問題がありますか。
# ver1.5
from sympy import *
var('W,Wi,Mi,l,a,x,E,I')
def my_Str_al2la(mySiki):
return str(mySiki).replace("-W*(a - l)","W*(l - a)").replace("a + 2*l","2*l + a")
MAC=W*(l-a-x)+Wi*(l-x)+Mi
MCB= +Wi*(l-x)+Mi
U = 1/(2*E*I)*integrate(MAC**2,(x,0 ,l-a)) \
+1/(2*E*I)*integrate(MCB**2,(x,l-a,l ))
rep={Wi:0,Mi:0}
print("#",my_Str_al2la(-diff(U,Mi).subs(rep).simplify()))
print("#",my_Str_al2la(-diff(U,Wi).subs(rep).simplify()))
# # W*(l - a)**2/(2*E*I)
# # W*(l - a)**2*(2*l + a)/(6*E*I)
問題ありません。結論を先に言うと：
・仮想モーメントMiとと仮想集中力Wiを同時に梁の曲げモーメント式に入れておき、エネルギーUを作ってから目的とする仮想量で偏微分し、最後に全ての仮想量を 0 に戻す(Wi=,Mi=0)
という手順は Castigliano の定法に対して正しいです。
・そのとき得られる結果は「偏微分した仮想量が物理的に表す変位（回転・たわみ）に対応する」点に注意してください。...
いつもの? sympyの実行環境と 参考のおすすめです。
いつもと違うおすすめです。
・(いつもと違うおすすめです。)Youtube
??? Qiitaの途中リンクの調子が悪いです。
sympyのdoc
参考文献
・p51
2019年度京都大学大学院工学系研究科機械工学専攻修士課程の入試問題 3-1
・(旧版) P130
表5.1 片持ちばりの公式 6)