Microsoft PowerPoint - 3D-FEM-1.ppt [互換モード]

Size: px
Start display at page:

Download "Microsoft PowerPoint - 3D-FEM-1.ppt [互換モード]"

Transcription

1 三次元弾性解析コード (/3) プログラムの概要 22 年夏学期 中島研吾 科学技術計算 Ⅰ(82-27) コンピュータ科学特別講義 Ⅰ(8-2)

2 FEM3D-Part 2 対象とする問題 弾性体 Z ヤング率 E ポアソン比 ν U Z Z Z Y Y 直方体 一辺長さの立方体 ( 六面体 ) 要素 各方向に Y Z 個 境界条件 対称条件 U U 強制変位 U ma move

3 FEM3D-Part 3 構成式 : 応力 ~ひずみ関係 ヤング率 E 応力とひずみは比例 比例定数をヤング率 Eとする ( 各物質に固有の値 ) σ σ Eε ε E ポアソン比 ν 方向に荷重をかけると 横方向 (YZ) (YZ) にも変形 縮み割合をポアソン比 ν とする 各物質に固有の値 ε 金属では.3 程度 水 :.5 ゴム : ほぼ.5 非圧縮 νε σ ν E ε νε ε σ

4 FEM3D-Part YZν.3 とすると Z U Z Z ε ε. ε νε.3 Z u u ε. 3 Y Y Y Y

5 FEM3D-Part 5 有限要素法の処理 支配方程式 ガラーキン法 : 弱形式 要素単位の積分 要素マトリクス生成 全体マトリクス生成 境界条件適用 連立一次方程式

6 FEM3D-Part 6 初期化 有限要素法の処理 : プログラム 制御変数読み込み 座標読み込み 要素生成 (: 節点数 ICELO: 要素数 ) 配列初期化 ( 全体マトリクス 要素マトリクス ) 要素 全体マトリクスマッピング (IndeItem) マトリクス生成 要素単位の処理 (do cel ICELO) 要素マトリクス計算 全体マトリクスへの重ね合わせ 境界条件の処理 連立一次方程式 共役勾配法 (CG) 応力計算

7 FEM3D-Part 7 三次元要素の定式化 三次元弾性力学方程式 ガラーキン法 要素マトリクス生成 宿題 プログラムの実行 データ構造 プログラムの構成 計算効率について Paraew による可視化

8 FEM3D-Part 8 二次元への拡張 : 三角形要素 任意の形状を扱うことができる 特に一次要素は精度が悪く 一部の問題を除いてあまり使用されない

9 FEM3D-Part 9 二次元への拡張 : 四角形要素 一次元要素と同じ形状関数を 軸に適用することによって 四角形要素の定式化は可能である 三角形と比較して特に低次要素の精度はよい しかしながら 各辺が座標軸に平行な長方形でなければならない 差分法と変わらない 3 このような形状を扱うことができない 2

10 FEM3D-Part アイソパラメトリック要素 (/3) 各要素を 自然座標系 (η) の正方形要素 [± ±] に変換する 3 2 η 各要素の全体座標系 (global coordnate)() における座標成分を 自然座標系における形状関数 [] ( 従属変数の内挿に使うのと同じ []) を使用して変換する場合 このような要素をアイソパラメトリック要素 (soparametrc element) という

11 FEM3D-Part アイソパラメトリック要素 (2/3) 3 η 各節点の座標 :( ) ( 2 2 ) ( 3 3 ) ( ) 各節点における 方向変位 :u u 2 u 3 u u ( η) u η) ( ( η)

12 FEM3D-Part 2 アイソパラメトリック要素 (3/3) η - - 高次の補間関数を使えば 曲線 曲面も扱うことが可能となる そういう意味で 自然座標系 と呼んでいる

13 FEM3D-Part 3 2D 自然座標系の形状関数 (/2) 自然座標系における正方形上の内挿多項式は下式で与えられる : u α α 2 α3η α η 各節点での条件より : η u u2 u3 u α α 2 u u2 u3 u α3 α u u u2 u u2 u 3 3 u u

14 FEM3D-Part 2D 自然座標系の形状関数 (2/2) 元の式に代入して u について整理すると以下のようになる : η 3 u u u u u u 形状関数 は以下のようになる : ( η) ( )( η) 2( η) ( )( η) 3 ( η) ( )( η) ( η) ( )( η) 双一次 (b-lnear) 要素とも呼ばれる 各節点における の値を計算してみよ Y 方向変位 v についても同様 v v 2v2 3v3 v

15 FEM3D-Part 5 三次元への拡張 四面体要素 : 二次元における三角形要素 任意の形状を扱うことができる 特に一次要素は精度が悪く 一部の問題を除いてあまり使用されない 高次の四面体要素は広く使用されている 本講義では低次六面体要素 ( アイソパラメトリック要素 ) を使用する tr-lnear 変位法 自由度 : 変位 各節点上で 3 成分 (uvw) が定義される

16 FEM3D-Part 6 3D 自然座標系の形状関数 ( η ζ ) 8 2 ( η ζ ) η ζ 8 3 ( η ζ ) η ζ 8 ( η ζ ) η ζ 8 8 u ( η ζ ) u 8 v ( η ζ ) v 8 ( )( η)( ζ ) w ( η ζ ) w ( )( )( ) ( )( )( ) ( )( )( ) 5 ( η ζ ) 8 6 ( η ζ ) η ζ 8 7 ( η ζ ) η ζ 8 8 ( η ζ ) η ζ 8 ( ) 8 7 ( ) 5 6 ( ) ( η ζ ) ( ) ( ) 3 2 ( )( η)( ζ ) ( )( )( ) ( )( )( ) ( )( )( ) ( ) ( ) ( )

17 FEM3D-Part 7 三次元要素の定式化 三次元弾性力学方程式 ガラーキン法 要素マトリクス生成 宿題 プログラムの実行 データ構造 プログラムの構成 計算効率について Paraew による可視化

18 FEM3D-Part 8 弾性力学の支配方程式 つりあい式 ひずみ ~ 変位関係式 ひずみ ~ 応力関係式

19 FEM3D-Part 9 三次元のつりあい式な成応力の独立な成分は 6 つ τ τ { τ σ τ τ τ σ σ τ τ τ τ σ τ τ τ τ τ τ σ Y τ σ τ Z σ τ τ Z

20 FEM3D-Part 2 垂直ひずみ~ 変位の関係 PQ P Q ε u d u d ( u ) d d u d R P u / P d R v / Q Q ε ε ε u u v w

21 FEM3D-Part 2 せん断ひずみ~ 変位の関係 R R γ u v d u / P v / Q γ γ v w w u P d Q

22 FEM3D-Part 22 応力 ひずみ関係ず関係 σ σ ν ν ν ν ε ε ( ) E τ σ σ ν ν ν ν ν γ ε ε 2 ( ) ( ) ( ) E τ τ ν ν γ γ ( ) τ ν γ 2

23 FEM3D-Part 23 ひずみ 応力関係ず関係 ν ν ν ε ε ν ν ν ν ν ν σ σ ( )( ) ( ) ( ) E γ ε ν ν ν ν τ σ ( ) ( ) γ γ ν ν τ τ [ ] D { [ ]{ ε σ D 非圧縮性材料 (ν~.5) の場合 特別な扱い必要

24 FEM3D-Part 2 - 方向のつりあい式に注目 σ τ τ σ τ τ ガラーキン法 [ ] { σ τ τ d [ ] { [ ] { [ ] { σ d τ d τ d [ ] { d

25 FEM3D-Part 25 - 方向のつりあい式 (/3) [ ] { σ [ ] { d σ d [ ] { σ Q [ ] { ] σ d [ ] { σ S n S ds n ds [ ] { ] [ ] { [ ] σ d σ d { σ 同様にして以下の弱形式が得られる ( 表面積分省略 ): 一次元ガウスの公式 d 部分積分の公式 [ ] { [ ] { [ ] { σ d τ d τ d [ ] { d [ ] { [ ] { [ ] σ d τ d { τ d [ ] { d

26 FEM3D-Part 26 - 方向のつりあい式 (2/3) [ ] { [ ] { [ ] d τ d { τ d [ ] { d (*) σ σ E ν ε νε ν 2ν νε ( )( ) ( ) [ ] E ν u νv w ν 2ν ν ( )( ) ( ) [ ] τ τ E 2 γ E [ ] u v ( ν ) 2( ν ) E E 2 γ [ ] u w ( ν ) 2( ν )

27 FEM3D-Part 27 要素内の変位分布形状関数要素内の変位分布 形状関数 [ ] [ ] [ ] { { { W w v U u [ ] [ ] [ ] { { { W w v U u [ ] [ ] [ ] { { { U u U u U u [ ] [ ] [ ] [ ] [ ] { { { { { { { v v U u U u U u [ ] [ ] { { W w W w E ( )( ) ( )[ ] [ ] [ ] ( ) { { { 2 E E W U E ν ν ν ν ν σ ( ) ( ) [ ] [ ] ( ) { { 2 2 U E E ν γ ν τ ( ) ( ) [ ] [ ] ( ) { { 2 2 W U E E ν γ ν τ

28 FEM3D-Part 28 - 方向のつりあい式 (3/3) 方向のつりあい式 (3/3) ( ) ( )( ) ( )( ) ( ) ν ν E b E a E D とすると [ ] [ ] [ ] { { { W a a U D σ ( )( ) ( )( ) ( ) ν ν ν ν ν [ ] [ ] [ ] [ ] [ ] [ ] [ ] { { { { { { { W b U b b U b τ τ [ ] [ ] { { W b U b τ () [ ] { [ ] { [ ] { [ ] { * d d d d τ τ σ [ ] [ ] [ ] [ ] [ ] [ ] ( ) { { U d b D [ ] [ ] [ ] [ ] ( ) { { [ ] [ ] [ ] [ ] { { W d b a d b a [ ] { d

29 FEM3D-Part 29 Y- 方向のつりあい式 [ ][ ] [ ] [ ] [ ] [ ] ( ) { { d b D [ ] [ ] [ ] [ ] ( ) { { [ ] [ ] [ ] [ ] { { [ ] { d Y W d b a U d b a [ ] { d Y Z- 方向のつりあい式 [ ][ ] [ ] [ ] [ ] [ ] ( ) { { d b U d b W d b D [ ] [ ] [ ] [ ] ( ) { { [ ] [ ] [ ] [ ] { { [ ] { d Z d b a U d b a [ ] { d Z

30 FEM3D-Part 3 3 つのつりあい式 変位の 3 成分を未知数 33 つの方程式 3つの方程式はカップルしている ( 独立では無い ) 形状関数ベクトル []:8 行列 [] [] [ ] [ ] 等 :8 8 行列 σ τ τ τ τ σ τ τ σ Y Z

31 FEM3D-Part 3 要素マトリクス :2 2 行列 各節点上の (uvw) 成分が物理的にも強くカップルしているので3 自由度をまとめて扱う :8 8 行列 j ( ) 8 7 ( ) ( ) 5 6 ) ( ) 3 ( ) 2 ( η ζ ) ( ) ( ) j j 2 j3 j2 j 22 j 32 j3 j 23 j 33 ( j K8)

32 FEM3D-Part 32 要素マトリクス :2 2 行列 各節点上の (uvw) 成分を独立に扱うことも可能であるが 8 8とする利点については来週以降も説明する j ( ) 8 7 ( ) ( ) 5 6 ) ( ) 3 ( ) 2 ( η ζ ) ( ) ( ) j j 2 j3 j2 j 22 j 32 j3 j 23 j 33 ( j K8)

33 FEM3D-Part 33 要素マトリクス :-j j 成分 方向 (/3) j j2 j j 2 j3 j 22 j 32 j 23 j 33 ( j K 8 ) { D[ ][ ] b( [ ] [ ] [ ] [ ]) d{ U { a [ ] [ ] b ( [ ] [ ] ) d { { a[ ] [ ] b[ ] [ ] d{ W j ( ) ( ) 3 ( ) ( ) 5 6 ( ) 2 ( η ζ ) ( ) { D ( ) j b j j { a j b j j 2 { a j b j j3 d d d ( ) ( )

34 FEM3D-Part 3 要素マトリクス :-j j 成分 Y 方向 (2/3) j j 2 j3 j2 j 22 j 32 j3 j 23 j 33 ( j K 8 ) { D[ ][ ] b( [ ] [ ] [ ] [ ]) d{ { a [ ] [ ] b ( [ ] [ ] ) d { U { a[ ] [ ] b[ ] [ ] d{ W { a j b j j 2 j 22 ( ) 8 ( ) 7 3 ( ) ( ) 5 6 ( ) 2 ( η ζ ) ( ) d { D ( ) j b j j { a j b j j 23 d d ( ) ( )

35 FEM3D-Part 35 要素マトリクス :-j j 成分 Z 方向 (3/3) j j 2 j3 j2 j 22 j 32 j3 j 23 j 33 ( j K 8 ) { D[ ][ ] b( [ ] [ ] [ ] [ ]) d{ W { a [ ] [ ] b ( [ ] [ ] ) d { U { a[ ] [ ] b[ ] [ ] d{ { a j b j j 3 { a j b j j 32 j33 ( ) 8 ( ) 7 3 ( ) ( ) 5 6 ( ) 2 ( η ζ ) ( ) d d { D ( ) j b j j d ( ) ( )

36 FEM3D-Part 36 あとは積分を求めれば良いあ積分を求良 自然座標系 (ηζ) 上で定義 ガウス積分公式が ηζ 使える ( 三次元 ) しかし 微分が M L d d d f I ) ( ζ η ζ η [ ] j j M j L f W W W ) ( ζ η ) ( ζ LM: ηζ 方向の積分点数積分点座標値 j W W W ) ( j ζ η : 積分点での重み係数 : 積分点の座標値

37 FEM3D-Part 37 自然座標系における偏微分 (/) 自然座標系における偏微分 (/) 偏微分の公式より以下のようになる : ζ η ) ( η η η η ζ η ) ( ζ ζ ζ ζ ζ η η η η η ) ( ζ ζ ζ ζ は定義より簡単に求められるが ζ η は定義より簡単に求められるが を実際の計算で使用する

38 FEM3D-Part 38 自然座標系における偏微分 (2/) 自然座標系における偏微分 (2/) マトリックス表示すると : [ ] J [ ] ζ ζ ζ η η η ζ η ζ ζ ζ ζ [ ] J : ヤコビのマトリクス [ ] J : ヤコビのマトリクス (Jacob matr Jacoban)

39 FEM3D-Part 39 自然座標系における偏微分 (3/) 自然座標系における偏微分 (3/) の定義より簡単に求められる J J J J J η η η η η η J η η η J J ζ ζ ζ ζ ζ ζ J ζ ζ ζ

40 FEM3D-Part 自然座標系における偏微分 (/) 自然座標系における偏微分 (/) 従って下記のように偏微分を計算できる ヤコビアン (3 3 行列 ) の逆行列を求める [ ] η η η η η J ζ ζ ζ ζ ζ

41 FEM3D-Part 要素単位での積分 j j 2 j3 j2 j 22 j 32 j3 j 23 j 33 ( j K8) ( ) 8 7 ( ) 3 ( ) ( ) 5 6 ( ) 2 ( η ζ ) ( ) ( ) ( ) { ( ) j D j b j j d j j j D b d

42 FEM3D-Part 2 自然座標系での積分 j j 2 j3 j2 j 22 j 32 j3 j 23 j 33 ( j K8) ( ) 8 7 ( ) 3 ( ) ( ) 5 6 ( ) 2 ( η ζ ) ( ) ( ) ( ) j j j D b d j j j D b ddd D j j j b dt det J d dηdζζ

43 FEM3D-Part 3 ガウスの積分公式 j j 2 j3 j2 j 22 j 32 j3 j 23 j 33 ( j K8) ( ) 8 7 ( ) 3 ( ) ( ) 5 6 ( ) 2 ( η ζ ) ( ) ( ) ( ) D j b j j det J ddηdζ I L M j f ( η ζ ) ddηdζ [ W ] W j W f ( η j ζ )

44 FEM3D-Part 残りの手順 ここまでで 要素ごとの積分が可能となる あとは : 全体マトリクスへの重ね合わせ 境界条件処理 連立一次方程式を解く 詳細は来週以降の講義でプログラムの内容を解説 詳細は来週以降の講義でプログラムの内容を解説しながら説明する

45 FEM3D-Part 5 要素 全体マトリクス重ね合わせ要素全体トリク重ね合わ (2) (2) ) ( ) ( ) ( ) ( (2) (2) (2) { ]{ [ f φ (2) (2) 3 (2) 2 (2) (2) (2) 3 (2) 2 (2) (2) (2) (2) (2) (2) 3 (2) 33 (2) 32 (2) 3 (2) 2 (2) 23 (2) 22 (2) 2 (2) (2) 3 (2) 2 (2) f f f φ φ φ () () () () () () () () () { ]{ [ f f φ φ (2) (2) (2) (2) 3 (2) 2 (2) f φ () () 3 () 2 () () 3 () 2 () () 3 () 2 () () 3 () 33 () 32 () 3 () 2 () 23 () 22 () f f f f φ φ φ φ f φ Φ Φ { ]{ [ B D F K Φ Φ Φ Φ B B B B D D D D Φ Φ Φ B B B D D D

46 FEM3D-Part 6 三次元要素の定式化 三次元弾性力学方程式 ガラーキン法 要素マトリクス生成 宿題 プログラムの実行 データ構造 プログラムの構成 計算効率について Paraew による可視化

47 FEM3D-Part 7 宿題 ガウスの積分公式を使用して以下の四角形の面積を求めよ ( プログラムを作って 計算機で計算してください ) 3 2 :(..) 2:(. 2.) 3:(3. 5.) :(2..) I d det J ddζ ζ

48 FEM3D-Part 8 ヒント (/2) 座標値によってヤコビアン ( ヤコビの行列 ) を計算 ガウスの積分公式 (n2) に代入する I m n f j η ) d d η [ W ] W j f ( j ) ( η mplct REAL*8 (A-HO-Z) real*8 W(2) real*8 POI(2) W().d W(2).d POI() d POI(2) d SUM.d do jp 2 do p 2 FC F(POI(p)POI(jp)) SUM SUM W (p)*w (jp)*fc enddo enddo

49 FEM3D-Part 9 ヒント (2/2) [ ] η η η η J J det η η η η η η η η ( )( ) ( )( ) η η η η ) ( ) ( 2 ( )( ) ( )( ) η η η η ) ( ) ( 3

50 FEM3D-Part 5 三次元要素の定式化 三次元弾性力学方程式 ガラーキン法 要素マトリクス生成 宿題 プログラムの実行 データ構造 プログラムの構成 計算効率について Paraew による可視化

51 FEM3D-Part 5 対象とする問題 弾性体 Z ヤング率 E ポアソン比 ν U Z Z Z Y Y 直方体 一辺長さの立方体 ( 六面体 ) 要素 各方向に Y Z 個 境界条件 対称条件 U U 強制変位 U ma move

52 FEM3D-Part 52 YZν.3 とすると Z U Z Z ε ε. ε νε.3 Z u u ε. 3 Y Y Y Y

53 FEM3D-Part 53 ファイルコピー インストール (/2) 三次元弾性解析コード >$ cd <$fem> >$ cp /home3/sengon/documents/class/fem/fem3d.tar. >$ tar vf fem3d.tar >$ cd fem3d >$ ls c f run FEM インストール (C) >$ cd <$fem>/fem3d/c >$ mae >$ ls../run/sol sol FEM インストール (FORRA) >$ cd <$fem>/fem3d/f >$ mae >$ ls../run/sol sol

54 FEM3D-Part 5 ファイルコピー インストール (2/2) メッシュジェネレータインストール >$ cd <$fem>/fem3d/run >$ g95 O3 mgcube.f o mgcube

55 FEM3D-Part 55 計算の流れメッシュ生成 計算 ファイル名称固定 mgcube メッシュジェネレータ cube. メッシュファイル sol FEM 本体 IPU.DA 制御データ test.np 可視化用出力変位 3 成分 σ

56 FEM3D-Part 56 メッシュ生成 Z U Z >$ cd <$fem>/fem3d/run / >$./mgcube Y Z 各辺長さを訊いてくる このように 入れてみる Z >$ ls cube. 生成を確認 Y Y cube.

57 FEM3D-Part 57 制御ファイル :IPU.DA IPU.DA cube. fname MEHOD PRECOD terprema( 不使用 ) 2 IER..3 ELAS POISSO fname: メッシュファイル名 MEHOD: 反復解法 ( に固定 ) PRECOD: 前処理手法 :Bloc-LU-GS:Bloc 対角スケーリング terprema: ( 不使用 ) IER: ELAS: 反復回数上限 ヤング率 POISSO: ポアソン比 あとで.999などの場合を試してみよ

58 FEM3D-Part 58 >$ cd <$fem>/fem3d/run >$./sol ( 反復法の収束履歴 ) E E E-9 実行 ### DISPLACEME at (maymazma)) 33-3.E- -3.2E-.E >$ ls test.np 生成を確認 test.np Z U Z この点での変位が表示されている 33 番節点 ( 3 ) Z Y Y

59 FEM3D-Part 59 三次元要素の定式化 三次元弾性力学方程式 ガラーキン法 要素マトリクス生成 宿題 プログラムの実行 データ構造 プログラムの構成 計算効率について Paraew による可視化

60 FEM3D-Part 6 メッシュファイル構成 :cube. 番号は から始まっている 節点データ 節点数 節点番号 座標 要素データ 要素数 要素タイプ 要素番号 材料番号 コネクティビティ 節点グループデータ グループ数 グループ内節点数 グループ名 グループ内節点

61 FEM3D-Part 6 メッシュ生成コード :mgcube.f(/5) mplct REAL*8 (A-HO-Z) real(nd8) dmenson(::) allocatable :: Y real(nd8) dmenson(::) allocatable :: Y character(len8) :: GRIDFILE HHH nteger dmenson(::) allocatable :: IW!C!C !C II.!C !C wrte (**) 'YZ' read (**) YZ P YP Y ZP Z D.d IODO P*YP*ZP ICELO *Y *Z IBODO P*YP allocate (IW(IODO)) IW 方向節点数 Y 方向節点数 Z 方向節点数 総節点数総要素数 Y 平面の節点数

62 FEM3D-Part 62 メッシュ生成コード :mgcube.f(2/5)!c cou b do ZP do j YP cou cou (-)*IBODO (j-)*p IW(coub) enddo enddo cou b 2 do ZP j do P cou cou (-)*IBODO (j-)*p IW(coub) enddo enddo cou b 3 do j YP do P cou cou (-)*IBODO (j-)*p IW(coub) enddo enddo cou b ZP do j YP do P cou cou (-)*IBODO (j-)*p IW(coub) enddo enddo mn の節点を IW(b) に格納 (byp*zp) Z Y U Z Z Y

63 FEM3D-Part 63 メッシュ生成コード :mgcube.f(2/5)!c cou b do ZP do j YP cou cou (-)*IBODO (j-)*p IW(coub) enddo enddo cou b 2 do ZP j do P cou cou (-)*IBODO (j-)*p IW(coub) enddo enddo cou b 3 do j YP do P cou cou (-)*IBODO (j-)*p IW(coub) enddo enddo cou b ZP do j YP do P cou cou (-)*IBODO (j-)*p IW(coub) enddo enddo YYmnの節点をIW(b2) に格納 (bp*zp) ZP) Z U Z Z Y Y

64 FEM3D-Part 6 メッシュ生成コード :mgcube.f(2/5)!c cou b do ZP do j YP cou cou (-)*IBODO (j-)*p IW(coub) enddo enddo Z U Z cou b 2 do ZP Z j do P cou cou (-)*IBODO (j-)*p IW(coub) Y enddo enddo cou Y b 3 do j YP do P cou cou (-)*IBODO (j-)*p IW(coub) enddo enddo cou b ZP do j YP do P cou cou (-)*IBODO (j-)*p IW(coub) enddo enddo ZZmn の節点を IW(b3) に格納 (bp*yp) ZZmaの節点をIW(b) に格納 (bp*yp) YP)

65 FEM3D-Part 65 メッシュ生成コード :mgcube.f(3/5)!c!c !C GeoFEM data!c !C wrte (**) 'GeoFEM grdfle name?' GRIDFILE 'cube.' open (2 fle GRIDFILE status'unnown'form'formatted') wrte(2'()') IODO cou do ZP do j YP do P dfloat(-)*d YY dfloat(j-)*d ZZ dfloat(-)*d cou cou wrte (2'(3(pe6.6))') cou YY ZZ enddo enddo enddo wrte(2'()') ICELO IELMYPL 36 wrte(2'()') (IELMYPL ICELO) 節点数節点番号 座標 要素数要素タイプ ( 使わないデータ )

66 FEM3D-Part 66 cube.: 節点データ 要素数 要素タイプ (YZ) 25 5*5*5.E.E.E 2.E.E.E 3 2.E.E.E 3.E.E.E 5.E.E.E 6.E.E.E 7.E.E.E 8 2.E.E.E 9 3.E.E.E ( 途中省略 ) 2.EE.E.E 22.E.E.E 23 2.E.E.E 2 3.E.E.E 25.E.E.E 6 ** move

67 FEM3D-Part 67 メッシュ生成コード :mgcube.f(/5) cou mat do Z do j Y do cou cou n (-)*IBODO (j-)*p n2 n n3 n2 P n n3 - n5 n IBODO n6 n2 IBODO n7 n3 IBODO n8 n IBODO wrte (2'()') cou mat n n2 n3 n & & n5 n6 n7 n8 enddo enddo enddo mat: 材料番号 ()

68 FEM3D-Part 68 cube.: 要素データ ( 途中省略 )

69 FEM3D-Part 69 メッシュ生成コード :mgcube.f(5/5) IGO IB YP*ZP IB2 P*ZP IB IB3 P*YP IB2 IB P*YP IB3 wrte (2'()') IGO wrte (2'()') IB IB2 IB3 IB HHH 'mn' wrte (2'(a8)') HHH wrte (2'()') (IW() YP*ZP) HHH 'Ymn' wrte (2'(a8)') HHH wrte (2'()') ( ) (IW(2) ( P*ZP) HHH 'Zmn' wrte (2'(a8)') HHH wrte (2'()') (IW(3) P*YP) HHH 'Zma' wrte (2'(a8)') HHH wrte (2'()') (IW() P*YP) IGO IB グループ総数 (mnymnzmnzma) Zma) 累積数 ( 以下略 )!C stop end deallocate (IW) close (2)

70 FEM3D-Part 7 cube.: 節点グループデータ mn Ymn Zmn Zma ( 以下使用せず )

71 FEM3D-Part 7 メッシュ生成 実は技術的には大きな課題 複雑形状 大規模メッシュ 並列化が難しい 市販のメッシュ生成アプリケーション FEMAP CAD データとのインタフェース move

72 FEM3D-Part 72 三次元要素の定式化 三次元弾性力学方程式 ガラーキン法 要素マトリクス生成 宿題 プログラムの実行 データ構造 プログラムの構成 計算効率について Paraew による可視化

73 FEM3D-Part 73 有限要素法の処理 : プログラム 初期化 制御変数読み込み 座標読み込み 要素生成 (: 節点数 E: 要素数 ) 配列初期化 ( 全体マトリクス 要素マトリクス ) 要素 全体マトリクスマッピング (IndeItem) マトリクス生成 要素単位の処理 (do cel E) 要素マトリクス計算 全体マトリクスへの重ね合わせ 境界条件の処理 連立一次方程式 共役勾配法 (CG) 応力計算

74 FEM3D-Part 7 test メインプログラム nput_cntl 制御データ入力 nput_grd メッシュファイル入力 fnd_node 節点探索 mat_con 行列コネクティビティ生成 msor ソート mat_con 行列コネクティビティ生成 mat_ass_manass man 係数行列生成 jacob ヤコビアン計算 mat_ass_bc 境界条件処理 三次元弾性解析コードfem3d の構成 solve33 線形ソルバー制御 recover_stress 応力計算 output_ucd 可視化処理 cg_3 CG 法計算 jacob ヤコビアン計算

75 FEM3D-Part 75 全体処理 #nclude <stdo.h> #nclude <stdlb.h> FILE* fp_log; #defne GLOBAL_ALUE_DEFIE #nclude "pfem_utl.h" etern vod IPU_CL(); etern vod IPU_GRID(); etern vod MA_CO(); etern vod MA_CO(); etern vod MA_ASS_MAI(); etern vod MA_ASS_BC(); ASS etern vod SOLE33(); etern vod RECOER_SRESS(); etern vod OUPU_UCD(); nt man() { /** Logfle for debug **/ f( (fp_logfopen("log.log""w")) ULL){ fprntf(stdout"nput fle cannot be opened! n"); et(); IPU_CL(); IPU_GRID(); MA_CO(); MA_CO(); MA_ASS_MAI(); MA_ASS_BC() ; SOLE33(); RECOER_SRESS(); OUPU_UCD() ;

76 FEM3D-Part 76 Global 変数表 :pfem_utl.h/f(/3) 変数名種別サイズ I/O 内容 fname C [8] I メッシュファイル名 P I I 節点数 ICELO I I 要素数 ODGRPtot I I 節点グループ数 YZ R [][3] I 節点座標 ICELOD I [ICELO][8] I 要素コネクティビティ ODGRP_IDE I [ODGRPtot] I 各節点グループに含まれる節点数 ( 累積 ) ODGRP_IEM I [ODGRP_IDE[ODG RPO]] I 節点グループに含まれる節点 ODGRP_AME C8 [ODGRP_IDE[ODG RPO]] I 節点グループ名 L U I O 各節点非対角成分数 ( 上三角 下三角 ) PL PU I O 非対角成分総数 ( 上三角 下三角 ) D R [9*] O 全体行列 : 対角ブロック B R [3*] O 右辺ベクトル 未知数ベクトル

77 FEM3D-Part 77 Global 変数表 :pfem_utl.h/f(2/3) 変数名種別サイズ I/O 内容 ALUG R [9*] O 全体行列 : 対角ブロックの完全 LU 分解 ALAU R [9*PL][9*PU] O 全体行列 : 上 下三角成分 ndel ndeu I [] O 全体行列 : 非零非対角ブロック数 teml temu I [PL] [PU] O 全体行列 : 上 下三角ブロック ( 列番号 ) IL IU I [] O 各節点の上 下三角ブロック数 IALIAUIAU I [][L][][U][][U] O 各節点の上 下三角ブロック ( 列番号 ) IWK I [][2] O ワーク用配列 MEHOD I I 反復解法 ( に固定 ) PRECOD I I 前処理手法 (: ブロック SSOR: ブロック対角スケーリング ) IER IERactual I I 反復回数の上限 実際の反復回数 RESID R I 打ち切り誤差 (.e-8 に設定 ) SIGMA_DIAG R I LU 分解時の対角成分係数 (. に設定 ) pfemiarra I [] O 諸定数 ( 整数 ) pfemrarra R [] O 諸定数 ( 実数 )

78 FEM3D-Part 78 Global 変数表 :pfem_utl.h/f(3/3) 変数名種別サイズ I/O 内容 O8th R I.25 PQ PE P R [2][2][8] O 各ガウス積分点における η ζ ( ~ 8) POS WEI R [2] O 各ガウス積分点の座標 重み係数 COL COL2 I [] O ソート用ワーク配列 SHAPE R [2][2][2][8] O 各ガウス積分点における形状関数 (~8) P PY PZ R [2][2][2][8] O 各ガウス積分点における ( ~ 8 ) DEJ R [2][2][2] O 各ガウス積分点におけるヤコビアン行列式 ELAS POISSO R I ヤング率 ポアソン比 SIGMA_ AU_ R [][3] O 節点における垂直 せん断応力成分

79 FEM3D-Part 79 制御ファイル入力 :IPU_CL #nclude <stdo.h> #nclude <stdlb.h> #nclude "pfem_utl.h" /** **/ vod IPU_CL() { FILE *fp; f( (fpfopen("ipu.da""r")) ULL){ fprntf(stdout"nput fle cannot be opened! n"); et(); fscanf(fp"%s"fname); f(f " f fscanf(fp"%d %d"&mehod&precod); fscanf(fp"%d"&terprema); fscanf(fp"%d"&ier); fscanf(fp "%lf %lf" &ELAS &POISSO); fclose(fp); f( ( terprema < ) ){ terprema ; f( ( terprema > ) ){ terprema ; SIGMA_DIAG.; SIGMA.; RESID.e-8; SE ; pfemrarra[] RESID; pfemrarra[] SIGMA_DIAG; pfemrarra[2] SIGMA; IPU.DA cube. fname MEHOD PRECOD 2 terprema( 不使用 ) IER..3 ELAS POISSO pfemiarra[] IER; pfemiarra[] MEHOD; pfemiarra[2] PRECOD; pfemiarra[3] SE; pfemiarra[] terprema;

80 FEM3D-Part 8 メッシュ入力 :IPU_GRID(/2) #nclude <stdo.h> #nclude <stdlb.h> #nclude "pfem_utl.h" #nclude "allocate.h" vod IPU_GRID() { FILE *fp; nt jnncelse; nt YPEIMA; /** **/ f( (fpfopen(fname"r")) ULL){ fprntf(stdout"nput f( fle cannot be opened! n"); et(); ODE fscanf(fp"%d"&); P; YZ(KREAL**)allocate_matr(seof(KREAL)3); for(;<;){ for(j;j<3;j){ YZ[][j].; YZ[][3] for(;<;){ for( ;<;){ fscanf(fp"%d %lf %lf %lf"&&yz[][]&yz[][]&yz[][2]);

81 FEM3D-Part 8 allocate deallocate 関数 #nclude <stdo.h> #nclude <stdlb.h> vod* allocate_vector(nt sent m) { vod d*; *a; f ( ( a(vod * )malloc( m * se ) ) ULL ) { fprntf(stdout"error:memor does not enough! n vector n"); et(); return a; allocate を FORRA 並みに簡単にやるための関数 vod deallocate_vector(vod *a) { free( a ); vod** allocate_matr(nt sent mnt n) { vod **aa; nt ; f((aa(vod ( (vod ** )malloc( m * seof(vod*) )) ) ULL ) { fprntf(stdout"error:memor does not enough! aa n matr n"); et(); f ( ( aa[](vod * )malloc( m * n * se ) ) ULL ) { fprntf(stdout"error:memor does not enough! n matr n"); et(); for(;<m;) aa[](char*)aa[-]se*n; return aa; vod deallocate matr(vod **aa) vod deallocate_matr(vod **aa) { free( aa );

82 FEM3D-Part 82 メッシュ入力 :IPU_GRID(2/2) /** **/ ELEME fscanf(fp"%d"&icelo); ICELOD(KI**)allocate_matr(seof(KI)ICELO8); t( ICELO for(;<icelo;) fscanf(fp"%d"&ype); ICELOD[][j] の中身としては から始まる通し節点番号がそのまま読み込まれている 要素番号は から番号付け /** **/ for(cel;cel<icelo;cel){ fscanf(fp"%d %d %d %d %d %d %d %d %d %d"&&ima &ICELOD[cel][]&ICELOD[cel][]&ICELOD[cel][2]&ICELOD[cel][3] l][] l][2] l][3] &ICELOD[cel][]&ICELOD[cel][5]&ICELOD[cel][6]&ICELOD[cel][7]); ODE grp. nfo. fscanf(fp"%d"&odgrptot); ODGRP_IDE(KI* )allocate_vector(seof(ki)odgrptot); ODGRP_AME (CHAR8*)allocate_vector(seof(CHAR8)ODGRPtot); for(;<odgrptot;) ODGRP_IDE[]; for(;<odgrptot;) fscanf(fp"%d"&odgrp_ide[]); nnodgrp_ide[odgrptot]; ODGRP_IEM(KI*)allocate_vector(seof(KI)nn); for(;<odgrptot;){ S ODGRP_IDE[]; E ODGRP_IDE[]; fscanf(fp"%s"odgrp_ame[].name); nn E - S; f( nn! ){ for(s;<e;) fscanf(fp"%d"&odgrp_iem[]); fclose(fp);

83 FEM3D-Part 83 メッシュ入力 :IPU_GRID(2/2) /** **/ ELEME fscanf(fp"%d"&icelo); ICELOD(KI**)allocate_matr(seof(KI)ICELO8); t( ICELO for(;<icelo;) fscanf(fp"%d"&ype); /** **/ for(cel;cel<icelo;cel){ fscanf(fp"%d %d %d %d %d %d %d %d %d %d"&&ima &ICELOD[cel][]&ICELOD[cel][]&ICELOD[cel][2]&ICELOD[cel][3] l][] l][2] l][3] &ICELOD[cel][]&ICELOD[cel][5]&ICELOD[cel][6]&ICELOD[cel][7]); ODE grp. nfo. fscanf(fp"%d"&odgrptot); ODGRP_IDE(KI* )allocate_vector(seof(ki)odgrptot); ODGRP_AME (CHAR8*)allocate_vector(seof(CHAR8)ODGRPtot); for(;<odgrptot;) ODGRP_IDE[]; for(;<odgrptot;) fscanf(fp"%d"&odgrp_ide[]); nnodgrp_ide[odgrptot]; ODGRP_IEM(KI*)allocate_vector(seof(KI)nn); for(;<odgrptot;){ S ODGRP_IDE[]; E ODGRP_IDE[]; fscanf(fp"%s"odgrp_ame[].name); nn E - S; f( nn! ){ for(s;<e;) fscanf(fp"%d"&odgrp_iem[]); 節点グループの中身も から始まる通し節点番号がそのまま読み込まれている fclose(fp);

84 FEM3D-Part 8 fem3d: いくつかの特徴 非対角成分 上三角 下三角を分けて記憶 U ブロックとして記憶 ベクトル : 節点 3 成分 L 行列 : 各ブロック 9 成分 行列の各成分ではなく 節点上の3 変数に基づくブロックとして処理する

85 有限要素法で得られるマトリクス FEM3D-Part 85 疎行列 Φ Φ Φ F F F D D D 疎行列 が多い A( j) のように正方行列の Φ Φ Φ Φ F F F F D D D D A(j) のように正方行列の全成分を記憶することは疎行列では非効率的 Φ Φ Φ Φ Φ F F F F F D D D D D 行列では非効率的 密 行列向け Φ Φ Φ Φ Φ F F F F F D D D D D 有限要素法 : 非零非対角成分の数は高々 数百 規模 例えば未知数が 8 個あるとすると記憶容量 ( ワード数 ) は例えば未知数が 個あるとすると記憶容量 ( ワド数 ) は 正方行列 :O( 6 ) 非零非対角成分数 :O( ) 非零成分のみ記憶するのが効率的

86 FEM3D-Part 86 d.fd.c におけるマトリクス関連変数 変数名型サイズ内容 I - 未知数総数 PLU I - 連立一次方程式係数マトリクス非対角成分総数 Dag(:) R 連立一次方程式係数マトリクス対角成分 U(:) R 連立一次方程式未知数ベクトル Rhs(:) R 連立一次方程式右辺ベクトル Inde(:) I : Item(:) I PLU AMat(:) R PLU 係数マトリクス非対角成分要素番号用一次元圧縮配列 ( 非対角成分数 ) 係数マトリクス非対角成分要素番号用一次元圧縮配列 ( 非対角成分要素 ( 列 ) 番号 ) 係数マトリクス非対角成分要素番号用一次元圧縮配列 ( 非対角成分 ) 非零非対角成分のみを格納する Compressed Row Storage 法を使用している

87 行列ベクトル積への適用 ( 非零 ) 非対角成分のみを格納疎行列向け方法 FEM3D-Part 87 ( 非零 ) 非対角成分のみを格納 疎行列向け方法 Compressed Row Storage (CRS) 対角成分 ( 実数 ) Dag () 対角成分 ( 実数 ) Inde() 非対角成分数に関する一次元配列 ( 通し番号 ) ( 整数 ) () 非対角成分の要素 ( 列 ) 番号 Item() 非対角成分の要素 ( 列 ) 番号 ( 整数 nde()) AMat() 非対角成分 ( 実数 d ()) ( 実数 nde()) {Y [A]{ Φ Φ 2 2 F F D D do Y() Dag()*() d I d ( ) I d () Φ Φ Φ Φ Φ Φ F F F F F F D D D D D D do Inde(-) Inde() Y() Y() Amat()*(Item()) enddo Φ Φ Φ Φ Φ Φ F F F F F F D D D D D D enddo Φ Φ F F D D

88 FEM3D-Part 88 行列ベクトル積への適用 ( 非零 ) 非対角成分のみを格納 疎行列向け方法 Compressed Row Storage (CRS) {Q[A]{P for(;<;){ W[Q][] Dag[] * W[P][]; for(inde[];<inde[];){ W[Q][] AMat[]*W[P][Item[]]; []]

89 FEM3D-Part 89 行列ベクトル積 : 密行列 とても簡単 a a2... a a a a a a M M a a a a 2 a a2... a a {Y [A]{ do j Y(j).d do Y(j) Y(j) A(j)*() enddo enddo

90 FEM3D-Part 9 Compressed Row Storage (CRS)

91 FEM3D-Part 9 Compressed Row Storage (CRS):C 対角成分 Dag[] Dag[] Dag[2] Dag[3] Dag[].5 Dag[5] Dag[6] Dag[7]

92 FEM3D-Part 92 Compressed Row Storage (CRS):C

93 FEM3D-Part 93 Compressed Row Storage (CRS):C 非対角成分数 Inde[] Inde[] 2 Inde[2] Inde[3] Inde[] Inde[5] Inde[6] Inde[7] Inde[8] 25 PLU (Inde[]) Inde[]~Inde[]- 番目が 行目の非対角成分

94 FEM3D-Part 9 Compressed Row Storage (CRS):C 非対角成分数 Inde[] Inde[] Inde[2] Inde[3] Inde[] Inde[5] Inde[6] Inde[7] Inde[8] 25 Inde[]~Inde[]- 番目が 行目の非対角成分 PLU 25 (Inde[])

95 FEM3D-Part 95 Compressed Row Storage (CRS):C 例 : Item[ 6] AMat[ 6].5 Item[8] 2 AMat[8] 2.5

96 FEM3D-Part 96 Compressed Row Storage (CRS):C Dag [] 対角成分 ( 実数 []) Inde[] 非対角成分数に関する一次元配列 ( 通し番号 )( 整数 []) Item[] 非対角成分の要素 ( 列 ) 番号 ( 整数 [Inde[]]) AMat[] 非対角成分 ( 実数 [Inde[]]) {Y[A]{ for(;<;){ Y[] Dag[] * []; for(inde[];<inde[];){ Y[] AMat[]*[Item[]];

97 FEM3D-Part 97 ブロックとして記憶 (/3) 記憶容量が減る Inde Item に関する記憶容量を削減できる j j j j

98 FEM3D-Part 98 ブロックとして記憶 (2/3) 計算効率 間接参照 ( メモリに負担 ) と計算の比が大きくなる ベクトル スカラー共に効く :2 倍以上の性能 連続領域 キャッシュに載る ループあたりの計算量増加 do 3* do Y() D()*() (3*-2) do nde(-) nde() 2 (3*-) tem() 3 (3*) enddo enddo Y() Y() AMA()*() Y(3*-2) D(9*-8)*D(9*-7)*2D(9*-6)*3 Y(3*-) D(9*-5)*D(9*-)*2D(9*-3)*3 Y(3*I ) D(9*-2)*D(9*-)*2D(9*I )*3 do nde(-) nde() enddo enddo tem() (3*-2) 2 (3*-) 3 (3*) Y(3*-2) Y(3*-2)AMA(9*-8)*AMA(9*-7)*2 & AMA(9*-6)*3 Y(3*-) Y(3*-)AMA(9*-5)*AMA(9*-)*2 & AMA(9*-3)*3 Y(3*I ) Y(3*I )AMA(9*-2)*AMA(9*-)*2 & AMA(9* )*3

99 FEM3D-Part 99 三次元要素の定式化 三次元弾性力学方程式 ガラーキン法 要素マトリクス生成 宿題 プログラムの実行 データ構造 プログラムの構成 計算効率について Paraew による可視化

100 FEM3D-Part マイクロプロセッサの動向 CPU 性能 メモリバンド幅のギャップ

101 FEM3D-Part スカラープロセッサ CPU- キャッシュ - メモリの階層構造 FAS CPU Regster Cache 小容量 (MB): 一時置き場高価大きい ( 億以上のトランジスタ ) SLOW Man Memor 大容量 (GB) 廉価

102 FEM3D-Part 2 ベクトルプロセッサ ベクトルレジスタと高速メモリ ector Processor 単純構造のDOループの並列処理 単純 大規模な演算に適している ector Regster do A() B() C() enddo er FAS Man Memor

103 FEM3D-Part 3 典型的な挙動 : CG 法 ( 疎行列 FEM) 3..E2 2.5 % of pea GFL LOPS % of pea GFL LOPS.E.E..E.E5.E6.E7 IBM-SP3: DOF: Problem Se 問題サイズが小さい場合はキャッシュの影響のため性能が良い.E-.E.E5.E6.E7 DOF: Problem Se Earth Smulator: 大規模な問題ほどベクトル長が長くなり 性能が高い

104 FEM3D-Part プロセッサに応じたチューニング メモリ参照の最適化 に尽きる j A(j)

105 FEM3D-Part 5 プロセッサに応じたチューニング ( 続き ) ベクトルプロセッサ ループ長を大きくとる スカラープロセッサ キャッシュを有効利用 細切れにしたデータを扱う PC クラスタなどでは キャッシュサイズを大きくできない一方で メモリレイテンシ ( 立ち上がりオーバーヘッド ) の減少 メモリバンド幅の増加の傾向 しかしマルチコア化で帳消し 共通事項 メモリアクセスの連続性 局所性 計算順序の変更によって計算結果が変わる可能性について注意すること

106 FEM3D-Part 6 Scalar スカラープロセッサの代表的なチューニング ループアンローリング ループのオーバーヘッドの削減 ロード ストアの削減 ブロック化 キャッシュミスの削減

107 FEM3D-Part 7 Scalar ループアンローリング ロード ストアの削減 (/) ループ処理に対する演算の割合が増す do j do A() A() B()*C(j) enddo enddo do j - 2 do A() A() B()*C(j) A() A() B()*C(j) enddo enddo 2K 東大での計算時間.2338E E E- do j -3 do A() A() B()*C(j) A() A() B()*C(j) A() A() B()*C(j2) A() A() B()*C(j3) enddo enddo

108 FEM3D-Part 8 Scalar ループアンローリング ロード ストアの削減 (2/) ロード : メモリ キャッシュ レジスタ ストア : ロードの逆 ロード ストアが少ないほど効率良い FAS CPU Regster Cache SLOW Man Memor

109 FEM3D-Part 9 Scalar ループアンローリング ロード ストアの削減 (3/) do j do A() A() B()*C(j) ストア ロード ロードロード enddo enddo A()B()C(j) に対して 各ループでロード ストアが () ()C( j) に対して 各ルプでドストアが発生する

110 FEM3D-Part Scalar ループアンローリング ロード ストアの削減 (/) do j -3 do A() A() B()*C(j) ロードロードロード A() A() B()*C(j) A() A() B()*C(j2) A() A() B()*C(j3) ストア enddo enddo 同一ループ内で複数回表れている場合には 最初にロード 最後にストアが実施され その間レジスターに保持される 計算順序に注意

111 FEM3D-Part Scalar ループ入替によるメモリ参照最適化 (/2) YPE-A do do j A(j) A(j) B(j) enddo enddo YPE-B do j do A(j) A(j) B(j) enddo enddo j A(j) FORRA では A(j) のアドレスは A() ) A(2) A(3) A() A(2) A(22) A() A(2) A() のように並んでいる C は逆 :A[][] A[][] A[][2] A[-][] A[-][] A[-][-] この順番にアクセスしないと効率悪い ( ベクトル スカラーに共通 )

112 FEM3D-Part 2 Scalar ループ入替によるメモリ参照最適化 (2/2) YPE-A do do j A(j) A(j) B(j) enddo enddo YPE-B do j do enddo enddo A(j) A(j) B(j) YPE-A for (j; j<; j){ for (; <; ){ A[][j] A[][j] B[][j]; YPE-B for (; <; ){ for (j; j<; j){ A[][j] A[][j] B[][j]; 2K 東大での計算時間 (FORRA) ### ### 52 A 3.25E-2 B E-3 ### ### 2 A E- B 7.825E-3 ### ### 536 A E- B.5625E-2 2 ### ### 28 A E- B 3.25E-2 ### ### 256 A E- B.6875E-2 ### ### 372 A 2.523E B E-2 ### ### 358 A.92875E B E-2 ### ### 96 A E B.25E-

113 FEM3D-Part 3 Scalar ブロック化によるキャッシュミス削減 (/7) do do j A(j) A(j) B(j) enddo enddo for (j;j<;j){ for (;<;){ A[j][] A[j][] B[][j]; 計算の処理の都合でこのような順番で計算せざるを得ない場合

114 FEM3D-Part Scalar キャッシュの有効利用 キャッシュ 実際は6bte~28bteの細かいキャッシュラインに分かれている キャッシュライン単位でメモリへのリクエストが行われる FAS SLOW CPU Regster Cache Man Memor 小容量 (MB) 高価大きい ( 億以上のトランジスタ ) 大容量 (GB) 廉価 LB(ranslate Looasde Buffer) アドレス変換バッファ 仮想アドレスから実アドレスへの変換機能 LB 用のキャッシュ 通常 28 8 bte 程度 : リンク時に可変 キャッシュ が有効に使われていれば LB ミスも起こりにくい

115 FEM3D-Part 5 Scalar ブロック化によるキャッシュミス削減 (2/7) AB でメモリアクセスパターンが相反 特に B は効率が悪い j j A[j][] B[][j]

116 FEM3D-Part 6 Scalar ブロック化によるキャッシュミス削減 (3/7) 例えば キャッシュラインサイズを ワードとすると配列の値は以下のようにキャッシュに転送される j A[j][]

117 FEM3D-Part 7 Scalar ブロック化によるキャッシュミス削減 (/7) したがって A[][][ ] をアクセスしたら A[][][ ] A[][] ] A[][2] ] A[][3] が A[][9] をアクセスしたらA[][9] A[2][] A[2][] A[2][2] がそれぞれキャッシュ上にあるということになる 2 j A[j][] 9

118 FEM3D-Part 8 Scalar ブロック化によるキャッシュミス削減 (5/7) したがって 以下のようなブロック型のパターンでアクセスすれば キャッシュを有効利用できる j j A[j][] B[][j]

119 FEM3D-Part 9 Scalar ブロック化によるキャッシュミス削減 (6/7) で囲んだ部分はキャッシュに載っている j j A[j][] B[][j]

120 FEM3D-Part 2 Scalar ブロック化によるキャッシュミス削減 (7/7) 2 2 ブロック for (j;j<;j){ for (;<;){ A[j][] A[j][] B[][j]; for (j;j<-;j2){ for (;<-;2){ A[j ][ ] A[j ][ ] B[ ][j ]; A[j ][] A[j ][] B[ ][j]; A[j][ ] A[j][ ] B[][j ]; A[j][] A[j][] B[][j]; 2K 東大での実行時間 (FORRA) ### ### 2 BASIC E E-2 ### ### 536 BASIC 6.25E E-2 ### ### 372 BASIC E E- ### ### 358 BASIC E E- ### ### 96 BASIC 8.375E E-

121 FEM3D-Part 2 Scalar チューニング : まとめ スカラープロセッサ 密行列 疎行列の場合はもっと難しい ( 研究途上の課題 ) しかし 基本的な考え方は変わらない メモリアクセスの効率化

122 FEM3D-Part 22 ブロックとして記憶 (3/3) 計算の安定化 対角成分で割るのではなく 対角ブロックの完全 LU 分解を求めて解く 特に悪条件問題で有効 j j j j

123 FEM3D-Part 23 三次元要素の定式化 三次元弾性力学方程式 ガラーキン法 要素マトリクス生成 宿題 プログラムの実行 データ構造 プログラムの構成 計算効率について Paraew による可視化

124 FEM3D-Part 2 Paraew ファイルを開く 図の表示 イメージファイルの保存

125 25 UCD Format (/3) Unstructured Cell Data 要素の種類キーワード 点 線 三角形 四角形 四面体 角錐 三角柱 六面体 二次要素 線 2 三角形 2 四角形 2 四面体 2 pt lne tr quad tet pr prsm he lne2 tr2 quad2 tet2 角錐 2 pr2 三角柱 2 六面体 2 prsm2 he2

126 26 UCD Format (2/3) Orgnall for AS mcroas Etenson of the UCD fle s np here are two tpes of formats. Onl old tpe can be read b Paraew.

127 27 UCD Format (3/3): Old Format ( 全節点数 ) ( 全要素数 ) ( 各節点のデータ数 ) ( 各要素のデータ数 ) ( モデルのデータ数 ) ( 節点番号 ) ( 座標 ) (Y 座標 ) (Z 座標 ) ( 節点番号 2) ( 座標 ) (Y 座標 ) (Z 座標 ) ( 要素番号 ) ( 材料番号 ) ( 要素の種類 ) ( 要素を構成する節点のつながり ) ( 要素番号 2) ( 材料番号 ) ( 要素の種類 ) ( 要素を構成する節点のつながり ) ( 要素のデータ成分数 ) ( 成分 の構成数 ) ( 成分 2の構成数 ) ( 各成分の構成数 ) ( 要素データ成分 のラベル )( 単位 ) ( 要素データ成分 2のラベル )( 単位 ) ( 各要素データ成分のラベル )( 単位 ) ( 要素番号 ) ( 要素データ ) ( 要素データ 2) ( 要素番号 2) ( 要素データ ) ( 要素データ 2) ( 節点のデータ成分数 ) ( 成分 の構成数 ) ( 成分 2 の構成数 ) ( 各成分の構成数 ) ( 節点データ成分 のラベル )( 単位 ) ( 節点データ成分 2 のラベル )( 単位 ) ( 各節点データ成分のラベル )( 単位 ) ( 節点番号 ) ( 節点データ) ( 節点データ2) ( 節点番号 2) ( 節点データ ) ( 節点データ 2)

Microsoft PowerPoint - 3D-FEM-1.ppt [互換モード]

Microsoft PowerPoint - 3D-FEM-1.ppt [互換モード] 三次元弾性解析コード (/3) プログラムの概要 2 年夏学期 中島研吾 科学技術計算 Ⅰ(482-27) コンピュータ科学特別講義 Ⅰ(48-24) FEM3D-Part 2 対象とする問題 弾性体 Z ヤング率 E ポアソン比 ν U Z Z Z Y Y 直方体 一辺長さの立方体 ( 六面体 ) 要素 各方向に Y Z 個 境界条件 対称条件 U @ U Y @Y U Z @Z 強制変位 U Z

More information

Microsoft PowerPoint - FEM3D-C [互換モード]

Microsoft PowerPoint - FEM3D-C [互換モード] 有限要素法による 三次元定常熱伝導解析プログラム C 言語編 202 年夏季集中講義中島研吾 並列計算プログラミング (66-2057) 先端計算機演習 (66-4009) FEM3D 2 対象とする問題 : 三次元定常熱伝導 Z T T=0@Z= ma Z T 定常熱伝導 + 発熱 一様な熱伝導率 直方体 T 一辺長さの立方体 ( 六面体 ) 要素 各方向にX Y Z 個 境界条件 Q 0 X Y

More information

Microsoft PowerPoint - elast.ppt [互換モード]

Microsoft PowerPoint - elast.ppt [互換モード] 弾性力学入門 年夏学期 中島研吾 科学技術計算 Ⅰ(48-7) コンピュータ科学特別講義 Ⅰ(48-4) elast 弾性力学 弾性力学の対象 応力 弾性力学の支配方程式 elast 3 弾性力学 連続体力学 (Continuum Mechanics) 固体力学 (Solid Mechanics) の一部 弾性体 (lastic Material) を対象 弾性論 (Theor of lasticit)

More information

Microsoft PowerPoint - 03-FEM3D-C.ppt [互換モード]

Microsoft PowerPoint - 03-FEM3D-C.ppt [互換モード] 有限要素法による 三次元定常熱伝導解析プログラム C 言語編 中島研吾東京大学情報基盤センター FEM3D 2 対象とする問題 : 三次元定常熱伝導 Z T T=0@Z= ma Z T 定常熱伝導 + 発熱 一様な熱伝導率 直方体 T 一辺長さの立方体 ( 六面体 ) 要素 各方向にX Y Z 個 境界条件 Q 0 X Y X Y T=0@Z= ma 体積当たり発熱量は位置 ( メッシュの中心の座標

More information

FEM原理講座 (サンプルテキスト)

FEM原理講座 (サンプルテキスト) サンプルテキスト FEM 原理講座 サイバネットシステム株式会社 8 年 月 9 日作成 サンプルテキストについて 各講師が 講義の内容が伝わりやすいページ を選びました テキストのページは必ずしも連続していません 一部を抜粋しています 幾何光学講座については 実物のテキストではなくガイダンスを掲載いたします 対象とする構造系 物理モデル 連続体 固体 弾性体 / 弾塑性体 / 粘弾性体 / 固体

More information

Microsoft PowerPoint - FEM3D-F [互換モード]

Microsoft PowerPoint - FEM3D-F [互換モード] 有限要素法による 三次元定常熱伝導解析プログラム Fortran 編 202 年夏季集中講義中島研吾 並列計算プログラミング (66-2057) 先端計算機演習 (66-4009) FEM3D 2 対象とする問題 : 三次元定常熱伝導 Z T T=0@Z=z ma Z T z 定常熱伝導 + 発熱 一様な熱伝導率 直方体 T z 一辺長さの立方体 ( 六面体 ) 要素 各方向にX Y Z 個 境界条件

More information

(Microsoft PowerPoint - \221\34613\211\361)

(Microsoft PowerPoint - \221\34613\211\361) 計算力学 ~ 第 回弾性問題の有限要素解析 (Ⅱ)~ 修士 年後期 ( 選択科目 ) 担当 : 岩佐貴史 講義の概要 全 5 講義. 計算力学概論, ガイダンス. 自然現象の数理モデル化. 行列 場とその演算. 数値計算法 (Ⅰ) 5. 数値計算法 (Ⅱ) 6. 初期値 境界値問題 (Ⅰ) 7. 初期値 境界値問題 (Ⅱ) 8. マトリックス変位法による構造解析 9. トラス構造の有限要素解析. 重み付き残差法と古典的近似解法.

More information

Microsoft PowerPoint - 03-FEM3D-F.ppt [互換モード]

Microsoft PowerPoint - 03-FEM3D-F.ppt [互換モード] 有限要素法による 三次元定常熱伝導解析プログラム Fortran 編 中島研吾東京大学情報基盤センター FEM3D 2 対象とする問題 : 三次元定常熱伝導 Z T T=0@Z=z ma Z T z 定常熱伝導 + 発熱 一様な熱伝導率 直方体 T z 一辺長さの立方体 ( 六面体 ) 要素 各方向にX Y Z 個 境界条件 Q z 0 X Y X Y T=0@Z=z ma 体積当たり発熱量は位置

More information

パソコンシミュレータの現状

パソコンシミュレータの現状 第 2 章微分 偏微分, 写像 豊橋技術科学大学森謙一郎 2. 連続関数と微分 工学において物理現象を支配する方程式は微分方程式で表されていることが多く, 有限要素法も微分方程式を解く数値解析法であり, 定式化においては微分 積分が一般的に用いられており. 数学の基礎知識が必要になる. 図 2. に示すように, 微分は連続な関数 f() の傾きを求めることであり, 微小な に対して傾きを表し, を無限に

More information

Microsoft PowerPoint - 3D-FEM-2.ppt [互換モード]

Microsoft PowerPoint - 3D-FEM-2.ppt [互換モード] 三次元弾性解析コード (2/3) マトリクス生成 22 年夏学期 中島研吾 科学技術計算 Ⅰ(482-27) コンピュータ科学特別講義 Ⅰ(48-24) FEM3D-Prt2 2 初期化 有限要素法の処理 : プログラム 制御変数読み込み 座標読み込み 要素生成 (: 節点数 ICELTOT: 要素数 ) 配列初期化 ( 全体マトリクス 要素マトリクス ) 要素 全体マトリクスマッピング (IndeItem)

More information

<4D F736F F F696E74202D20906C8D488AC28BAB90DD8C7689F090CD8D488A D91E F1>

<4D F736F F F696E74202D20906C8D488AC28BAB90DD8C7689F090CD8D488A D91E F1> 人工環境設計解析工学構造力学と有限要素法 ( 第 回 ) 東京大学新領域創成科学研究科 鈴木克幸 固体力学の基礎方程式 変位 - ひずみの関係 適合条件式 ひずみ - 応力の関係 構成方程式 応力 - 外力の関係 平衡方程式 境界条件 変位規定境界 反力規定境界 境界条件 荷重応力ひずみ変形 場の方程式 Γ t Γ t 平衡方程式構成方程式適合条件式 構造力学の基礎式 ひずみ 一軸 荷重応力ひずみ変形

More information

Microsoft PowerPoint - シミュレーション工学-2010-第1回.ppt

Microsoft PowerPoint - シミュレーション工学-2010-第1回.ppt シミュレーション工学 ( 後半 ) 東京大学人工物工学研究センター 鈴木克幸 CA( Compter Aded geerg ) r. Jaso Lemo (SC, 98) 設計者が解析ツールを使いこなすことにより 設計の評価 設計の質の向上を図る geerg の本質の 計算機による支援 (CA CAM などより広い名前 ) 様々な汎用ソフトの登場 工業製品の設計に不可欠のツール 構造解析 流体解析

More information

Microsoft PowerPoint - 2_FrontISTRと利用可能なソフトウェア.pptx

Microsoft PowerPoint - 2_FrontISTRと利用可能なソフトウェア.pptx 東京大学本郷キャンパス 工学部8号館2階222中会議室 13:30-14:00 FrontISTRと利用可能なソフトウェア 2017年4月28日 第35回FrontISTR研究会 FrontISTRの並列計算ハンズオン 精度検証から並列性能評価まで 観測された物理現象 物理モデル ( 支配方程式 ) 連続体の運動を支配する偏微分方程式 離散化手法 ( 有限要素法, 差分法など ) 代数的な数理モデル

More information

静的弾性問題の有限要素法解析アルゴリズム

静的弾性問題の有限要素法解析アルゴリズム 概要 基礎理論. 応力とひずみおよび平衡方程式. 降伏条件式. 構成式 ( 応力 - ひずみ関係式 ) 有限要素法. 有限要素法の概要. 仮想仕事の原理式と変分原理. 平面ひずみ弾性有限要素法定式化 FEM の基礎方程式平衡方程式. G G G ひずみ - 変位関係式 w w w. kl jkl j D 構成式応力 - ひずみ関係式 ) (. 変位の境界条件力の境界条件境界条件式 t S on V

More information

Microsoft PowerPoint - 10.pptx

Microsoft PowerPoint - 10.pptx m u. 固有値とその応用 8/7/( 水 ). 固有値とその応用 固有値と固有ベクトル 行列による写像から固有ベクトルへ m m 行列 によって線形写像 f : R R が表せることを見てきた ここでは 次元平面の行列による写像を調べる とし 写像 f : を考える R R まず 単位ベクトルの像 u y y f : R R u u, u この事から 線形写像の性質を用いると 次の格子上の点全ての写像先が求まる

More information

Microsoft PowerPoint - 10.pptx

Microsoft PowerPoint - 10.pptx 0. 固有値とその応用 固有値と固有ベクトル 2 行列による写像から固有ベクトルへ m n A : m n n m 行列によって線形写像 f R R A が表せることを見てきた ここでは 2 次元平面の行列による写像を調べる 2 = 2 A 2 2 とし 写像 まず 単位ベクトルの像を求める u 2 x = v 2 y f : R A R を考える u 2 2 u, 2 2 0 = = v 2 0

More information

NS NS Scalar turbulence 5 6 FEM NS Mesh (A )

NS NS Scalar turbulence 5 6 FEM NS Mesh (A ) 22 3 2 1 2 2 2 3 3 4 NS 4 4.1 NS............ 5 5 Scalar turbulence 5 6 FEM 5 6.1 NS.................................... 6 6.2 Mes A )................................... 6 6.3.....................................

More information

情報処理概論(第二日目)

情報処理概論(第二日目) 情報処理概論 工学部物質科学工学科応用化学コース機能物質化学クラス 第 8 回 2005 年 6 月 9 日 前回の演習の解答例 多項式の計算 ( 前半 ): program poly implicit none integer, parameter :: number = 5 real(8), dimension(0:number) :: a real(8) :: x, total integer

More information

Microsoft PowerPoint - H21生物計算化学2.ppt

Microsoft PowerPoint - H21生物計算化学2.ppt 演算子の行列表現 > L いま 次元ベクトル空間の基底をケットと書くことにする この基底は完全系を成すとすると 空間内の任意のケットベクトルは > > > これより 一度基底を与えてしまえば 任意のベクトルはその基底についての成分で完全に記述することができる これらの成分を列行列の形に書くと M これをベクトル の基底 { >} による行列表現という ところで 行列 A の共役 dont 行列は A

More information

行列、ベクトル

行列、ベクトル 行列 (Mtri) と行列式 (Determinnt). 行列 (Mtri) の演算. 和 差 積.. 行列とは.. 行列の和差 ( 加減算 ).. 行列の積 ( 乗算 ). 転置行列 対称行列 正方行列. 単位行列. 行列式 (Determinnt) と逆行列. 行列式. 逆行列. 多元一次連立方程式のコンピュータによる解法. コンピュータによる逆行列の計算.. 定数項の異なる複数の方程式.. 逆行列の計算

More information

テンソル ( その ) テンソル ( その ) スカラー ( 階のテンソル ) スカラー ( 階のテンソル ) 階数 ベクトル ( 階のテンソル ) ベクトル ( 階のテンソル ) 行列表現 シンボリック表現 [ ]

テンソル ( その ) テンソル ( その ) スカラー ( 階のテンソル ) スカラー ( 階のテンソル ) 階数 ベクトル ( 階のテンソル ) ベクトル ( 階のテンソル ) 行列表現 シンボリック表現 [ ] Tsor th-ordr tsor by dcl xprsso m m Lm m k m k L mk kk quott rul by symbolc xprsso Lk X thrd-ordr tsor cotrcto j j Copyrght s rsrvd. No prt of ths documt my b rproducd for proft. テンソル ( その ) テンソル ( その

More information

Microsoft PowerPoint - H22制御工学I-10回.ppt

Microsoft PowerPoint - H22制御工学I-10回.ppt 制御工学 I 第 回 安定性 ラウス, フルビッツの安定判別 平成 年 6 月 日 /6/ 授業の予定 制御工学概論 ( 回 ) 制御技術は現在様々な工学分野において重要な基本技術となっている 工学における制御工学の位置づけと歴史について説明する さらに 制御システムの基本構成と種類を紹介する ラプラス変換 ( 回 ) 制御工学 特に古典制御ではラプラス変換が重要な役割を果たしている ラプラス変換と逆ラプラス変換の定義を紹介し

More information

Microsoft PowerPoint - FEM3D-C.ppt [互換モード]

Microsoft PowerPoint - FEM3D-C.ppt [互換モード] 3D-FEM n C: Stead State Heat Conducton Kengo aajma Informaton echnolog Center Programmng for Parallel Computng (66-2057) Semnar on Advanced Computng (66-009) FEM3D 2 3D Stead-State Heat Conducton Z =0@Z=

More information

構造力学Ⅰ第12回

構造力学Ⅰ第12回 第 回材の座屈 (0 章 ) p.5~ ( 復習 ) モールの定理 ( 手順 ) 座屈とは 荷重により梁に生じた曲げモーメントをで除して仮想荷重と考える 座屈荷重 偏心荷重 ( 曲げと軸力 ) 断面の核 この仮想荷重に対するある点でのせん断力 たわみ角に相当する曲げモーメント たわみに相当する ( 例 ) 単純梁の支点のたわみ角 : は 図 を仮想荷重と考えたときの 点の支点反力 B は 図 を仮想荷重と考えたときのB

More information

OpenFOAM(R) ソースコード入門 pt1 熱伝導方程式の解法から有限体積法の実装について考える 前編 : 有限体積法の基礎確認 2013/11/17 オープンCAE 富山富山県立大学中川慎二

OpenFOAM(R) ソースコード入門 pt1 熱伝導方程式の解法から有限体積法の実装について考える 前編 : 有限体積法の基礎確認 2013/11/17 オープンCAE 富山富山県立大学中川慎二 OpenFOAM(R) ソースコード入門 pt1 熱伝導方程式の解法から有限体積法の実装について考える 前編 : 有限体積法の基礎確認 2013/11/17 オープンCAE 勉強会 @ 富山富山県立大学中川慎二 * OpenFOAM のソースコードでは, 基礎式を偏微分方程式の形で記述する.OpenFOAM 内部では, 有限体積法を使ってこの微分方程式を解いている. どのようにして, 有限体積法に基づく離散化が実現されているのか,

More information

< B795FB8C6094C28F6F97CD97E12E786477>

< B795FB8C6094C28F6F97CD97E12E786477> 長方形板の計算システム Ver3.0 適用基準 級数解法 ( 理論解析 ) 構造力学公式集( 土木学会発行 /S61.6) 板とシェルの理論( チモシェンコ ヴォアノフスキークリ ガー共著 / 長谷川節訳 ) 有限要素法解析 参考文献 マトリックス構造解析法(J.L. ミーク著, 奥村敏恵, 西野文雄, 西岡隆訳 /S50.8) 薄板構造解析( 川井忠彦, 川島矩郎, 三本木茂夫 / 培風館 S48.6)

More information

スライド 1

スライド 1 CAE 演習 有限要素法のノウハウ ( 基礎編 ) 1. はじめに 有限要素法はポピュラーなツールである一方 解析で苦労している人が多い 高度な利用技術が必要 ( 解析の流れに沿って説明 ) 2. モデル化 要素の選択 3. メッシュ分割の工夫 4. 境界条件の設定 5. 材料物性の入力 6.7. 解析の結果の検証と分析 2. モデル化 要素の選択 モデルを単純化していかに解析を効率的 高精度に行うか?

More information

memo

memo 数理情報工学特論第一 機械学習とデータマイニング 4 章 : 教師なし学習 3 かしまひさし 鹿島久嗣 ( 数理 6 研 ) kashima@mist.i.~ DEPARTMENT OF MATHEMATICAL INFORMATICS 1 グラフィカルモデルについて学びます グラフィカルモデル グラフィカルラッソ グラフィカルラッソの推定アルゴリズム 2 グラフィカルモデル 3 教師なし学習の主要タスクは

More information

<4D F736F F F696E74202D E94D58B9393AE82F AC82B782E982BD82DF82CC8AEE E707074>

<4D F736F F F696E74202D E94D58B9393AE82F AC82B782E982BD82DF82CC8AEE E707074> 地盤数値解析学特論 防災環境地盤工学研究室村上哲 Mrakam, Satoh. 地盤挙動を把握するための基礎. 変位とひずみ. 力と応力. 地盤の変形と応力. 変位とひずみ 変形勾配テンソルひずみテンソル ひずみテンソル : 材料線素の長さの 乗の変化量の尺度 Green-Lagrange のひずみテンソルと Alman のひずみテンソル 微小変形状態でのひずみテンソル ひずみテンソルの物理的な意味

More information

破壊の予測

破壊の予測 本日の講義内容 前提 : 微分積分 線形代数が何をしているかはうろ覚え 材料力学は勉強したけど ちょっと 弾性および塑性学は勉強したことが無い ー > ですので 解らないときは質問してください モールの応力円を理解するとともに 応力を 3 次元的に考える FM( 有限要素法 の概略 内部では何を計算しているのか? 3 物が壊れる条件を考える 特に 変形 ( 塑性変形 が発生する条件としてのミーゼス応力とはどのような応力か?

More information

memo

memo 数理情報工学演習第一 C プログラミング演習 ( 第 5 回 ) 2015/05/11 DEPARTMENT OF MATHEMATICAL INFORMATICS 1 今日の内容 : プロトタイプ宣言 ヘッダーファイル, プログラムの分割 課題 : 疎行列 2 プロトタイプ宣言 3 C 言語では, 関数や変数は使用する前 ( ソースの上のほう ) に定義されている必要がある. double sub(int

More information

技術者のための構造力学 2014/06/11 1. はじめに 資料 2 節点座標系による傾斜支持節点節点の処理 三好崇夫加藤久人 従来, マトリックス変位法に基づく骨組解析を紹介する教科書においては, 全体座標系に対して傾斜 した斜面上の支持条件を考慮する処理方法として, 一旦, 傾斜支持を無視した

技術者のための構造力学 2014/06/11 1. はじめに 資料 2 節点座標系による傾斜支持節点節点の処理 三好崇夫加藤久人 従来, マトリックス変位法に基づく骨組解析を紹介する教科書においては, 全体座標系に対して傾斜 した斜面上の支持条件を考慮する処理方法として, 一旦, 傾斜支持を無視した . はじめに 資料 節点座標系による傾斜支持節点節点の処理 三好崇夫加藤久人 従来, マトリックス変位法に基づく骨組解析を紹介する教科書においては, 全体座標系に対して傾斜 した斜面上の支持条件を考慮する処理方法として, 一旦, 傾斜支持を無視した全体座標系に関する構造 全体の剛性マトリックスを組み立てた後に, 傾斜支持する節点に関して対応する剛性成分を座標変換に よって傾斜方向に回転処理し, その後は通常の全体座標系に対して傾斜していない支持点に対するのと

More information

例 e 指数関数的に減衰する信号を h( a < + a a すると, それらのラプラス変換は, H ( ) { e } e インパルス応答が h( a < ( ただし a >, U( ) { } となるシステムにステップ信号 ( y( のラプラス変換 Y () は, Y ( ) H ( ) X (

例 e 指数関数的に減衰する信号を h( a < + a a すると, それらのラプラス変換は, H ( ) { e } e インパルス応答が h( a < ( ただし a >, U( ) { } となるシステムにステップ信号 ( y( のラプラス変換 Y () は, Y ( ) H ( ) X ( 第 週ラプラス変換 教科書 p.34~ 目標ラプラス変換の定義と意味を理解する フーリエ変換や Z 変換と並ぶ 信号解析やシステム設計における重要なツール ラプラス変換は波動現象や電気回路など様々な分野で 微分方程式を解くために利用されてきた ラプラス変換を用いることで微分方程式は代数方程式に変換される また 工学上使われる主要な関数のラプラス変換は簡単な形の関数で表されるので これを ラプラス変換表

More information

以下 変数の上のドットは時間に関する微分を表わしている (ex. 2 dx d x x, x 2 dt dt ) 付録 E 非線形微分方程式の平衡点の安定性解析 E-1) 非線形方程式の線形近似特に言及してこなかったが これまでは線形微分方程式 ( x や x, x などがすべて 1 次で なおかつ

以下 変数の上のドットは時間に関する微分を表わしている (ex. 2 dx d x x, x 2 dt dt ) 付録 E 非線形微分方程式の平衡点の安定性解析 E-1) 非線形方程式の線形近似特に言及してこなかったが これまでは線形微分方程式 ( x や x, x などがすべて 1 次で なおかつ 以下 変数の上のドットは時間に関する微分を表わしている (e. d d, dt dt ) 付録 E 非線形微分方程式の平衡点の安定性解析 E-) 非線形方程式の線形近似特に言及してこなかったが これまでは線形微分方程式 ( や, などがすべて 次で なおかつそれらの係数が定数であるような微分方程式 ) に対して安定性の解析を行ってきた しかしながら 実際には非線形の微分方程式で記述される現象も多く存在する

More information

cp-7. 配列

cp-7. 配列 cp-7. 配列 (C プログラムの書き方を, パソコン演習で学ぶシリーズ ) https://www.kkaneko.jp/cc/adp/index.html 金子邦彦 1 本日の内容 例題 1. 月の日数配列とは. 配列の宣言. 配列の添え字. 例題 2. ベクトルの内積例題 3. 合計点と平均点例題 4. 棒グラフを描く配列と繰り返し計算の関係例題 5. 行列の和 2 次元配列 2 今日の到達目標

More information

Microsoft PowerPoint - H22制御工学I-2回.ppt

Microsoft PowerPoint - H22制御工学I-2回.ppt 制御工学 I 第二回ラプラス変換 平成 年 4 月 9 日 /4/9 授業の予定 制御工学概論 ( 回 ) 制御技術は現在様々な工学分野において重要な基本技術となっている 工学における制御工学の位置づけと歴史について説明する さらに 制御システムの基本構成と種類を紹介する ラプラス変換 ( 回 ) 制御工学 特に古典制御ではラプラス変換が重要な役割を果たしている ラプラス変換と逆ラプラス変換の定義を紹介し

More information

09.pptx

09.pptx 講義内容 数値解析 第 9 回 5 年 6 月 7 日 水 理学部物理学科情報理学コース. 非線形方程式の数値解法. はじめに. 分法. 補間法.4 ニュートン法.4. 多変数問題への応用.4. ニュートン法の収束性. 連立 次方程式の解法. 序論と行列計算の基礎. ガウスの消去法. 重対角行列の場合の解法項目を変更しました.4 LU 分解法.5 特異値分解法.6 共役勾配法.7 反復法.7. ヤコビ法.7.

More information

OpenCAE勉強会 公開用_pptx

OpenCAE勉強会 公開用_pptx OpenCAE 勉強会岐阜 2013/06/15 ABAQUS Student Edition を用い た XFEM き裂進展解析事例報告 OpenCAE 学会員 SH 発表内容 ABAQUS Student Edition とは? ABAQUS Student Edition 入手方法など - 入手方法 / インストール - 解析 Sample ファイルの入手方法 etc. XFEM について -XFEM

More information

FrontISTR による熱応力解析 東京大学新領域創成科学研究科人間環境学専攻橋本学 2014 年 10 月 31 日第 15 回 FrontISTR 研究会 < 機能 例題 定式化 プログラム解説編 熱応力解析 / 弾塑性解析 >

FrontISTR による熱応力解析 東京大学新領域創成科学研究科人間環境学専攻橋本学 2014 年 10 月 31 日第 15 回 FrontISTR 研究会 < 機能 例題 定式化 プログラム解説編 熱応力解析 / 弾塑性解析 > FronISR による熱応力解析 東京大学新領域創成科学研究科人間環境学専攻橋本学 214 年 1 月 31 日第 15 回 FronISR 研究会 < 機能 例題 定式化 プログラム解説編 熱応力解析 / 弾塑性解析 > FronISR に実装されている定式化を十分に理解し, 解きたい問題に対してソースコードを自由にカスタマイズ ( 要素タイプを追加, 材料の種類を追加, ユーザサブルーチンを追加

More information

슬라이드 1

슬라이드 1 SoilWorks for FLIP 主な機能特徴 1 / 13 SoilWorks for FLIP Pre-Processing 1. CADのような形状作成 修正機能 AutoCAD感覚の使いやすいモデリングや修正機能 1 CADで形状をレイヤー整理したりDXFに変換しなくても Ctrl+C でコピーしてSoilWorks上で Ctrl+V で読込む 2. AutoCAD同様のコマンドキー入力による形状作成

More information

FORTRAN( と C) によるプログラミング 5 ファイル入出力 ここではファイルからデータを読みこんだり ファイルにデータを書き出したりするプログラムを作成してみます はじめに テキスト形式で書かれたデータファイルに書かれているデータを読みこんで配列に代入し 標準出力に書き出すプログラムを作り

FORTRAN( と C) によるプログラミング 5 ファイル入出力 ここではファイルからデータを読みこんだり ファイルにデータを書き出したりするプログラムを作成してみます はじめに テキスト形式で書かれたデータファイルに書かれているデータを読みこんで配列に代入し 標準出力に書き出すプログラムを作り FORTRAN( と C) によるプログラミング 5 ファイル入出力 ここではファイルからデータを読みこんだり ファイルにデータを書き出したりするプログラムを作成してみます はじめに テキスト形式で書かれたデータファイルに書かれているデータを読みこんで配列に代入し 標準出力に書き出すプログラムを作ります FORTRAN の場合 OPEN 文でファイルを開いた後 標準入力の場合と同様に READ 文でデータを読みこみます

More information

Microsoft PowerPoint - 夏の学校(CFD).pptx

Microsoft PowerPoint - 夏の学校(CFD).pptx /9/5 FD( 計算流体力学 ) の基礎理論 性能 運動分野 夏の学校 神戸大学大学院海事科学研究科勝井辰博 流体の質量保存 流体要素内の質量の増加率 [ 単位時間当たりの増加量 ] 単位時間に流体要素に流入する質量 流体要素 Fl lm (orol olm) v ( ) ガウスの定理 v( ) /9/5 = =( ) b=b =(b b b ) b= b = b + b + b アインシュタイン表記

More information

連立方程式の解法

連立方程式の解法 連立方程式の解法連立方程式をエクセルを用いて解く方法は以下の 2 種類が考えられます 1) エクセルの行列関数を用いる 2) VBA でヤコビ法やガウスザイデル法を用いる ここでは両方について説明します 1) エクセルの行列関数を用いる方法エクセルは表計算ですから行と列に並んだ数値を扱うのは得意です 連立方程式は次のように行列を用いて表すことができます 連立方程式が行列形式で表されることを考慮して解法を考えてみます

More information

耳桁の剛性の考慮分配係数の計算条件は 主桁本数 n 格子剛度 zです 通常の並列鋼桁橋では 主桁はすべて同じ断面を使います しかし 分配の効率を上げる場合 耳桁 ( 幅員端側の桁 ) の断面を大きくすることがあります 最近の桁橋では 上下線を別橋梁とすることがあり また 防音壁などの敷設が片側に有る

耳桁の剛性の考慮分配係数の計算条件は 主桁本数 n 格子剛度 zです 通常の並列鋼桁橋では 主桁はすべて同じ断面を使います しかし 分配の効率を上げる場合 耳桁 ( 幅員端側の桁 ) の断面を大きくすることがあります 最近の桁橋では 上下線を別橋梁とすることがあり また 防音壁などの敷設が片側に有る 格子桁の分配係数の計算 ( デモ版 ) 理論と解析の背景主桁を並列した鋼単純桁の設計では 幅員方向の横桁の剛性を考えて 複数の主桁が協力して活荷重を分担する効果を計算します これを 単純な (1,0) 分配に対して格子分配と言います レオンハルト (F.Leonhardt,1909-1999) が 1950 年初頭に発表した論文が元になっていて 理論仮定 記号などの使い方は その論文を踏襲して設計に応用しています

More information

位相最適化?

位相最適化? 均質化設計法 藤井大地 ( 東京大学 ) 位相最適化? 従来の考え方 境界形状を変化させて最適な形状 位相を求める Γ t Ω b Γ D 境界形状を変化させる問題点 解析が進むにつれて, 有限要素メッシュが異形になり, 再メッシュが必要になる 位相が変化する問題への適応が難しい Γ Γ t t Ω b Ω b Γ D Γ D 領域の拡張と特性関数の導入 χ Ω ( x) = f 0 f x Ω x

More information

memo

memo 計数工学プログラミング演習 ( 第 4 回 ) 2016/05/10 DEPARTMENT OF MATHEMATICA INFORMATICS 1 内容 リスト 疎行列 2 連結リスト (inked ists) オブジェクトをある線形順序に並べて格納するデータ構造 単方向連結リスト (signly linked list) の要素 x キーフィールド key ポインタフィールド next x->next:

More information

スライド 1

スライド 1 H25 創造設計演習 ~ 振動設計演習 1~ 1 ゆれない片持ち梁の設計 振動設計演習全体 HP(2011 年度まで使用 今は閲覧のみ ): http://hockey.t.u-tokyo.ac.jp/shindousekkei/index.html M4 取付ネジ 2 Xin 加振器 50mm 幅 30mm 材料 :A2017または ABS 樹脂 計測点 :Xout 2mm? Hz CAD 所望の特性になるまで繰り返す?

More information

Microsoft PowerPoint - ロボットの運動学forUpload'C5Q [互換モード]

Microsoft PowerPoint - ロボットの運動学forUpload'C5Q [互換モード] ロボットの運動学 順運動学とは 座標系の回転と並進 同次座標変換行列 Denavit-Hartenberg の表記法 多関節ロボットの順運動学 レポート課題 & 中間試験について 逆運動学とは ヤコビアン行列 運動方程式 ( 微分方程式 ) ロボットの運動学 動力学 Equation of motion f ( ( t), ( t), ( t)) τ( t) 姿勢 ( 関節角の組合せ ) Posture

More information

PowerPoint Presentation

PowerPoint Presentation 付録 2 2 次元アフィン変換 直交変換 たたみ込み 1.2 次元のアフィン変換 座標 (x,y ) を (x,y) に移すことを 2 次元での変換. 特に, 変換が と書けるとき, アフィン変換, アフィン変換は, その 1 次の項による変換 と 0 次の項による変換 アフィン変換 0 次の項は平行移動 1 次の項は座標 (x, y ) をベクトルと考えて とすれば このようなもの 2 次元ベクトルの線形写像

More information

4 月 東京都立蔵前工業高等学校平成 30 年度教科 ( 工業 ) 科目 ( プログラミング技術 ) 年間授業計画 教科 :( 工業 ) 科目 :( プログラミング技術 ) 単位数 : 2 単位 対象学年組 :( 第 3 学年電気科 ) 教科担当者 :( 高橋寛 三枝明夫 ) 使用教科書 :( プロ

4 月 東京都立蔵前工業高等学校平成 30 年度教科 ( 工業 ) 科目 ( プログラミング技術 ) 年間授業計画 教科 :( 工業 ) 科目 :( プログラミング技術 ) 単位数 : 2 単位 対象学年組 :( 第 3 学年電気科 ) 教科担当者 :( 高橋寛 三枝明夫 ) 使用教科書 :( プロ 4 東京都立蔵前工業高等学校平成 30 年度教科 ( 工業 ) 科目 ( プログラミング技術 ) 年間授業計画 教科 :( 工業 ) 科目 :( プログラミング技術 ) 単位数 : 2 単位 対象学年組 :( 第 3 学年電気科 ) 教科担当者 :( 高橋寛 三枝明夫 ) 使用教科書 :( プログラミング技術 工業 333 実教出版 ) 共通 : 科目 プログラミング技術 のオリエンテーション プログラミング技術は

More information

有限要素法法による弾弾性変形解析 (Gmsh+Calculix)) 海洋エネルギギー研究センター今井 問題断面が1mmx1mm 長さ 20mmm の鋼の一端端を固定 他他端に点荷重重をかけた場場合の先端変変位および最大応力を求求める P Equation Chapter 1 Section 1 l

有限要素法法による弾弾性変形解析 (Gmsh+Calculix)) 海洋エネルギギー研究センター今井 問題断面が1mmx1mm 長さ 20mmm の鋼の一端端を固定 他他端に点荷重重をかけた場場合の先端変変位および最大応力を求求める P Equation Chapter 1 Section 1 l 有限要素法法による弾弾性変形解析 (Gmsh+Calculix)) 海洋エネルギギー研究センター今井 問題断面が1mmx1mm 長さ 20mmm の鋼の一端端を固定 他他端に点荷重重をかけた場場合の先端変変位および最大応力を求求める P Equation Chapter 1 Section 1 l δ 1 形状の作作成 (Gmsh) c: gmsh test1 フォルダを作る http://geuz.org/gmsh/#

More information

2018/6/12 表面の電子状態 表面に局在する電子状態 表面電子状態表面準位 1. ショックレー状態 ( 準位 ) 2. タム状態 ( 準位 ) 3. 鏡像状態 ( 準位 ) 4. 表面バンドのナローイング 5. 吸着子の状態密度 鏡像力によるポテンシャル 表面からzの位置の電子に働く力とポテン

2018/6/12 表面の電子状態 表面に局在する電子状態 表面電子状態表面準位 1. ショックレー状態 ( 準位 ) 2. タム状態 ( 準位 ) 3. 鏡像状態 ( 準位 ) 4. 表面バンドのナローイング 5. 吸着子の状態密度 鏡像力によるポテンシャル 表面からzの位置の電子に働く力とポテン 表面の電子状態 表面に局在する電子状態 表面電子状態表面準位. ショックレー状態 ( 準位. タム状態 ( 準位 3. 鏡像状態 ( 準位 4. 表面バンドのナローイング 5. 吸着子の状態密度 鏡像力によるポテンシャル 表面からzの位置の電子に働く力とポテンシャル e F z ( z z e V ( z ( Fz dz 4z e V ( z 4z ( z > ( z < のときの電子の運動を考える

More information

Microsoft PowerPoint - NA03-09black.ppt

Microsoft PowerPoint - NA03-09black.ppt きょうの講義 数値 記号処理 2003.2.6 櫻井彰人 NumSymbol@soft.ae.keo.ac.jp http://www.sakura.comp.ae.keo.ac.jp/ 数値計算手法の定石 多項式近似 ( 復習 )» 誤差と手間の解析も 漸化式» 非線型方程式の求解 数値演算上の誤差 数値計算上の誤差 打ち切り誤差 (truncaton error)» 使う公式を有限項で打ち切る

More information

<4D F736F F D2094F795AA95FB92F68EAE82CC89F082AB95FB E646F63>

<4D F736F F D2094F795AA95FB92F68EAE82CC89F082AB95FB E646F63> 力学 A 金曜 限 : 松田 微分方程式の解き方 微分方程式の解き方のところが分からなかったという声が多いので プリントにまとめます 数学的に厳密な話はしていないので 詳しくは数学の常微分方程式を扱っているテキストを参照してください また os s は既知とします. 微分方程式の分類 常微分方程式とは 独立変数 と その関数 その有限次の導関数 がみたす方程式 F,,, = のことです 次までの導関数を含む方程式を

More information

スライド 1

スライド 1 数値解析 2019 年度前期第 13 週 [7 月 11 日 ] 静岡大学創造科学技術大学院情報科学専攻工学部機械工学科計測情報講座 三浦憲二郎 講義アウトライン [7 月 11 日 ] 関数近似と補間 最小 2 乗近似による関数近似 ラグランジュ補間 T.Kanai, U.Tokyo 関数近似 p.116 複雑な関数を簡単な関数で近似する 関数近似 閉区間 [a,b] で定義された関数 f(x)

More information

JSMECM教育認定

JSMECM教育認定 一般社団法人日本機械学会 018/09/6 計算力学技術者 級問題集 ( 固体力学分野 )018 年度版 ( 第 9 版 3 刷 ) P 項目誤正 175 問 -6/ 上 8 行 1 1 sin cos sin cos rs y y xy rs y x xy i 計算力学技術者 級 ( 固体力学分野の有限要素法解析技術者 ) の認定の範囲 認定技術者の技術レベル本認定を取得した技術者は, 基本的な固体力学の問題に対して,

More information

PowerPoint Presentation

PowerPoint Presentation 応用数学 Ⅱ (7) 7 連立微分方程式の立て方と解法. 高階微分方程式による解法. ベクトル微分方程式による解法 3. 演算子による解法 連立微分方程式 未知数が複数個あり, 未知数の数だけ微分方程式が与えられている場合, これらを連立微分方程式という. d d 解法 () 高階微分方程式化による解法 つの方程式から つの未知数を消去して, 未知数が つの方程式に変換 のみの方程式にするために,

More information

モデリングとは

モデリングとは コンピュータグラフィックス基礎 第 5 回曲線 曲面の表現 ベジェ曲線 金森由博 学習の目標 滑らかな曲線を扱う方法を学習する パラメトリック曲線について理解する 広く一般的に使われているベジェ曲線を理解する 制御点を入力することで ベジェ曲線を描画するアプリケーションの開発を行えるようになる C++ 言語の便利な機能を使えるようになる 要素数が可変な配列としての std::vector の活用 計算機による曲線の表現

More information

Microsoft PowerPoint - 第3回目.ppt [互換モード]

Microsoft PowerPoint - 第3回目.ppt [互換モード] 第 3 回プログラミング応用 目的ファイル入出力 1. ファイルの概念 2. ファイルの読み込み 3. ファイルの書き込み CPU 演算 判断 ファイルの概念 内部記憶装置 OS 機械語プログラム 入力装置 キーボード 出力装置 ディスプレイ ファイル 外部記憶装置ハードディスク CD-ROM CPU が外部とデータをやり取りするための媒介 printf 関数や scanf 関数でもうすでにファイルのやり取りの基本は学んでいる

More information

<4D F736F F D208D5C91A297CD8A7793FC96E591E631308FCD2E646F63>

<4D F736F F D208D5C91A297CD8A7793FC96E591E631308FCD2E646F63> 第 1 章モールの定理による静定梁のたわみ 1-1 第 1 章モールの定理による静定梁のたわみ ポイント : モールの定理を用いて 静定梁のたわみを求める 断面力の釣合と梁の微分方程式は良く似ている 前章では 梁の微分方程式を直接積分する方法で 静定梁の断面力と変形状態を求めた 本章では 梁の微分方程式と断面力による力の釣合式が類似していることを利用して 微分方程式を直接解析的に解くのではなく 力の釣合より梁のたわみを求める方法を学ぶ

More information

Matrix and summation convention Kronecker delta δ ij 1 = 0 ( i = j) ( i j) permutation symbol e ijk = (even permutation) (odd permutation) (othe

Matrix and summation convention Kronecker delta δ ij 1 = 0 ( i = j) ( i j) permutation symbol e ijk = (even permutation) (odd permutation) (othe Matr ad summato covto Krockr dlta δ ( ) ( ) prmutato symbol k (v prmutato) (odd prmutato) (othrs) gvalu dtrmat dt 6 k rst r s kt opyrght s rsrvd. No part of ths documt may b rproducd for proft. 行列 行 正方行列

More information

Microsoft Word - 補論3.2

Microsoft Word - 補論3.2 補論 3. 多変量 GARC モデル 07//6 新谷元嗣 藪友良 対数尤度関数 3 章 7 節では 変量の対数尤度を求めた ここでは多変量の場合 とくに 変量について対数尤度を求める 誤差項 は平均 0 で 次元の正規分布に従うとする 単純化のため 分散と共分散は時間を通じて一定としよう ( この仮定は後で変更される ) したがって ij から添え字 を除くことができる このとき と の尤度関数は

More information

...Y..FEM.pm5

...Y..FEM.pm5 . 剛塑性有限要素法 名古屋大学大学院工学研究科. はじめに. 剛塑性体の構成式.. 降伏条件.. 構成方程式 ([D] マトリックス ). 節点速度 ひずみ速度関係..[B] マトリックス.. 四角形一次要素の [B] マトリックス.4 4 仮想仕事の原理 ( 剛性マトリックス ([K] マトリックス )).5 非線形方程式の解法.5. 直接代入法.5.wto-Raphso 法.6 非圧縮性の拘束と数値積分.7

More information

今後の予定 6/29 パターン形成第 11 回 7/6 データ解析第 12 回 7/13 群れ行動 ( 久保先生 ) 第 13 回 7/17 ( 金 ) 休講 7/20 まとめ第 14 回 7/27 休講?

今後の予定 6/29 パターン形成第 11 回 7/6 データ解析第 12 回 7/13 群れ行動 ( 久保先生 ) 第 13 回 7/17 ( 金 ) 休講 7/20 まとめ第 14 回 7/27 休講? 今後の予定 6/29 パターン形成第 11 回 7/6 データ解析第 12 回 7/13 群れ行動 ( 久保先生 ) 第 13 回 7/17 ( 金 ) 休講 7/20 まとめ第 14 回 7/27 休講? 数理生物学演習 第 11 回パターン形成 本日の目標 2 次元配列 分子の拡散 反応拡散モデル チューリングパタン 拡散方程式 拡散方程式 u t = D 2 u 拡散が生じる分子などの挙動を記述する.

More information

Microsoft PowerPoint - zairiki_3

Microsoft PowerPoint - zairiki_3 材料力学講義 (3) 応力と変形 Ⅲ ( 曲げモーメント, 垂直応力度, 曲率 ) 今回は, 曲げモーメントに関する, 断面力 - 応力度 - 変形 - 変位の関係について学びます 1 曲げモーメント 曲げモーメント M 静定力学で求めた曲げモーメントも, 仮想的に断面を切ることによって現れる内力です 軸方向力は断面に働く力 曲げモーメント M は断面力 曲げモーメントも, 一つのモーメントとして表しますが,

More information

数学 t t t t t 加法定理 t t t 倍角公式加法定理で α=β と置く. 三角関数

数学 t t t t t 加法定理 t t t 倍角公式加法定理で α=β と置く. 三角関数 . 三角関数 基本関係 t cot c sc c cot sc t 還元公式 t t t t t t cot t cot t 数学 数学 t t t t t 加法定理 t t t 倍角公式加法定理で α=β と置く. 三角関数 数学. 三角関数 5 積和公式 6 和積公式 数学. 三角関数 7 合成 t V v t V v t V V V V VV V V V t V v v 8 べき乗 5 6 6

More information

5-仮想仕事式と種々の応力.ppt

5-仮想仕事式と種々の応力.ppt 1 以上, 運動の変数についての話を終える. 次は再び力の変数に戻る. その前に, まず次の話が唐突と思われないように 以下は前置き. 先に, 力の変数と運動の変数には対応関係があって, 適当な内積演算によって仕事量を表す ことを述べた. 実は,Cauchy 応力と速度勾配テンソル ( あるいは変位勾配テンソル ) を用いると, それらの内積は内部仮想仕事を表していて, そして, それは外力がなす仮想仕事に等しいという

More information

Microsoft Word - 1B2011.doc

Microsoft Word - 1B2011.doc 第 14 回モールの定理 ( 単純梁の場合 ) ( モールの定理とは何か?p.11) 例題 下記に示す単純梁の C 点のたわみ角 θ C と, たわみ δ C を求めよ ただし, 部材の曲げ 剛性は材軸に沿って一様で とする C D kn B 1.5m 0.5m 1.0m 解答 1 曲げモーメント図を描く,B 点の反力を求める kn kn 4 kn 曲げモーメント図を描く knm 先に得られた曲げモーメントの値を

More information

Microsoft Word ã‡»ã…«ã‡ªã…¼ã…‹ã…žã…‹ã…³ã†¨åłºæœ›å•¤(佒芤喋çfl�)

Microsoft Word ã‡»ã…«ã‡ªã…¼ã…‹ã…žã…‹ã…³ã†¨åłºæœ›å•¤(佒芤喋çfl�) Cellulr uo nd heir eigenlues 東洋大学総合情報学部 佐藤忠一 Tdzu So Depren o Inorion Siene nd rs Toyo Uniersiy. まえがき 一次元セルオ-トマトンは数学的には記号列上の行列の固有値問題である 固有値問題の行列はふつう複素数体上の行列である 量子力学における固有値問題も無限次元ではあるが関数環上の行列でその成分は可換環である

More information

航空機の運動方程式

航空機の運動方程式 オブザーバ 状態フィードバックにはすべての状態変数の値が必要であった. しかしながら, システムの外部から観測できるのは出力だけであり, すべての状態変数が観測できるとは限らない. そこで, 制御対象システムの状態変数を, システムのモデルに基づいてその入出力信号から推定する方法を考える.. オブザーバとは 次元 m 入力 r 出力線形時不変システム x Ax Bu y Cx () の状態変数ベクトル

More information

<4D F736F F D2097CD8A7793FC96E582BD82ED82DD8A E6318FCD2E646F63>

<4D F736F F D2097CD8A7793FC96E582BD82ED82DD8A E6318FCD2E646F63> - 第 章たわみ角法の基本式 ポイント : たわみ角法の基本式を理解する たわみ角法の基本式を梁の微分方程式より求める 本章では たわみ角法の基本式を導くことにする 基本式の誘導法は各種あるが ここでは 梁の微分方程式を解いて基本式を求める方法を採用する この本で使用する座標系は 右手 右ネジの法則に従った座標を用いる また ひとつの部材では 図 - に示すように部材の左端の 点を原点とし 軸線を

More information

Microsoft Word ●IntelクアッドコアCPUでのベンチマーク_吉岡_ _更新__ doc

Microsoft Word ●IntelクアッドコアCPUでのベンチマーク_吉岡_ _更新__ doc 2.3. アプリ性能 2.3.1. Intel クアッドコア CPU でのベンチマーク 東京海洋大学吉岡諭 1. はじめにこの数年でマルチコア CPU の普及が進んできた x86 系の CPU でも Intel と AD がデュアルコア クアッドコアの CPU を次々と市場に送り出していて それらが PC クラスタの CPU として採用され HPC に活用されている ここでは Intel クアッドコア

More information

多次元レーザー分光で探る凝縮分子系の超高速動力学

多次元レーザー分光で探る凝縮分子系の超高速動力学 波動方程式と量子力学 谷村吉隆 京都大学理学研究科化学専攻 http:theochem.kuchem.kyoto-u.ac.jp TA: 岩元佑樹 iwamoto.y@kuchem.kyoto-u.ac.jp ベクトルと行列の作法 A 列ベクトル c = c c 行ベクトル A = [ c c c ] 転置ベクトル T A = [ c c c ] AA 内積 c AA = [ c c c ] c =

More information

Microsoft PowerPoint - 三次元座標測定 ppt

Microsoft PowerPoint - 三次元座標測定 ppt 冗長座標測定機 ()( 三次元座標計測 ( 第 9 回 ) 5 年度大学院講義 6 年 月 7 日 冗長性を持つ 次元座標測定機 次元 辺測量 : 冗長性を出すために つのレーザトラッカを配置し, キャッツアイまでの距離から座標を測定する つのカメラ ( 次元的なカメラ ) とレーザスキャナ : つの角度測定システムによる座標測定 つの回転関節による 次元 自由度多関節機構 高増潔東京大学工学系研究科精密機械工学専攻

More information

Microsoft PowerPoint - KN-2006NOV16.ppt

Microsoft PowerPoint - KN-2006NOV16.ppt 局所細分化メッシュに基づく並列有限 要素法における前処理付き反復法 Preconditioned Iterative Methods for Parallel Finite-Element Applications with Adaptive Mesh Refinement 中島研吾 (1) 兵藤守 (2) (1) 東京大学大学院理学系研究科地球惑星科学専攻 (2) 地球シミュレータセンター固体地球シミュレーション研究グループ

More information

<8D828D5A838A817C A77425F91E6318FCD2E6D6364>

<8D828D5A838A817C A77425F91E6318FCD2E6D6364> 4 1 平面上のベクトル 1 ベクトルとその演算 例題 1 ベクトルの相等 次の問いに答えよ. ⑴ 右の図 1 は平行四辺形 である., と等しいベクトルをいえ. ⑵ 右の図 2 の中で互いに等しいベクトルをいえ. ただし, すべてのマス目は正方形である. 解 ⑴,= より, =,= より, = ⑵ 大きさと向きの等しいものを調べる. a =d, c = f d e f 1 右の図の長方形 において,

More information

NAPRA

NAPRA 研究の動機 圧縮応力下の破壊現象 主要な亀裂の破壊に支配される引張応力下の破壊現象と異なり, 亀裂群の破壊パターンが多様 物理亀裂の進展条件 数理多様な破壊パターン 理論解析と数値解析 現状理論解析が主, 数値解析を従 将来 数値解析が主, 理論解析を従 数値解析を前提とした数理問題の設定が必要 背景 解析が困難な破壊現象 亀裂の三次元的進展 圧縮応力下での破壊 破壊問題を解くために FEM に導入される技巧

More information

PowerPoint Presentation

PowerPoint Presentation Non-linea factue mechanics き裂先端付近の塑性変形 塑性域 R 破壊進行領域応カ特異場 Ω R R Hutchinson, Rice and Rosengen 全ひずみ塑性理論に基づいた解析 現段階のひずみは 除荷がないとすると現段階の応力で一義的に決まる 単純引張り時の応カーひずみ関係 ( 構成方程式 ): ( ) ( ) n () y y y ここで α,n 定数, /

More information

<4D F736F F D E4F8E9F82C982A882AF82E98D7397F1>

<4D F736F F D E4F8E9F82C982A882AF82E98D7397F1> 3 三次における行列 要旨高校では ほとんど 2 2 の正方行列しか扱ってなく 三次の正方行列について考えてみたかったため 数 C で学んだ定理を三次の正方行列に応用して 自分たちで仮説を立てて求めていったら 空間における回転移動を表す行列 三次のケーリー ハミルトンの定理 三次における逆行列を求めたり 仮説をたてることができた. 目的 数 C で学んだ定理を三次の正方行列に応用する 2. 概要目的の到達点として

More information

Autodesk Inventor Skill Builders Autodesk Inventor 2010 構造解析の精度改良 メッシュリファインメントによる収束計算 予想作業時間:15 分 対象のバージョン:Inventor 2010 もしくはそれ以降のバージョン シミュレーションを設定する際

Autodesk Inventor Skill Builders Autodesk Inventor 2010 構造解析の精度改良 メッシュリファインメントによる収束計算 予想作業時間:15 分 対象のバージョン:Inventor 2010 もしくはそれ以降のバージョン シミュレーションを設定する際 Autodesk Inventor Skill Builders Autodesk Inventor 2010 構造解析の精度改良 メッシュリファインメントによる収束計算 予想作業時間:15 分 対象のバージョン:Inventor 2010 もしくはそれ以降のバージョン シミュレーションを設定する際に 収束判定に関するデフォルトの設定をそのまま使うか 修正をします 応力解析ソルバーでは計算の終了を判断するときにこの設定を使います

More information

数学 ⅡB < 公理 > 公理を論拠に定義を用いて定理を証明する 1 大小関係の公理 順序 (a > b, a = b, a > b 1 つ成立 a > b, b > c a > c 成立 ) 順序と演算 (a > b a + c > b + c (a > b, c > 0 ac > bc) 2 図

数学 ⅡB < 公理 > 公理を論拠に定義を用いて定理を証明する 1 大小関係の公理 順序 (a > b, a = b, a > b 1 つ成立 a > b, b > c a > c 成立 ) 順序と演算 (a > b a + c > b + c (a > b, c > 0 ac > bc) 2 図 数学 Ⅱ < 公理 > 公理を論拠に定義を用いて定理を証明する 大小関係の公理 順序 >, =, > つ成立 >, > > 成立 順序と演算 > + > + >, > > 図形の公理 平行線の性質 錯角 同位角 三角形の合同条件 三角形の合同相似 量の公理 角の大きさ 線分の長さ < 空間における座漂とベクトル > ベクトルの演算 和 差 実数倍については 文字の計算と同様 ベクトルの成分表示 平面ベクトル

More information

中学 1 年生 e ライブラリ数学教材一覧 学校図書 ( 株 ) 中学 1 年 数学 文字式式の計算 項と係数 中学 1 年 数学 次式 中学 1 年 数学 項のまとめ方 中学 1 年 数学 次式の加法 中学 1 年 数学 77

中学 1 年生 e ライブラリ数学教材一覧 学校図書 ( 株 ) 中学 1 年 数学 文字式式の計算 項と係数 中学 1 年 数学 次式 中学 1 年 数学 項のまとめ方 中学 1 年 数学 次式の加法 中学 1 年 数学 77 中学 1 年生 e ライブラリ数学教材一覧 学校図書 ( 株 ) 中学 1 年 数学 1 14-20 正の数 負の数正の数 負の数 14- ある基準から考えた量の表現 中学 1 年 数学 14- 正の数 中学 1 年 数学 14- 負の数 中学 1 年 数学 14- 量の基準を表す数 中学 1 年 数学 15- 反対の性質をもつ量の表現 中学 1 年 数学 17- 数直線 中学 1 年 数学 18-19

More information

PowerPoint Presentation

PowerPoint Presentation 工学部 6 7 8 9 10 組 ( 奇数学籍番号 ) 担当 : 長谷川英之 情報処理演習 第 7 回 2010 年 11 月 18 日 1 今回のテーマ 1: ポインタ 変数に値を代入 = 記憶プログラムの記憶領域として使用されるものがメモリ ( パソコンの仕様書における 512 MB RAM などの記述はこのメモリの量 ) RAM は多数のコンデンサの集合体 : 電荷がたまっている (1)/ いない

More information

PowerPoint Presentation

PowerPoint Presentation CAE 演習 :Eas-σ lite に よる応力解析 目標 : 機械工学実験 はりの曲げと応力集中 の有限要素法による応力解析を行う 用語 CAD: Computer Aided Design CAE: Computer Aided Engineering コンピュータシミュレーション CAM: Computer Aided Manufacturing スケジュール. 有限要素法の基礎と応用例 2.

More information

Microsoft Word - NumericalComputation.docx

Microsoft Word - NumericalComputation.docx 数値計算入門 武尾英哉. 離散数学と数値計算 数学的解法の中には理論計算では求められないものもある. 例えば, 定積分は, まずは積分 ( 被積分関数の原始関数をみつけること できなければ値を得ることはできない. また, ある関数の所定の値における微分値を得るには, まずその関数の微分ができなければならない. さらに代数方程式の解を得るためには, 解析的に代数方程式を解く必要がある. ところが, これらは必ずしも解析的に導けるとは限らない.

More information

微分方程式による現象記述と解きかた

微分方程式による現象記述と解きかた 微分方程式による現象記述と解きかた 土木工学 : 公共諸施設 構造物の有用目的にむけた合理的な実現をはかる方法 ( 技術 ) に関する学 橋梁 トンネル ダム 道路 港湾 治水利水施設 安全化 利便化 快適化 合法則的 経済的 自然および人口素材によって作られた 質量保存則 構造物の自然的な性質 作用 ( 外力による応答 ) エネルギー則 の解明 社会的諸現象のうち マスとしての移動 流通 運動量則

More information

学習指導要領

学習指導要領 (1) いろいろな式 学習指導要領紅葉川高校学力スタンダードア式と証明展開の公式を用いて 3 乗に関わる式を展開すること ( ア ) 整式の乗法 除法 分数式の計算ができるようにする 三次の乗法公式及び因数分解の公式を理解し そ 3 次の因数分解の公式を理解し それらを用いて因数れらを用いて式の展開や因数分解をすること また 分解することができるようにする 整式の除法や分数式の四則計算について理解し

More information

様々なミクロ計量モデル†

様々なミクロ計量モデル† 担当 : 長倉大輔 ( ながくらだいすけ ) この資料は私の講義において使用するために作成した資料です WEB ページ上で公開しており 自由に参照して頂いて構いません ただし 内容について 一応検証してありますが もし間違いがあった場合でもそれによって生じるいかなる損害 不利益について責任を負いかねますのでご了承ください 間違いは発見次第 継続的に直していますが まだ存在する可能性があります 1 カウントデータモデル

More information

Microsoft PowerPoint - sales2.ppt

Microsoft PowerPoint - sales2.ppt 最適化とは何? CPU アーキテクチャに沿った形で最適な性能を抽出できるようにする技法 ( 性能向上技法 ) コンパイラによるプログラム最適化 コンパイラメーカの技量 経験量に依存 最適化ツールによるプログラム最適化 KAP (Kuck & Associates, Inc. ) 人によるプログラム最適化 アーキテクチャのボトルネックを知ること 3 使用コンパイラによる性能の違い MFLOPS 90

More information

Microsoft Word - 201hyouka-tangen-1.doc

Microsoft Word - 201hyouka-tangen-1.doc 数学 Ⅰ 評価規準の作成 ( 単元ごと ) 数学 Ⅰ の目標及び図形と計量について理解させ 基礎的な知識の習得と技能の習熟を図り それらを的確に活用する機能を伸ばすとともに 数学的な見方や考え方のよさを認識できるようにする 評価の観点の趣旨 式と不等式 二次関数及び図形と計量における考え方に関 心をもつとともに 数学的な見方や考え方のよさを認識し それらを事象の考察に活用しようとする 式と不等式 二次関数及び図形と計量における数学的な見

More information

memo

memo 計数工学プログラミング演習 ( 第 3 回 ) 2016/04/26 DEPARTMENT OF MATHEMATICAL INFORMATICS 1 内容 ポインタ malloc 構造体 2 ポインタ あるメモリ領域 ( アドレス ) を代入できる変数 型は一致している必要がある 定義時には値は不定 ( 何も指していない ) 実際にはどこかのメモリを指しているので, #include

More information

演習2

演習2 神戸市立工業高等専門学校電気工学科 / 電子工学科専門科目 数値解析 2017.6.2 演習 2 山浦剛 (tyamaura@riken.jp) 講義資料ページ h t t p://clim ate.aic s. riken. jp/m embers/yamaura/num erical_analysis. html 曲線の推定 N 次多項式ラグランジュ補間 y = p N x = σ N x x

More information

変 位 変位とは 物体中のある点が変形後に 別の点に異動したときの位置の変化で あり ベクトル量である 変位には 物体の変形の他に剛体運動 剛体変位 が含まれている 剛体変位 P(x, y, z) 平行移動と回転 P! (x + u, y + v, z + w) Q(x + d x, y + dy,

変 位 変位とは 物体中のある点が変形後に 別の点に異動したときの位置の変化で あり ベクトル量である 変位には 物体の変形の他に剛体運動 剛体変位 が含まれている 剛体変位 P(x, y, z) 平行移動と回転 P! (x + u, y + v, z + w) Q(x + d x, y + dy, 変 位 変位とは 物体中のある点が変形後に 別の点に異動したときの位置の変化で あり ベクトル量である 変位には 物体の変形の他に剛体運動 剛体変位 が含まれている 剛体変位 P(x, y, z) 平行移動と回転 P! (x + u, y + v, z + w) Q(x + d x, y + dy, z + dz) Q! (x + d x + u + du, y + dy + v + dv, z +

More information

Microsoft PowerPoint - 講義PPT2019.ppt [互換モード]

Microsoft PowerPoint - 講義PPT2019.ppt [互換モード] . CA 演習 :as σ lite による応力解析 目標 : 機械工学実験 はりの曲げと応力集中 の有限要素法による応力解析を行う CAD: Computer Aided Design CA: Computer Aided ngineering コンピュータシミュレーション CAM: Computer Aided Manufacturing スケジュール. 有限要素法の基礎と応用例. as σの使い方の説明.

More information

喨微勃挹稉弑

喨微勃挹稉弑 == 全微分方程式 == 全微分とは 変数の関数 z=f(, ) について,, の増分を Δ, Δ とするとき, z の増分 Δz は Δz z Δ+ z Δ で表されます. この式において, Δ 0, Δ 0 となる極限を形式的に dz= z d+ z d (1) で表し, dz を z の全微分といいます. z は z の に関する偏導関数で, を定数と見なし て, で微分したものを表し, 方向の傾きに対応します.

More information

Microsoft Word - elastostatic_analysis_ docx

Microsoft Word - elastostatic_analysis_ docx 静弾性解析 1. 定式化と離散化の概要 1.1 線形弾性体の定式化 Fig.1 に示される線形弾性体の境界値問題を考える. ただし, 微小変形を仮定する.Fig.1 N において,N を次元数とすると, は有界領域であり, はその境界である. ここで, d は変位境界条件が与えられる境界, t は応力境界条件が与えられる境界である. d と t の間には, d および t d の関係が成り立つとする.

More information

Stage 並列プログラミングを習得するためには : 1 計算機リテラシ, プログラミング言語 2 基本的な数値解析 3 実アプリケーション ( 例えば有限要素法, 分子動力学 ) のプログラミング 4 その並列化 という 4 つの段階 (stage) が必要である 本人材育成プログラムでは1~4を

Stage 並列プログラミングを習得するためには : 1 計算機リテラシ, プログラミング言語 2 基本的な数値解析 3 実アプリケーション ( 例えば有限要素法, 分子動力学 ) のプログラミング 4 その並列化 という 4 つの段階 (stage) が必要である 本人材育成プログラムでは1~4を コンピュータ科学特別講義 科学技術計算プログラミング I ( 有限要素法 ) 中島研吾 東京大学情報基盤センター 1. はじめに本稿では,2008 年度冬学期に実施した, コンピュータ科学特別講義 I 科学技術計算プログラミング ( 有限要素法 ) について紹介する 計算科学 工学, ハードウェアの急速な進歩, 発達を背景に, 第 3 の科学 としての大規模並列シミュレーションへの期待は, 産学において一層高まっている

More information

<4D F736F F F696E74202D F A282BD94BD959C89F A4C E682528D652E707074>

<4D F736F F F696E74202D F A282BD94BD959C89F A4C E682528D652E707074> 発表の流れ SSE を用いた反復解法ライブラリ Lis 4 倍精度版の高速化 小武守恒 (JST 東京大学 ) 藤井昭宏 ( 工学院大学 ) 長谷川秀彦 ( 筑波大学 ) 西田晃 ( 中央大学 JST) はじめに 4 倍精度演算について Lisへの実装 SSEによる高速化 性能評価 スピード 収束 まとめ はじめに クリロフ部分空間法たとえば CG 法は, 理論的には高々 n 回 (n は係数行列の次元数

More information

に対して 例 2: に対して 逆行列は常に存在するとは限らない 逆行列が存在する行列を正則行列 (regular matrix) という 正則である 逆行列が存在する 一般に 正則行列 A の逆行列 A -1 も正則であり (A -1 ) -1 =A が成り立つ また 2 つの正則行列 A B の積

に対して 例 2: に対して 逆行列は常に存在するとは限らない 逆行列が存在する行列を正則行列 (regular matrix) という 正則である 逆行列が存在する 一般に 正則行列 A の逆行列 A -1 も正則であり (A -1 ) -1 =A が成り立つ また 2 つの正則行列 A B の積 2 逆行列 逆行列の計算は 連立一次方程式を数値的に解くために利用される 気象学の分野では線形系の応答問題を数値的に解くときに用いられることも多い ここでは計算機を用いて逆行列を求める方法を学ぶ 2.1 はじめにたとえば 次のような連立一次方程式を解くことを考える このような 2 元連立一次方程式は 代入法や消去法によって容易に解くことができる 解法をプログラミング言語によって記述することも困難ではない

More information