JP3924701B2 - 連続的材料分布を用いた位相最適化手法 - Google Patents

連続的材料分布を用いた位相最適化手法 Download PDF

Info

Publication number
JP3924701B2
JP3924701B2 JP2003175543A JP2003175543A JP3924701B2 JP 3924701 B2 JP3924701 B2 JP 3924701B2 JP 2003175543 A JP2003175543 A JP 2003175543A JP 2003175543 A JP2003175543 A JP 2003175543A JP 3924701 B2 JP3924701 B2 JP 3924701B2
Authority
JP
Japan
Prior art keywords
finite element
convergence
calculated
material density
stress
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Expired - Fee Related
Application number
JP2003175543A
Other languages
English (en)
Other versions
JP2004348691A (ja
Inventor
惠三 石井
和己 松井
賢二郎 寺田
眞二 西脇
Original Assignee
株式会社くいんと
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by 株式会社くいんと filed Critical 株式会社くいんと
Priority to JP2003175543A priority Critical patent/JP3924701B2/ja
Publication of JP2004348691A publication Critical patent/JP2004348691A/ja
Application granted granted Critical
Publication of JP3924701B2 publication Critical patent/JP3924701B2/ja
Anticipated expiration legal-status Critical
Expired - Fee Related legal-status Critical Current

Links

Images

Description

【0001】
【発明の属する技術分野】
本発明は、有限要素法を用いた最適構造設計に使用される手法および装置である。
【0002】
【従来の技術】
連続体の位相最適化手法は1988年デンマーク工科大学のMartin Philip Bendsoeとミシガン大学のNoboru Kikuchiによって提示された(Generating Optimal Topologies in Structural Design using a Homogenization Method.,Martin Philip Bendsoe and Noboru Kikuchi,Computer Methods in Applied Mechanics and Engineering71(1988)197−224.:文献1)。この方法は対象構造物を含む大きな設計領域を無数の小さな穴が無限に空いた多孔質体と見做し、穴の空き方に周期性を仮定する。これに均質化法という応用数学の理論を適用して、ユニットセルと呼ばれる周期の最小単位の形状に周期境界条件を課して解析することにより、各有限要素でそれぞれ穴の状態に応じて等価になる均質な応力−歪マトリックスを算出する。そして、総量を制約された材料で最大の剛性を持つ形状を求める最適設計問題を、多孔質体の穴の大きさを制御することで、領域内に材料を最適に配置し、最大の剛性を持つ形状を求めた。彼らが用いた穴の形状は2次元問題で正方形であった。この手法は1991年Katsuyuki SuzukiとNoboru Kikuchiによって更に改良され、ユニットセルの穴の形状を、正方形を含む長方形にすることと、長方形の穴の長い方の辺を力のかかる方向に回転することにより、より効率の良い異方性材料の配置が可能となり、少ない繰り返し計算でより質の高い位相形状を算出することに成功した。(A homogenization method for shape and topology optimization.,Katsuyuki Suzuki and Noboru Kikuchi,Computer Methods in Applied Mechanics and Engineering 93(1991)291−318.:文献2)。しかしながら、複雑な問題に適用した場合に、必ずしも穴の大きさは空隙、または充填状態に分かれるとは限らず、中間の状態が多く残る場合や、チェッカーボードと呼ばれる空隙、充填の状態を交互に繰り返す位相形状を得ることもしばしばであった。図1に位相最適化によく用いられる梁の問題を示す。図1上の丸のマーク(13)は対称拘束条件を表し、このモデルが上下に1/2対称であることを示す。図1左下の三角形に丸が2つついているマーク(14)は、このモデルの左下端を垂直方向にはすべりを許容し、水平方向には拘束することを表している。図1の右上の矢印(12)は右から水平方向に荷重が負荷されることを示す。この条件で、体積を領域全体の30%に制約して最適化した結果が図2である。図2の下部に示すスケールは、ユニットセル内の材料の充填状態を空隙状態が0、充填状態を1としたときに、両者間を10等分し、それを色の濃さで示したものである。薄い色の状態や、チェッカーボード状の模様が残っているのが見て取れる。この状態を改善するために、手法を踏襲しながらユニットセルの穴の形状を図4に示すような3次元の格子型にして、より鮮明な位相形状を得る手法も発明者の一人により提示されている(フレームベース・ユニットセルを用いた位相最適化の研究、石井惠三、青村茂、菊池昇、日本機械学会論文集(C編)67巻654号、499−506:文献3)が、多くの場合、算出された結果に画像処理フィルターを適用しなければ、結果をそのまま利用することは困難である。
【0003】
【発明が解決しようとする課題および目的】
前記技術で算出した位相形状を各有限要素で残った材料と、充填された状態との比で密度表示すると、図2に示すように空隙状態と充填状態との中間的な密度や、チェッカーボードと呼ばれる充填状態と空隙状態が交互に配置される分布が頻繁に現れ、設計者が算出した位相形状を基に構造形状を決定するためには、複雑な画像処理フィルターを最適化の繰り返し処理の中に組み込むなどの高度の技量を必要とした。本発明は一切のフィルター処理を用いることなく従来の位相最適化の計算結果に現れる中間密度およびチェッカーボードの分布を極力なくし、この手法を利用する技術者に形状の認識がよりしやすいシステムを提供することを目的とする。
【0004】
【課題を解決するための手段】
上記目的を達成するために、本発明による位相最適化システムは、従来の技術では有限要素分割されたモデルの各有限要素で、要素内では一定となる均質化された材料定数を定義するのに対し、本発明によるシステムは各有限要素の頂点である節点において材料密度を定義し、各有限要素内では節点の密度値と、要素で定義された補間近似関数を用いて設計変数(密度値)が要素内で連続に分布することを仮定する。これは、要素内で一定となる設計変数を定義していた従来の手法は、要素間での設計変数の連続性を阻害し、これがチェッカーボードを生む原因になると予測したことによる。
【0005】
本発明による位相最適化システムは、構造物を含む設計領域内に最適な連続的材料分布を形成して構造物の最適位相形状をコンピュータによって生成するシステムであって、対象構造の有限要素モデルデータおよび位相最適化制御データの入力手段と、各有限要素で定義された補間近似関数を用いて各有限要素の頂点である節点の材料密度から各有限要素内の任意の点の材料密度を計算し、前記計算した各有限要素内の任意の点の材料密度を用いて応力一歪マトリックスデータを作成する手段と、前記算出された応力−歪マトリックスを用いて有限要素解析を行う手段と、目的関数、制約関数およびそれぞれの感度を計算する感度解析の手段と、前記有限要素解析および感度解析により算出された目的関数と予め定められた閾値との関係を評価し収束判定する手段と、収束判定の結果に応じて設計変数である節点における材料密度を更新する手段と、前記材料密度を更新した後、前記応力一歪マトリックスデータを作成する手段と、前記有限要素解析を行う手段と、前記感度解析を行う手段と、前記収束判定する手段と、前記設計変数を更新する手段による処理を前記収束判定する手段により収束判定がなされるまで繰り返す手段を具備することを特徴とする。
【0008】
【発明の実施形態】
本発明の実施形態を図面により詳細に説明する。図5に本発明の最適位相形状を生成するプログラムのフローチャートを示す。
【0009】
図5の(1.1)において、マイクロソフト社のオペレーティングシステム、Windowsを採用したパーソナルコンピュータで稼動する、Electronic Data Systems社(アメリカ合衆国)の有限要素法プリポストプロセッサFEMAPで作成した有限要素データ、即ち節点座標情報、要素結合情報、要素と材料定数の対応情報、要素と物理定数の対応情報、材料定数情報、物理定数情報、拘束条件情報、荷重条件情報を図5に示す本発明の最適位相形状を生成するプログラムに入力する。
【0010】
図5の(1.2)において、要素毎の設計/非設計情報(設計する部分と設計しない部分を異なる番号で区別する)、制約体積(設計領域内に残す材料の体積が設計領域全体積に占める割合)、繰り返し計算を始めるステップ番号、繰り返し計算を終了させるステップ番号、収束判定に用いる閾値、設計変数の変動幅を図5に示す本発明の最適位相形状を生成するプログラムに入力する。尚、これらの値は設計の条件や設計者が解の精度を考慮して任意に決定できる。
【0011】
図5の(2)において、設計変数として定義した要素内の任意の節点iにおける材料密度ρを用いて、要素内の任意の場所における材料密度ρ(x)を以下のように計算する。
Figure 0003924701
ここでmは要素を構成する節点の数、Nは要素で定義された補間近似関数とする。
次に、密度に依存する応力−歪マトリックスは、材料が充填された状態の応力−歪マトリックスDと材料密度ρ(x)のp乗の積と仮定して算出する。
Figure 0003924701
ここで、Dは後の有限要素解析に用いる応力−歪マトリックス、xは領域内の任意の場所を示す座標、ρは設計変数として仮定した仮想的な材料密度で値は0(ゼロ)より大きく1より小さい。pは任意の正の実数で本発明では数値実験により2.5を用いた。
【0012】
図5の(3)において、前項において算出された各節点における応力−歪マトリックスを用いて対象モデルを有限要素解析し、平均コンプライアンス(領域内の総歪みエネルギーの2倍の値)を計算する。
有限要素モデルは、節点の番号と座標値、要素を構成する節点番号のデーブル、モデルの材料定数、物理定数、拘束条件、節点および要素に負荷される荷重によって構成され、未知数である節点の変位u(Ux,Uy,Uz)と要素内の歪εの関係を次のように定義する。
Figure 0003924701
ここで、uは要素を構成する節点における仮想変位ベクトル、εは要素内の歪ベクトル、Bは要素で仮定した補間近似関数の微分で構成する歪−変位関係マトリックスである。また要素内における仮想変位uによって起こる歪εと応力σの関係は次のように決定される。
Figure 0003924701
ここで、Dは応力−歪関係マトリックスと呼ばれ、線形弾性問題では等方性材料の場合はヤング率とポアソン比のみで決定される。便宜上式(4)の右辺の歪εに掛かる係数を D=ρと定義する。σは要素内応力ベクトルである。本発明では従来の手法が要素内で一定なDを仮定していたのに対し、式(2)で表される各節点で計算される応力−歪マトリックスに、要素で定義した補間近似関数を適用して、要素内で連続したDを仮定する。
仮想仕事の原理を適用して、ポテンシャルエネルギーを停留させる条件から次の関係式が導かれる。
Figure 0003924701
ここでKは各要素の剛性マトリックスを重ね合わせた全体剛性マトリックス、fは同様に要素の節点荷重成分を重ね合わせた全体荷重ベクトルとする。またΩは設計領域を示す。この平衡方程式を解いて得られた節点変位ベクトルuから前述した関係に従って歪および応力を計算し、これらの値から目的関数である次式の平均コンプライアンスを計算する。
Figure 0003924701
【0013】
図5の(4.1)において、前項で算出された目的関数の感度は、各要素で次のように計算される。
目的関数である平均コンライアンスの節点iにおける設計変数ρに関する感度は次式のようになる。
Figure 0003924701
【0014】
図5の(4.2)において、同様に制約関数g(体積)の節点iにおける設計変数ρに関する感度は次式のように計算される。
Figure 0003924701
【0015】
図5の(5)において、上記有限要素解析を行う手段により算出された目的関数とあらかじめ定められた閾値との関係を評価し収束判定する。
第n回目に計算された目的関数値(平均コンプライアンス)をΠ、第n−1回目に計算された目的関数値をΠn−1とし、最適化制御データとして入力した収束判定の閾値をrとすると、収束判定は次のように行う。
Figure 0003924701
式(9)が満足されるか、最適化制御データとして入力した繰り返し計算を終了させるステップ番号に早く到達したほうの結果を採用する。
【0016】
図5の(6)は、図5の(5)において収束しなかった場合の処理を示す。即ち、収束しなかった場合は、最適化規準法、凸線形化最適化法、逐次線形化法などの数値最適化手法により節点の密度値を更新する。
本発明では最適化規準法を使用した。最適化規準法は次のような手順で行う。1)ラグランジュ未定乗数Λを導入してラグランジュ関数を次のように定義する。
Figure 0003924701
ここで、ρは設計変数(仮想的な材料密度)、Πは目的関数(ここでは平均コンプライアンス)、V(ρ)は仮想的な材料密度による設計領域の体積、Vcは最適化制御データとして入力した制約体積を示す。
2)このラグランジュ関数が停留することが最適化の必要条件になるから、式(10)の変分を取ると、
Figure 0003924701
となる。δρ,δΛの任意性から次式が成り立つ。
Figure 0003924701
これを利用して、設計変数aおよびラグランジュ未定乗数Λに関し以下の更新ルールを適用する。
Figure 0003924701
ここでηは適当な減速パラメータで、数値実験により0.75とした。kは計算繰り返し回数を示す。
式(13)および式(14)中の右辺の括弧内の値が繰り返し計算の中で常に1になるように設計変数、ラグランジュ未定乗数を更新する。但し、算出した設計変数は上限値が1.0、下限値は0.0と決められているので、どちらかに達した場合は、上下限値を採用する。さらに、設計変数の変動幅が大きくなると収束が不安定になることから、変動幅にも最適化制御データで入力をした値(moving bound)で制約を与える。本発明では0.15を用いた。
【0017】
図5の(6)の処理、即ち設計変数を更新した後に、図5の(7)に示すように、(2)、(3)、(4)、(5)および(6)を収束判定が収束となるか、繰り返し計算を終了させるステップ番号に達するまで繰り返し、最適な位相形状を得る。
【0018】
本発明の有効性を示すために2つの具体例をあげる。
(実施例1:)
図1に示す垂直方向150.0mm、水平方向25.0mm、厚さ1.0mmの梁(11)の対称性を考慮した1/2モデルの対称境界上(13)に荷重(12)p=1.0Kgfが負荷された場合の位相最適化を行なう。拘束条件は、対称境界上に垂直方向に拘束を与え、左端下(14)に水平方向の拘束を与える。材料定数はヤング率E=21000.0Kgf/mm、ポアソン比v=0.3を使用した。モデルの有限要素分割は横×縦×奥行き:75×25×1、即ち要素辺長1mmの立方体要素を用いてモデルを構築した。
最適化の条件は、目的関数として平均コンプライアンスの最小化、制約条件は領域に残す材料を全領域の30%とする。また、設計/非設計情報は全領域を設計の対象とした。繰り返し計算を終了させるステップ番号は50とし収束判定誤差は10−6とした。
図2に従来の要素毎に設計変数を定義し、応力−歪マトリックスの作成にSuzuki、Kikuchi(文献2)の長方形の穴を3次元の直方体に拡張したユニットセルを用い、均質化法を適用して計算した結果を示す。図2の下部に示すスケールは穴を除いた材料の占める体積を、元の要素の体積で除した体積比、即ち要素の材料充填率として、空隙から充填までを0から1の間に正規化した値を10個のランクに色分けしたものである。
そして図3に本発明によるシステムで計算した結果を示す。図3の下部に示すスケールは要素を構成する各節点において計算された設計変数(密度値)の要素内平均値を10個のランクに色分けしたものである。
両者を見比べると本発明によるシステムで求めた位相形状は明らかに従来のものより簡素な形状で、しかも中間値やチェッカーボード状の位相形状も残っていない。
【0019】
(実施例2:)
図6に示すミッチェルトラスと呼ばれる、横40.0mm、縦55.0mm、厚さ1.0mmの直方体の板(61)の上端から垂直方向に15.0mm、左端から横方向に20.0mmの位置に中心を持つ半径10.0mmの穴をあけ、穴の周上(62)(図6中の太い輪で示す部分)を拘束する。そして下部左端から20.0mmの板中心に1.0Kgfの荷重(63)を左向きに(図6中矢印の方向に)載荷する。材料定数はヤング率E=21000.0Kgf/mm、ポアソン比v=0.3を使用した。
最適化の条件は、目的関数は平均コンプライアンスの最小化、制約条件は領域に残す材料を全領域の15%とする。また、設計/非設計情報は全領域を設計の対象とした。繰り返し計算を終了させるステップ番号は50とし収束判定誤差は10−6とした。
図7に従来の手法により、要素で設計変数を定義し、応力−歪マトリックスの作成にSuzuki、Kikuchi(文献2)の長方形の穴を3次元の直方体に拡張したユニットセルを用い、均質化法を適用して計算した結果を示す。図7の下部に示すスケールは穴を除いた材料の占める体積を、元の要素の体積で除した要素の材料充填率である。そして図8に本発明の手法で計算した結果を示す。図8の下部に示すスケールは要素を構成する各節点において計算された設計変数(密度値)の要素内平均値を10個のランクに色分けしたものである。
【発明の効果】
前記2種類の実施例において、図2および図7を見ても明らかなように、チェッカーボード状の位相形状が算出され、これらの結果を設計に利用するためには複雑なフィルター処理が避けられない。しかし,本発明によればフィルター処理を全く用いなくとも、従来の手法を使用した場合に比べ、位相最適化の結果に現れる中間密度やチェッカーボードを作らず、認識しやすい形状を得ることができる。
【図面の簡単な説明】
【図1】図1は梁の1/2対称モデルの有限要素分割と境界条件である。
【図2】図2は図1に示すモデルを従来の手法で計算して得られた位相形状である。
【図3】図3は図1に示すモデルを本発明による手法で計算して得られた位相形状である。
【図4】図4は文献3で使用されたユニットセル形状である。
【図5】図5は本発明によるシステムのフローチャートである。
【図6】図6はミッチェルトラスと呼ばれるモデルの有限要素分割と境界条件である。
【図7】図7は図6に示すモデルを従来の手法で計算して得られた位相形状である。
【図8】図8は図6に示すモデルを本発明による手法で計算して得られた位相形状である。
【符号の説明】
11 梁モデル
12 梁モデルに載荷される荷重
13 梁モデルの1/2対称境界の拘束条件
14 梁端部の拘束条件
61 ミッチェルトラスモデル
62 ミッチェルトラス円周部の拘束条件
63 ミッチェルトラスに載荷される荷重

Claims (1)

  1. 構造物を含む設計領域内に最適な連続的材料分布を形成して構造物の最適位相形状をコンピュータによって生成するシステムであって、対象構造の有限要素モデルデータおよび位相最適化制御データの入力手段と、各有限要素で定義された補間近似関数を用いて各有限要素の頂点である節点の材料密度から各有限要素内の任意の点の材料密度を計算し、前記計算した各有限要素内の任意の点の材料密度を用いて応力一歪マトリックスデータを作成する手段と、前記算出された応力−歪マトリックスを用いて有限要素解析を行う手段と、目的関数、制約関数およびそれぞれの感度を計算する感度解析の手段と、前記有限要素解析および感度解析により算出された目的関数と予め定められた閾値との関係を評価し収束判定する手段と、収束判定の結果に応じて設計変数である節点における材料密度を更新する手段と、前記材料密度を更新した後、前記応力一歪マトリックスデータを作成する手段と、前記有限要素解析を行う手段と、前記感度解析を行う手段と、前記収束判定する手段と、前記設計変数を更新する手段による処理を前記収束判定する手段により収束判定がなされるまで繰り返す手段を具備することを特徴とする最適位相形状を生成するシステム。
JP2003175543A 2003-05-19 2003-05-19 連続的材料分布を用いた位相最適化手法 Expired - Fee Related JP3924701B2 (ja)

Priority Applications (1)

Application Number Priority Date Filing Date Title
JP2003175543A JP3924701B2 (ja) 2003-05-19 2003-05-19 連続的材料分布を用いた位相最適化手法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
JP2003175543A JP3924701B2 (ja) 2003-05-19 2003-05-19 連続的材料分布を用いた位相最適化手法

Publications (2)

Publication Number Publication Date
JP2004348691A JP2004348691A (ja) 2004-12-09
JP3924701B2 true JP3924701B2 (ja) 2007-06-06

Family

ID=33534835

Family Applications (1)

Application Number Title Priority Date Filing Date
JP2003175543A Expired - Fee Related JP3924701B2 (ja) 2003-05-19 2003-05-19 連続的材料分布を用いた位相最適化手法

Country Status (1)

Country Link
JP (1) JP3924701B2 (ja)

Families Citing this family (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US8126684B2 (en) * 2009-04-10 2012-02-28 Livermore Software Technology Corporation Topology optimization for designing engineering product
JP6537216B2 (ja) * 2015-12-03 2019-07-03 ザ・リージェンツ・オブ・ザ・ユニバーシティ・オブ・ミシガン 異なる厚さの構造セグメントについて短縮された長さの境界を使用する形態最適化
CN109472052A (zh) * 2018-10-12 2019-03-15 中国航天空气动力技术研究院 一种高速机械手支座结构拓扑优化设计方法
EP3647973A1 (en) 2018-11-04 2020-05-06 Dassault Systèmes Designing a mechanical part with topology optimization

Also Published As

Publication number Publication date
JP2004348691A (ja) 2004-12-09

Similar Documents

Publication Publication Date Title
CN109313670A (zh) 在计算机辅助设计应用中生成晶格建议的方法和系统
CN111737835B (zh) 基于三周期极小曲面的三维多孔散热结构的设计与优化方法
Sabiston et al. 3D topology optimization for cost and time minimization in additive manufacturing
Jansen et al. A hybrid density/level set formulation for topology optimization of functionally graded lattice structures
Van Marcke et al. An improved pore network model for the computation of the saturated permeability of porous rock
US10503149B2 (en) System, method, and computer program for creating united cellular lattice structure
Fuhrimann et al. Data-driven design: Exploring new structural forms using machine learning and graphic statics
US10121279B1 (en) Systems and methods for generating a mesh
JP5854067B2 (ja) 不均質材料のシミュレーションモデルの作成方法、不均質材料のシミュレーション方法、及びプログラム
CN111523270B (zh) 一种改进的连续体结构拓扑优化后处理方法
EP1242978B1 (en) Mesh generator for and method of generating meshes in an extrusion process
JP2023525792A (ja) 製造及び構造的パフォーマンスを促進する全体的な厚さの制御によるコンピュータ支援のジェネレーティブデザイン
JP3988925B2 (ja) 混合格子型解適合格子法を用いた数値解析装置
Oishi et al. Finite elements using neural networks and a posteriori error
JP3924701B2 (ja) 連続的材料分布を用いた位相最適化手法
Abdi Evolutionary topology optimization of continuum structures using X-FEM and isovalues of structural performance
Tan et al. CFD-Micromesh: A fast geometric modeling and mesh generation tool for 3D microsystem simulations
CN111814383B (zh) 一种基于b样条密度法的自支撑结构拓扑优化设计方法
JPH1125293A (ja) メッシュ生成方法
JP3861259B2 (ja) 位相最適化システムおよび位相最適化方法
JP5774404B2 (ja) 解析装置、その方法及びそのプログラム
Ma et al. An automated approach to quadrilateral mesh generation with complex geometric feature constraints
CN110765506A (zh) 实体模型的多分辨率等几何拓扑优化方法
JP2005352818A (ja) 不均質材料のシミュレーションモデル作成方法
JP2003178100A (ja) 流体解析方法、及び、その流体解析方法を用いた流体解析装置

Legal Events

Date Code Title Description
A131 Notification of reasons for refusal

Free format text: JAPANESE INTERMEDIATE CODE: A131

Effective date: 20060523

A521 Written amendment

Free format text: JAPANESE INTERMEDIATE CODE: A523

Effective date: 20060718

A521 Written amendment

Free format text: JAPANESE INTERMEDIATE CODE: A523

Effective date: 20060921

A521 Written amendment

Free format text: JAPANESE INTERMEDIATE CODE: A523

Effective date: 20061002

TRDD Decision of grant or rejection written
A01 Written decision to grant a patent or to grant a registration (utility model)

Free format text: JAPANESE INTERMEDIATE CODE: A01

Effective date: 20070130

A61 First payment of annual fees (during grant procedure)

Free format text: JAPANESE INTERMEDIATE CODE: A61

Effective date: 20070214

R150 Certificate of patent or registration of utility model

Free format text: JAPANESE INTERMEDIATE CODE: R150

LAPS Cancellation because of no payment of annual fees