★ FEATURED ARTICLE 这个系列写到第六篇了。前几期我们一直在打磨一件事怎么让GPT真正成为会写计算代码的助手而不是一个只会复读示例的聊天机器人。从提示词基础、代码骨架拆解到让GPT自我检查一步步来。今天我想聊一个更贴近实战的场景——单晶材料的弹性各向异性分析。手里只有一组弹性刚度常数比如铜的C11、C12、C44要算出单晶沿不同方向的杨氏模量并把三维各向异性曲面画出来。这事用Python加NumPy、Matplotlib做并不难难点在于物理公式的约定、Voigt记号的坑、数值验证的锚点都藏在一行行代码里。如果提示词写得不够清楚GPT生成的代码大概率能跑跑出来的数却可能错得离谱。这篇会用完整实例把需求拆解→提示词设计→GPT生成→人工审查→物理验证整条链路走一遍适合正在用GPT辅助编程的计算材料、力学方向学生也适合所有想知道提示工程到底怎么落地的Python用户。1. 为什么选弹性张量可视化当这一期的实战案例1.1 能跑的代码不等于算得对的代码用GPT生成代码最迷惑人的地方就是它给你的程序往往能直接运行结果也打印得整整齐齐甚至图表都画出来了。但如果你没有一套物理上的验证方法这些漂亮的输出可能全是错的。这种可运行但不可靠的问题在纯Web开发、爬虫脚本里可能影响不大但在计算材料科学和力学编程里是致命的——你拿错的数据可能会写进论文、进报告、影响产品选型。弹性常数的计算就是一个典型场景。它看起来只是矩阵求逆、方向向量点乘、画一个三维曲面但中间藏着好几种行业约定。GPT训练数据里相关代码很多但来自不同教科书、不同软件、不同版本的脚本约定并不统一。你不把约定钉死它就会自由发挥它一自由发挥结果就可能在某个细节上全盘出错。所以这一期的核心并不是让GPT写一个能画图的脚本而是怎么通过提示词让GPT写一个物理上站得住脚的脚本。1.2 一个看似简单但信息量很大的材料问题单晶材料的弹性各向异性在材料设计里非常有用。比如半导体衬底要选晶向以减少热应力高温合金叶片要理解不同晶粒取向的刚度差异锂电池正极材料颗粒在充放电过程中的应变也和各向异性弹性密切相关。描述这个问题的标准方式是画一张杨氏模量方向曲面图从原点出发朝每个晶向画一条射线射线的长度等于该方向的杨氏模量。球面的凸起越明显说明材料弹性各向异性越强。计算过程本身不复杂。只要有了材料的弹性刚度矩阵C6x6Voigt记号取逆得到柔度矩阵S然后用方向余弦n去缩并S就能得到任意方向的杨氏模量。难就难在缩并的写法、单位、归一化和坐标约定。这恰恰是提示工程能发挥最大价值的地方你不一定要把公式背得滚瓜烂熟但你必须知道有哪些地方容易出错然后把这些不确定性通过提示词交给GPT处理并让它输出可验证的结果。1.3 这一篇做完你能拿到什么按我的经验看完这篇并且跟着跑一遍你至少能带走四样东西。第一一套可复用的弹性常数分析提示词模板以后换材料、换晶系都能改改就上。第二一条需求拆解→对话迭代→物理验证的GPT协作思路这套思路不限于弹性计算做其他物理或力学计算也通用。第三几个典型坑的识别方法尤其围绕Voigt记号、采样方式和数值锚点。第四一个从临时脚本升级成小工具的工作流让GPT生成的代码不只跑一次而是沉淀下来反复用。2. 动手之前先把物理和数据的约定钉死2.1 弹性常数到底怎么进程序先说输入数据。这一篇我们用铜Cu做例子用的是经典单晶弹性常数单位GPaC11 168.4C12 121.4C44 75.4。这三个数组成立方晶系的6x6刚度矩阵Voigt记号排布如下import numpy as np C np.array([ [168.4, 121.4, 121.4, 0.0, 0.0, 0.0], [121.4, 168.4, 121.4, 0.0, 0.0, 0.0], [121.4, 121.4, 168.4, 0.0, 0.0, 0.0], [0.0, 0.0, 0.0, 75.4, 0.0, 0.0], [0.0, 0.0, 0.0, 0.0, 75.4, 0.0], [0.0, 0.0, 0.0, 0.0, 0.0, 75.4], ], dtypefloat)注意这里的Voigt索引映射11→022→133→223→313→412→5。这个映射必须和后续计算保持一致。GPT如果没被告知很可能默认一套它自己熟悉的映射或者干脆在缩并时把剪切分量搞混。我在提示词里会明确写出映射规则并要求它不要改动。关于单位还有一个细节有人习惯把C除以GPa有人习惯转成Pa再算。我建议全程用GPa因为柔度矩阵的数值在GPa量级下比较温和输出的杨氏模量也是GPa不容易产生数量级错误。GPT如果中途偷偷把单位换掉你最后的数字会莫名其妙多个10的9次方。后面验证环节这个单位问题要重点盯。2.2 杨氏模量的方向依赖一个公式理清所有编码逻辑对任意晶系杨氏模量的通用定义是1 / E(n) n_i · n_j · n_k · n_l · S_ijkl其中n是笛卡尔坐标系中的单位方向向量S_ijkl是四阶柔度张量。这个式子看着抽象但物理含义很直白朝某个方向施加单位名义应力测出该方向的应变倒数就是模量。GPT要写这一段最稳妥的方案不是让它去实现四阶张量缩并而是用立方晶系的解析公式1 / E S11 - 2·(S11 - S12 - S44/2)·(n1²n2² n2²n3² n3²n1²)这个公式对立方晶系完全准确而且只需要柔度矩阵S里的三个独立分量S11、S12、S44。为什么强调这个公式因为它能大幅降低Voigt记号带来的实现风险。你不需要把6x6矩阵展开成3x3x3x3的四阶张量不需要处理工程剪应变和物理剪应变的换算一个等式就能算完所有方向。对于初学者这是最不容易翻车的路线。拿到S后三个分量这样取S11 S[0,0]S12 S[0,1]S44 S[3,3]。然后方向向量必须归一化。比如[111]方向Miller指数对应笛卡尔向量(1,1,1)归一化要除以√3。这个归一化看着简单但GPT写循环时很容易漏掉或者只在某个分支里做了最后极值点找错位置。2.3 五要素法拆需求喂给GPT之前你自己先得想清楚我每次让GPT写这类物理计算脚本都会先在草稿纸上按五个要素过一遍需求。这一步花五分钟能省后面两个小时调试。输入数据弹性矩阵C是什么、单位是什么、是硬编码还是读文件。物理模型用什么公式、什么约定、适用于什么晶系。计算操作求逆、采样、遍历、找极值。输出形式表格还是三维图极值要不要标注。验证锚点哪些方向的数值有解析解程序跑完得和这些值对得上。五个要素里验证锚点最容易被忽略。很多人只会要求画一张漂亮的图但不会要求GPT输出[100]、[110]、[111]三个方向的参考值。导致的结果是图看起来很正常实际模量差了30%也毫无察觉。这一篇后面所有验证手段本质上都是在补强这个环节。3. 三轮与GPT的协作从第一版脚本到带校验的完整工具3.1 第一轮提示词把任务说明书写清楚我不会一上来就让GPT给完整代码而是先给一段结构化需求明确告诉它物理背景、公式、输入输出和要避免的坑。你可以直接复制下面这版提示词试试你是一位计算材料学方向的资深Python工程师。我需要你帮我写一个Python脚本用于分析立方晶系单晶材料的弹性各向异性。输入弹性刚度矩阵CVoigt记号6x6单位GPa。我会在脚本中硬编码一个铜的示例但代码要能方便地替换成其他立方晶系材料。任务计算柔度矩阵S inv(C)。使用立方晶系杨氏模量公式1/E S11 - 2*(S11 - S12 - S44/2)(n1^2n2^2 n2^2n3^2 n3^2n1^2)。在球面上均匀采样方向计算每个方向的E。求出E的最大值、最小值以及对应的晶向并输出这组结果。用Matplotlib绘制三维曲面曲面半径表示该方向的E用颜色映射E的大小标注极值点。技术约束Voigt索引映射11-022-133-223-313-412-5。方向向量必须先归一化。全程单位使用GPa。请用中文注释关键步骤。最后请在代码里加入一个自检函数计算[100]、[110]、[111]三个方向的杨氏模量并打印。请注意我把计算哪些方向作为自检直接写进了需求。这不是为了让程序变复杂而是为了让GPT在写代码时心里有一个正确答案的约束。它知道后面要打印这三个值就不太敢在公式里乱来。3.2 第一轮产出基础脚本和人工检查按我的经验上面这份提示词生成的代码骨架通常已经能用。GPT给的核心计算部分清理之后大致长这样def young_modulus_cubic(S11, S12, S44, direction): n direction / np.linalg.norm(direction) term (n[0]**2 * n[1]**2 n[1]**2 * n[2]**2 n[2]**2 * n[0]**2) return 1.0 / (S11 - 2.0 * (S11 - S12 - S44 / 2.0) * term) S np.linalg.inv(C) S11, S12, S44 S[0, 0], S[0, 1], S[3, 3] for label, vec in [([100], [1, 0, 0]), ([110], [1, 1, 0]), ([111], [1, 1, 1])]: E young_modulus_cubic(S11, S12, S44, np.array(vec, dtypefloat)) print(label, round(E, 3))人工检查时我最先看三处。第一S44取的是S[3,3]不是C[3,3]也不是1/C[3,3]——因为我们已经取了逆了。第二公式里有没有S44/2这个半系数GPT经常把这个半系数丢掉那样算出来的E会整体偏大。第三方向向量有没有归一化。比如[110]向量长度是√2如果忘了归一化算出来的E会直接差2倍而且方向模依然正常很难一眼看出来。3.3 第二轮让GPT补上Born稳定性与输入合法性检查第一版能跑但我还会追加一轮提示词让GPT加入物理合理性检查。对弹性材料来说刚度矩阵必须满足Born机械稳定性条件。立方晶系的判据有三条C11 2*C12 0C11 - C12 0C44 0。用物理语言说分别是体积模量、正方剪切抗力、菱形剪切抗力都必须是正的。这轮对话的提示词不用太长直接追加请在上面的脚本里加入Born稳定性检查C112*C12必须大于0C11-C12必须大于0C44必须大于0。如果不满足程序要抛出带明确信息的异常并停止计算。为什么值得加这一步因为实际工作中我们会批量处理很多材料的数据可能来自文献、数据库或者AI估算数据不一定可靠。如果有一个负的C44混进去图照样能画方向也能算但物理上完全不成立。程序能自动拦截比人工一个个盯要靠谱。实测下来GPT生成的断言代码基本不需要改它会写成直接抛ValueError的形式。3.4 第三轮三维可视化与极值标注图形部分是这类任务的重头戏。球坐标采样时我习惯用规则网格θ从0到πφ从0到2π各取几十到上百个点得到一个经纬网格。这样用plot_surface画出来的曲面是规整的颜色映射也连续。核心绘图代码大致如下theta np.linspace(0, np.pi, 100) phi np.linspace(0, 2 * np.pi, 200) TH, PH np.meshgrid(theta, phi) Nx np.sin(TH) * np.cos(PH) Ny np.sin(TH) * np.sin(PH) Nz np.cos(TH) E np.empty_like(TH) for i in range(TH.shape[0]): for j in range(TH.shape[1]): n np.array([Nx[i, j], Ny[i, j], Nz[i, j]]) E[i, j] young_modulus_cubic(S11, S12, S44, n) X E * Nx Y E * Ny Z E * Nz fig plt.figure(figsize(10, 10)) ax fig.add_subplot(111, projection3d) norm plt.Normalize(E.min(), E.max()) surf ax.plot_surface(X, Y, Z, facecolorsplt.cm.viridis(norm(E)), rstride1, cstride1, linewidth0, antialiasedTrue)极值标注我一般建议单独做先在整个E矩阵里找最大最小值把对应方向的单位向量乘上模量得到三维坐标再用ax.scatter和ax.text标注。这一步GPT通常写得比我快但它找极值时容易搞错索引把np.argmax(E)的位置直接当成方向向量然后从坐标数组里拿错轴。我会在提示词里明确要求极值方向请在Nx、Ny、Nz数组里通过索引同步获取。4. 验证数字对不对不能靠看起来像4.1 三个方向解析值锚点我拿到任何GPT生成的弹性计算脚本第一件事永远是跑三个晶向的解析结果。以铜的参数为例我用解析式计算了三个轴向的杨氏模量结果如下晶向解析式数值GPa[100]1 / S1166.7[110]1 / (S11 - (S11 - S12 - S44/2))130.3[111]1 / (S11 - (2/3)·(S11 - S12 - S44/2))191.2这三个数字不是随便抄的。铜的多晶杨氏模量通常在110到130 GPa之间单晶在[100]方向软、约66.7 GPa在[111]方向硬、约191 GPa各向异性比大约2.9倍。如果你跑出来的[111]方向不是190附近而是100左右那几乎可以肯定公式里丢了S44/2或者方向的归一化出了问题。把这些锚点直接贴在代码注释旁边是我个人很推荐的做法。因为三个月后你再翻这个脚本不需要重新翻教科书一眼就能确认程序是否正常。4.2 拿着锚点反过来审GPT代码我做一个额外的强烈建议把锚点直接写进GPT的自检函数并用断言或误差容限来判断。例如E100 young_modulus_cubic(S11, S12, S44, np.array([1, 0, 0])) E111 young_modulus_cubic(S11, S12, S44, np.array([1, 1, 1])) assert abs(E100 - 66.7) 1.0, f[100]方向模量异常: {E100:.2f} GPa assert abs(E111 - 191.2) 1.0, f[111]方向模量异常: {E111:.2f} GPa这个方法妙在它不依赖人眼而是让GPT生成的代码在每次运行前先自证清白。如果数据源被替换成一种不满足这些锚点的材料程序会主动报错而不是悄悄给你一张错误但精致的图。这也是提示工程里把验证内建到任务里的典型思路。你不需要在对话里反复追问GPT你确定吗——你让代码自己校验。4.3 参考值之外的合理性质疑解析值之外图形本身的对称性也是一个验证手段。立方晶系具有面心立方或者体心立方对应的对称性杨氏模量曲面应该在三个坐标轴平面上呈现两两对称。具体来说[100]、[010]、[001]三个方向模量相等[110]、[011]、[101]三个方向也相等。跑完图之后我会肉眼确认曲面是否在x轴、y轴、z轴方向上保持相同的膨胀程度。如果某个轴方向的曲面明显更大或者更瘪多半是采样坐标写错了。另外还可以做一个各向异性比的自检E_max / E_min。对铜来说大约是191.2 / 66.7 ≈ 2.87。如果你算出来是一个接近1的数值说明你的采样可能只覆盖了一个小锥角或者公式被简化成了各向同性。GPT在这种地方常常不会主动提醒但你自己心里要有数。5. GPT在这个案例里反复踩的三个坑我也跟着踩过5.1 坑一把Voigt记号的逆用错这是最大的一个坑而且它藏得很深。很多人包括我自己早期都干过直接S np.linalg.inv(C)然后拿着6x6的S矩阵把S[0,1]当成S12把S[3,3]当成S44然后代入公式算方向模量。严格说这样得到的是工程剪应变约定下的柔度矩阵和四阶柔度张量之间存在剪切分量换算关系。为什么这里又没问题呢因为我们用的立方晶系解析公式里面的S11、S12、S44恰好就是工程约定下6x6柔度矩阵里的对应元素。也就是说inv(C)之后直接取S[0,0]、S[0,1]、S[3,3]是可行的。但GPT如果好心帮你把柔度矩阵展开成3x3x3x3的张量再缩并就可能在剪切分量的系数上出问题要么忘了除以4要么除错了。所以我在这类题里会明确要求除非你有百分百把握处理四阶张量的分量换算否则直接用立方晶系解析公式。这个限制条件看起来很啰嗦实际上是让GPT避开自己最容易犯错的领域。5.2 坑二把Miller指数当普通坐标但忘了归一化计算[110]方向时方向向量是(1,1,0)归一化之后才是(0.707,0.707,0)。GPT有时会在脚本里顺手写n direction / np.linalg.norm(direction)看起来很安全。但如果你在函数外部已经归一化过一次函数内部又归一化一次按理说没问题问题往往出在循环里把单位向量和坐标数组混用最后画的曲面某些方向被拉长极值点也偏了。排查技巧其实很简单全程序只保留一份方向向量生成的逻辑不要既在采样函数里归一化又在球坐标生成时手动除模。只要看到np.linalg.norm出现在多个地方我就会多留个心眼。GPT生成的代码最常见的bug不是没有归一化而是归一化了两次但其中一次是把向量分量平方后当模长用这种错误肉眼极难发现只能靠锚点数值兜底。5.3 坑三让GPT自由发挥采样方式导致曲面误导如果提示词只写在球面上均匀采样GPT可能会用随机点采样比如再利用随机方向向量然后试图做散点图或插值。问题在于随机采样画出来的点云很难看出连续曲面而如果它用规则网格但θ和φ的范围搞反了比如θ从0到2π、φ从0到π曲面就会被翻转看起来像一个上下颠倒的气球。这种图用来发朋友圈还行用来做技术判断就不可靠了。我的要求是网格采样必须用θ ∈ [0, π]、φ ∈ [0, 2π]并且生成点后把Nx、Ny、Nz显式参与计算不要内部再偷偷转置。为了让GPT老老实实画图我还会加一句不要使用随机采样除非你同时提供三维插值方案。这样它就只能在规则网格这条路上走出错概率大大降低。5.4 对付沉默的假设让GPT亮出它的假设GPT出错的根源很多时候是它默默做了假设却不告诉你。比如它可能假设你用的是多晶平均弹性常数或者假设晶格坐标系与笛卡尔坐标系完全重合或者假设单位是Pa不是GPa。这些假设在代码里可能根本没有注释但你看到结果时已经完全被带入它的节奏里。我在提示词末尾会固定加一句请用三句话说明你这个脚本里最重要的三个假设。实测下来这句话非常有用。GPT常常会写出类似我假设输入矩阵已经满足Voigt排列这样的话帮你确认它有没有在关键处犯错。这也是提示工程里一个被低估的技巧让模型输出它认为的边界条件比单纯让它解释代码更有效。6. 把临时脚本变成长期能用的材料计算小工具6.1 设计上尽量让计算与展示分离一次性的脚本跑完就扔当然没问题。但如果你和我一样经常要处理多组材料数据我会建议让GPT帮你把代码函数化。第一版脚本通常把C矩阵写死在文件里计算和画图混在一起。升级时我会提这样的需求请把代码重构为两个文件elastic_calc.py和plot_elastic.py。elastic_calc.py里放所有计算函数包括刚度矩阵求逆、杨氏模量计算、Born稳定性检查、方向自检plot_elastic.py负责读入一个JSON文件并画图。JSON格式为{name: Cu, C11: 168.4, C12: 121.4, C44: 75.4}。这样做的好处是显而易见的以后拿到一种新材料只需要往JSON文件里加一段不用改任何逻辑。我还会让GPT在elastic_calc.py里加一个if __name__ __main__的入口专门打印自检结果。这样一来命令行直接python elastic_calc.py就能确认基础数值是否正确。6.2 一次跑多个材料的对比当脚本函数化之后批量对比就顺理成章了。让GPT生成一个循环读取一个包含多种材料的JSON数组分别计算锚点模量和各向异性比并输出一张汇总表。我会额外要一列各向异性比E_max/E_min因为这是选材时非常关注的指标比值大的材料定向方向性明显适合需要各向异性的设计比值接近1的材料更接近各向同性适合通用结构件。铜的比值大约2.87钨可能只有1.0左右铝大约1.2这个范围差异很大。批量打印出来你根本不需要看三维曲面就能快速筛选目标材料。三维图则可以用来做报告展示或者排查异常材料时使用。6.3 用GPT做代码维护时的提示词模式脚本沉淀下来之后我很少自己从头写维护性改动而是继续让GPT干活。但维护性对话和从零生成代码不同我的提示词会换成另一种结构这是我现在正在使用的弹性计算脚本。请你在不改动核心公式的前提下做三件事把C矩阵读取改成从JSON加载在输出表格中增加E_min对应的晶向增加一个命令行参数--material可以指定只计算某一种材料。这种约束性修改提示词的关键是不改动核心公式这六个字。你不说这句GPT可能顺手帮你优化公式排版、改变变量名结果代码改得面目全非你还得花时间重新核对物理含义。反过来你说清楚边界它就会老老实实做局部改动。6.4 和更复杂的计算场景衔接弹性常数分析只是起点。同一套柔度矩阵还经常用来算体积模量K (C11 2*C12)/3、剪切模量G的相关组合也可以进一步算德拜温度、弹性波速度或者给有限元模拟提供材料卡片。GPT在这些环节依然是很好的助手但每走一步都需要你亲手把物理锚点铺设好。我的经验法则是每引入一个公式就要求GPT同时给出至少两个极端情况的解析结果每引入一个新材料就多跑一次已知材料做对照。这套方法论正是从这一篇的三维杨氏模量曲面案例里提炼出来的。最后说点个人实操中的体会。用GPT辅助物理编程最忌讳的心理是它跑起来了应该就对了。我在铜这个例子上踩过很多次坑后来养成的习惯是每个脚本必须至少有一个物理锚点没有锚点之前不看图。你可以把[100]方向的66.7 GPa、[111]方向的191.2 GPa打印在控制台最显眼的位置也可以把它们写成断言。对我来说这比任何代码规范都管用。以后你再让GPT帮你写任何材料计算脚本时不妨也先花两分钟想想这个任务里哪几个数字是我闭着眼都能猜出来的想清楚这个问题你的提示词会立刻上一个档次。
阅读完成 · 觉得有帮助？
