虚拟器官插件开发教程7从细胞到组织——传导、折返与 svFSI 的 EP 算例版本声明块工具/软件svFSI / svMultiPhysicsBSD 系许可EP 内建openCARP v19.02026-04本体 Academic Public License v1.1Chaste 2026.1BSD-3-Clause语言/环境Python 3.11 numpy行波演示Fortran/C 求解器官方算例本文目标把第 6 篇的单细胞 qNet推进到组织尺度——理解传导方程、摸清 svFSI 内建电生理的能与不能、建立波长-折返的基质判据一句话结论单域方程dV/dt D·∂²V/∂x² − I_ion I_stim的最小 numpy 实现即可复现传导行波实测 D 每翻倍 CV 约乘 √2本文实跑 71.5→106.7→154.3→219.9 cm/ssvFSI/svMultiPhysics 内建cepModel_AP/FN/BO/TTP电生理注册表C 类CepModTtp含I_Na、I_bNa官方svFSI-Tests的08-cep/03-benchmark_tTP对应 Niederer 2011 N-version benchmarkPhil Trans R Soc A 369:43313×7×20 mm 组织块openCARP 本体是非商业Academic PL v1.1商用找 NumeriCor而Chaste 是 GPL属于过时误传2026.1 为 BSD-3。〇、本篇要解决的认知问题单域monodomain与双域bidomain模型到底差在哪插件开发者该选哪个SimVascular 的 svFSI/svMultiPhysics 究竟有没有心脏电生理内建哪些模型、官方给了哪些算例openCARP、carputils、meshalyzer、Chaste 各自的许可证是什么商用插件哪些能直接链接Niederer 2011 N-version benchmark 对电生理插件验证有什么用动作电位波速怎么变成传导速度→激动波长→折返基质的心律失常风险语言一、机制解析1.1 从细胞到组织只差一项扩散语义全变第 6 篇里一个细胞的状态是 ODE 系统铺成组织后膜电位多了空间耦合成为反应-扩散reaction-diffusion方程单域模型monodomain只解一个跨膜电位场 V∂V/∂t D∇²V − I_ion/C_m I_stim。各向异性写进张量 D。计算量小适合插件迭代与批量筛选。双域模型bidomain胞内、胞外两个电位场 coupled能给出胞外电位/伪 ECG 与刺激阈值场。文献与代码中统称单域/双域中文首次出现给英文原文即 monodomain/bidomain。两者的工程分水岭要不要细胞外空间。做传导速度与折返机制单域够用要往第 8 篇的伪 ECG 或临床电极对比走才需要双域或单域离线积分的近似。尺度链插件视角 通道阻滞(第5篇) → 单细胞AP/qNet(第6篇) → 组织传导/波长(本篇) → 器官ECG(第8篇) Hill/Markov ODE 积分 反应-扩散 PDE 容积导体叠加1.2 波速为什么是风险语言行波速度即传导速度conduction velocityCV。两个直接推论激动波长 λ CV × ERP有效不应期一个波长内组织来不及恢复兴奋。折返基质折返环周长 ≥ λ 时波前追上的组织已脱不应期折返可持续——这是室速/室颤维持的经典判据电生理标测EP study中的程序刺激EPS本质就是在测不应期与传导判断有没有这种可被点燃的环路。药物插件在这里的输出语义就变了第 6 篇给APD90/qNet 延长了本篇要回答传导慢了没、不应期变了没、λ 小到能让哪块组织折返。智源虚拟生理心脏官方叙事中的组织传导速度/波长/兴奋易损区/折返波层级见第 8 篇讨论其器官层说的就是这一段。1.3 工具箱三条组织级仿真路线与许可证分层铁律 3工具定位版本基线许可证商用插件可用性svFSI / svMultiPhysicsSimVascular 体系多物理求解器EP 内建svMultiPhysics 2026-09 仍活跃原名 svFSIplus 已并入2024-11-20 官方改名提交svMultiPhysics/svZeroDSolverBSD-3svFSI “MIT-like”可BSD/MIT 系openCARP离子通道→器官组织级专用v19.02026-04-02支持 CellML 导入carputils 编排本体 Academic Public License v1.1非商业carputilsApache-2.0meshalyzerGPL-3.0商用必须向 NumeriCor 购买授权子组件分层核对ChasteC 多尺度框架heart/mesh/pde/linalg…CellML→C 代码生成2026.13-clause BSD网上Chaste 是 GPL为过时说法勿信旧页面可openCARP 论文口径Comput Methods Programs Biomed 2021;208:106223doi:10.1016/j.cmpb.2021.106223。它继承 CARP/acCELLerate 脉络组织-器官级功能最全——但立项前先过许可证矩阵把 openCARP 当BSD 同名软件商用是本系列见过的高危误判。一个诚实的替代姿势把 openCARP 用于研发期对照商用交付切到 svMultiPhysics/Chaste 路线插件层用统一契约隔离后端。1.4 svFSI/svMultiPhysics 内建 EP有什么、缺什么官方口径先纠正一个流传很广的错误“SimVascular 只管血流”——不成立。svFSI 简介原文含 “blood flow simulation including fluid-structure interaction and cardiac electrophysiology”svMultiPhysics 覆盖 “solid and fluid mechanics, diffusion, and electrophysiology … whole heart dynamics”。内建电生理模型注册表svFSI README Features 表Doxygen 可查实现类注册名模型实现类性质cepModel_APAliev-PanfilovCepModAp简化两变量cepModel_FNFitzhugh-NagumoCepModFn概念模型cepModel_BOBueno-Orovio-Cherry-Fenton—螺旋波常用cepModel_TTPten Tusscher-Panfilov 全离子模型CepModTtp成员含I_Na、I_bNa人室肌三型endo/M/epi 参数在 TP04/TP06 谱系TP06 出处 Am J Physiol Heart Circ Physiol 291(5):H2396–H2411doi:10.1152/ajpheart.00109.2006。注意CepModTtp已经显式拆出I_Na/I_bNa晚钠——这给第 9 节药物阻滞层自建留了接口把浓度→阻滞率映射接到电流项缩放就是插件要写的代码铁律 4官方没有 hERG 阻滞→I_Kr 缩放的药理学模块GUI 也没有 EP 面板EP 全部在求解器输入文件层配置。官方算例库svFSI-TestsEP 相关部分算例路径内容08-cep/01-2Dsqr_AP2D 方块 Aliev-Panfilov 传播08-cep/03-benchmark_tTPNiederer et al. 2011 N-version benchmarkPhil Trans R Soc A 369:43313×7×20 mm 立方块模板svFSI_master.inp08-cep/04-2Dspiral_BOBO 模型螺旋波08-cep/05-Purkinje浦肯野-心肌偶联传导06-ustruct/03-LV-Guccione-active兴奋-收缩耦合EC couplingGuccione 本构 active 收缩衔接第 10 篇svFSI 的 JOSS 论文2022明确口径 “simulating the complex excitation-contraction coupling…”。N-version benchmark 的价值多个独立团队在同一几何/同一模型上的互测基准——你的 EP 求解器接进插件后先复跑它对齐官方数值再谈药物场景第 15 篇验证方法学会展开。浦肯野刺激通路用 SimTK 官方分发的Purkinje Plugin跨平台安装包 2022-07-20“used to create a Purkinje network on a surface model of the heart”——SimVascular 体系里唯一官方电生理刺激插件患者特异全管线看 SimCardiodocsSimCardio.html。二、完整代码与逐行剖析2.1 工程动作定位官方算例并跑通bash# 获取算例库与求解器命令形态按各仓库 README以官方文档为准gitclone https://github.com/SimVascular/svFSI-Tests# 官方算例库gitclone https://github.com/SimVascular/svFSIcdsvFSI# svFSI 为 Fortran 多物理求解器官方调用形态调研核实mpiexec-npNbuild/svFSI-build/bin/svFSIinput.inp# 输入 .inp模板为 svFSI_master.inpcd../svFSI-Tests/08-cep/03-benchmark_tTP# Niederer 2011 N-version 基准# 该算例几何为 3×7×20 mm 组织块——正是 1.2 节环路周长 vs 波长判据能直接套用的尺寸要点svFSI 的 EP 场景全部在.inp输入文件层配置GUI 不面板化svMultiPhysics现行默认后端走 XML 输入且无原生 Windows需 WSL。2.2 svFSI.inp 关键段解读示意结构字段以官方仓库为准# ---- svFSI_master.inp 电生理相关关键段教学示意逐段语义解读 ---- [EP] cepModel cepModel_TTP ; 内建注册表名AP/FN/BO/TTP 四选一 ; 选 TTP ten Tusscher-Panfilov 全离子模型CepModTtp ; 药物插件的挂载点把 Hill 阻滞率接到 I_Na/I_bNa/IKr 电流缩放 ;官方无药理学层这段接到的代码就是第 17 篇 cardiotox 插件本体 cellType 2 ; 透壁细胞类型 endo/M/epi 选档 → 第 8 篇伪ECG的梯度来源 stimNodesFile stim.dat ; 起搏点文件浦肯野方案改用 Purkinje Plugin 生成的网络文件 [Diffusion] D tensor 1 ; 各向异性扩散张量需要纤维方向数据无纤维各向同性退化 ; 纤维方向缺失时 CV 会被系统性高估——传导结论先过此检查读这段的正确姿势不求逐字背诵记住三件事——模型名在注册表里选、起搏靠输入文件、各向异性靠纤维数据。2.3 1D 单域行波CV 标度 ERP 波长Python已实跑# -*- coding: utf-8 -*-第7篇实跑验证1D 单域(monodomain)反应-扩散行波 CV 标度 ERP 波长/折返基质 反应项 第6篇同款最小平台期细胞玩具演示模型非 TP04/TP06/ORd 本体 真实传导速度以 svFSI-Tests 08-cep/03-benchmark_tTPNiederer 2011 3×7×20mm等官方算例为准。 importnumpyasnp E_NA,E_K,E_CA60.0,-85.0,60.0defxinf(V,Vh,k):return1.0/(1.0np.exp(np.clip(-(V-Vh)/k,-50,50)))# 激活门defhinf(V,Vh,k):return1.0/(1.0np.exp(np.clip((V-Vh)/k,-50,50)))# 失活门defionic(V,h,s,f,w,gNa4.0,gKr0.06,gCa0.05,gNal0.002,gK10.09):mxinf(V,-40.,4.)f11.0/(1.0np.exp(np.clip((V40.)/10.,-50,50)))# IK1 内向整流因子return(gNa*m*h*(V-E_NA)gCa*s*f*(V-E_CA)gNal*(V-E_NA)gKr*w*(V-E_K)gK1*f1*(V-E_K))# 合电流外向为正defgate_step(V,h,s,f,w,dt,tau_h2.,tau_s30.,tau_f150.,tau_w80.):hdt*((hinf(V,-60.,3.)-h)/tau_h)sdt*((xinf(V,-25.,7.)-s)/tau_s)fdt*((hinf(V,-45.,5.)-f)/tau_f)wdt*((xinf(V,-15.,8.)-w)/tau_w)returnh,s,f,wdefmonodomain_1d(D0.002,L8.0,dx0.01,T120.0,dtNone):dV/dt D·∂²V/∂x² − I_ion I_stim显式 Euler。 dt 自动满足扩散稳定条件 dt ≤ 0.4·dx²/D——单域方程的第一号数值坑。ifdtisNone:dtmin(0.05,0.4*dx*dx/D)nxint(L/dx);ntint(T/dt)Vnp.full(nx,-85.);hhinf(V,-60.,3.);sxinf(V,-25.,7.)fhinf(V,-45.,5.);wxinf(V,-15.,8.)t_upnp.full(nx,np.nan);xnp.arange(nx)*dxforitinrange(nt):tit*dt lapnp.zeros(nx)lap[1:-1](V[2:]-2*V[1:-1]V[:-2])/dx**2# 中心差分两端保持 0Neumann 零流Inp.where(xdx*2,40.,0.)ift1.0elsenp.zeros(nx)# 左端 1ms 局域刺激Vdt*(D*lap-ionic(V,h,s,f,w)I)# 扩散项反应项行波的全部来源h,s,f,wgate_step(V,h,s,f,w,dt)freshnp.isnan(t_up)(V-10.)# 阈值上时记录各点激活时刻t_up[fresh]t finnp.isfinite(t_up);xs,tsx[fin],t_up[fin]keepxsdx*20# 去掉刺激附近非线性启动段pnp.polyfit(xs[keep],ts[keep],1)# t_up≈x/CVb → 斜率倒数CVresidnp.max(np.abs(ts[keep]-np.polyval(p,xs[keep])))return1.0/p[0],resid,int(keep.sum()),dtprint( CV 随扩散系数细胞间耦合强度/纤维化程度代理的标度 )prevNoneforDin[0.001,0.002,0.004,0.008]:cv,resid,n,dtmonodomain_1d(DD)linef D{D:g}dt{dt:.3g}→ CV{cv*1000:6.1f}cm/s (拟合残差{resid:.2f}ms, n{n})ifprev:linef CV 比{cv/prev:.2f}vs √2{np.sqrt(2):.2f}prevcv;print(line)deferp_cardiac_like(dt0.01,SI_listNone):向量化 S1-S2所有间期一次积分。返回能再次激发的最短 S1-S2 间期。ifSI_listisNone:SI_listnp.arange(40,520,10).astype(float)Nlen(SI_list)Vnp.full(N,-85.);hhinf(V,-60.,3.);sxinf(V,-25.,7.)fhinf(V,-45.,5.);wxinf(V,-15.,8.)firednp.zeros(N,bool);NTint((SI_list[-1]200)/dt)foritinrange(NT):tit*dt Inp.where(t1.,40.,0.)np.where((tSI_list)(tSI_list1.),40.,0.)Vdt*(-ionic(V,h,s,f,w)I)h,s,f,wgate_step(V,h,s,f,w,dt)fired|(tSI_list1.)~fired(V-30.)# S2 后出现可传播兴奋已脱敏idxnp.where(fired)[0]returnfloat(SI_list[idx[0]])ifidx.sizeelsefloat(nan)erperp_cardiac_like()cv,*_monodomain_1d(D0.002)lam(cv*10)*erp# cm/ms→mm/ms(×10)再×ERP(ms)波长(mm)print(f\nERP(再激发最短间期) ≈{erp:.0f}msD0.002 时 CV{cv*1000:.0f}cm/s)print(f激动波长 λ CV×ERP {lam:.1f}mm)print( 折返基质判据环周长 ≥ 波长 → 折返可维持以 3×7×20 mm 组织块为例 )for(a,b)in[(3,7),(7,20),(3,20)]:loop2*(ab)print(f{a}×{b}mm 环路 周长{loop:3d}mm vs λ{lam:.1f}mm → (可维持折返iflooplamelse折返自灭波前撞上不应期波尾))print( 敏感性ERP 缩短缺血/晚钠抑制类比让小块组织也能折返 )forkin[1.0,0.5,0.25]:l2lam*kprint(f ERP×{k:4.2f}→ λ{l2:5.1f}mm : 3×7环(20mm) (可折返if20l2else不能)f; 7×20环(54mm) (可折返if54l2else不能))实跑输出Python 3.10 numpy 2.2.6 复验一致D0.001 dt0.04 → CV 71.5 cm/s (拟合残差 0.26 ms, n777) D0.002 dt0.02 → CV 106.7 cm/s (拟合残差 0.24 ms, n776) CV 比1.49 vs √21.41 D0.004 dt0.01 → CV 154.3 cm/s (拟合残差 0.20 ms, n776) CV 比1.45 vs √21.41 D0.008 dt0.005 → CV 219.9 cm/s (拟合残差 0.14 ms, n774) CV 比1.42 vs √21.41 ERP(再激发最短间期) ≈ 40 msD0.002 时 CV107 cm/s 激动波长 λ CV×ERP 42.7 mm 3×7 mm 环路 周长 20 mm vs λ42.7 mm → 折返自灭 7×20 mm 环路 周长 54 mm vs λ42.7 mm → 可维持折返 ERP×0.25 → λ10.7 mm : 3×7环(20mm) 可折返逐行要害dt min(0.05, 0.4*dx*dx/D)显式 Euler 解扩散的稳定性条件是dt ≤ dx²/(2D)留 0.4 安全系数。CV 比接近理论标度 √2 本身就说明网格/步长进入了波速已收敛区——这行是插件里最便宜的自检。t_up记录每点首次越过 −10 mV的时刻表polyfit一次线性拟合出 1/CV残差即行波质量指标本文 0.14–0.26 ms残差大波还没进入稳定传播段或被边界反射污染。ERP 用向量化 S1-S2 程序刺激所有间期一次积分与电生理实验室的期前刺激协议同构——把实验范式写成代码比背公式更重要。玩具参数下 CV 数量级几十~两百 cm/s不代表人心室肌真实值把它当CV 对 D 的标度关系与测量方法的验证。真实定量结论交给cepModel_TTP官方算例传导各向异性还需要纤维方向数据缺失时结果系统性偏快——svFSI 官方边界。相似 API 对比svFSI 的 CEP 是求解器内建、输入文件配置openCARP 的对应物是.inp carputils编排并可经CellML 导入外部细胞模型许可证提醒openCARP 本体非商业Chaste 则是代码生成路线CellML→Cchaste_codegenMyokit 能导出 C/CUDA/OpenCL 代码嵌入任意 PDE 框架——插件层选择多契约层不变进来的是阻滞率数组出去的是 APD/CV/λ。三、常见报错与排查ValueError: The truth value of an array with more than one element is ambiguous。现象第 4 篇能跑的标量 HH 速率函数搬进 1D 网格立刻崩。根因def f(u): return 10.0 if abs(u)1e-7 else ...是标量写法numpy 数组不能进if。解法np.where向量化奇点替换本文调试实录改为np.where(np.abs(u)1e-6, 10.0, u/(1.0-np.exp(-u/10.0)))并用np.errstate屏蔽 0/0 警告。全域瞬间一起兴奋CV 测出来是天文数字。现象激活时刻表平坦、拟合斜率≈0、CV≈1e15 cm/s。根因两种扩散稳定条件被违反导致数值爆掉后被截断或反应项激活门瞬时化强耦合让传播退化为非再生性扩散。解法先dt ≤ 0.4dx²/D再让激活门有有限时间常数本文最终用与第 6 篇同源的最小平台期细胞作反应项HH 瞬时 m 门不适合粗网格。RuntimeWarning: overflow encountered in exp/invalid value in add随后 NaN 一片。根因exp参数没有 clipV 一旦过冲 ±500 mV 即溢出污染全场。解法所有 Boltzmann 写np.clip(..., -50, 50)批量作业加 NaN 自检行铁律 8/10。“Chaste 是 GPL 所以商用要慎”。这是旧网页/旧课件的过时说法Chaste 2026.1 为 3-clause BSD可商用链接。反例要小心FDA/CiPA 的 R 仓库确是GPL-3.0闭源分发插件须做进程级隔离第 6 篇细讲——别把两件事混着记。把 openCARP 直接接进商用插件。许可证是Academic Public License v1.1非商业商用须向 NumeriCor 购买其子组件再分层carputils 是 Apache-2.0、meshalyzer 是 GPL-3.0——链接策略逐个核铁律 3。四、动手练习把D0.002改0.008确认 CV 比在 1.4±0.1判定输出两行 CV 相除再把dx加密一倍、D 不变确认 CV 变化 10%网格收敛检查。令SI_list步长从 10 ms 改 2 ms观察 ERP 估值的抖动是否 5 ms判定两次输出差值解释为什么 S1-S2 协议是测出来的不应期而不是模型参数。用官方svFSI-Tests复跑08-cep/03-benchmark_tTPWSL 或 Linux与你所在小组的插件 CV/APD 输出按 Niederer N-version 方式列表对比判定产出一张三行对照表即可交付。五、小结与下一篇预告本篇把第 6 篇的单细胞输出推到组织单域方程 最小平台期细胞就能实测 CV∝√D 与 λCV×ERP并套在 3×7×20 mm 的 benchmark 几何上判折返基质svFSI/svMultiPhysics 的cepModel_TTP/CepModTtp与官方算例库是插件的现成宿主openCARP/Chaste 各有分工但许可证先核再用。细胞层参数APD/ERP来自第 6 篇的 qNet 体系组织层的 ERP/CV 又将成为第 8 篇器官级伪 ECG 的输入。下一篇08跨到器官尺度三透壁 AP 组装伪 ECG提取 J-Tpeak 与 Tpeak-Tend接通 CiPA 第四工作流的风险评分。本篇认知问题回显FAQQ1心脏组织仿真的单域模型和双域模型区别是什么插件该选哪个A单域monodomain只解跨膜电位 V方程∂V/∂tD∇²V−I_ion/C_mI_stim双域bidomain解胞内胞外两场可出胞外电位/伪ECG。做传导速度、波长、折返筛机制选单域计算省一半以上要与真实体表/电极信号比对才上双域。Q2svFSI 和 svMultiPhysics 内建了哪些心脏电生理模型官方算例在哪里A模型注册表为cepModel_APAliev-Panfilov、cepModel_FNFitzhugh-Nagumo、cepModel_BOBueno-Orovio、cepModel_TTPten Tusscher-Panfilov 全离子类CepModTtp含I_Na/I_bNa。算例在 GitHub 组织仓库svFSI-Tests08-cep/01-2Dsqr_AP、03-benchmark_tTP、04-2Dspiral_BO、05-Purkinje刺激网络用 SimTK 官方 Purkinje Plugin。Q3openCARP 和 Chaste 的许可证分别是什么商用插件能直接用吗AopenCARP 本体是 Academic Public License v1.1非商业商用向 NumeriCor 购买carputils 为 Apache-2.0、meshalyzer 为 GPL-3.0需分层核查Chaste 2026.1 为 3-clause BSD 可商用Chaste 是 GPL是过时误传。svMultiPhysics/svZeroDSolver 为 BSD-3、svMorph 为 MIT。Q4Niederer 2011 N-version benchmark 在电生理插件验证中的作用是什么A它是多团队在相同 3×7×20 mm 几何、相同 TP06 类模型上的互测基准svFSI-Tests 提供08-cep/03-benchmark_tTP算例模板svFSI_master.inp。插件接入新求解器后端后先复跑该算例对齐参考数值再进药物场景能把实现错误与科学分歧分开属于模型验证的第一道闸。Q5传导速度和有效不应期怎么组合成折返风险判据A激动波长 λCV×ERP折返环周长 ≥ λ 时波前遇到的是已脱不应期的组织折返可维持。例本文实跑 CV107 cm/s、ERP40 ms → λ42.7 mm3×7 mm 小环20 mm不足以折返7×20 mm 环54 mm可以ERP 缩到 1/4 时小环也变危险——这就是插件在组织层要产出的风险语言。