正交基底气泡函数要素数值分析方法、正交基底气泡函数要素数值分析程序以及正交基...的制作方法

文档序号:6553115阅读:189来源:国知局

专利名称::正交基底气泡函数要素数值分析方法、正交基底气泡函数要素数值分析程序以及正交基...的制作方法
技术领域
:本发明涉及一种正交基底气泡函数要素数值分析方法、正交基底气泡函数要素数值分析程序、以及正交基底气泡函数要素数值分析装置,其用来对基于使用气泡函数要素的有限元法的分析(有限元分析),使用只成为计算效率优秀的对角项的质量行列,进行可靠性高的数值仿真。
背景技术
:对以往的气泡函数要素进行说明。图47为表示以往的二维气泡函数要素的说明图,图48为表示以往的三维气泡函数要素的说明图。如图47以及图48所示,使用三角形(四面体)要素的气泡函数要素,使用在各个要素中形成三角形(四面体)的3(4)个点以及重心点这4(5)个节点,在等参数座标系[r,s]({r,s,t})中通过下述公式(1)来表示(例如参照下述非专利文献1、非专利文献2、非专利文献3)。uh|Ωe=Σα=1N+1Φαuα+φBuB=ΦTu=uTΦ]]>…式(1)Φα=Ψα-1N+1φB,α=1···N+1]]>…式(2)式(1)中的Φα、φB是表示气泡函数要素的形状函数,uα、uB表示三角形(四面体)的各个顶点的值(分析物理量)、重心点的值(分析物理量)、N表示空间维数。如果通过矢量形式来描述,则形状函数变为下面的式(3)~式(6)。二维ФT=[Ф1Ф2Ф3φB]…式(3)uT=[u1u2u3uB]…式(4)三维ФT=[Ф1Ф2Ф3Ф4φB]…式(5)uT=[u1u2u3u4uB]…式(6)式(2)中的ψa是使用二维或三维一次要素的形状函数,通过下述公式(7)、(8)来表示。二维ψ1=1-r-s,ψ2=r,ψ3=s…式(7)三维ψ1=1-r-s-t,ψ2=r,ψ3=s,ψ4=t…式(8)形状函数φB称作气泡函数。气泡函数位于要素边界上,将每一个要素定义为,其值为0,重心点的值为1。在不稳定(非恒定)问题中,空间方向的离散化中使用气泡函数要素的有限元方程式,能够表示为下述式(9)。Mu·+F(u)=0]]>…式(9)式(9)的u是应当求出的未知分析物理量(污染物质浓度、温度、流量、水深、流速、压力、位移等),M是质量行列式,F(u)是时间微分项以外归纳而成的项。作为式(9)的时间方向的离散化,基于泰勒(テイラ一)扩展的4段解法,表示为下述式(10)~式(13)(例如参照下述非专利文献4)。4段解法<第1步>un+1/4=un-M-1Δt4F(un)]]>…式(10)<第2步>un+2/4=un-M-1Δt3F(un+1/4)]]>…式(11)<第3步>un+3/4=un-M-1Δt2F(un+2/4)]]>…式(12)<第4步>un+1=un-M-1ΔtF(un+3/4)…式(13)式(10)~式(13)的上述标n表示当前时刻n时已知的分析物理量,n+1表示从时刻n经过了微小时间Δt之后的未知分析物理量。非专利文献1D.N.Arnold,F.BrezziandM.Fortin,“AStableFiniteElementfortheStokesEquations”,Calcolo,Vol.23,1984,pp.337-pp.344非专利文献2J.C.Simo,F.ArmeroandC.A.Taylor,“StableandTime-DessipativeFiniteElementMethodsfortheIncompressibleNavier-StrokesEquationsinAdvectionDominatedFlows”,InternationalJournalforNumericalMethodsinEngineering,Vol.38,1995,pp.1475-pp.1506非专利文献3松本純一,“気泡関数を用いた非圧縮性粘性流れ分析のための2レベル-3レベル有限元法”,応用力学論文集(土木学会),第7卷,2004年8月,339页-346页非专利文献4畑中勝守,“多段階有限元法による非圧縮粘性流体の順·逆分析に関する計算力学の研究”,中央大学博士论文,1993年3月这里如上述式(10)~式(13)所示,为了使用4段解法,根据该已知的分析物理量un求出未知分析物理量un+1,需要质量行列式的逆行列式。图49为表示以往的RotatingCone问题的分析中所使用的分析模型的说明图,图50为表示以前的RotatingCone问题的分析中所使用的初始条件的等高线的说明图,图51为表示通过以前的匹配质量行列式所计算出的RotatingCone问题的分析结果的等高线的说明图(参照非专利文献4)。设图49的分析模型4900为初始状态(0周时)。之后,如果使用质量行列式的逆行列式,从初始状态(0周时)开始旋转给定周,例如5周,图50中所示的分析模型4900的初始状态中的等高线模型5000,就变为图51所示的分析结果(等高线模型5100)。上述质量行列式,在时间方向的离散化中使用气泡函数要素,因此一般成为较疏的分布行列(匹配质量行列)。因此,为了求出该分布行列式(匹配质量行列式)的逆行列式,在数值分析上需要较多的存储容量与计算时间,存在装置本体的成本增加,以及分析处理延迟这一问题。为了解决该问题,通常使用将质量行列式的各行的成分相加(集中化),而只在对角项中保留成分的近似行列式(集中质量行列式)。在使用集中质量行列式的情况下,行列式的成分仅仅是对角项,因此逆行列式变为取各个对角成分的倒数的行列式,与在数值分析上使用没有进行近似的匹配质量行列式及其逆行列式的情况相比,仅仅需要很少的存储容量与计算时间就能够执行分析。但是,在使用上述集中质量行列式的情况下,集中质量行列式不是与原来的质量行列式相同的行列式,而是近似行列式。因此,如果从图49中所示的分析模型4900的初始状态(0周时)开始旋转给定周,例如5周,图50中所示的分析模型4900的初始状态中的等高线模型5000,就变为图52所示的分析结果(等高线模型5200)。这样,由于分析模型4900旋转5周后的等高线模型5200,相对初始状态中的等高线模型5000或图51中所示的等高线模型5100大幅变形,因此存在计算精度恶化,分析结果的可靠性较低这一问题。
发明内容本发明为了解决上述以往技术的问题点,目的在于提供一种能够实现简单且可靠性高的有限元分析的正交基底气泡函数要素数值分析方法、正交基底气泡函数要素数值分析程序、以及正交基底气泡函数要素数值分析装置。为解决上述问题,实现发明目的,本发明的正交基底气泡函数要素数值分析方法、正交基底气泡函数要素数值分析程序、以及正交基底气泡函数要素数值分析装置,取得分析对象的各个要素的匹配质量行列式,根据上述分析对象的各个要素的气泡函数,生成将通过上述第2取得工序所取得的各个要素的匹配质量行列式对角化了的各个要素的对角质量行列式,并根据分析对象的已知分析物理量,与所生成的各个要素的对角质量行列式,对上述分析对象的动作进行分析。另外,还可以通过将各个要素中的气泡函数的积分值,代入到各个要素的匹配质量行列式中,来生成各个要素的对角质量行列式。进而,还可以根据所生成的要素的对角质量行列式,计算出分析对象全区域的对角质量行列式,计算出上述分析对象全区域的对角质量行列式的逆行列式,并根据分析对象的已知分析物理量,与上述分析对象全区域的对角质量行列式,以及上述逆行列式,对上述分析对象的动作进行分析。通过本发明,在使用气泡函数要素的有限元分析中,不进行质量行列式的近似,使用质量行列式满足成为对角行列式的下述式(14)的气泡函数。⟨φB,1⟩Ωe=||φB||Ωe2=CAe]]>…式(14)二维C=34]]>三维C=45]]>⟨φB,1⟩Ωe=∫ΩeφBdΩ]]>…式(15-1)||φB||Ωe2=⟨φB,φB⟩Ωe=⟨φB2,1⟩Ωe=∫ΩeφB2dΩ]]>…式(15-2)Ae=∫ΩedΩ]]>…式(15-3)上述式(15-1)~式(15-3),是进行气泡函数要素的公式化所必需的积分,<·,·>Ωe表示要素区域Ωe的积分,Ae表示要素区域Ωe的面积(体积)。通过这样,不需要使用通过将各个要素的匹配质量行列式内的成分相加而近似得到的集中质量行列式,而能够导入气泡函数要素的基底(形状函数)正交的条件,使用各个要素的高精度的对角质量行列式,进行分析。另外,还能够简单地计算出分析对象全区域的对角质量行列式及其逆行列式,能够使用未知的分析对象物理量对分析对象的行为进行高精度的分析。本发明的相关正交基底函数要素数值分析方法、正交基底函数要素数值分析程序、以及正交基底函数要素数值分析装置,通过简单且高精度地进行分析处理,起到了能够实现分析对象的行为分析的可靠性的这一效果。另外,还起到了能够实现存储容量的降低以及分析时间的缩短等有效的分析方法(正交基底函数要素数值分析方法)这一效果。图1为表示本发明的实施方式的相关数值分析装置的硬件构成的方框图。图2为表示本发明的实施方式的相关数值分析装置的功能构成的方框图。图3为表示本发明的实施方式的相关数值分析装置的数据存储部所存储的分析物理量信息表的说明图。图4为表示本发明的实施方式的相关数值分析装置中的数值分析处理顺序的流程图。图5为表示二维气泡函数要素的说明图。图6为表示三维气泡函数要素的说明图。图7为表示形成正交基底的三角形气泡函数的形状φB的说明图。图8为表示形成正交基底的三角形气泡函数的形状φB2的说明图。图9为表示形成正交基底的三角形气泡函数的形状φB的说明图。图10为表示形成正交基底的三角形气泡函数的形状φB2的说明图。图11为表示二维气泡函数的说明图。图12为表示三维气泡函数的说明图。图13为表示RotatingCone问题的分析中的计算区域的说明图。图14为表示RotatingCone问题的分析中的分析模型之网格的说明图。图15为表示RotatingCone问题的分析中的分析模型的流速(流况)的说明图。图16为表示RotatingCone问题的分析中的初始条件的概观图。图17为表示RotatingCone问题的分析中的初始条件的等高线的说明图。图18为表示使用气泡函数(式(47))中的匹配质量行列式的计算结果的概观图。图19为表示使用气泡函数(式(47))中的匹配质量行列式的计算结果的等高线的说明图。图20为表示使用气泡函数(式(47))中的集中质量行列式的计算结果的概观图。图21为表示使用气泡函数(式(47))中的集中质量行列式的计算结果的等高线的说明图。图22为表示使用气泡函数(式(48))中的匹配质量行列式的计算结果的概观图。图23为表示使用气泡函数(式(48))中的匹配质量行列式的计算结果的等高线的说明图。图24为表示使用气泡函数(式(48))中的集中质量行列式的计算结果的概观图。图25为表示使用气泡函数(式(48))中的集中质量行列式的计算结果的等高线的说明图。图26为表示使用气泡函数(式(49))中的对角质量行列式的计算结果的概观图。图27为表示使用气泡函数(式(49))中的对角质量行列式的计算结果的等高线的说明图。图28为表示RotatingCone问题的分析中的分析模型的网格(节点数438413,要素数2539200)的说明图。图29为表示RotatingCone问题的分析中的垂直位置=0的剖面的初始条件的概观图。图30为表示RotatingCone问题的分析中的垂直位置=0部分的初始条件的等高线的说明图。图31为表示气使用泡函数(式(50))中的匹配质量行列式的垂直位置=0部分的计算结果的概观图。图32为表示使用气泡函数(式(50))中的匹配质量行列式的垂直位置=0部分的计算结果的等高线的说明图。图33为表示使用气泡函数(式(50))中的集中质量行列式的垂直位置=0部分的计算结果的概观图。图34为表示使用气泡函数(式(50))中的集中质量行列式的垂直位置=0部分的计算结果的等高线的说明图。图35为表示使用气泡函数(式(51))中的匹配质量行列式的垂直位置=0部分的计算结果的概观图。图36为表示使用气泡函数(式(51))中的匹配质量行列式的垂直位置=0部分的计算结果的等高线的说明图。图37为表示使用气泡函数(式(51))中的集中质量行列式的垂直位置=0部分的计算结果的概观图。图38为表示使用气泡函数(式(51))中的集中质量行列式的垂直位置=0部分的计算结果的等高线的说明图。图39为表示气泡函数(式(52))中的使用对角质量行列式的垂直位置=0部分的计算结果的概观图。图40为表示气泡函数(式(52))中的使用对角质量行列式的垂直位置=0部分的计算结果的等高线的说明图。图41为表示等参数座标系r的说明图(之1)。图42为表示等参数座标系r的说明图(之2)。图43为表示使用式(81)、(93)的三角形气泡函数的形状φB的说明图。图44为表示使用式(81)、(93)的三角形气泡函数的形状φB2的说明图。图45为表示使用式(132)的三角形3级气泡函数的形状的说明图。图46为表示使用式(81)、(93)的三角形气泡函数与使用式(132)的三角形3级气泡函数之积的形状的说明图。图47为表示以往的二维气泡函数要素的说明图。图48为表示以往的三维气泡函数要素的说明图。图49为表示以往的RotatingCone问题的分析中所使用的分析模型的说明图。图50为表示以往的RotatingCone问题的分析中所使用的初始条件的等高线的说明图。图51为表示通过以往的匹配质量行列式所计算出来的RotatingCone问题的分析结果的等高线的说明图。图52为表示通过以前的集中质量行列式所计算出来的RotatingCone问题的分析结果的等高线的说明图。图中200-数值分析装置,201-数据存储部,202-第1取得部,203-第2取得部,204-生成部,205-第1计算部,206-第2计算部,207-分析部。具体实施例方式下面对照附图,对本发明的相关正交基底气泡函数要素数值分析方法、正交基底气泡函数要素数值分析程序、以及正交基底气泡函数要素数值分析装置(以下简称作“数值分析装置”)的理想实施方式进行详细说明。(数值分析装置的硬件结构)下面对本实施方式的相关数值分析装置的硬件结构进行说明。图1为表示本发明的实施方式的相关数值分析装置的硬件结构的方框图。数值分析装置具有CPU101、ROM102、RAM103、HDD(硬盘驱动器)104、HD(硬盘)105、FDD(软盘驱动器)106、作为可插拔的存储介质之一例的FD(软盘)107、显示器108、I/F(接口)109、键盘110、鼠标111、扫描仪112、以及打印机113。另外,各个构成部通过总线100分别连接起来。这里,CPU101负责数值分析装置全体的控制。ROM102存储引导程序等程序。RAM103用作CPU101的工作区域。HDD104在CPU101的控制下,控制对HD105的数据读/写。HD105在HDD104的控制下存储所写入的数据。FDD106在CPU101的控制下,控制对FD107的数据读/写。FD107在FDD106的控制下存储所写入的数据,或将FD107中所存储的数据读取到数值分析装置中。可拔插的存储介质,除了FD107之外,还可以是CD-ROM(CD-R、CD-RW)、MO、DVD(DigitalVersatileDisk)、存储卡等。显示器108显示出以光标、图标或工具箱等为代表的文档、图像、功能信息等数据。该显示器108例如可以采用CRT、TFT液晶显示器、等离子显示器等。I/F109通过通信线路与互联网等网络相连接,并经由该互联网与其他装置相连接。另外,I/F109负责网络与内部的接口,控制来自外部装置的数据输入输出。I/F109例如可以采用调制解调器或LAN适配器等。键盘110中具有用来输入文字、数字、各种指示等的按键,进行数据的输入。另外也可以是触摸屏式的输入板或十字键等。鼠标111进行光标的移动、选择范围或窗口的移动、以及大小的变更等。另外,还可以是轨迹球或操纵杆等,只要作为点击设备具有同样的功能就可以。扫描仪112光学读取图像,并将图像数据获取到数值分析装置内。另外,打印机113印刷图像数据或文件数据。打印机113例如可以采用激光打印机或喷墨打印机。(数值分析装置的功能构成)对本发明的实施方式的相关数值分析装置的功能构成进行说明。图2为表示本发明的实施方式的相关数值分析装置的功能构成的方框图。如图2所示,数值分析装置200由数据存储部201、第1取得部202、第2取得部203、生成部204、第1计算部205、第2计算部206、以及分析部207构成。数据存储部201具有分析对象的数值分析中所使用的分析物理量信息表。这里对分析物理量信息表进行说明。图3为表示本发明的实施方式的相关分析物理量信息表的方框图。图3中,分析物理量信息表中,存储有关于分析物理量ID、分析物理量名(污染物质、热、流体、构造等)、以及已知分析物理量。另外,作为已知分析物理量存储有已知物性值、边界分析物理量、以及初始分析物理量等。分析物理量ID与分析物理量名,是通过用户的输入操作而指定的信息,一旦指定了分析物理量ID或分析物理量名,就能够提取相关联的已知物性值、边界分析物理量以及初始分析物理量。已知物性值,是决定由用户的输入操作指定了分析物理量ID或分析物理量名的分析对象,是什么样的物质的信息。作为已知物性值可以列举出例如扩散系数、密度、比热、热传导系数、粘性系数、水底面摩擦系数、杨氏模量、泊松比等。边界分析物理量,是指由分析对象的各个要素所构成的分析区域(分割成三角形或四面体的网格)在什么样的条件下被进行分析的分析条件的信息。作为边界分析物理量,可以列举出例如物质浓度、温度、流速、压力、流量、水深、位移等。关于边界分析物理量,在分析对象例如是污染物质或热的情况下,在污染物质或热的发生(提供)源的部分中,物质浓度=发生(提供)源值,温度=发生(提供)源值。另外,如果分析对象是流体且边界是壁,则流速=0,流量=0,如果分析对象是构造且边界能够固定物体,则能够设为位移=0。如果是恒定问题(分析物理量不依赖于时间,时时刻刻都不变的问题),通过上述已知物性值与边界分析物理量,就能够进行分析。在非恒定问题(分析物理量依赖于时间,时刻分析变化的问题)的情况下,还需要作为初始条件的初始分析物理量。初始分析物理量是分析开始时的分析物理量的值。另外,数据存储部201中,虽然未图示,但作为初始设定信息,为了对时间方向分析逐次分析计算,而将时间增加量、计算步骤数等信息,也存储在分析物理量信息表中。进而,还存储有分割成了分析对象范围的网格的分割宽度、节点数、要素数、对应于节点编号的节点坐标值、以及对应于要素编号的要素结合信息等网格信息。该网格在分析对象范围是二维空间的情况下,为三角形的要素,在分析对象范围是三维空间的情况下,是四面体的有限形状的要素。网格可以由图1中所示的扫描仪112读取并取得。该数据存储部201,具体的说,例如能够通过图1中所示的RAM103、HD105或FD107来实现其功能。第1取得部202取得分析对象的已知分析物理量。具体的说,通过操作输入装置(例如图1中所示的键盘110或鼠标111),选择成为分析对象的物理量的分析物理量ID,从数据存储部201提取已知分析物理量。该第1取得部202,具体的说,例如通过由CPU101执行存储在图1中所示的ROM102、RAM103、HD105、FD107等中的程序,或通过图1中所示的I/F109来实现其功能。第2取得部203,取得分析对象的各个要素(网格)的匹配质量行列式。第2取得部203,具体的说,例如通过由CPU101执行存储在图1中所示的ROM102、RAM103、HD105、FD107等中的程序,或通过图1中所示的I/F109来实现其功能。生成部204根据分析对象的各个要素的气泡函数,生成将通过第2取得部203所取得的各个要素的匹配质量行列式对角化而成的各个要素的对角质量行列式。生成部204具体的说,例如通过由CPU101执行存储在图1中所示的ROM102、RAM103、HD105、FD107等中的程序,或通过图1中所示的I/F109来实现其功能。第1计算部205,根据由生成部204所生成的各个要素的对角质量行列式,计算出分析对象全区域的对角质量行列式。第1计算部205,具体的说,例如通过由CPU101执行存储在图1中所示的ROM102、RAM103、HD105、FD107等中的程序,或通过图1中所示的I/F109来实现其功能。第2计算部206,计算出由第1计算部205所计算出的分析对象全区域的对角质量行列式的逆行列式。第2计算部206,具体的说,例如通过由CPU101执行存储在图1中所示的ROM102、RAM103、HD105、FD107等中的程序,或通过图1中所示的I/F109来实现其功能。分析部207,根据第1取得部202所取得的已知分析物理量,与生成部204所生成的各个要素的对角质量行列式,对分析对象的行为进行分析。具体的说,根据第1取得部202所取得的已知分析物理量、第1计算部205所计算出的分析对象全区域的对角质量行列式、以及第2计算部206所计算出的逆行列式,对分析对象的行为进行分析。另外,还有不需要分析对象全区域的对角质量行列式,而只使用各个要素的对角质量行列式,对分析对象的动作进行分析的解法(例如参照O.C.Zienkiewicz,R.L.Taylor,S.J.SherwinandJ.Perio,“OnDiscontinuousGalerkinMethods”,InternationalJournalforNumericalMethodsinEngineering,Vol.58,2003,pp.1119-1148),这种情况下,将生成部204与第1计算部205看作同一个,能够进入第2计算部206。分析部207具体的说,例如通过由CPU101执行存储在图1中所示的ROM102、RAM103、HD105、FD107等中的程序,或通过图1中所示的I/F109来实现其功能。对本实施方式的相关数值分析装置200的数值分析处理顺序进行说明。图4为表示本发明的实施方式的相关数值分析装置200的数值分析处理顺序的流程图。首先,通过第1取得部202取得分析对象的已知分析物理量(步骤S401)。接下来,通过第2取得部203,取得各个要素(要素等级)的匹配质量行列式(步骤S402)。接下来,对每一个要素积分气泡函数(步骤S403),并将步骤S403中积分得到的值代入到各个要素(要素等级)的匹配质量行列式中,通过这样计算出各个要素(要素等级)的对角质量行列式(步骤S404)。接下来,通过各个要素(要素等级)的对角质量行列式的总和(叠加重ね合せ),计算出分析对象全区域的对角质量行列式(步骤S405)。之后,计算出该对角质量行列式的逆行列式(步骤S406),根据分析对象的已知分析物理量、分析对象全区域的对角质量行列式、及其逆行列式,对分析对象的行为进行分析。具体的说,例如能够通过代入到上述式(10)~式(13)中,进行分析。对数值分析装置200的具体分析处理内容进行说明。<1.正交基底气泡函数要素>质量行列式M,在对基于有限元法的分析(有限元分析),给对某个进行物理量进行过插补的函数乘以权重函数,并在对象区域中进行了积分的情况下出现。在使用气泡函数要素的情况下,通过在全体系重合要素群的作业而形成,所述要素群是对成为分析对象的任意区域有区分地用三角形(四面体)形状以要素数Ne分割而成的。并使用要素区域中的积分,通过下式(16)来表示。M=Σe=1Ne⟨Φ,ΦT⟩Ωe=Σe=1NeMij(e),i=1···N+2,j=1···N+2]]>…式(16)为了让气泡函数要素的基底(形状函数)垂直于公式(16),必需满足下述式(18)、式(20)。另外,在式(18)的导入时,使用式(17)。另外,在式(20)的导入时,使用式(19)。⟨Φα,φB⟩Ωe=1N+1⟨1,φB⟩Ωe-1N+1||φB||Ωe2=0,α=1···N+1]]>…式(17)[公式9]根据式(17)⟨1,φB⟩Ωe=||φB||Ωe2]]>…式(18)[公式10]⟨Φα,Φβ⟩Ωe=⟨Ψα,Ψβ⟩Ωe-1(N+1)2⟨1,φB⟩Ωe=0,α≠β,β=1···N+1]]>…式(19)[公式11]根据式(19)⟨1,φB⟩Ωe=(N+1)2⟨Ψα,Ψβ⟩Ωe]]>…式(20)通过式(18)、(20)得到以下的关系式(21)。⟨φB,1⟩Ωe=||φB||Ωe2=CAe]]>…式(21)C=1Ae(N+1)2⟨Ψα,Ψβ⟩Ωe,α≠β]]>二维C=34]]>三维C=45]]><2.满足正交条件的气泡函数>提出通过下式(22)所表示的气泡函数,作为能够满足式(18)、(20)的气泡函数。φB=α1φ^B+α2φ~B+φ‾Bα1+α2+1]]>…式(22)α1,α2是未知量,φB分别是不同的气泡函数。通过将式(22)代入到式(18)中,得到下述式(23)。β1α12+β2α22+β3+β4α1α2+β5α1+β6α2=0]]>…式(23)β1=⟨φ^B,1⟩Ωe-||φ^B||Ωe2,β2=⟨φ~B,1⟩Ω2-||φ~B||Ωe2,β3=⟨φ‾B,1⟩Ωe-||φ‾B||Ωe2,]]>β4=⟨φ~B,1⟩Ωe+⟨φ^B,1⟩Ωe-2⟨φ^B,φ~B⟩Ωe,]]>β5=⟨φ‾B,1⟩Ωe+⟨φ^B,1⟩Ωe-2⟨φ^B,φ‾B⟩Ωe,]]>β6=⟨φ‾B,1⟩Ωe+⟨φ~B,1⟩Ωe-2⟨φ~B,φ‾B⟩Ωe]]>或者β1=(⟨φ^B,1⟩Ωe-||φ^B||Ωe2)Ae-1,]]>β2=(⟨φ~B,1⟩Ωe-||φ~B||Ωe2)Ae-1,]]>β3=(⟨φ‾B,1⟩Ωe-||φ‾B||Ωe2)Ae-1,]]>β4=(⟨φ~B,1⟩Ωe+⟨φ^B,1⟩Ωe-2⟨φ^B,φ~B⟩Ωe)Ae-1,]]>β5=(⟨φ‾B,1⟩Ωe+⟨φ^B,1⟩Ωe-2⟨φ^B,φ‾B⟩Ωe)Ae-1,]]>β6=(⟨φ‾B,1⟩Ωe+⟨φ~B,1⟩Ωe-2⟨φ~B,φ‾B⟩Ωe)Ae-1]]>通过将式(22)代入到式(20)中,得到下式(24)。α2=γ1α1+γ2…式(24)γ1=-⟨φ^B,1⟩Ωe-CAe⟨φ~B,1⟩Ωe-CAe,γ2=-⟨φ‾B,1⟩Ωe-CAe⟨φ~B,1⟩Ωe-CAe]]>或者γ1=-⟨φ^B,1⟩ΩeAe-1-C⟨φ~B,1⟩ΩeAe-1-C,γ2=-⟨φ‾B,1⟩ΩeAe-1-C⟨φ~B,1⟩ΩeAe-1-C]]>通过将式(24)代入到式(23)中,最终得到以α1为未知量的下式(25)。αα12+bα1+c=0]]>…式(25)α=β1+β2γ12+β4γ1,]]>b=2β2γ1γ2+β4γ2+β5+β6γ1,c=β2γ22+β3+β6γ2]]>式(25)是关于α1的二次方程式,因此通过解的公式,能够由下式(26)求出α1。α1=-b±b2-4ac2a]]>…式(26)<3.形成正交基底的具体气泡函数>对形成正交基底的具体气泡函数进行说明。图5为表示二维气泡函数要素的说明图,图6为表示三维气泡函数要素的说明图。为了导出形成正交基底的具体气泡函数,使用将图5与图6中所示的三角形(四面体)的要素区域,用重心点分割成3(4)个小三角形(小四面体)wi,i=1…N+1而得到的ξ次方气泡函数(0<ξ<∞)。ξ次方气泡函数对每一个小三角形(小四面体)使用等参数座标系{r,s}({r,s,t}),定义为如下述式(27)与式(28)那样。二维φBξ=3ξ(1-r-s)ξinw13ξrξinw23ξsξinw3]]>…式(27)三维φBξ=4ξ(1-r-s-t)ξinw14ξrξinw24ξsξinw34ξtξinw4]]>…式(28)通过下式(29),能够求出该气泡函数的要素区域中的积分值。⟨φBξ,1⟩Ωe=N!Πγ=1N(ξ+γ)Ae]]>…式(29)将式(22)使用ξ次方气泡函数定义为下式(30)。φ^B=φBξ^,]]>φ~B=φBξ~,]]>φ‾B=φBξ‾]]>…式(30)ξ为各个幂的值如果将式(30)的各个幂的值如下式(31)进行设定,则对于二维(三角形)、三维(四面体)的气泡函数,α1、α2为就变为下述式(32)~式(35)。ξ^=110,]]>ξ~=15,]]>ξ‾=34]]>…式(31)[公式22]二维(三角形)α1=-b+b2-4ac2a=3.1871016···,α2=-4.5742685···]]>…式(32)α1=-b-b2-4ac2a=3.1457828···,α2=-3.9426812···]]>…式(33)三维(四面体)α1=-b+b2-4ac2a=3.1277974···,α2=-3.8884541···]]>…式(34)α1=-b-b2-4ac2a=2.6785040···,α2=-4.0779944···]]>…式(35)如果将式(30)的各个幂的值如下式(36)进行设定(与式(31)不同的值),则对于二维(三角形)、三维(四面体)的气泡函数,α1、α2就变为下述式(37)~式(40)。ξ^=3,]]>ξ~=2,]]>ξ‾=65]]>…式(36)[公式24]二维(三角形)α1=-b+b2-4ac2a=-0.1393869···,α2=-0.6433844···]]>…式(37)α1=-b-b2-4ac2a=1.2696544···,α2=-2.2134591···]]>…式(38)三维(四面体)α1=-b+b2-4ac2a=0.7333511···,α2=-1.6387018···]]>…式(39)α1=-b-b2-4ac2a=1.2578218···,α2=-2.2006346···]]>…式(40)图7为表示形成使用了式(32)的α1、α2的正交基底的三角形气泡函数的形状φB的说明图。图8为表示形成使用了式(32)的α1、α2的正交基底的三角形气泡函数的形状φB2的说明图。另外,图9为表示形成使用了式(37)的α1、α2的正交基底的三角形气泡函数的形状φB的说明图。图10为表示形成使用了式(37)的α1、α2的正交基底的三角形气泡函数的形状φB2的说明图。如式(32)、式(37)以及图7~图10所示,满足式(21)的气泡函数(形状)存在有无数个,因此形成正交基底的具体气泡函数及其形状几乎没有意义,正交基底气泡函数要素的本质在于发明了与气泡函数要素的基底相垂直的条件式(21)的部分。<4.正交基底气泡函数要素的要素等级的质量行列式>下述式(41)~式(46)中,示出了二维(三角形)、三维(四面体)中的通常的气泡函数要素的要素等级的匹配质量行列式、集中质量行列式、以及正交基底气泡函数要素的要素等级的质量行列式(对角质量行列式)。三角形气泡函数要素要素等级(level)的匹配质量行列式⟨Φ,ΦT⟩Ωe=Mij(e)=Ae6Ae12Ae120Ae12Ae6Ae120Ae12Ae12Ae600000]]>+19⟨1,φB⟩Ωe-2-2-23-2-2-23-2-2-233330+19||φB||Ωe2111-3111-3111-3-3-3-39]]>…式(41)要素等级的集中质量行列式(集中了的情况)⟨Φ,ΦT⟩Ωe=Mij(e)≈diag(Σj=1N+2Mij(e))]]>=Ae30000Ae30000Ae300000+13⟨1,φB⟩Ωe-10000-10000-100003]]>…式(42)要素等级的对角质量行列式(成为正交基底的气泡函数)⟨Φ,ΦT⟩Ωe=Mij(e)=Ae120000Ae120000Ae13000034Ae]]>…式(43)式(43)的要素等级的对角质量行列式M(e)ij,由生成部204将气泡函数φB的积分值CAe代入到式(41)的要素等级的匹配质量行列式M(e)ij的<1,φB>Ωe,以及‖φB‖2Ωe中而生成。之后,通过由第1计算部205取各个要素的要素等级的对角质量行列式M(e)ij的总和(合并),能够计算出分析对象全区域的对角质量行列式M。另外,通过第2计算部206能够计算出对角质量行列式M的逆行列式M-1。这样,由分析部207将已知的分析物理量u、对角质量行列式M以及逆行列式M-1,代入到上述式(10)~式(13)中,就能够对分析对象的行为进行分析。四面体气泡函数要素要素等级的匹配质量行列式⟨Φ,ΦT⟩Ωe=Mij(e)=Ae10Ae20Ae20Ae200Ae20Ae10Ae20Ae200Ae20Ae20Ae10Ae200Ae20Ae20Ae20Ae10000000]]>+1416⟨1,φB⟩Ωe-2-2-2-24-2-2-2-24-2-2-2-24-2-2-2-2444440+116||φB||Ωe21111-41111-41111-41111-4-4-4-4-416]]>…式(44)要素等级的集中质量行列式(集中化后的情况)⟨Φ,ΦT⟩Ωe=Mij(e)≈diag(Σj=1N+2Mij(e))]]>=Ae400000Ae400000Ae400000Ae4000000+14⟨1,φB⟩Ωe-100000-100000-100000-1000004]]>…式(45)要素等级的对角质量行列式(成为正交基底的气泡函数)⟨Φ,ΦT⟩Ωe=Mij(e)=Ae2000000Ae2000000Ae2000000Ae200000045Ae]]>…式(46)式(46)的要素等级的对角质量行列式,由生成部204将气泡函数φB的积分值CAe代入到式(44)的要素等级的匹配质量行列式的<1,φB>Ωe,以及‖φB‖2Ωe中而生成。之后,通过由第1计算部205取各个要素的要素等级的对角质量行列式M(e)ij的总和(合并),能够计算出分析对象全区域的对角质量行列式M。另外,通过第2计算部206能够计算出对角质量行列式M的逆行列式M-1。这样,由分析部207将已知的分析物理量u、对角质量行列式M以及逆行列式M-1,代入到上述式(10)~式(13)中,就能够对分析对象的行为进行分析。<5.关于气泡函数>对气泡函数进行说明。图11为表示二维气泡函数的说明图,图13为表示三维气泡函数的说明图。气泡函数能够选择任意的函数,只要满足形状函数的定义(重心点为1,其他点为0的函数)就可以。作为检验的比较计算,采用使用的最多的以下所示的两个气泡函数以及形成正交基底的气泡函数。使用图11以及图12中所示的等参数座标系{r,s}({r,s,t}),如式(47)~式(52)进行定义。二维(三角形)较多地使用的气泡函数1φB=27(1-r-s)rs,⟨φB,1⟩Ωe=920Ae,]]>||φB||Ωe2=81280Ae]]>…式(47)较多地使用的气泡函数2φB=3(1-r-s)inw13rinw23sinw3,⟨φB,1⟩Ωe=13Ae,||φB||Ωe2=16Ae]]>…式(48)成为正交基底的气泡函数φB=α1φ^B+α2φ~B+φ‾Bα1+α2+1,]]>⟨φB,1⟩Ωe=34Ae,]]>||φB||Ωe2=34Ae]]>…式(49)[公式28]三维(四面体)较多地使用的气泡函数1φB=256(1-r-s-t)rst,⟨φB,1⟩Ωe=32105Ae,]]>||φB||Ωe2=819251975Ae]]>…式(50)较多地使用的气泡函数2φB=4(1-r-s-t)inw14rinw24sinw34tinw4,⟨φB,1⟩Ωe=14Ae,||φB||Ωe2=110Ae]]>…式(51)成为正交基底的气泡函数φB=α1φ^B+α2φ~B+φ‾Bα1+α2+1,]]>⟨φB,1⟩Ωe=45Ae,]]>||φB||Ωe2=45Ae]]>…式(52)<6.RotatingCone问题>作为发明方法的检验计算,进行用来对数值分析方法的性能进行比较讨论的、有名的基准点(benchmark)问题即RotatingCone问题(物质浓度的移流问题)的分析(参照非专利文献4)。图13为表示RotatingCone问题的分析中的计算区域的说明图,图14为表示RotatingCone问题的分析中的分析模型的网格(节点数4921,要素数9600)的说明图,图15为表示RotatingCone问题的分析中的分析模型的流速(流况)的说明图。图16为表示RotatingCone问题的分析中的初始条件(初始浓度分布)的概观图,图17为表示RotatingCone问题的分析中的初始条件(初始浓度分布)的等高线的说明图。该RotatingCone问题,是在非恒定移流方程式中的图13的计算区域中,在图14中假设图15的稳定的一定流速(流况),对如图16以及图17所示的物理量(污染物质浓度等)的移流状况进行分析的问题。实际的物理量(污染物质浓度等),浓度一边扩散且随着流况而被传送,但在消除了扩散的效果进行计算的情况下(计算移流方程式的情况下),在将浓度乘以流量(流れ)并包围了计算区域时,存在如下问题即由于只有使物理量移流的作用,因此物理量的分布理想的是必须与原来的初始分布相一致。如果通过该问题进行数值分析,就能够对作为对象的方法所具有的近似误差进行评价,检验计算精度。图18为表示使用气泡函数(式(47))中的匹配质量行列式的计算结果的概观图,图19为表示使用气泡函数(式(47))中的匹配质量行列式的计算结果的等高线的说明图,另外,图20为表示使用气泡函数(式(47))中的集中质量行列式的计算结果的概观图,图21为表示使用气泡函数(式(47))中的集中质量行列式的计算结果的等高线的说明图。图22为表示使用气泡函数(式(48))中的匹配质量行列式的计算结果的概观图,图23为使用表示气泡函数(式(48))中的匹配质量行列式的计算结果的等高线的说明图。另外,图24为表示使用气泡函数(式(48))中的集中质量行列式的计算结果的概观图,图25为表示使用气泡函数(式(48))中的集中质量行列式的计算结果的等高线的说明图。另外,图26为表示使用气泡函数(式(49))中的对角质量行列式的计算结果的概观图,图27为表示使用气泡函数(式(49))中的对角质量行列式的计算结果的等高线的说明图。作为分析方法,使用对成为检验对象的所有气泡函数,在空间方向的离散化中考虑了数值分析的稳定化的气泡函数要素稳定化法(参照非专利文献3),在时间方向的离散化中使用4段解法,微小时间Δt使用π/400。图18、图19、图22以及图23,是使用气泡函数(式(47)、式(48))中的匹配质量行列式的计算结果。对该计算结果,使用气泡函数(式(47)、式(48))中的集中质量行列式所得到的计算结果是,图20、图21、图24以及图25中,物理量(圆锥)的前进方向的后方中产生振动,圆锥的分布状况也被破坏。这是浓度的物理上无意义的衰减、扩散变得显著,浓度应当达到最大的地点也大幅偏差,且变成了数值分析上可靠性较低的计算结果。另外,气泡函数(式(49))中的使用对角质量行列式的计算结果是,圆锥的前进方向的后方没有发生振动,其分布也没有破坏,确保了与使用气泡函数(式(47))、式(48))的匹配质量行列式的计算结果同等的计算精度,得到了数值分析上的可靠性高的计算结果。<7.基于RotatingCone问题的三维分析的检验>作为发明方法的三维分析的检验,与二维分析一样,进行RotatingCone问题(物质浓度的移流问题)的分析(参照非专利文献4)。图28为表示RotatingCone问题的分析中的分析模型的网格(节点数438413,要素数2539200)的说明图。图28为在图13的二维分析模型中具有垂直长度=2(-1~1)的模型。分析模型的流速(流况)以及浓度分布,在垂直方向上没有变化,与二维分析模型的情况下一样。图29为表示RotatingCone问题的分析中垂直位置=0的剖面的初始条件(初始浓度分布)的概观图,图30为表示RotatingCone问题的分析中垂直位置=0部分的初始条件(初始浓度分布)的等高线的说明图。图31为表示使用气泡函数(式(50))中的匹配质量行列式的垂直位置=0部分的计算结果的概观图,图32为表示使用气泡函数(式(50))中的匹配质量行列式的垂直位置=0部分的计算结果的等高线的说明图。另外,图33为表示使用气泡函数(式(50))中的集中质量行列式的垂直位置=0部分的计算结果的概观图,图34为表示使用气泡函数(式(50))中的集中质量行列式的垂直位置=0部分的计算结果的等高线的说明图。图35为表示使用气泡函数(式(51))中的匹配质量行列式的垂直位置=0部分的计算结果的概观图,图36为表示使用气泡函数(式(51))中的匹配质量行列式的垂直位置=0部分的计算结果的等高线的说明图。另外,图37为表示使用气泡函数(式(51))中的集中质量行列式的垂直位置=0部分的计算结果的概观图,图38为表示使用气泡函数(式(51))中的集中质量行列式的垂直位置=0部分的计算结果的等高线的说明图。另外,图39为表示使用气泡函数(式(52))中的对角质量行列式的垂直位置=0部分的计算结果的概观图,图40为表示使用气泡函数(式(52))中的对角质量行列式的垂直位置=0部分的计算结果的等高线的说明图。作为分析方法,使用对成为检验对象的所有气泡函数,在空间方向的离散化中考虑了数值分析的稳定化的气泡函数要素稳定化法(参照非专利文献3),在时间方向的离散化中使用4段解法,微小时间Δt使用π/600。图31、图32、图35以及图36,是使用气泡函数(式(50)、式(51))中的匹配质量行列式的计算结果。对该计算结果,使用气泡函数(式(50)、式(51))中的集中质量行列式所得到的计算结果是,图33、图34、图37以及图38中,物理量(圆锥)的前进方向的后方中产生振动,圆锥的分布状况也被破坏。这是浓度的物理上无意义的衰减、扩散变得显著,浓度应当达到最大的地点也大幅偏差,且变成了数值分析上可靠性较低的计算结果。另一方面,气泡函数(式(52))中的使用对角质量行列式的计算结果是,圆锥的前进方向的后方没有发生振动,其分布也没有破坏,确保了与使用气泡函数(式(50))、式(51))的匹配质量行列式的计算结果同等的计算精度,得到了数值分析上的可靠性高的计算结果。<8.一维正交基底气泡函数要素与要素等级的质量行列式>以上发明了二维、三维正交基底气泡函数要素并示出了其效果。但有限元法的数值计算中,虽然成为问题的分析几乎都是二维或三维的,但正交基底气泡函数要素也能够定义为一维的。这里,示出一维气泡函数要素,与ξ次方气泡函数,使用空间维数N对正交条件式(14)(或式(21))统一进行描述。另外,示出了一维(线性气泡函数要素)的要素等级的质量行列式。图41为表示等参数座标系r的说明图(之1)。图41中所示的等参数座标系r,如果通过式(1)、(2)的表现来描述,则一维的ΦT、uT、ψα变为下述式(53)~(55)。一维ΦT=[Φ1Φ2φB]…式(53)uT=[u1u2uB]…式(54)ψ1=1-r,ψ2=r…式(55)图42为表示等参数座标系r的说明图(之2)。图42中所示的等参数座标系r中,使用中间点将线性要素区域分割为w1、w2,气泡函数的要素区域中的积分值变为式(29)的一维ξ次方气泡函数,如下式(56)进行定义。一维φBξ=2ξ(1-r)ξinw12ξrξinw2]]>…式(56)式(56)、(27)、(28)的ξ(η)次方气泡函数(0<η<∞),如下述式(57)所示,对气泡函数在空间中进行偏微分所得到的积分值也能够统一描述,保持下述式(58)、(59)所示的积分值的关系。⟨∂φBξ∂xk,∂φBη∂xl⟩Ωe=(N+1)N!ξηΠγ=-1N-2(ξ+η+γ)Ae(Σα=1N+1∂Ψα∂xk∂Ψα∂xl)]]>…式(57)⟨Ψα,φBϵ⟩Ωe=1N+1⟨1,φBϵ⟩Ωe,α=1···N+1]]>…式(58)⟨Ψα,Σkd=1KDαkdφBϵkd⟩Ωe=Σkd=1KDαkd⟨Ψα,φBϵkd⟩Ωe=1N+1(Σkd=1KDαkd⟨1,φBϵkd⟩Ωe)]]>=1N+1⟨1,Σkd=1KDαkdφBϵkd⟩Ωe,α=1···N+1,0<ϵkd<∞]]>…式(59)式(57)的xk,x1=x1……xN,分别表示空间方向,式(59)的αkd表示实数,KD表示正的整数。式(59)通过式(58)导出。如果设C=(N+1)/(N+2),则如果使用空间维数N,将正交基底气泡函数要素的条件式(14)(或式(21))统一描述为一维~三维,就变为下述式(60)。⟨φB,1⟩Ωe=||φB||Ωe2=N+1N+2Ae]]>…式(60)下述式(61)~(63)中,示出了一维(线性)中的通常的气泡函数要素的要素等级的匹配质量行列式、集中质量行列式、以及正交基底气泡函数要素的要素等级的质量行列式(对角质量行列式)。线性气泡函数要素要素等级的匹配质量行列式⟨Φ,ΦT⟩Ωe=Mij(e)=Ae3Ae60Ae6Ae30000]]>+14⟨1,φB⟩Ωe-2-22-2-22220+14||φB||Ωe211-211-2-2-24]]>…式(61)要素等级的集中质量行列式(集中化后的情况)⟨Φ,ΦT⟩Ωe=Mij(e)≈diag(Σj=1N+2Mij(e))]]>=Ae2000Ae20000+12⟨1,φB⟩Ωe-1000-10002]]>…式(62)要素等级的对角质量行列式(成为正交基底的气泡函数)⟨Φ,ΦT⟩Ωe=Mij(e)=Ae6000Ae600023Ae]]>…式(63)<9.扩展后的正交基底气泡函数要素>使用气泡函数要素的流体分析、构造分析、固有值分析等中,一般来说,除了一次要素的形状函数的ψα的积分值之外,还需要下述式(64)~(66)的气泡函数的积分值。<ψα,φB>Ωe,α=1…N+1…式(64)‖φB‖Ωe2…式(65)⟨∂φB∂xk,∂φB∂xl⟩Ωe]]>…式(66)通过下述式(67)与式(64),得到下述关系式(68)。Σα=1N+1Ψα=1]]>…式(67)⟨1,φB⟩Ωe=⟨Σα=1N+1Ψα,φB⟩Ωe=Σα=1N+1⟨Ψα,φB⟩Ωe]]>…式(68)在气泡函数为具有下式(69)的关系的气泡函数的情况下,能够将必要的气泡函数的积分值(64)~(66)的式(64),代入到式(70)中。⟨Ψα,φB⟩Ωe=1N+1⟨1,φB⟩Ωe,α=1···N+1]]>…式(69)<φB,1>Ωe…式(70)式(69)的关系式,只要没有导入非专利文献3中所出现的特殊的气泡函数,由于例如在式(47)、(48)、(50)、(51)等一般广泛使用的气泡函数中成立,因此对于在正交条件式(21)的导出中出现的式(17)、(19)也使用并,假定式(22)是满足式(69)的气泡函数。在式(22)右边的气泡。在使用正交基底气泡函数要素的情况下,式(70)、(65)、(66)中,式(70)、(65)能够根据正交条件式(60)来得知,式(66)无法通过正交条件来决定。在使用式(56)、(27)、(28)的ξ次方气泡函数,形成了正交基底气泡函数的情况下,实际上式(66)由下述式(71)来表示。⟨∂φB∂xk,∂φB∂xl⟩Ωe=D(ξ1,ξ2,ξ3)Ae(Σα=1N+1∂Ψα∂xk∂Ψα∂xl)]]>…式(71)D(ξ1,ξ2,ξ3)=(N+1)N!(Σld=13αld)2Σmd=13Σnd=13αmdαndξmdξndΠγ=-1N-2(ξmd+ξnd+γ)]]>…式(72)ξ1=ξ^,]]>ξ2=ξ~,]]>ξ3=ξ,α1=α1(ξ1,ξ2,ξ3),α2=α2(ξ1,ξ2,ξ3),α3=1,α1(ξ1,ξ2,ξ3)=α1+(ξ1,ξ2,ξ3)=-b(ξ1,ξ2,ξ3)+b2(ξ1,ξ2,ξ3)-4a(ξ1,ξ2)c(ξ2,ξ3)2a(ξ1,ξ2)α1-(ξ1,ξ2,ξ3)=-b(ξ1,ξ2,ξ3)-b2(ξ1,ξ2,ξ3)-4a(ξ1,ξ2)c(ξ2,ξ3)2a(ξ1,ξ2)]]>从式(71)可以得知,D(ξ1、ξ2、ξ3)的值成为变量ξ1、ξ2、ξ3的函数。实际上,如果使用式(31)、(36)的值,计算式(72)(α1(ξ1、ξ2、ξ3)=α1+(ξ1、ξ2、ξ3)),则如下述式(73)~(75)所示,所计算出的值因所采用的变量而不同。一维D(110,15,34)=325.0049920···,]]>D(3,2,65)=5.9124806···]]>…式(73)二维D(110,15,34)=982.4790415···,]]>D(3,2,65)=14.3484025···]]>…式(74)三维D(110,15,34)=4347.7594268···,]]>D(3,2,65)=46.8184447···]]>…式(75)本发明中,作为正交基底气泡函数要素的扩展,新定义了能够利用通过正交条件式(60)所得到的积分值,决定式(66)的积分值的下述式(76)。⟨∂φB∂xk,∂φB∂xl⟩Ωe:=⟨∂φB∂xk,ΦT⟩Ωe⟨Φ,ΦT⟩Ωe-1⟨Φ,∂φB∂xl⟩Ωe]]>=(-⟨φB,1⟩Ωe∂ΨT∂xk0)⟨Φ,ΦT⟩Ωe-1(∂Ψ∂xl0⟨1,φB⟩Ωe)]]>=(N+1)3N+2Ae(Σα=1N+1∂Ψα∂xk∂Ψα∂xl)]]>…式(76)一维ψT=[ψ1ψ2]…式(77)二维ψT=[ψ1ψ2ψ3]…式(78)三维ψT=[ψ1ψ2ψ3ψ4]…式(79)将在与基底正交的条件式(60)中添加了式(76)的条件之后的下述式(80),作为扩展后的正交基底气泡函数要素的条件式。⟨φB,1⟩Ωe=||φB||Ωe2=⟨∂φB∂xk,∂φB∂xl⟩Ωe{(N+1)2(Σα=1N+1∂ψα∂xk∂ψα∂xl)}-1=N+1N+2Ae]]>…式(80)<10.满足扩展后的正交基底气泡函数要素的气泡函数>从式(71)可以得知,使用式(56)、(27)、(28)的ξ次方气泡函数而形成了正交基底气泡函数的式(66)的积分值,由变量ξ1、ξ2、ξ3的值来决定。本申请中,考虑满足扩展后的正交基底气泡函数要素的条件式(80)的变量ξ1、ξ2、ξ3,而示出了这样的变量实际存在。虽然变量为ξ1、ξ2、ξ3这三个,但考虑将其中的两个变量如下式(81)进行固定,作为ξ1的1变量问题,求出满足扩展后的正交基底气泡函数要素的条件的气泡函数。ξ2=ξ~=85,]]>ξ3=ξ-=25]]>…式(81)作为ξ1的1变量问题,定义了下式(82)。f(ξ1)=(α1(ξ1)+α2(ξ1)+1)2{D(ξ1)-(N+1)3(N+2)}]]>…式(82)α1(ξ1)=α1+(ξ1,85,25),]]>α2(ξ1)=α2(ξ1,85,25),]]>D(ξ1)=D(ξ1,85,25)]]>存在满足扩展后的正交基底气泡函数要素的条件的气泡函数,是指实数的集合R中,存在f(ξ1)=0的变量。现在,考虑下述式(83)~(85)所示的实数的集合R的部分集合S。一维1.8310≤S≤1.8315…式(83)二维5.2678≤S≤5.2683…式(84)三维1.8061≤S≤1.8066…式(85)在ξ1∈S中,一维~三维中变为下式(86)。b2(ξ1)-4a(ξ1)c>0…式(86)因此,α1、α2在实数的范围中具有值。另外,ξ1∈S中的一维~三维中,分别成立下述式(87)以及(88),因此式(72)取有界限的值,式(82)变为实变量的函数。α1(ξ1)+α2(ξ1)+1≠0…式(87)Πγ=-1N-2(ξmd+ξnd+γ)≠0,md,nd=1,2,3]]>…式(88)如果使用式(82)与部分集合S,就能够得到下述式(89)~(91)。一维f(1.8310)<0,f(1.8315)>0,df(ξ1c)dξ1c>0]]>…式(89)二维f(5.2678)>0,f(5.2683)<0,df(ξ1c)dξ1c<0]]>…式(90)三维f(1.8061)>0,f(1.8066)<0,df(ξ1c)dξ1c<0]]>…式(91)对于任意的ξ1c∈S,一维中式(89)的第1式为负,第2式为正,第3式为正(单调增加),二维、三维中,式(90)、(91)的第1式为正,第2式为负,第3式为负(单调减少),因此一维~三维中,部分集合S的范围内,f(ξ1)=0的解分别唯一存在。具体的说,部分集合S的范围内,如果使用二分法求出f(ξ1)的解,则得到下述式(92)~(94)的值。一维ξ1=ξ^=1.8312889···,]]>α1=20.4651928…,α2=-22.5761007……式(92)二维ξ1=ξ^=5.2680399···,]]>α1=-2.6201232…,α2=3.1609628……式(93)三维ξ1=ξ^=1.8063444···,]]>α1=-17.3739109…,α2=17.5493578……式(94)图43、图44中,示出了使用式(81)、(93)的三角形气泡函数的形状φB、三角形气泡函数的形状φB2。或者,使用下述式(95)、(96),能够示出存在f(ξ1)=0的解。g(ξ1)=ξ1-f(ξ1)f′(ξ1)]]>…式(95)|g′(ξ1)|=|f(ξ1)f′′(ξ1)(f′(ξ1))2|]]>…式(96)f′(ξ1)=df(ξ1)dξ1,]]>f′′(ξ1)=d2f(ξ1)dξ12,]]>g′(ξ1)=dg(ξ1)dξ1]]>一维~三维中,h为正的常数,对任意的ξ1c∈S,得到下述式(97)的结果。|g′(ξ1c)|<h<1…式(97)根据平均值的定理,对任意的ξ1a,ξ1b∈S,得到下式(98)。|g(ξ1a)-g(ξ1b)|=|g′(ξ1c)||ξ1a-ξ1b|<h|ξ1a-ξ1b|,0<h<1…式(98)因此,根据缩小映射原理(定理),由让式(95)反复逐次进行计算的下式(99)所得到的数列{ξ1m},集中为在部分集合S的范围内满足f(ξ1)=0的唯一的解。ξ1m+1=ξ1m-f(ξ1m)f′(ξ1m),m=0,1,2,3,···]]>…式(99)另外,为了表示出扩展后的正交基底气泡函数要素的存在,作为满足式(80)的三个条件的气泡函数,虽然将α1、α2、ξ1作为三个未知量,采用适当的气泡函数φB1、φB2、φB3,但是,其他的如下式(100)的气泡函数φB所示,如果将α1、α2、α3作为三个未知量,并采用适当的气泡函数φB1、φB2、φB3、φB4,也能够得到满足式(80)的条件的气泡函数。这样,满足扩展后的正交基底气泡函数要素的条件式的气泡函数,并不是唯一的,而存在有好几个,具体的气泡函数及其形状几乎没有意义,扩展后的正交基底气泡函数要素中的重要意义的,是变为式(80)的条件的气泡函数的积分值。φB=α1φB1+α2φB2+α3φB3+φB4α1+α2+α3+1]]>…式(100)<11.扩展后的正交基底气泡函数要素的要素等级的粘性项(扩散项)的行列式>通常的气泡函数、正交基底气泡函数要素的数值分析中,式(66)所必需的要素等级的粘性项(扩散项)的行列式如下述式(101)~(103)所示。一维⟨∂Φ∂xk,∂ΦT∂xl⟩Ωe=Ae∂Ψ1∂xk∂Ψ1∂xl∂Ψ1∂xk∂Ψ2∂xl0∂Ψ2∂xk∂Ψ1∂xl∂Ψ2∂xk∂Ψ2∂xl0000]]>+⟨∂φB∂xk,∂φB∂xl⟩Ωe1411-211-2-2-24]]>…式(101)二维⟨∂Φ∂xk,∂ΦT∂xl⟩Ωe=Ae∂Ψ1∂xk∂Ψ1∂xl∂Ψ1∂xk∂Ψ2∂xl∂Ψ1∂xk∂Ψ3∂xl0∂Ψ2∂xk∂Ψ1∂xl∂Ψ2∂xk∂Ψ2∂xl∂Ψ2∂xk∂Ψ3∂xl0∂Ψ3∂xk∂Ψ1∂xl∂Ψ3∂xk∂Ψ2∂xl∂Ψ3∂xk∂Ψ3∂xl00000]]>+⟨∂φB∂xk,∂φB∂xl⟩Ωe19111-3111-3111-3-3-3-39]]>…式(102)三维⟨∂Φ∂xk,∂ΦT∂xl⟩Ωe=Ae∂Ψ1∂xk∂Ψ1∂xl∂Ψ1∂xk∂Ψ2∂xl∂Ψ1∂xk∂Ψ3∂xl∂Ψ1∂xk∂Ψ4∂xl0∂Ψ2∂xk∂Ψ1∂xl∂Ψ2∂xk∂Ψ2∂xl∂Ψ2∂xk∂Ψ3∂xl∂Ψ2∂xk∂Ψ4∂xl0∂Ψ3∂xk∂Ψ1∂xl∂Ψ3∂xk∂Ψ2∂xl∂Ψ3∂xk∂Ψ3∂xl∂Ψ3∂xk∂Ψ4∂xl0∂Ψ4∂xk∂Ψ1∂xl∂Ψ4∂xk∂Ψ2∂xl∂Ψ4∂xk∂Ψ3∂xl∂Ψ4∂Ψk∂Ψ4∂xl000000]]>+⟨∂φB∂xk,∂φB∂xl⟩Ωe1161111-41111-41111-41111-4-4-4-4-416]]>…式(103)将式(101)~(103)代入到式(76)中并进行整理之后的扩展过后正交基底气泡函数要素的要素等级的粘性项(扩散项)的行列式,如下述式(104)~(106)所示。一维⟨∂Φ∂xk,∂ΦT∂xl⟩Ωe=Ae∂Ψ1∂xk∂Ψ1∂xl∂Ψ1∂xk∂Ψ2∂xl0∂Ψ2∂xk∂Ψ1∂xl∂Ψ2∂xk∂Ψ2∂xl0000]]>+23Ae(Σα=1N+1∂Ψα∂xk∂Ψα∂xl)11-211-2-2-24]]>…式(104)二维⟨∂Φ∂xk,∂ΦT∂xl⟩Ωe=Ae∂Ψ1∂xk∂Ψ1∂xl∂Ψ1∂xk∂Ψ2∂xl∂Ψ1∂xk∂Ψ3∂xl0∂Ψ2∂xk∂Ψ1∂xl∂Ψ2∂xk∂Ψ2∂xl∂Ψ2∂xk∂Ψ3∂xl0∂Ψ3∂xk∂Ψ1∂xl∂Ψ3∂xk∂Ψ2∂xl∂Ψ3∂xk∂Ψ3∂xl00000]]>+34Ae(Σα=1N+1∂Ψα∂xk∂Ψα∂xl)111-3111-3111-3-3-3-39]]>…式(105)三维⟨∂Φ∂xk,∂ΦT∂xl⟩Ωe=Ae∂Ψ1∂xk∂Ψ1∂xl∂Ψ1∂xk∂Ψ2∂xl∂Ψ1∂xk∂Ψ3∂xl∂Ψ1∂xk∂Ψ4∂xl0∂Ψ2∂xk∂Ψ1∂xl∂Ψ2∂xk∂Ψ2∂xl∂Ψ2∂zk∂Ψ3∂xl∂Ψ2∂xk∂Ψ4∂xl0∂Ψ3∂xk∂Ψ1∂xl∂Ψ3∂xk∂Ψ2∂xl∂Ψ3∂xk∂Ψ3∂xl∂Ψ3∂xk∂Ψ4∂xl0∂Ψ4∂xk∂Ψ1∂xl∂Ψ4∂xk∂Ψ2∂xl∂Ψ4∂xk∂Ψ3∂xl∂Ψ4∂xk∂Ψ4∂xl000000]]>+45Ae(Σα=1N+1∂Ψα∂xk∂Ψα∂xl)1111-41111-41111-41111-4-4-4-4-416]]>…式(106)以上,示出了扩展后的正交基底气泡函数要素的要素等级的粘性项(扩散项)的行列式(104)~(106),但本发明的相关公式表述或程序,也可以不将式(76)代入到行列式(101)~(103)中进行整理,而是直接使用行列式(101)~(103)。<12.正交基底气泡函数要素与正规化气泡函数>式(1)能够如下式(107)、(108)进行变形。uh|Ωe=Σα=1N+1Ψαuα+φBuB′]]>…式(107)=Σα=1N+1Ψαuα+φSuS′]]>…式(108)uB′=uB-1N+1(Σα=1N+1uα),]]>uS′=⟨1,φB⟩ΩeAeuB′=N+1N+2uB′]]>…式(109)φS=Ae⟨1,φB⟩ΩeφB=N+2N+1φB]]>…式(110)式(108)是使用气泡函数自身的积分值变为Ae(一维中为各个要素e的线,二维中为各个要素e的面积,三维中为各个要素e的体积)的气泡函数,进行公式变形而得到的,式(110)称作正规化气泡函数。正交基底气泡函数要素中,除了式(1)的表现形式之外,还能够使用式(107)的表现形式或式(108)的表现形式,进行定型化或编程、执行计算。另外,在使用式(1)、(107)、(108)的表现形式所得到的有限元方程式(变分方程式)中,能够通过静态压缩操作来将每一个要素删除重心点自由度uB、uB’、uS’,但也可以使用删除了重心点的自由度的有限元方程式(变分方程式),执行定型化或编程、以及计算执行。使用正规化气泡函数作为正交基底气泡函数要素的气泡函数所得到的积分值,如下述式(111)、(112)所示。⟨1,φS⟩Ωe=N+2N+1⟨1,φB⟩Ωe=Ae]]>…式(111)||φS||Ωe2=(N+2)2(N+1)2||φB||Ωe2=N+2N+1Ae]]>…式(112)进而,扩展过的正交基底气泡函数要素中,得到下述式(113)。⟨∂φS∂xk,∂φS∂xl⟩Ωe=(N+2)2(N+1)2⟨∂φB∂xk,∂φB∂xl⟩Ωe=(N+1)(N+2)Ae(Σα=1N+1∂Ψα∂xk∂Ψα∂xl)]]>…式(113)<13.正交基底气泡函数要素与稳定化作用控制项>为了缓和气泡函数要素所具有的数值不稳定性,存在使用下述式(114)的气泡函数的粘性项(稳定化作用控制项)而提高分析精度的方法。v′⟨∂φB∂xk,∂φB∂xl⟩ΩeuB′]]>…式(114)式(114)的稳定化作用控制参数v’(-∞<v’<∞),是具有粘性系数的维数(在将相当于粘性系数的变量无维数化了的情况下,无维数)的参数,通过适当决定其值,能够实现高精度解法。该方法以往使用气泡函数要素中通常的Galerkin法,只给所得到的有限元方程式添加重心点的粘性项(稳定化作用控制项)(参照松本純一,梅津剛,“気泡関数を用ぃた非压縮性流体の有限元法分析に関する一考察”,日本応用数理学会平成8年度年会,1996年,46页~47页、松本純一,梅津剛,川原睦人,“線形型気泡関数を用ぃた非压縮性粘性流体分析と適応型有限元法”,応用力学論文集(土木学会),2卷,1999年8月,223页-232页),当前,关于权重函数也采用1个气泡函数(3级气泡函数),只导入重心点的粘性项(稳定化作用控制项)(参照J.MatumotoandM.Kawahara,“ShapeIdentificationforFluid-StructureInteractionProblemUsingImprovedBubbleElement”,InternationalJournalofComputationalFluidDynamics,Vol.15,2001,pp.33-45、非专利文献3)。3级气泡函数是满足下式(115)~(117)的气泡函数。…式(115)_B为3级气泡函数…式(116)…式(117)通过式(67)、(115)得到下式(118)。…式(118)另外,在为具有下述式(119)的关系的3级气泡函数的情况下,条件式(115)~(117)中式(115)能够由式(120)代用。…式(119)…式(120)假设具有式(119)的关系的3级气泡函数,作为满足上述条件式(120)、(116)、(117)的具体的3级气泡函数,使用下式(121)。…式(121)式中的v是粘性系数(扩散系数)。将式(121)代入到式(120)、(116)、(117)中并进行整理,便得到下式(122)。δ1δ2δ3=a11a12a13a21a22a23a31a32a33-1b1b2b3]]>…式(122)或者,式(121)右边的4个气泡函数,具体定义为下式(123)。…式(123)根据式(123)、(59),可以得知式(121)满足式(119)。式(123)的变量定义为下式(124)。η1=4,η2=3,η3=2,η4=1…式(124)如果使用式(31)的变量作为α1(ξ1、ξ2、ξ3)=α1+(ξ1、ξ2、ξ3),则通过式(122)得到下式(125)~(127)。一维δ1=-7.4533261…,δ2=13.9037038…,δ3=-7.4557822……式(125)二维δ1=-12.8855927…,δ2=21.8201295…,δ3=-9.9378406……式(126)三维δ1=73.7388138…,δ2=-102.3508563…,δ3=27.6071956……式(127)如果使用式(36)的变量作为α1(ξ1、ξ2、ξ3)=α1+(ξ1、ξ2、ξ3),则通过式(122)得到下式(128)~(130)。一维δ1=-5.2897450…,δ2=10.5097066…,δ3=-6.2084329……式(128)二维δ1=-9.2351821…,δ2=16.2678560…,δ3=-8.0666407……式(129)三维δ1=-24.4550047…,δ2=37.6758731…,δ3=-14.3507923……式(130)如果使用式(92)~(94)、式(81)的变量作为α1(ξ1、ξ2、ξ3)=α1+(ξ1、ξ2、ξ3),则由式(122)得到下式(131)~(133)。一维δ1=-6.8577197…,δ2=12.6283403…,δ3=-6.8566234……式(131)二维δ1=-7.0364016…,δ2=13.3033720…,δ3=-7.1674625……式(132)三维δ1=-12.0454679…,δ2=20.3290654…,δ3=-9.2229704……式(133)图45、图46中示出了使用了式(132)的三角形3级气泡函数的形状(v’=v的情况),以及使用了式(93)、(81)的三角形3级气泡函数与使用了式(132)的三角形3级气泡函数的乘积的形状(v’=v的情况)。根据以上,即使在使用正交基底气泡函数要素的情况下,也存在满足式(120)(或式(115)、(116)、(117))的3级气泡函数要素,能够导入稳定化作用控制项(114)。另外,由于示出了3级气泡函数的存在,而采用了具有满足式(120)、(116)、(117)这3个条件的3个未知量δ1、δ2、δ3与4个适当的气泡函数的3级气泡函数(121),但此外还可以如下式(134)所示,导入具有与稳定化作用控制参数v’相同维数的参数,在式(134)的右边使用适当的3个气泡函数作为未知量δ1、δ2,求出满足式(120)、(116)的条件的3级气泡函数。之后,如果如下式(135)所示定义稳定化作用控制参数v’,则还能导出式(117),从而能够示出成为式(120)(或式(115))、(116)、(117)的3级气泡函数的存在。这样,满足3级气泡函数的条件式的气泡函数并不是唯一的,而是存在有好几个,因此具体的气泡函数及其形状基本上没有意义,3级气泡函数中具有重要意义的是,为了导入稳定化作用控制项而必需的式(120)或(式(115))、(116)、(117)这些条件。…式(134)其中v_(-∞<v_<∞)是与v′具有相同维数的参数…式(135)<14.基于正交基底气泡函数要素与正规化气泡函数的稳定化作用控制项>在使用正规化气泡函数的情况下,也能够导入稳定化作用控制项,并通过下式(136)来表示。v′⟨∂φS∂xk,∂φS∂xl⟩ΩeuS′]]>…式(136)在采用正规化气泡函数作为正交基底气泡函数要素中所使用的气泡函数的情况下,3级正规化气泡函数能够如下式(137)进行定义,其积分值可描述为下式(138)~(142),[公式71]…式(137)…式(138)…式(139)=v′⟨∂φS∂xk,∂φS∂xl⟩Ωe=Ae2⟨1,φB⟩Ωe2v′⟨∂φB∂xk,∂φB∂xl⟩Ωe]]>…式(140)…式(141)…式(142)根据上述式(138)~(142),如果存在满足式(115)~(117)、(119)、(120)的3级气泡函数,则式(138)~(142)也自动满足,因此可以得知,3级正规化气泡函数也同时存在。<15.正交基底气泡函数要素与气泡函数的稳定化作用行列式>采用式(107)的表现形式的气泡函数要素的重心点的自由度uB’,能够通过称作静态压缩的操作来对每一个要素删除。结果所得到的有限元方程式,与通过稳定化有限元法所导出的有限元方程式之间,存在密切的关系(参照非专利文献2)。通过静态压缩所得到的气泡函数的稳定化作用行列式,一般通过下式(143)的形式来表示。τ~eB=K~(v′)B-1⟨1,φB⟩2Ae]]>…式(143)…式(144)L是应当求出的位置分析物理量的数目。式(144),是在使用了满足式(69)的气泡函数的情况下,由式(70)、(65)、(66)的气泡函数的积分值所构成的行列式。稳定化作用控制参数v’,被决定为使得下式(145)的关系满足。τ~eBu′B≈τ~eRu′B]]>…式(145)…式(146)uB′T=uB1′···uBL′]]>…式(147)式(146)是能够期待提高分析精度的稳定化行列式。为了求出使得式(145)的关系满足的稳定化作用控制参数v’,而考虑决定使得下式(148)的评价函数最小的稳定化作用控制参数v’这一问题。J=12υ(v′)BTυ(v′)B]]>…式(148)υ(v′)B=[K~(v′)Bτ~eRuB′-⟨1,φB⟩2AeuB′]]]>…式(149)为了求出使得式(148)最小的稳定化作用控制参数v’,而采用作为式(148)的停留条件的下述关系式(150),求出稳定化作用控制参数v’。∂J∂v′=[∂υ(v′)B∂v′]Tυ(v′)B=0]]>…式(150)由于3级气泡函数在权重函数中导入,因此能够使用对于每一个应当求出的未知分析物理量的偏微分方程式而不同的3级气泡函数。因此,稳定化作用控制参数v’如果使用下式(151)的3级气泡函数,则能够只定义应当求出的未知分析物理量的数目L个的稳定化作用控制参数vm’。…式(151)式(151)的系数δ1、δ2、δ3能够由式(122)来决定,式(151)满足式(120)、(116)。另外,代替式(117)得到下述式(152)的稳定化作用控制项。…式(152)在采用式(152)的稳定化作用控制项的情况下,式(144)变为下式(153)。…式(153)为了求出满足式(145)的关系的稳定化作用控制参数vm’,而考虑对使得下式(154)的评价函数最小的稳定化作用控制参数vm’进行确定的问题。J=12υ(vm′)BTυ(vm′)B,m=1···L]]>…式(154)υ(vm′)B=[K~(vm′)Bτ~eRuB′-⟨1,φB⟩2AeuB′]]]>…式(155)为了求出使得式(154)最小的稳定化作用控制参数vm’,而采用作为式(154)的停留条件的下述关系式(156),求出稳定化作用控制参数vm’。∂J∂vm′=[∂υ(vm′)B∂vm′]Tυ(vm′)B=0]]>…式(156)式(155)是矢量,将矢量的成分表示为下式(157)。υ(vm′)BT=υ(v1′)B1···υ(vL′)BL]]>…式(157)从式(157)可以得知,由于稳定化作用控制参数vm’在各个矢量成分中独立定义,因此作为对使得式(145)的关系满足的稳定化作用控制参数vm’进行确定的这一问题,可以考虑对式(157)的每一个成分m,使用下式(158)的评价函数的最小化问题。J=12{υ(vm′)Bm}2,m=1···L]]>…式(158)为了求出使得式(158)最小的稳定化作用控制参数vm’,而采用作为式(158)的停留条件的下述关系式(159),求出稳定化作用控制参数vm’。∂J∂vm′=[∂υ(vm′)Bm∂vBm′]υ(vm′)m=0]]>…式(159)式(157)中,稳定化作用控制参数vm’在各个矢量成分中独立定义,因此用来求出从式(154)、(158)的停留条件所导出的稳定化作用控制参数vm’的式(156)、(159),结果变为同一个式子。作为其他方法,还可以考虑以使得下式(160)的关系满足的方式而对稳定化作用控制参数v’进行确定的问题。τ~eB-1uB′≈τ~eR-1uB′]]>…式(160)为了求出使得式(160)的关系满足的稳定化作用控制参数v’,而考虑决定使得下式(161)的评价函数最小的稳定化作用控制参数v’这一问题。J=12w(v′)BTw(v′)B]]>…式(161)w(v′)B=[K~(v′)BuB′-⟨1,φB⟩2Aeτ~eR-1uB′]]]>…式(162)为了求出使得式(161)最小的稳定化作用控制参数v’,而采用作为式(161)的停留条件的下述关系式(163),求出稳定化作用控制参数v’。∂J∂v′=[∂w(v′)B∂v′]Tw(v′)B=0]]>…式(163)对于稳定化作用控制参数v’,如果使用式(151)的3级气泡函数,则能够导出式(152)的稳定化作用控制项,只定义应当求出的未知分析物理量的数目L个的稳定化作用控制参数vm’。为了求出满足式(160)的关系的稳定化作用控制参数vm’,而考虑对使得下式(164)的评价函数最小的稳定化作用控制参数vm’进行确定的这一问题。J=12w(vm′)BTw(vm′)B,m=1···L]]>…式(164)w(vm′)B=[K~(vm′)BuB′-⟨1,φB⟩2Aeτ~eR-1uB′]]]>…式(165)为了求出使得式(164)最小的稳定化作用控制参数vm’,而采用作为式(164)的停留条件的下述关系式(166),求出稳定化作用控制参数vm’。∂J∂vm′=[∂w(vm′)B∂vm′]Tw(vm′)B=0]]>…式(166)式(165)是矢量,矢量的成分表示为下式(167)。w(vm′)BT=w(v1′)B1···w(vL′)BL]]>…式(167)从式(167)可以得知,由于稳定化作用控制参数vm’在各个矢量成分中被独立定义,因此作为对使得式(160)的关系满足的稳定化作用控制参数vm’进行确定的这一问题,可以考虑对式(167)的每一个成分m,使用下式(168)的评价函数的最小化问题。J=12{w(vm′)Bm}2,m=1···L]]>…式(168)为了求出使得式(168)最小的稳定化作用控制参数vm’,而采用作为式(168)的停留条件的下述关系式(169),求出稳定化作用控制参数vm’。∂J∂vm′=[∂w(vm′)Bm∂vm′]w(vm′)Bm=0]]>…式(169)式(167)中,稳定化作用控制参数vm’在各个矢量成分中独立定义,因此用来求出从式(164)、(168)的停留条件所导出的稳定化作用控制参数vm’的式(166)、(169),结果变为同一个式子。<16.正交基底气泡函数要素与正规化气泡函数的稳定化作用行列式>采用式(108)的表现形式的气泡函数要素的重心点的自由度uS’,能够通过称作静态压缩的操作来对每一个要素删除。通过静态压缩所得到的正规化气泡函数的稳定化作用行列式,一般通过下式(170)的形式来表示。τ~eS=K~(v′)S-1⟨1,φS⟩Ωe2Ae=K~(v′)B-1⟨1,φB⟩Ωe2Ae=τ~eB]]>…式(170)K~(v′)S=K~(v′)BAe2⟨1,φB⟩Ωe2]]>…式(171)式(170)的第二项,是在使用了满足式(141)的关系的气泡函数的情况下,由式(111)~(113)各自的第一项的正规化气泡函数的积分值所构成的行列式,在通过式(70)、(65)、(66)来表示的情况下,能够像第三项那样进行描述。稳定化作用控制参数v’,被确定为使得下式(172)的关系满足。τ~eSuS′≈τ~eRuS′]]>…式(172)uS′=⟨1,φB⟩ΩeAeuB′]]>…式(173)为了求出使得式(172)的关系满足的稳定化作用控制参数v’,而考虑对使得下式(174)的评价函数最小的稳定化作用控制参数v’进行确定这一问题。J=12υ(v′)STυ(v′)S]]>…式(174)υ(v′)S=[K~(v′)Sτ~eRuS′-⟨1,φS⟩2AeuS′]]]>…式(175)为了求出使得式(174)最小的稳定化作用控制参数v’,而采用作为式(174)的停留条件的下述关系式(176),求出稳定化作用控制参数v’。∂J∂v′=[∂υ(v′)S∂v′]Tυ(v′)S=0]]>…式(176)由于3级正规化气泡函数也在权重函数中导入,因此能够使用在每一个应当求出的未知分析物理量的偏微分方程式中不同的3级气泡函数。因此,稳定化作用控制参数v’如果使用下式(177)的3级气泡函数,则能够只定义应当求出的未知分析物理量的数目L个的稳定化作用控制参数vm’。…式(177)式(177)的系数δ1、δ2、δ3能够由式(122)来决定,式(177)满足式(142)、(139)。另外,代替式(140)得到下式(178)的稳定化作用控制项。…式(178)在采用式(178)的稳定化作用控制项的情况下,式(171)变为下式(179)。K~(vm′)S=K~(vm′)BAe2⟨1,φB⟩Ωe2]]>…式(179)式(179)的第一项,是由式(111)~(113)各自的第一项的正规化气泡函数的积分值所构成的行列式,在通过式(70)、(65)、(66)来表示的情况下,能够像第二项那样进行描述。为了求出满足式(172)的关系的稳定化作用控制参数vm’,而考虑对使得下式(180)的评价函数最小的稳定化作用控制参数vm’进行扩展这一问题。J=12υ(vm′)STυ(vm′)S,m=1···L]]>…式(180)υ(vm′)S=[K~(vm′)Sτ~eRuS′-⟨1,φS⟩2AeuS′]]]>…式(181)为了求出使得式(180)最小的稳定化作用控制参数vm’,而采用作为式(180)的停留条件的下述关系式(182),求出稳定化作用控制参数vm’。∂J∂vm′=[∂υ(vm′)S∂vm′]Tυ(vm′)S=0]]>…式(182)式(181)是矢量,矢量的成分表示为下式(183)。υ(vm′)ST=u(v1′)S1···u(vL′)SL]]>…式(183)从式(183)可以得知,由于稳定化作用控制参数vm’在各个矢量成分中独立定义,因此作为对使得式(172)的关系满足的稳定化作用控制参数vm’进行确定的这一问题,可以考虑如下问题即对式(183)的每一个成分m,使用下式(158)的评价函数的进行最小化。J=12{υ(vm′)Sm}2,m=1···L]]>…式(184)为了求出使得式(184)最小的稳定化作用控制参数vm’,而采用作为式(184)的停留条件的下述关系式(185),求出稳定化作用控制参数vm’。∂J∂vm′=[∂υ(vm′)Sm∂vm′]υ(vm′)Sm=0]]>…式(185)式(183)中,稳定化作用控制参数vm’在各个矢量成分中独立定义,因此用来求出从式(180)、(184)的停留条件所导出的稳定化作用控制参数vm’的式(182)、(185),结果变为同一个式子。作为其他方法,还可以考虑以使得下式(186)的关系满足满足的方式而确定稳定化作用控制参数v’的问题。τ~eS-1uS′≈τ~eR-1uS′]]>…式(186)为了求出使得式(186)的关系满足的稳定化作用控制参数v’,而考虑对使得下式(187)的评价函数最小的稳定化作用控制参数v’进行确定的这一问题。J=12ω(v′)STω(v′)S]]>…式(187)ω(v′)S=[K~(v′)SuS′-⟨1,φS⟩2Aeτ~eR-1uS′]]]>…式(188)为了求出使得式(187)最小的稳定化作用控制参数v’,而采用作为式(187)的停留条件的下述关系式(189),求出稳定化作用控制参数v’。∂J∂v′=[∂w(v′)S∂v′]Tw(v′)S=0]]>…式(189)对于稳定化作用控制参数v’,如果使用式(177)的3级气泡函数,则能够导出式(178)的稳定化作用控制项,并只定义应当求出的未知分析物理量的数目L个的稳定化作用控制参数vm’。为了求出满足式(186)的关系的稳定化作用控制参数vm’,而考虑对使得下式(190)的评价函数最小的稳定化作用控制参数vm’进行确定的这一问题。J=12w(vm′)STw(vm′)S]]>…式(190)w(vm′)S=[K~(vm′)SuS′-⟨1,φS⟩2Aeτ~eR-1uS′]]]>…式(191)为了求出使得式(190)最小的稳定化作用控制参数vm’,而采用作为式(190)的停留条件的下述关系式(192),求出稳定化作用控制参数vm’。∂J∂vm′=[∂w(vm′)S∂vm′]Tw(vm′)S=0]]>…式(192)式(191)是矢量,矢量的成分表示为下式(193)。w(vm′)ST=w(v1′)S1···w(vL′)SL]]>…式(193)从式(193)可以得知,由于稳定化作用控制参数vm’在各个矢量成分中独立定义,因此作为对使得式(186)的关系满足的稳定化作用控制参数vm’进行确定的这一问题,可以考虑对式(193)的每一个成分m使用下式(194)的评价函数的最小化问题。J=12{w(vm′)Sm}2,m=1···L]]>…式(194)为了求出使得式(194)最小的稳定化作用控制参数vm’,而采用作为式(194)的停留条件的下述关系式(195),求出稳定化作用控制参数vm’。∂J∂vm′=[∂w(vm′)Sm∂vm′]w(vm′)Sm=0]]>…式(195)式(193)中,稳定化作用控制参数vm’在各个矢量成分中独立定义,因此用来求出从式(190)、(194)的停留条件所导出的稳定化作用控制参数vm’的式(192)、(195),结果变为同一个式子。<17.气泡函数与正规化气泡函数的稳定化作用控制参数的关系>正交基底气泡函数要素中,在采用了式(1)、(107)的表现形式的情况下,如果使用式(70)、(65)、(66)的积分,并采用式(150)、(156)、(159)、(163)、(166)、(169)中的任一个,在采用了式(108)的表现形式的情况下,如果使用式(111)~(113)各自的第一项积分,并采用式(176)、(182)、(185)、(189)、(192)、(195)中的任一个,就能够计算出稳定化作用控制参数。现在,如果使用式(70)、(65)、(66)的积分,对式(176)、(182)、(185)、(189)、(192)、(195)进行描述,则变为下述式(196)~(201)那样。[∂u(v′)S∂v′]Tυ(v′)S=Ae2⟨1,φB⟩Ωe2[∂υ(v′)B∂v′]Tυ(v′)B=0]]>…式(196)[∂u(vm′)S∂vm′]Tυ(vm′)S=Ae2⟨1,φB⟩Ωe2[∂υ(vm′)B∂vm′]Tυ(vm′)B=0]]>…式(197)[∂u(vm′)Sm∂vm′]Tυ(vm′)Sm=Ae2⟨1,φB⟩Ωe2[∂υ(vm′)B∂vm′]Tυ(vm′)Bm=0]]>…式(198)[∂w(v′)S∂v′]Tw(v′)S=Ae2⟨1,φB⟩Ωe2[∂w(v′)B∂v′]Tw(v′)B=0]]>…式(199)[∂w(vm′)S∂vm′]Tw(vm′)S=Ae2⟨1,φB⟩Ωe2[∂w(vm′)B∂vm′]Tw(vm′)B=0]]>…式(200)[∂w(vm′)Sm∂vm′]Tw(vm′)Sm=Ae2⟨1,φB⟩Ωe2[∂w(vm′)Bm∂vm′]Tw(vm′)Bm=0]]>…式(201)通过使用了式(70)、(65)、(66)的积分的式子(196)~(201)的第二项,求出稳定化作用控制参数的关系式,变得与用来求出式(150)、(156)、(159)、(163)、(166)、(169)的稳定化作用控制参数的关系式相同(相等)。因此可以得知根据式(150)、(156)、(159)、(163)、(166)、(169),使用式(70)、(65)、(66)的气泡函数的积分所求出的稳定化作用控制参数,与根据式(176)、(182)、(185)、(189)、(192)、(195),使用式(111)~(113)各自的第一项的正规化气泡函数的积分而求出的稳定化作用控制参数,变为相同的值(相等的值)。如上所述,在1)气泡函数、使得气泡函数自身的积分值变为Ae的正规化气泡函数、2)为了缓和气泡函数所具有的数值不稳定而导入的稳定化作用控制项、3级气泡函数、3级标准化气泡函数、以及3)稳定化作用控制参数的计算中所使用的稳定化作用行列式的公式化或编程中,需要式(70)、(65)、(66)的气泡函数的积分值,或式(111)~(113)各自的第一项的正规化气泡函数的积分值。扩展后的正交基底气泡函数要素,其所需要的积分值都已事先(不进行复杂的积分)由式(80)给出,因此在具有与匹配质量行列式同等的分析精度的对角质量行列式的生成,或对气泡函数所具有的数值不稳定性进行改良的稳定化作用控制项的导入过程中所出现的一系列的作业中,能够非常合理地进行处理。关于正交基底气泡函数要素,从二维、三维的分析结果可以得知,采用存储容量、计算时间的效率较好的解法(例如4段解法),能够得到具有可靠性的结果。进而,根据扩展过的正交基底气泡函数要素的条件式(80),不但是喷墨的液滴动作分析、都市或河流的泛滥分析、地震等所产生的海啸的分析,在构造物的应力分析、汽车发动机的固有振动分析等中,对于为了改良所需要的气泡函数要素的积分值或气泡函数所具有的数值的不稳定性而必要的处理,能够合理且统一地进行处理,因此还具有能够容易地进行在上述多种多样的数值分析中的应用的这一效果。如上所述,通过本实施方式,在使用气泡函数要素的数值分析中,能够实现存储容量、计算时间等计算效率较好的分析方法(正交基底气泡函数要素数值分析方法)。另外,本实施方式中所说明的正交基底气泡函数要素数值分析方法,能够通过在个人计算机、工作站、PC机群、以及超级计算机等计算机中执行预先存储的程序来实现。该程序记录在硬盘、软盘、CD-ROM、MO、DVD等计算机可读存储媒体中,由计算机从存储媒体中读出并执行。另外,该程序还可以是能够经由互联网等网络进行传输的传输媒体。如上所述,本发明的相关正交基底气泡函数要素数值分析方法、正交基底气泡函数要素数值分析程序、以及正交基底气泡函数要素数值分析装置,在基于使用线、三角形、以及四面体形状的网格的有限元法的数值分析中,非常有用。权利要求1.一种正交基底气泡函数要素数值分析方法,其特征在于,包括取得工序,其中取得分析对象的各个要素的匹配质量行列式;生成工序,其中基于所述分析对象的各个要素的气泡函数,生成将通过所述取得工序所取得的各个要素的匹配质量行列式对角化了的各个要素的对角质量行列式;以及分析工序,其中根据所述分析对象的已知分析物理量,与通过所述生成工序所生成的各个要素的对角质量行列式,对所述分析对象的动作进行分析。2.如权利要求1所述的正交基底气泡函数要素数值分析方法,其特征在于所述生成工序,通过将所述各个要素中的所述气泡函数的积分值,代入到所述各个要素的匹配质量行列式中,来生成所述各个要素的对角质量行列式。3.如权利要求1或2所述的正交基底气泡函数要素数值分析方法,其特征在于,包括第1计算工序,其中根据所述生成工序所生成的各个要素的对角质量行列式,计算出所述分析对象全区域的对角质量行列式;以及第2计算工序,其中计算出由所述第1计算工序所计算出的所述分析对象全区域的对角质量行列式的逆行列式,所述分析工序中,根据所述分析对象的已知分析物理量,与所述分析对象全区域的对角质量行列式,以及所述第2计算工序所计算出的逆行列式,对所述分析对象的行为进行分析。4.一种正交基底气泡函数要素数值分析程序,其特征在于,在计算机中执行取得工序,其中所述取得分析对象的各个要素的匹配质量行列式;生成工序,其中根据所述分析对象的各个要素的气泡函数,生成将通过所述取得工序所取得的各个要素的匹配质量行列式对角化了的各个要素的对角质量行列式;以及分析工序,其中根据所述分析对象的已知分析物理量,与通过所述生成工序所生成的各个要素的对角质量行列式,对所述分析对象的动作进行分析。5.如权利要求4所述的正交基底气泡函数要素数值分析程序,其特征在于所述生成工序中,通过将所述各个要素中的所述气泡函数的积分值,代入到所述各个要素的匹配质量行列式中,来生成所述各个要素的对角质量行列式。6.如权利要求4或5所述的正交基底气泡函数要素数值分析程序,其特征在于,在计算机中执行第1计算工序,其中根据所述生成工序所生成的各个要素的对角质量行列式,计算出所述分析对象全区域的对角质量行列式;以及第2计算工序,其中计算出由所述第1计算工序所计算出的所述分析对象全区域的对角质量行列式的逆行列式,所述分析工序中,根据所述分析对象的已知分析物理量,与所述分析对象全区域的对角质量行列式,以及所述第2计算工序所计算出的逆行列式,对所述分析对象的行为进行分析。7.一种正交基底气泡函数要素数值分析装置,其特征在于,包括取得机构,其取得所述分析对象的各个要素的匹配质量行列式;生成机构,其基于所述分析对象的各个要素的气泡函数,生成将通过所述取得机构所取得的各个要素的匹配质量行列式对角化了的各个要素的对角质量行列式;以及分析机构,其基于所述分析对象的已知分析物理量,与通过所述生成机构所生成的各个要素的对角质量行列式,对所述分析对象的动作进行分析。8.如权利要求7所述的正交基底气泡函数要素数值分析程序,其特征在于所述生成机构,通过将所述各个要素中的所述气泡函数的积分值,代入到所述各个要素的匹配质量行列式中,来生成所述各个要素的对角质量行列式。9.如权利要求7或8所述的正交基底气泡函数要素数值分析装置,其特征在于,包括第1计算机构,其根据所述生成机构所生成的各个要素的对角质量行列式,计算出所述分析对象全区域的对角质量行列式;以及第2计算机构,其计算出由所述第1计算机构所计算出的所述分析对象全区域的对角质量行列式的逆行列式,所述分析机构,根据所述分析对象的已知分析物理量,与所述分析对象全区域的对角质量行列式,以及所述第2计算机构所计算出的逆行列式,对所述分析对象的行为进行分析。全文摘要本发明公开一种正交基底气泡函数要素数值分析方法及其装置,首先,通过第1取得部(202)取得分析对象的已知分析物理量(S401)。接下来,通过第2取得部(203),取得各个要素的要素等级的匹配质量行列式(S402)。接下来对每一个要素积分气泡函数(S403),并将步骤S403中积分得到的值代入到各个要素的要素等级的匹配质量行列式中,由此得出各个要素的要素等级的对角质量行列式(S404)。接下来,通过各个要素的要素等级的对角质量行列式的总和(合并),计算出分析对象全区域的对角质量行列式(S405)。之后,计算出该对角质量行列式的逆行列式(S406),根据分析对象的已知分析物理量、分析对象全区域的对角质量行列式、及其逆行列式,对分析对象的动作进行分析(S407)。文档编号G06F17/13GK101065751SQ20058004049公开日2007年10月31日申请日期2005年11月25日优先权日2004年11月26日发明者松本纯一申请人:独立行政法人科学技术振兴机构,独立行政法人产业技术综合研究所
网友询问留言 已有0条留言
  • 还没有人留言评论。精彩留言会获得点赞!
1