はじめに
立体要素関係（線形、弾性、03/12/2更新）
CalculiX用データ作成(CGX色スケール変更、グレースケール化)
Blender-CalculiX入力データ作成
(シェル要素、ダイヤカット缶、PCCPシェル)
Salome-CalculiX入力データ作成
3Dプリンター用データ作成用
梁要素関係（幾何学非線形、弾性、
弾塑性）
小物(せん断力図、モーメント図、たわみ図描画など)
板要素関係（幾何学非線形、弾性）
今後の予定
構造景観関係？？（番外編）
線形回帰
代数計算等サブルーチン
その他(学生の番号をランダムに並び替え)
成績処理関係?？？（更に番外編）
BASIC(15個の振り子)
15個の振り子の動きの主題によるミニマル・ミュージック
その他の他（私以外、構造研以外の学生等によるプログラム等）
手引き書関係のリンク
解説関係のリンク
FEMソフト関係のリンク
私が自分の研究用や学生の教材用に作った Fortran77* の プログラム （一部、学生を指導して作ってもらった、 より正確には、命令して作らせたものも含む） を少しずつ整理して、 ここに公開していこうかと考えています。 以下のプログラムは、一応、無責任／無保証のフリーソフトですので （著作権参照）、 以下のプログラムを使って何か問題が起きても （例えば、計算結果を信じて構造物を設計したら壊れたとか その下敷きになった人が死んだとか、 バグでおかしな結果が出たのを真に受けて 「常識を覆す現象が数値的に確認された」とトンデモ論文を 書いて恥を書いたとか）、 一切、責任は取りません。 プログラムのオリジナルをここから取ってきたことが 分かるように明示されているなら、 改造や再配布は自由ですが （仮にもそんなことをしたい人がいらっしゃるなら* ）、 「改造した新たなプログラムの改造や再配布を禁止すること」 は禁止します（基本的に コピーレフトの考えに準拠します。著作権参照）。
*注意： Cygwin や Linux の gcc (GNU Compiler Collection) のフォートラン77コンパイラー g77 では end do を使える ことが分かってからは、 end do の doループも（行番号付きのdoループと）混在させて 使ってたりします。 過去に作った使えそうなプログラムを流用しながら、 入り用のものを（取り敢えずは動けばいいから）手っ取り早く作るという 思想で作っているので、 変数は宣言せずに後からどんどん付け加えているし、（改造に便利な） goto文（計算型も含む）で気兼ねなくあちこちに飛びまくってるし、 メインはやたらと長いし……といった具合に、時代の流れとは逆行する 非教育的な(典型的な) スパゲッティープログラム群です（念のため）。
バグの報告は、 こちらまで。
今後の方針： Linux上（または Windows の Cygwin 上）の gcc の フォートランコンパイラー（現行ではg77）を使い続けたいと考えているが、 g77はGCC 3.4.xまでで，GCC 3.5以降はg95 (ではなく、厳密には そこから分岐?したgfortran)をサポートするとのこと なので、ゆくゆくは gfortran への転向を視野に入れておきたい。 と言っているうちに Vine Linux4.0でg77の他にgfortranも 入れられるようになったので、少しずつ調査開始 (g77メモ 参照。日本語はutf-8でないとだめみたい)。
tyokus.f
（tyokus.html02/8/1版）：
直交異方性材料の直方体構造物を、
8節点24自由度の
直方体要素（直方六面体 立体要素）の
線形／弾性の剛性行列で解くプログラム。
剛性マトリクスはスカイライン法で解いている。
将来的に非線形解析＋弧長法とかに移行することもふまえて、
敢えて対称行列を非対称行列と見なしたまま解いてる
（というか、単に、対称行列の解法に特化する気力がなかった。これを
やれば更に速くなり、最大分割数も増やせるだろう）。
tyokum.f
（tyokum.html02/8/1版）：
tyokus.f の剛性方程式を解くサブルーチンをメインプログラムに入れることで、
配列サイズをもう少し大きく切れるようにした
（もっと、いいやり方がある筈だけど、取り敢えず）。
cygwin の g77.exe の最新版では（あるいは XP との組み合わせでは？）、
こんなにでかい配列は切れなくなってたりするようです
（コンパイルはできるけど、実行してもエラーになるとか）。
そのときは、適宜 配列を小さくして下さい。
tyoku.f
（tyoku.html表示）：
tyokus.f をスカイライン法にする前のプログラム。
剛性マトリクスをフルマトリクスで解いているので、
tyokus.f の数十倍ぐらい計算速度が遅い（改造、教材、比較用）。
tyokus.f, tyokum.f, tyoku.f への入力データの与え方の例は、
hipp.f
（hipp.html, 02/8/1版）：
片持ち梁の先端を引っ張る計算をtyokus.fにさせるための入力データを
作るプログラム。
mage.f
（mage.html, 02/8/1版）：
片持ち梁の先端に集中荷重を与える計算をtyokus.fにさせるための入力
データを作るプログラム。
mage5p.f（06/11/27版）：
tyokum.f用の5点曲げ試験の入力データを作成(調整中)。
tyoku2関係
tyoku2.f：tyokum.f を２種類の材料で
解けるように改造（調整中）。これは鋼材部分の要素数を固定させた試験例なので、
下記のtyoku2ou.fとmage2ou.fを使うこと。
mage2.f：
集成材梁の上下面に軸方向に溝を掘って鋼リブを挿入して補強した剛性梁の
曲げをtyoku2.fに半解析させるための入力データを作るプログラム（調整中）。
tyoku2ou.f
tyokum.f(というより下記のtyokuon.f)を2種類の材料に対応させ、
ouryogu2.fでgnuplot用の応力表示の
出力に対応。
mage2ou.f
上下に同じ桁高の鋼板を挿入した集成材梁の片持ち梁の曲げを半解析で
tyoku2ou.fに解かせるデータを作成するプログラム。
tyokuon.f：tyoku2.f の全節点に
温度荷重などを与える場合を想定して、荷重ベクトルは全節点を一括して
入力するように変更。
gnuplot
の3次元プロット用の出力を追加（03/5/6）
cygwin の g77.exe の最新版では（あるいは XP との組み合わせでは？）、
こんなにでかい配列は切れなくなってたりするようです
（コンパイルはできるけど、実行してもエラーになるとか）。
そのときは、適宜 配列を小さくして下さい。
hippou.f：
集成材梁の上下面に軸方向に溝を掘って鋼リブを挿入して補強した剛性梁
（片持ち梁）の全ての箇所の温度が均等にΔT℃変化し、
集成材の全ての箇所の含水率が均等にΔH％だけ変化した場合の
各節点の変位をtyokuon.fに1/4解析させるための入力データを作るプログラム。
鋼材のみ、木材のみで解いた限りにおいては、
与えた線膨張係数から予想される軸方向変位、軸直角方向変位が求まることは
確認した（調整中02/11/28）。ouryoku.fで応力を求めるのに必要な
行列を出力するように変更し名前をhippon.fからhippou.fに改めた（03/1/9）。
ouryoku.f：
tyokuon.fで求めた変位に、ひずみ−変位行列と応力−ひずみ行列をかけて
各要素の応力を求める（調整中03/1/9）。
gnuplot
の３次元プロット用のデータを出力するように書き換えた
ouryogu.f。
tyokuon3.f：
tyokuon.f を3種類の材料で解けるようにした。
hippon3.f：
鋼板と集成材の間に接着剤要素を入れる場合。
接着剤は木材に比べて十分に堅いので、接着剤なしの場合（hippou.f）と
特に差がない。
hipponx.f：
接着剤要素を入れるよりは、鋼板と隣接する集成材の要素分割を
細かくした方が、変形の大きい鋼板隣接部付近の乱れが小さくなるかも知れないと、
hippon3.f の接着剤要素の材料定数を集成材要素の材料定数と
置き換えてその部分の集成材だけ要素分割を細かくできるようにして
みた。しかし、やはり鋼板にじかに隣接している要素の変形が大きくなり、
変位や応力に乱れが入る。
hippon0.f：
という訳で、hipponx.f の鋼板隣接部付近の集成材要素には、
乾燥収縮荷重を与えないで解いてみることにした。
すると、多少は鋼板隣接部付近の集成材の変形が、
段階的に変化するようにはなった。でも今一つ。
ouryogu3.f：
tyokuon3 < hippon3.d（または hipponx.d または hippon0.d）で
出力されるファイルから応力の３次元プロット（gnuplot の splot）用の
データファイルを作る。
ccxkataz.f
ccxkataz.f90(gfortran用)：
FEMソフト関係のリンク
のとこでも紹介している
CalculiX
で、せん断を受ける片持ち梁を直方体要素で解析するための例題(beam8p.inp)
の要素分割を任意に設定できるようにし、節点番号も、梁断面の左上から
右下に順番に振り直すようにしたデータ作成プログラム。
./ccxkataz>ccxbeam8p.inp みたいにしてデータを作成し、
ccx_1.4 ccxbeam8p のように実行。
ちなみに生成されるデータ例は、ccxbeam8p.inp。
wirte文とformat文をfortran77的に書いた例は
ccxbeam8p.f。
ccxkataza.f: 上記を線形座屈解析用に修正中。
ccxkataiza.f: ccxkataza.fを2軸対称I型梁用に改造中。
ccxkata.f： 上記のccxkataz.fを左右対称条件から半解析とし、 固定端は軸方向変位の他には中立軸の鉛直変位のみを固定し、それ以外は、 固定端面に沿って動けるようにした。 自由端の荷重は、beam8p.inpと同様に、すべての節点に同じ荷重が分配されているので、 1要素辺りの荷重の総計が同じになるように重みを付けて分配する方式に 変更する予定。
gyaku.f： ccxで解いた片持ち梁の荷重とたわみを入力して、せん断補正係数を逆算する プログラム。 ただ、ccxの変位出力が有効数字5桁だから、逆算精度は落ちるかも知れない。 ccxの変位出力の桁数を増やすオプションってあるのかな。
ccxonsaito.f
ccxonsaito.f90(gfortran用utf-8, 片持ち梁対応)
プレストレス床版を想定した木部材2段を左右の鋼板で固定した箱桁を
半解析で解くためのinpファイルを生成。自由形式。
ccxbunpu.f： 木材の3点曲げ試験を想定して、 支点部分にめり込み防止の鋼板を敷いて、 載荷部分を等分布載荷するinpファイル生成。自由形式。
その他、CalculiX関係の入力データ作成プログラム各種は、 千田さん作のものが、 この辺にあります。
tasudat.f (自由形式utf-8はtasudatw.f) ： ccxで計算した変位や座屈モードを初期座標に足し算して、 変形図やモード図をcgxに描かせるためのデータを作るプログラム。 ccx1.5とか古いバージョンでは、frdファイルの変位にモードを足す tasu.fで対処できていたが、 ccx1.7とかでは、frdファイルに(シェル要素だと?)拡張された節点等が 書き込まれるようになってしまったので、inpファイルの節点座標に モードを足すようにした。
cgxppm.f：
cgxで描いた応力図などは、応力値の最大値から最小値までが
(赤、黄、緑、青、紫みたいな感じの)21段階の色スケールで表示される。
二つ以上の応力図を比較する場合、
それぞれの色スケールが(それぞれの最大、最小に対応して)違うので、
色スケールを同じにして比較したいということがある。
そんな場合に、最大値から最小値の幅が広い方の図の最大値の赤を+10,
中央値の緑を0, 最小値の紫を-10として、
その色スケール内で比較したい幅の狭い方の図の最大値、最小値を
整数で与えると、色スケールを(最大、最小の幅の広い方に合わせて)
変更したカラー図とグレースケール図を出力する。
例えば、最大、最小の広い図
tau1.pngと
最大、最小の狭い図
tau2.pngがあった場合、
まず、
tau1.pngをgimp等でアスキーppm化して
cgxinp.ppmというファイル名で保存し、./cgxppm を実行すると、
最大、最小を聞かれるから、これのグレースケールが欲しい場合には、
10,-10と入力する。
すると、cgxout.ppmとcgxout.pgmが得られる。
最大、最小を変化させない場合、
cgxout.ppmは入力ファイルと同じである。
cgxout.pgmはグレースケールファイルで、これをgimpなどでpng化
したものがtau1g.pngである。
次に、最大、最小の狭い図
tau2.png
を同様にgimp等でアスキーppm化して
cgxinp.ppmというファイル名で保存し、./cgxppm を実行し、
最大、最小を聞かれたら、
tau2.pngの応力の最大値、最小値が、
tau1.pngの色スケール上では、
最大：赤(+10), 中央：緑(0), 最小：紫(-10)の21段階のどの辺に対応するかを
読み取って入力する。
たぶん、10,-1ぐらいだと思うので、そう入力する。
そうして得られたcgxout.ppmをpng化したのが
tau2c.png,
cgxout.pgmをpng化したのが
tau2g.png.
という訳で、色スケールを(大体)合わせたカラーでの比較：
tau1.pngと
tau2c.png
(元の色スケールtau2.png).
それをグレースケールにして比較：
tau1g.pngと
tau2g.png.
しかし、cgxの色スケールは、カラーの図を白黒印刷すると
グレースケールにならないという大きな問題がある
(私の色感性だと、赤-黄-緑-青-紫というグラデーションも
今一つ分かりにくい)。
そこで、白黒印刷してもある程度はグレースケールになるような
(白-黄-緑-黒の)色スケールの出力cgxout2.ppmも出すようにした。
という訳で、白黒印刷してもグレースケールがある程度は成立するカラー図での比較：
tau1c2.pngと
tau2c2.png
(元の色スケールtau2.png).
Blender2.40でVideoscapeでエクスポートしたobjファイルから、 ccxの8節点シェル要素と6節点シェル要素のinpファイルを作るプログラム群。
kakutyuu.f(EUC-JP)
kakutyuu.f(utf-8)
角柱や円筒のobjファイルを吐き出す。
ble6ccx.f90を通せばinpファイルを作れる。
ダイヤカット缶(PCCPシェル)をccxの6節点シェル要素で 解くinpファイルを作るプログラム群。
ble6ccx.f90(gfortran用utf-8):
tubure6.f:
n角形の折り畳み円筒1段ぶんのobjファイルを吐き出すプログラム
tubure.f:
n角形の折り畳み円筒n段ぶんのobjファイルを吐き出すプログラム
作成の経緯など
以上のプログラムを利用しながら、 工藤さんが作成した円筒折り紙構造各種を解析するプログラム群はこちら
unvc3d4.f90(Salomeの四面体メッシュunvファイルから
inpファイルを生成する滝田さんのプログラムを修正して、体積を表示するようにした
)
以下は、Gmshを経由する方式。
Gmshを経由せずにunvから直接inpを作る滝田さんのプログラムがオススメ。
ccxc3d8.f90(直方体要素C3D8用)
ccxs6.f90(シェル要素S6用、以下のccxc3d8はccxs6に読み替えて)
Salomeで直方体要素のメッシュを切り、それをUNV形式で保存したものを
Gmshでマージで読み込んでからAbaqus inp形式で保存すると、
salome.inp
みたいなファイルができる。
salome.inp
の冒頭の節点番号と座標の書かれた部分を
setten.inp
にコピーして冒頭に節点数を書き込み、
salome.inp
の末尾の要素番号と要素を囲む節点番号の書かれた部分を
youso.inp
にコピーして冒頭に要素数を書き込み、
ccxc3d8を実行すると、
ccxc3d8.inp
みたいなファイルができてるので、後は、
cgxを利用するなりプログラムを組むなりして拘束節点と載荷節点を設定する。
シェル要素用のccxs6.f90は、
Salomeが吐き出した3頂点の三角形要素に対してble6ccx.fの要領で
中間節点を振っている。
ccxs6s3.f90
ccxs6で生成した中間節点を利用して、入力データの三角形を
4分割した三角形メッシュデータを出力する。
s3setten.inpとs3youso.inpに下記のccxs6の要領でメッシュデータを
書きこんでccxs6s3を実行すると、4倍の分割のメッシュを
s3setten.inpとs3youso.inpに上書きする。
ccxs6s3をもう一回実行すると、4倍メッシュのs6要素用inpファイルを
画面出力し、16倍メッシュのs3データを
s3setten.inpとs3youso.inpに上書きする。
実行するたびにメッシュを4倍ずつ細かくできる。
例えば、
柴田さんの
ハニカムメッシュ生成プログラムで生成したs3ハニカムメッシュを
ccxs6s3に2回か3回通して、
16倍メッシュにした例。
zya3dinp.f90
daiya3dinp.f90
シェル要素のような厚さのないデータの場合、
inpファイルをSalomeやMentatに読み込んでから厚さを指定しても、
厚さのある充実物体のモデルにはなってくれないので、
3Dプリンターでモデル生成できない。
そこで、蛇腹折り円筒やダイヤカット円筒を
シェル要素ではなく、厚さを指定して4面体要素で分割したinpファイルを
吐き出すプログラム。調整中。このinpファイルをMentatに読み込んでエクスポートした
stlファイルはSalomeに読み込むとエラーが出るそうなので、stlを直接 出力するように修正中。
daiya3dstl.f90
zya3dstl.f90
zya3dstlt.f90(implicit*4)
hari.f
（hari.html）：
弾性（曲がり）梁の座屈や大変位問題を、直線梁要素の
有限変位・有限要素法で解くプログラム
（スカイライン法＋弧長法）。
入力データの与え方の例は、
oiraa.f:
柱のオイラー座屈をhari.f で解かせるための
入力データを作るプログラム。
oiraa2.f:
柱のオイラー座屈を二次元問題としてhari.f で解かせるための
入力データを作るプログラム
(正方形断面や円断面などの二軸対称断面を座屈させるには、
面外変位を拘束して二次元問題にしないと解けない、たぶん)。
enko.f
（enko.html）：
等曲げを受ける円弧アーチの横ねじれ座屈をhari.f で解かせるための
入力データを作るプログラム。
kata.f
（kata.html）：
端部荷重を受ける片持ち梁の横ねじれ座屈をhari.f で解かせるための
入力データを作るプログラム。
katakukei.f
上記kata.fのねじり定数とそりねじり定数をちゃんと求めてから、
座屈荷重を計算するように修正。
katai.f
kata.fのI型梁用。
学生時代に作った冗長な（メインのやたらと長い）プログラムを、 その後、様々な目的に応じて書き加えたりしているうちに、 つぎはぎだらけの分かりにくいプログラムになってしまった （一応、動くことは動きますが）。 せめて、もう少し整理しようかと、プログラムを解読していると、 「過去の自分は一体 何を考えてこういう不可解な書き方をしたのだろう？」 と思わされる箇所が随所に満ちている （だから、私が書き込んだ解説は過去の私の意図？を誤解している 可能性もあります）。 ボツにならなかった計算例はこの辺とか 。
harisen.f
(SJIS版はharisen.f)：
弾性（曲がり）梁の座屈や大変位問題を、直線梁要素の
有限変位・有限要素法で解くプログラム hari.f を
二軸曲げに対してせん断変形を考慮できるように修正。
hari.f では初等梁理論の剛性行列を用いていたのを、
Timoshenko梁の剛性行列に変えた（調整中）。
oirasen.f:
柱のオイラー座屈をharisen.f で解かせるための
入力データを作るプログラム。
katasen.f
(SJIS版はkatasen.f)：
端部荷重を受ける片持ち梁の横ねじれ座屈を harisen.f に解かせるための
入力データを作るプログラム。
enkosen.f：
等曲げを受ける円弧アーチの横ねじれ座屈（修正Vlasovの解析解）を
harisen.f に解かせるための入力データを作るプログラム。
これを等曲げを受ける直線梁の横ねじれ座屈（Trahairの解析解）に
特化したものがenkosen2.f。
hariso.f：
harisen.f を集成材梁の弾塑性解析用に改造中。
oiraso.f:
柱のオイラー座屈をhariso.f で解かせるための
入力データを作るプログラム。
kataso.f：
端部荷重を受ける片持ち梁の横ねじれ座屈を hariso.f に解かせるための
入力データを作るプログラム。
enkoso.f：
等曲げを受ける円弧アーチの横ねじれ座屈（修正Vlasovの解析解）を
hariso.f に解かせるための入力データを作るプログラム。
harisori.f：
hariso.fの軸ひずみ増分の計算にそりの影響を考慮。
harisoj.f：
harisori.fを、
座屈や破断が生じなくても、最初に除荷が生じた時点で計算を打ち切るようにした。
kukei.f
（kukei.html, 02/6/18）：
kukeiw.f(Linux用utf-8版)：
（薄肉とは見なせない）長方形断面のねじり定数、そりねじり定数
（ついでに断面積、二軸回りの断面二次モーメント）を
断面（の面積）積分で求めるプログラム。
梁の弾塑性解析をファイバーモデルとかでやるときに、
断面の分割要素ごとに（降伏判定して）弾性係数を掛けながら全断面に面積積分して
各種断面定数を
求める際のサブルーチンとしても利用できるかも知れない
（その時は断面一次モーメントとかも必要だろうけど）。
但し、プログラム内にも書いたけど、そりねじり定数については、
正解との比較をしていないので、この求め方で当たっているかどうかは不明
（薄肉と見なせない長方形断面のそりねじり定数をちゃんと求めてる文献を知ってる
人がいたら教えて下さい）。
! 長方形の核がひし形になるのかの確認 ! x軸方向の格点は、ex=Ix/(A*y)=b/6 ! y軸方向の格点は、ey=Iy/(A*x)=h/6 となるのはいいとして、 ! 斜め方向が、それらの格点を結んだ直線でいいという説明がしっくりこない。 ! exとeyを結んだ線分上に載荷された荷重Pは、ex上とey上に載荷された ! PxとPyとに分解できるから、それぞれの応力分布を重ね合わせると、 ! 作用点から最も遠い角の応力は 0 のままだから、この線分上に格点がある ! みたいな説明がなされるけど、この重ね合わせができるとうことは、 ! σ =(Px*ex)*y/Ix + (Py*ey)*x/Ix がx軸からθ 傾いた軸回りの曲げについても ! 成り立つと仮定していることになる。つまり、この式が、 ! 任意の角度で計算したP*e*r/I に等しいと。 ! ほんとうだろうか。 ! 簡単のため、正方形(b=h=1)について、頂点がx軸にある状態(θ =45度）から ! 底辺が水平になる普通の状態（θ =90度）まで、1度ずつ回転させながら、 ! 格点eを計算して、そのx座標(e*cosθ )とy座標(e*sinθ )をプロットしてみたら、 ! どうなるかと。 ! θ =45度のときの正解は、ひし形の核で、√ 2*b/12=0.11785で、これは合ってそう。 ! θ =90度のときは、分母のcosθ =0になって、おかしい答えが出るので、 ! θ =89度を見ると、e=0.161, 90度のときの正解は、長方形の核だからh/6=0.1666.. ! 両端は、まあ合ってそう。 ! で、途中の座標を gnuplot でプロットしてみたら、みごとに直線だ。
smv.f90(UTF-8版)
smv.f(EUC-JP版)
gnusmvp(gnuplotでpng画像を出力)
gnusmvf(gnuplotでXfigファイルを出力)
構造力学の演習で解かされるような簡単な梁の
せん断力図、曲げモーメント図、たわみ図をgnuplotで描かせるための
データを作るプログラムを準備中。
ちなみに、Windowsを使ってる人は岩熊先生の
cbeam
の方が遥かに有用で高機能ですので(念のため)。
smv.f の中で、境界条件と荷重条件を設定して実行すると、
sendan.d, mage.d, tawami.d, ziku.d が出力される。
gnuplot<gnusmvpとしてgnuplotを実行すれば、
こんな感じのpng画像(smv.png)が、
gnuplot<gnusmvfとしてgnuplotを実行すれば、
Xfigファイル(smv.fig)が得られる
(それをLaTeXに取り込む場合は、この辺参照)。
tokuron.f90
構造力学特論用。滝田さんのプログラムを修正中。
tora2.f90
構造情報学用。トラス要素2要素の例
torasu.f90
torasu_sj.f90(Windows用Shift JIS版)
構造情報学用。2次元トラス9要素まで。
torasu2.f90
(UTF-8版。日本語版以外のWindowsの場合は、こちらの方が文字化けしにくいかもしれません)
torasu2sj.f90(Windows日本語版用Shift JIS版)
近藤さんが上記のtorasu.f90を節点力を出せるように改造してくれたもの。
hone2.f90(UTF-8版。日本語版以外のWindowsの場合は、こちらの方が文字化けしにくいかもしれません)
hone2sj.f90(Windows日本語版用Shift JIS版)
近藤さんが下記のhone.f90を節点力を出せるように改造してくれたもの。
hone.f90
hone_sj.f90(Windows用Shift JIS版)
構造情報学用。
2次元骨組み9要素まで。
材料諸元を与えてTrahairの座屈公式で
片持ち梁の横ねじれ座屈荷重を求める。
各種実験値との比較に使ったらしい残骸は、
tuzino.f
yokonezire.f(SHFT-JIS)(たぶんTrahairの片持用)
yone.f90(utf-8でf90)(たぶんTrahairの片持用)
上記を田村さんが改造して、座屈設計ガイドラインに載っている公式で、
各種の境界条件や荷重条件に対応できるようにしたものが以下。
yonekyou.f90(utf-8でf90)(各種境界条件・荷重条件用。これは両端固定、中央集中荷重の例)
sounyuukg.f
鋼板挿入集成材梁のせん断補正係数を、
西野文雄・長谷川彰夫『新体系土木工学7構造物の弾性解析』(土木学会)に
書かれている細長い長方形断面に対する方法を参考にして、
薄木先生の方法(
薄木征三，後藤文彦，キッシュ ラヨシュ：
挿入リブ鋼板で補剛した集成材の曲げ耐荷力，
構造工学論文集，Vol. 49A, 889-894, 2003.)
およびその変種で求めてみる。
上記の面積積分dA=木材幅(b-ts)dyの方法で求めると、
挿入鋼板の深さを0にすればkは5/6になる。が、
鋼板がある状態で、集成材と鋼板のヤング率を
同じにしても、5/6より(鋼板が深いほど)大きめになる。
一方、面積積分dA=換算幅(b-ts+(Es/Ew)ts)dyの方法で求めると、
鋼板がどんなに深くても鋼板と集成材のヤング率を同じにすれば、
k=5/6になる。
なお、せん断補正係数に材料定数比を反映させるとすれば、
EIが関係するIに対してはヤング率比で、
GAが関係するAに対してはせん断弾性係数比で、
といったことを試しにやってみるとかなり小さめのkになってしまう。
後藤の今のところの見積りは、
FEMから逆算したkから察するに
材料定数比に応じてkが極端に小さくなったりはしないような気がしているが、
これについては今後も検討したい。
gyaku4.f: 逆対称4点曲げで、鋼板挿入集成材梁が 曲げより先にせん断で壊れる短さのLを算出。
hari.f と同じアルゴリズムで、 板の大変位解析のプログラムを長方形要素で 学生さん （ここに時々 出没するteruさんとかいう人、今はもう学生じゃないけど） に作ってもらったやつを そのままここに置いておきます。 kasou.f は、hari.f と同じように 局所系の荷重ベクトルと全体系の荷重ベクトルとの間の変換行列を仮想仕事式から 求めたもの、kani.f は、 それを座標変換行列で代用した簡易解法です （尤も解析結果に有意差はないようですが）。 いずれも行列のスカイライン化はしてません。 データ入力の例は、input.txt にあります。 kani.f の定式化と、その数値計算例は、 ita.pdfに書いてます （ ボツにならなかったやつはこの辺とか ）。 プログラムの作者から、こんなものは公開しないでくれという要請があれば、 ここの板関係は削除するかも。
hari.f の剛性行列を 3次元Timoshenko梁のものに取り換える（済）。
hari.f の出力から座屈モードを計算させるプログラムの例（いずれ）。
hari.f を弾塑性解析用に改造（改悪）したプログラム（矩形断面集成材梁用を今年中）。
tyokum.f で二種類の材料（木材と鋼材とか）の合成構造を解けるようにする。
更に、各要素の温度歪が与えられた時に、
それを節点温度荷重に変換して解けるようにする（済）。
tyokus.f の直方六面体要素を使って、幾何学非線形解析の
プログラムを作ってみたい（予定立たず）。
直方体の8節点の座標値から各要素の向きを
方向余弦（9成分）で表して、それで座標変換することにすれば、
回転パラメーターを使わずに、
u,v,wの節点自由度だけを変数として定式化できそうだ。
ただ、方向余弦を変位自由度で表す際に分母に2乗和のルートとかがある
項が入ってくるから、増分式はそれなりに煩雑になりそう。
ppm.f
（ppm.html, 02/12/6）
spe2dg.f
（spe2dg.html, 02/12/6）
アスキーフォーマットのppm画像ファイルの R値、G値、B値それぞれの
平均と標準偏差を求める。また、
ppm画像データの各行の輝度信号(Y=0.31R+0.59G+0.11B)の横方向の平均を
1列に並べたデータ（yokosima.d）と、
ppm画像データの各列の輝度信号(Y=0.31R+0.59G+0.11B)の縦方向の平均を
1列に並べたデータ（tatesima.d）を出力する。ここまでは
ppm.fがやる。
次に、
貴家仁志
『よくわかるディジタル画像処理』(CQ出版社)付属の
1次元FFTプログラム（fftvec）に、
yokosima.d とtatesima.d を
入力して、画像ファイルの横縞成分の振幅スペクトルyokospe.dと
縦縞成分の振幅スペクトルtatespe.dを求める。
それらをspe2dg.fに入力すると、
パワーと周波数で両対数プロット
した場合の横縞パワーと縦縞パワーの回帰直線の
傾きの平均、標準偏差の平均、傾きの差などを
定量化のための数値指標として出力する。
やろうとしてることは、例えば、官能検査？の目的に応じて こんな感じ の木材（とか構造景観とかの）画像のRGB値の平均や、フーリエ解析で 得られる縦縞や横縞のスペクトルを両対数でプロットした 回帰直線の傾きと標準偏差などなどを説明変数として、 この手のアンケートで得られる 大衆（なり特定の個人なり）の評価を目的変数として 重回帰分析したりして、 「どの辺から木材らしく見られなくなるかの閾値」とか、 大衆が好む景観の図形的特徴を数値指標で表せたりしないだろうかという あまり科学的でない趣味の一環。
senyuug.f(GIMP用),
senyuu.f（Paint Shop用）
適当な３色だけで描かれた画像の中の、それぞれの色の部分の面積比
などを求めるプログラム。
例えば木橋の写真で、木部材の橋全体に占める割合や
木橋の全景に占める割合を求める場合は、まず、
橋の部分を切り出して背景白の画像に貼り付け、
木部材の部分を赤一色で塗りつぶし、その他の部材を黒一色で
塗りつぶした画像を ppm のアスキー形式で保存したものを入力する。
構造景観などの解析の際に必要になりそうな画像の合成
（写真から橋だけを切り取って他の背景に貼り付けるとか）
の方法については
この辺。
ppm18.f（04/12/14）
ppm.f をPaint Shop Pro 8
のppm 形式に対応させた例。
ppmzeng.f(GIMP用),
ppmzen.f（Paint Shop用）
ppm18.f の横縞や縦縞を出力する部分を取り除いたもの。
RGB関係の数値だけがほしい場合用。
景観の全景に対して R値、G値、B値それぞれの平均と標準偏差を求める場合など用。
ppmhasig.f(GIMP用),
ppmhasi.f（Paint Shop用）
白い背景に橋の部分の写真を切り取って貼り付けた画像の橋の部分のRGB等を求める。
ppmzen.fの全景に対する結果と合わせて、
全景に対する橋の部分のRGBの比などを求めたいとき用。
hensa.f
これはたぶん、ウェブで集計した景観アンケートとかの点数データを
naga.dみたいな形式に整理して、
hensa<naga.d で平均とか標準偏差とかを出すためのものみたい。
一応、置いておく。
tirabari.f(GIMP用)
二値化またはグレスケール化したppm画像の中の黒い部分の面積比、
散らばり具合を求める。
耐候性鋼材の錆び具合を調べるセロテープ試験をスキャナーに
読み込んで二値化したデータから、セロテープについた錆びの部分の面積比や
錆びの散らばり具合を定量化して、目視による外観評点とどれくらい相関が
あるかを調べてみようかと。
tira.f(GIMP用)
akikoma.sh
100枚以上とかの大量のデータをtirabari.fで1枚ずつ処理するのは大変だから、
まとめて処理する時は、
akikoma.shみたいなスクリプトを
書いておいて、このスクリプトの中で
tira.fを走らせると、
akikoma.csv
のように、ファイル名、測定箇所記号、サビ面積比、散らばり度、正規化Y値、その標準偏差の順に並ぶコンマ区切りデータが得られる。
kaikisen.f
kaikisen.f90(utf-8)
x,yデータを入力して、回帰直線y=ax+bの係数a,bや相関係数、標準偏差などを
求める。
maku.f
耐候性鋼材の膜厚の測定などで、測定箇所ごとに複数のデータがある場合、
それぞれの測定箇所ごとのデータの平均値、標準偏差、変動係数を出力する。
kaikir.f90(utf-8)
荷重ー変位関係(b.csvみたいな)
の初期勾配に対して、
線形回帰したときに最も決定係数の大きくなる範囲を調べる。
toumen.f90(utf-8)
等値面積置換で極限支持力を見積もる用。
塑性域で荷重が落ちていく
荷重ー変位曲線(dpみたいな)
に水平線を引きながら、
荷重ー変位曲線が水平線より上に来る領域と
下に来る領域の面積が近くなる水平線の高さを見積もる。
mxv.f90(マトリクス×ベクトルの掛け)算
mxmx.f90(マトリクス×マトリクスの掛け)算
gausu.f90(代数方程式を解く)
gau.f
gau.f90(gfortran用utf-8)
gautai.f(対称行列用)
まずは、ガウスの消去法ぐらいは自分でも作っておこうかと。
1次方程式AX=BのAとBを与えてXを求める。
hairetu.f90(gfortranで整合配列にするとうまくいかない例)
hairetud.f90(gfortranでもサブルーチンからサブルーチンを呼び出せばうまくいく例)
g77では、メインプログラムでa(100,100)などと大きめの配列を切っておいて、 サブルーチン内では、引数として受け渡されたnを使って、 a(n,n)などと小さめの配列を切り直しても問題なかったが、 gfortranでは、サブルーチン内でもメインプログラムと同じ大きさの 配列を切らないとうまく引き渡せなくなってしまった。 ただし、サブルーチン内からサブルーチンを呼び出す場合は、 特に問題は起きないので、 整合配列を使いたい場合は、一端サブルーチン内に入ってからサブルーチンを 呼び出すようにすればよいのかもしれない。
hensati.f90(平均、標準偏差、偏差値の計算)
yotei20.f90(出力例)
yotei.f90(古いバージョン)
yotei.f(漢字が文字変数として扱えない場合用)
youbi.sed(漢字を後から置換用)
私は長年、
予定表は縦開きの小型スケジュール帳(DAIGO Appoint Diary)を使っている。
メモはいつでも即座に取り出せて即座に書き込めることが重要なので、
この縦開きスケジュール帳をお尻のポケットに、
小型シャープペンを胸ポケットの免許証入れに入れて常に携帯している。
なぜ縦開きかというと、横開きのものは、
お尻のポケットに入れて座っているうちに、
綴じ部分が劣化してボロボロになってしまうので、縦開きであることが
重要なのだ。
まあ、縦開きでもそれなりに劣化してきたり、
トイレなどで落とすことがたまにあるので、
事前に来年用のスケジュール帳を3冊ぐらい買っている。
以前、私が縦開きのスケジュール帳を使っているところを見た学生に
「コロンボみたいですねえ」と言われて、
刑事コロンボの好きな私はますます気を良くして、
縦開きスケジュール帳を使い続けている。
さて、私は最近(2011)、 Andoroidのスマホを使うようになった (2011/3/11から使い始めたので、最初に受信したのは緊急地震速報で、 設定を間違えた誤受信かなぐらいに思ってしまった) のだが、 スケジュール帳を電子ファイルに移行できれば、 検索性や記録性など、様々な利点がある。 と思い、様々なスケジュール帳アプリを試してみたのだが、 どのアプリのインターフェースも 「お尻から手帖を出して、パラパラとめくって確認し、 必要があれば、 胸からシャープペンを取って書き込む」 という実に使いやすいインターフェースと比べると、 あまりに使いにくい。 ただ、日付の書いてある手帳に書き込むということが したいだけなのに、 「新規登録」だの「タイトル」だの 「開始日」「終了日」の選択だの、 メモを書きこもうとするまでに余計なことに頭を使わされて とてもストレスである。 ただ、メモを書くということ以外に頭を使いたくない私にとっては、 テキストエディタでテキストファイルに書き込むのが一番だ。 ただ、日付と曜日はあらかじめ記入されててほしいので、 そのための道具を作ってみた。 春分の日とかの有効期限は2025年あたりまでかな。 祝日だの振替の考え方が間違ってたりするかもしれないが、 こんなふうなテキストファイルが 得られる。 これをスマートフォンのテキストエディタで閲覧・編集していくのだが、 フリック入力になれていない私には、 紙の手帳に書き込むのに比べると、まだまだストレスは高い。 でもまあ、これくらいなら、 適応できるかもしれない。
webhtml.f
webhtml.f90
webhtmlhuna.f90
（船川データベース化用に600dpi画像と300dpi画像を区別してリンク）
(gfortran用)
まだ、調整の余地はあるけども、
ディレクトリ内の画像ファイルのサムネイルを
ImageMagickを利用して生成してから、
ウェブアルバムを作成する。
gThumbが2.12.3当たりにバージョンアップしてから、
各種ウェブツール類との連携やインポート機能みたいな
(私には余計な)機能ばかりが増えて、
最も重宝していたウェブアルバムの作成機能とかが見当たらなくなってしまった
(gThumb 2.14.3では、
ファイル→エクスポート→ウェブアルバム)。
というわけで、自力で対処することにした。
zaseki.f
座席などをランダムに決めたいときなど、
番号をランダムに並び替えて出力。
例えば、学生数が57だったら、1から57をランダムに並び替えて
それに100とかを足して101から157がランダムに並んだファイルを生成。
あとは、
101から157までを57人の学籍番号と置換する置換ファイルの
101から157の列をランダム順のものに取り換えて、
htmlやTeXの座席表に埋め込んだ101から157の文字列を学籍番号に置換するとか。
omomi.f
omomi.f90
(gfortran用utf-8, doループの変数を整数化)
工学資源学研究科では、教員の個人評価書というのを毎年 書かされるのだが、
自己評価した教育領域、研究領域、国際・社会領域、管理・運営領域の
それぞれに対して、例えば50才未満の准教授だったら、
0.3, 0.5, 0.1, 0.1の重みをかけて合計点を求めるのだが、
重みは、合計が1になるようにすれば、
0.2の増減で調節していいことになっている。
それだったら、点数が最大になるように最適に重みを設定した方がいい。
その計算用。
というか、
こういうのは、雛形ファイルをワードで書かせるのはやめにして、
せめて表計算ツールで最適化された重みと点数が
自動計算されるような雛形ファイルを用意してほしい
（せめて、各領域の評価点数を算出するめんどくさい計算式は、
表計算ツールに最初から入れておいてほしい。
だからこういう書式をワープロ書式で書かせるのはつくづくやめてほしい）。
しかも、点数の自己評価は、評価委員会の評価とずれるので、
評価委員会の点数が出てから、改めて最適な重みを設定し直していいことに
してほしい。というか、評価委員会で各人に最適の重みを自動計算してくれるのが
望ましいのだが。
omomi19.f90
2018年度からは、教員活動評価というのを全学で実施するようになり、
教育、研究、診療、社会貢献、産学連携、国際、管理運営といった領域の
重みを多少 増減できるのだが、この重みは、
年度始めの活動計画を出す段階で決めてしまわなければならず、
達成状況の点数を見て最適な重みを設定するということができない。
2019年度も同様なのだが（なぜか年度始めの活動計画を報告する時期が、
ありがたいことに11月までずれ込んだのだが）、
昨年度の審査会評価の点数を見て、
今年度も昨年度とほぼ同様の活動と評価になると想定されるなら、
昨年度の点数が最適化されるように重みを設定すればいいのかもしれない。
ちなみに、私の昨年度の審査会による評価点は3.28だが、
重みを最適化してみたら3.42になった。
結構 変わるのでは？
sadokuirai.f90
委員会の論文集等で、査読システムを使う経済的余裕がない場合に、
メールで査読依頼を出したりする際に、
査読者ごとに担当する論文名とかを流し込んでコピペ用テキストを生成。
kyouyou.csv(教養基礎教育科目コード)
senmon.csv(土木の専門科目コード)
mokuhyou0.csv(学習・教育目標マトリックス)
moku0.csv(上記の科目名なし入力データ用)
kamoku.sed
工学資源学部学部の成績ファイルの科目コードを科目名に変換する。
コマンドプロンプトで、
$ sed -f kamoku.sed seiseki.mae > seiseki.ato
のように実行。
pre05.f
進級判定のための前処理をするプログラム。まず、
工学資源学部学部の成績ファイルをcsv形式でもらってきて、
このうち、
学籍番号、科目コード、単位、修得フラグのみを下記のように
書き出して hoge.csv とかに保存する。
但し、最終行に、学籍番号が9999999の適当なデータを1行加えておく
(最終行を1行コピーして学籍番号を9999999に書き換える)。
7500001,8012345,1,1 7500001,8012346,2,1 7500001,8012347,1,1 7500001,8012348,2,0 7500001,8012349,1,1 7500002,8012345,1,1 7500002,8012346,2,1 7500002,8012347,1,0 7500002,8012348,2,1 7500002,8012349,1,1 7500002,8012350,3,1 9999999,8012350,3,1
pre05<hoge.csv>hoge.d のように実行すると、下記のように 冒頭の1行目に学生数、2行目から、 一人につき 学籍番号を1行目、修得科目数を2行目、 科目コードと単位と修得フラグの組を3行目から修得科目数ぶん、 という具合に並ぶデータに並び替えてくれる。 これでエクセル流のcsvよりはだいぶ処理がしやすくなった。
2 7500001 5 8012345,1,1 8012346,2,1 8012347,1,1 8012348,2,0 8012349,1,1 7500002 6 8012345,1,1 8012346,2,1 8012347,1,0 8012348,2,1 8012349,1,1 8012350,3,1
han34.f, sotu.f
という訳で、
まずは、3年から4年に進学できるかどうかの判定プログラムhan34.f
(プログラム名は実はハンサングンをもじってある)
と、卒業判定プログラムsotu.f。
これが更にポートフォリオの自動生成だの
JABEEツールだのに
発展する可能性は皆無。
例えば3年生の成績ファイルの中から、 3年生に開講されている科目の成績だけを抜き出したい場合、 成績ファイル(csv)の科目コードの右となりにデータの入ってない一列を 挿入する。 そして、 senmon3.sedを使って、
sed -f senmon3.sed seisekimae.csv > seisekiato.csv
のように実行すれば、3年生に開講されている科目コードの右列にだけ 1が入るので、ソートして、それだけ抜き出せば良い。
siitex.f
学習・教育目標マトリックスデータ
moku0.csv
から、達成度チェックシートを自動生成するプログラム。
実行すると、siito.texというTeXファイルを出力するので、
$ sed -f
kamoku.sed
siito.tex > siiton.tex
を実行して科目コードを科目名に変換する
(出力例)。
soukatex01.f(01/平成13年度入学生用)
soukatex02.f(02/平成14年度入学生用)
soukatex03.f(03/平成15年度入学生用)
soukatex03h.f(03/平成15年度入学生編入生用)
soukatu.f(コンマ区切りで出力)
学務の成績データを上記の
pre05.f
で整形したデータを入力して、
総括表と
各学生の達成度チェック一覧表を自動生成。
JABEEでは、一人一人の学生が
学習・教育目標を
達成していることを示さなければならないのだけど、
どうやってそれを示すかというと、まず、学生の一人一人に、
それぞれの学習・教育目標に対応する授業科目を何単位ぶん修得したかという
表(総括表)を書かせて、それを教官が学生一人一人に対してチェックしていく。
この作業は学生にとっても教官にとっても時間と労力を消耗するので、
学務の成績データからなんとか自動的に総括表をTeX形式で作成したいと考えた
(まだ、調整段階。下の
moku0.csv
さえ、適切に作れれば、それなりにいけそうな感じ。
シラバスデータからmoku0.csvが自動生成されるのが理想だけどなあ
)。
moku0.csv
まずは、各科目ごとに
学習・教育目標を
設定したファイル(moku0.csv)を作る。
1行ごとに、
'科目コード', '単位数', '必修だったら20/選択必修だったら15/選択だったら10',
'学習・教育目標のA-1からF-2について、◎だったら2/○だったら1/なしだったら0'
を書き込む。
ちなみに、科目名と科目コードとの対応は、
mokuhyou0.csv参照
(科目名の変換は上記のkamoku.sed参照)。
上記のpre05.fを用いて前処理した成績データ hoge.d を
作成し、./soukatu<hoge.d>hoge.csv のように実行。
一人一人の総括表が順番に並ぶcsvファイルが出力される。
kentiku.f
kentiku90.f
建築デザインコースの科目を修得してない学生の人数確認用。
happ.f
卒論・修論発表会のときに、それぞれの卒論・修論の
新奇性、信頼性、発表のうまさ、
質疑に答えられたかなどについて教官が5段階評価をつけるのだけど、
発表者ごとに評価する教官の数が一定しないようなデータを集計して、
各発表者ごとに評価者数で平均した点数をコンマ区切りデータで出力する
プログラム。データの並びは、
1行目が発表者数、2行目が学籍番号、3行目がその発表者の評価者人数、
4行目以降は、評価者名(半角英数)と各項目の評価点を評価者の数だけ、
で、次の発表者の学籍番号と続く。
こういう入力作業を発表会の最中にノートパソコンで行う訳だけど、
私の場合、表計算ツール(OpenOfficeなりエクセルなり)のセルを
マウスだのタブだので移動しながら打ち込むというのがとても苦手なので、
学籍番号だけが一列に並んだテキストデータを用意しておいて、
それにviで以下のように各採点者の点数を打ち込んでいくならマウスが
不要だしラクかなと。
2 7502401 3 got,3,4,2,1 ogi,5,5,5,5 ham,1,1,1,1 7502402 2 tak,1,2,3,4 tok,4,3,2,1
happy.f
しかし、こういう採点をやるとほぼ必ず、採点項目(新奇性、信頼性、発表、質疑)の
一部を空欄にする人が現れるものだろう。
そのとき、一部が空欄の採点者の評価はすべて無効にするというやり方が簡単では
あるけれど、空欄の部分を平均化の人数から差し引いて、
採点項目ごとに違う分母で平均を取るという方法もあり得る。
という訳で、採点項目が空欄だった場合は、そこに8より大きい数(例えば9)と
書き込むと、その項目は点数に加算されず、そのぶん平均化の採点者数
にも加算されないようにした。
happa.f
ABCDの4段階評価をする05年度用。
データは、
2 7502422 3 got,a,z,b,c ogi,b,a,b,c oik,c,z,a,z 7502421 2 tak,z,c,a,a tok,b,a,a,c
みたいに。空欄のところはzを入れておくと平均化の分母の 評価者数には入れない。
happc.f
セッションごとに採点者の人数がバラバラの場合、
採点者の人数をいちいち数えて入れるのはめんどくさいので、
採点者名に「zzz」という行が現れるまで、各採点者の評点を読み込むことにした。
なお、学籍番号に加えてEUCコードで氏名を添える場合は、""でくくること
(学籍番号,"氏名"とならんでるテキストデータは表計算ソフトで
簡単に吐き出せると思うが)。
出力では、トラブルを避けるため氏名は一番最後の列に出す。
2 7502422,"上田次郎" got,a,z,b,c ogi,b,a,b,c oik,c,z,a,z zzz,z,z,z,z 7502421,"山田奈緒子" tak,z,c,a,a tok,b,a,a,c zzz,z,z,z,z
長さを変えた15個の振り子を一斉に揺らす実験（動画）
を見てなかなか面白いと思い、
これを(ちょっと今やらされてる高大接続テキスト云々でも
使えないかなと思って)、
実際に作って学生に実験させたりする教材にできないかなあと思い立った。
既に、
前野昌弘さんが、
ActionScriptで作成した3次元Flash動画を公開しており、
私なんかが出る幕ではないのだが、
具体的に振り子1本1本の長さを何cmずつにすればいいのか、
長さの設定具合によって球の動きがそろうまでどれくらいの時間がかかるのか
ということが知りたかった。
というか、もしかして、
そのなんとかテキストに、
「手計算で調べようとしてもなかなかよくわからない現象を可視化したりして
調べるのに、数値シミュレーションという手法がある」
といったことを書いてみてもいいかもしれないかなと思って、
高校生でも数学Bの授業でやってる
BASICで、
15個の振り子の簡単なシミュレーションプログラムを作ってみようかなと。
ちなみに、15個の振り子の動きが周期的にそろうようにするには、
角振動数が一定の増分で変化するように振り子の長さを変える。
振り子の長さを1cmずつ短くするとかでは、動きはそろわない。
だから、その角振動数の増分を一定にしたときの具体的な長さも出力したい。
hurikoc.bas
Linuxターミナル内での実行に適した
Chipmunk Basicでの
利用を想定している。
実行例は、
こんな感じ。
huriko10.bas(Linux用utf-8コード)
huriko10sj.bas(Windows用sjisコード)
高校の数学Bの教科書でも使われている
(仮称)十進BASICで、
Microsoft互換モードのグラフィック機能
(circleで丸を描ける)を利用して、15個の振り子の動きを可視化表示できるようにした
(Microsoft互換モードはLinux版でも使える)。
高校の教科書ではグラフィック機能は使ってないと思うけど。
実行画面例。
Vine Linux 5.2だと、編集画面で日本語を入力することはできない。
ただ、viなど他のエディターでUTF-8の日本語を入力して、
print文で日本語を画面表示するぶんには、問題なく表示する。
input文のキーボード入力で、文字が2つずつ入力されたりする不具合がある。
ウィンドウマネージャーをGnomeではなくtwmとかにすると、その不具合はなくなる。
15個の振り子の長さのいくつかの計算例 59.07 秒周期 1 cos(ωt)= 0.999473 L= 60 最長 2 cos(ωt)= 0.999336 L= 56.963977 3 cos(ωt)= 0.999183 L= 54.152692 4 cos(ωt)= 0.999014 L= 51.5445 5 cos(ωt)= 0.998829 L= 49.120298 6 cos(ωt)= 0.998629 L= 46.863179 7 cos(ωt)= 0.998413 L= 44.758135 8 cos(ωt)= 0.998181 L= 42.791804 9 cos(ωt)= 0.997933 L= 40.952261 10 cos(ωt)= 0.99767 L= 39.228837 11 cos(ωt)= 0.99739 L= 37.611958 12 cos(ωt)= 0.997095 L= 36.093021 13 cos(ωt)= 0.996784 L= 34.664272 14 cos(ωt)= 0.996458 L= 33.31871 15 cos(ωt)= 0.996115 L= 32.05 最短 43.4 秒周期 1 cos(ωt)= -0.994538 L= 72 最長 2 cos(ωt)= -0.994315 L= 66.666973 3 cos(ωt)= -0.994088 L= 61.90525 4 cos(ωt)= -0.993856 L= 57.63605 5 cos(ωt)= -0.99362 L= 53.793717 6 cos(ωt)= -0.993379 L= 50.323185 7 cos(ωt)= -0.993134 L= 47.177986 8 cos(ωt)= -0.992885 L= 44.318692 9 cos(ωt)= -0.992631 L= 41.711673 10 cos(ωt)= -0.992372 L= 39.3281 11 cos(ωt)= -0.99211 L= 37.14315 12 cos(ωt)= -0.991842 L= 35.135351 13 cos(ωt)= -0.991571 L= 33.286059 14 cos(ωt)= -0.991295 L= 31.579016 15 cos(ωt)= -0.991014 L= 30 最短 49.1秒周期 1 cos(ωt)= -0.591728 L= 80 最長 2 cos(ωt)= -0.591623 L= 74.456055 3 cos(ωt)= -0.591519 L= 69.469047 4 cos(ωt)= -0.591414 L= 64.966799 5 cos(ωt)= -0.591309 L= 60.888461 6 cos(ωt)= -0.591205 L= 57.182437 7 cos(ωt)= -0.5911 L= 53.804752 8 cos(ωt)= -0.590996 L= 50.717735 9 cos(ωt)= -0.590891 L= 47.888965 10 cos(ωt)= -0.590786 L= 45.290422 11 cos(ωt)= -0.590682 L= 42.897782 12 cos(ωt)= -0.590577 L= 40.689854 13 cos(ωt)= -0.590473 L= 38.648102 14 cos(ωt)= -0.590368 L= 36.756258 15 cos(ωt)= -0.590263 L= 35 最短 59.35秒周期 1 cos(ωt)= 0.99923 L= 120 最長 2 cos(ωt)= 0.999175 L= 111.581531 3 cos(ωt)= 0.999118 L= 104.018842 4 cos(ωt)= 0.99906 L= 97.199753 5 cos(ωt)= 0.998999 L= 91.029878 6 cos(ωt)= 0.998937 L= 85.429348 7 cos(ωt)= 0.998872 L= 80.330207 8 cos(ωt)= 0.998806 L= 75.674344 9 cos(ωt)= 0.998738 L= 71.41183 10 cos(ωt)= 0.998668 L= 67.499571 11 cos(ωt)= 0.998596 L= 63.900218 12 cos(ωt)= 0.998523 L= 60.581271 13 cos(ωt)= 0.998447 L= 57.514343 14 cos(ωt)= 0.99837 L= 54.674548 15 cos(ωt)= 0.99829 L= 52.04 最短
ホームセンターで売っているキャスター付きスチール棚で作る際のこつ：
振り子の長さを随時 微調節できるように、上端の固定部の一方は
クリップで固定位置を動かせるようにしておく。
クリップは揺らしているうちに少しずつずれてくるので、
調整が終了したら、
養生テープなどでひもを固定する。
ひもの上端をスチールの棒に固定する際、
振り子の根元がスチール棒の下部から回転するように、
数センチの短いひもでスチール棒の下部に結び目を作って固定する。
最長振り子と最短振り子の差をなるべく大きくとった方が、 誤差の影響は少なくなるが、振り子が長くなりすぎると ねじれやすくなったりの問題が生じる。 V字につるすのは、隣同士ぶつかったり、回転運動にならないようにするためだが、 長い方の振り子のV字は、幅を広くするといった工夫も有効かもしれない。
2012/7/27: docomo携帯の有料番組?の BeeTV「おやすみチャンネル」 （を制作している千代田ラフト）の 取材を受けた。 「お母さんの子守歌をママさん合唱団が歌う編」「世界的に有名な物理学者による超絶難解授業編」など、眠気を誘う映像を冗談感覚で紹介する番組に、 15個の振り子を撮らせてほしいとのこと。 「振り子の映像に誘眠やリラックスの効果がある」といった 科学的に確認されていない主張にお墨付きを与えることに利用したりはしないことに 念を押した上で了承したものの、有料番組なので、実際のところはどうかわからない。 というわけで、 私自身は 「振り子の映像に誘眠やリラックスの効果がある」といった主張は 一切していないことを念のためここで明言しておく。 なお、15個の振り子の方はスチール棚の柱がじゃまだったので、 12個の振り子の方を主に撮影していたようだった。
hurikomml.bas
以下はmidi形式音楽ファイル
hurikoc.mid(ハ長調音階, 25cm-12.5cm)
huriko3.mid(3度音階, 25cm-12.5cm)
hurikoo.mid(沖縄音階, 25cm-12.5cm)
hurikoj.mid(CDEbEGABbCペンタトニック系音階, 25cm-12.5cm)
huriko2.mid(長いけど振り子の動きをそれなりに追えるハ長調音階, 40cm-20cmだったかな?)
huriko.mid(最初の試作。短く収まるように時間間隔をはしょりすぎたら独特の旋律的?音階になった, 40cm-20cmだったかな?)
これはまったくのおまけというか遊びだが、
振り子の運動の
計算結果を眺めていたら、
例によって私にはこの独特の周期性が、
どうもミニマルミュージックみたいだなあと思えてきた。
で、各振り子の振幅をドレミファソラシドレミファソラシドの
15段階の音階にして
(ちなみに半音間隔もやったけど、これは調性感が失われるので
私の好みではなく却下)、
一つの振り子に一つずつの音色を与えて音楽にしてみた。
手順は、
hurikoc.basを修正した
hurikomml.basで
こんな感じのMML形式ファイルを吐き出して
微修正してちゃんとした形のMML形式にし、
mml2midを使ってmidi化。
なんかメシアン風というか、
ペトルシュカ風というか、まあまあ面白い。
まあ、ここには出さないが後は趣味の範疇で色々やってみたい。
が、こういうところから物理やプログラミングに興味を持つ人もいるかもしれないので、
少し置いておく。
梅原さんのプログラムここから
2017年3月修了予定の梅原さん（コンクリ研）のプログラムに、
リンクを張っておきます。文字コードはsjisかな。
ftoughness.f90
このプログラムは角柱の曲げ試験における応力と変位で囲まれた部分の面積を求め、
そこから曲げタフネス係数を求めるものであり、ファイルは2列の数字のみのファイル
(CSV形式)に対応しており、1列目に荷重(kN)、2列目に変位(mm)となっていること。
また、荷重は正の値、変位に関しては負の値であること。
プログラムに従い、ファイル名、試験時のスパンや供試体の寸法を入力すると応力-変位曲線
のスパン1/150の変位までの囲まれた面積、曲げタフネス係数を算出し結果を出力する。
ra0.f90
このプログラムはある一定のピッチにて計測した粗さのデータから
算術平均粗さRaを算出するものであり、ファイルは1列の数字のみの
ファイル(CSV形式)に対応している。
プログラムに従い算術平均粗さを求めたいファイル名を入力すると
計算し、結果を出力する。
staticmodulus.f90
このプログラムでは2つのひずみゲージを用い測定したコンクリートの静弾性係数を
求めるものであり、ファイルは3列の数字のみのファイル(CSV形式)に対応しており、
1列目は荷重(MN)、2，3列目にひずみ(μ)であること。また圧縮は負の値を、ひずみは
引張が正の値、圧縮が負の値であること。
プログラムに従い、ファイル名、供試体の寸法等を入力すると、静弾性係数を算出し
結果を出力する。
梅原さんのプログラムここまで
ここから下に紹介しているのは他人の頁への「リンク」であって、 私が作ったものではありません （ですので、「はじめに」に書いたことは、 ここから下には適用されません（念のため）。
FORTRANプログラミング初級編（秋元博路さん）
Ｆｏｒｔｒａｎ９０プログラミング
岩熊先生による
構造力学講義ノート案と有限要素法プログラム
ホクト・システム-ちょっといい話（有限要素法・よもやま話）
『もっと知りたいコンクリート講座』（
武蔵工業大学都市基盤工学科構造材料工学研究室、
だれでも質問コーナー）
技術情報Q&A掲示板（環境海洋技術官ネットワーク）
Civil Engineering Excercise III（「コンピューターを微分方程式を解くための道具として使う方法」
斉木 功さん）
有限要素法の基礎
（Shinsuke SAKAIさん）
構造解析特論のページ
(岩崎英治さん。
アイソパラメトリック要素の数値積分や積分点あたりの話も詳しく解説されていて有用)
有限要素法プログラミング演習（fortran, C のソースコードつき）（開発版）
(渡辺浩志さん。
非線形の有限要素解析プログラムを実際にゼロから作ろうとしている人にとって、
とても有用)
The Internet-College of FEM
（福森さん）
ＣＡＥ大事典
（ＮＣネットワーク）
検索用鍵語：
copyleft,
立体要素, ソリッド要素, 6面体要素, アイソパラメトリック要素ではないけど,
集成材構造, ハイブリット集成材梁, 鋼板挿入集成材梁,
数値シミュレーション, 数値シュミレーション,
数値計算, 数値解析, 有限要素モデル, 梁要素モデル, 立体要素モデル,
木構造, ハイブリット梁, ハイブリット構造
