JP4599301B2 - 設計支援装置及びコンピュータプログラム - Google Patents
設計支援装置及びコンピュータプログラム Download PDFInfo
- Publication number
- JP4599301B2 JP4599301B2 JP2006008057A JP2006008057A JP4599301B2 JP 4599301 B2 JP4599301 B2 JP 4599301B2 JP 2006008057 A JP2006008057 A JP 2006008057A JP 2006008057 A JP2006008057 A JP 2006008057A JP 4599301 B2 JP4599301 B2 JP 4599301B2
- Authority
- JP
- Japan
- Prior art keywords
- equation
- error
- condition
- procedure
- solution
- 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
Links
Images
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
- G06F17/10—Complex mathematical operations
- G06F17/11—Complex mathematical operations for solving equations, e.g. nonlinear equations, general mathematical optimization problems
- G06F17/13—Differential equations
Description
本発明は、設計支援装置に関し、特に、複雑な境界条件や初期条件を伴う偏微分方程式や、積分方程式、漸化式写像差要素で表現される方程式、又は非線型方程式などの近似解を求める設計支援装置に関する。
従来、複雑な境界条件や初期条件を伴う３次元物理空間における偏微分方程式、積分方程式、変分原理を支配方程式とする問題などは、ＦＤＭ，ＦＥＭ、ＢＥＭ、Ｒａｙｌｅｉｇｈ−Ｒｉｔｚ法、重み付き残差法、グリーン関数法によってコンピュータ解析されてきた。ここに、設計支援で問題となる支配方程式は、具体的には、弾性の基礎式（大変形構造問題）、放物型偏微分方程式（熱や物質の拡散）、オイラー方程式（ポテンシャル流れ）、ナヴィエ・ストークス方程式（非線形偏微分方程式、非圧縮性及び圧縮性流れ、化学反応解析、燃焼解析）波動方程式、制御、最適化（目的関数の極大化・極小化）など、工学のあらゆる場面に及んでいる。
しかし、偏微分項の係数が、独立変数の関数となったとき（支配方程式が半線形化、準線形化したとき）、ＦＥＭ，ＢＥＭでも扱いが困難となる。
また、支配方程式が非線形であれば、解の存在、一意性が不明であったり、解の収束性が保障されないことが多い。例えば、流体方程式において、レイノズル数が大きいときや、境界形状、境界条件が単純でないときは、解の収束が得られないことが多い。
ところで、制御・最適化分野では、解の安定性解析、ファンデアポル方程式、マシュー方程式、カオス、アトラクタ、リミットサイクルなどについては、高度な解析が可能である。しかし、特性根を所望の領域に閉じ込めたり、リミットサイクルなどの軌跡を自由に設計するところまでは到達しておらず、試行錯誤で軌跡を作画するレベルにとどまっている。また、ポアンカレ切断面のコンセプトに基づき、漸化式写像作用素方程式が研究され、カオス的挙動が盛んになっているものの、漸化式Ｕｋ+１ ＝ ｆ（Ｕｋ）の究極収束点であるＵ∞なる固定点（不動点ということもある）は、カオティック挙動ゆえに求めることができないことが多い。（本発明は、漸化式写像作用素方程式の解が安定固定点或いは不安定固定点でアルトにかかわらず求解の展望を与える。）
また、例えば、Ｒｏｕｔｈ−Ｈｕｒｗｉｔｚの安定判別理論は、制御システムの特性を代数方程式で表し、代数方程式の全ての根の実数部が左半平面にある条件を表現したものがＲｏｕｔｈ−Ｈｕｒｗｉｔｚの安定規範である。この安定規範は、ｎ次の代数方程式の（ｎ−１）次以下、０次までの係数項Ａ１〜Ａｎ-１が満足する複数の拘束条件を与えることであり、ｎ次元空間において拘束条件を満足する領域を限定することと等価である。しかし、拘束条件により限定された拘束領域自体を表現することは行われず、もっぱら制御システムが安定か否かの判定に主眼が置かれていた。
しかしながら、特許文献１に開示された「多気筒内燃機関、そのクランク軸のクランク配置決定方法」では、（ｎ+１）個のクランクをクランク軸に不等角で配置し、不平衡力をゼロにする手法において、ｎ次元ベクトル空間を導入し、拘束領域（以下、特許文献１に従って、「許容領域」）という。）を可視化し、解法を容易にしている。具体的には、基準クランク角度０°とし、クランク角θｉ（ｉは０〜ｎ）を、ｚ＝ｅｘｐ（ｊθｉ）として単位円上の点ｚとして表現する。すると、クランク角全体は、（ｎ+１）成分のｚを有するベクトルＶとなる。さらに、不平衡力をゼロとする拘束条件のひとつは、単位ベクトルＩとクランク角ベクトルＶの内積がゼロであると表現されることが判明した。なお、ベクトルＶの成分が、単位円上にあることも拘束条件のひとつである。こうして、「許容領域」の幾何学的形状を描くことが可能となり、求める解は、例えば、許容領域内の谷底の点にあることが予測されるなど、特許文献１に開示された手法は非線形最適化問題の強力な一般的手法のひとつである。
特開２００５−９３２３号公報（図６、段落００３９）
しかし、特許文献１記載の発明では、内燃機関のクランク軸に関し、１次や２次の不平衡力、高次のＸ型起振力を扱えるが、取り扱い対象は限定的である。そのため、より一般的に、例えば、機械構造物のバランシング最適化、クランク軸の内部モーメントや捩じり振動の抑制などについて取り扱うことはできない。
また、特許文献１記載の発明では、円周等分方程式を拡張し、目的関数（最小化すべき関数）の根が単位円上で不等分の位置に存在する問題しか取り扱えない。
そこで、本発明の課題は、複雑な境界条件や初期条件を伴う３次元物理空間における偏微分方程式、積分方程式、漸化式方程式、非線型方程式、変分原理、そのほか、制御理論、最適化問題、などの支配方程式を代数化して高次元ベクトル空間を導入し、ベクトル成分によって表現される拘束条件が許容する部分空間（「許容領域」）を可視化し、近似誤差を最小とさせる許容領域の点を方程式の解として求める具体的手法を提供することである。
上述した課題を解決するための第１の手段は、工学分野での設計に際し、支配方程式及びその初期条件及び／又は境界条件からなる拘束条件が設定され、未知関数uに関する該支配方程式を、コンピュータを用いて解く設計支援装置であって、前記支配方程式は作用素Ｄ及び定数Ｃとして（Ｄu−Ｃ=０）なる形に表現され、
該未知関数uは、実n次元のベクトル変数の実関数とし、前記拘束条件の数をn個より少ないｍ個とする支配方程式の解を求める設計支援装置において、
前記コンピュータは、境界条件や初期条件を伴う３次元物理空間における偏微分方程式や積分方程式、さらには漸化式方程式や非線形連立方程式を扱う設計支援装置に組み込まれており、
前記コンピュータは、
前記支配方程式、及び、その初期条件及び／又は境界条件、若しくは前記方程式に伴う拘束条件を入力する方程式入力部と、
前記未知関数uを、前記実n次元に対応する有限個数ｎの規格直交基底の線形結合で近似するときの前記ｎを入力する近似次数入力部と
前記近似により生じた前記支配方程式の誤差であって、ｍａｘノルム、１−ノルム、２−ノルム、ミンコフスキーノルム、ヘルダーノルム、又は距離の公理を満足するノルムで表わされる誤差‖Ｄu−Ｃ‖を計算するＣＰＵと、
前記拘束条件によって関係付けられた前記線形結合の結合係数の数をｍ個とするとき、（ｎ−ｍ）次元実空間の領域を許容領域として表示する許容領域表示部と、
前記方程式の解を出力する解出力部とを備え、
前記ＣＰＵに前記方程式入力部より、初期条件及び／又は境界条件、若しくは前記方程式に伴う拘束条件を入力された前記支配方程式を入力して、前記近似次数入力部より入力された前記許容領域内の任意の点を候補点として、前記ＣＰＵにおいて前記候補点における前記誤差を計算して、
前記各候補点における前記誤差‖Ｄu−Ｃ‖が最小ならしめる前記許容領域の候補点を前記方程式の解として求め、
該方程式の解を前記解出力部を介して記拘束条件によって関係付けられた前記線形結合の結合係数の数ｍを、（ｎ−ｍ）次元実空間の領域が、３次元以下の許容領域となるように設定し、
前記許容領域表示部は、前記許容領域が２次元以下のときは２次元以下のグラフィック表示を行い、許容領域が３次元以上のときは３次元のグラフィック表示を行い、又は、２次元領域に階層化された複数のグラフとして表示し、許容領域の曲面が３次元以上のときは、３次元のグラフィック表示を行い、又は、２次元領域に階層化された立体グラフ又は等高線図として表示することを特徴とする設計支援装置にある。
該未知関数uは、実n次元のベクトル変数の実関数とし、前記拘束条件の数をn個より少ないｍ個とする支配方程式の解を求める設計支援装置において、
前記コンピュータは、境界条件や初期条件を伴う３次元物理空間における偏微分方程式や積分方程式、さらには漸化式方程式や非線形連立方程式を扱う設計支援装置に組み込まれており、
前記コンピュータは、
前記支配方程式、及び、その初期条件及び／又は境界条件、若しくは前記方程式に伴う拘束条件を入力する方程式入力部と、
前記未知関数uを、前記実n次元に対応する有限個数ｎの規格直交基底の線形結合で近似するときの前記ｎを入力する近似次数入力部と
前記近似により生じた前記支配方程式の誤差であって、ｍａｘノルム、１−ノルム、２−ノルム、ミンコフスキーノルム、ヘルダーノルム、又は距離の公理を満足するノルムで表わされる誤差‖Ｄu−Ｃ‖を計算するＣＰＵと、
前記拘束条件によって関係付けられた前記線形結合の結合係数の数をｍ個とするとき、（ｎ−ｍ）次元実空間の領域を許容領域として表示する許容領域表示部と、
前記方程式の解を出力する解出力部とを備え、
前記ＣＰＵに前記方程式入力部より、初期条件及び／又は境界条件、若しくは前記方程式に伴う拘束条件を入力された前記支配方程式を入力して、前記近似次数入力部より入力された前記許容領域内の任意の点を候補点として、前記ＣＰＵにおいて前記候補点における前記誤差を計算して、
前記各候補点における前記誤差‖Ｄu−Ｃ‖が最小ならしめる前記許容領域の候補点を前記方程式の解として求め、
該方程式の解を前記解出力部を介して記拘束条件によって関係付けられた前記線形結合の結合係数の数ｍを、（ｎ−ｍ）次元実空間の領域が、３次元以下の許容領域となるように設定し、
前記許容領域表示部は、前記許容領域が２次元以下のときは２次元以下のグラフィック表示を行い、許容領域が３次元以上のときは３次元のグラフィック表示を行い、又は、２次元領域に階層化された複数のグラフとして表示し、許容領域の曲面が３次元以上のときは、３次元のグラフィック表示を行い、又は、２次元領域に階層化された立体グラフ又は等高線図として表示することを特徴とする設計支援装置にある。
第２手段は、第１手段において、前記支配方程式は、微分作用素、積分作用要素、漸化式の写像作用素、又は非線形方程式における非線形作用素をＤ、定数をＣ、未知関数をｕとするとき、Ｄｕ＝Ｃであり、ｕを線形結合で近似し且つ方程式全体を代数化した後のＤｕ＝Ｃを使って、誤差であるノルム||Ｄｕ−Ｃ||を計算する。
第３手段は、第１手段において、ＣＰＵは、線形結合の係数ｂｉ（ｉは１からｎまで）を根とする代数方程式において、根のすべてを実数に限定する。
第４手段は、工学分野での設計に用い、支配方程式及びその初期条件及び／又は境界条件からなる拘束条件が設定され、未知関数uに関する該支配方程式を解く為に、前記支配方程式は作用素Ｄ及び定数Ｃとして（Ｄu−Ｃ＝０）なる形に表現され、
該未知関数uは、実n次元のベクトル変数の実関数とし、前記拘束条件の数をn個より少ないｍ個とする支配方程式の解を求めるＣＰＵを具えたコンピュータにおいて、
前記コンピュータは、前記支配方程式、及び、その初期条件及び／又は境界条件を入力する方程式入力部と、
前記未知関数を、前記実n次元に対応する有限個数ｎの規格直交基底の線形結合で近似するときの前記ｎを入力する近似次数入力部と、
前記近似により生じた前記支配方程式の誤差であって、ｍａｘノルム、１−ノルム、２−ノルム、ミンコフスキーノルム、ヘルダーノルム、又は距離の公理を満足するノルムで表わされる誤差‖Ｄu−Ｃ‖を計算するＣＰＵと、
前記初期条件及び／又は境界条件によって関係付けられた前記線形結合の結合係数の数をｍ個とするとき、（ｎ−ｍ）次元実空間の領域を許容領域として表示する許容領域表示部とを具えるとともに、
前記ＣＰＵに組み込まれた前記支配方程式の解を出力する解出力部として機能させるコンピュータプログラムに基づいて
前記ＣＰＵで前記方程式入力部より、初期条件及び／又は境界条件、若しくは前記方程式に伴う拘束条件を入力された前記支配方程式を入力して、前記近似次数入力部より入力された前記許容領域内の任意の点を候補点として、前記候補点における前記誤差を計算して、
前記各候補点における前記誤差‖Ｄu−Ｃ‖が最小ならしめる前記許容領域の候補点を前記方程式の解として求め
該方程式の解を前記解出力部を介して記拘束条件によって関係付けられた前記線形結合の結合係数の数ｍを、（ｎ−ｍ）次元実空間の領域が、３次元以下の許容領域となるように設定し、
前記許容領域表示部は、前記許容領域が２次元以下のときは２次元以下のグラフィック表示を行い、許容領域が３次元以上のときは３次元のグラフィック表示を行い、又は、２次元領域に階層化された複数のグラフとして表示し、許容領域の曲面が３次元以上のときは、３次元のグラフィック表示を行い、又は、２次元領域に階層化された立体グラフ又は等高線図として表示することを特徴とするコンピュータにある。
該未知関数uは、実n次元のベクトル変数の実関数とし、前記拘束条件の数をn個より少ないｍ個とする支配方程式の解を求めるＣＰＵを具えたコンピュータにおいて、
前記コンピュータは、前記支配方程式、及び、その初期条件及び／又は境界条件を入力する方程式入力部と、
前記未知関数を、前記実n次元に対応する有限個数ｎの規格直交基底の線形結合で近似するときの前記ｎを入力する近似次数入力部と、
前記近似により生じた前記支配方程式の誤差であって、ｍａｘノルム、１−ノルム、２−ノルム、ミンコフスキーノルム、ヘルダーノルム、又は距離の公理を満足するノルムで表わされる誤差‖Ｄu−Ｃ‖を計算するＣＰＵと、
前記初期条件及び／又は境界条件によって関係付けられた前記線形結合の結合係数の数をｍ個とするとき、（ｎ−ｍ）次元実空間の領域を許容領域として表示する許容領域表示部とを具えるとともに、
前記ＣＰＵに組み込まれた前記支配方程式の解を出力する解出力部として機能させるコンピュータプログラムに基づいて
前記ＣＰＵで前記方程式入力部より、初期条件及び／又は境界条件、若しくは前記方程式に伴う拘束条件を入力された前記支配方程式を入力して、前記近似次数入力部より入力された前記許容領域内の任意の点を候補点として、前記候補点における前記誤差を計算して、
前記各候補点における前記誤差‖Ｄu−Ｃ‖が最小ならしめる前記許容領域の候補点を前記方程式の解として求め
該方程式の解を前記解出力部を介して記拘束条件によって関係付けられた前記線形結合の結合係数の数ｍを、（ｎ−ｍ）次元実空間の領域が、３次元以下の許容領域となるように設定し、
前記許容領域表示部は、前記許容領域が２次元以下のときは２次元以下のグラフィック表示を行い、許容領域が３次元以上のときは３次元のグラフィック表示を行い、又は、２次元領域に階層化された複数のグラフとして表示し、許容領域の曲面が３次元以上のときは、３次元のグラフィック表示を行い、又は、２次元領域に階層化された立体グラフ又は等高線図として表示することを特徴とするコンピュータにある。
第５手段は、第４手段において、前記支配方程式は、微分作用素をＤ、定数をＣ、未知関数をｕとするとき、Ｄｕ＝Ｃであり、ｕを線形結合で近似し且つ方程式全体を代数化した後のＤｕ＝Ｃを使って、方程式の誤差であるノルム||Ｄｕ−Ｃ||を計算する。
第６手段は、第４手段において、ＣＰＵは、線形結合の係数ｂｉ（ｉは１からｎまで）を根とする代数方程式において、根のすべてを実数に制限するか、又は、根が１対をなす共役根であるときは複素平面内の一部の領域に存在するものに限定する。
第７手段は、第４手段において、
前記ＣＰＵに組み込まれたコンピュータプログラムは
未知関数ｕ、初期条件、境界条件を入力する第１手順と、
前記ｕを変数分離し、各変数の関数を、各々異なる正規直交基底の線形結合で近似する第２手順と、
前記方程式を前記正規直交基底で表わし、その表現に含まれる線形結合定数を根とする代数方程式を生成する第３手順と、
前記初期条件、前記境界条件を前記正規直交基底で表わし、１以上の拘束条件を得る第４手順と、
前記代数方程式の根を前記拘束条件の下で求める第５手順と、
前記根が前記代数方程式の前記誤差を最小ならしめるか否かを判定する第６手順とを含み、
前記誤差が最小であると判定されれば、前記根を前記方程式の解とし、
前記誤差が最小でないと判定されれば、第５手順に戻る。
前記ＣＰＵに組み込まれたコンピュータプログラムは
未知関数ｕ、初期条件、境界条件を入力する第１手順と、
前記ｕを変数分離し、各変数の関数を、各々異なる正規直交基底の線形結合で近似する第２手順と、
前記方程式を前記正規直交基底で表わし、その表現に含まれる線形結合定数を根とする代数方程式を生成する第３手順と、
前記初期条件、前記境界条件を前記正規直交基底で表わし、１以上の拘束条件を得る第４手順と、
前記代数方程式の根を前記拘束条件の下で求める第５手順と、
前記根が前記代数方程式の前記誤差を最小ならしめるか否かを判定する第６手順とを含み、
前記誤差が最小であると判定されれば、前記根を前記方程式の解とし、
前記誤差が最小でないと判定されれば、第５手順に戻る。
第８手段は、第４手段において、
前記ＣＰＵに組み込まれたコンピュータプログラムは
未知関数ｕ、初期条件、境界条件を入力する第１手順と、
前記ｕを正規直交基底の線形結合で近似する第２手順と、
前記方程式を定積分し、前記定積分に含まれる線形結合定数を根とする代数方程式を生成する第３手順と、
前記初期条件、前記境界条件を前記正規直交基底で表わし、１以上の拘束条件を得る第４手順と、
前記代数方程式の根を前記拘束条件の下で求める第５手順と、
前記根が前記代数方程式の前記誤差を最小ならしめるか否かを判定する第６手順とを含み、
前記誤差が最小であると判定されれば、前記根を前記方程式の解とし、
前記誤差が最小でないと判定されれば、第５手順に戻る。
前記ＣＰＵに組み込まれたコンピュータプログラムは
未知関数ｕ、初期条件、境界条件を入力する第１手順と、
前記ｕを正規直交基底の線形結合で近似する第２手順と、
前記方程式を定積分し、前記定積分に含まれる線形結合定数を根とする代数方程式を生成する第３手順と、
前記初期条件、前記境界条件を前記正規直交基底で表わし、１以上の拘束条件を得る第４手順と、
前記代数方程式の根を前記拘束条件の下で求める第５手順と、
前記根が前記代数方程式の前記誤差を最小ならしめるか否かを判定する第６手順とを含み、
前記誤差が最小であると判定されれば、前記根を前記方程式の解とし、
前記誤差が最小でないと判定されれば、第５手順に戻る。
第９手段は、第４手段において、
前記ＣＰＵに組み込まれたコンピュータプログラムは
未知関数ｕ、初期条件、境界条件を入力する第１手順と、
前記ｕを正規直交基底の線形結合で近似する第２手順と、
前記線形結合に含まれる線形結合定数を根とする代数方程式を生成する第３手順と、
前記初期条件、前記境界条件を前記正規直交基底で表わし、１以上の拘束条件を得る第４手順と、
前記代数方程式の根を前記拘束条件の下で求める第５手順と、
前記根が前記代数方程式の前記誤差を最小ならしめるか否かを判定する第６手順とを含み、
前記誤差が最小であると判定されれば、前記根を前記方程式の解とし、
前記誤差が最小でないと判定されれば、第５手順に戻る。
前記ＣＰＵに組み込まれたコンピュータプログラムは
未知関数ｕ、初期条件、境界条件を入力する第１手順と、
前記ｕを正規直交基底の線形結合で近似する第２手順と、
前記線形結合に含まれる線形結合定数を根とする代数方程式を生成する第３手順と、
前記初期条件、前記境界条件を前記正規直交基底で表わし、１以上の拘束条件を得る第４手順と、
前記代数方程式の根を前記拘束条件の下で求める第５手順と、
前記根が前記代数方程式の前記誤差を最小ならしめるか否かを判定する第６手順とを含み、
前記誤差が最小であると判定されれば、前記根を前記方程式の解とし、
前記誤差が最小でないと判定されれば、第５手順に戻る。
第１及び第４手段によれば、未知関数を正規直交基底で展開したために方程式に生じた誤差を、最小にすることが出来る。
第２及び第５手段によれば、方程式を代数化した後の近似方程式の誤差を使用して解の精度を検証することが出来る。
第１及び第４手段によれば、許容領域の谷底や山頂をリアルに注目しながら、求解を進めることが出来る。
第３及び第６手段によれば、線形結合係数の値を絞り込むことが出来る。
第７手段によれば、変数分離を必要とする方程式を簡単に解くことが出来る。
第８手段によれば、定積分を必要とする方程式を簡単に解くことが出来る。
第９手段によれば、変数分離も定積分も必要としない方程式を簡単に解くことが出来る。
本発明によれば、複雑な境界条件や初期条件を伴う３次元物理空間における偏微分方程式、積分方程式、漸化式方程式、非線型方程式、変分原理、その他、制御理論、最適化問題、などの支配方程式を代数化して高次元ベクトル空間を導入し、ベクトル成分によって表現される拘束条件が許容する部分空間（「許容領域」）を可視化し、近似により生じた方程式の誤差を最小とさせる許容領域の点を方程式の解として求める具体的手法を提供することができる。
特に、本発明によれば、従来のように境界条件が方程式に付随して取り扱われるのではなく、方程式も境界条件も代数化されるため両者が対等に扱われる。代数化された方程式の解の許容領域及び解の曲面は、３次元表示などで容易に可視化できるので、解法の見通しが極めてよい。
本発明の実施形態について説明する前に、まず、本発明の数学的原理について説明する。
［本発明の数学的原理］
独立変数（ｘ、ｙ、ｚ、ｔ）及びこれらの微分演算子を含み、未知関数ｆに関する方程式は、例えば、ナビエ・ストークス方程式の場合、密度ρ、圧力ｐはスカラー量であるが、速度（ｕ、ｖ、ｗ）はベクトル量である。このように、未知関数も独立変数もベクトルが使用され且つ微積分作用素を含む関数方程式では、独立変数空間の領域Ωでその方程式が成立し、領域Ωの境界Γで拘束条件（初期条件、境界条件）が与えられる。しかし、境界条件が複雑になれば、関数方程式を解析的に解くことは困難となってくる。そこで、関数方程式を近似的に解くことが考えられるが、本発明は、その手法のひとつ実現する装置及びコンピュータプログラムに関する。
本発明では、未知関数ｆ（ｘ）を関数基底で展開することを考える。具体的には、計算の簡易性のため、正規直交基底ＯＮＢ（ｘ）を有限個使用する線形結合により未知関数ｆ（ｘ）を近似する。ＯＮＢ（ｘ）には、方程式に応じて、フーリエ級数、ルジャンドル多項式、ラゲール多項式、エルミート多項式などが選択できる。
この場合、方程式に微積分作用素が含まれていても、微積分の結果は、同一のＯＮＢ（ｘ）の組で表現されなければならない。したがって、高階微分、２重積分３重積分などが実行されても、その結果は、必ず、同一のＯＮＢ（ｘ）の線形結合となる。同様に、拘束条件（初期条件、境界条件）も同一のＯＮＢ（ｘ）の線形結合となる。このように、拘束条件も含めて方程式が代数化されるのでコンピュータによる近似解の求解に適している。
今、ＯＮＢ（ｘ）の数をｎ個とし、拘束条件（初期条件、境界条件）によって関係付けられる線形結合係数の数をｍ個とすると、方程式は（ｎ−ｍ）次元実空間での領域で成立することになる。この（ｎ−ｍ）次元実空間の領域Ωを改めて「許容空間」と呼ぶことにする。
ところで、未知関数のＯＮＢ（ｘ）が意味を持つためには、代数方程式
（ｔ−ｂ１）（ｔ−ｂ２）．．．．．．．．（ｔ−ｂｎ）
の根がすべて実でなければならない。そこで、ｎが小さいときは根の判別式を利用することも出来る。ここで、ｍ個の線形結合係数が、拘束条件（初期条件、境界条件）によって関係付けられているため、代数方程式については、（ｎーｍ）の実根を求めればよいことになる。代数方程式が解けるということは、元の関数方程式が解けるということと同じであるが、代数方程式は、未知方程式ｆ（ｘ）のＯＮＢ（ｘ）展開近似に基づくため、解の正確さまでが保障されるものではない。
（ｔ−ｂ１）（ｔ−ｂ２）．．．．．．．．（ｔ−ｂｎ）
の根がすべて実でなければならない。そこで、ｎが小さいときは根の判別式を利用することも出来る。ここで、ｍ個の線形結合係数が、拘束条件（初期条件、境界条件）によって関係付けられているため、代数方程式については、（ｎーｍ）の実根を求めればよいことになる。代数方程式が解けるということは、元の関数方程式が解けるということと同じであるが、代数方程式は、未知方程式ｆ（ｘ）のＯＮＢ（ｘ）展開近似に基づくため、解の正確さまでが保障されるものではない。
そこで、近似により発生した方程式の誤差を、許容領域の中で最小とすることが考えられる。ここに、微分作用素をＤ、定数をＣ、未知関数ｕ、例えば、方程式をＤｕ＝Ｃとすると、方程式の誤差とは（Ｄｕ−Ｃ）のことであり、Ｄｕはｕの正規直交基底展開により完全に代数化されたものを使用する。本発明では、単なる誤差ではなく二乗誤差を最小とすることを考える。言い換えると、二乗誤差を最小とさせる線形結合係数の組である（ｎ−ｍ）次元実空間の点が求める解となる。なお、二乗誤差に限らず、より高次の誤差を最小としてもよい。そのため、解の精度を検証する基準は、複数の中から選択し得る。そこで、これ以降は、本発明における求解のことを、「（ｎ−ｍ）次元実空間の中で、目的関数Ｌを最小化する解を求めること」と表現することがある。目的関数Ｌは、場合によっては、二乗誤差などとなる。
以下、図面を参照して本発明の実施の形態について説明する。但し、本実施形態に記載されている構成部品の寸法、材質、形状、その相対的配置等に特定的な記載があっても、本発明をそれに限定する趣旨ではない。
図１は、本発明の設計支援装置のブロック図である。この設計支援装置は、方程式入力部１、近似次数入力部２、ＣＰＵ５、許容領域表示部３、解出力部４からなる。
方程式入力部１は、未知関数ｕの方程式、初期条件ＩＣ、境界条件ＢＣを格納し、ＣＰＵ５に出力する。
扱える方程式は、物理学上の微分方程式その他を含め制限がない。その理由は、方程式そのものと境界条件、初期条件、若しくは方程式に伴う拘束条件が、対等に代数化され、解法が容易だからである。
境界条件としては、従来型の境界条件、すなわち、境界の法線方向の未知関数の値を与えるもの（ノイマン型）、境界上の点で未知関数の値を与えるもの（ディリクレ型）、ノイマン型とディリクレ型を併せて値を与えるロバン型がある。常微分方程式の２点境界値問題に関していえば、両端のそれぞれで値を指定する場合は両端ディリクレ型であり、一方端で値を指定する場合は単にディリクレ型と呼ばれ、一方端で傾きを指定する場合は、ノイマン型に相当する。
しかし、本発明で取り扱える境界条件は、例えば、ロバン型にさらに２階の微分項を負荷した境界条件であってもよい。本発明では、いくら高階の境界条件は、選択された正規直交基底の線形結合として代数化されるからである。また、境界上の線分上で境界条件を与えたいときは、線分上で境界条件を積分すればよい。また、初期条件と境界条件は、本発明では、対等に代数化されているので、方程式に対して対等に拘束条件を与えることが特徴である。
近似次数入力部２は、未知関数ｕ（ｘ）の近似次数ｎ（初期条件ＩＣ、境界条件ＢＣにふさわしい正規直交系を選択し、使用する正規直交系基底ＯＮＢ（ｘ）の数をｎとする）を、ユーザが決定し、近似次数ｎをＣＰＵ５に出力する。
許容領域表示部３は、初期条件ＩＣ、境界条件ＢＣを満足する領域（以下、「許容領域」という）を、ＣＰＵ５から受け取り、グラフィック表示する。初期条件・境界条件により関連図けられるＯＮＢ（ｘ）の数がｍであるとき、許容領域は、（ｎ−ｍ）次元の実空間である。
解出力部４は、ＣＰＵ５が計算した解を、プリンタやディスプレイなどに出力する。なお、ユーザは、許容領域を見ながら、方程式ｕの求解を進めていくため、許容領域表示部３とは別に、解出力部４にもディスプレイを備えるのがよい。
こうして、ＣＰＵ５は、許容領域内の停留点（峰の頂上や谷底）やその付近の点などを候補点として、目的関数の値を計算し、多数の候補点での計算結果の中から、目的関数が最小となる候補点を解とする。ここに、目的関数は、近似方程式に生じた誤差、誤差の二乗などであり、許容領域の各点で定美される。
図２は、本実施形態の設計支援装置の動作のフローチャートである。まず、Ｓ１において、未知関数ｕの方程式、初期条件ＩＣ，境界条件ＢＣが方程式入力部１に入力される。
次に、Ｓ２に先立ち、Ｓ２１において、正規直交系が選択され、Ｓ２に進む。ここに、正規直基底ＯＮＢとしては、初期条件ＩＣや境界条件ＢＣに応じて、例えば、フーリエ級数の基底ｓｉｎ（ｊｗｔ）、ｃｏｓ（ｊｗｔ）などを選択すればよい。また、Ｓ２に先立ち、Ｓ２２において、未知関数ｕの近似次数ｎが決定され、Ｓ２に進む。
Ｓ２では、ＣＰＵ５において、未知関数ｕをｎ個の正規直交基底ＯＮＢの線形結合で近似する。線形結合係数の値を求めることが、方程式の解法そのものである。
次に、Ｓ３において、初期条件ＩＣや境界条件ＢＣが、正規直基底ＯＮＢによって表現される。
Ｓ３に続いて、Ｓ３１において、初期条件ＩＣや境界条件ＢＣを満たす許容領域の局面がグラフィック表示される。しかし、グラフィック表示は最大３次元表示であるため、種々の断面が表示される必要がある。
次に、Ｓ３２において、許容領域の中から未知関数ｕの解を与えそうな候補点Ｐｉ（後述するように、ｉは、１から始まり、インクリメントされる）を選択する。
Ｓ４では、候補点Ｐｉにおける目的関数Ｌの値Ｌ（Ｐｉ）が計算される。
次に、Ｓ５において、目的関数値Ｌ（Ｐｉ）が出力される。そして、Ｓ３２に戻り、次の候補点が選択され、次いで、Ｓ４に進んで、次の候補点での目的関数値が計算される。
こうして、複数の（多数の）候補点で、目的関数値を計算し、Ｓ６で、計算結果を集積する。
最後に、Ｓ７において、目的関数の値を最小値たらしめる容量領域の点が決定され、これが方程式の解となる。
［２次元ポアソン方程式］
図３には、ポアソン方程式と境界条件を示す。図４は、実施例１の解法のフローチャートである。図３、図４を参照して、実施例１について説明する。
まず、Ｓ１１において、未知関数ｕ（ｘ、ｙ）が満たす２次元ポアソン方程式と矩形領域内境界での境界条件が入力される。
次に、Ｓ１２において、未知関数ｕを正規直交基底ＯＮＢ（ψｉ、φｉ）で展開する。ここでは、説明の簡単のため正規直交関数を５つ使用する。
ｕ（ｘ、ｙ）＝ Ｘ（ｘ）・Ｙ（ｙ）として、
Ｘ（ｘ）＝（ｐ１ψ１（ｘ）+ｐ２ψ２（ｘ）+ｐ３ψ３（ｘ）） （１）
Ｙ（ｙ）＝（ｋ１φ（ｙ）+ｋ2φ（ｙ）） （２）
ｕ（ｘ、ｙ）は、（ψｉφｊ）で表現できる。ここで、τｋ（ｘ、ｙ）＝（ψｉ（ｘ）・φｊ（ｙ））（ｋ＝１ｔｏ ６）も正規直交系をなすことは簡単に証明できる。なぜなら、（ψｉφｊ）を積分すると、ψｉψｊφｋφｍ（ｉ，ｊ＝１ ｔｏ３；ｋ、ｍ＝１ ｔｏ２）の積分となるので、ψｉ、φｉの直交関数性を使えば、τｋも正規直交性を有するからである。したがって、ｕ（ｘ、ｙ）は、正規直交基底ＯＣＮ（τｋ（ｘ、ｙ））を使って、
ｕ（ｘ、ｙ）＝ａ１τ１（ｘ、ｙ）
+ａ２τ２（ｘ、ｙ）
+ａ３τ３（ｘ、ｙ）
+ａ４τ４（ｘ、ｙ）
+ａ５τ５（ｘ、ｙ）
+ａ６τ６（ｘ、ｙ） （３）
ここで 、
τ１＝ψ１（ｘ）・φ１（ｙ）
τ２＝ψ１（ｘ）・φ２（ｙ）
τ３＝ψ２（ｘ）・φ１（ｙ）
τ４＝ψ２（ｘ）・φ２（ｙ）
τ５＝ψ３（ｘ）・φ１（ｙ）
τ６＝ψ３（ｘ）・φ１（ｙ） （４）
ｕ（ｘ、ｙ）＝ Ｘ（ｘ）・Ｙ（ｙ）として、
Ｘ（ｘ）＝（ｐ１ψ１（ｘ）+ｐ２ψ２（ｘ）+ｐ３ψ３（ｘ）） （１）
Ｙ（ｙ）＝（ｋ１φ（ｙ）+ｋ2φ（ｙ）） （２）
ｕ（ｘ、ｙ）は、（ψｉφｊ）で表現できる。ここで、τｋ（ｘ、ｙ）＝（ψｉ（ｘ）・φｊ（ｙ））（ｋ＝１ｔｏ ６）も正規直交系をなすことは簡単に証明できる。なぜなら、（ψｉφｊ）を積分すると、ψｉψｊφｋφｍ（ｉ，ｊ＝１ ｔｏ３；ｋ、ｍ＝１ ｔｏ２）の積分となるので、ψｉ、φｉの直交関数性を使えば、τｋも正規直交性を有するからである。したがって、ｕ（ｘ、ｙ）は、正規直交基底ＯＣＮ（τｋ（ｘ、ｙ））を使って、
ｕ（ｘ、ｙ）＝ａ１τ１（ｘ、ｙ）
+ａ２τ２（ｘ、ｙ）
+ａ３τ３（ｘ、ｙ）
+ａ４τ４（ｘ、ｙ）
+ａ５τ５（ｘ、ｙ）
+ａ６τ６（ｘ、ｙ） （３）
ここで 、
τ１＝ψ１（ｘ）・φ１（ｙ）
τ２＝ψ１（ｘ）・φ２（ｙ）
τ３＝ψ２（ｘ）・φ１（ｙ）
τ４＝ψ２（ｘ）・φ２（ｙ）
τ５＝ψ３（ｘ）・φ１（ｙ）
τ６＝ψ３（ｘ）・φ１（ｙ） （４）
次に、Ｓ１３において、ポアソン方程式をＯＮＢで代数化する。そのため、Ｘ（ｘ）、Ｙ（ｙ）を微分した結果を、やはり、正規直交関数ψｉ、φｉによって表現すると、
ｄＸ（ｘ）／ｄｘ＝ｑ１ψ１（ｘ）+ｑ２ψ２（ｘ）+ｑ３ψ３（Ｘ）
ｄＹ（ｙ）／ｄｙ＝ｒ１φ１（ｙ）+ｒ２φ２（ｙ） （５）
同様に、
ｄ２Ｘ（ｘ）／ｄｘ２＝ｌ１ψ１（ｘ）+ｌ２ψ２（ｘ）+ｌ３ψ３（Ｘ）
ｄ２Ｙ（ｙ）／ｄｙ２＝ｍ１φ１（ｙ）+ｍ２φ２（ｙ） （６）
ｄＸ（ｘ）／ｄｘ＝ｑ１ψ１（ｘ）+ｑ２ψ２（ｘ）+ｑ３ψ３（Ｘ）
ｄＹ（ｙ）／ｄｙ＝ｒ１φ１（ｙ）+ｒ２φ２（ｙ） （５）
同様に、
ｄ２Ｘ（ｘ）／ｄｘ２＝ｌ１ψ１（ｘ）+ｌ２ψ２（ｘ）+ｌ３ψ３（Ｘ）
ｄ２Ｙ（ｙ）／ｄｙ２＝ｍ１φ１（ｙ）+ｍ２φ２（ｙ） （６）
図３（Ｃ）に示すように、未知関数ｕ（ｘ、ｙ）のｘに関する２次変微分項はＸ“Ｙであり、未知関数ｕ（ｘ、ｙ）のｙに関する２次変微分項はＸＹ”である。したがって、方程式Ｐｏｉｓｓｏｎ Ｅｑ．（ＰＥ）も、τ１〜τ６の線形結合で表される。
ＰＥ＝ｂ１τ１（ｘ、ｙ）
+ｂ２τ２（ｘ、ｙ）
+ｂ３τ３（ｘ、ｙ）
+ｂ４τ４（ｘ、ｙ）
+ｂ５τ５（ｘ、ｙ）
+ｂ６τ６（ｘ、ｙ） （７）
以上のように、ポアソン方程式を正規直交基底ＯＮＢの線形結合で代数化できた。
ＰＥ＝ｂ１τ１（ｘ、ｙ）
+ｂ２τ２（ｘ、ｙ）
+ｂ３τ３（ｘ、ｙ）
+ｂ４τ４（ｘ、ｙ）
+ｂ５τ５（ｘ、ｙ）
+ｂ６τ６（ｘ、ｙ） （７）
以上のように、ポアソン方程式を正規直交基底ＯＮＢの線形結合で代数化できた。
次に、Ｓ１４において、境界条件を代数化する。図３（Ｂ）に示す領域境界Γ上の線分ＡＢ上で、ｕ（ｘ、ｙ）＝ｇ１（ｘ）が与えられたとき、線形結合近似と境界条件との関係を見ていく。
ところで、その準備として、線分ＡＢ上の境界値問題ではなく、１次元２点境界値問題を振り返る。未知関数Ｘ（ｘ）が、ｘ＝０のとき４であり、ｘ＝１のとき１であるとする。
ここに、Ｘ（ｘ）を正規直交基底の線形結合で近似する。
Ｘ（ｘ）≒ｋ１ψ１（ｘ）+ｋ２ψ２（ｘ） （８）
この近似を使って、原点での境界条件を表現すると、
ｋ１ψ１（０）+ｋ２ψ２（０）＝４ （９）
となるが、これは（Ｘ（０）-４）とディラックのデルタ関数δ（ｘ）の内積に等しいことを下式に示す。
ここに、Ｘ（ｘ）を正規直交基底の線形結合で近似する。
Ｘ（ｘ）≒ｋ１ψ１（ｘ）+ｋ２ψ２（ｘ） （８）
この近似を使って、原点での境界条件を表現すると、
ｋ１ψ１（０）+ｋ２ψ２（０）＝４ （９）
となるが、これは（Ｘ（０）-４）とディラックのデルタ関数δ（ｘ）の内積に等しいことを下式に示す。
次に、２次元境界値問題について考える。今、線素ＡＢ上のみで無限大で、線素ＡＢ上以外ではゼロであるデルタ関数であるδＡＢを定義する。２次元境界値問題は、近似を張る領域Ωの部分領域である線素ＡＢ上へのＸ（ｘ）の正射影がｇ１（ｘ）となるということと等価である。線素ＡＢ上の境界条件式ｇ１（ｘ）は境界Γをベースとする関数空間ＩΓひとつのベクトル（関数の点）である。なお、ΩとΓの空間次元はそれぞれ２と１であるが、関数空間としての次元ｄｉｍＩΩ、ｄｉｍＩΓはいずれも無限大である。すなわち、領域Ωで定義されたベクトルＸ（ｘ）はＩΓからははみ出しているが、そのＩΓへの正射影がＩΓの中では最もよい近似であると考えられる。言い換えると、近似ベクトルＸ（ｘ）からｇ１（ｘ）を差し引いたベクトルは、ＩΓをベースとする関数空間のいかなるベクトルとも直交する。
下式は、線素ＡＢ上の関数空間の基底ベクトルδＡＢを表す式である。
これを使用すると、線分ＡＢ上への正射影は、下式となる。
この正射影をＯＰ（ｘ、ｙ）をｇ１（ｘ）に等置したものが線分ＡＢ上の境界条件である。
近似関数ｕ（ｘ、γ）と境界条件式ｇ１（ｘ）との関係は、
ｕ（ｘ、γ）＝ａ１τ１（ｘ、γ）
+ａ２τ２（ｘ、γ）
+ａ３τ３（ｘ、γ）
+ａ４τ４（ｘ、γ）
+ａ５τ５（ｘ、γ）
+ａ６τ６（ｘ、γ）＝ｇ１（ｘ） （１０）
を上述した基底ベクトルδＡＢを使用して積分したものである。
同様に、境界線素ＢＣ，ＣＤ，ＡＤ上からも定数１次式が得られる。
このように、ベクトルＡの成分が満たす４つの拘束条件が得られた。したがって、拘束条件が表す領域をビジュアルに表示することにより、領域の谷底や山頂などに、安定解、不安定解などの存在を予測することが出来る。このように、本発明のポアソン方程式の解法は、極めて見とおしのよいものである。
下式は、線素ＡＢ上の関数空間の基底ベクトルδＡＢを表す式である。
近似関数ｕ（ｘ、γ）と境界条件式ｇ１（ｘ）との関係は、
ｕ（ｘ、γ）＝ａ１τ１（ｘ、γ）
+ａ２τ２（ｘ、γ）
+ａ３τ３（ｘ、γ）
+ａ４τ４（ｘ、γ）
+ａ５τ５（ｘ、γ）
+ａ６τ６（ｘ、γ）＝ｇ１（ｘ） （１０）
を上述した基底ベクトルδＡＢを使用して積分したものである。
同様に、境界線素ＢＣ，ＣＤ，ＡＤ上からも定数１次式が得られる。
このように、ベクトルＡの成分が満たす４つの拘束条件が得られた。したがって、拘束条件が表す領域をビジュアルに表示することにより、領域の谷底や山頂などに、安定解、不安定解などの存在を予測することが出来る。このように、本発明のポアソン方程式の解法は、極めて見とおしのよいものである。
次に、Ｓ１５において、２次元許容領域（ｘｙ平面）において、６次代数方程式の２実根を求める。６次の代数方程式は、
（ｔ−ａ１）（ｔ−ａ２）（ｔ−ａ３）（ｔ−ａ４）（ｔ−ａ５）（ｔ−ａ６）＝０ （１１）
２実根を求めればよい理由は、６つの線形結合定数のうち４つは境界条件により既に関連付けられているためである。
（ｔ−ａ１）（ｔ−ａ２）（ｔ−ａ３）（ｔ−ａ４）（ｔ−ａ５）（ｔ−ａ６）＝０ （１１）
２実根を求めればよい理由は、６つの線形結合定数のうち４つは境界条件により既に関連付けられているためである。
次に、Ｓ１６において、Ｓ１５で求めた２実根により目的関数（方程式に生じた二乗誤差）が最小化されるかどうか判定する。実施例１のポアソン方程式Δｕ＝０の二乗誤差とは（Δｕ−０）２＝（Δｕ）２である。ここに、Δｕは代数化表現を使用する。最小化される場合は既に求めた２実根がポアソン方程式の近似解とするが（Ｓ１７）、最小化されない場合はＳ１５に戻る。
以上、実施例１について説明した。
なお、正規直交基底ＯＮＢの数を増やせば、ベクトルＡ（未知関数ｕ（ｘ、ｙ）の線形結合近似の定数係数を成分とするベクトル）の成分数も同様に増加するが、本発明の解法は、ベクトルＡを決定することに尽きる。すなわち、４つの拘束条件の下で、６次代数方程式（９）の根を決定すれば、方程式の解が求められたことになる。
しかし、方程式の解を確定するためには、ポアソン方程式の誤差又は２乗誤差などを最小化しなければならない。
しかし、方程式の解を確定するためには、ポアソン方程式の誤差又は２乗誤差などを最小化しなければならない。
以上、矩形境界条件について説明したが、境界が矩形でない場合でもかまわない。例えば三角形状、台形状の境界でも簡単に取り扱える。また、矩形との変換則があれば、矩形でない境界を対照とする境界値問題にも対応できる。さらに、例えば、円領域のジューコフスキー変換などが利用できる。また、正規直交基底は、境界の形に即して選択するとよい。また、直線線分ＡＢ上で、ｕ（ｘ、γ）＝ｇ１（ｘ）のような境界条件を与えたが、これに限らず、ｄｕ／ｄｙ＋ｓｕ+ｔ＝０のような境界条件であっても、ｄｕ／ｄｙを代数化できるため、代数化した境界条件を線分ＡＢに沿って積分すればよい。また、境界条件を形成する式の作用素が線形ならば、高階の式（例えば、ｙ“+ｉｙ‘+ｊｙ＝ｋ）であっても代数化できる。なお、非線形作用素であっても扱える可能性がある。
また、定数項ゼロのポアソン方程式について説明したが、定数項があっても手順は全く同様である。また、同様の手順で、ラプラス方程式、波動方程式などの近似解を求めることが出来る。
本発明では、複雑な境界条件を代数化して扱う。本来、境界条件は方程式そのものとは独立したものである。従来の境界条件は、境界での変位指定（ディリクレ型）、対称面(軸)での対称条件（ノイマン型）などが通例であるが、今後は、より複雑な境界条件の指定が必要となる。本発明は、代数化手法により、それに対応することが出来る。
［非線形２階微分方程式］
図５は、実施例２の解法のフローチャートである。
まず、Ｓ２１において、方程式と１次元２点境界条件が入力される（方程式ｙ“（ｘ）＝（２／３）ｙ２（ｘ）を、ｙ（０）＝４、ｙ（１）＝０の下で求解する）。
次に、Ｓ２２において、未知関数ｙ（ｘ）を４つ程度の正規直交基底で近似する。すなわち、
ｙ（ｘ）＝ｂ１ψ１（ｘ）+ｂ２ψ２（ｘ）+ｂ３ψ３（ｘ）+ｂ４ψ４（ｘ）（１２）
と近似する。
ψとして、フーリエ級数、ラゲール多項式、ルジャンドル多項式などが境界条件に対応して使用される。
ｙ（ｘ）＝ｂ１ψ１（ｘ）+ｂ２ψ２（ｘ）+ｂ３ψ３（ｘ）+ｂ４ψ４（ｘ）（１２）
と近似する。
ψとして、フーリエ級数、ラゲール多項式、ルジャンドル多項式などが境界条件に対応して使用される。
次に、Ｓ２３において、方程式を０から１まで積分する。そのため、まず、（ｙ（ｘ））２、を０から１まで積分すると、積分結果Ｉは、
Ｉ＝（２／３）（ｂ１２+ｂ２２+ｂ３２+ｂ４２） （１３）
ｙ“は、正規直交系で表現して、ｙ“＝ｄ１ψ１+ｄ２ψ２+ｄ３ψ３+ｄ４ψ４としてもよいが、ｙ”を０から１まで積分する都合上、
ｙ‘＝ｃ１ψ１（ｘ）+ｃ２ψ２（ｘ）+ｃ３ψ３（ｘ）+ｃ４ψ４（ｘ）
を使用する。
ｙ”の０から１までの定積分は、ｙ’（１）とｙ’（０）の差分である。
ところで、定数係数ｃ１〜ｃ３は、ｂ１〜ｂ３に基づいて生成されたものであり、ｂ１、ｂ２、ｂ３を成分とする列ベクトルに、３行３列の変換行列（Ｒ）を乗算して生成される。したがって、ｙ”の０から１までの定積分Ｉ‘は、
Ｉ‘＝ｃ１（ψ１（１）−ψ１（０））
+ｃ２（ψ２（１）-ψ２（０））
+ｃ３（ψ３（１）-ψ３（０））
+ｃ４（ψ４（１）-ψ４（０）） （１４）
微分方程式は、式（１３）と式（１４）を等置したものである。
Ｉ＝（２／３）（ｂ１２+ｂ２２+ｂ３２+ｂ４２） （１３）
ｙ“は、正規直交系で表現して、ｙ“＝ｄ１ψ１+ｄ２ψ２+ｄ３ψ３+ｄ４ψ４としてもよいが、ｙ”を０から１まで積分する都合上、
ｙ‘＝ｃ１ψ１（ｘ）+ｃ２ψ２（ｘ）+ｃ３ψ３（ｘ）+ｃ４ψ４（ｘ）
を使用する。
ｙ”の０から１までの定積分は、ｙ’（１）とｙ’（０）の差分である。
ところで、定数係数ｃ１〜ｃ３は、ｂ１〜ｂ３に基づいて生成されたものであり、ｂ１、ｂ２、ｂ３を成分とする列ベクトルに、３行３列の変換行列（Ｒ）を乗算して生成される。したがって、ｙ”の０から１までの定積分Ｉ‘は、
Ｉ‘＝ｃ１（ψ１（１）−ψ１（０））
+ｃ２（ψ２（１）-ψ２（０））
+ｃ３（ψ３（１）-ψ３（０））
+ｃ４（ψ４（１）-ψ４（０）） （１４）
微分方程式は、式（１３）と式（１４）を等置したものである。
次に、Ｓ２４において、境界条件を代数化する。代数化された境界条件は、
ｂ１ψ１（１）+ｂ２ψ２（１）+ｂ３ψ３（１）+ｂ４ψ４（１）＝１ （１５）
ｂ１ψ１（０）+ｂ２ψ２（０）+ｂ３ψ３（０）+ｂ４ψ４（０）＝４ （１６）
である。これを０から１まで積分して、下式の境界条件を得る。
ｂ１ψ１（１）+ｂ２ψ２（１）+ｂ３ψ３（１）+ｂ４ψ４（１）＝１ （１５）
ｂ１ψ１（０）+ｂ２ψ２（０）+ｂ３ψ３（０）+ｂ４ψ４（０）＝４ （１６）
である。これを０から１まで積分して、下式の境界条件を得る。
次に、Ｓ２５において、４次代数方程式の２実根を２つの境界条件（１５）（１６）の下で求める。代数方程式は、
（ｔ−ｂ１）（ｔ−ｂ２）（ｔ−ｂ３）（ｔ−ｂ４）＝０
境界条件（１５）、（１６）の下で未知の根は２つに減少する。
（ｔ−ｂ１）（ｔ−ｂ２）（ｔ−ｂ３）（ｔ−ｂ４）＝０
境界条件（１５）、（１６）の下で未知の根は２つに減少する。
なお、境界条件（１１）、（１２）を幾何学的に可視化すれば、その領域の谷底や山頂の存在を手がかりに求解を行ない易くなるのが本発明の特徴である。
次に、Ｓ２６において、残りの２つの線形結合定数係数を求めるには、２乗誤差最小化が使用される。これは、線形結合定数係数で表現されたｙ"とｙ２が満たす方程式ｙ”＝（３／２）ｙ２の２乗誤差Ｌ＝（ｙ”-（３／２）ｙ２）２を最小ならしめる線形結合定数係数が求められる（Ｓ２７）。
以上、実施例２について説明した。
なお、方程式に非線形項（３／２）ｙ２がなければＬは、ｂ1〜ｂ４の２次形式となって、Ｌを滞留させる条件は、ｂ１〜ｂ４の１次方程式系を誘導するから、手順は簡単であるが、Ｌには非線形項があるため、ｂ１〜ｂ４の4次式となる。
また、本方程式で最適界が定まった後は、式（１１）を参照すると、
ｂ１２+ｂ２２+ｂ３２+ｂ４２＝定数
となる（パーセヴアル等式近似）。この意味で、線形結合係数ベクトルＢの先端は、高次元球面上にある。解Ｂの先端は、目的関数Ｌを最小にしており、高次元曲面上で原点を中心とする最も半径の小さな超球と接すると見ても良い。しかし、クランク軸を分析したＭＣＰ手法（特開２００５−９３２３号公報）の場合のように、半径ＳＱＲＴ（ｎ）（クランク軸数は（ｎ+１））の球上にはない。ＭＣＰ手法のときは、球面を固定して拘束面を動かすが、本発明では、拘束面を固定して、球の半径を動かしている。
ｂ１２+ｂ２２+ｂ３２+ｂ４２＝定数
となる（パーセヴアル等式近似）。この意味で、線形結合係数ベクトルＢの先端は、高次元球面上にある。解Ｂの先端は、目的関数Ｌを最小にしており、高次元曲面上で原点を中心とする最も半径の小さな超球と接すると見ても良い。しかし、クランク軸を分析したＭＣＰ手法（特開２００５−９３２３号公報）の場合のように、半径ＳＱＲＴ（ｎ）（クランク軸数は（ｎ+１））の球上にはない。ＭＣＰ手法のときは、球面を固定して拘束面を動かすが、本発明では、拘束面を固定して、球の半径を動かしている。
[Ｍａｔｈｉｅｕ方程式]
本方程式は、Ｈｅｌｍｈｏｌｔｚ方程式を楕円筒座標、回転楕円体座標で解くときに現れる。ここでは、
方程式ｄ２ｙ／ｄｔ２+（1-２ｑｃｏｓ（２ｔ））ｙ＝０
を取り上げる。図６は、実施例３の解法のフローチャートである。
境界条件は、出発時刻ｔ＝０として、
ｙ（０）＝ｙ（π）
ｙ’（０）＝ｙ’（π）
である。
方程式ｄ２ｙ／ｄｔ２+（1-２ｑｃｏｓ（２ｔ））ｙ＝０
を取り上げる。図６は、実施例３の解法のフローチャートである。
境界条件は、出発時刻ｔ＝０として、
ｙ（０）＝ｙ（π）
ｙ’（０）＝ｙ’（π）
である。
この方程式を解くために、まず、Ｓ３１において、方程式、拘束条件が入力される。
次に、Ｓ３２において、未知関数を４つの基底ＯＮＢで近似する。
未知関数ｙの正規直交基底展開は、
ｙ＝ｂ１ψ１（ｔ）+ｂ２ψ２（ｔ）+ｂ３ψ３（ｔ）+ｂ４ψ４（ｔ） （１７）
未知関数ｙの正規直交基底展開は、
ｙ＝ｂ１ψ１（ｔ）+ｂ２ψ２（ｔ）+ｂ３ψ３（ｔ）+ｂ４ψ４（ｔ） （１７）
次に、Ｓ３３において、方程式を代数化する。ｙ‘、ｙ“は、式（１５）と同様に、
ｙ’＝ｃ１ψ１（ｔ）+ｃ２ψ２（ｔ）+ｃ３ψ３（ｔ）+ｃ４ψ４（ｔ） （１８）
ｙ”＝ｄ１ψ１（ｔ）+ｄ２ψ２（ｔ）+ｄ３ψ３（ｔ）+ｄ４ψ４（ｔ） （１９）
ここに、ｂ１〜ｂ４を成分とする列ベクトルを（Ｂ）、変換行列を（Ｋ１）とすると、
ｃ１〜ｃ４を成分とする列ベクトル（Ｃ）＝（Ｋ１）（Ｂ）であり、
ｄ１〜ｄ４を成分とする列ベクトル（Ｄ）＝（Ｋ２）（Ｂ）である。
式（１５）及び式（１７）を用いると、Ｍａｔｈｉｅｕ方程式が代数化される。
ｙ’＝ｃ１ψ１（ｔ）+ｃ２ψ２（ｔ）+ｃ３ψ３（ｔ）+ｃ４ψ４（ｔ） （１８）
ｙ”＝ｄ１ψ１（ｔ）+ｄ２ψ２（ｔ）+ｄ３ψ３（ｔ）+ｄ４ψ４（ｔ） （１９）
ここに、ｂ１〜ｂ４を成分とする列ベクトルを（Ｂ）、変換行列を（Ｋ１）とすると、
ｃ１〜ｃ４を成分とする列ベクトル（Ｃ）＝（Ｋ１）（Ｂ）であり、
ｄ１〜ｄ４を成分とする列ベクトル（Ｄ）＝（Ｋ２）（Ｂ）である。
式（１５）及び式（１７）を用いると、Ｍａｔｈｉｅｕ方程式が代数化される。
次に、Ｓ３４において、境界条件が代数化される。境界条件は、式（１７）、式（１８）にｔ＝０、πを代入すれば得られるが、ベクトル形式では、
（ψＴ（０））（Ｂ）＝（ψＴ（π））（Ｂ） （２０）
（ψＴ（０））（Ｋ１）（Ｂ）＝（ψＴ（π））（Ｋ１）（Ｂ） （２１）
これを０からπまで積分したものが境界条件である。
（ψＴ（０））（Ｂ）＝（ψＴ（π））（Ｂ） （２０）
（ψＴ（０））（Ｋ１）（Ｂ）＝（ψＴ（π））（Ｋ１）（Ｂ） （２１）
これを０からπまで積分したものが境界条件である。
次に、Ｓ３５，Ｓ３６，Ｓ３７において方程式の二乗誤差を最小にする代数方程式の解を求める。代数方程式は、
（ｔ-ｂ１）（ｔ-ｂ２）（ｔ−ｂ３）（ｔ−ｂ４）
となる。この４次代数方程式が、拘束条件（２０）、（２１）の下で解かれる。したがって、２つの実根を求めればよい。
目的関数Ｌは、ベクトル形式で、
Ｌ＝∫（ｙ”+（1-２ｑｃｏｓ２ｔ）ｙ）２ｄｔ
＝（Ｂ）Ｔ（Ｋ１）Ｔ・（∫（（ψ（ｔ））・（ψ（ｔ））Ｔｄｔ）・（Ｋ２）（Ｂ）
+２（Ｂ）Ｔ（Ｋ２）Ｔ・（∫（（ψ（ｔ））・（１-２ｑｃｏｓ２ｔ）・（ψ（ｔ））Ｔｄｔ）・（Ｂ）+（Ｂ）Ｔ・（∫（（ψ（ｔ））・（１−２ｑｃｏｓ２ｔ）・（ψ（ｔ））Ｔｄｔ）・（Ｂ）
となる。積分は、０からπまでである。目的関数Ｌは、パラメトリック発振であろうが、時間ｔの関数であろうが、（Ｂ）の２次曲面形式になる。
（ｔ-ｂ１）（ｔ-ｂ２）（ｔ−ｂ３）（ｔ−ｂ４）
となる。この４次代数方程式が、拘束条件（２０）、（２１）の下で解かれる。したがって、２つの実根を求めればよい。
目的関数Ｌは、ベクトル形式で、
Ｌ＝∫（ｙ”+（1-２ｑｃｏｓ２ｔ）ｙ）２ｄｔ
＝（Ｂ）Ｔ（Ｋ１）Ｔ・（∫（（ψ（ｔ））・（ψ（ｔ））Ｔｄｔ）・（Ｋ２）（Ｂ）
+２（Ｂ）Ｔ（Ｋ２）Ｔ・（∫（（ψ（ｔ））・（１-２ｑｃｏｓ２ｔ）・（ψ（ｔ））Ｔｄｔ）・（Ｂ）+（Ｂ）Ｔ・（∫（（ψ（ｔ））・（１−２ｑｃｏｓ２ｔ）・（ψ（ｔ））Ｔｄｔ）・（Ｂ）
となる。積分は、０からπまでである。目的関数Ｌは、パラメトリック発振であろうが、時間ｔの関数であろうが、（Ｂ）の２次曲面形式になる。
以上、実施例３について説明した。
なお、境界条件は、
ｙ（ｔ＋π）＝ｅｘｐ（ｉπν）ｙ（ｔ）
で与えられることがあり（Ｆｌｏｑｕｅｔの解）、この場合の境界条件のベクトル表現は、
（ψ（π））Ｔ・（Ｂ）＝ｅｘｐ（ｉπν）・（ψ（０））Ｔ・（Ｂ）
（ψ（π））Ｔ・（Ｋ１）・（Ｂ）＝ｅｘｐ（ｉπν）・（ψ（０））Ｔ・（Ｂ）
となる。ｂ１〜ｂ４の間に成立する２つの１次形式の係数は一般に複素数となるが、方程式（１８）の根がすべて実根であることに変わりはない。
ｙ（ｔ＋π）＝ｅｘｐ（ｉπν）ｙ（ｔ）
で与えられることがあり（Ｆｌｏｑｕｅｔの解）、この場合の境界条件のベクトル表現は、
（ψ（π））Ｔ・（Ｂ）＝ｅｘｐ（ｉπν）・（ψ（０））Ｔ・（Ｂ）
（ψ（π））Ｔ・（Ｋ１）・（Ｂ）＝ｅｘｐ（ｉπν）・（ψ（０））Ｔ・（Ｂ）
となる。ｂ１〜ｂ４の間に成立する２つの１次形式の係数は一般に複素数となるが、方程式（１８）の根がすべて実根であることに変わりはない。
[Ｈｉｌｌ方程式]
ｍｄ２ｙ／ｄｔ２+ｋ（ｔ）ｙ＝０のように、ばね定数が時間的に変化し、ｋ（ｔ+Ｔ）＝ｋ（ｔ）である場合にＨＩｌｌ方程式が登場する。図７は、実施例４の解法のフローチャートである。
まず、Ｓ４１において、Ｈｉｌｌ方程式及び拘束条件を入力する。
次に、Ｓ４２において、方程式系を生成する。すなわち、ｙ＝ｙ１、ｄｙ／ｄｔ＝ｙ２とおき、ｙ１、ｙ２を成分とする列ベクトル（Ｘ）に対して、Ｈｉｌｌ方程式のベクトル表現は、ｙ＝ｙ１、ｄｙ／ｄｘ＝ｙ２とおくと、
次に、Ｓ４３において、方程式系を、それぞれ２基底で近似する。
ｙ１＝ｂ１φ１（ｔ）+ｂ２φ２（ｔ） （２２）
ｙ２＝ｃ１φ１（ｔ）+ｃ２φ２（ｔ） （２３）
ｙ１＝ｂ１φ１（ｔ）+ｂ２φ２（ｔ） （２２）
ｙ２＝ｃ１φ１（ｔ）+ｃ２φ２（ｔ） （２３）
次に、Ｓ４４において、境界条件を代数化する。境界条件は、
ｙ１（Ｔ）＝ｙ１（０）
ｙ２（Ｔ）＝ｙ２（０）
ｄｙ１／ｄｔ（０）＝ｄｙ２／ｄｔ（０）
ｄｙ１／ｄｔ（Ｔ）＝ｄｙ２／ｄｔ（Ｔ） （２４）
式（２１）、式（２２）を、４つの境界条件（２４）に代入すれば、境界条件が代数化され、それを０からＴまで積分して境界条件を得る。
ｙ１（Ｔ）＝ｙ１（０）
ｙ２（Ｔ）＝ｙ２（０）
ｄｙ１／ｄｔ（０）＝ｄｙ２／ｄｔ（０）
ｄｙ１／ｄｔ（Ｔ）＝ｄｙ２／ｄｔ（Ｔ） （２４）
式（２１）、式（２２）を、４つの境界条件（２４）に代入すれば、境界条件が代数化され、それを０からＴまで積分して境界条件を得る。
次に、Ｓ４５からＳ４７において、４次代数方程式
（ｔ-ｂ１）（ｔ-ｂ２）（ｔ−ｃ１）（ｔ−ｃ２）＝０
の４つの実根が求められる。なお、目的関数は、下式で表わされる。
（ｔ-ｂ１）（ｔ-ｂ２）（ｔ−ｃ１）（ｔ−ｃ２）＝０
の４つの実根が求められる。なお、目的関数は、下式で表わされる。
以上、実施例４について説明した。
なお、実施例４では、４つの線形結合係数（ｂ１、ｂ２、ｃ１、ｃ２）が、４つの拘束関係（２４）で関連付けられている。これでは、拘束が強すぎるため、３基底展開、４基底展開を考えるとよい。また、滑らかな接続条件
（ｄｙ１／ｄｔ（Ｔ）＝ｄｙ１／ｄｔ（０）、ｄｙ２ｄｔ（Ｔ）＝ｄｙ２／ｄｔ（０））
を課せば、基底数をそれだけ増やしやすくなる。
（ｄｙ１／ｄｔ（Ｔ）＝ｄｙ１／ｄｔ（０）、ｄｙ２ｄｔ（Ｔ）＝ｄｙ２／ｄｔ（０））
を課せば、基底数をそれだけ増やしやすくなる。
［積分問題］
ｙ‘＝ｆ（ｘ）、ｙ（０）＝Ｃ の解はｙ=∫ｆ（ｔ）ｄｔ+Ｃ（積分は０からｘまでであるが、これを本発明にしたがって求解する。図８は、実施例５の解法のフローチャートである。
まず、Ｓ５１において、方程式及び境界条件（拘束条件）が入力される。
次に、Ｓ５２において、未知関数ｙを３基底で近似する。
ｙ＝ｂ１ψ１（ｘ）+ｂ２ψ２（ｘ）+ｂ３ψ３（ｘ） （２５）
ｙ＝ｂ１ψ１（ｘ）+ｂ２ψ２（ｘ）+ｂ３ψ３（ｘ） （２５）
次に、Ｓ５２において、式（２５）を方程式に代入して代数化し方程式を積分する。
次に、Ｓ５４において、拘束条件を代数化する。代数化結果は、
ｂ１ψ１（０）+ｂ２ψ２（０）+ｂ３ψ３（０）＝Ｃ （２６）
これが境界条件である。
ｂ１ψ１（０）+ｂ２ψ２（０）+ｂ３ψ３（０）＝Ｃ （２６）
これが境界条件である。
次に、Ｓ５５〜５７において、ｂ１〜ｂ３が満たす代数方程式の実根を求める。
代数方程式は、
（ｔ-ｂ１）（ｔ−ｂ２）（ｔ−ｂ３）
＝ｔ３−（ｂ１+ｂ２+ｂ３）ｔ２+（ｂ１ｂ２+ｂ２ｂ３+ｂ３ｂ１）ｔ+ｂ１ｂ２ｂ３
（２７）
この代数方程式に、拘束条件式（２６）から、例えば、
ｂ３＝（Ｃ−ｂ１ψ１（０）−ｂ２ψ２（０））／ψ３
を代入すれば、求める根はｂ１、ｂ２だけとなる。したがって、２次元実空間である許容領域の中で、２つの実根を探せばよい。
Ｓ５６では、目的関数Ｌは、
２乗誤差（ｂ１ψ１（ｘ）+ｂ２ψ２（ｘ）+ｂ３ψ３（ｘ）−ｆ（ｘ））２
を１〜αまで積分したものである。αは有限であるが、必要に応じて充分大なる値とする。
代数方程式は、
（ｔ-ｂ１）（ｔ−ｂ２）（ｔ−ｂ３）
＝ｔ３−（ｂ１+ｂ２+ｂ３）ｔ２+（ｂ１ｂ２+ｂ２ｂ３+ｂ３ｂ１）ｔ+ｂ１ｂ２ｂ３
（２７）
この代数方程式に、拘束条件式（２６）から、例えば、
ｂ３＝（Ｃ−ｂ１ψ１（０）−ｂ２ψ２（０））／ψ３
を代入すれば、求める根はｂ１、ｂ２だけとなる。したがって、２次元実空間である許容領域の中で、２つの実根を探せばよい。
Ｓ５６では、目的関数Ｌは、
２乗誤差（ｂ１ψ１（ｘ）+ｂ２ψ２（ｘ）+ｂ３ψ３（ｘ）−ｆ（ｘ））２
を１〜αまで積分したものである。αは有限であるが、必要に応じて充分大なる値とする。
以上、実施例５について説明した。
［初期値問題］
ここでは、次の方程式を解く。図９は、実施例６の解法のフローチャートである。方程式と初期条件は、下式である。
ここでは、次の方程式を解く。図９は、実施例６の解法のフローチャートである。方程式と初期条件は、下式である。
まず、Ｓ６１において、方程式及び初期条件（拘束条件）が入力される。
次に、Ｓ６２において、ｕ＝Ｔ（ｔ）Ｘ（ｘ）と変数分離し、Ｔ（ｘ）及びＸ（ｘ）をそれぞれ２項の基底で近似する。
Ｔ（ｔ）＝ａ１φ１（ｔ）+ａ２φ２（ｔ）
Ｘ（ｘ）＝ｃ１ζ１（ｘ）+ｃ２ζ２（ｘ）
すると、ｕ＝ａ１ｃ１φ１（ｔ）ζ１（ｘ）
+ａ１ｃ２φ１（ｔ）ζ２（ｘ）
+ａ２ｃ１φ２（ｔ）ζ１（ｘ）
+ａ２ｃ２φ２（ｔ）ζ２（ｘ）
＝ｂ１τ１（ｔ、ｘ）
+ｂ２τ１（ｔ、ｘ）
+ｂ３τ１（ｔ、ｘ）
+ｂ４τ４（ｔ、ｘ） （２８）
これは、ポアソン方程式の場合と同様の手順である。
Ｔ（ｔ）＝ａ１φ１（ｔ）+ａ２φ２（ｔ）
Ｘ（ｘ）＝ｃ１ζ１（ｘ）+ｃ２ζ２（ｘ）
すると、ｕ＝ａ１ｃ１φ１（ｔ）ζ１（ｘ）
+ａ１ｃ２φ１（ｔ）ζ２（ｘ）
+ａ２ｃ１φ２（ｔ）ζ１（ｘ）
+ａ２ｃ２φ２（ｔ）ζ２（ｘ）
＝ｂ１τ１（ｔ、ｘ）
+ｂ２τ１（ｔ、ｘ）
+ｂ３τ１（ｔ、ｘ）
+ｂ４τ４（ｔ、ｘ） （２８）
これは、ポアソン方程式の場合と同様の手順である。
次に、Ｓ６３において、初期条件を代数化し、０からαまで積分し、初期条件が下式のとおり得られる。
αは有限であるが、必要に応じて充分大なる値とする。
次に、Ｓ６４において、式（２７）を使って、ポアソン方程式の場合と同様に方程式を代数化する。
次に、Ｓ６５〜Ｓ６７において、実施例１から５と同様に、方程式の二乗誤差を最小とさせる、４次代数方程式の３根が求められる。拘束条件がひとつあるからである。
目的関数は、下式のとおりである。
αは有限であるが、必要に応じて充分大なる値とする。
目的関数は、下式のとおりである。
［アトラクタの設計］
実施例５のＨｉｌｌ方程式の取り扱いを応用してアトラクタを設計を行うことが出来る。アトラクタは、２度と同じ場所を通過しないが、有界領域の中で周期的軌跡の運動を継続する。図１０は許容領域Σ内の閉曲線Ｃｐｃから１つのアトラクタを生成する概念図、図１１はアトラクタを生成するためのフローチャートである。
まず、Ｓ７１において、Ｈｉｌｌ方程式の目的関数Ｌが入力される。
次に、Ｓ７２において、ｂ１、ｃ１、を一定として、許容領域中の閉曲線を求める。
次に、Ｓ７３において、ｂ１（ｔ）、ｃ１（ｔ）として、目的関数Ｌが時間ｔに依存しないようにｋ（ｔ）を定めると、アトラクタが生成される（Ｓ７４）。他のアトラクタを生成するには、Ｓ７３に戻る。
本発明は、物理学、化学、機械工学などにおいて解を求めることが必要な、複雑な境界条件や初期条件を伴う３次元物理空間における偏微分方程式や積分方程式、さらには漸化式方程式や非線形連立方程式を扱う設計支援装置に利用可能である。ここに、設計支援で問題となる支配方程式としては、弾性の基礎式（大変形構造問題）、放物型偏微分方程式（熱や物質の拡散）、オイラー方程式（ポテンシャル流れ）、ナヴィエストークス方程式（非線形偏微分方程式、非圧縮性及び圧縮性流れ、化学反応解析、燃焼解析の波動方程式、制御、最適化（目的関数の極大化・極小化）問題などがある。
また、ナビエ・ストークス方程式などの３次元非線形方程式にいたるまでもなく、２点境界値問題の範囲でも、Ｍａｔｈｉｅｕ方程式、Ｈｉｌｌ方程式などの振動問題の解法に利用可能である。例えば、Ｍａｔｈｉｅｕ方程式を適用する平歯車などの歯車伝達機構の設計、その他、高圧送電線の風力振動対策、高層ビルクレーン吊鎖の安定化・作業環境改善策、平歯車などの歯車伝達機構の設計、鉄道・電車・電線の張力変動対策、に応用可能である。
産業分野を具体的に特定すると、極限設計部品に関しては、長大橋、タンク、容器、エンジン部品（ピストン、フレーム、カバー、クランク、テンションボルト、タービン、ブレード、クリスマスツリー、ロータ、ハブ）、ノッキングデトネーションを回避したガスエンジン、車室、ボイラ配管、ドラム、炉壁、ダブルハル、プロペラ、エンジン燃焼室廻り部品（熱応力・材料耐性の耐火面、給配気弁、ピストンクラウン）に利用可能である。
また、排出ガス成分を充分に考慮したストーカー炉・回転炉、炉壁溶損予測を施した焼却炉、設置状況・周囲地形を配慮しリサーキュレーションを回避したエアフィンクーラー、白煙発生対策を施した冷却塔・ラジエータ、スタックレイン回避方法を適用した煙突、密度差・噴流混合を勘案したＬＮＧ・ＬＰＧタンクや液化ガスタンク、スロッシングを抑制したタンク構造、応力解析を施す航空機部品（翼、ブレード）、貫流ファン、扇風機翼、トンネル排気ファン、ロータリーキルン、溶鉱炉、石炭サイロに利用可能である。
観点を変えると、ステファン問題（境界条件が時間変化する問題、例えば、物体間の接触非線形問題、地下ＬＮＧタンク稼動後の土壌冷結・タンク浮上）解析ビジネス、土壌浸透圧問題を含む土木解析ビジネス、電磁波障害解析ビジネス、雷雲・落雷解析ビジネス、気象予測解析ビジネス、海洋・湖沼・プールの汚染防止解析ビジネスなどに利用可能である。
また、カオス・アトラクタ・リミットサイクルの軌跡をコンピュータで創生することに利用可能である。
１ 方程式入力部
２ 近似次数入力部
３ 許容領域表示部
４ 解出力部
５ ＣＰＵ
２ 近似次数入力部
３ 許容領域表示部
４ 解出力部
５ ＣＰＵ
Claims (9)
- 工学分野での設計に際し、支配方程式及びその初期条件及び／又は境界条件からなる拘束条件が設定され、未知関数uに関する該支配方程式を、コンピュータを用いて解く設計支援装置であって、前記支配方程式は作用素Ｄ及び定数Ｃとして（Ｄu−Ｃ=０）なる形に表現され、
該未知関数uは、実ｎ次元のベクトル変数の実関数とし、前記拘束条件の数をn個より少ないｍ個とする支配方程式の解を求める設計支援装置において、
前記コンピュータは、境界条件や初期条件を伴う３次元物理空間における偏微分方程式や積分方程式、さらには漸化式方程式や非線形連立方程式を扱う設計支援装置に組み込まれており、
前記コンピュータは、
前記支配方程式、及び、その初期条件及び／又は境界条件、若しくは前記方程式に伴う拘束条件を入力する方程式入力部と、
前記未知関数uを、前記実n次元に対応する有限個数ｎの規格直交基底の線形結合で近似するときの前記ｎを入力する近似次数入力部と
前記近似により生じた前記支配方程式の誤差であって、ｍａｘノルム、１−ノルム、２−ノルム、ミンコフスキーノルム、ヘルダーノルム、又は距離の公理を満足するノルムで表わされる誤差‖Ｄu−Ｃ‖を計算するＣＰＵと、
前記拘束条件によって関係付けられた前記線形結合の結合係数の数をｍ個とするとき、（ｎ−ｍ）次元実空間の領域を許容領域として表示する許容領域表示部と、
前記方程式の解を出力する解出力部とを備え、
前記ＣＰＵに前記方程式入力部より、初期条件及び／又は境界条件、若しくは前記方程式に伴う拘束条件を入力された前記支配方程式を入力して、前記近似次数入力部より入力された前記許容領域内の任意の点を候補点として、前記ＣＰＵにおいて前記候補点における前記誤差を計算して、
前記各候補点における前記誤差‖Ｄu−Ｃ‖が最小ならしめる前記許容領域の候補点を前記方程式の解として求め、
該方程式の解を前記解出力部を介して記拘束条件によって関係付けられた前記線形結合の結合係数の数ｍを、（ｎ−ｍ）次元実空間の領域が、３次元以下の許容領域となるように設定し、
前記許容領域表示部は、前記許容領域が２次元以下のときは２次元以下のグラフィック表示を行い、許容領域が３次元以上のときは３次元のグラフィック表示を行い、又は、２次元領域に階層化された複数のグラフとして表示し、許容領域の曲面が３次元以上のときは、３次元のグラフィック表示を行い、又は、２次元領域に階層化された立体グラフ又は等高線図として表示することを特徴とする設計支援装置。 - 請求項１において、
前記支配方程式は、微分作用素、積分作用要素、漸化式の写像作用素、又は非線形方程式における非線形作用素をＤ、定数をＣ、未知関数を前記ｕとするとき、Ｄｕ＝Ｃであり、前記ｕを前記線形結合で近似し且つ前記方程式全体を代数化した後のＤｕ＝Ｃを使って、前記誤差である前記ノルム||Ｄｕ−Ｃ||を計算することを特徴とする設計支援装置。 - 請求項１において、
前記ＣＰＵは、前記線形結合の係数ｂｉ（ｉは１から前記ｎまで）を根とする代数方程式において、
前記根のすべてを実数に限定することを特徴とする設計支援装置。 - 工学分野での設計に用い、支配方程式及びその初期条件及び／又は境界条件からなる拘束条件が設定され、未知関数uに関する該支配方程式を解く為に、前記支配方程式は作用素Ｄ及び定数Ｃとして（Ｄu−Ｃ＝０）なる形に表現され、
該未知関数ｕは、実ｎ次元のベクトル変数の実関数とし、前記拘束条件の数をn個より少ないｍ個とする支配方程式の解を求めるＣＰＵを具えたコンピュータにおいて、
前記コンピュータは、前記支配方程式、及び、その初期条件及び／又は境界条件を入力する方程式入力部と、
前記未知関数を、前記実n次元に対応する有限個数ｎの規格直交基底の線形結合で近似するときの前記ｎを入力する近似次数入力部と、
前記近似により生じた前記支配方程式の誤差であって、ｍａｘノルム、１−ノルム、２−ノルム、ミンコフスキーノルム、ヘルダーノルム、又は距離の公理を満足するノルムで表わされる誤差‖Ｄu−Ｃ‖を計算するＣＰＵと、
前記初期条件及び／又は境界条件によって関係付けられた前記線形結合の結合係数の数をｍ個とするとき、（ｎ−ｍ）次元実空間の領域を許容領域として表示する許容領域表示部とを具えるとともに、
前記ＣＰＵに組み込まれた前記支配方程式の解を出力する解出力部として機能させるコンピュータプログラムに基づいて
前記ＣＰＵで前記方程式入力部より、初期条件及び／又は境界条件、若しくは前記方程式に伴う拘束条件を入力された前記支配方程式を入力して、前記近似次数入力部より入力された前記許容領域内の任意の点を候補点として、前記候補点における前記誤差を計算して、
前記各候補点における前記誤差‖Ｄu−Ｃ‖が最小ならしめる前記許容領域の候補点を前記方程式の解として求め
該方程式の解を前記解出力部を介して記拘束条件によって関係付けられた前記線形結合の結合係数の数ｍを、（ｎ−ｍ）次元実空間の領域が、３次元以下の許容領域となるように設定し、
前記許容領域表示部は、前記許容領域が２次元以下のときは２次元以下のグラフィック表示を行い、許容領域が３次元以上のときは３次元のグラフィック表示を行い、又は、２次元領域に階層化された複数のグラフとして表示し、許容領域の曲面が３次元以上のときは、３次元のグラフィック表示を行い、又は、２次元領域に階層化された立体グラフ又は等高線図として表示することを特徴とするコンピュータ。 - 請求項４において、
前記コンピュータプログラムの支配方程式は、微分作用素、積分作用要素、漸化式の写像作用素、又は非線形方程式における非線形作用素をＤ、定数をＣ、未知関数を前記ｕとするとき、Ｄｕ＝Ｃであり、前記ｕを前記線形結合で近似し且つ前記方程式全体を代数化した後のＤｕ＝Ｃを使って、前記誤差である前記ノルム||Ｄｕ−Ｃ||を計算することを特徴とするコンピュータ。 - 請求項４において、
前記ＣＰＵに組み込まれたコンピュータプログラムは、前記線形結合の係数ｂｉ（ｉは１から前記ｎまで）を根とする代数方程式において、
前記根のすべてを実数に限定することを特徴とするコンピュータ。 - 請求項４において、
前記ＣＰＵに組み込まれたコンピュータプログラムは
未知関数ｕ、初期条件、境界条件を入力する第１手順と、
前記ｕを変数分離し、各変数の関数を、各々異なる正規直交基底の線形結合で近似する第２手順と、
前記方程式を前記正規直交基底で表わし、その表現に含まれる線形結合定数を根とする代数方程式を生成する第３手順と、
前記初期条件、前記境界条件を前記正規直交基底で表わし、１以上の拘束条件を得る第４手順と、
前記代数方程式の根を前記拘束条件の下で求める第５手順と、
前記根が前記代数方程式の前記誤差を最小ならしめるか否かを判定する第６手順とを含み、
前記誤差が最小であると判定されれば、前記根を前記方程式の解とし、
前記誤差が最小でないと判定されれば、第５手順に戻ることを特徴とするコンピュータ。 - 請求項４において、
前記ＣＰＵに組み込まれたコンピュータプログラムは
未知関数ｕ、初期条件、境界条件を入力する第１手順と、
前記ｕを正規直交基底の線形結合で近似する第２手順と、
前記方程式を定積分し、前記定積分に含まれる線形結合定数を根とする代数方程式を生成する第３手順と、
前記初期条件、前記境界条件を前記正規直交基底で表わし、１以上の拘束条件を得る第４手順と、
前記代数方程式の根を前記拘束条件の下で求める第５手順と、
前記根が前記代数方程式の前記誤差を最小ならしめるか否かを判定する第６手順とを含み、
前記誤差が最小であると判定されれば、前記根を前記方程式の解とし、
前記誤差が最小でないと判定されれば、第５手順に戻ることを特徴とするコンピュータ。 - 請求項４において、
前記ＣＰＵに組み込まれたコンピュータプログラムは
未知関数ｕ、初期条件、境界条件を入力する第１手順と、
前記ｕを正規直交基底の線形結合で近似する第２手順と、
前記線形結合に含まれる線形結合定数を根とする代数方程式を生成する第３手順と、
前記初期条件、前記境界条件を前記正規直交基底で表わし、１以上の拘束条件を得る第４手順と、
前記代数方程式の根を前記拘束条件の下で求める第５手順と、
前記根が前記代数方程式の前記誤差を最小ならしめるか否かを判定する第６手順とを含み、
前記誤差が最小であると判定されれば、前記根を前記方程式の解とし、
前記誤差が最小でないと判定されれば、第５手順に戻ることを特徴とするコンピュータ。
Priority Applications (3)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| JP2006008057A JP4599301B2 (ja) | 2006-01-16 | 2006-01-16 | 設計支援装置及びコンピュータプログラム |
| US11/646,285 US20070179762A1 (en) | 2006-01-16 | 2006-12-28 | Design aiding apparatus and computer program |
| EP06027001.4A EP1808775A3 (en) | 2006-01-16 | 2006-12-28 | Design aiding apparatus and computer program |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| JP2006008057A JP4599301B2 (ja) | 2006-01-16 | 2006-01-16 | 設計支援装置及びコンピュータプログラム |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| JP2007188454A JP2007188454A (ja) | 2007-07-26 |
| JP4599301B2 true JP4599301B2 (ja) | 2010-12-15 |
Family
ID=37909778
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| JP2006008057A Expired - Fee Related JP4599301B2 (ja) | 2006-01-16 | 2006-01-16 | 設計支援装置及びコンピュータプログラム |
Country Status (3)
| Country | Link |
|---|---|
| US (1) | US20070179762A1 (ja) |
| EP (1) | EP1808775A3 (ja) |
| JP (1) | JP4599301B2 (ja) |
Families Citing this family (11)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| JP4778486B2 (ja) * | 2007-07-10 | 2011-09-21 | 日本電信電話株式会社 | 多項式間距離算出装置、方法、プログラムおよび記録媒体、一変数最近実多項式算出装置、方法、プログラムおよび記録媒体 |
| US8238675B2 (en) * | 2008-03-24 | 2012-08-07 | Microsoft Corporation | Spectral information recovery for compressed image restoration with nonlinear partial differential equation regularization |
| JP5282649B2 (ja) * | 2008-09-25 | 2013-09-04 | 富士通株式会社 | レイアウト評価装置、レイアウト評価プログラム、ダミールール生成装置及びダミールール生成プログラム |
| JP5645463B2 (ja) * | 2010-04-30 | 2014-12-24 | 修三 前田 | ３次元形状を３次元形状へ変換する方法 |
| JP5227384B2 (ja) * | 2010-10-12 | 2013-07-03 | 住友ゴム工業株式会社 | 構造格子を用いたシミュレーション方法 |
| CN103615739B (zh) * | 2013-11-27 | 2015-11-18 | 广东电网公司电力科学研究院 | 燃烧锅炉运行控制方法与系统 |
| CN103646143B (zh) * | 2013-12-17 | 2016-04-06 | 辽宁石油化工大学 | 一种火灾环境下lpg储罐涂层防护系统性能分析的方法 |
| US9256709B2 (en) * | 2014-02-13 | 2016-02-09 | Taiwan Semiconductor Manufacturing Company, Ltd. | Method for integrated circuit mask patterning |
| CN108491646A (zh) * | 2018-03-28 | 2018-09-04 | 青岛理工大学 | 一种隧道近距穿越重要建筑物爆破引起振动损伤的鉴定方法 |
| WO2020098704A1 (zh) * | 2018-11-13 | 2020-05-22 | 苏州润迈德医疗科技有限公司 | 基于造影图像获取血管评定参数的方法、装置及系统 |
| CN112036039B (zh) * | 2020-09-01 | 2022-12-23 | 内蒙古科技大学 | 一种共轭传热材料热性能的高精度数值分析方法 |
Citations (5)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| JPH01145723A (ja) * | 1987-08-28 | 1989-06-07 | Hitachi Ltd | プログラム生成方法 |
| JPH06148038A (ja) * | 1992-11-09 | 1994-05-27 | Hitachi Ltd | 振動解析装置 |
| JP2005009323A (ja) * | 2003-06-16 | 2005-01-13 | Mitsubishi Heavy Ind Ltd | 多気筒内燃機関、そのクランク軸のクランク配置決定方法 |
| JP2005202634A (ja) * | 2004-01-15 | 2005-07-28 | Keio Gijuku | 六面体有限要素数値解析方法 |
| JP2005267214A (ja) * | 2004-03-18 | 2005-09-29 | Canon Inc | 有限要素法による制御プログラムおよび記憶媒体 |
Family Cites Families (6)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US5313406A (en) * | 1992-06-10 | 1994-05-17 | Temet Instruments Oy | Procedures for analyzing multicomponent FT-IR spectra for unknown mixtures of gases |
| US5966524A (en) * | 1997-07-24 | 1999-10-12 | Lucent Technologies Inc. | 3-D electromagnetic infinite element |
| US6963824B1 (en) * | 1999-02-19 | 2005-11-08 | Davidson Joseph K | Method and apparatus for geometric variations to integrate parametric computer aided design with tolerance analyses and optimization |
| JP2004527860A (ja) * | 2001-05-25 | 2004-09-09 | パラメトリック・オプティミゼーション・ソリューションズ・リミテッド | 改善されたプロセス制御 |
| JP3861012B2 (ja) * | 2002-01-30 | 2006-12-20 | 三菱重工業株式会社 | 多気筒内燃機関 |
| US7283140B2 (en) * | 2005-06-21 | 2007-10-16 | Microsoft Corporation | Texture montage |
-
2006
- 2006-01-16 JP JP2006008057A patent/JP4599301B2/ja not_active Expired - Fee Related
- 2006-12-28 EP EP06027001.4A patent/EP1808775A3/en not_active Withdrawn
- 2006-12-28 US US11/646,285 patent/US20070179762A1/en not_active Abandoned
Patent Citations (5)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| JPH01145723A (ja) * | 1987-08-28 | 1989-06-07 | Hitachi Ltd | プログラム生成方法 |
| JPH06148038A (ja) * | 1992-11-09 | 1994-05-27 | Hitachi Ltd | 振動解析装置 |
| JP2005009323A (ja) * | 2003-06-16 | 2005-01-13 | Mitsubishi Heavy Ind Ltd | 多気筒内燃機関、そのクランク軸のクランク配置決定方法 |
| JP2005202634A (ja) * | 2004-01-15 | 2005-07-28 | Keio Gijuku | 六面体有限要素数値解析方法 |
| JP2005267214A (ja) * | 2004-03-18 | 2005-09-29 | Canon Inc | 有限要素法による制御プログラムおよび記憶媒体 |
Also Published As
| Publication number | Publication date |
|---|---|
| US20070179762A1 (en) | 2007-08-02 |
| EP1808775A3 (en) | 2013-05-01 |
| JP2007188454A (ja) | 2007-07-26 |
| EP1808775A2 (en) | 2007-07-18 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| JP4599301B2 (ja) | 設計支援装置及びコンピュータプログラム | |
| Boukharfane et al. | A combined ghost-point-forcing/direct-forcing immersed boundary method (IBM) for compressible flow simulations | |
| US7647211B2 (en) | Constraint-based method of designing a route for a transport element | |
| Toal | Some considerations regarding the use of multi-fidelity Kriging in the construction of surrogate models | |
| Le Beau et al. | SUPG finite element computation of compressible flows with the entropy and conservation variables formulations | |
| Masters et al. | Review of aerofoil parameterisation methods for aerodynamic shape optimisation | |
| Colonna et al. | Numerical simulation of dense gas flows on unstructured grids with an implicit high resolution upwind Euler solver | |
| Firrone et al. | Passive control of vibration of thin-walled gears: advanced modelling of ring dampers | |
| Brahimi et al. | Blade dynamical response based on aeroelastic analysis of fluid structure interaction in turbomachinery | |
| Shershnev et al. | HyCFS, a high-resolution shock capturing code for numerical simulation on hybrid computational clusters | |
| US9361413B1 (en) | Systems and methods for simulating contact between physical objects | |
| Xia et al. | Particle swarm optimization of aerodynamic shapes with nonuniform shape parameter–based radial basis function | |
| Barrett et al. | Integrated free‐form method for aerostructural optimization of wind turbine blades | |
| Madadi et al. | Aerodynamic design of S-shaped diffusers using ball–spine inverse design method | |
| Scholten et al. | Uncoupled method for static aeroelastic analysis | |
| Ghoreyshi et al. | Numerical simulation and reduced-order aerodynamic modeling of a lambda wing configuration | |
| Deng et al. | Reduced-order multiscale modeling of plastic deformations in 3D alloys with spatially varying porosity by deflated clustering analysis | |
| Chen et al. | Towards high-accuracy deep learning inference of compressible turbulent flows over aerofoils | |
| Tissera et al. | Computational fluid dynamics methods for hypersonic flow around blunted-cone-cylinder-flare | |
| Vandevoorde et al. | Comparison of algorithms for unsteady flow calculations in inlet and exhaust systems of IC engines | |
| Shan et al. | Turbulence Modeling via Data Assimilation and Machine Learning for Separated Flows over Airfoils | |
| Pereira | Generalized finite element methods for three-dimensional crack growth simulations | |
| Housman et al. | Time-derivative preconditioning methods for multicomponent flows—part ii: Two-dimensional applications | |
| Cao et al. | Nacelle Inlet Optimization at High Angles of Attack Based on the Ensemble Indicator Method | |
| Watanabe et al. | Application of High Fidelity Simulation to LE-X Engine Development |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| A621 | Written request for application examination |
Free format text: JAPANESE INTERMEDIATE CODE: A621 Effective date: 20080819 |
|
| A977 | Report on retrieval |
Free format text: JAPANESE INTERMEDIATE CODE: A971007 Effective date: 20090417 |
|
| A131 | Notification of reasons for refusal |
Free format text: JAPANESE INTERMEDIATE CODE: A131 Effective date: 20090626 |
|
| A521 | Request for written amendment filed |
Free format text: JAPANESE INTERMEDIATE CODE: A523 Effective date: 20090824 |
|
| A131 | Notification of reasons for refusal |
Free format text: JAPANESE INTERMEDIATE CODE: A131 Effective date: 20100115 |
|
| RD02 | Notification of acceptance of power of attorney |
Free format text: JAPANESE INTERMEDIATE CODE: A7422 Effective date: 20100201 |
|
| A521 | Request for written amendment filed |
Free format text: JAPANESE INTERMEDIATE CODE: A523 Effective date: 20100316 |
|
| 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: 20100910 |
|
| A01 | Written decision to grant a patent or to grant a registration (utility model) |
Free format text: JAPANESE INTERMEDIATE CODE: A01 |
|
| A61 | First payment of annual fees (during grant procedure) |
Free format text: JAPANESE INTERMEDIATE CODE: A61 Effective date: 20100927 |
|
| FPAY | Renewal fee payment (event date is renewal date of database) |
Free format text: PAYMENT UNTIL: 20131001 Year of fee payment: 3 |
|
| LAPS | Cancellation because of no payment of annual fees |