企业官网建设流程全解析

PharmaPy 干法制粒二次开发教程（3）：对象模型与 flowsheet 机制——SimulationExec、流程拓扑串与单元装配

版本声明块

• 工具/软件：CryPTSys/PharmaPy（SimulationExec执行器）
• 语言/环境：Python 3.9 + numpy<2 + scipy<1.14 + assimulo==3.4.3
• 本文目标：让你把"一个单元 = 一个类实例"再升格为"一条流程 = 一串拓扑 + 一套属性注入"，跑通官方 guide 的 R01→HOLD01→CR01→F01 四单元流程，并看清今后如何把干法制粒单元插进同一套插座。

一句话结论：PharmaPy 用SimulationExec(path_phys, flowsheet='R01 --> HOLD01 --> CR01 --> F01')的-->字符串定义流程拓扑，单元以类实例挂到sim并注入Utility/Kinetics/Inlet/Phases属性完成装配，SolveFlowsheet(kwargs_run=...)按每单元runtime/sundials_opts求解，结果经DynamicResult.time/result.字段与GetStreamTable(basis='mass'|'mole')读取。

〇、本篇要解决的认知问题

学完前两篇，你已经会搭建环境、跑通官方测试，能确认一个批式反应器能求解。但"跑单元"和"跑流程"之间，隔着一层框架的装配心智。本篇要回答四个问题：

1. PharmaPy 为什么用"对象 + 拓扑字符串"来描述流程，而不是一个写死的大模型输入文件？
2. 一个单元在sim上到底是怎么被"喂饱"的——Utility、Kinetics、Inlet、Phases四个属性分别装什么、各自的赋值时机？
3. flowsheet='R01 --> HOLD01 --> CR01 --> F01'这条串在SimulationExec里会被解释成什么？缓冲罐、结晶器、过滤器在串中的角色有何不同？
4. 求解器和结果到底长什么样——SolveFlowsheet(kwargs_run=...)怎么传每单元的求解参数，DynamicResult与GetStreamTable如何取数？

这四个问题一旦通透，篇 5 你要新增的RollerCompactor/Comminution/Classifier就只是"再造一个类 + 改一行拓扑串"的工作。

一、机制解析

为什么是"对象 + 拓扑串"

PharmaPy 的定位是"pythonic"——官方仓库 README 原话是让用户能simulate the dynamics of standalone, drug substance unit operations，也能setting up and simulating pharmaceutical flowsheets, i.e., interconnected unit operations。它的立身哲学不是"喂一个总装脚本"，而是：每个单元操作就是一个普通 Python 类，挂在你创建的仿真执行器对象上，再用一个字符串描述它们之间的连接。

这套设计有三个现实收益：

1. 增量构建。你可以先单独构造一个单元、只喂必要属性，求解看行为；满意后再把下一个单元挂上去。单元的"独立性"和"连接性"被解耦。二次开发时，我们新写的RollerCompactor可以先脱离流程单独跑，再接入拓扑串——这正是调试的最小闭环。

2. 属性注入赶走构造函数地狱。一个BatchCryst的全参数签名有十几项（target_comp、method、scale、vol_tank、controls、adiabatic、h_conv……），但真正随时间变化的边界条件（冷却水、进料流股、动力学、初始相）并不宜全塞进构造函数。PharmaPy 用sim.单元名.Utility = ...这种赋值即注入的方式，把"过程对象"和"配置属性"分开：类型是sim上挂的类，配置是后面逐项点的属性。想让单元换一种水浴，改一行赋值即可，不碰构造函数。

3. 拓扑可读、可审计。'R01 --> HOLD01 --> CR01 --> F01'一行就把链条画完，设备和物流的上下游关系一目了然，且F01这类带别名（It's Process）的单元还能在同一流程里做替代对比。对干法制粒的三段式（辊压→粉碎→分级），这样一行拓扑就是自文档。

上图中的"Feeder → Roller Compactor → Mill → Classifier → 颗粒产品"正是我们要在 PharmaPy 上用同类拓扑串FEED --> RC --> MILL --> CLASS --> OUT二次开发出来的流程——先把三段式画在拓扑上，再逐段填动力学模型（篇 6~12）。

价值锚点：本篇顺完官方 guide 的"反应器→缓冲罐→结晶器→过滤器"链条，就等于把框架默认的"插座"摸清了。干法制粒所有单元——无论新增多少——最终都插进这同一个插座：先建模成类，再sim.Xxx = 类(...)，再点属性，再进拓扑串，再求解。形式不变，变的只是动力学。

属性注入模式的价值

看官方 guide 的写法，sim.R01建好后并不是一条config字典，而是连续四行"挂属性"：

sim.R01=PlugFlowReactor(diam_in,num_discr=50,isothermal=False)sim.R01.Utility=cw# 冷媒：冷却水sim.R01.Kinetics=kinetics# 反应动力学sim.R01.Inlet=liquid_in# 进料流股（动态边界）sim.R01.Phases=liquid_init# 初始相（罐内初始状态）

这四条属性对应单元操作的四类"输入来源"，它们合起来才是单元的"完整装配"：

• Utility：公用工程，把过程参数（温度设定、流速）从物料动力学里剥离出来。干法制粒中，这正是辊压机辊温/冷却回路的位置。
• Kinetics：微观速率，决定反应或结晶动力学。干法制粒里对应压缩性指数、破碎选择函数/破碎函数——我们篇 6~11 落地的东西最终都从这里塞进去。
• Inlet：动态进料边界。
• Phases：初始相（含固相 PSD）。这是唯一与相有关的输入，也是干法制粒"数字孪生"最需要精耕的字段——SolidPhase的distrib/x_distrib就是 PBM 状态向量。

判定法则：凡是"这个单元在求解开始时需要先知道的状态"塞Phases；"求解过程中持续驱动的外部量"塞Utility/Inlet/Kinetics。理解这一点，篇 4 讲SolidPhase的矩/离散双表示时你就不会迷路。

二、完整代码与逐行剖析

下面给出官方 guide（PharmaPy_guide.html）风格的最小可运行 flowsheet：一列管式反应器（R01）、一个动态缓冲罐（HOLD01）、一个1D-FVM结晶器（CR01）、一个过滤器（F01）。path_phys指向官方 data 目录中的物性数据。API 签名逐字摘自官方 guide 与 readthedocs。

importnumpyasnpfrompharmapy.SimExecimportSimulationExecfrompharmapy.ReactorsimportPlugFlowReactorfrompharmapy.CrystallizersimportBatchCrystfrompharmapy.PhasesimportLiquidPhase,SolidPhasefrompharmapy.StreamsimportLiquidStreamfrompharmapy.KineticsimportRxnKineticsfrompharmapy.UtilitiesimportCoolingWaterfrompharmapy.SolidLiquidSepimportDynamicCollector,Filter# ---------- 0. 物性数据目录（data/ 下 JSON） ----------path_phys='path/to/PharmaPy/data'# ---------- 1. 定义流程拓扑串 ----------graph='R01 --> HOLD01 --> CR01 --> F01'# 先建执行器，再往 sim 上挂单元sim=SimulationExec(path_phys,flowsheet=graph)# ---------- 2. 通用工程 / 动力学 / 流股 / 相 ----------# 冷却水（R01 恒温外套水源）cw=CoolingWater(mass_flow=5.0,temp_in=298.15)# 一级不可逆反应动力学（A -> C）kinetics=RxnKinetics(rate_constants={'k1':1e-4})# 反应器进料：以摩尔浓度给出liquid_in=LiquidStream(path_phys,temp_in=298.15,mole_conc={'A':2.0,'C':0.0})# 反应器初始相（批充料）liquid_init=LiquidPhase(path_phys,init_temp=298.15,mass_frac={'A':1.0},vol=1e-3)# ---------- 3. 逐单元装配 ----------sim.R01=PlugFlowReactor(diam_in=0.05,num_discr=50,isothermal=False)sim.R01.Utility=cw sim.R01.Kinetics=kinetics sim.R01.Inlet=liquid_in sim.R01.Phases=liquid_init sim.HOLD01=DynamicCollector()# 缓冲罐，无需额外属性# 结晶器：1D-FVM PBM，scale=1e-9 配合 µm 级固相网格x_gr=np.geomspace(1,1500,num=35)# 对数等距 35 bin，[1, 1500] µmdistrib_init=np.zeros(35);distrib_init[10]=1e9# 一个种子 binsolid_cry=SolidPhase(path_phys,x_distrib=x_gr,distrib=distrib_init,mass_frac=[0,0,1,0,0])sim.CR01=BatchCryst(target_comp='C',method='1D-FVM',scale=1e-9,controls={'temp':298.15})sim.CR01.Phases=solid_cry# 过滤器：直径 / 阻挡系数 / 膜阻sim.F01=Filter(diam=0.1,alpha=1.0,Rm=1e12)# ---------- 4. 逐单元求解参数 ----------runtime_reactor=3600.0# sruntime_cryst=7200.0# srun_kwargs={'R01':{'runtime':runtime_reactor},'HOLD01':{'runtime':runtime_reactor},'CR01':{'runtime':runtime_cryst,'sundials_opts':{'maxh':60}},'F01':{'runtime':runtime_cryst},}sim.SolveFlowsheet(kwargs_run=run_kwargs)# ---------- 5. 取结果 ----------res_cry=sim.CR01# 该单元自带 resultplot=sim.CR01.plot_profiles()# 矩 μ0..μ3、温度、浓度/过饱和度图t=res_cry.result.time# 时间向量mu0=res_cry.result.moments[:,0]# 0 阶矩（总粒子数密度）演化stream_table=sim.GetStreamTable(basis='mass')# 全流程流股表（pandas）stream_table.head()

逐段剖析：

• 步骤 0-1：SimulationExec(path_phys, flowsheet=graph)是流程的总指挥。graph里的-->即"物流指向"，SimulationExec会在求解时按此构建连接关系。
• 步骤 2：所有"输入材料"在挂单元之前先构造好。注意CoolingWater(mass_flow=..., temp_in=...)、LiquidStream(path_phys, ...)都带temp级参数——这些是实时边界，和步骤 3 里"谁用它们"完全解耦。
• 步骤 3：装配点。sim.R01 = PlugFlowReactor(diam_in=..., num_discr=50, isothermal=False)先建"空槽"，再用Utility/Kinetics/Inlet/Phases逐项点属性填满。sim.HOLD01 = DynamicCollector()是缓冲罐，只需要挂上即可。结晶器BatchCryst(target_comp='C', method='1D-FVM', scale=1e-9, ...)，其Phases装的是SolidPhase——这个scale=1e-9与x_gr = np.geomspace(1, 1500, num=35)的 µm 级网格配合（PharmaPy 内部把 µm 基准矩换算成米基准体积分数时乘以(1e-6)**conv_exp）。sim.F01 = Filter(diam, alpha, Rm)则是管线式单元。
• 步骤 4：run_kwargs的键名 = 拓扑串里的单元名，每个值为该单元的求解参数字典，runtime是仿真时长，sundials_opts={'maxh': 60}限定了吃紧单元的最大步长（结晶这类刚性系统常需要）。
• 步骤 5：剖面型结果sim.CR01.plot_profiles()直接出图；result.time给时间向量，result.moments之类的字段按绝对值取状态；GetStreamTable(basis='mass')汇总全流程所有流股的 pandas 表格（temp/pres/mass/moles/vol/mass_frac/vol_flow，basis可换'mole'）。

关键：sim.CR01既是"执行器上的一个槽"，求解后它自己同时携带result对象。结果对象在框架里由DynamicResult.__dict__.update(**results)填充，result.time与result.<field>都能直取。所以你要取某个单元的剖面，就sim.单元名.result.字段。

三、常见报错与排查

报错 1：拓扑串格式写错——sim = SimulationExec(path_phys, flowsheet='R01-->CR01')
现象：解析期即报流程解读失败。根因：flowsheet需要把节点名用空格包好的箭头-->分隔，裸-->或不带空格往往无法正确切分出节点。解法：统一写'R01 --> HOLD01 --> CR01 --> F01'，节点名与箭头之间存在空格；并核对每个节点名与后续run_kwargs键、sim.Xxx赋值名三者完全一致（大小写敏感）。

报错 2：某单元漏注入Kinetics/Phases——sim.R01.Kinetics未赋值就求解
现象：SolveFlowsheet执行到该单元时抛 “attribute not found” 或盒模型无初值报错。根因：官方引导下RxnKinetics、初始相必须由用户显式点属性，漏一个单元在运行时才开始报告。解法：装配后将官方 guide 的四条属性逐一核对——反应器要Utility/Kinetics/Inlet/Phases四件套齐，结晶器至少要喂Phases里的SolidPhase。排查时可先只解单单元，缩小范围。

报错 3：scale/distrib维度或量级不匹配——固相网格与scale冲突
现象：矩或体积累计出现 NaN / 负质量 / 数量级离谱。根因：scale=1e-9与x_distrib的 µm 量级网格必须是配套的同一约定；SolidPhase的distrib长度必须与x_distrib的 bin 数一致（如上面都是num=35）。解法：网格np.geomspace(1, 1500, num=35)与scale成组绑定，改其中一个就要同步另一个；先跑官方结晶示例自测plot_profiles()的 μ0 是否单调不炸再改参数。

报错 4：sundials_opts缺失导致刚性系统解不动
现象：结晶/反应系统求解步长不断缩短、报积分器步长过小或收敛失败。根因：1D-FVMPBM 与多相耦合是典型刚性，默认 SUNDIALS 设置不足以追踪快速反应。解法：在run_kwargs对应单元补'sundials_opts': {'maxh': 60}（限制最大步长），必要时调小maxh或调整容差；确认重刚assimulo==3.4.3已随 conda-forge 装好。

报错 5：plot_profiles()出空图 / “no data”
现象：图框架出来了但曲线为空。根因：该单元没解到终点或result未填充（往往是因为run_kwargs漏了它的键，导致该单元被跳过）。解法：确认run_kwargs里确实给该单元传了runtime；再确认sim.单元名.result非空后再plot_profiles()。

四、动手练习

练习 1（必做）：逐字抄下 §二 的完整代码并运行（把path_phys指向你本机PharmaPy/data目录，必要时借官方reactor_tests.py里的物性文件名）。判定标准：脚本无报错，sim.CR01.plot_profiles()出图，stream_table = sim.GetStreamTable(basis='mass')返回的是一个 pandas DataFrame（可print(stream_table.columns)）。

练习 2：把basis从'mass'换成'mole'，对比GetStreamTable返回的流股表在质量项 vs 摩尔项上的差异。再在run_kwargs中给CR01加'maxh'从 60 改到 20，观察积分是否更稳。判定标准：数值量纲正确、变更前后物理趋势不变。

练习 3（进阶，预告）：新增一个纯"占位"单元sim.M01 = Xxx(...)并把它插进拓扑串末尾，但先不往run_kwargs里加它的键。运行观察错误现象，对照 §三「报错 5」。判定标准：你能说出报错与"漏装配/漏求解参数"的因果关系，并补上正确键使其通过。

五、小结与下一篇预告

本篇完成了一次认知跃迁：单元（类）+ 装配（属性注入）+ 拓扑（字符串）三件套构成 PharmaPy flowsheet 的完整装配。SimulationExec(path_phys, flowsheet=graph)定义框架，Utility/Kinetics/Inlet/Phases决定每个单元吃什么、怎么动、初始什么状态，SolveFlowsheet(kwargs_run=...)按每单元runtime/sundials_opts推进，DynamicResult.time/result.字段与GetStreamTable(basis=...)承接结果。今后所有干法制粒单元——RollerCompactor、Comminution、Classifier——都会套用这套插座，只是把Kinetics换成压缩性/破碎动力学，把Phases换成SolidPhase承载 PSD。

下一篇（篇 4固相/PBM 基础设施）将深入本篇反复出现的SolidPhase：x_distrib/distrib的离散物流 vs 矩表示（moments/num_mom/getMoments）、1D-FVM与moments双方法的选择依据。理解了 PSD 怎么在SolidPhase里被表示、被取矩，你才真正握住干法制粒数字孪生的"状态中枢"。

本篇认知问题回显（FAQ）

Q1：PharmaPy 为什么用"对象 + 拓扑字符串"，而不是一个大配置文件？
因为它要同时支持"独立单单元"和"互联 flowsheet"两种用法。每个单元操作是一个普通类（可脱离流程单独构造调试）；流程连接用一个可读的-->字符串描述。增量构建、属性注入、拓扑可审计三方面的收益，让二次开发时"先单测单元、再进流程"成为可能。

Q2：Utility、Kinetics、Inlet、Phases四个属性各自装什么、何时赋值？
Utility装公用工程（如CoolingWater，过程设定温度/流速）；Kinetics装微观速率（反应动力学，干法处对应压缩性/破碎动力学）；Inlet装动态进料流股；Phases装初始相（含SolidPhase的 PSD）。赋值时机在单元类挂到sim之后、SolveFlowsheet之前，用"赋值即注入"完成装配，与构造函数解耦。

Q3：flowsheet='R01 --> HOLD01 --> CR01 --> F01'会被怎么解释？缓冲罐、结晶器、过滤器角色有何不同？
SimulationExec读字符串，把-->理解为物流指向，据此构建上下游连接。HOLD01 = DynamicCollector()只需挂上、角色存中继；CR01 = BatchCryst(...)是执行 PBM 的核心单元，需喂Phases；F01 = Filter(diam, alpha, Rm)是管线式分离单元。三者都在拓扑串中占一个节点，但所需装配的复杂度不同。

Q4：SolveFlowsheet(kwargs_run=...)怎么传参？DynamicResult与GetStreamTable如何取数？
BatchCryst(…)是执行 PBM 的核心单元，需喂Phases；F01 = Filter(diam, alpha, Rm)` 是管线式分离单元。三者都在拓扑串中占一个节点，但所需装配的复杂度不同。

Q4：SolveFlowsheet(kwargs_run=...)怎么传参？DynamicResult与GetStreamTable如何取数？
run_kwargs的键 = 拓扑串中的单元名，每值是一个求解参数字典，常用runtime（时长）与sundials_opts={'maxh': 60}（求解器限步）。结果经sim.单元名.result.time取时间向量、result.字段取状态量；sim.CR01.plot_profiles()出剖面图；sim.GetStreamTable(basis='mass'|'mole')返回全流程流股的 pandas 表格。