feat(components): 新增蒸汽盘管阀门模型及优化轮转设备计算逻辑

新增 `SteamCoilVal` 类用于建模蒸汽盘管阀门，支持根据输入参数计算开度、
流量和_sigma值。同时为 `WheelS2V2`、`WheelS3` 和 `WheelS3V3` 增加了新参数
`beta_P6`、`beta_P7` 等以增强建模精度，并更新了相关传热传质计算公式。

此外，将 `BaseComponents.get_func_by_engine` 改为类方法以提升可维护性，
并修正数据清洗采样频率由15分钟改为60分钟，以及模型训练中对缺失列的处理方式。

apps/DHU/train.py

@@ -74,6 +74,10 @@ def train(*inputs,config=None):
                 DHU_type   = equp_type,
                 exist_Fa_H = config_reader.get_equp_info(each_eaup_name,'存在回风口','bool'),
                 exist_Fa_B = config_reader.get_equp_info(each_eaup_name,'存在补风口','bool'),
+                other_info={
+                    'heatingcoil_1_Fs_rated': config_reader.get_equp_info(each_eaup_name,'前蒸汽盘管额定流量','float'),
+                    'heatingcoil_2_Fs_rated': config_reader.get_equp_info(each_eaup_name,'后蒸汽盘管额定流量','float'),
+                }
             )
             
             # 清洗数据
@@ -85,7 +89,7 @@ def train(*inputs,config=None):
                 fill_zero     = False,
                 save_log      = f'{config_reader_path}/data/train/clean_log/{each_eaup_name}.txt',
             )
-            equp_data = equp_data.resample('15min').mean().dropna()
+            equp_data = equp_data.resample('60min').mean().dropna()
             save_data(f'{config_reader_path}/data/train/data_his_clean',f'{each_eaup_name}.pkl',equp_data)
             
             if not config_reader.get_app_info(each_eaup_name,'模型训练','训练模型','bool'):

components/_base_components.py

@@ -22,7 +22,8 @@ class BaseComponents(BaseModel):
         super().__init__()
         self.name = name
     
-    def get_func_by_engine(engine:str) -> dict:
+    @classmethod
+    def get_func_by_engine(cls,engine:str) -> dict:
         if engine == 'pymc':
             EXP   = pm.math.exp
             WHERE = pm.math.switch

components/coil_steam.py

@@ -90,3 +90,85 @@ class SteamCoil(BaseComponents):
 
 def get_root(b0,b1,FP):
     return (np.sqrt(b0**2+4*b1*FP)-b0)/(2*b1)
+
+
+class SteamCoilVal(BaseComponents):
+    def __init__(self, name,Fs_rated):
+        super().__init__(name)
+        self.Fs_rated = Fs_rated
+        if self.Fs_rated is None:
+            raise Exception('Fs_rated must be provided')
+        
+    def model(
+        self,
+        TinA,ToutA,FA,
+        param,
+        engine
+    ):
+        FUNC = self.get_func_by_engine(engine)
+        EXP  = FUNC['EXP']
+        Val = smooth_exp_transition(
+            x          = (ToutA - TinA) / 100,
+            s          = param['s'] - param['s_offset'] * TinA / 100,
+            a          = param['a'],
+            b          = param['b'],
+            d          = param['d'],
+            FUNC_where = FUNC['WHERE']
+        ) * 100
+        sigma = EXP((ToutA - TinA) / 100 * param['alpha2'] + param['alpha1'])
+        Fs    = self.get_Fs_by_Val_and_Tin(Val,TinA,engine=engine)
+        return {'Val':Val,'sigma':sigma,'Fs':Fs}
+    
+    def prior(self):
+        param = {
+            's'       : pm.Normal(f'{self.name}_s', mu=0, sigma=10,initval=0),
+            'a'       : pm.HalfNormal(f'{self.name}_a', sigma=10,initval=1),
+            'b'       : pm.Normal(f'{self.name}_b', mu=0, sigma=10,initval=0),
+            'd'       : pm.HalfNormal(f'{self.name}_d', sigma=10,initval=1),
+            's_offset': pm.HalfNormal(f'{self.name}_s_offset', sigma=1,initval=0.1),
+            'alpha1'  : pm.Normal(f'{self.name}_alpha1', mu=0, sigma=10,initval=-5),
+            'alpha2'  : pm.HalfNormal(f'{self.name}_alpha2', sigma=10,initval=5)
+        }
+        return param
+
+    def get_Fs_by_Val_and_Tin(self,Val,Tin,engine):
+        WHERE = self.get_func_by_engine(engine)['WHERE']
+        
+        def piecewise_function(x, a, b, s, k):
+            y_max = a * s + b + a / k
+            y = WHERE(
+                x < s, 
+                a * x + b, 
+                y_max - (a / k) * np.exp(-k * (x - s))
+            )
+            return y
+        
+        Fs = piecewise_function(
+            Val / 100,
+            a = 4.48157826,
+            b = -0.62511713,
+            s = 0.47809967 + Tin / 100 * (-0.37605304),
+            k = 16.03141911 + Tin / 100 * 1.74054886
+        )
+        return Fs * self.Fs_rated
+
+def smooth_exp_transition(x, s, a, b, d,FUNC_where):
+    """
+    光滑连接的线性-指数分段函数
+    
+    参数：
+    - x: 输入变量
+    - s: 分段点
+    - a: 线性部分斜率
+    - b: 线性部分截距
+    - d: 指数增长率
+    """
+    # 计算常数项
+    c = a / (d * np.exp(d * s))
+    e = a * s + b - a / d
+    y = FUNC_where(
+        x <= s,
+        a * x + b,
+        c * np.exp(d * x) + e
+    )
+    return y

components/wheel2.py

@@ -88,16 +88,25 @@ class WheelS2V2(BaseComponents):
         beta_P3 = param['beta_P3']
         beta_P4 = param['beta_P4']
         beta_P5 = param['beta_P5']
+        beta_P6 = param['beta_P6']
         
         RinP = get_RH_from_Tdb_and_Hr(TinP,HinP,engine)
         
         # 处理侧
-        HdiffP  = (beta_P1 * RinP**beta_P4 * HinP * TinR * EXP(-beta_P5 * FP) + beta_P2)/1000
+        HdiffP  = (
+            beta_P1 * RinP**beta_P4 
+            * HinP 
+            * (
+                TinR * (1 - EXP(-beta_P5 * FR)) 
+                + beta_P6 * TinP * (1 - EXP(-beta_P5 * FP))
+            )
+            + beta_P2
+        ) / 1000
         WdiffP  = HdiffP * FP
         HoutP   = HinP - HdiffP
         DoutP   = get_Dew_from_HumRatio(HoutP,engine)
         Q_lat_P = WdiffP * cls.CONSTANT['h_ads']
-        Q_sen_P = beta_P3 * (TinR - TinP) * FP #TODO
+        Q_sen_P = beta_P3 * (TinR - TinP) * FP
         TdiffP  = (Q_lat_P + Q_sen_P) / (FP * cls.CONSTANT['c_p_air'])
         ToutP   = TinP + TdiffP
         
@@ -120,6 +129,7 @@ class WheelS2V2(BaseComponents):
             'beta_P3': pm.TruncatedNormal(f'{self.name}_beta_P3',mu=1,sigma=2,initval=1.5,lower=0),
             'beta_P4': pm.TruncatedNormal(f'{self.name}_beta_P4',mu=1,sigma=0.3,initval=1,lower=0),
             'beta_P5': pm.TruncatedNormal(f'{self.name}_beta_P5',mu=5,sigma=2,initval=5,lower=0),
+            'beta_P6': pm.Normal(f'{self.name}_beta_P6',mu=0,sigma=1,initval=0),
         }
         return param
 

components/wheel3.py

@@ -101,6 +101,95 @@ class WheelS3(BaseComponents):
         }
         return param
 
+class WheelS3V2(BaseComponents):
+    # 适用于不能够准确测量预冷侧入口的空气温湿度
+    def __init__(self, name):
+        super().__init__(name)
+    
+    @classmethod
+    def model(
+        cls,
+        TinP,HinP,FP,
+        TinR,HinR,FR,
+        TinC,HinC,FC,
+        engine  : str,
+        param   : dict,
+    ):
+        FUNC = cls.get_func_by_engine(engine)
+        EXP  = FUNC['EXP']
+        
+        beta_P1 = param['beta_P1']
+        beta_P2 = param['beta_P2']
+        beta_P3 = param['beta_P3']
+        beta_P4 = param['beta_P4']
+        beta_P5 = param['beta_P5']
+        beta_P6 = param['beta_P6']
+        beta_P7 = param['beta_P7']
+        beta_C1 = param['beta_C1']
+        beta_C2 = param['beta_C2']
+        beta_C3 = param['beta_C3']
+        beta_C4 = param['beta_C4']
+        
+        RinP = get_RH_from_Tdb_and_Hr(TinP,HinP,engine)
+        # RinC = get_RH_from_Tdb_and_Hr(TinC,HinC,engine)
+        
+        # 处理侧
+        HdiffP  = (
+            beta_P1 * RinP ** beta_P4  
+            * HinP
+            * (
+                TinR * (1 - EXP(-beta_P5 * FR)) 
+                + beta_P6 * TinC * (1 - EXP(-beta_P5 * FC))
+                + beta_P7 * TinP * (1 - EXP(-beta_P5 * FP))
+            ) 
+            + beta_P2 
+        ) / 1000
+        WdiffP  = HdiffP * FP
+        HoutP   = HinP - HdiffP
+        DoutP   = get_Dew_from_HumRatio(HoutP,engine)
+        Q_lat_P = WdiffP * cls.CONSTANT['h_ads']
+        Q_sen_P = beta_P3 * (TinR - TinP) * FP
+        TdiffP  = (Q_lat_P + Q_sen_P) / (FP * cls.CONSTANT['c_p_air'])
+        ToutP   = TinP + TdiffP
+        
+        # 冷却侧
+        TdiffC    = beta_C1 * EXP(-beta_C2 * EXP(-beta_C3 * (TinR - TinC))) * EXP(-beta_C4 * FC)
+        ToutC     = TinC + TdiffC
+        HdiffC    = 0
+        WdiffC    = HdiffC * FC
+        HoutC     = HinC - HdiffC
+        DoutC     = get_Dew_from_HumRatio(HoutC,engine)
+        Q_total_C = TdiffC * FC * cls.CONSTANT['c_p_air']
+        
+        # 再生侧
+        WdiffR    = WdiffP + WdiffC
+        HoutR     = (HinR * FR + WdiffR) / FR
+        DoutR     = get_Dew_from_HumRatio(HoutR,engine)
+        Q_total_R = Q_lat_P + Q_sen_P + Q_total_C
+        TdiffR    = Q_total_R / (FR * cls.CONSTANT['c_p_air'])
+        ToutR     = TinR - TdiffR
+        
+        return {
+            'ToutP':ToutP,'HoutP':HoutP,'DoutP':DoutP,'FP':FP,
+            'ToutR':ToutR,'HoutR':HoutR,'DoutR':DoutR,'FR':FR,
+            'ToutC':ToutC,'HoutC':HoutC,'DoutC':DoutC,'FC':FC,
+        } 
+        
+    def prior(self):
+        param = {
+            'beta_P1': pm.TruncatedNormal(f'{self.name}_beta_P1',mu=5,sigma=10,initval=5,lower=0),
+            'beta_P2': pm.Normal(f'{self.name}_beta_P2',sigma=1,initval=0),
+            'beta_P3': pm.TruncatedNormal(f'{self.name}_beta_P3',mu=1,sigma=2,initval=1.5,lower=0),
+            'beta_P4': pm.HalfNormal(f'{self.name}_beta_P4',sigma=1,initval=0.1),
+            'beta_P5': pm.HalfNormal(f'{self.name}_beta_P5',sigma=10,initval=1),
+            'beta_P6': pm.Normal(f'{self.name}_beta_P6',mu=0,sigma=1,initval=0),
+            'beta_P7': pm.Normal(f'{self.name}_beta_P7',mu=0,sigma=1,initval=0),
+            'beta_C1': pm.TruncatedNormal(f'{self.name}_beta_C1',mu=60,sigma=10,initval=60,lower=10),
+            'beta_C2': pm.TruncatedNormal(f'{self.name}_beta_C2',mu=30,sigma=10,initval=30,lower=1),
+            'beta_C3': pm.TruncatedNormal(f'{self.name}_beta_C3',mu=0.05,sigma=0.1,initval=0.05,lower=0),
+            'beta_C4': pm.TruncatedNormal(f'{self.name}_beta_C4',mu=1,sigma=1,initval=1,lower=0),
+        }
+        return param
 
 class WheelS3V3(BaseComponents):
     def __init__(self, name):
@@ -135,9 +224,12 @@ class WheelS3V3(BaseComponents):
         # 处理侧
         HdiffP  = (
             beta_P1 * RinP ** beta_P4  
-            * HinP 
-            * (TinR + beta_P6 * TinC + beta_P7 * TinP) 
-            * EXP(-beta_P5 * FP) 
+            * HinP
+            * (
+                TinR * (1 - EXP(-beta_P5 * FR)) 
+                + beta_P6 * TinC * (1 - EXP(-beta_P5 * FC))
+                + beta_P7 * TinP * (1 - EXP(-beta_P5 * FP))
+            ) 
             + beta_P2 
         ) / 1000
         WdiffP  = HdiffP * FP
@@ -177,7 +269,7 @@ class WheelS3V3(BaseComponents):
             'beta_P2': pm.Normal(f'{self.name}_beta_P2',sigma=1,initval=0),
             'beta_P3': pm.TruncatedNormal(f'{self.name}_beta_P3',mu=1,sigma=2,initval=1.5,lower=0),
             'beta_P4': pm.HalfNormal(f'{self.name}_beta_P4',sigma=1,initval=0.1),
-            'beta_P5': pm.TruncatedNormal(f'{self.name}_beta_P5',mu=5,sigma=2,initval=5,lower=0),
+            'beta_P5': pm.HalfNormal(f'{self.name}_beta_P5',sigma=10,initval=1),
             'beta_P6': pm.Normal(f'{self.name}_beta_P6',mu=0,sigma=1,initval=0),
             'beta_P7': pm.Normal(f'{self.name}_beta_P7',mu=0,sigma=1,initval=0),
             'beta_C1': pm.TruncatedNormal(f'{self.name}_beta_C1',mu=0.5,sigma=1,initval=0.5,lower=0),

model/DHU/DHU_AB.py

@@ -30,10 +30,11 @@ class DHU_AB(BaseDevice):
         wheel_2       = None,
         coolingcoil_2 = 'CoolingCoil2',
         coolingcoil_3 = 'CoolingCoil2',
-        heatingcoil_1 = 'SteamCoil',
-        heatingcoil_2 = 'SteamCoil',
+        heatingcoil_1 = 'SteamCoilVal',
+        heatingcoil_2 = 'SteamCoilVal',
         mixed_1       = 'Mixed',
         mixed_2       = 'Mixed',
+        other_info    = None
     ) -> None:
         super().__init__()
         self.DHU_type = DHU_type.replace('DHU_','')
@@ -58,11 +59,13 @@ class DHU_AB(BaseDevice):
         }
         self.exist_Fa_H = exist_Fa_H
         self.exist_Fa_B = exist_Fa_B
+        self.other_info = other_info if other_info is not None else {}
         self.record_load_info(
             components_str = self.components_str,
             DHU_type       = self.DHU_type,
             exist_Fa_H     = self.exist_Fa_H,
-            exist_Fa_B     = self.exist_Fa_B
+            exist_Fa_B     = self.exist_Fa_B,
+            other_info     = self.other_info
         )
     
     @property
@@ -73,9 +76,15 @@ class DHU_AB(BaseDevice):
         }
         output ={}
         for comp_name,comp_model in self.components_str.items():
+            if comp_model == 'SteamCoilVal':
+                output[comp_name] = coil_steam.SteamCoilVal(
+                    name     = comp_name,
+                    Fs_rated = self.other_info[f'{comp_name}_Fs_rated']
+                )
+                continue
             for comp_map_k,comp_map_v in comp_map.items():
                 if comp_model.startswith(comp_map_k):
-                    output[comp_name] = getattr(comp_map_v,comp_model)(comp_name)
+                    output[comp_name] = getattr(comp_map_v,comp_model)(name = comp_name)
         return output
     
     @property
@@ -125,6 +134,10 @@ class DHU_AB(BaseDevice):
                 columns[f'{heatingcoil_idx}_Fs'] = f'{heatingcoil_idx}_Fs'
             elif isinstance(self.components[heatingcoil_idx],coil_steam.SteamCoil):
                 columns['wheel_2_ToutC'] = 'wheel_2_ToutC'
+                columns['wheel_1_ToutR'] = 'wheel_1_ToutR'
+                columns['mixed_2_ToutA'] = 'mixed_2_ToutA'
+            elif isinstance(self.components[heatingcoil_idx],coil_steam.SteamCoilVal):
+                columns[f'{heatingcoil_idx}_Val'] = f'{heatingcoil_idx}_Val'
             else:
                 raise Exception('WRONG')
         return columns
@@ -156,9 +169,9 @@ class DHU_AB(BaseDevice):
             )
             for std_name,name in self.model_observe_data_columns.items():
                 if name not in observed_data.columns:
-                    continue
+                    raise Exception(f'Missing column: {name}')
                 observed_data = observed_data.rename(columns={name:std_name})
-
+                
                 std_name_equp,std_name_point = std_name.rsplit('_',1)
                 sigma = {
                     'wheel_2_DoutP'   : 0.3,
@@ -167,11 +180,22 @@ class DHU_AB(BaseDevice):
                     'heatingcoil_1_FP': 10000,
                     'heatingcoil_2_FP': 10000,
                 }
+                if std_name in ['heatingcoil_1_Val','heatingcoil_2_Val']:
+                    sigma = res[std_name_equp]['sigma']
+                else:
+                    sigma = {
+                        'wheel_2_DoutP'   : 0.3,
+                        'heatingcoil_1_Fs': 20,
+                        'heatingcoil_2_Fs': 20,
+                        'heatingcoil_1_FP': 10000,
+                        'heatingcoil_2_FP': 10000,
+                    }.get(std_name,1)
+                
                 observe(
                     name     = std_name,
                     var      = res[std_name_equp][std_name_point],
                     observed = observed_data,
-                    sigma    = sigma.get(std_name,1)
+                    sigma    = sigma
                 )
             
             self.param_posterior = pm.find_MAP(maxeval=50000,include_transformed=False)
@@ -233,20 +257,26 @@ class DHU_AB(BaseDevice):
     
     def optimize(
         self,
-        cur_input_data : pd.DataFrame,
-        wheel_1_TinR_ub: float = 120,
-        wheel_1_TinR_lb: float = 70,
-        wheel_2_TinR_ub: float = 120,
-        wheel_2_TinR_lb: float = 70,
-        constrains     : list  = None,
-        logging        : bool  = True
+        cur_input_data: pd.DataFrame,
+        wheel_1_TinR  : tuple = (70,120),
+        wheel_2_TinR  : tuple = (70,120),
+        fan_2_Hz      : tuple = (30,50),
+        constrains    : list  = None,
+        logging       : bool  = True,
+        target        : str   = 'summary_Fs',
+        target_min    : bool  = True
     ) -> list:
         constrains = [] if constrains is None else constrains
         cur_input_data = cur_input_data.iloc[[0],:]
-        opt_var_boundary  = {
-            'wheel_1_TinR':{'lb':wheel_1_TinR_lb,'ub':wheel_1_TinR_ub},
-            'wheel_2_TinR':{'lb':wheel_2_TinR_lb,'ub':wheel_2_TinR_ub},
-        }
+        
+        opt_var_boundary = {}
+        if wheel_1_TinR is not None:
+            opt_var_boundary['wheel_1_TinR'] = {'lb':min(wheel_1_TinR),'ub':max(wheel_1_TinR)}
+        if wheel_2_TinR is not None:
+            opt_var_boundary['wheel_2_TinR'] = {'lb':min(wheel_2_TinR),'ub':max(wheel_2_TinR)}
+        if fan_2_Hz is not None:
+            opt_var_boundary['fan_2_Hz'] = {'lb':min(fan_2_Hz),'ub':max(fan_2_Hz)}
+            
         opt_var_value = cur_input_data.loc[:,list(opt_var_boundary.keys())]
         oth_var_value = (
             cur_input_data
@@ -258,8 +288,8 @@ class DHU_AB(BaseDevice):
             opt_var_boundary = opt_var_boundary,
             opt_var_value    = opt_var_value,
             oth_var_value    = oth_var_value,
-            target           = 'summary_Fs',
-            target_min       = True,
+            target           = target,
+            target_min       = target_min,
             constrains       = constrains,
             logging          = logging,
             other_kwargs     = {'NIND':2000,'MAXGEN':50}
@@ -649,7 +679,7 @@ class AirFlow_DHU_AB:
         
         elif type == 'DHU_B':
             wheel_1_FaP    = Fa_X
-            wheel_1_FaC    = None
+            wheel_1_FaC    = np.nan
             wheel_1_FaR    = Fa_P
             wheel_2_FaP    = Fa_S
             wheel_2_FaC    = Fa_X + Fa_H - Fa_S

model/utils/fit_utils.py

@@ -10,8 +10,11 @@ def observe(name,var,observed,sigma=1):
         if name not in observed.columns:
             raise Exception(f'observed data中未找到{name}')
         observed = observed.loc[:,name].values
-    mu    = pm.Deterministic(f'{name}_mu',var)
-    sigma = pm.HalfNormal(f'{name}_sigma',sigma=sigma)
+    mu = pm.Deterministic(f'{name}_mu',var)
+    if isinstance(sigma,(int,float)):
+        sigma = pm.HalfNormal(f'{name}_sigma',sigma=sigma)
+    else:
+        sigma = sigma
     pm.Normal(name,mu=mu,sigma=sigma,observed=observed)
 
 def reorder_posterior(prior:dict,posterior:dict):