CN100348403C - 用于模拟流体向模具内腔中的注入的方法和装置 - Google Patents

Abstract

本发明涉及一种用于在考虑传热作用的情况下分析流体流动的方法和装置，具体地说，涉及从模制状态到固化状态的一种相态变化。特别是该方法和装置可以应用于分析注入模制处理，以由热塑性或热固性聚合物产生模制的聚合物部件，在一个实施例中，该方法可以用于确定填充模具内腔所需的压力和在填充和压实注入模具内腔期间所引起的压力变化梯度。这些分析结果可以用于确定浇口的数量和定位、确定该部件的最佳材料以及优化在模制过程中所应用的处理条件。

Description

用于模拟流体向模具内腔中的注入的方法和装置

技术领域

本发明涉及流体向内腔中的流入的三维模拟(modeling)的领域，更具体地说，在一个实施例中，涉及用于生产模制(molded)的聚合物部件的注入模制处理的模拟。

背景技术

近些年来，在很多工业中注入模制的塑料部件的使用已经急剧增加。电子设备、消耗商品、医药设备以及汽车部件的制造商正以前所未有的速度、越来越多地生产他们的产品和用在他们的产品中的塑料部件。同时，竞争压力正迫使注入模制业制造商找寻新的实现最佳设计的方法，以便和生产过程的设计更好地配合。当在设计开发过程后期才发现对于部件或模具结构需要改进时，实现所需的变化就会滞后和成本迅速上升。想保证它们的部件可生产并且实现最优化的一些公司，必须开始利用计算机辅助设计的工程技术来仿真或模拟注入模具中复杂的流动，以便更好地理解制造过程和在设计阶段的早期将这种知识结合到部件设计中。

当设计注入模具和要在其中制造的部件时，有很多因素应当考虑。如考虑一些参数例如部件总的几何尺寸、最小和最大壁厚、模具中的浇口的数目和位置(通过该浇口注入液态聚合物)、模具中的使内腔中的气体通过逸出的排气孔的数目和位置、聚合物成分和特性以及收缩容许量。由于密切的相关关系，单纯根据最终部件的形状和功能，部件和模具设计是不能可靠为依据的，而且还应该考虑制造过程的影响。

最好利用计算机辅助设计的工程仿真，以便为工程师提供关于在注入模制过程中在模具内腔内很可能发生的虚拟和数字反馈，使他们能更好地理解和预见所设想的部件设计中的动态过程，以便可以基本消除传统的高成本的接近制造的反复试验。计算机辅助设计的工程仿真的使用有利于在设计阶段优化部件设计、模具设计和制造过程参数，其中可以方便地以最低成本和紧凑的日程表实现所需改变。

对于注入模制方法和与高生产率、注模部件质量相关联的难题的主要讨论见诸于主要读物“Moldflow Design Principles：Quality and Productivityby Design”，由本专利申请的受让人Moldflow Pty.Ltd.，Kilsynh，Victoia分销，这里整个引用其公开文件，可供参考。

简而言之，该注模方法为复杂的两步法。在称为填充(filling)阶段的第一步，加压迫使聚合物材料进入模具内腔，直到内腔容积被充满。之后，在称为压实(packing)阶段的第二步，维持对聚合物的压力，以使聚合物材料进一步流入内腔，以补偿由于材料固化和收缩造成的收缩量。当部件充分固化时，可以从模具中脱出该部件。可以注入模制热塑性和热固性聚合物。

当模制热塑性聚合物时，将模具的内腔表面或腔壁温度维持在一低于注入的材料的熔点的温度。随着材料流入内腔，沿内腔腔壁液体材料形成固化层。这一层可以称为冷凝层，在填充过程中，取决于过程的条件和所使用的材料，该层厚度可变化。冷凝层的厚度是重要的，因为冷凝层会降低内腔中流动的有效的通道宽度，这是由于热塑性聚合物的热流变特性通常影响流动材料的粘度。

早期的分析仿真技术依靠二维有限元法，已发现该方法对于仿真相对简单的薄壁的部件的注入模制是有利的。例如在PCT/AU98/00130国际申请(转让给本申请的受让人)中讨论了更先进的仿真技术，这里整个引用其公开文件，可供参考。然而，在厚或复杂的部件中，其中熔化的塑料可以沿所有方向流动，传统的依靠规定厚度的平面区域的薄壁分析假设，通常不能预计这种类型的流动。为了实现高精确度和预见性，例如在其中将形成熔接流线，产生气阱并且流动将会超前或滞后的情况下，希望建立一种全三维仿真。

为了分析三维注入模制部件的设计，通常希望利用一计算机化的固体模拟组件例如Pro-Engtneer^TM、CATIA^TM、I-DEAS^TM、Solid Edge^TM或其它组件着手，该组件通常用在机械设计和制图(drafting)应用中。模拟组件可以用于产生部件几何特征的光仿真三维描述，称之为固体模型。目前，直接基于固体模型的有限元分析代码利用一些节点来定义有限元，例如四面体和六面体。为了说明在注入模制中涉及的物理过程，通常希望计算有限元模型中的每个节点的五个量，即压力、3个正交的矢量分量和温度。假定一个适宜的模型可包含几百个节点，求解这个复杂的数字问题是困难的，需要大量的计算机资源。

授予Nakano的5835379号美国专利和相关的EP0689467A1号欧洲专利中请，这里整个引用其公开文件，可供参考，其建议一种降低填充阶段在注入模制有限元模型中需确定的变量的数目的方法，以便使能利用与其它方法所需相比较少的计算机资源进行计算。Nakano讨论了流导(flowconductance)的概念，以将变量的数目降低到2，即压力和流导。按照这里所用的术语，可互换地使用“流导”、“流体传导”和“流度”并认为是同义的。利用流导变量将材料粘度的影响纳入该计算中。这就必然涉及粘度数据的外推，可以预料这将导致粘度计算中的明显的误差。此外在填充和压实阶段，接近内腔腔壁的材料开始冷凝，该层的厚度随时间增加直到该部分被脱出，被认为当采用Nakano的方法时，会导致附加误差。

经常测量作为温度和切变速率(shear rate)的函数的熔化的聚合物粘度η。当进行测量时不可能测量接近材料固化温度的低温度下的粘度，因为在这些低温度下，粘度相对高在合理的切变速率下，我们已经发现热粘滞耗散(viscous dissipation)是明显的。相应地，通常在这一温度范围内进行粘度测量，在该温度下，熔化的聚合物明显容易流动。在测量之后，形成一函数使之适合测量的数据。可以得到适合数据的函数通常在高温下是合理的，在很宽的切变范围内保持，然而，当考虑低于或接近固化温度的腔壁处的材料时，需要在实验范围之外适当外推粘度。这将导致粘度值的误差，又使得流导产生误差，因为粘度和流导是按照如下关系式相关的：

${\▿}^2\mathrm{\κ}=-\frac{1}{\mathrm{\η}}$ (1)

其也可以按如下形式表示：

*(ηκ)＝-1    (2)

此外，正如Nakano所指出的，流导在腔壁处具有小的数值，随着离腔壁的距离增加而增加。在腔壁处，流导的数值接近或等于0，这是因为有很高的粘度。在流体流动中的通常的粘度边界条件是在腔壁处粘度是0。当发生滑动时，在腔壁处粘度不是0，可以按几种方法处理。利用流导近似，可以将流导设置到一很小的非0的数值。但是，其结果实际上在冷凝层中流导的数值变化其数量量级小于熔化的流体。这一宽范围的数值变化导致在计算流导时产生Nakano方案和其它方案的数值误差。

在压实阶段，还可以通过首先对流导求解得到压力数值。在压实时，我们已确定，冷凝层变得十分明显，其程度使最终大部分节点低于固化温度。因此，在压实阶段，以上提到的问题是夸大的。

发明内容

常规的三维模拟技术主要依靠基本的质量守衡方程和动量守衡方程，未预测在注入模具内腔的流体流动。然而，为了实现高精度和高预测能力，在注入模制过程的完全三维仿真中，我们已经确定还建议应用能量守衡原理。为了同时满足3个要求，而不耗费非一般的计算机资源，已经根据本发明基于各种新的模拟假设和方法，开发了一些创建方法。

因此，本发明的一个目的是模拟在将腔壁中的流体流动时采用能量守衡原理。该流体可以是粘性液体例如注入模制过程中的熔化的聚合物，或在铸造过程中使用的熔化的或半固化的材料，已经发现在两个起始流动和固化过程中考虑基本的能量守衡原理的热影响是重要的。

例如，已确认在填充和压实阶段期间，有3个主要的传热机理要考虑：来自进入的熔化的材料对流、向模具腔壁的传导以及粘滞耗散，它们均与由在流动的聚合物内部的切变产生的热能有关。另外，还有其它的机理，例如压缩发热效应，是由于压缩产生的热量以及减压导致的冷却。所有3个主要机理均显著有助于在填充期间的能量均衡。因此，精确地仿真注入模制通常当熔化的聚合物流入模具时，需要对熔化的聚合物的非等温分析。然而，在压实阶段期间，材料的流动显著降低以及主要传热机理是向模具腔壁的传导和从高压区到低压区的熔化材料的对流。当最初加压时在最早压实阶段的对流作用十分强。还应指出，压实压力可以低于或高于填充压力，并且可以为一分布曲线，即可以使它们随时间变化。随之材料可以流入模具和经过一馈送系统从模具中流出，这取决于内腔中的压力与馈送系统处施加的压力之间的差。压实阶段的分析通常也需要将熔化材料看成可压缩的流体。

虽然，常规的方法确定在将该流体流动进入的区域所分割的多个小的元的流导，但根据本发明的一个实施例，计算在每个元的每个节点的流导。由于从一个节点到另一个节点，流导可以按各种各样所希望的方式中任何一种方式变化，这一技术使得能够确定更精确的流导数值。这些方式的差异取决于用于离散该问题的有限元的类型。与使用元的数值相比，更希望使用节点的数值，因为，作为点的数值，可更精确地确定其它量例如温度、压力和粘度，再者，一个节点到另一个节点可从外推。这种外推能力还能够更精确地确定固化的和熔化的聚合物之间的界面位置，这已经发现与仿真预测精度是特别相关的。

因此，本发明的一个目的是利用一维解析函数(function)，来描述在一节点的局部温度分布。本发明的另一个目的是确定一维解析函数随时间的变化，以说明热对流。本发明的再一个目的是确定一维解析函数随时间的变化以说明粘滞发热。本发明的再一个目的是利用一维解析函数确定明确的温度的对流模式。本发明的再一个目的是提供一种各向异性的有限元网格精细化，其包括一种到内腔腔壁的距离的计算算法。

此外，根据一个实施例，本发明纳入了冷凝层的概念，以便改善上述常规的仿真技术存在的某些问题。利用一判据来确定模具内腔腔壁处在固化的和熔化的聚合物之间的界面。这一判据可以是温度、最小粘度值、两者的组合，或者某些其它物理量。在冷凝区，设有明显的聚合物移动，这样就不需要计算该处的流导。可以从对温度、压力和流导的求解过程中除去在冷凝区中确定的节点。仅需要计算在冷凝层中的温度。然而，应指出，也可以计算在冷凝层中除温度以外的特性。例如，材料结构将随状态变化而变化。半结晶聚合物的结晶度随温度变化，结晶仅可说是一旦熔化的聚合物已从熔化状态充分冷却的一种形式。应力状态也取决于温度和液相与固相之间的差。聚合物可以是粘弹性的，这样即使在固化时它们的应力状态也会随时间和温度变化。

实际上，通过仅计算温度，这就减少了人们对于在冷凝层中的要确定的变量的数目，并因此显著地提高了计算速度，同时减少了所需的存储器。通过除去冷凝区的节点，消除了粘度的外推并使得还消除了与其相关的误差。此外，通过从分析中除去冷凝层的节点，消除了很小数值的流导，需求解的方程比常规的方法能更好地适应。此外，通过计算在这些区域中的流导，已经发现，不仅问题的规模明显地降低了，而且数字方法的稳定性也提高了。

因此，本发明的再一个目的是确定模具内腔中的以及在网格中定义的元的内部的固体/液体界面的位置。本发明的再一个目的是普遍改进按照材料特性的线性变化的元的公式化(formulation)。本发明的再一个目的是确定包括固体/液体界面的元中的有效粘度函数。本发明的再一个目的是除去冷凝层的节点和求解域中的元。本发明的再一个目的是计算已经固化的区域中的有效压力。

本发明的一个实施例涉及模拟流体向一限定三维内腔的模具的注入的方法，该方法包括步骤：提供一个三维固体的计算机模型，根据该固体模型离散一求解域，规定边界条件，利用质量、动量和能量守衡方程对填充阶段过程变量求解，以便提供对于至少该部分求解域各对应的填充阶段，然后，以类似的方式对压实阶段过程变量求解，以便对于至少该部分求解域提供各对应的压实阶段的解，并确定填充和压实阶段的解是否可接受。填充和压实阶段过程变量可以包括：密度、流度、模具内腔填充时间、模具内腔压实时间、压力、切变速率、切变应力、温度、速度、粘度和体积收缩。在填充和压实阶段的解不可接受的情况下，可以修改边界条件和/或离散求解域，并重复进行分析，直到得到可接受的结果。为了便于用户确定该结果是否可接受，可以图形显示填充和压实阶段的各种解，例如填充时间、压实时间、密度、压力、切变速率、切变应力、温度、速度、粘度和体积收缩。

本发明的一个实施例涉及模拟将向一限定三维内腔的模具的流体注入的方法，该方法包括步骤：提供一个三维固体的计算机模型，根据该固体模型离散一求解域，规定边界条件，利用利用质量、动量和能量守衡方程对填充阶段过程变量求解，以便提供对于至少某些部分求解域的填充阶段的解，并确定填充模具内腔期间流体的注入的解是否可接受。

离散步骤可以包括如下子步骤：根据固体模型，通过将该模型再分为由多个节点限定的多个相连的元产生有限元网格。网格产生子步骤包括：在该模型中的各厚和薄区域内产生一各向异性的网格，这样网格精细化沿厚度方向提供了增高的分辨率，而不明显增加沿纵向的网格精细化。

边界条件可以包括这样一些参数如：流体成分、流体注入位置、流体注入温度、流体注入压力、流体注入体积流量、模具温度、内腔尺寸、内腔构形和模具分割平面，以及它们的变化。

对于填充阶段过程变量求解可以包括如下子步骤：在其中取决于温度的至少某些求解域，对于流度、压力、速度和粘度求解。温度又可取决于对流、传导和/或粘滞耗散传热的所起作用。可以迭代地计算速度和/或粘度，直到压力收敛。

这一方法可以包括如下子步骤：根据按指定时间增量的速度，确定在内腔中的流体自由表面渐屈(evolution)，其中迭代确定表面渐屈，直到充满内腔。

一旦已经模拟了填充内腔的仿真，该方法就可以利用质量、动量和能量守衡方程在至少部分求解域，在填充终止时根据过程变量的各对应状态对于压实阶段过程变量求解，以便提供在至少部分求解域的压实阶段各对应的解。然后可以确定对于模具内腔压实期间的流体注入，该压实阶段的解是否可接受。

此外，对于压实阶段过程变量的求解可以包括：在至少部分求解域，对于流度、压力、速度和粘度的求解，其中粘度取决于温度。温度又可取决于对流、传导和/或粘滞耗散传热的所起作用。可以迭代地计算速度和/或粘度，直到压力收敛。还可以确定按照边界条件产生的部件的质量特性。该质量特性可以包括：部件密度、体积收缩、部件质量和部件体积，并可以与速度和粘度一起迭代地计算，直到完成预定的压力为分布曲线。

附图简要说明

描述据信为本发明的特征的新颖特征并在所提出的权利要求中作为区别特征。根据优选和示范性的的实施例，在结合附图的如下的详细说明中更具体地介绍本发明连同其优点。其中：

图1示意表示根据本发明的一个实施例的、适用于所公开的模拟向模具内腔的流体注入的方法的计算机硬件装置；

图2示意表示根据本发明的一个实施例的、概括某些操作步骤的上层系统的流程；

图3A和3B示意表示在模具内腔中的一内部点以及沿一条将内部点和模具腔壁结合的直线相关的温度分布；

图4示意表示对于在模具内腔中的每一点i的温度分布的一维解析函数；

图5是根据本发明的一个实施例的按全局(global)时间增量在总体温度求解中所包括的步骤的流程图形式的示意图；

图6示意表示在一节点的速度和其与周围的元的关系；

图7表示在一元的表面上的速度的示意图；

图8示意表示在一节点的速度和其与一相连的上游四面体的元的关系；

图9示意表示在一元的表面上的速度和其与一相连的上游四面体的元的关系；

图10A和10B示意表示在模具内腔中的一个内部点以及利用在上游侧元节点处的一维解析函数的上游侧温度内插；

图11示意表示四面体的元表面的分割组合；

图12示意表示一用于在一个边缘分割四面体的元表面的模板；

图13A示意表示用于在两个相邻的边缘分割四面体的元表面的模板；

图13B示意表示用于分割四棱锥的模板；

图13C示意表示用于在两个相对的边缘分割四面体的元表面的模板；

图14A示意表示用于在三个边缘按一公用共享表面分割四面体的元的模板；

图14B示意表示用于在三个边缘按一公用节点分割四面体的元模板；

图14C示意示意表示用于在三棱锥的模板；

图14D示意表示用于在接连的三个边缘分割四面体的元的模板；

图15A示意表示用于在四个相对的边缘分割四面体的元的模板；

图15B示意表示用于在四个相邻的边缘分割四面体的元的模板；

图16示意表示用于在五个边缘分割四面体的元的模板；

图17示意表示用于在六个边缘分割四面体的元的模板；

图18A和18B示意表示节点层编号和在确定到内腔腔壁的距离时使用的连通性树(tree)，以限制可能的搜索空间；

图19示意表示用于利用固化温度和对于在每个液体节点处的温度分布的一维解析函数确定元的内部的固体/液体界面的位置的方法；

图20示意表示包含固体/液体界面的一个元内部的粘度变化；

图21示意表示在固体/液体界面附近的流体流动；

图22示意表示用于确定一个与每一腔壁节点和每一内部节点相关联的核心节点的方法；

图23A和23B示意表示腔壁、内部和核心节点之间的关系，和用于确定冷凝点处的压力的方法；

图24示意表示根据本发明的一个实施例的填充阶段的流程图，其中概括某些操作步骤；

图25示意表示根据本发明的一个实施例的压实阶段的流程图，其中概括某些操作步骤；

本发明的详细介绍

可以分别通过质量、动量和能量守衡方程将在三维仿真中的注入模制过程的模拟描述如下：

$\frac{\∂\mathrm{\ρ}}{\∂t}+\▿\·\left(\mathrm{\ρ}\stackrel{\→}{v}\right)=0---\left(3\right)$

其中ρ是密度，v是速度，t是时间，代表关于一位置矢量的梯度，→代表矢量的量值。

$\mathrm{\ρ}(\frac{\∂\stackrel{\→}{v}}{\∂t}+\stackrel{\→}{v}\·\▿\stackrel{\→}{v})+\▿P=\mathrm{\ρ}\stackrel{\→}{g}---\left(4\right)$

其中是P是动量通量，g是重力。

$\mathrm{\ρ}(\frac{\∂\hat{U}}{\∂t}+\stackrel{\→}{v}\·\▿\hat{U})+P:D+\▿\stackrel{\→}{q}=\mathrm{\ρQ}---\left(5\right)$

其中是U是比内能，D是变形速率张量，q是热通量，Q是表示其它热源或热沉的比热速率。

公式4一般按照温度和压力表达，而对于其中Q＝0，即没有附加热源或热沉的(例如热塑性材料)情况下，取如下形式：

$\mathrm{\ρ}{C}_P(\frac{\∂T}{\∂t}+\stackrel{\→}{v}\·\▿T)+\frac{T}{\mathrm{\ρ}}\frac{\∂p}{\∂T}(\frac{\∂p}{\∂t}+\stackrel{\→}{v}\·\▿p)=\mathrm{\η}{\stackrel{\·}{\mathrm{\γ}}}^2+\▿\·(\mathrm{\κ}\▿T)---\left(6\right)$

其中，C_p是比热，T是温度，p是压力，η是粘度，

是切变速率，κ是热传导率。这一κ应不受流导干扰。

通常这些方程是难分析处理的，从而利用离散数字程序将它们求解，例如利用边界元法、有限元法、有限差法、有限体积法和无网格法。根据有限元法，求解需要将求解域划分为一组更小的子域，并将各方程变换为描述更小的子域的集合的离散方程。将有限元法应用于三维域导致形成称为三维子域的元。根据一种方法，可以使用四个节点线性四面体的元。其中利用每个四面体的元的顶点处的节点限定各元，木语线性是指元内插的多项式阶或形状函数。

最通用的有限元实现导致域的离散，这样计算在节点的各点处的所关注的区域变量，和在一个元内的分析所需的材料特性是恒定的。为了提高精度，特别是在高梯度区，已经确定希望将各方程公式化，使得这样限定在各节点的区域变量(例如压力、温度和速度)和材料特性(例如粘度)，除非另有说明，根据该元形状函数内插在一个元内。已经确定按照这种方法的离散化，不仅保证附加的计算精度，而且还便于使用更高阶的先进的内插，下面将更完整地进行介绍。

归根结蒂，物理过程仿真的精度明显地取决于仿真中使用的模型怎样能够充分地描述作为该过程的特征的重要的现象。对于注入模制，有很多反映过程的复杂本质的具体的现象，为了精确地预测仿真结果，需要将它们正确地模拟。这些现象是：接近固体/液体界面的很大的温度变化梯度；固体/液体界面边界位置随时间和模具内腔中的位置的变化；在模制期间所有位置和所有时间的聚合物材料特性的大的变化；以及在模制期间随着注入的液体材料的固化造成的固相材料体积随时间的增加。

可以用很多方法捕捉这些特定的现象。一种方法是精细化网格。这等于要求充分较少元的尺寸，以保证在邻近元内部的线性内插，可以充分捕捉区域变量的快速空间变化。然而，这种方法的缺点是，元和节点的数目明显增加，这经常使计算费用过高。特别是它应用于如下这样的注入模制的仿真时尤其如此，即在其中部分尺寸例如厚度和长度可以差两个或更多个数量级，以及在其中区域变量的最大的空间梯度是横越最小的空间尺寸即厚度求出的。

一种替换方法是使用按最高阶形状函数例如多项式函数的元。更高阶形状函数能够以预期所需的元较少的方式，在一个元内对区域变量进行非线性内插。然而用于限定在这些元内的更高阶形状函数所需的节点的数目增加了，因而，该方程组仍然计算费用过高。

因此，已经确定，可以单独或组合采用某些模拟函数和假设，以便明显地提高仿真精确度，而不明显增加节点的数目，不会伴随对计算时间和资源的不利的影响。

迄今，如图1中所示，注入模制的仿真引起明显的计算工作量是适用于所公开的根据本发明的用于模拟向模具内腔的流体注入的方法的计算机硬件装置10。装置10可以是具有所需计算速度和容量的便携式计算机、小型计算机或其它适合的计算机，以便支持下面将更充分讨论的功能。计算机10通常包括一个或多个中央处理单元12，用于执行包含在有关模型仿真的软件代码中的指令。提供存储装置14例如随机存取存储器和只读存储器，用于临时性或永久性地保存模型代码以及计算机10所需的其它操作软件。永久性的非易失性读/写存储器例如硬盘，通常用于在其使用和闲置时间期间存储代码，以及存储由软件产生的数据。计算机10还可包括一个或多个输入装置16例如键盘和读盘器，用于接收输入例如来自用户的数据和信息，还可包括一个或多个输出装置16例如监视器和打印机，用于以图形和其它形式提供仿真结果。另外，提供通信总线和I/O端口，以便链接所有的组成部分，并按照需要能与其它计算机和计算机网络通信。

图2高度简化地示意表示根据本发明的一个实施例的、概括某些操作步骤的上层系统的流程。在第一步骤20，如上所述，产生和提供模具内腔的三维固体计算机辅助设计模型。然后，在步骤30，利用各种各样的方法中任意一种定义该模型的求解域并将其离散，例如通过有限元分析，其中通过根据固体模型产生一有限元网格来产生有限元模型。该网格由多个接连的由共享的节点限定的固体元组成。为了进行分析，在步骤40，用户利用所形成的有限元模型或所限定的其它离散的求解域来规定边界条件。这些边界条件例如可以包括：几何限制条件和与流体成分相关的数值、流体注入位置、流体注入温度、流体注入压力、流体注入体积流量、模具温度、内腔尺寸、模具分开平面以及其它初始设计数据。此外，这些条件可以不是固定值而是与时间有关或具有分布为分布曲线。例如，注入压力可以随时间变化，或模具温度可以在不同的区域是不同的，这是由于在先冷却分析或实际冷却的结果。

一旦已经输入边界条件，在步骤50，计算机10执行根据仿真模型的指令，首先计算或求解与该节点有关的填充阶段过程变量。正如下面将对于在图24中表示的填充阶段流程图更详细讨论的一样，这些变量可以包括：流度、模具内腔填充时间、压力、切变速率、应力、速度、粘度和温度。此外，可以计算的并不限于这些变量，然而，有一些基本变量可以用于求解在计算这些例如结晶动力和纤维取向分布课题时所包括其它变量。此外，可以求解作为可压缩流体的填充，在这种情况下，在压实阶段计算中包括的很多质量项目(例如密度、质量和体积收缩率)也可以在填充阶段中计算。根据一个实施例，可以根据流体在填充阶段是不可压缩的而在压实阶段是可压缩的假设进行仿真。根据另一个实施例，可以根据流体在填充阶段和在压实阶段都是可压缩的假设。此外，参照图24和25正如下面将更详细讨论的，为了在压力、速度和粘度之前求解流度并不是必须的，没有一个对求解流度是完全必须的。这通常仅是当将一种例如Nakano所讨论的方法用于求解质量和动量守衡方程的方法时所需要的。也可以利用例如Stokes和Navier-Stokes方程求解质量和动量守衡。因此，这里所提出的原理可应用于例如Nakano、Stokes和Navier-Stokes的物理过程的所有说明。

一旦仿真在分析中达到该其中确定模具内腔已充满的状态，在步骤60，计算机10执行根据仿真模型的指令，以便接着计算或求解与该节点有关的压实阶段过程变量。正如下面对于参照图25中所表示的压实阶段流程图将更详细讨论的，这些变量除了流度、压实时间、压力、切变速率、应力、速度、粘度和温度以外，还可以包括：根据仿真模型产生的部件的质量特性例如密度和体积收缩率。

当完成分析时，在步骤70，可以将分析结果以各种方式输出。例如，为了使用户能通过视觉检查，可以将相关变量以图形形式显示变暗的固体模型，或者可以电方式输出以便进一步处理或分析。在步骤80，如果填充阶段和压实阶段的结果看来是可接受的，则在步骤100仿真终止，用户可以继续将该设计交付制造。由于规定的边界条件包括与注入模具的构形和过程参数相关的信息，用户将该设计交付注入模具的制造和直接产生的注入模制过程的图纸。

然而，如果用户确定在步骤70的仿真结果是不可接受的或不是优选的，则在步骤90用户任选修改一个或多个边界条件和/或将仿真求解域离散，之后迭代地重复仿真步骤50到70，直到用户满意该结果之时为止。不可接受的结果的示例包括：模型分析的不稳定性或操作故障例如造型空洞(shortshot)，其中模具内腔未完全充满、或者在填充期间产生过高的温度、温度、速度和压力，它们可能降低部件聚合物材料的特性，或者在部件中引起过大的残余应力，这将对产品生产率产生不利的影响，并可能导致部件永久性的缺陷。通过在设计过程之前提供这种高精度分析仿真能力，可以避免在初始生产运行期间的费用过高和下游侧延迟。

在了解图24和25中的填充阶段和在压实阶段的流程图的操作步骤之前，最好研究一下基础理论和其中采用的某些方法的基础。

如上所述，常规的三维注入模制的模拟技术的某些局限性与在固化期间很大的温度梯度和固体/液体界面边界移动、很大材料特性变化以及在模制期间熔化材料的固化有关。已经确定粗略分组的五种改进的方法与能量守衡原理的应用有关，其在预测精确度和/或计算时间方面有明显的提高。首先使用一维解析函数来描述在一节点的局部温度分布。接下来限定一维解析函数随时间的变化，以计入热对流。第三涉及与描述一维解析函数的变化描述，以计入粘滞发热。再接下来利用一维解析函数对显式(explicit)温度对流模式进行描述。第五提出有限元网格精细化包括一到内腔腔壁的距离的计算算法。虽然公认并将能量守衡原理应用于注入模制仿真的模拟，这里被分为5类的一般区，但这些并不构成限制，相反，这种方式表现仅是为了更清晰。

首先了解使用一维解析函数来在模具填充和压实期间的任一时间点描述在一节点的局部温度分布，存在于在内腔中的材料之内的指定内部点和模具腔壁的一部分之间的温度分布，该部分最接近该内部点。通常，可以取沿一条直线的温度分布，该条直线为该内部点和模具腔壁之间的最短距离。例如参见图3A和3B，它们示意表示在流体中沿一条将模具腔壁的内部点和模具腔壁结合的直线的温度分布曲线，其作为到模具腔壁的距离d的函数。

可以将各种解析函数或离散函数用于描述温度分布。根据一个实施例，可以使用一种误差函数公式，因为其物理意义是在半无限固体中热传导的解。根据本发明，任何一个和每一个内部点像半无限固体的域的一部分一样进行处理。每一个点与其它每一点无关处理。使用改进的关于在半无限固体中的热传导的方程，以便提供内腔中的所有位置上的点的、在注入模制期间所有时间及局部温度分布的更准确的估值，以及考虑所有有意义的传热机理。

在半无限固体中的温度分布随时间和距有限边界的距离而变化。在本应用中，可以将有限边界限定为模具内腔腔壁。当在半无限固体中的有限元处的温度T_∞(起初为均匀温度)瞬间变为另一温度T_w时，在固体中的温度分布随距腔壁的距离及时间而根据如下函数变化：

$T=T_w+(T_{\∞}-T_w)\mathrm{erf}(d/2\sqrt{\mathrm{at}})---\left(7\right)$

其中：

erf为误差函数；

d为距腔壁的距离；

t为时间；

T_∞为距腔壁的距离无限大之处的温度；

T_w为腔壁温度；

α为固体的热扩散率。

例如参见由R.B.Bird，W.E..Stewart，E.N.Lightfoot所著由JohnWiley&Sons1960年出版的 Transport Penomena。

已经发现这一函数在应用于注入模制的仿真时具有几个有用的特性。例如，其可以捕捉在短和长距离方面和时间的温度的空间和时间变化。其形式简单并易于计算。有几个可调整的参数。其具有物理含义，因为热传导是主要的传热机理，并且在短的接触时间，该方程紧密近似对于有限块(slab)的解，即对于一个在两侧由腔壁界定的固体，与在模具内腔中的情况相似。

这个函数可以便于地用于提高的温度的精确度，因为在短的时间区间内，在注入模制期间，温度分布在相对大的距离内是高度非线性的。由于通常希望保持相对低数量的节点和元，以便降低计算费用，与要模拟的局部的内腔几何尺寸相比较，元的尺寸相对大。此外，这些元具有线性形状的函数。因而，它们不能适当地描述在温度变化为高度非线性的情况下的温度变化。这一误差函数公式可以描述非线性的温度分布。

在应用中，可以对于在三维域中的每一个和任何一个内部点将这一函数定义。分别对于每一个节点将用于每一个节点的方程参数定义。例如参见图4，图4示意表示对于在模具内腔中的每一点i的温度分布的一维解析函数。图中表示了3个点。

下一个考虑了热对流的改进的方法与一维解析函数随时间的的变化的定义相关。在某些情况下，使用在以上提出的误差函数方程可能不足以将注入模制期间的流体温度的特有变化适当地模拟。在严格意义上，误差函数方程打算应用于流体，即范围为半无限和静态。然而，在实际注入模制应用中，范围为无限的，以内腔腔壁为界限，对于很多过程，流体处于移动中。特别是流体移动引入另一种显著的传热机理，即对流，其分布应当包括在总的能量守衡分析中，以提高仿真预测精确度。

考虑注入模制的物理过程。材料接触能进行热传导，在没有其它传热机理的情况下，利用误差函数方程可适当地描述在短时间内的热传导。然而，对流传热叠加在这种机理上。对流传热具有由于阻碍由传导引起的热损耗的作用，因为在填充阶段期间，连续地注入模制将来自上游侧的材料带到所考虑的任何一点。

应注意，误差函数方程包括关于时间的一项。因此，已经确定，一种方法建议要考虑热输入，对流传热对延迟对于在任何一点所考虑的传导通过时间有影响。这样，当全局实时时间按恒定速率继续时，所描述的从填充到压实的注入通过时间，可以考虑每个节点，以便带有其本身的节点热时钟，该时钟局部变化并且通常可能不像全局时间一样按相同速率推进。

例如，如果在内腔中的一内部点处的流体是静止的，则只要这是真实的，节点热时钟像全局时间一样按相同速率标记时间。然而，如果在内腔中的一内部点处的流体是十分快地移动的，使得对流传热占优势，则由于传导引起的任何局部热损耗整体上由于对流由热输入所抵消。在这种情况下，该节点热时钟是静止的，完全不前进，即便当全局时钟已运行也是如此。

在流动居中以使得传导和对流都不占优势的情况下，对于任一节点i的热时钟可以根据如下方程递增：

$t_i^k=t_i^{k-1}+\exp (-F_{\mathrm{pe}}{\mathrm{Pe}}_i)\mathrm{\Δt}---\left(8\right)$

其中：

t_i ^k是在时间步k在关于节点i的误差函数方程中使用的节点的热时钟；

t_i ^k-1是在时间步k-1在关于节点i的误差函数方程中使用的节点的热时钟；

Pe_i是对于节点i的Peclet数；

F_pe是调谐参数；

Δt是全局时钟时间增量。以及其中还有：

Pe_i＝d_i·v_i/α_i    (9)

其中：

d_i是距节点i处的腔壁的距离；

v_i是节点i的速度；

α_i是节点i的热扩散率。

多个方程中的任一个可以用于描述节点热时钟时间推进，选择这一经验方程是因为其具有所希望的特性。特别是：

当v→0时，exp(-F_pePe)Δt →Δt；

当v→∞时，exp(-F_pePe)Δt→0；

当α→0时，exp(-F_pePe)Δt→0；以及

当α→0时，exp(-F_pePe)Δt→Δt。

利用对于经验数据的调谐仿真结果，求出调谐参数常数F_pe。该数据可以是根据注入模制实验的显式的温度分布数据，或者间接数据例如注入压力测量值。已经确定通常可以将常数F_pe设置在约0到1的范围内。

再一种改进的方法建议了一维解析函数的变化以计入粘滞发热。由于流体的运动需要机械功，热动力原理还需要通过这种作功将热量添加到系统中。在流体运动的情况下，这种作功采取粘滞耗散或摩擦热的形式。在注入模制期间遍及模具内腔的内部的摩擦热随时间变化。在任何时间区间，这种热能的所起作用由在该时间区间内的局部切变应力和局部切变速率的乘积给定。已经确定，这种热能所起的作用往往会在固体/液体边界附近最重要，在其中速度和应力两者都高，速度高是由于接近冷凝的材料的相对低的温度，应力高是由于这种粘性凭体的运动。

可以将摩擦热认为是一种在内腔中的每一点的局部的随时间变化的热源。这种局部添加的热，在其标准形式下，不能由以上参照的误差函数方程说明。对这一点，在由T_∞定义的最大值以上，没有任何机制增加流体中的温度；而在注入试验中，已经确定固体/液体边界附近的流体温度不仅可能超过T_∞，而且还可以超过注入模具内腔之点的流体温度。

因此，对于精确的模型，在仿真中的内腔中的流体温度，在分析中应当考虑摩擦热。当这样做时，即使在对流和传导之后，也可能发生局部温度升高在T_∞之上。通过在每一时间步调节远区域的温度T_∞可以说明这种温度漂移。通过在该仿真中替换在每一时间步终点的温度T_∞可以实现调节，例如通过将在能量方程求解之后所计算的节点温度代入误差函数方程和按照如下方程得出在每一节点i处的T_∞的数值来实现：

$T_{{\∞}_i}=T_{w_i}+\frac{(T_i-T_{W_i})}{\mathrm{erf}(d_i/2\sqrt{{\mathrm{\α}}_i{t}_i})}---\left(10\right)$

因此，对于任何一个内部节点，有一机理来描述局部的一维温度分布和其随空间和时间的变化。这就使得比简单的线性内插利用更好的近似，甚至由元形状函数负担的更高阶的多项式内插。该新的函数可以用于显著地提高对流传热计算的精确度，和计算注入模制过程期间内在任何时间的固体/液体界面的位置。

图5以流程图示意表示根据本发明的一个实施例的按全局时间增量在总体温度求解中所包括的步骤，下面将更充分地讨论其中的细节。在第一步骤15，计算在该时间增量内的对流传热的的所起作用。然后，在步骤25，计算在粘滞耗散传热的的所起作用。之后，在步骤35，考虑在该时间增量内的热传导，求解总能量方程。然后更新对于每个节点的温度T_∞，其结果可用于求解动量方程。确定在填充阶段期间的流动前部的渐屈，以及量化材料特性变化，下面将更详细地进行讨论。

已经由本发明提出的其它区域涉及利用一维解析函数对显式温度对流模式的描述。对用于注入模制的能量守衡方程的求解，部分地需要计算对流对于总传热过程的的所起作用。同时根据本发明的一个实施例，已经预期一些机理，以便处理对流传热，采用一个显式模式，其中单独计算对于每个节点的温度变化的所起作用。

由于显式模式能够找到一在一指定时间增量内将流入目标节点的上游侧点，及在该上游侧点处的温度的插值。可以将先前介绍一维解析函数用于提供在该上游侧点而不是其它点处的更精确的内插值。采用拉格朗日内插来确定上游侧的各项，换句话说，称为平流项。可以利用相同的机理来求解该能量方程求解所需的温度的平流(对流)和求解该在填充阶段期间模具内腔内的自由表面和流前位置所需的物质的平流(浓度)。

用于能量和浓度函数方程的Peclet数分别是十分大和无穷大。用于仿真这种特征的常规的有限元法或者是不稳定的或者是往往会包含破坏精确度的程度不能允许的扩散。这是因为常规的最小平方计算法处理该问题是非通常的方法。已经确定处理混合不同方程的通常(nutural)方法连同其特性都遵循该函数。在注入模制范围中的流体流动的情况下，这解释为跟踪在连续介质中的任意“质点”的路径。

该算法的一个目的是对指定的三维域中的任意指定点利用速度场在最后的时间步内跟踪该现在占据节点的质点的原点。换句话说，该概念是查看上游侧，以了解现占据该节点的质点在一个时间增量之前在何处。一旦原点被确定，可以进行内插，以确定在该上游侧点处的浓度函数或温度的数值。然后，对于这一指定的数值，可以在当前的时间增量精确地更新所考虑的节点处的数值。

下面参照图6，从一节点开始，以节点A表示速度矢量，不失普遍性，在这些图中的元表示为三角形，而不是四面体形，以便于说明，然而，如上所述，该元还可以包括六面体形，或其它适宜的元的形状，这取决于所采用的方法和模型的构形。

粗略查看该图可看出当前占据节点A的质点来自该具有相对面L的三角形。然而，一种算法通过测试以了解哪一个三角形/四面体形满足所有所需数据来确认这一点。测试应当通过计算机快速实现。在下面将讨论在子程序upwelm1.F中执行这部分计算，以确定由环绕的三角形/四面体形中间的哪一个是正确的。

一旦记录了正确的三角形/四面体形，可以认为该质点来自表面L以外。通过将当前的时间步长度与对于一质点从表面L开始横越三角形/四面体形所需的时间比较，可以证实这一点。如果需要搜索以便延伸超出L进入一接连的三角形/四面体形，由也在下面讨论的子程序contig.F确定该正确的邻近三角形/四面体形并相应地重排列该局部的形状函数值。在一个三角形实例中，假设按L排列形状函数是〔0.2，0.8，0.0〕。则如果在接连的元中的关于L的节点号数的次序相反，形状函数排列变为〔0.8，0.2，0.0〕，假设第三节点和A两者是在元中的最后的节点，因此在表面L处的的数值为0.0。

在图7中示意地表示下一个出口(issue)。计算关于L的相关位置处的速度，并用于确定具有该速度的质点可能来源于哪一个表面。在下面讨论的子程序upwelm2.F中对此进行计算。应注意，upwelm1.F和upnem2.F之间的主要差别在于，在前者中，搜索从相对于该相连的元的上游侧元的节点开始，而在后者中，搜索从相对于同一元的内部的上游侧表面的一个表面开始。此外，应注意，如果质点来自该相对面之外，则在重新排列在congtig.F中编号之后，可在upwelm2.F中完成进一步的搜索。

upwelm1.F和upnem2.F两个子程序从全局和局部坐标之间的基本关系开始：

x-x₄＝J·L    (11)

这些字母以黑体字显示，以代表矢量(x.L)或矩阵(J)，其中：x是任一全局位置矢量(坐标)；L是局部坐标矢量，J是两个坐标系统中转换的导数行列式(雅可比)(Jacobian)，x₄是在全局坐标系统中的该元的第四节点的位置。

参照图8，可以将upwelm1.F子程序描述如下。在节点A开始，希望求出与一个元的网格面(facet)的第一迎流向相交，或者如果时间步对达到该网格面是不够的，则希望求出在该元内的迎流向位置。如在图中所示的，A处的速度是已知的。

在第一种情况下，假设A是该元中的第四节点，在图中表示四面体中的节点次序。这样产生：

L＝J^-1·(x-x₄)    (12)

下一步骤是求出在何处该质点轨迹相交于由与节点A相对的网格面限定的(无限)平面。注意，除非速度平行于该平面，必然在某处有一相交。如果

x＝x₄-kv    (13)

则位置x在该质点的轨迹上。这里，k是x₄起的上游侧时间。将方程(13)代入方程(12)，这一位置的局部坐标是：

L＝J^-1·(-kv)    (14)

该轨迹相交于其中L₄＝0的平面，因此

$\underset{i=1}{\overset{3}{\mathrm{\Σ}}}{L}_i=1---\left(15\right)$

展开这一式产生：

$-k\underset{i=1}{\overset{3}{\mathrm{\Σ}}}\underset{j=1}{\overset{3}{\mathrm{\Σ}}}{J}_{\mathrm{ij}}^{-1}{v}_j=1---\left(16\right)$

使k成为该公式的左端项：

$k=\frac{-1}{\underset{i=1}{\overset{3}{\mathrm{\Σ}}}\underset{j=1}{\overset{3}{\mathrm{\Σ}}}{J}_{\mathrm{ij}}^{-1}{v}_j}---\left(17\right)$

对于这一要作为候选的网格面，k必须是正的(即节点A的网格面上游侧)和有限的。一个无限大数值意指该速度平行于该相对的表面L。如果k是正的和有限的，下一步骤是检查与网格面平面相交是否在网格面的边界之内。当且仅当网格面节点相关的L的三个分量是非负的，这点才真。然后可以通过将来自方程(11)的k代入方程(14)来计算这些分量。

在第二种情况下，假设A不是第四节点，而是某一节点m，其中m是1、2或3中的任一个。则轨迹与该相对A的网格面的平面的相交发生在L_m＝0之处。假设这是在某位置x，因而：

x＝x_m-kv    (18)

将方程(18)代入方程(12)，并考虑第m个节点，(该符号利用m以描述局部的和全局的节点号数)产生：

$L=\underset{j=1}{\overset{3}{\mathrm{\Σ}}}{J}_{\mathrm{mj}}^{-1}(x_m^j-k{v}_j-x_4^j)=0---\left(19\right)$

该关于x的上标指示坐标分量(例如，x²＝Y，和x³＝Z)。使k成为分式的左端项，产生：

$k=\frac{\underset{j=1}{\overset{3}{\mathrm{\Σ}}}{J}_{\mathrm{mj}}^{-1}(x_m^j-x_4^j)}{\underset{j=1}{\overset{3}{\mathrm{\Σ}}}{J}_{\mathrm{mj}}^{-1}{v}_j}---\left(20\right)$

对于这一需作为候选的网格面，k必须是正的(即节点A的网格面上游侧)和有限的。一个无限大的数值意指该速度平行于相对的表面L。如果k是正的和有限的，则下一步骤是检查与该网格面平面的相交是否在该网格面的边界内，如果仅与该网格面节点相关的所有三个分量的是非负的，则这点是成立的。然后，通过将根据方程(20)的k代入方程(14)可以计算这些分量。

当求出正确的网格面时，下一步骤是通过比较在相交点的k，按照时间步Δt，确定该迎流向的点是否仍在上游侧。如果需要进一步迎流向，调用子程序contig.f以求出与新相交的网格面邻近的元，并将逻辑opface设置为正确，准备调用子程序wpwelm2.F。然而，如果迎流向的点在当前的元内部，通过在方程(4)中利用Δt而不是k，可以确定局部坐标，在这一情况中听任opface是不成立的。

下面可以参照图9介绍upwelm2.F子程序。首先，通过在网格面L进行内插，计算在当前相交点x₁处的局部速度V。由于在从原始节点A到达L时已利用了某些时间，原始时间步Δt仅余下某些时间Δt_rem。有三个候选表面通过测试轨迹与它们的每个平面的相交来检查。通过对于当前的网格面中的每个节点即图9中的节点a，b，c设置L_i＝0，可以实现这一点。

如果：

x＝x₁-kv    (21)

则位置x在轨迹上。因此，在轨迹上的各点的局部坐标是：

L＝J^-1·(x₁-kv-x₄)    (22)

如下方程用于每个节点a，b，c。为了简化，假设如下的情况是针对仅测试相对节点a的网格面的。在第一种情况下，假设节点是在该元中的第四节点，在图F中表示四面体形的节点次序。设L₄＝0，以便求出网格面平面相交点，这意指方程(15)能提供：(利用方程(22)展开)

$\underset{i=1}{\overset{3}{\mathrm{\Σ}}}\underset{j=1}{\overset{3}{\mathrm{\Σ}}}{J}_{\mathrm{ij}}^{-1}(x_1^j-k{v}_j-x_4^j)=1---\left(23\right)$

对此式求解k：

$k=\frac{\underset{i=1}{\overset{3}{\mathrm{\Σ}}}\underset{j=1}{\overset{3}{\mathrm{\Σ}}}{J}_{\mathrm{ij}}^{-1}(x_1^j-x_4^j)-1}{\underset{i=1}{\overset{3}{\mathrm{\Σ}}}\underset{j=1}{\overset{3}{\mathrm{\Σ}}}{J}_{\mathrm{ij}}^{-1}{v}_j}---\left(24\right)$

在第二种情况下，假设a不是第四节点，而是某一节点m，其中m是1，2或3中的任一个，则在L_m＝0之处，产生轨迹与相对节点的网格面平面的交点。将其代入方程(22)产生：

$\underset{j=1}{\overset{3}{\mathrm{\Σ}}}{J}_{\mathrm{mj}}^{-1}(x_1^j-k{v}_j-x_4^j)=0---\left(25\right)$

对此式求解k，产生：

$k=\frac{\underset{j=1}{\overset{3}{\mathrm{\Σ}}}{J}_{\mathrm{mj}}^{-1}(x_1^j-x_4^j)}{\underset{j=1}{\overset{3}{\mathrm{\Σ}}}{J}_{\mathrm{mj}}^{-1}{v}_j}---\left(26\right)$

为了确定这一网格面是否有效，可以进行对于upwelm1.F上述的相同的检查。当求出正确的网格面时，下一步骤是确定迎流向的点是否仍在上游侧，通过按照Δt_rem比较在该交点处的k。如果需要进一步迎流向的，可以调用子程序congtig.f，以便求出该邻近新的相交网格面的元，并可以重复上述操作程序。在这种情况下，Δt_rem适当地递降。然而，如果该迎流向的点在当前的元内部，通过在方程(12)中利用Δt_rem而不是k可以确定局部坐标，并在其后退出upwelm2.F。

一旦按照上述方法已识别在指定时间增量内的将流动到目标节点的该上游侧点，需要确定该上游侧点的温度，一旦用于精确确定该上游侧点温度的方法采用一维分析温度函数，如在图10A和10B中所示。正如以前显示三角形的元以便易于形象化，不过该概念应用于三维元而不会失去普遍性。

通常，对于每个节点，确定将“流动”或对流到在所要求时间增量内的节点的对应上游侧点。其可以是在该节点构成为一部分的元中的一个节点，或可以是在分开某一距离的另一个元中的一点，上游侧点与上游侧节点可以也可以不一致。如果不一致，需要内插以估计该上游侧点的温度。根据一种方法，可以采用两个步骤。

在第一个步骤中，估计到上游侧点的内膛腔壁的距离。虽然，在这一实施例中，将腔壁用作为基准表面，但在另一个实施例中，可以利用局部其它数据，例如冷凝层的表面。在确定到腔壁的距离时，可以使用各种方法。由于通常上游侧点将不与节点一致，通过对到在该元中包围该上游侧点的所有节点的腔壁的距离进行线性内插，可以确定其到腔壁的距离，这种内插利用元形状函数。

在第二步骤中，使用在该元的每个节点的一维解析函数，以便通过代换到上游侧点的腔壁的距离，提供在上游侧点的温度的估计值。通过根据元形状函数对四个节点估计值的内插，提供上游侧点的最终温度。则这是按指定时间增量的温度值，其对流到所考虑的节点。

在另一个实施例中，将到腔壁的距离和描述上游例点温度分布的一维解析函数取作为在上游侧元中的最接近的节点。

虽然有需处理的特定情况，例如上游侧点是注入位置，但这些情况并不与基本概念相背离，或有损于使用一维解析函数来估计温度。

再一种改进的方法涉及各相异性的网格精细化，其包括到模具内腔的距离的计算算法。各向异性的网格精细化是一种强有力的措施，其增加模具内腔的分辨率，而不会因精细化不必要地增加计算机的负担。精细化所体现的过程细分，已有的四面体的元为更小的四面体。然而，该过程仅适用于呈现函数不足之处。如果呈现充足的层数，找准变量求解的精确度，则不必进行改造。根据本发明的一个实施例执行该算法如下。

已有的四面体网格用层号标注。首先，标注在内腔腔壁上的所有节点为0。然后，将与层0相邻的所有节点标注为1。接着，将与层1相邻的所有节点标注为2。按照顺序递增的层号重复进行，直到将所有的节点指定层号数。其目的是确定多少层存在于模具内腔的每个区域中的以及定义仅在需要的那些区域中的所需的层号数。

在层数不足的每个区域中，将四面体的元分割。分割过程导致元的形状比(aspect ratio)大于原有的四面体。因而，精细化的网格变为各向异性在于该形状比明显偏离单位值。由于分割过程发生在接连的各层之间，精细化朝沿厚度方向增加元和节点密度进行偏置。由于正如前面介绍的，为了增加仿真精度贯穿该部分厚度区域变量变化十分明显，这是精确的方法。图11到174表示可以将四面体分割的各种可能的方法。

具体地说，图11示意表示四面体元分割的各种组合。在如下的图中所有的元分割模板关于在图11中所示的元的表面产生三种图形之一。可以精细化一、二或三个边缘。如果这点成立，则在精细化期间网格表面连通性将是不可能的。图12示意表示用在一个边缘分割四面体的模板。图13A示意表示用于在两个相邻的边缘分割四面体元的模板，其包含首先分割为两个元，其中之一是四棱锥，然后，继续将该四棱锥分割为一些四棱锥，根据图13B中所示的模板。作为当分割四棱锥时的一个通用规则，选择其底中的对角线，按包含该节点的对角线以四个底中节点的最低数目进行分割。图13c是一用于在两个相对的边缘上分割四面体元的模板，第一次分割产生两个四面体元，每个利用适用的模板需要在一个边缘上进行分割。

图14A示意表示用于在具有一共同的面的三个边缘上分割四面体元的模板，图14B示意表示在具有一共同的节点的三个边缘分割四面体元。后者形成一个四面体元和一个三棱锥元，它们可以根据图14C中所示的模板进行分割。应注意，节点A将是出现的六个节点中具有最低节点号数的节点。这就保证了将由第一步骤所分割的两个四边形按照包括该最低节点的对角线进行分割。此外，第一步骤生成一个四面体元和一个以四边形为底的棱锥，它们可以被继续分割。图14D是用于在接连的三个边缘上分割四面体元的模板，根据哪一个节点具有最低的节点号数，形成其中的两个可能的组合。

图15A是用于在四个相对的边缘分割四面体元的模板。第一对称分割结果是形成两个三棱锥，又可再进行分割，如图14C中所示。图15B是用于在四个相邻的边缘上分割四面体元的模板。该分割产生两个四面体元和两个四边形为底的棱锥，它们可以利用以上讨论的适当的模板和规则进行分割。

图16是一用于在五个边缘上分割四面体元的模板，其产生两个四面体元、一个三棱锥和一个以四边形为底的棱锥，其中的两个可以利用上面讨论的适当的模板进行分割。最后，图17是一用于在六个边缘上分割四面体元的模板。切去四面体元的每个角，产生一具有八个三角形表面的中心主体，在一个四边形平面上将中心主体分割为二，产生两个四边形棱锥，根据上面讨论的模板可以进一步将其分割。

再次涉及分层，已经确定最佳的层数是6，最佳的形状比是5，可以采用在这一范围之外的数值。但是当层数降低时，最大形状比增加，精确度通常会降低。因而，当层数增加而形状比降低时，元的数目很快地增加，计算费用势必变得高费用的。

为了计算对于每个内部节点的到腔壁的距离，根据节点所在的层和为了帮助找出最接近的腔壁位置而构成的树，将节点编号。按照本发明的一个实施例过程算法执行如下。

首先，将在内腔腔壁上的所有节点标注为层0。然后将与该层0相邻的所有节点标注为层1，接着将与该层1相邻的所有节点标注为层2，按照顺序递增地层号数重复进行，直到所有的节点被指定层号数。在此以后，对于每个节点，构成一用于所有连接节点的节点关联数，其中这些节点具有一可下降到层号数0的逐渐降低的层号数，即腔壁节点。然后，突出目标节点，在该树中由层0的各节点形成的所有表面上将到腔壁的距离记录为在该组可能的距离中的最短的距离，这种方法具有的优点是将对于最短距离的可能的搜索空间变窄。在利用具有层号数0的节点形成的一个表面上，不能突出该节点的情况下，则将该距离取作到最近层0的节点的距离。

在图18A和图18B中示意地表示了这一算法。在图18A中，将节点层编号，并且该连续性树限制了用于突出目标节点到表面“A”和“B”的搜索空间，其是远小于在内腔中所有可能的表面的集合(set)。因此是高效的，在图18中所示的情况下，其中目标节点不能在其连接性树中的表面“A”或“B”上突出，将距离d取作为最接近的腔壁节点。正向以上所讨论的，不仅通过考虑和应用能量守恒原理，而且还通过建议以该形成的冷凝或固化层近似模具内腔腔壁，已经实现明显提高仿真预测精度。常规的三维注入模制的模拟技术的局限性经常与大的温度梯度、固体/液体界面边界的移动、在固化期间的大的材料特性变化以及在模制过程中熔化材料的固化有关。

粗略分组的以上五种改进的方法涉及对固化的材料的便于处理，这些方法已经确认为在预测精度和/或计算时间方面比常规的仿真方法有明显的改进。首先与确定固体/液体界面的位置有关。其次是遍及其材料特性线性变化的元的公式化。第三与确定在包含固体/液体界面的元内的有效粘度函数有关。第四涉及从求解域中消除冷凝点和元。第五是计算在已经固化的区域中的有效压力。

首先注意确定在一个元中的固体/液体界面位置，通过该元存在一固体/液体界面，就是说，在一个包含固体区域和液体区域的元内，材料特性通常按照横越界面边界的量值的几个数量级变化。特别是，在固体区域中的粘度(如果适合完全使用术语粘度)实际上是无限大，而在液体区域，粘度是有限的并且相对小。因此，在有限元方法中利用的简单常规元组合程序将导致对该元对流动的阻力的实际过高估计，因为对该元的刚度朝着高粘度值严重偏置。无论实际冷凝多少的元，这点都成立。最终结果是对在流动过程中的压力下降过高预测。

根据本发明的一个实施例，通过利用关于温度的一维解析函数，可以解决这一问题，以便更精确地估计固体/液体界面的位置，然后限定与在一个元中的冷凝材料数量相关的有效的节点粘度。通过检查温度分布曲线与聚合物材料的固化温度的相交，选择界面的位置。固化温度T_s限定了固体和液体状态之间的转变。

在图19中示意表示了该过程。为了易于描述，再次利用二维的元，但是该改进的概念可以应用于三维的元而不会失去普遍性。该图表示通过利用固化温度T_s和在每一个液体节点的关于温度分布的一维解析函数，未确定在一个元内部的固体/液体界面。

下一种改进的方法关于遍及每个元的具有材料特性线性变化的元的以公式表示(formulation)。很多共同的有限元实现导致域的离散，这样计算在该节点的各点处的区域变量，分析所需要的材料特性是一个元内部是不变的。根据本发明虽然这一技术可以在某些较少区域的应用中提供可接受的结论，但在将有效元方法应用于注入换压中的液体流动，各方程以式表示，这样区域变量(例如压力、温度、流导和速度)和材料特性(例如粘度)限定在各节点处的，并且除非另作说明，根据元件形状函数将它们内插在一个元内部。按这种方法的离散提供了计算的附加精度，并已经发现对于更高阶的先进的内插是特别有利的。

341例如，对于由Nakano给出的对于流导的方程(该方程利用材料粘度特性)可以重新改写如下：

$\frac{\∂}{\∂x_k}\left(\mathrm{\η}\frac{\∂k}{\∂x_k}\right)=-1---\left(27\right)$

其中下标k是分别标示每个坐标方向x，y和z的指数1、2和3。对该元进行Galerkin加权积分，给出：

${\∫}_{{\mathrm{\Ω}}_e}\frac{\∂}{\∂x_k}\left[{\mathrm{\η}}^e\frac{\∂k^e}{\∂x_k}\right]N_{i}d{\mathrm{\Ω}}_e=-{{\∫}_{\mathrm{\Ω}}}_e{N}_{i}d{\mathrm{\Ω}}_e---\left(28\right)$

其中该元流导和粘度分别由下式给出：

$k^e=\underset{j=1}{\overset{4}{\mathrm{\Σ}}}{N}_j{N}_j---\left(29\right)$

${\mathrm{\η}}^e=\underset{j=1}{\overset{4}{\mathrm{\Σ}}}{N}_j{N}_j---\left(30\right)$

其中Ω_e指该元的体积，Ni是在元i中的每一个节点的元形状函数，下标e是元的指数。应用收敛定理给出：

${\∫}_e{\mathrm{\η}}^e{N}_1\frac{\∂k^e}{\∂x_k}{n}_k{N}_{i}d{\mathrm{\Γ}}_e-{{\∫}_{\mathrm{\Ω}}}_e{\mathrm{\η}}^e\frac{\∂k^e}{\∂x_k}\frac{\∂N_i}{\∂x_k}d{\mathrm{\Ω}}_e=-{{\∫}_{\mathrm{\Ω}}}_e{N}_{i}d{\mathrm{\Ω}}_e---\left(31\right)$

其中n_k是在元边界Г处的正交矢量。

在元组合在内部元表面期间，由于由相邻各元的相等和相反的作用，假设第一项变为0。这一点成立，是由于基于节点值内插的元的粘度沿着两个相邻元之间的边界是一致的。此外，将流导的梯度假设为元边界的每一侧是相等的。仅在相邻的元中正交矢量具有相反的符号，其导致因来自相邻各元件的作用相反而变为0。在填充的区域内的内腔边界，当施加一已知的数值边界条件(k＝0)时，无需该边界积分。然而在需要腔壁滑动的情况下，该边界积分项也是需要的。在任何一种情况下，需要估计在自由表面边界和在注入边界处的边界积分项。

可以将证明在一个元内部的线性特性变化的方程(29)和(30)代入方程(31)，给出：

${-{{\∫}_{\mathrm{\Ω}}}_e\mathrm{\η}}_1{N}_1\frac{\∂N_j}{\∂x_k}{n}_k{N}_i{k}_{j}d{\mathrm{\Γ}}_e+{{\∫}_{\mathrm{\Ω}}}_e{\mathrm{\η}}_1{N}_1\frac{\∂N_j}{\∂x_k}\frac{\∂N_i}{\∂x_k}{k}_{j}d{\mathrm{\Ω}}_e={{\∫}_{\mathrm{\Ω}}}_e{N}_{1}d{\mathrm{\Ω}}_e---\left(32\right)$

已知利用自然坐标形状函数对于一个四面体的精确积分是：

${{{\∫}_{\mathrm{\Ω}}}_{e}N}_i^{\mathrm{\α}}{N}_j^{\mathrm{\β}}{N}_k^{\mathrm{\γ}}{N}_l^{\mathrm{\δ}}=\frac{\mathrm{\α}!\mathrm{\β}!\mathrm{\γ}!}{(\mathrm{\α}+\mathrm{\β}+\mathrm{\γ}+\mathrm{\δ}+3)!}---\left(33\right)$

其中：

N_i是形状函数，

V是元的体积，

α、β和γ是积分指数。

利用按照方程(32)的积分关系式，给出：

$\left\{\underset{f=1}{\overset{4}{\mathrm{\Σ}}}\left({\mathrm{\δ}}_{\mathrm{fi}}\right)\right[\frac{A\underset{i=1}{\overset{3}{\mathrm{\Σ}}}(1+{\mathrm{\δ}}_{i1}){\mathrm{\η}}_1}{12}\frac{\∂N_j}{\∂x_k}{n}_k]+[\frac{\mathrm{\Σ}{\mathrm{\η}}_{1}V}{4}\frac{\∂N_j}{\∂x_k}\frac{\∂N_i}{\∂x_k}\left]\right\}{k}_j=\frac{V}{4}---\left(34\right)$

其中仅在形成每个表面的三个节点的范围内，在边界积分的粘度的求和。为了完整，这里在每一表面f上表示该积分。该表面的编号与相对的顶点的号数相一致。

这样表明对于需要材料粘度特性的流导的方程，怎样可以根据元形状函数作为遍及每个元的随粘度变化的有限元系统以式表示。

改进的再一个方面涉及确定在包含固体/液体界面的元中的有效粘度函数。最一般的部分冷凝的元发生在模具腔壁附近，在该处它们包含“冷”的腔壁节点和“热”的腔内节点。根据通过各个元的粘度的线性变化的假设并由于利用按照线性形状函数的元而升高的粘度线性变化，通常会过高估计这些元的有效粘度，在很多情况下产生明显的误差，最一般和最严重的是部分地冷凝的元。例如参阅图20，该图将线性粘度分布曲线和更实际的高度非限性的分布曲线相比较。这一线性假设导致在仿真中的人为高的流动阻力，其反映可能比实际值多于100％的注入压力。因此该线性粘度变化的近似在部分冷凝的元中是十分差的。同时对于充分细小的网格该误差渐近到零，这一求解通常通过计算机实现是不实际的或有效的。因此，不按有效粘度函数技术实施，很大程度上可能影响仿真精度。

根据本发明的一个实施例，大为提高在部分冷凝的元中的粘度仿真精度涉及计算有效粘度。更具体地说，可以在元组合程序中计算和局部使用在包含固体/液体界面的元中的每个冷凝点处的有效粘度。正因这样，这些节点粘度是瞬态的，并且不替换利用在每个节点处的温度和切变速率的状态以及材料粘度函数的资料所计算的全局的数值。

通过分别将节点上方和下方表示预定的固化温度为“热”和“冷”，将部分冷凝的元定义为未完全“热”(所有节点指定为“热”)或未完全“冷”(即所有节点指定为“冷”)。冷的节点可以指定任意大的粘度，即粘度值充分大，保证在冷凝部分不发生流动。固体状态的粘度通常为10⁶兆帕的数量级。这种技术提供了在部分冷却的元中的更实际的有效粘度。完全冷凝的或完全热的元的粘度无需改动。

在对有效粘度法的一个实施例的说明之前，提供对其某些主要调整的概述是有所帮助的。可以对非等温流动求解的任何一种方法包括Nakano法以及更多的常规的Stoke和Navier Stokes解法进行粘度的改进。当能够进行切变加热时，可以通过相关的例行程序进行适当的粘度改进。

此外，在元组合期间进行改进。对于每一个元，改进四个局部的节点粘度。全局节点粘度是不改变的。在其节点全是冷的或全是热的元的中的粘度是不改进的。对于这些元仍然使用常规的或线性粘度变化的近似。

将在该元中的冷凝的体积百分率用于计算在每个节点体积内的冷凝的体积百分率。在这个特定的范围，“冷凝“是指低于T_s的体积部分，节点体积仅指在当前的元内的那部分。分别计算在该元中的冷和热节点的平均粘度vis_h和vis_c，每个节点体积可以包含冷凝和熔化部分，而与节点冷热无关。如果节点是冷的，则在冷凝的部分内的粘度是节点粘度，如果该节点是热的则为vis_c。与之相似，在熔化部分内的粘度是节点粘度，如果节点是热的；如果该节点是冷的则为vis_h。最后则是确定整个节点体积的有效粘度并局部用于该节点。根据流动平衡于局部冷凝的表面的假设，通过确定该节点体积中的冷凝和熔化的部分的平行流动阻力，可以导出关于有效粘度的表达式。

表1列出在这一方法中所使用的术语，以易于进行参考。

表1

变量	定义
		Ts	固化温度
nc	在低于Ts的元上的节点数目
		nh	在高于Ts的元上的节点数目
tfroz	在元中冷凝的体积的百分率(意味着低于Ts的冷凝)
		frzc	每个冷节点的节点体积中的冷凝的体积的百分率
frzh	每个热节点的节点体积中的冷凝的体积的百分率
		visi	元的各节点的粘度(i＝1到4)-全局数值
eteti	元的各节点的粘度(i＝1到4)-局部用于元的组合的数值

visc	冷节点的全局粘度nc的平均值
		vish	热节点的全局粘度nh的平均值

对于每个元，执行如下的步骤，以确定局部的节点粘度etet_i。下面介绍各方程的导出。根据该方法的这一实施例，首先确定nc和nh。如果nc＝0或nh＝0，这就不是一个部分冷凝的元。设etet_i＝vis_i并继续到下一个元。接着计算vis_c和vis_h。然后根据如下的方程计算frzc和frzh：

$\mathrm{frzc}=\mathrm{Min}[\frac{4\·\mathrm{tfroz}}{\mathrm{nc}},1]---\left(35\right)$

$\mathrm{frzh}=\frac{4\·\mathrm{tfroz}-\mathrm{nc}\·\mathrm{frzc}}{\mathrm{nh}}---\left(36\right)$

最后，检查四个节点；并计算每个的粘度etet_i。如果节点是冷的，利用方程(37)，如果是热的，利用方程(38)

${\mathrm{etet}}_i=\frac{{\mathrm{vis}}_i{\mathrm{vis}}_h}{(1-w){\mathrm{vis}}_i+w\·{\mathrm{vis}}_h},$ $w=\frac{{\mathrm{frzh}}^{1/3}}{{\mathrm{frzc}}^{1/3}+{(1-\mathrm{frzc})}^{1/3}}---\left(37\right)$

${\mathrm{etet}}_i=\frac{{\mathrm{vis}}_i{\mathrm{vis}}_c}{(1-w){\mathrm{vis}}_c+w{\·\mathrm{vis}}_i},$ $w=\frac{{\mathrm{frzh}}^{1/3}}{{\mathrm{frzh}}^{1/3}+{(1-\mathrm{frzh})}^{1/3}}---\left(38\right)$

按照下述导出量frzc和frzh。首先令当前的元的体积为V，则四个节点的体积每个是V/4。在该元中的冷凝的体积是V*tfroz.通过最终在各冷节点nc之间分布冷凝体积进行近似，在每个节点体积中具有冷凝的体积V*tfroz/nc。由于这一冷凝体积不可能大于节点体积，实际的冷凝体积是V*tfroz/nc和V/4中的最小值，即Min(V*tfroz/nc，V/4)。将其转换为体积百分率，得出方程(35)，即：

$\mathrm{frzc}=\mathrm{Min}[\frac{4\·\mathrm{tfroz}}{\mathrm{nc}},1]---\left(35\right)$

在分布到冷节点之后任何冷凝“剩下的”或其余的，最终分布到热节点nh。指定到冷节点冷凝体积是V×nc*frzc/4，因而剩余的冷凝体积是V*tfroz-V*nc*fizc/4。在各热节点nh之间对此进行分割，并转换为体积百分率，形成方程(36)，即

$\mathrm{frzh}=\frac{4\·\mathrm{tfroz}-\mathrm{nc}\·\mathrm{frzc}}{\mathrm{nh}}---\left(36\right)$

按照下述导出当地的节点粘度，令节点体积为V_t。假设该体积部分冷凝，将总体分割为冷体积V_c和热体积V_h，具有的对应的粘度为η_c和η_h，很明显，V_t＝V_c+V_h。此外，该局部的流动应平行于在节点体积中的局部冷却表面。图21示意地表示这种状况。虽然该图是二维的，但该可代表三维。

下一步骤是确定在流动方向上的整个节点体积V_t的有效粘度。对于该流动的有效粘度并不影响该流动求解。调用这一有效粘度η_t。η_t的严格计算对于在模型中的奇特形状的体积往往是不实际的，其需要复杂的体积积分。代之以可以利用简化的方案，当每个节点体积冷凝时，其提供正确的状态。

现在，具有均匀粘度η的某些体积的流阻R与ηL/A成比例，其中L是沿流动方向通过该体积的路径长度，A是与流动方向垂直的横截面积。由1/R_t＝1/R_c+1/R_h确定节点体积的流阻。对于每个体积计算L/A通常是不易的。相反，应注意，L/A的尺寸为1/长度，因而近似为1/A随体积变化。因此，应用如下的关系式：

$\frac{V_t^{1/3}}{{\mathrm{\η}}_t}=\frac{V_c^{1/3}}{{\mathrm{\η}}_c}+\frac{V_h^{1/3}}{{\mathrm{\η}}_h}---\left(39\right)$

由该公式可以得到：

${\mathrm{\η}}_t=\frac{{\mathrm{\η}}_c{\mathrm{\η}}_h{V}_t^{1/3}}{{\mathrm{\η}}_h{V}_c^{1/3}+{\mathrm{\η}}_c{V}_h^{1/3}}---\left(40\right)$

节点体积中的冷凝百分率为F，其中F是frzc或frzh，这取决于该节点是冷的还是热的。于是V_c＝FV_t，V_h＝(1-F)V_t。将它们代入方程(40)，得出：

${\mathrm{\η}}_t=\frac{{\mathrm{\η}}_c{\mathrm{\η}}_h}{{\mathrm{\η}}_h{F}^{1/3}+{\mathrm{\η}}_c{(1-F)}^{1/3}}---\left(41\right)$

可以将方程(41)归一化。例如，如果η_c＝η_h，则方程(41)可还原到η_t＝η_c。因此，将方程(41)改进到：

${\mathrm{\η}}_t=\frac{{\mathrm{\η}}_c{\mathrm{\η}}_h[F^{1/3}+{(1-F)]}^{1/3}}{{\mathrm{\η}}_h{F}^{1/3}+{\mathrm{\η}}_c{(1-F)}^{1/3}}---\left(42\right)$

最终的问题是在节点体积的两个部分中使用什么粘度。如果与节点体积相关的节点i是冷的，利用η_c＝vis_i并近似为η_h＝vis_h，与之相似，如果节点i是热的，利用η_h＝vis_i并且近似为η_c＝vis_c。最终的方程对于etet_i则是方程(37)或方程(38)，这取决于节点；是冷的还是热的。

方程(42)的特性在有限的情况下是正确的。例如，如果F＝0，则η_t＝η_h，按照要求，与之相似，如果F＝1，则按照要求，η_t＝η_c。此外，η_t随F的变化是平滑的并且是单调的。如果η_c＝η_h，则按照要求η_t＝η_c。最后，当η_c＞＞η_b时，η_t＝η_h[1+(F/1-F)^1/3]。这为正确的形式，因为η_t从在F＝0时的η_h的数值随F增加。

如果需要，可以利用两个较简单的表达式；然而，它们可能导致预测的状态的精确度下降。一个表达式如下：

η_t＝η_cF+η_h(1-F)    (43)

414按照方程(43)中的一个结果是在通常的情况下η_c＞＞η_h，该方程形成为η_t＝η_cF。实际上由于η_c非常大，这个方程预计为没有通过该节点的流动，即使当大部分的节点的体积被熔化时。相应地，预测的注入压力往往大于预期值。

第二表达式如下：

${\mathrm{\η}}_t=\frac{{\mathrm{\η}}_c{\mathrm{\η}}_h}{{\mathrm{\η}}_{h}F+{\mathrm{\η}}_c(1-F)}---\left(44\right)$

这一方程像方程(42)一样，在有限的情况下提供正确的状态。即，如果F＝0，按照要求η_t＝η_h。如果F＝1，按照要求，η_t＝η_c。η_t随F的变化是平滑的和单调的。如果η_c＝η_h，按照要求η_t＝η_c。最后，当η_c＞＞η_h时，η_t＝η_h[1+(F/1-F)]。当η_t随F从F＝0时的η_h的数值增加。

然而，有两个结果，利用方程(44)。首先，通过假设流阻正比于η/V，而不是像对于方程(42)的η/V^1/3，导出方程(44)。因此，该基本假设尺寸上是不正确的。其次，实际上已经证明方程(42)预测的注入压力接近正确的数值。对于一管状模型测试实例，预期的注入压力为0.62兆帕，利用方程(44)6层Darcy模型提供的压力为0.47兆帕，与之对比利用方程(42)提供的0.51兆帕。16层提供的压力分别为0.56兆帕和0.55兆帕。按照6和16层精细化，线性粘度公式提供的压力分别为1.59兆帕和0.63兆帕。

有效粘度方法的测试已经表明，tfroz的数值是相当精确的。此外，虽然，本方法采用对于“冷”节点的均匀分布，可替换的近似为应利用元节点温度对分布进行加权。另外一个是通过对通过每个节点体积的T_s表面进行内插和计算体积直接计算frzc。可以利用一种类似的技术用于将剩余的冷凝均匀分布到“热”节点，虽然这些替换方案可能耗用额外的计算资源。

虽然，对于冷和热区域粘度使用F^1/3和(1-F)^1/3的权重是基于经验的权重，提供附加的定量的支持，预期该权重可以基于在部分节点体积范围内的积分。另外，对于冷节点的热区域的平均热节点粘度的使用可能会过低估计热区域粘度，这是由于近似为粘度随冷凝表面而增加。这又可能导致轻微过低估计注入压力，正如上面所指出的。然而，已经确定对于热节点中的冷部分的平均冷节点粘度的使用产生可忽略的误差，这是由于仅需要vis_c充分大，以消除通过冷凝部分的流动。

现在，关于从求解域中消除冷凝的节点和元的技术，在填充和/或压实期间的某一阶段，根据在元中的所有节点的温度已下降低于固化温度T_s这一事实，使全部的元基本上是不可移动的。问题是怎样处理在流体流动分析的范围内的这些元。一种方案是通过将固体区域考虑为具有十分高粘度的区域，作为流体连续地处理这些元。通过对流体粘度进行外推，或通过指定一任意的大的数值，可以得到“固体”粘度。然而，这可能导致方程组病态。可能使收敛延迟或完全中止，并因此计算费用增加。

根据本发明的一个实施例，将对于流导和压力的求解域还原为仅包括具有固化温度以上的某些节点的那些元。通过从求解域中完全消去冷凝的元，易于实现对主要流体的收敛问题，同时还达到节省计算时间。将对于非滑动流动的常用的边界条件应用在冷凝区域的边界，而不是应用在模具内腔。在模制过程的早期阶段中，这两个边界是一致的。通过将在边界处的流导函数设为0应用这些边界条件。应注意，通常将该完全的内腔域用于温度场计算。

再一种改进方法涉及计算在已经固化的区域中的有效压力。虽已选择在模具内腔中的作为固体或液体的材料特性，通常还希望计算固体中的压力。实际上，冷凝的聚合物仍要经历压力波动，这是由于相距很小距离的邻近熔化聚合物，可能由固体聚合物传递冷凝压缩应力和小的偏转。为了向冷凝结定指定压力，通常从熔化的节点传播压力是不充分的，因为这种方法不能针对冷凝过程。其中所填充的内腔的某些部分可能仍然是熔化的并受到注入压力，而因为它们处于距熔化的区域为这样一个距离，没有实现有意义的压力传递。

此外，希望知道压力，以便确定材料的形态。利用压力以及温度资料确定材料密度以接连地确定在每个节点的质量，在内腔中的固体或液体。

此外，对于制造商通常合理地将压力传感器齐平安装在实际的模具内腔中，以便提供某些信息，利用这些信息可以控制模制机械和过程，以及试图确定与模制部件质量相关的参数。考虑该过程的实际物理质量，从而推断精确的仿真很可能需要十分复杂和费用高的利用计算机实现的多种物理量的仿真，以及将流体动力学和固体力学的结合。这时，然而，这些技术是不易产生的。

根据本发明的一种解决方案已经开发了用于通过设计对于内腔的外冷凝层的核心压力来估计固体压力的技术。更具体地说，在冷凝和熔化区域中，设计在冷凝的外侧节点上的核心压力。例如，可以根据通常的流体有限元法计算来确定核心压力，或者可以根据冷凝的固体材料的压力衰减，由于冷却该衰减随聚合物温度的降低而发展。通过按几何特性的分解的过程，将所有的未在部件核心的所有节点指定一核心节点，一旦冷凝它们的压力将是最小的。然后仅实际计算在冷却过程中其压力衰减的节点是给出核心属性的那些节点，可以对每个外侧节点指定核心节点，与其最接近的作为其核心节点。然而，用于确定哪些节点是核心节点的算法实际上取决于外侧节点。即，每个核心节点具有至少一个其所依赖的外侧节点；否则，核心节点必不给出节点的属性。

下面介绍用于产生核心节点的算法，首先将在内腔腔壁上的所有节点标作为层0。然后，将与层0相邻的所有节点标记为层1。接着将与层1相邻的所有节点标记为层2，并按顺序递增的层号数重复进行，直到将所有的节点指定层号数。此后，计算对于每个节点到最接近的内腔腔壁的距离。此外，记录哪一个是最接近对于每个内部节点的内腔腔壁节点。将其层号大于或等于该贯通厚度的要保证的层号的所有节点标记为潜在的核心节点。

该贯通厚度的要保证的层数是网格精细化处理的结果，如前所述。下面从在内腔腔壁上的每个节点即层号数0开始。越过依次与每个节点相连的节点接连地移动到更高的层号数，直到达到一潜在的核心节点之时为止。通常，对于每个内腔腔壁节点将有几个潜在的核心节点。该选作为与腔壁节点相关的核心节点的节点是具有最大比值d_cw/d_cn的节点，其中：d_cw是核心节点到内腔腔壁的距离，d_cn是核心节点到腔壁节点的距离，对核心节点，是在前一层编号步骤过程中作为最接近的内腔腔壁节点确定的节点。

这种算法其特征是当薄和厚的区域彼此相邻时，在薄区域中的冷凝点具有将来自厚层的核心压力指定到其上的。对于最小的实际物理过程已正确地发现这一点。例如，在将一狭窄的浇口连接到一较厚的转轮进给系统，在浇口已冷凝之后，其上的压力报告关于其模仿较厚的进给系统的压力。这就表示了通过固体聚合物的压应力的传递。

图22表示对于核心节点、腔壁节点和内部节点进行按层编号以及按该比值采用的距离的定义，以便确定对于每个腔壁节点指定哪一个核心节点。

图22A和23C表示核心节点、腔壁节点和内部节点之间的关系，以及它们随着时间的增加和冷凝的压力变化。

已经详细地介绍了各种精细化和算法，现在注意力转到图24和25，以便理解根据本发明怎样在填充和压实阶段实现改进。图24示意表示填充阶段分析，其概括了根图2中的步骤50中的某些过程和步骤。进行回顾，根据图2中的顶层流程图，一旦已经离散三维固体模型，以便提供模型求解域，并已设置边界条件，在步骤50中该仿真首先对填充阶段过程变量进行求解。

在步骤110中，步骤50中的填充阶段分析从假设内腔是空的开始。在步骤120中，将所有的区域变量例如压力、温度和速度进行初始化。此后，在步骤130将质量守恒和动量守恒方程用于求解流度k，在步骤140求解压力，在步骤150计算速度，对于一部分存解域，在使用有限元方法的情况下，确定在网格中的每个节点处的流度、压力和速度。

通过应用Daroy流动方程可以求解流度，然后利用Laplace方程求解压力，因为压力是流度的函数。另外，可以利用Navier-Stokes技术根据动量守恒求解压力。一旦已经解出流度和压力，可以直接计算速度。接着，在步骤160，根据能量守恒原理，利用温度计算粘度。然后，在步骤170，该模型检查以了解压力是否已经收敛，其很可能与初期情况不同。如为否，迭代地重复步骤130到160，直到压力收敛。

一旦压已经收敛，在自由表面渐屈步骤180中，该仿真以递增方式将流体流动波前向前进。接着在步骤190中考虑对流传热所起作用，其中实现显式(explicit)拉格朗日温度分析，以说明对在内腔中的熔化聚合物的连续流入的热能作用。然后，在步骤200中对在温度方程中的其余的项求解，以说明传导、粘滞耗散和任何其它预期的热作用，例如在热塑性材料中的固化热和在热固性材料中的反应热。

一旦所有的温度作用已经量化并结合到仿真中，在步骤210中，该算法进行检查以了解模具内腔是否已到内腔充满之时为止。通常该算法可以循环量级达100到200次，以仿真模具内腔的填充，当完成时，然后仿真进行到步骤60，正如在图2中的流程图中所描述的，其中执行压实阶段分析。

下面参照图25，其示意表示根据图2概括了步骤60中的某些过程子步骤的压实阶段。步骤60中的压实阶段分析从在步骤50中的填充阶段形成的所有变量的初始状态开始，再次将质量守恒和动量守恒方程应用于步骤220中，对于至少一部分求解域，求解流度k，在步骤230中求解压力，在步骤240中计算速度。一旦已经解出流度和压力，可以直接计算速度。接着根据能量守恒原理，在步骤250中利用温度计算粘度。然后在步骤260中，模型进行检查，以了解压力是否已经收敛，其很可能不同初期的情况，如果为否，跌代重复步骤220到150，直到压力收敛。

一旦压力已经收敛，接着在步骤220中考虑对流传热的所起作用，其中实现显式拉格朗日温度分析，已说明在内腔中该减少的但仍继续流入熔化的聚合物。然后，在步骤280中求解温度方程中的其余的项，以说明传导，粘滞耗散和其它热作用。

一旦已经将所有的温度作用量化并结合到仿真中，在步骤290中计算和更新部件特性例如密度、体积收缩、质量和冷凝体积。之后，在步骤300该算法进行检查，以了解是否完成PREC的压力分布曲线。当完成时，仿真则前进到步骤20，如在图2的流程图中所示的，其中为了由设计工程师考虑分析输出仿真结果。

虽然在填充阶段和压实阶段中的很多步骤是相似的，对于总的仿真的预计精确度分别所起的作用是变化的，特别是对于应用能量守恒原理而言。例如，在填充阶段期间，其中存在相对高的流量，粘滞耗散热作用可能是很明显的，而在压实阶段期间，其中存在相对低的流量，传导和对流所起的作用决定了传热分析。

根据本发明的一个实施例，关于与能量守恒原理相关的以上讨论的很多方法，可以如下的技术实施。最好在填充阶段流程图中的步骤190中：利用一维解析函数来描述在一节点的局部温度分布；限定一维解析函数的随时间的变化，以说明热对流；利用一维解析函数描述显式温度对流模式和有限元精细化，包括到内腔腔壁的距离的计算算法。在步骤200可以实现说明粘滞发热的一维解析函数的描述。这些技术也可以在压实阶段流程图中的对应步骤中实施，然而，由于不同的条件在仿真结果中它们的影响程度较小。

与之相似，关于与固化材料的处理有关的一些精细化，在填充阶段的步骤220中在对流度求解时可以实施如下的技术：确定固体，液体界面的位置；将所有按照材料特性线性变化的元以式表示；确定在包含固体、液体界面的元中的有限粘度函数，从求解域中消去冷凝的节点和元，在步骤200中可以有效地实现在已固化的区域中的有效压力的计算，以便为求解压力。这些技术也可以在填充阶段的对应步骤中实现；然而，在仿真结果中它们的影响程度较小。

虽然这里已经对被认为是示范性的和优选的实施例进行了介绍，对于本领域的技术人员来说很明显，根据这里作的介绍，本发明可以有其它改进和替换方案。所有这些改进和替换方案均认为在本发明的范围内，例如虽然该公开总的说是针对热塑性材料的透入模制的模拟，但本发明的技术可应用于热固性聚合物的注入模制。如上所述，在这些仿真中，在能量守恒分析中可以考虑热反应的所起作用，正如在热塑性材料模制中的固化热一样。将关于其它传热和热能的考虑可以结合到分析中，例如模具的热容量及模具的实际的冷却或加热。

作为其中作为基础的一般原理的技术也可以应用于在同一模具中，同时或顺序地将不同材料的流动进行模拟，以及其它材料和流体的模拟。例如，这里提出的本发明的某些或全部原理可应用于其它模制过程。这可以包括聚合物处理过程，例如挤压、吹模、加压模制、热成塑和注入模制各种变型例如气体辅助注入模制、加压注入模制和共通入模制。此外，这些原理也可以应用于气体工业例如金属浇铸以及食品处理即到处为流体或其行为类似于流体的物质的流动。

因此，希望按照如下权利要求和其等效物限定和区分的本发明中由书面的专利保护。

Claims (36)

1.一种用于模拟流体向一限定三维内腔的模具的注入的方法，该方法包含步骤：

(a)提供一个限定该内腔的三维固体计算机模型；

(b)根据该固体模型将求解域离散；

(c)指定边界条件；

(d)对求解域中的至少第一部分对于填充阶段过程变量求解，以便对求解域中的至少第一部分提供对应的填充阶段的解；

(e)部分地基于在填充终止时的过程变量的对应状态，在求解域的至少某些第一部分上对压实阶段过程变量求解，以便因此而提供对于求解域的该至少某些第一部分的对应的压实阶段解；以及

(f)确定对应的填充阶段解和压实阶段解的至少其中之一是否可接受，其中步骤(d)和(e)的至少其中之一包括子步骤：

使用关于在对应的求解域部分内的第一多个内部节点或元中的每一个的过程变量的分布的第一描述；以及

使用对在求解域中的第二部分的过程变量的分布的第二描述，该求解域包括第一多个内部节点或元，该第二描述使用质量守恒、动量守恒和能量守恒方程，

其中，该方法执行步骤(d)或步骤(e)，或步骤(d)和步骤(e)两者。

2.根据权利要求1所述的方法，其中在指定的时间步，第一描述论述独立于其它节点或元的第一多个内部节点或元中的每一个。

3.根据权利要求1所述的方法，其中由关于第一多个内部节点或元的过程变量的第一描述提供的在第一点的过程变量的数值不必等于由关于第二多个内部节点或元的过程变量的第一描述提供的在第一点的过程变量的数值。

4.根据权利要求1所述的方法，其中由关于一节点或元的过程变量的第一描述提供的在第一点的过程变量的数值用在第二描述。

5.根据权利要求1所述的方法，其中的第一描述是一维分析函数或离散函数。

6.根据权利要求1所述的方法，其中第一描述为描述关于一个节点或元的温度或内能的分布。

7.根据权利要求6所述的方法，其中第一描述为或者近似为对在固体中的一维热传导的解。

8.根据权利要求7所述的方法，其中第一描述包含一误差函数或者由误差函数导出。

9.根据权利要求6所述的方法，其中第一描述为关于一个节点或元的温度分布的一维描述并包含一误差函数。

10.根据权利要求1所述的方法，其中填充阶段过程变量和压实阶段过程变量是从该由密度、流度、模具内腔填充时间、模具内腔压实时间、压力、变形速率、切应力、温度、内能、速度、速度斜率、流量、粘度和体积收缩组成的组合中选择的。

11.根据权利要求1所述的方法，还包含步骤：

(g)在对应的填充阶段的解和压实阶段的解的至少其中之一被确定为不可接受的情况下，改进经离散的求解域和边界条件的至少其中之一；以及

(h)重复步骤(d)到(g)，直到对应的填充阶段的解或压实阶段的解被确定是可接受的为止。

12.根据权利要求1所述的方法，还包含步骤：

以图形格式显示从由填充时间、压力、变形速率、切应力、温度、速度和粘度组成的组中选择的填充阶段的解。

13.根据权利要求1所述的方法，还包含步骤：

以图形格式显示由密度、压实时间、压力、变形速率、温度、速度、粘度和体积收缩组成的组中所选择的压实阶段的解。

14.根据权利要求1所述的方法，其中离散步骤(b)包括根据模型通过将模型再分割为多个由多个节点限定的相连的元来产生有限元网格的子步骤。

15.根据权利要求1所述的方法，其中所述边界条件是由液体成分、流体注入位置、流体注入温度、流体注入压力、流体注入体积流量、模具温度、内腔尺寸、内腔构形和模具分开平面，以及它们的变化组成的组中所选择的。

16.根据权利要求15所述的方法，其中利用质量守恒和动量守恒方程的求解步骤(d)和(e)中至少之一包含子步骤：

(i)求解对于求解域中的至少某一部分的流度(130，220)。

17.根据权利要求16所述的方法，其中利用质量守恒和动量守恒方程的求解步骤(d)和(e)中至少之一包含子步骤：

(ii)求解对于求解域中至少某一部分的压力(140，230)；以及

(iii)求解对于求解域中至少某一部分的速度(150，240)。

18.根据权利要求1所述的方法，其中该利用能量守恒方程的求解步骤(d)和(e)中至少之一包括子步骤：计算对于求解域中至少某一部分的粘度。

19.根据权利要求18所述的方法，其中该粘度计算子步骤是基于温度的。

20.根据权利要求1所述的方法，其中求解步骤(d)和(e)中至少之一包括子步骤：迭代地计算速度和粘度中至少之一，直到压力收敛为止。

21.根据权利要求15所述的方法，其中求解步骤(d)包括子步骤：根据速度确定(180)内腔中流体的自由表面渐屈。

22.根据权利要求14所述的方法，其中求解步骤(d)和(e)中至少之一包括子步骤：根据对流传热作用、传导传热作用和粘滞耗散作用中至少之一计算温度。

23.根据权利要求21所述的方法，其中确定自由表面渐屈，直到内腔被充满为止。

24.根据权利要求1所述的方法，还包括步骤：

(g)计算一部件的质量特性。

25.根据权利要求24所述的方法，其中该质量特性是从部件密度、体积收缩、部件质量和部件体积组成的组中所选择的。

26.根据权利要求24所述的方法，其中迭代地计算速度、粘度和质量特性中至少之一，直到完成预定的压力分布曲线。

27.根据权利要求24所述的方法，其中该网格产生子步骤包含在厚和薄的区域中产生各向异性的网格，使得在基本上沿纵向不增加网格精细化的情况下，该网格精细化提高沿厚度方向的分辨率。

28.根据权利要求1所述的方法，其中步骤(d)或(e)包括子步骤：在求解能量守恒方程时利用显式模式。

29.一种用于模拟流体向一限定三维内腔的模具中的注入的方法，该方法包含步骤：

(a)提供限定该内腔的三维固体计算机模型；

(b)根据该固体模型将求解域离散；

(c)规定边界条件；以及

(d)利用关于求解域中的至少一部分的质量守恒、动量守恒和能量守恒方程对于过程变量求解，其中步骤(d)包括子步骤：在求解能量守恒方程时利用显式模式。

30.根据权利要求29所述的方法，其中步骤(d)的显式模式包括子步骤：利用一维分析函数、从一维分析函数导出的数据或一离散函数来描述关于一节点周围的温度分布。

31.根据权利要求29所述的方法，其中步骤(d)的显式模式包括子步骤：利用局部变化并且可以与全局时钟不同地行进的热时钟。

32.根据权利要求31所述的方法，其中该热时钟是一节点热时钟。

33.根据权利要求31所述的方法，其中热时钟的使用包括子步骤：计算Peclet数，以便估计在节点处的对流和传导的相对分布。

34.根据权利要求29所述的方法，其中步骤(d)包括确定在对应于在前一时间步的质点位置的上游侧位置处的温度。

35.根据权利要求29所述的方法，其中步骤(d)包括子步骤：计算由于粘滞耗散、压缩热、减压热、固化热以及反应热中至少之一引起的传热作用。

36.一种用于模拟流体向一限定三维内腔的模具中的注入的计算机硬件设备，该设备包括：

存储器，用于储存限定指令集的程序代码；

处理器，用于操作来运行由所述程序代码提供的指令，其中，所述程序代码促使处理器：

(a)根据限定内腔的三维固体计算机模型离散化求解域；

(b)规定边界条件；

(c)在求解域的至少某些第一部分上对填充阶段过程变量求解，以由此而为求解域的至少第一部分提供相应的填充阶段解；

(d)部分地基于在填充终止时的过程变量的相应状态，在求解域的至少某些第一部分上对压实阶段过程变量求解，以由此而为求解域的该至少某些第一部分提供对于相应的压实阶段解；以及

(e)确定相应的填充阶段解和压实阶段解中至少之一是否可接受，其中步骤(c)和(d)中至少之一包括子步骤：

使用关于在求解域的相应部分内的第一多个内部节点或元中的每一个的过程变量的分布的第一描述；以及

使用对在包括第一多个内部节点或元的求解域的至少第二部分中的过程变量的分布的第二描述，该第二描述使用质量守恒、动量守恒和能量守恒方程，

其中，所述程序代码促使处理器执行步骤(c)或步骤(d)，或步骤(c)和步骤(d)两者。

Images