航空装备腐蚀数字孪生-数字模型验证及确认（上）

确认实验与传统实验的主要区别在于：传统实验非常强调在受控环境中对过程进行测量。只有在受控环境中，实验人员才能可重复的测量物理过程。而在确认实验中，实验的受控程度和可重复性并非处于重要位置，取而代之的是要准确地测量非受控的实验条件。如果只要准确测量出装备服役的环境条件，则由于气候和天气等造成的环境参数变化就变得不再模糊。对于非受控的确认实验，则需要合理的实验设计及大量测试，细致地识别系统激励、输入和响应，同时充分分析和量化服役环境的动态特性。

航空装备模型确认的分层思路：由于航空装备等复杂系统中包含大量的不确定性因素，因此模型确认工作需要通过将复杂系统层层分解来降低其复杂性，这种思想类似于IT算法的分治理念。为了针对性地进行理论分析，建立数值模型，以及严谨的实验设计，复杂系统通常可以分为全系统、子系统、组件和单元四个层次或更多，如图4所示，系统层次越低，影响因素越少，耦合程度越低。

     

一般的分层实施思路：将全系统按照功能划分为不同的子系统，子系统按照不同的学科知识划分为组件或基准，最终按照学科专业的解耦划分到单元层级。图5给出了航空装备典型腐蚀模型分层示意图。

     

图5中由上往下第一层为全系统层，该层应包含航空装备系统所有真实组件，实际物理层应包含几何结构信息、材料属性数据和组合而成的诸多不同的功能系统等；数值模型层包含所有组件的耦合关系、真实边界条件、初始值和各类激励源。对于全局层的腐蚀整体性能退化，一般在设计阶段只能依靠少量的样机实验获取数据，大多数的数据来源依然是装备服役期间的运行或实验数据。

图5中由上往下第二层为子系统层，根据常规的航空装备功能划分，每个子系统应独立承担一个功能，如结构、动力、电子和武器等子系统，各个子系统中仍然包含了较多的耦合关系，如动力系统的腐蚀应包含应力腐蚀、高温氧化和环境腐蚀等。但相对于全系统层，子系统层相应的测试实验数据相对比较容易获取。

图5中由上往下第三层为基准层，即图5为了简化分层将单元层与基准层进行了合并处理，因为整体分层是面向腐蚀学科的，并不包含信号工程、能源动力工程等学科，所以这种简化操作是相对合理的。按照严格意义上的基准层定义，每个基准案例都是相对独立的（各个基准案例之间解耦处理），且只考虑通用的2～3个物理现象。单元层在基准层的框架下进一步解耦处理，每个单元只包含1个物理或化学现象。对于腐蚀仿真的分层构想图，基准层的每个案例以某种具体的腐蚀案例为主导，考虑一定的环境或工况不确定性。处于这一层面的组件，其几何结构、物性参数及环境条件容易获取，因此可以得到数量较多、质量较高的测试数据，是三个层级中不确定性最低的层级。

对仿真系统模型的验证确认技术研究最早开始于对仿真模型的校验，至今已有40多年的历史了。20世纪90年代以来大型复杂仿真系统研究取得了很大的发展。国外从20世纪80年代中期开始对分布交互仿真技术进行研究和开发，1995年美国国防部又提出了高层体系结构。对模型V&V研究的重点从仿真模型的校验方法研究为主转向如何更加全面系统地对仿真系统进行模型V&V上来。

进入20世纪70年代，随着工业和航天业的迫切需要，美国电气与电子工程师学会（IEEE）对V&V理论进行了新的定义，验证：评估开发极端的软件产品能否满足前一阶段对它们定义的要求的过程。确认：测试计算机程序并评估结果，以确保符合具体要求的过程。可以看到，IEEE对V&V的定义都突出和强调了“过程”，表明V&V应用于M&S的全寿命周期且持续进行，“前一阶段”也体现了“验证”过程的多样复杂的特点。

到20世纪90年代，美国大力发展国防和推进军事现代化发展，美国国防部成立了国防部建模与仿真实验室（DMSO），研究并判定IEEE针对V&V的定义是否符合他们的需求，公布了密切结合他们需求的定义，验证：确定模型实现是否准确表示了开发人员对模型的概念描述的过程；确认：从模型预期用途的角度出发，确定模型对现实世界的表述准确程度的过程。通过对比不难看出，美国国防部较清晰地描述了V&V在建模和仿真中的重要性，保留了其“过程”的特性，同时，比IEEE所定义的V&V更全面和侧重工业应用。1992年NASA对仿真的投资范围包括仿真算法、应用、网格剖分、数据展示、流体动力学模型、验证和确认，总投资1399.8万美元左右，而模型的验证和确认这一项就投资了713.5万美元，占总投资的50.97%，可见NASA认定的模型验证和确认在仿真技术中的重要地位。

由于DMSO对V&V的定义十分全面和宏观，因此也基本适用于各个行业。但定义只能解释“是什么”，并不能够告诉我们“怎么做”。因此，为落实V&V在建模和仿真中的作用，计算流体力学（CFD）学界和固体力学学界分别制定了适用于本行业的指南，即著名的《AIAA指南》和《ASME指南》。

《AIAA指南》全称为《计算流体力学中的验证与确认》，定义了许多基本概念和专业术语，并指出了如何在CFD中执行V&V的流程，为其他领域展开V&V提供了思路和方法。《ASME指南》是在AIAA指南提出的8年后，由固体力学界的计算固体力学确认与验证委员会提出的。《ASME指南》首次提出了在M&S过程中可将V&V分开的想法，并将它们分别称为数学建模和物理建模，并指出虽然V&V是持续在M&S的全过程，但可以采用分层、分级的方式逐步进行验证与确认，且层与层之间可以存在不同的V&V要求。

模型验证主要包含两方面工作：第一是代码验证，第二是解验证。

《ASME指南》（ASME，2006）对“代码验证”的定义如下：确定在计算代码中正确实现数值算法并识别软件中错误的过程。

代码验证可以分为两种活动：软件质量保证（SQA）和数值算法验证。数值算法验证用于确定会影响计算结果的数值精度的所有数值算法的软件实现中的数学准确性，其主要目的是搜集证据，证明代码中的数值算法得到了正确实现，且能够如预期般正常使用。软件质量保证（SQA）着重确定是否正确实现了构成软件系统一部分的代码，以及代码在指定的计算机硬件和软件环境中是否得出可重复使用的结果。此类环境包括计算机操作系统、编译程序、功能库等。建模与仿真中代码验证的综合视图如图6所示。

解验证，也称为计算验证，定义如下：确定输入数据的准确性、解的数值精度及某个特定仿具中输出数据的准确性的过程。

解验证尝试确定并量化执行计算机仿真代码时出现的三类误差源。第一类为计算分析师在准备计算机模拟代码输入时造成的错误、误差或者差错；第二类为在数字计算机上计算数学模型的离散解时出现的数值误差；第三类为计算分析师处理由仿真代码生成的输出数据时造成的错误、误差或者差错。第一类和第三类错误来源与第二类有很大区别。第二类误差源不包括构建或者组建数学模型时出现的误差或者近似值，而第一类和第三类误差源包括了人为误差，但不包括任何其他误差源。在复杂系统的大型计算分析中，人为误差是很难发现的。即便在较小规模的分析中，如果未刻意采用程序性或者数据校验法来检测可能的误差，那么仍旧无法发现人为误差。解验证中处理的误差源如图7所示。

     

2.1.1  软件质量保证

软件质量保证确保了计算机代码真实贴近原始的数学模型，主要包括开发过程确立、版本控制、验证与确认、质量和可靠性评估。

第一步，开发过程确立，主要包括确认编程语言、架构设计和软件开发过程模型。

（1）确认编程语言，科学计算所用的大多是现代、高级编程语言都支持程序编程和面向对象的编程，选择编程语言时需要考虑多重因素，如涉及的数学计算及数学关系、软件缺陷率和后期维护成本的影响等，需要根据预计的代码规模、后期可维护性等综合因素选取适合的编程语言。

（2）架构设计，是在开发工作开始前确定软件子系统及其接口的过程，其输出通常是包含流程图、伪代码等形式在内的描述软件子系统及其结构的文档。

（3）软件开发过程模型，又称为软件生命周期模型，主要有4个：瀑布模型、螺旋模型、原型化模型及迭代模型，这些模型都包含了需求、设计、执行、验证确认的行为，其区别在于管理者组织和实施这些行为的方式。对于传统瀑布模型，各行为阶段间相对独立，只有前一阶段完成后才开始下一个阶段，如图8所示；螺旋模型是综合了瀑布模型的优点，同时引入了风险分析和成本管理的思想，在此基础上各阶段行为具有迭代特性，如图9所示；原型化模型是适用于用户不能对软件提出十分明确需求的情况，此时通过快速原型开发与用户共同确认需求、概念模型和最终产品；迭代模型适用于开发新的模型，也可用于需要对现有模型进行大范围再利用的场合，迭代模型的每轮开发中都包含了瀑布模型的四大阶段，但并非所有的迭代都要重新执行完整的四个阶段，如有些新一轮的开发可以使用已有的需求分析和概念模型。

     

     

第二步，版本控制，指对软件开发过程中各种程序代码、配置文件及说明文档等文件变更的管理。一个良好的版本控制系统可以让用户清楚地了解到更改的内容及更改时间等信息，同时允许软件开发人员撤销对代码的任何更改，并追溯至之前的版本。使用版本控制工具的主要概念及定义如表1所示。

     

第三步，验证与确认，科学计算中的软件确认与验证和模型V&V有所区别，此处的验证指的是确保仿真软件满足各类软件工程规范及设计文档，而确认指的是确保软件满足用户的实际使用需求。表2简单介绍了科学软件正确性的典型评估方法。

     

通常软件工程中的测试案例可以采用用户需求、标准规范、详细设计和源代码，而测试层级分为验收测试、系统测试、集成测试、模块测试和单元测试。图10常被称为软件开发与测试的V字模型。需要强调指出，无论使用哪种测试方法，都不可能使用所有可能的输入值组合进行测试，因为这意味着海量的测试工作量，因此如何采用尽可能少的测试案例来暴露最多的代码问题，就必须引入测试覆盖率原则。

     

第四步，质量和可靠性评估，通常的软件质量定义为是否满足客户需求和要求，以及是否符合规范，可靠性即在确定的条件下软件无故障运行的概率。评估软件可靠性的定量方法有缺陷密度分析和复杂性分析。

缺陷密度分析指评估软件中的缺陷数量，其缺点在于比较依赖于测试程序的完善度，如当测试程序不够完善时，可能导致某些缺陷不能被发现。

复杂性分析是将代码质量属性转换为代码可靠性，一般可以认为，当一个子程序越复杂，代码的可靠性越低，反之则可靠性越高。对于代码可靠性的评估可以从代码行数、NPATH复杂度、循环复杂度、条件嵌套深度4点入手。代码行数，即统计每个子程序源代码的数量，这是最简单直接的衡量复杂度的方法，建议子程序的可执行源代码数控制在250以内；NPATH复杂度指经过某个子程序的可能执行路径的数量，一般来说，此数量应该控制在200以下；循环复杂度指子程序中的循环语句或逻辑语句的数量加1，一般上限值建议为10；条件嵌套深度指IF语句的嵌套深度，一般更高程度的嵌套更容易出错，且更难以理解和追溯。

2.1.2  数值算法验证

在模型验证过程中，数值算法验证的目的在于确保代码能够真实地表示出数理模型（也称为控制方程），其形式以带有若干初始值、边界条件和耦合关系的微积分方程组成。因此，数值算法验证解决的问题是数值算法的正确性和将该算法转化为源代码实施的正确性。

在Oberkampf编写的《科学计算中的验证和确认》中，重复使用了代码验证的概念，即在代码验证过程中将数值算法验证再一次称为代码验证。本文为了避免读者混淆，直接采用数值算法验证描述相关内容。表4-3按照验证的标准严格程度，由低到高排列给出了5种方法。

     

表3中，简单测试和代码比较两种方法的优点在于可以用于没有精确解的数学模型，但它们的严格程度也最低。其余三种验证方法的使用要求数学模型存在解析解或可经证实的替代解。

第一级，简单测试，具体包含三种测试方法：

（1）对称性检验，即当研究对象具有对称的几何结构，且初始值和边界条件也存在对称关系时，该对象的数值解也是对称的，这种对称性可以作为检验条件；

（2）守恒测试，大多数模型遵循物理量的守恒，如电化学腐蚀类模型遵循电量守恒、质量守恒，应力腐蚀模型遵循能量守恒和质量守恒，则守恒性可以作为检验条件；

（3）伽利略不变性测试，即认为物理定律对所有惯性系都有效，因此，数学模型和离散方程的解需要遵守伽利略不变性原则。需要指出该测试是基于经典力学或牛顿体系的，对于原子/分子尺度的模型极有可能并不适用。

第二级，代码比较，是评估算法 正确性最常用的方法之一，本质是将两种不同算法的输出进行对比分析。适用于两种算法代码的是相同的模型或者模型存在经过专家/权威评估验证的基准代码。需要指出，存在两种算法代码均存在谬误且其计算结果一致的可能，因此代码比较不能取代严格的代码验证。

第三级，离散误差，适用于存在解析解的模型，是较为常用的验证方法。它使用单一的网格间距和时间步长，定量评估数值解与解析解之间的误差。该方法的主要缺点在于离散误差评估完毕后，还需要进一步给出误差大小的评价方法。

对于瞬态模型，如研究随着时间变化的腐蚀速率。航空装备构件的腐蚀模型一般同时存在空间网络离散和时间步长离散（少有稳态模型），因此可以先验证空间网格离散，其次验证时间步长离散。这种方法虽然简单易行，但是存在无法发现空间网格离散和时间步长离散潜在交互性风险，尤其是当对流传质占优的时候，Pe数较大，网格尺寸和时间步长的交互作用必须同时进行耦合验证，这种耦合验证一般较为耗时费力。空间和时间数值算法验证中使用的网格级如表4所示。

     

由于时间只有一个方向，所以时间步长细化较为简单，本文暂且不进行讨论。而网格的细化设计空间网格离散，随着研究对象几何结构复杂程度的变化，网格细化对网格整体的分辨率和质量影响较大，而且在实际仿真工作中，网格的剖分往往也占据50%左右的时间，故本文从网格的均匀性和一致性方面对网格细化进行讨论。

网格质量的量化指标一般由单元长宽比、偏斜率和拉伸率（网格有粗到细过度的尺寸变化速率）组成。均匀性的网格细化指的是在各个方向上，单元长宽比较为一致的加密网格，如图11所示。

     

需要指出，网格细化的均匀性是在整个区域内整体加密的，加密前后几何区域内网格的整体分布趋势与细化前相同。

网格细化的一致性指的是网格加密时，网格质量不变或者有所提高。图12所示网格一致性加密示意图。若采用Knupp提出的全局网格质量度量s，即当s =1时，网格质量最好，如正方形、正方体和等边三角形等；当s =0时，网格质量最差，会出现网格反转或交叉等现象，一般s 较小时，网格的性能指标如偏斜度、拉伸率和曲率都比较大，因此可以用s 来表征网格一致性。无法回避的事实是，随着网格的加密及s 趋近于1，网格各项同性的趋势也会越来越明显，网格剖分的工作量也会大大提高，并进一步影响求解效率。尤其是对于航空装备中复杂的几何问题，网格一致性和求解效率是十分矛盾的。

     

第四级，收敛测试，评估相对于数学模型的解析解，数值解中的误差（离散误差）是否随着网格尺寸和时间步长的细化而减小。它除了要评估离散误差的大小，还需要确保该误差是否会随着网格尺寸和时间步长细化而减小，因此收敛测试是严格代码验证所遵循的最低标准。需要指出，此处的收敛与日常所说的迭代计算收敛有一定区别。此处的收敛指的是当网格间距与时间步长减少到极限时，离散方程的解逼近数学模型的解析解，迭代计算收敛往往指的是计算结果是否满足求解算法的容差或迭代次数限制（一致性，即离散方程需要随着离散参数接近零而使数学模型方程接近极限。稳定性，指的是对于初值敏感问题，数值误差在迭代过程中是否趋于可控，若可控则认为离散方法是稳定的。需要注意的是，数值算法稳定性只适用于离散方程，与数学模型本身的可导性或连续性无关）。

第五级，精度阶验证，是最严格的数值算法验证标准，它不仅检查数值解的收敛与否，还关注离散误差是否随着网格尺寸和时间步长的进一步细化以理论速率（形式精度阶）减少。与理论速率相对的是离散误差的实际减少速率（观测精度阶），它需要两个级别的网格或时间步长进行计算。精度阶测试是模型验证过程中最难满足的测试，它对算法中细微的缺陷都非常敏感，如电流密度边界条件，坐标变换或者浮点数精度等。由于其严苛的特性，故图13给出了稳态问题阶验证程序的详细步骤。

     

目前，模型精度阶验证的方法有4种，Thomas等人认为不同验证方法之间的对比如表5所示。

     

残值法实质上是将模型的解析解直接代入离散方程获取截断误差，由于不需要任何迭代，所以误差评估简单易行。但其缺点显而易见，如无法评估常见的第一类边界条件（对于材料力学而言就是直接以应力和位移作为边界条件，对于电化学而言就是直接以电压或电流作为边界条件，对于传热学就是温度边界条件等），无法适用于解耦方式求解的偏微分方程组。

统计法由Heber等人提出，只采用一套网格，通过在关注区域上不断连续缩小尺寸比例，然后随机抽样进行获取数值离散误差。其优点在于不需要系统性的细化网格，故消耗时间较少；其缺点在于统计方法的收敛性需要得到保证，且随着尺寸比例的缩小关注点会向边界位置靠拢等。

降尺度方法类似于统计法的升级版，它具有诸多统计法的优势，但解决了统计法收敛性的潜在风险。进一步通过选择网格缩放的焦点，使得我们可以考虑航空装备复杂几何结构情况下的腐蚀边界条件验证。不过，由于降尺度方法忽略了离散误差传递至缩小域中的可能性，因此对实际收敛速度的判断偏向乐观。

如图13第4步所述，在阶验证的过程中，需要构造模型的精确解，同时在前述的离散误差评估和收敛测试两个验证过程中，也需要模型的精确解进行支撑。由于数值计算通常涉及复杂的不同物理场耦合的偏微分方程组（PDEs），如应力腐蚀涉及弹塑性方程、损伤动力学方程，甚至包含热传递方程，冲刷腐蚀涉及N-S方程、粒子追踪方程及界面损伤方程等，这些偏微分方程组往往难以获取解析解，因此研究人员发展出不同的方法，如虚构解和替代解方法用于模型验证工作。模型验证中高精度解构造方法如表6所示。

     

PDEs一般有两种形式，微分形式和积分形式，两者依次又称为强形式和弱形式。强形式的PDEs要求解必须是可微分，而弱形式的解可以包容不连续、不可微的情况。因此，弱形式的解只能有限的满足PDEs的要求，本文考虑的均是强形式问题，但需要指出强形式的解必然满足弱形式要求，故上述精确解验证的内容同样适于常见的有限元方法。

在科学计算中，为了确认特定数学模型的仿真精度是否达到预期要求，就需要进行解验证。它包括代码输入、代码结果的精度、所关注情况的仿真精度，以及涉及的任何后处理工作的精度。解验证主要包括三个方面，输入数据验证、后处理工具验证和数值误差估算，所含内容如表7所示。

     

输入数据就是运行科学计算代码所需的信息。例如，计算模型中的几何结构、边界条件、材料物性等。

后处理工具指的是作用于科学计算代码输出的任何软件。如果后处理涉及任何类型的数值近似，如离散、积分、内插等，那么应验证这些工具的精度阶（如通过阶验证），否则需遵照标准软件工程实践。常用后处理软件有origin、tecplot、ParaView等，origin是科学绘图、数据分析软件，功能包括数据统计、曲线拟合及峰值分析；tecplot（见图14左）是绘图和数据分析的通用软件，具有丰富的数据可视化功能，如曲线绘制、标量场展示、2D和3D矢量场可视化等；ParaView（见图14右）是一个开源跨平台的应用程序框架，不仅可以对2D和3D数据进行分析和展示，还支持并行数据处理，它采用C++编写，可以进行二次开发。

     

数值误差包括四种：舍入误差、统计抽样误差、迭代误差、离散误差。

舍入误差：在仿真中，计算机将实数存储在有限的内存中，为了使内存足够存储这些数字，需要限制数字的精度和范围，这样造成的误差称为舍入误差。各种编程语言指定的浮点精度如表8所示。

     

注：a表示与编译程序相关，b表示可通过“kind”属性访问。摘自《科学计算中的确认与验证》P193。

统计抽样误差：如果仿真稳定且确定，将不会有统计抽样误差，仿真计算中的系统响应量大多为平均量，统计时会产生抽样误差，一些科学计算方法本身就具有随机性（如有限元法，直接模拟蒙特卡洛、格子玻尔兹曼），通过增加实现次数、迭代次数或时间步长，通过对系统响应量的收敛情况进行评估，从而估算统计抽样误差。

迭代误差：控制方程中近似解与精确解之间的差，使用迭代法（或者松弛迭代法）求解代数方程，任何时候都可能出现这类误差。迭代法按效率升序，稳健性降序可以归纳为试位法、正割法、牛顿法。迭代收敛准则可以依靠迭代之间的差或者迭代残值进行确定。迭代误差的估算方法有机器零位法及局部收敛速度法。

离散误差：数学模型解析解与离散方程数值解之间的差。所有的仿真过程中几乎都会产生离散误差。数值求解过程中对数学模型的离散会造成离散误差。离散过程主要有两种，数学模型的离散和域的离散。

数学模型的离散：数学模型中的控制方程基本由偏微分或者积分方程组成。离散方法有：有限差分法，用有限差分代替方程中的偏导数；有限体积法，用单元格的数值平均值代替体积积分；有限单元法，用加权残值法或者变分原理将偏微分方程化为代数方程。

域的离散：仿真需要对数学模型的定义域进行离散处理。对于与时间无关的稳态问题，将空间划分为节点、单元格或者元素，应用于离散方程。对于与时间有关的非稳态问题，仿真也通过把时间域分割为合适的时间步长Dt来实现。对于复杂的三维几何体，需要投入大量的精力与时间才能生成合适的空间网格。离散误差也极大地受到网格质量和粗细程度的影响。有两种基本的网格划分方式：排序规则的结构化网格，这种网格在复杂几何体上应用时的难度也更大；无规则的非结构化网格，这种网格生成时较为简单。结构化网格和非结构化网格如图15所示。

     

网格质量由网格收敛指标（GCI）确定，由Roache提出，利用绝对值将离散误差估值转换为不确定度估值，即

式中，Fs为安全系数，f1为细网格解，f2为粗网格解，rp为网格细化系数。

对于复杂的模型，使用网格细化来缩小离散误差的计算成本十分高昂。在三维模型中，使用网格对分法进行网格细化时，单元数量会成8倍增加，计算成本也会提高8倍。因此，一种成本较低的缩小离散误差的方法是针对性的局部解自适应。解自适应主要包含确定自适应的区域及完成自适应的方法两个方面。

推动自适应过程中找到一个合理的标准是解自适应最大的困难之一，经常使用基于解值或者梯度的高阶重构方法，或者基于对局部截断误差或者残差的估算。如果使用推动自适应过程的方法，一般最终会在空间域上得到一个变化的权重函数，利用这个权重函数，即可运用不同的方法实现局部解自适应。一般从粗网格开始，使用网格生成工具以递归的方式优化网格，完成自适应网格的重划分。自适应的方法主要有四种：基于解梯度的自适应、基于解曲率的自适应、基于离散误差的自适应和基于截断误差的自适应。

有两种方法实现不基于网格重划分的网格自适应：网格细化和网格移动。网格细化方法细分权重函数较大的单元和粗化权重函数较小的单元。网格移动要保持单元数量和网格的连通性，可采用线性弹性方程、变分法、拉普拉斯法等，也可以同时使用网格细化和移动方法。

根据工程实践经验，大多数解验证活动重点关注仿真中的数值误差估算，上文讨论了舍入误差、统计抽样误差、迭代误差和离散误差，前3个数值误差必须足够小才能不影响离散误差的估计。仿真中舍入误差和离散误差是始终存在的，而迭代误差和统计抽样误差取决于应用和所选算法。科学计算中数值误差的量化分析还涉及不确定性估算和模型确认的计算结果精度问题，这些将在下篇展开叙述。（由于全文较长，本文将分为上下两篇）