How to use property_matrix method in Kiwi

Best Python code snippet using Kiwi_python

solverFunctions.py

Source:solverFunctions.py

1# -*- coding: utf-8 -*-2"""3Created on Tue Mar 23 11:29:26 202145@author: nicol6"""789import GeometryFunctions as geoF10from math import pi, exp11import numpy as np12import propertyFunctions as propF13import processFunctions as procF14import pickle15import os.path1617# property_matrix = loadmat('R1233zd.mat') 18# m_dot_in = 019# m_dot_out = 020# x_l_in = None21# x_l_out = None22# x_l = 023# h_in = 024# h_out = 025# Q_dot = 026# T = 36027# rho = 3028# m_cv = 0.129# dVdTheta=-10/1e630# omega = 3003132# A 'controlVolume' class is used to represent every chamber within the scroll machine. The object is used to 33# keep track of the evolving properties of the chamber as the crank rotation goes on. The following input are34# required to generate the chamber object in first place:35# - name: identifies the kind of chamber that is generated.36# - geo: object that contains all the geometric properties of the machine.37# - property_matrix: table in which all the thermodynamic properties of the fluid are stored.38# - dict_V: dictionaries in which the coefficients of the polynamial that approximate chambers' volume evolution are specified.39# - dict_dV: dictionaries in which the coefficients of the polynamial that approximate volume derivtives evolution are specified.40# - dict_area: dictionaries in which the coefficients of the polynamial that approximate surfaces evolution are specified.41# - T_in, rho_in, x_l_in: known thermodynamic properties.42class controlVolume:43    44    # Constructor function for the 'controlVolume' class45    def __init__(self, name, geo, property_matrix, dict_V, dict_dV, dict_area, T_in, rho_in, x_l_in, P_out = None, alpha = 1):46        47        self.name = name48        self.alpha = alpha49        self.dict_V = dict_V  50        self.dict_dV = dict_dV 51        self.dict_area = dict_area52        self.theta = 053        54        # The inlet values of temperature and density are directly assigned as initial values of temperature and55        # density for suction chambers and first pair of compression chambers. The volume is evaluated through the56        # polyval function applied to the polynomial specified in 'dict_V'57        if name == "sa" or name =="s1" or name =="s2" or (name == "c1" and alpha == 1) or (name == "c2" and alpha == 1) : 58            self.theta_shift = 059            (self.T_init, self.rho_init) = (T_in, rho_in)  60            self.V_init = np.polyval(dict_V[name][alpha-1], self.theta + self.theta_shift)61            # The volume of the suction chamber 's1' is initialized as 5e-8 (which is a convential value to represent62            # a null volume)63            if name =='s1': 64                self.V_init = 5e-865            self.dict_init = {'T' :  [0] , 'rho' :[0] }66        67        # The initial values of temperature and density in the discharge chambers are evaluated by means of the68        # function 'discharge_is_guess'. This function use the properties at machine inlet, the pressure ratio69        # and a first-try value of the isentropic efficiency to guess the initial value of temperature and density70        # inside the discharge chambers.71        elif name == "ddd" or name == 'd1' or name == 'd2' or name == 'dd':72            # A first-attempt value of the isentropic efficiency is specified 73            eta_is = 0.674            (self.T_init, self.rho_init) = procF.discharge_is_guess(property_matrix, rho_in, T_in, P_out, x_l_in, eta_is)75            # ***???***76            if (geo.theta_d < 3/2*pi and (name == 'ddd' or name =='dd')) or (geo.theta_d >= 3/2*pi and ( name == 'd1' or name == 'd2')):77                self.theta_shift = 2*pi78                self.x_out = x_l_in79            else:80                self.theta_shift = 081            self.V_init = np.polyval(dict_V[name][alpha-1], self.theta + self.theta_shift)82            self.dict_init = {'T' :  [0] , 'rho' :[0] }83            84        # The initial values of temperature and density in the inner compression chambers are evaluated by means 85        # of the function 'compression_is_guess'. This function use the properties at machine inlet, the pressure ratio86        # and a first-try value of the isentropic efficiency to guess the initial value of temperature and density87        # inside the inner compression chambers.88        elif (name == 'c1' or name == 'c2') and alpha > 1:89            self.theta_shift = 090            eta_is = 0.691            (self.T_init, self.rho_init) = procF.compression_is_guess(property_matrix, dict_V[name], self.alpha, rho_in, T_in,  x_l_in, eta_is)92            self.V_init =  np.polyval(dict_V[name][alpha-1], self.theta + self.theta_shift)93            self.dict_init = {'T' :  [0] , 'rho' :[0] }94        else:95            raise Exception('Invalid inputs')96        97        # Other properties are computed by means of the initial values of temperature and density computed thus far.98        self.x_l_init = x_l_in99        self.P_init = propF.propfunction_mixture(property_matrix, self.x_l_init, 'P', T = self.T_init , rho = self.rho_init)   100        self.Q_init = propF.propfunction_mixture(property_matrix, self.x_l_init, 'Q', T = self.T_init , rho = self.rho_init)    101        self.s_init = propF.propfunction_mixture(property_matrix, self.x_l_init, 's', T = self.T_init , rho = self.rho_init) 102        self.h_init = propF.propfunction_mixture(property_matrix, self.x_l_init, 'h', T = self.T_init , rho = self.rho_init)103        self.m_init =  self.rho_init * self.V_init   104        105        # 106        self.h_out = propF.propfunction_mixture(property_matrix, self.x_l_init, 'h', T = self.T_init , rho = self.rho_init)   107108        # 109        self.m_dot_dict_init = {'su' : 0., 'ex' : 0., 'fl1' : 0., 'fl2' : 0., 'r1' : 0., 'r2' : 0., 'r3' : 0., 'r4' : 0.}110        self.h_dict_init = {'su' : 0., 'ex' : 0., 'fl1' : 0., 'fl2' : 0., 'r1' : 0., 'r2' : 0., 'r3' : 0., 'r4' : 0.}111        self.x_l_dict_init = {'su' : 0., 'ex' : 0., 'fl1' : 0., 'fl2' : 0., 'r1' : 0., 'r2' : 0., 'r3' : 0., 'r4' : 0.}112        113        # - The 'mh_sum' attribute is used to compute the total amount of energy flowing through the considered chamber114        # - The 'm_sum' attribute is used to compute the total amount of mass flowing through the considered chamber115        # - The 'h_disc' attribute is used to compute the discharge enthalpy of the chamber.116        self.mh_sum = 0117        self.m_sum = 1e-10118        self.h_disc = 0119        120        121    # The 'next_values' method takes as input the property matrix, the crank angle and the step. Within this method the122    # functions that compute mass and temperature derivatives are called. The output of this method are the updated values123    # of all the properties within the considered control volume.124    def next_values(self, property_matrix, omega, theta_step, x_l = 0, Q_dot = 0):125        126        # Only volume is updated if the considered control volume is the suction 'sa' chamber.127        if self.name == 'sa':128            self.theta = self.theta + theta_step129            self.V = np.polyval(self.dict_V[self.name][self.alpha-1], self.theta)  130        else:131            # Computation of the mass derivative132            dmCVdTheta = procF.dmCV_dTheta(omega, self.m_dot_dict)133            # Computation of oil mass fraction derivative134            dxoildTheta = procF.dxoil_dTheta(omega, self.m, np.polyval(self.dict_dV[self.name][self.alpha-1], self.theta + self.theta_shift), self.m_dot_dict, self.x_l_dict,  self.x_l)135            # Computation of temperature derivative.136            dTdTheta = procF.dT_dTheta(property_matrix, omega, self.m, self.T, self.rho, np.polyval(self.dict_dV[self.name][self.alpha-1], self.theta + self.theta_shift),  dmCVdTheta, dxoildTheta, self.m_dot_dict, self.h_dict, self.x_l, Q_dot)137            138            # If chamber's volume is lower than 5e-8 the derivatives are set equal to 0 in order to avoid any139            # numerical issue related with division by very small numbers.140            if self.V <=  5e-8: # avoid error suction chamber 141                dTdTheta = 0142                dmCVdTheta = 0143            144            # Attributes are updated accordingly with the derivatives just computed145            self.theta = self.theta + theta_step146            self.T = self.T + theta_step*dTdTheta147            self.m = self.m + theta_step*dmCVdTheta148            self.x_l = self.x_l + theta_step*dxoildTheta       149            150            # Evaluation of chamber volume151            self.V =  np.polyval(self.dict_V[self.name][self.alpha-1], self.theta + self.theta_shift) 152            # If volume is lower than 5e-8, its value is set to 5e-8 to avoid error in the suction chamber153            if self.V <=  5e-8:154                self.V = 5e-8155            156            # Computation of fluid density within the control volume.157            self.rho = self.m/self.V158            # if self.name == 'dd':159            160            #     print(self.name,self.theta)161162            # Computation of other properties based on temperature (computed thanks to energy equation) and163            # density (computed via the definition of density itself).164            self.Q = propF.propfunction_mixture(property_matrix, self.x_l, 'Q', T = self.T , rho = self.rho)165            self.P = propF.propfunction_mixture(property_matrix, self.x_l, 'P', T = self.T , rho = self.rho)166            self.s = propF.propfunction_mixture(property_matrix, self.x_l, 's', T = self.T , rho = self.rho) 167            self.h = propF.propfunction_mixture(property_matrix, self.x_l, 'h', T = self.T , rho = self.rho) 168            169    # The 'record_values' function is used to record the value of the properties within each chamber at every crank 170    # angle step during the last crank angle rotation (after stopping criteria are met).171    def record_values(self):172                173        self.P_vect_final.append(self.P)174        self.T_vect_final.append(self.T)175        self.Q_vect_final.append(self.Q)176        self.s_vect_final.append(self.s)177        self.m_vect_final.append(self.m)178        self.x_l_vect_final.append(self.x_l)179        self.rho_vect_final.append(self.rho)180        self.theta_vect_final.append(self.theta)181        self.V_vect_final.append(self.V)182        # 'self.m_dot_vect_final' is a dictionary composed of different lists, one for each mass flow that occurs within183        # the scroll machine.184        for k in self.m_dot_vect_final: 185            self.m_dot_vect_final[k].append(self.m_dot_dict[k])186    187    def discharge(self, old_CV):188        189        190        self.T = old_CV.T191        self.rho = old_CV.rho192        self.x_l = old_CV.x_l193        self.m = old_CV.m194        self.V = old_CV.V195        self.Q = old_CV.Q196        self.P = old_CV.P197        self.s = old_CV.s198        self.h = old_CV.h199        200        self.theta = old_CV.theta201        202        if old_CV.name == 'ddd':203            204            self.mh_sum = old_CV.mh_sum 205            self.m_sum = old_CV.m_sum206            self.h_disc = old_CV.h_disc207            self.h_out = old_CV.h_out208        209210    # Initialize is a method of the 'controlVolume' class. It does not require any input, 'CV' is an optional input whose default value211    # is set to 'None'212    def initialize(self, CV = None):213    214        # If a CV is passed as input to the function, its properties are attached as 215        # init properties to the 'controlVolume' object216        if CV != None:217            self.theta = 0218            self.T_init = CV.T219            self.rho_init = CV.rho220            self.x_l_init = CV.x_l221            self.P_init = CV.P  222            self.Q_init = CV.Q  223            self.s_init = CV.s224            self.h_init = CV.h225            226            # 'V_init' is evaluated thanks to 'dict_V', where all the information about volumes evolution are stored.227            self.V_init = np.polyval(self.dict_V[self.name][self.alpha-1], self.theta + self.theta_shift)228            229            # The volume of the suction chamber 's1' is always initialized as 5e-8230            if self.name == 's1': 231                self.V_init = 5e-8232                233            # Mass inside the control volume is evaluated as volume times density234            self.m_init = self.V_init*self.rho_init 235           236            # The dictionaries related to mass flows are re-initialized.237            self.m_dot_dict_init = {'su' : 0., 'ex' : 0., 'fl1' : 0., 'fl2' : 0., 'r1' : 0., 'r2' : 0., 'r3' : 0., 'r4' : 0.}238            self.h_dict_init = {'su' : 0., 'ex' : 0., 'fl1' : 0., 'fl2' : 0., 'r1' : 0., 'r2' : 0., 'r3' : 0., 'r4' : 0.}239            self.x_l_dict_init = {'su' : 0., 'ex' : 0., 'fl1' : 0., 'fl2' : 0., 'r1' : 0., 'r2' : 0., 'r3' : 0., 'r4' : 0.}240            241        # Initial properties are attached to the properties attributes.242        self.T = self.T_init243        self.rho = self.rho_init244        self.x_l = self.x_l_init245        self.P = self.P_init246        self.Q = self.Q_init247        self.s = self.s_init248        self.h = self.h_init249        self.m = self.m_init     250        self.V = self.V_init251        252        # Initial properties related to mass flows are attached to the properties attributes.253        self.m_dot_dict = self.m_dot_dict_init254        self.h_dict = self.h_dict_init255        self.x_l_dict = self.x_l_dict_init256        257        # The dictionary containing the initial properties is initialized.258        self.dict_init['T'].append(self.T_init)  259        self.dict_init['rho'].append(self.rho_init)260        261        # If discharge entalpy is different from 0 and the difference between outlet enthalpy and discharge enthalpy262        # is lower than 0.1%, outlet enthalpy is set equal to discharge enthalpy263        if self.h_disc != 0 and (self.h_out - self.h_disc)/self.h_disc > 0.001  : 264            self.h_out = self.h_disc265     266        # Initialization of the attribute in which energy flowing through the control volume is stored267        self.mh_sum = 0268        # Initialization of the attribute in which mass flowing through the control volume is stored269        self.m_sum = 1e-10270        # Initialization of the attribute in which discharge enthalpy is stored271        self.h_disc = 0272        273        # Attributes for the final properties are initialized.274        self.P_vect_final = []275        self.T_vect_final = []276        self.Q_vect_final = []277        self.s_vect_final = []278        self.m_vect_final = []279        self.x_l_vect_final = []280        self.rho_vect_final = []281        self.theta_vect_final = []282        self.V_vect_final = []283        self.m_dot_vect_final = {'su' : [], 'ex' : [], 'fl1' : [], 'fl2' : [], 'r1' : [], 'r2' : [], 'r3' : [], 'r4' : []}284    285        286    # The 'update' method of the 'controlVolume' class computes pressure, quality, entropy and enthalpy 287    # of a given control volume based on its pressure and density.288    def update(self, property_matrix):    289        self.P = propF.propfunction_mixture(property_matrix, self.x_l, 'P', T = self.T , rho = self.rho)290        self.Q = propF.propfunction_mixture(property_matrix, self.x_l, 'Q', T = self.T , rho = self.rho)291        self.s = propF.propfunction_mixture(property_matrix, self.x_l, 's', T = self.T , rho = self.rho) 292        self.h = propF.propfunction_mixture(property_matrix, self.x_l, 'h', T = self.T , rho = self.rho) 293        294    295    # The 'mass_flow_su' method computes the mass flow rate at suction based on the assumption that it behaves296    # as an isentropic nozzle 297    def mass_flow_su(self, property_matrix, geo, rho_in, T_in, x_l_in, theta_step):298299        # Suction flow is non-zero if and only if the considered chamber is a suction chamber, otherwise the mass300        # flow rate is set equal to 0.         301        if self.name == 's1' or self.name == 's2':302            # Evaluation of the suction port area at the current crank angle303            A_su = np.polyval(self.dict_area['s'], self.theta)304            # Mass flow rate is computed based on the isentropic nozzle hypothesis305            (m_dot_su, h_su, _)  = procF.isentropic_nozzle(geo, property_matrix, A_su,  self.rho,  self.T, self.x_l, rho2 = rho_in,  T2 = T_in,  x_l2 = x_l_in)306            # Chamber dctionaries are updted307            self.m_dot_dict['su'] = m_dot_su308            self.h_dict['su'] = h_su            309            # Energy flow, mass flow and discharge enthalpy of the current object are updated.310            self.mh_sum = self.mh_sum + m_dot_su*h_su311            self.m_sum = self.m_sum + m_dot_su/(2*pi)*theta_step312            self.h_disc = self.mh_sum/self.m_sum/(2*pi)*theta_step313        else:314            (self.m_dot_dict['su'], self.h_dict['su']) = (0., 0.)315          316        317    # The 'mass_flow_ex' method computes the mass flow rate at discharge based on the assumption that it behaves318    # as an isentropic nozzle     319    def mass_flow_ex(self, property_matrix, geo, P_out, theta_step, Reed_valve):    320 321        # Discharge flow is non-zero if and only if the considered chamber is a discharge chamber, otherwise the mass322        # flow rate is set equal to 0. A simple model of a discharge valve is implemented: the flow is non-zero only323        # if the pressure inside the discharge chamber is higher than the pressure in the discharge pipe.324        if self.name == 'ddd' or self.name == 'dd': 325            A_port = np.polyval(self.dict_area['port'], self.theta)326            if Reed_valve == True:327                delta_P = 0;328                if self.P > P_out + delta_P: 329                    (m_dot_ex, h_ex, _) = procF.isentropic_nozzle(geo, property_matrix, A_port,  self.rho,  self.T, self.x_l, h2 = self.h_out, P2 = P_out, x_l2 = self.x_out)  330                else :331                    (m_dot_ex, h_ex, _) = (0,0,0)332            if Reed_valve == False:333                (m_dot_ex, h_ex, _) = procF.isentropic_nozzle(geo, property_matrix, A_port,  self.rho,  self.T, self.x_l, h2 = self.h_out, P2 = P_out, x_l2 = self.x_out)334            335            self.m_dot_dict['ex'] = m_dot_ex336            self.h_dict['ex'] = h_ex337            338339            self.mh_sum = self.mh_sum + m_dot_ex*h_ex340            self.m_sum = self.m_sum + m_dot_ex/(2*pi)*theta_step 341            self.h_disc = self.mh_sum/self.m_sum/(2*pi)*theta_step342        else:343            (self.m_dot_dict['ex'], self.h_dict['ex']) = (0., 0.)344        345    346    # The 'mass_flow_dis_dd' method computes the mass flow rate from chambers 'd1' and 'd2' to chamber 'dd' based347    # on the assumption that it behaves as an isentropic nozzle348    def mass_flow_dis_dd(self, property_matrix, geo, CV_d1):349        350        # The flow from chambers 'd1' and 'd2' to chamber 'dd' is non-zero only if the chamber we are considering is351        # the 'dd' chamber.352        if self.name == 'dd':353            # Evaluation of the port area that is present between chambers 'dd' and 'd1' (or 'd2').354            A_dis = np.polyval(self.dict_area['d'], self.theta)355            # Mass flow rate is computed based on the isentropic nozzle hypothesis356            (m_dot_dis, h_dis, _) = procF.isentropic_nozzle(geo, property_matrix, A_dis,  self.rho,  self.T, self.x_l, rho2 = CV_d1.rho,  T2 = CV_d1.T,  x_l2 = CV_d1.x_l)357            # Mass flow rate is multiplied by 2 because the same flow rate comes from 'd2' and 'd1'358            self.m_dot_dict['fl2'] = 2*m_dot_dis 359            # Chamber dctionaries are updted 360            self.h_dict['fl2'] = h_dis        361            CV_d1.m_dot_dict['fl2'] = - m_dot_dis362            CV_d1.h_dict['fl2'] = h_dis 363    364    365    # The 'leakages' method computes the mass flow rate through internal chambers of the machine.366    # The 'leakages_computation' function passes to the current method the correct objects, based on the kind leakage367    # that has to be modeled.368    def leakages(self, property_matrix, geo, exist, CVs1 = None, CVc1 = None, CVd1 = None, CVddd = None):369        370        A_fl = flank_area(geo)371        372        # The number of existing compression chambers is computed accordingly to the geometry373        N_c = geoF.N_c_calculation(geo,self.theta)374        375        # The correct name of the discharge chamber is extract from the input variable 'exist'376        if exist =='d':377            CVdis = CVd1378            name_dis = 'd1'379            coef = 1380        elif exist == 'ddd':381            CVdis = CVddd382            name_dis = 'ddd'383            coef = 2384            385        # The radial leakage from chamber 's1' to chamber 'sa' is always present and for this reason is computed386        # ahead of everythink else.387        if self.name == 's1' :388                phi_ssa = geoF.get_phi_ssa(geo, self.theta)389                # Radial leakage area is computed through the 'radial_area' function from the 'GeometryFunction' module390                A_sa = radial_area(geo, geo.phi_ie , max( (geo.phi_ie - self.theta), phi_ssa))391                # Mass flow rate is computed using the isentropic nozzle assumption392                (self.m_dot_dict['r1'], self.h_dict['r1'], Re) = procF.isentropic_nozzle(geo, property_matrix,  A_sa,  self.rho,  self.T, self.x_l, rho2 = self.rho_init,  T2 = self.T_init,  x_l2 = self.x_l_init, leakage_type = 'r')393                # Leakage Mass flow rate is affected by Reynolds number. Therefore, a correction factor is computed 394                # in the 'corrected_mass_flow' function, accordingly with the empyrical model presented in Ian Bell's thesis. 395                self.m_dot_dict['r1'] = corrected_mass_flow(geo, Re, self.m_dot_dict['r1'], 'r')396397        # Different leakage paths exist depending on the number of compression chambers and therefore on the crank angle 'theta'.398        # If zero compression chambers exist, the only existing leakage paths are:399            # - radial leakage 's1' - 'd2'400            # - flank leakage 's1' - 'd2'401            # - radial leakage 's2' - 'discharge'402            # - flank leakage 's2' - 'discharge'403        if N_c == 0:404            if self.name == 's1' :405                # radial leakage 's1' - 'd2'406                A_d2 = radial_area(geo, geo.phi_ie - self.theta, geo.phi_ie - self.theta - pi )407                (self.m_dot_dict['r3'], self.h_dict['r3'], Re) = procF.isentropic_nozzle(geo, property_matrix,  A_d2,  self.rho,  self.T, self.x_l, rho2 = CVdis.rho,  T2 = CVdis.T,  x_l2 = CVdis.x_l, leakage_type = 'r')408                self.m_dot_dict['r3'] = corrected_mass_flow(geo, Re, self.m_dot_dict['r3'], 'r')409                (self.m_dot_dict['r2'], self.h_dict['r2']) = (0., 0.)410                (self.m_dot_dict['r4'], self.h_dict['r4']) = (0., 0.)411                # flank leakage 's1' - 'd2'412                (self.m_dot_dict['fl2'], self.h_dict['fl2'], Re) = procF.isentropic_nozzle(geo, property_matrix,  A_fl,  self.rho,  self.T, self.x_l, rho2 = CVdis.rho,  T2 = CVdis.T,  x_l2 = CVdis.x_l, leakage_type = 'fl')413                self.m_dot_dict['fl2'] = corrected_mass_flow(geo, Re, self.m_dot_dict['fl2'], 'fl')414                (self.m_dot_dict['fl1'], self.h_dict['fl1']) = (0., 0.)                415            416            elif self.name == name_dis :417                # radial leakage 'discharge' - 's2'418                A_s2 = radial_area(geo, geo.phi_ie - self.theta, geo.phi_ie - self.theta - pi )419                (m_dot, h, Re) = procF.isentropic_nozzle(geo, property_matrix,  A_s2,  self.rho,  self.T, self.x_l, rho2 = CVs1.rho,  T2 = CVs1.T,  x_l2 = CVs1.x_l, leakage_type = 'r')420                m_dot = corrected_mass_flow(geo, Re, m_dot, 'r')421                (self.m_dot_dict['r2'], self.h_dict['r2']) = (coef*m_dot, h)422                (self.m_dot_dict['r1'], self.h_dict['r1']) = (0., 0.)423                (self.m_dot_dict['r2'], self.h_dict['r2']) = (0., 0.)424                (self.m_dot_dict['r4'], self.h_dict['r4']) = (0., 0.)            425                # flank leakage 'discharge' - 's2' 426                (m_dot, h, Re) = procF.isentropic_nozzle(geo, property_matrix,  A_fl,  self.rho,  self.T, self.x_l, rho2 = CVs1.rho,  T2 = CVs1.T,  x_l2 = CVs1.x_l, leakage_type = 'fl')427                m_dot = corrected_mass_flow(geo, Re, m_dot, 'fl')428                (self.m_dot_dict['fl1'], self.h_dict['fl1']) = (coef*m_dot, h)429                (self.m_dot_dict['fl2'], self.h_dict['fl2']) = (0., 0.)               430            431            else:432                (self.m_dot_dict['r1'], self.h_dict['r1']) = (0., 0.)433                (self.m_dot_dict['r2'], self.h_dict['r2']) = (0., 0.)434                (self.m_dot_dict['r3'], self.h_dict['r3']) = (0., 0.)435                (self.m_dot_dict['r4'], self.h_dict['r4']) = (0., 0.)     436                (self.m_dot_dict['fl1'], self.h_dict['fl1']) = (0., 0.)437                (self.m_dot_dict['fl2'], self.h_dict['fl2']) = (0., 0.)    438         439        # If a pair compression chambers exist, the existing leakage paths are:440            # - radial leakage 's1' - 'c21' (same as c11) 441            # - flank leakage 's1' - 'c21' (same as c11) 442            # - radial leakage 'c21' (same as c11) - 's1'443            # - radial leakage 'c21' (same as c11) - 's1'444            # - radial leakage 'c21' (same as c11) - 'sa'445            # - flank leakage 'c21' - 's2'446            # - flank leakage 'c21' - 'discharge'447            # - radial leakage 'discharge' - 'c2' (alpha = N_c)448            # - flank leakages 'discharge' - 'c2' (alpha = N_c)449        elif N_c == 1:450            if self.name == 's1' :451                # radial leakage 's1' - 'c21' (same as c11) 452                phi_ssa = geoF.get_phi_ssa(geo, self.theta)453                A_c21 = radial_area(geo, max( geo.phi_ie - self.theta - pi , phi_ssa) , geo.phi_ie - self.theta - pi)454                (self.m_dot_dict['r3'], self.h_dict['r3'], Re) = procF.isentropic_nozzle(geo, property_matrix,  A_c21,  self.rho,  self.T, self.x_l, rho2 = CVc1[1].rho,  T2 = CVc1[1].T,  x_l2 = CVc1[1].x_l, leakage_type = 'r')455                self.m_dot_dict['r3'] = corrected_mass_flow(geo, Re, self.m_dot_dict['r3'], 'r')456                (self.m_dot_dict['r2'], self.h_dict['r2']) = (0., 0.)457                (self.m_dot_dict['r4'], self.h_dict['r4']) = (0., 0.)                458                # flank leakage 's1' - 'c21' (same as c11) 459                (self.m_dot_dict['fl2'], self.h_dict['fl2'], Re) = procF.isentropic_nozzle(geo, property_matrix,  A_fl,  self.rho,  self.T, self.x_l, rho2 = CVc1[1].rho,  T2 = CVc1[1].T,  x_l2 = CVc1[1].x_l, leakage_type = 'fl')460                self.m_dot_dict['fl2'] = corrected_mass_flow(geo, Re, self.m_dot_dict['fl2'], 'fl')461                (self.m_dot_dict['fl1'], self.h_dict['fl1']) = (0., 0.)                462                463            elif self.name == 'c1' :  464                # radial leakage 'c21' (same as c11) - 's1'465                phi_ssa = geoF.get_phi_ssa(geo, self.theta)466                A_s2 = radial_area(geo, max( geo.phi_ie - self.theta - pi , phi_ssa) , geo.phi_ie - self.theta - pi)467                (self.m_dot_dict['r2'], self.h_dict['r2'], Re) = procF.isentropic_nozzle(geo, property_matrix,  A_s2,  self.rho,  self.T, self.x_l, rho2 = CVs1.rho,  T2 = CVs1.T,  x_l2 = CVs1.x_l, leakage_type = 'r')468                self.m_dot_dict['r2'] = corrected_mass_flow(geo, Re, self.m_dot_dict['r2'], 'r') 469                # radial leakage 'c21' (same as c11) - 's1'470                A_d2 = radial_area(geo, geo.phi_ie - self.theta - 2*pi, geo.phi_is)471                (self.m_dot_dict['r4'], self.h_dict['r4'], Re) = procF.isentropic_nozzle(geo, property_matrix,  A_d2,  self.rho,  self.T, self.x_l, rho2 = CVdis.rho,  T2 = CVdis.T,  x_l2 = CVdis.x_l, leakage_type = 'r')472                self.m_dot_dict['r4'] = corrected_mass_flow(geo, Re, self.m_dot_dict['r4'], 'r') 473                # radial leakage 'c21' (same as c11) - 'sa'474                A_sa = radial_area(geo,  max( (geo.phi_ie - self.theta), phi_ssa), phi_ssa)475                (self.m_dot_dict['r1'], self.h_dict['r1'], Re) = procF.isentropic_nozzle(geo, property_matrix,  A_sa,  self.rho,  self.T, self.x_l, rho2 = self.rho_init,  T2 = self.T_init,  x_l2 = self.x_l_init, leakage_type = 'r')476                self.m_dot_dict['r1'] = corrected_mass_flow(geo, Re, self.m_dot_dict['r1'], 'r') 477                (self.m_dot_dict['r3'], self.h_dict['r3']) = (0., 0.)                  478                # flank leakage 'c21' - 's2'479                (self.m_dot_dict['fl1'], self.h_dict['fl1'], Re) = procF.isentropic_nozzle(geo, property_matrix,  A_fl,  self.rho,  self.T, self.x_l, rho2 = CVs1.rho,  T2 = CVs1.T,  x_l2 = CVs1.x_l, leakage_type = 'fl')480                self.m_dot_dict['fl1'] = corrected_mass_flow(geo, Re, self.m_dot_dict['fl1'], 'fl') 481                # flank leakage 'c21' - 'discharge'482                (self.m_dot_dict['fl2'], self.h_dict['fl2'], Re) = procF.isentropic_nozzle(geo, property_matrix,  A_fl,  self.rho,  self.T, self.x_l, rho2 = CVdis.rho,  T2 = CVdis.T,  x_l2 = CVdis.x_l, leakage_type = 'fl')  483                self.m_dot_dict['fl2'] = corrected_mass_flow(geo, Re, self.m_dot_dict['fl2'], 'fl') 484            485            elif self.name == name_dis :486                # radial leakage 'discharge' - 'c2' (alpha = N_c)487                A_c2 = radial_area(geo, geo.phi_ie - self.theta - 2*pi, geo.phi_is)488                (m_dot, h, Re) = procF.isentropic_nozzle(geo, property_matrix,  A_c2,  self.rho,  self.T, self.x_l, rho2 = CVc1[1].rho,  T2 = CVc1[1].T,  x_l2 = CVc1[1].x_l, leakage_type = 'r')489                m_dot = corrected_mass_flow(geo, Re, m_dot, 'r') 490                (self.m_dot_dict['r1'], self.h_dict['r1']) = (coef*m_dot, h)491                (self.m_dot_dict['r2'], self.h_dict['r2']) = (0., 0.)492                (self.m_dot_dict['r3'], self.h_dict['r3']) = (0., 0.)493                (self.m_dot_dict['r4'], self.h_dict['r4']) = (0., 0.) 494                # flank leakages 'discharge' - 'c2' (alpha = N_c)495                (m_dot, h, Re) = procF.isentropic_nozzle(geo, property_matrix,  A_fl,  self.rho,  self.T, self.x_l, rho2 = CVc1[1].rho,  T2 = CVc1[1].T,  x_l2 = CVc1[1].x_l, leakage_type = 'fl')496                m_dot = corrected_mass_flow(geo, Re, m_dot, 'fl') 497                (self.m_dot_dict['fl1'], self.h_dict['fl1']) = (coef*m_dot, h)498                (self.m_dot_dict['fl2'], self.h_dict['fl2']) = (0., 0.)                               499        500            else :501                (self.m_dot_dict['r1'], self.h_dict['r1']) = (0., 0.)502                (self.m_dot_dict['r2'], self.h_dict['r2']) = (0., 0.)503                (self.m_dot_dict['r3'], self.h_dict['r3']) = (0., 0.)504                (self.m_dot_dict['r4'], self.h_dict['r4']) = (0., 0.)505                (self.m_dot_dict['fl1'], self.h_dict['fl1']) = (0., 0.)506                (self.m_dot_dict['fl2'], self.h_dict['fl2']) = (0., 0.)               507                508                509        # If more than 1 pair compression chambers exist, the existing leakage paths are:510            # - radial leakage 's1' - 'c21' (same as c11) 511            # - flank leakage 's1' - 'c21' (same as c11) 512            # - radial leakage 'c1' - 's1'513            # - radial leakage 'c1' - 'sa'514            # - radial leakage 'c1' - 'c2' (alpha = self.alpha + 1)515            # - flank leakage with 'c1' - 's2'516            # - flank leakage 'c1' - 'c2' (alpha = self.alpha + 1)517            # - radial leakage 'c1' - 'c2' (alpha = self.alpha - 1)518            # - radial leakage 'c1' - 'd2' (alpha = N_c)519            # - flank leakages 'c1' - 'c2' (alpha = self.alpha - 1)520            # - flank leakages 'c1' - 'd2'521            # - radial leakages 'c1' - 'c2' (alpha = self.alpha + 1)522            # - radial leakages 'c1' - 'c2' (alpha = self.alpha - 1)523            # - flank leakages 'c1' - 'c2' (alpha = self.alpha - 1)524            # - flank leakages 'c1' - 'c2' (alpha = self.alpha + 1)525            # - flank leakage 'c21' - 'discharge'526            # - flank leakages 'discharge' - 'c2' (alpha = N_c)527        elif N_c > 1:528            if self.name == 's1' :529                # radial leakage 's1' - 'c21' (same as c11)530                phi_ssa = geoF.get_phi_ssa(geo, self.theta)531                A_c21 = radial_area(geo, max( geo.phi_ie - self.theta - pi , phi_ssa) , geo.phi_ie - self.theta - pi)532                (self.m_dot_dict['r3'], self.h_dict['r3'], Re) = procF.isentropic_nozzle(geo, property_matrix,  A_c21,  self.rho,  self.T, self.x_l, rho2 = CVc1[1].rho,  T2 = CVc1[1].T,  x_l2 = CVc1[1].x_l, leakage_type = 'r')533                self.m_dot_dict['r3'] = corrected_mass_flow(geo, Re, self.m_dot_dict['r3'], 'r') 534                (self.m_dot_dict['r2'], self.h_dict['r2']) = (0., 0.)535                (self.m_dot_dict['r4'], self.h_dict['r4']) = (0., 0.)    536                # flank leakage 's1' - 'c21' (same as c11)537                (self.m_dot_dict['fl2'], self.h_dict['fl2'], Re) = procF.isentropic_nozzle(geo, property_matrix,  A_fl,  self.rho,  self.T, self.x_l, rho2 = CVc1[1].rho,  T2 = CVc1[1].T,  x_l2 = CVc1[1].x_l, leakage_type = 'fl')538                self.m_dot_dict['fl2'] = corrected_mass_flow(geo, Re, self.m_dot_dict['fl2'], 'fl') 539                (self.m_dot_dict['fl1'], self.h_dict['fl1']) = (0., 0.) 540541            elif self.name == 'c1' :542                if self.alpha == 1:543                    # radial leakage 'c1' - 's1'544                    phi_ssa = geoF.get_phi_ssa(geo, self.theta)545                    A_s1 = radial_area(geo, max( geo.phi_ie - self.theta - pi , phi_ssa) , geo.phi_ie - self.theta - pi)546                    (self.m_dot_dict['r2'], self.h_dict['r2'], Re) = procF.isentropic_nozzle(geo, property_matrix,  A_s1,  self.rho,  self.T, self.x_l, rho2 = CVs1.rho,  T2 = CVs1.T,  x_l2 = CVs1.x_l, leakage_type = 'r')547                    self.m_dot_dict['r2'] = corrected_mass_flow(geo, Re, self.m_dot_dict['r2'], 'r') 548                    # radial leakage 'c1' - 'sa'549                    A_sa = radial_area(geo, max( geo.phi_ie - self.theta, phi_ssa) ,phi_ssa)550                    (self.m_dot_dict['r1'], self.h_dict['r1'], Re) = procF.isentropic_nozzle(geo, property_matrix,  A_sa,  self.rho,  self.T, self.x_l, rho2 = CVs1.rho_init,  T2 = CVs1.T_init,  x_l2 = CVs1.x_l_init, leakage_type = 'r')551                    self.m_dot_dict['r1'] = corrected_mass_flow(geo, Re, self.m_dot_dict['r1'], 'r') 552                    # radial leakage 'c1' - 'c2' (alpha = self.alpha + 1)553                    A_c2 = radial_area(geo, geo.phi_ie - self.theta - 2*pi*self.alpha, geo.phi_ie - self.theta - 2*pi*self.alpha - pi)554                    (self.m_dot_dict['r4'], self.h_dict['r4'], Re) = procF.isentropic_nozzle(geo, property_matrix,  A_c2,  self.rho,  self.T, self.x_l, rho2 = CVc1[self.alpha + 1].rho,  T2 = CVc1[self.alpha + 1].T,  x_l2 = CVc1[self.alpha + 1].x_l, leakage_type = 'r')555                    self.m_dot_dict['r4'] = corrected_mass_flow(geo, Re, self.m_dot_dict['r4'], 'r') 556                    (self.m_dot_dict['r3'], self.h_dict['r3']) = (0., 0.)               557                    # flank leakage with 'c1' - 's2'558                    (self.m_dot_dict['fl1'], self.h_dict['fl1'], Re) = procF.isentropic_nozzle(geo, property_matrix,  A_fl,  self.rho,  self.T, self.x_l, rho2 = CVs1.rho,  T2 = CVs1.T,  x_l2 = CVs1.x_l, leakage_type = 'fl')559                    self.m_dot_dict['fl1'] = corrected_mass_flow(geo, Re, self.m_dot_dict['fl1'], 'fl') 560                    # flank leakage 'c1' - 'c2' (alpha = self.alpha + 1)561                    (self.m_dot_dict['fl2'], self.h_dict['fl2'], Re) = procF.isentropic_nozzle(geo, property_matrix,  A_fl,  self.rho,  self.T, self.x_l, rho2 = CVc1[self.alpha + 1].rho,  T2 = CVc1[self.alpha + 1].T,  x_l2 = CVc1[self.alpha + 1].x_l, leakage_type = 'fl')562                    self.m_dot_dict['fl2'] = corrected_mass_flow(geo, Re, self.m_dot_dict['fl2'], 'fl') 563                    564                elif self.alpha == N_c: 565                    # radial leakage 'c1' - 'c2' (alpha = self.alpha - 1)566                    A_c2 = radial_area(geo, geo.phi_ie - self.theta - 2*pi*(self.alpha - 1), geo.phi_ie - self.theta - 2*pi*(self.alpha - 1) - pi)567                    (self.m_dot_dict['r1'], self.h_dict['r1'], Re) = procF.isentropic_nozzle(geo, property_matrix,  A_c2,  self.rho,  self.T, self.x_l, rho2 = CVc1[self.alpha - 1].rho,  T2 = CVc1[self.alpha - 1].T,  x_l2 = CVc1[self.alpha - 1].x_l, leakage_type = 'r')568                    self.m_dot_dict['r1'] = corrected_mass_flow(geo, Re, self.m_dot_dict['r1'], 'r') 569                    # radial leakage 'c1' - 'd2' (alpha = N_c)570                    A_d2 = radial_area(geo, geo.phi_ie - self.theta - 2*pi*N_c, geo.phi_is)571                    (self.m_dot_dict['r4'], self.h_dict['r4'], Re) = procF.isentropic_nozzle(geo, property_matrix,  A_d2,  self.rho,  self.T, self.x_l, rho2 = CVdis.rho,  T2 = CVdis.T,  x_l2 = CVdis.x_l, leakage_type = 'r')572                    self.m_dot_dict['r4'] = corrected_mass_flow(geo, Re, self.m_dot_dict['r4'], 'r') 573                    (self.m_dot_dict['r2'], self.h_dict['r2']) = (0., 0.)   574                    (self.m_dot_dict['r3'], self.h_dict['r3']) = (0., 0.)                    575                    # flank leakages 'c1' - 'c2' (alpha = self.alpha - 1)576                    (self.m_dot_dict['fl1'], self.h_dict['fl1'], Re) = procF.isentropic_nozzle(geo, property_matrix,  A_fl,  self.rho,  self.T, self.x_l, rho2 = CVc1[self.alpha - 1].rho,  T2 = CVc1[self.alpha - 1].T,  x_l2 = CVc1[self.alpha - 1].x_l, leakage_type = 'fl')577                    self.m_dot_dict['fl1'] = corrected_mass_flow(geo, Re, self.m_dot_dict['fl1'], 'fl') 578                    # flank leakages 'c1' - 'd2'579                    (self.m_dot_dict['fl2'], self.h_dict['fl2'], Re) = procF.isentropic_nozzle(geo, property_matrix,  A_fl,  self.rho,  self.T, self.x_l, rho2 = CVdis.rho,  T2 = CVdis.T,  x_l2 = CVdis.x_l, leakage_type = 'fl')580                    self.m_dot_dict['fl2'] = corrected_mass_flow(geo, Re, self.m_dot_dict['fl2'], 'fl') 581                582                else : 583                    # radial leakages 'c1' - 'c2' (alpha = self.alpha + 1)584                    A_c2p1 = radial_area(geo, geo.phi_ie - self.theta - 2*pi*self.alpha, geo.phi_ie - self.theta - 2*pi*self.alpha - pi)585                    (self.m_dot_dict['r4'], self.h_dict['r4'], Re) = procF.isentropic_nozzle(geo, property_matrix,  A_c2p1,  self.rho,  self.T, self.x_l, rho2 = CVc1[self.alpha + 1].rho,  T2 = CVc1[self.alpha + 1].T,  x_l2 = CVc1[self.alpha + 1].x_l, leakage_type = 'r')586                    self.m_dot_dict['r4'] = corrected_mass_flow(geo, Re, self.m_dot_dict['r4'], 'r') 587                    # radial leakages 'c1' - 'c2' (alpha = self.alpha - 1)588                    A_c2m1 = radial_area(geo, geo.phi_ie - self.theta - 2*pi*(self.alpha - 1), geo.phi_ie - self.theta - 2*pi*(self.alpha - 1) - pi)589                    (self.m_dot_dict['r1'], self.h_dict['r1'], Re) = procF.isentropic_nozzle(geo, property_matrix,  A_c2m1,  self.rho,  self.T, self.x_l, rho2 = CVc1[self.alpha - 1].rho,  T2 = CVc1[self.alpha - 1].T,  x_l2 = CVc1[self.alpha - 1].x_l, leakage_type = 'r')590                    self.m_dot_dict['r1'] = corrected_mass_flow(geo, Re, self.m_dot_dict['r1'], 'r') 591                    (self.m_dot_dict['r2'], self.h_dict['r2']) = (0., 0.)   592                    (self.m_dot_dict['r3'], self.h_dict['r3']) = (0., 0.)        593                    # flank leakages 'c1' - 'c2' (alpha = self.alpha - 1)594                    (self.m_dot_dict['fl1'], self.h_dict['fl1'], Re) = procF.isentropic_nozzle(geo, property_matrix,  A_fl,  self.rho,  self.T, self.x_l, rho2 = CVc1[self.alpha - 1].rho,  T2 = CVc1[self.alpha - 1].T,  x_l2 = CVc1[self.alpha - 1].x_l, leakage_type = 'fl')595                    self.m_dot_dict['fl1'] = corrected_mass_flow(geo, Re, self.m_dot_dict['fl1'], 'fl') 596                    # flank leakages 'c1' - 'c2' (alpha = self.alpha + 1)597                    (self.m_dot_dict['fl2'], self.h_dict['fl2'], Re) = procF.isentropic_nozzle(geo, property_matrix,  A_fl,  self.rho,  self.T, self.x_l, rho2 = CVc1[self.alpha + 1].rho,  T2 = CVc1[self.alpha + 1].T,  x_l2 = CVc1[self.alpha + 1].x_l, leakage_type = 'fl')598                    self.m_dot_dict['fl2'] = corrected_mass_flow(geo, Re, self.m_dot_dict['fl2'], 'fl') 599                    600            elif self.name == name_dis :601                # radial leakages 'discharge' - 'c2' (alpha = N_c)602                A_c2 = radial_area(geo, geo.phi_ie - self.theta - 2*pi*N_c, geo.phi_is)603                (m_dot, h, Re) = procF.isentropic_nozzle(geo, property_matrix,  A_c2,  self.rho,  self.T, self.x_l, rho2 = CVc1[N_c].rho,  T2 = CVc1[N_c].T,  x_l2 = CVc1[N_c].x_l, leakage_type = 'r')604                m_dot = corrected_mass_flow(geo, Re, m_dot, 'r') 605                (self.m_dot_dict['r1'], self.h_dict['r1']) = (coef*m_dot, h)606                (self.m_dot_dict['r2'], self.h_dict['r2']) = (0., 0.)607                (self.m_dot_dict['r3'], self.h_dict['r3']) = (0., 0.)608                (self.m_dot_dict['r4'], self.h_dict['r4']) = (0., 0.)              609                # flank leakages 'discharge' - 'c2' (alpha = N_c)610                (m_dot, h, Re) = procF.isentropic_nozzle(geo, property_matrix,  A_fl,  self.rho,  self.T, self.x_l, rho2 = CVc1[N_c].rho,  T2 = CVc1[N_c].T,  x_l2 = CVc1[N_c].x_l, leakage_type = 'fl')611                m_dot = corrected_mass_flow(geo, Re, m_dot, 'fl') 612                (self.m_dot_dict['fl1'], self.h_dict['fl1']) = (coef*m_dot, h)613                (self.m_dot_dict['fl2'], self.h_dict['fl2']) = (0., 0.)      614        615            else :616                (self.m_dot_dict['r2'], self.h_dict['r2']) = (0., 0.)617                (self.m_dot_dict['r2'], self.h_dict['r2']) = (0., 0.)618                (self.m_dot_dict['r3'], self.h_dict['r3']) = (0., 0.)619                (self.m_dot_dict['r4'], self.h_dict['r4']) = (0., 0.) 620                (self.m_dot_dict['fl1'], self.h_dict['fl1']) = (0., 0.)621                (self.m_dot_dict['fl2'], self.h_dict['fl2']) = (0., 0.)             622623624# The isentropic nozzle model effectively captures the compressibility effects, although it does not capture the frictional 625# effects. The ratio of the isentropic nozzle mass flow rate prediction to that of the detailed model can be correlated626# with the Reynolds number. The methodology and the equations used for this correlation are thouroughly explained627# in appendix C of Ian Bell's thesis.628# The 'corrected_mass_flow' function takes as input the machine geometry, Reynolds number and leakage type. It provides629# as output the ratio of the isentropic nozzle mass flow rate prediction to that of the more detailed model.630def corrected_mass_flow(geo, Re, m_dot, leakage_type):631    632    # Non-dimensionalization parameters633    delta_0 = 10e-6634    L_0 = 0.005635    636    # Different empyrical coefficient are used if the leakage flow is radial or flank637    if leakage_type == 'r':638        a0 =  2.59321070e+4639        a1 =  9.14825434e-1640        a2 = -1.77588568e+2641        a3 = -2.37052788e-1642        a4 = -1.72347611e+05643        a5 = -1.20687600e+01644        a6 = -1.28861161e-02645        a7 = -1.51202604e+02646        a8 = -9.99674458e-01647        a9 =  1.61435039e-02648        a10 = 8.25533457e-01649        Re_star = 5.24358195e+03650        L = geo.t_s*1e-3651        delta = geo.delta_r     652    elif leakage_type == 'fl':   653        a0 = -2.63970396e+00654        a1 = -5.67164431e-01655        a2 = 8.36554999e-01656        a3 = 8.10567168e-01657        a4 = 6.17402826e+03658        a5 = -7.60907962e+00659        a6 = -5.10200923e-01660        a7 = -1.20517483e+03661        a8 = -1.02938914e+00662        a9 = 6.89497786e-01663        a10 = 1.09607735e+00664        Re_star = 8.26167178e+02665        L = geo.r_o*1e-3666        delta = geo.delta_f667    else :668        raise Exception('Invalid inputs')669    670    # If Reynolds number is equal to 0 no correction is applied to the mass flow671    if Re == 0:672        M = 1673    else:674        zeta = 1/(1 + exp(-0.01*(Re - Re_star)))   675        M = a0*(L/L_0)**a1 / (a2*(delta/delta_0) + a3) * (zeta*(a4*Re**a5 + a6) + (1-zeta)*(a7*Re**a8 + a9)) + a10676    677    return m_dot/M    678        679def radial_area(geo, phi_max, phi_min): 680    681    phi_0 = (geo.phi_i0 + geo.phi_o0)/2682    A_radial = geo.delta_r*geo.r_b/1e3* (0.5*(phi_max**2 - phi_min**2) - phi_0*(phi_max - phi_min))683    return A_radial684685def flank_area(geo):686    A_flank = geo.h_s/1e3*geo.delta_f687    688    return A_flank689690# The 'merging_process' function takes as input the discharge chambers objects 'CV_ddd', 'CV_dd', 'CV_d1' and691# 'CV_d2'. Its function is to compute the properties of the resulting chamber after the discharge chambers 692# 'd1', 'd2' and 'dd' have merged into chamber 'ddd'.693def merging_process(property_matrix, CV_ddd, CV_dd, CV_d1, CV_d2):694    695    # Energy of the 'ddd' chamber is set equal to the energy of the 'dd' chamber.696    CV_ddd.mh_sum = CV_dd.mh_sum 697    # Mass of the 'ddd' chamber is set equal to the mass of the 'dd' chamber. 698    CV_ddd.m_sum = CV_dd.m_sum699    # Discharge enthalpy of chamber 'ddd' is sert equal to discharge enthalpy of chamber 'dd'700    CV_ddd.h_disc = CV_dd.h_disc701    # Crank angle of chamber 'ddd' is set equal to crank angle of chambr 'dd' 702    CV_ddd.theta = CV_dd.theta703    704    # Volume and mass of chamber 'ddd' is equal to the sum of volumes and masses of chambers 'd1', 'd2' and 'dd'705    CV_ddd.V = CV_dd.V + CV_d1.V + CV_d2.V706    CV_ddd.m = CV_dd.m + CV_d1.m + CV_d2.m707 708    # Pressure of chamber 'ddd' is computed as the weigthed mean of pressure in chambers 'd1', 'd2' and 'dd'709    CV_ddd.P = (CV_dd.V * CV_dd.P + CV_d1.V * CV_d1.P + CV_d2.V * CV_d2.P)/CV_ddd.V710    # Quality of chamber 'ddd' is computed as the weigthed mean of quality in chambers 'd1', 'd2' and 'dd'711    CV_ddd.x_l = (CV_dd.m * CV_dd.x_l + CV_d1.m  * CV_d1.x_l + CV_d2.m  * CV_d2.x_l)/CV_ddd.m712    # Enthalpy of chamber 'ddd' is computed as the weigthed mean of enthalpy in chambers 'd1', 'd2' and 'dd'713    CV_ddd.h = (CV_dd.m * CV_dd.h + CV_d1.m  * CV_d1.h + CV_d2.m  * CV_d2.h)/CV_ddd.m714    #CV_ddd.rho = CV_ddd.m/CV_ddd.V715    716    # Temperature of chamber 'ddd' is computed using enthalpy and pressure717    CV_ddd.T = propF.propfunction_mixture(property_matrix, CV_ddd.x_l, 'T',  h = CV_ddd.h , P = CV_ddd.P)718    #CV_ddd.T = propfunction_mixture(property_matrix, CV_ddd.x_l, 'T',  rho = CV_ddd.rho , P = CV_ddd.P)719    # Temperature of chamber 'ddd' is computed using enthalpy and pressure720    CV_ddd.rho = propF.propfunction_mixture(property_matrix, CV_ddd.x_l, 'rho', h = CV_ddd.h , P = CV_ddd.P)721    722    # Properties of chamber 'ddd' are computed (AGAIN???) using the 'update' method723    CV_ddd.update(property_matrix)724    725726# The 'stopping_criteria' function verifies whether convergence is attained at the end of a complete rotation of the727# crank angle. 728def stopping_criteria(geo, CV_c1, CV_ddd, tol)  :   729    730    L = len(CV_ddd.dict_init['T']) - 1 731    732    ddd_cond_T = abs(CV_ddd.dict_init['T'][L] \733        - CV_ddd.dict_init['T'][L-1])/CV_ddd.dict_init['T'][L] <= tol734    735    ddd_cond_rho = abs(CV_ddd.dict_init['rho'][L] \736        - CV_ddd.dict_init['rho'][L-1])/CV_ddd.dict_init['rho'][L] <= tol737    738    alpha = 1739    c1_cond_T = {}740    c1_cond_rho = {}741    742    stop = True and ddd_cond_T and ddd_cond_rho743    744    while alpha <= geo.N_c_max:745        746        c1_cond_T[alpha] = abs(CV_c1[alpha].dict_init['T'][L] \747            - CV_c1[alpha].dict_init['T'][L-1])/CV_c1[alpha].dict_init['T'][L] <= tol748        c1_cond_rho[alpha] = abs(CV_c1[alpha].dict_init['rho'][L] \749            - CV_c1[alpha].dict_init['rho'][L-1])/CV_c1[alpha].dict_init['rho'][L] <= tol    750751        stop = stop and  c1_cond_T[alpha] and c1_cond_rho[alpha]752        753        alpha += 1754        755    stop = stop == False756    return stop757    758 759    760# The function 'performance_computation' is responsible for computing mass flow rate through the machine, power761# absorbed, isentropic efficiency and volumetric efficiency. It only requires the properties of the 'ddd' and 's1'762# chambers, outlet pressure and mechanical power losses (expressed by means of the torque 'tau_loss').763def performance_computation(geo, property_matrix, CV_ddd, CV_s1, P_out, omega, tau_loss): 764    765    # The mass flow rate coming into the compressor is expressed as the mass flow rate passing through chamber 's1'766    # plus the mass flow rate passing through chamber 's2'. Since these 2 chambers are modeled as identical only the first one767    # is considered and then it is multiplied by 2.768    m_dot_in = CV_s1.m_sum*2769    770    # The mass flow rate coming out of the machine coincide with the mass flow rate through chamber 'ddd'771    m_dot_out = -CV_ddd.m_sum772    773    # Overall mass flow rate is computed as the arithmetic mean mass flow rates in and out of the machine.774    m_dot = (m_dot_in + m_dot_out)/2775    776    # Volumetric efficiency is computed as the ratio between the actual mass flow rate and maximum volumetric777    # efficiency that could be attained by the machine. The maximum mass flow rate is computed as the volume778    # displaced by the machine in one rotation (stored in the 'geo' object as 'V_disp') times the initial 779    # density of the fluid (density at inlet)780    eta_v = m_dot/(omega/2/pi*geo.V_disp*CV_s1.rho_init)781    h_out_is = propF.propfunction_mixture(property_matrix, CV_ddd.x_l, 'h', s = CV_s1.s_init , P = P_out)782    783    # Power loss is computed as torque loss times angular velocity784    W_dot_loss = omega*tau_loss785    # Useful power is computed as mass flow rate times enthalpy increase between inlet and outlet786    W_dot_i = m_dot*(CV_ddd.h_disc - CV_s1.h_init)787    W_dot = (W_dot_i + W_dot_loss)788    789    # Isentropic efficiency is computed as the ratio between the difference of inlet enthalpy and outlet790    # enthalpy under the assumption of isentropic process, and the actual enthalpy difference between inlet and outlet.791    # The actual enthalpy difference between inlet and outlet is computes as the power absorbed by the machine divided by792    # the mass flow rate processed.793    eta_is = ( h_out_is - CV_s1.h_init )*m_dot/W_dot794    795    print('----------------------------------')796    print('m_dot_in:' ,str(m_dot_in), 'kg/s')797    print('m_dot_out:' ,str(m_dot_out), 'kg/s')798    print('eta_is:' ,str(eta_is), '-')799    print('eta_v:' ,str(eta_v), '-')800    print('Mechanical losses:' ,str(W_dot_loss), 'W')801    print('Power consumed:' ,str(W_dot), 'W')802    803    return m_dot, eta_is, eta_v, W_dot_loss, W_dot   804805    806# The 'convert_recorded_values' function use the attributes from all the control volumes to generate dictionaries807# that are used within the 'diversePlots' module to generate the plots for the post process tool.808def convert_recorded_values(geo, CV_s1, CV_c1, CV_d1, CV_dd, CV_ddd):809    810    # The 'm_dot_vect_final' dictionary contains the recorded value of the mass flow rate in all the chambers.811    # These values are now converted from lists to numpy arrays.812    for k in CV_s1.m_dot_vect_final: 813        # Suction814        CV_s1.m_dot_vect_final[k] = np.array(CV_s1.m_dot_vect_final[k])  815        # Discharge 1816        CV_d1.m_dot_vect_final[k] = np.array(CV_d1.m_dot_vect_final[k]) 817        # Discharge 'dd'818        CV_dd.m_dot_vect_final[k] = np.array(CV_dd.m_dot_vect_final[k]) 819        # Discharge 'ddd'820        CV_ddd.m_dot_vect_final[k] = np.array(CV_ddd.m_dot_vect_final[k]) 821        # Compression822        alpha = 1823        while alpha <= geo.N_c_max:824            CV_c1[alpha].m_dot_vect_final[k] = np.array(CV_c1[alpha].m_dot_vect_final[k])825            alpha += 1826    827    # Every other attributes of the 's1' chamber is converted into a numpy array828    vect_s1 = {'P' : np.array(CV_s1.P_vect_final), 'T' : np.array(CV_s1.T_vect_final), 829               'Q': np.array(CV_s1.Q_vect_final) , 's' : np.array(CV_s1.s_vect_final), 830               'm' : np.array(CV_s1.m_vect_final), 'x_l' : np.array(CV_s1.x_l_vect_final),831               'rho' : np.array(CV_s1.rho_vect_final),'theta' : np.array(CV_s1.theta_vect_final),832               'm_dot' : CV_s1.m_dot_vect_final, 'V': np.array(CV_s1.V_vect_final)}833    834    # Every other attributes of the 'c1' chambers is converted into a numpy array835    alpha = 1836    vect_c1 = {}837    while alpha <= geo.N_c_max:838        vect_c1[alpha] = {'P' : np.array(CV_c1[alpha].P_vect_final), 'T' : np.array(CV_c1[alpha].T_vect_final), 839                   'Q': np.array(CV_c1[alpha].Q_vect_final) , 's' : np.array(CV_c1[alpha].s_vect_final), 840                   'm' : np.array(CV_c1[alpha].m_vect_final), 'x_l' : np.array(CV_c1[alpha].x_l_vect_final), 841                   'rho' : np.array(CV_c1[alpha].rho_vect_final),'theta' : np.array(CV_c1[alpha].theta_vect_final),842                   'm_dot' : CV_c1[alpha].m_dot_vect_final, 'V': np.array(CV_c1[alpha].V_vect_final)}843        alpha += 1844    845    # Every other attributes of the 'd1' chamber is converted into a numpy array846    vect_d1 = {'P' : np.array(CV_d1.P_vect_final), 'T' : np.array(CV_d1.T_vect_final), 847               'Q': np.array(CV_d1.Q_vect_final) , 's' : np.array(CV_d1.s_vect_final), 848               'm' : np.array(CV_d1.m_vect_final), 'x_l' : np.array(CV_d1.x_l_vect_final), 849               'rho' : np.array(CV_d1.rho_vect_final),'theta' : np.array(CV_d1.theta_vect_final),850               'm_dot' : CV_d1.m_dot_vect_final, 'V': np.array(CV_d1.V_vect_final)}851    852    # Every other attributes of the 'dd' chamber is converted into a numpy array853    vect_dd = {'P' : np.array(CV_dd.P_vect_final), 'T' : np.array(CV_dd.T_vect_final), 854               'Q': np.array(CV_dd.Q_vect_final) , 's' : np.array(CV_dd.s_vect_final), 855               'm' : np.array(CV_dd.m_vect_final), 'x_l' : np.array(CV_dd.x_l_vect_final), 856               'rho' : np.array(CV_dd.rho_vect_final),'theta' : np.array(CV_dd.theta_vect_final),857               'm_dot' : CV_dd.m_dot_vect_final, 'V': np.array(CV_dd.V_vect_final)}858    859    # Every other attributes of the 'ddd' chamber is converted into a numpy array860    vect_ddd = {'P' : np.array(CV_ddd.P_vect_final), 'T' : np.array(CV_ddd.T_vect_final), 861               'Q': np.array(CV_ddd.Q_vect_final) , 's' : np.array(CV_ddd.s_vect_final), 862               'm' : np.array(CV_ddd.m_vect_final), 'x_l' : np.array(CV_ddd.x_l_vect_final), 863               'rho' : np.array(CV_ddd.rho_vect_final),'theta' : np.array(CV_ddd.theta_vect_final),864               'm_dot' : CV_ddd.m_dot_vect_final, 'V': np.array(CV_ddd.V_vect_final)}865866    867    # The 'i' index iterates over the crank angle 'theta'868    for i in range(len(vect_s1['theta'])):869        # The items() method returns a view object. The view object contains the key-value pairs of the dictionary, 870        # as tuples in a list.871        for k, val in vect_s1.items():872            # This if statement checks return true at the angles in which the discharge chambers are merged, therefore873            # when the 'ddd' chamber is defined and the 'd1'-'dd' chambers are not.874            if ((geo.theta_d >= geo.theta_m and (vect_s1['theta'][i] > geo.theta_m and vect_s1['theta'][i] <= geo.theta_d))  875                    or (geo.theta_m >= geo.theta_d and (vect_s1['theta'][i] > geo.theta_m or vect_s1['theta'][i] <= geo.theta_d))) and k != 'theta' :876                877                # Since the 'dd'-'d1' chambers are not defined at this value of theta (chambers are merged into 'ddd') their878                # properties are substituted with 'nan' and will not be displayed in the final plot.879                if k == 'm_dot':880                    for n in vect_s1[k]:881                        vect_d1[k][n][i] = np.nan882                        vect_dd[k][n][i] = np.nan883                else:      884                    vect_d1[k][i] = np.nan885                    vect_dd[k][i] = np.nan886                887            # Since the if statements returns false, it means that at this value of theta the 'ddd' chamber is not defined888            # and the 'dd'-'d1' chambers are not merged. The 'ddd' chamber properties are substituted with 'nan'.889            elif  k != 'theta':890                if k == 'm_dot':891                     for n  in vect_d1[k]:892                         vect_ddd[k][n][i] = np.nan893                else:894                    vect_ddd[k][i] = np.nan895            896            # The innermost compression chamber becomes the discharge chamber right after the discharge angle. The 897            # properties of this chamber are substituted with 'nan' and will not be displayed in the final plot898            if vect_s1['theta'][i] > geo.theta_d and k != 'theta' :899                 if k == 'm_dot':900                     for n in vect_s1[k]:901                         vect_c1[geo.N_c_max][k][n][i] = np.nan902                 else:903                     vect_c1[geo.N_c_max][k][i] = np.nan904                     905    return vect_s1, vect_c1, vect_d1, vect_dd, vect_ddd   906907908# The 'leakage_computation' function wraps up the leakage flow that has to be modeled and passes the correct 909# inputs to the correct function.910def leakage_computation(property_matrix, geo, exist, CV_s1, CV_c1, CV_d1, CV_ddd, leakage):911    912    if leakage == True:913        CV_s1.leakages(property_matrix, geo, exist, CVc1 = CV_c1, CVd1 = CV_d1, CVddd = CV_ddd)914        CV_d1.leakages(property_matrix, geo, exist, CVs1 = CV_s1, CVc1 = CV_c1, CVddd = CV_ddd)915        CV_ddd.leakages(property_matrix, geo, exist, CVs1 = CV_s1, CVc1 = CV_c1, CVddd = CV_ddd)916        alpha = 1917        while alpha <= geoF.N_c_calculation(geo,CV_s1.theta):918            CV_c1[alpha].leakages(property_matrix, geo, exist, CVs1 = CV_s1, CVc1 = CV_c1, CVd1 = CV_d1, CVddd = CV_ddd)919            alpha += 1920921    else:922        pass923924925926# The 'full_resolution' function is called by the main script. It uses as input the geometry of the machine ('geo', 'dictV',927# 'dict_dV', 'dict_area') and the known properties ('T_in', 'rho_in', 'x_l_in', 'P_out') to compute mass flow rate,928# isentropic efficiency, volumetric efficiency and absorbed power.929def full_resolution(geo, property_matrix, dict_V, dict_dV, dict_area, T_in, rho_in, x_l_in, P_out, tol_glob, 930                    tol_pressure, theta_step_base, omega, tau_loss, leakage, heat_transfer, Reed_valve):931    932    # Suction chambers control volumes are generated933    print('----------------------------------')934    print('Starting new simulation')935    CV_s1 = controlVolume("s1", geo, property_matrix, dict_V, dict_dV, dict_area, T_in, rho_in, x_l_in)936    CV_sa = controlVolume("sa", geo, property_matrix, dict_V, dict_dV, dict_area, T_in, rho_in, x_l_in)937    938    # Initialization of the properties within the suction chambers control volumes939    CV_s1.initialize()940    CV_sa.initialize()941    942    # Generation of the compression chambers control volumes. In general, there are more than a pair of compression943    # chambers in a scroll machine, depending on the geometry parameters that are specfied. The maximum number944    # of compression chamber is specified in the 'geo' object as the attribute 'N_c_max'.945    alpha = 1946    CV_c1 = {}947    # Compression chambers are generated until 'alpha' becomes higher that 'N_c_max'. These chambers are placed in948    # a list and number thanks to the index 'alpha'949    while alpha <= geo.N_c_max:950        # generation of the alpha-th compression chamber951        CV_c1[alpha] = controlVolume("c1", geo, property_matrix, dict_V, dict_dV, dict_area, T_in, rho_in, x_l_in, alpha = alpha) 952        # Initialization of the properties within the alpha-th compression chambers control volumes953        CV_c1[alpha].initialize()954        alpha += 1955    956    # Discharge chambers control volumes are generated957    CV_ddd = controlVolume("ddd", geo, property_matrix, dict_V, dict_dV, dict_area, T_in, rho_in, x_l_in, P_out) 958    CV_d1 = controlVolume("d1", geo, property_matrix, dict_V, dict_dV, dict_area, T_in, rho_in, x_l_in, P_out) 959    CV_dd = controlVolume("dd", geo, property_matrix, dict_V, dict_dV, dict_area, T_in, rho_in, x_l_in, P_out)960    961    # Initialization of properties within the discharge chambers control columes962    CV_ddd.initialize()963    CV_d1.initialize()964    CV_dd.initialize()965    966    # The iterative procedure begins here. The stop criterion for the while loop are defined within the 967    # 'stopping_criteria' function.968    while stopping_criteria(geo, CV_c1, CV_ddd, tol_glob):969        970        # Beginning of the rotation971        theta = 0972        # Definition of the step for the crank angles. Initially the 'theta_step_base' is used. It is updated if necessary973        theta_step = theta_step_base974        975        # Inside this while loop the process is simulated from theta = 0 to theta = discharge_angle976        while theta < geo.theta_d - theta_step:977        978            # If leakage computation is turn on, leakages are evaluated. 979            leakage_computation(property_matrix, geo, 'ddd', CV_s1, CV_c1, CV_d1, CV_ddd, leakage)980            # print('s1',CV_s1.h_dict['r3'],CV_c1[1].h_dict['r2'])981            # print('c1',CV_c1[1].h_dict['r4'], CV_c1[2].h_dict['r1'])982            # print('ddd',CV_c1[2].h_dict['r4'], CV_ddd.h_dict['r1']/2)983            984            # Mass flow rates through 'ddd' chamber and 's1' chambers are evaluated using their respective method.985            CV_s1.mass_flow_su(property_matrix, geo, rho_in, T_in, x_l_in, theta_step)986            CV_ddd.mass_flow_ex(property_matrix, geo, P_out, theta_step, Reed_valve) 987988            # The 'next_values' method is called for the 's1' chamber. The method updates property attributes based989            # on mass and temperature derivatives computed through mass and energy balances.990            CV_s1.next_values(property_matrix, omega, theta_step)991            alpha = 1992            # The 'next_values' method is called for the compression chambers (in number of N_c_max). The method updates property attributes based993            # on mass and temperature derivatives computed through mass and energy balances.994            while alpha <= geo.N_c_max:995                CV_c1[alpha].next_values(property_matrix, omega, theta_step)996                alpha += 1997            # The 'next_values' method is called for the 'ddd' chamber. The method updates property attributes based998            # on mass and temperature derivatives computed through mass and energy balances.999            CV_ddd.next_values(property_matrix, omega, theta_step)1000            1001            # Crank angle is updated.1002            theta += theta_step1003            1004        # The 'discharge' method is applied to the discharge chambers 'd1' and 'dd' after the discharge angle is overtaken.1005        CV_d1.discharge(CV_c1[geo.N_c_max])1006        CV_dd.discharge(CV_ddd)1007        1008        # The merging process of the discharge chambers is simulated within this while loop. This part of the 1009        # simulation goes on until pressure difference of chambers 'd1' and 'dd' is lower than 'tol_pressure',1010        # the crank angle must also reach the value (discharge_angle + 0.66pi).1011        while abs((CV_d1.P - CV_dd.P)/CV_d1.P) > tol_pressure and theta < geo.theta_d + 2*pi/3:1012        1013            # If leakage computation is turned on, leakages are evaluated.1014            leakage_computation(property_matrix, geo, 'd', CV_s1, CV_c1, CV_d1, CV_ddd, leakage)1015            # print('s1',CV_s1.h_dict['r3'],CV_c1[1].h_dict['r2'])1016            # print('d1',CV_c1[1].h_dict['r4'], CV_d1.h_dict['r1'])1017          1018            # Mass flow rates through 'ddd' chamber and 's1' chambers are evaluated using their respective method.1019            CV_s1.mass_flow_su(property_matrix, geo, rho_in, T_in, x_l_in, theta_step)1020            CV_dd.mass_flow_ex(property_matrix, geo, P_out, theta_step, Reed_valve)1021            1022            # Mass flow rate through the discharge port is evaluated1023            CV_dd.mass_flow_dis_dd(property_matrix, geo, CV_d1)1024            1025            # The 'next_values' method is called for the 's1' chamber. The method updates property attributes based1026            # on mass and temperature derivatives computed through mass and energy balances.1027            CV_s1.next_values(property_matrix, omega, theta_step)1028            alpha = 11029            # The 'next_values' method is called for the compression chambers (in number of N_c_max). The method updates property attributes based1030            # on mass and temperature derivatives computed through mass and energy balances.1031            while alpha <= geo.N_c_max-1:1032                CV_c1[alpha].next_values(property_matrix, omega, theta_step)1033                alpha += 11034            # The 'next_values' method is called for the 'dd' and 'd1' chambers. The method updates property attributes based1035            # on mass and temperature derivatives computed through mass and energy balances.1036            CV_dd.next_values(property_matrix, omega, theta_step)1037            CV_d1.next_values(property_matrix, omega, theta_step)1038            1039            # Crank angle is updated1040            theta += theta_step1041            1042        # The 'merging_process' function is called right after merging conditions are met. This function computes the 1043        # properties of the fluid within the 'ddd' chambers that is created after chambers 'd1', 'd2' and 'dd' have 1044        # merged.1045        merging_process(property_matrix, CV_ddd, CV_dd, CV_d1, CV_d1)1046        # After the merging angle, the shift angle can be set to 0.1047        CV_ddd.theta_shift = 01048        1049        # In this while loop the remaining of the crank rotation is simulated1050        while theta < 2*pi:1051            1052            # If leakage computation is turn on, leakages are evaluated.1053            leakage_computation(property_matrix, geo, 'ddd', CV_s1, CV_c1, CV_d1, CV_ddd, leakage)1054            # print('s1',CV_s1.h_dict['r3'], CV_c1[1].h_dict['r2'])1055            # print('ddd',CV_c1[1].h_dict['r4'], CV_ddd.h_dict['r1']/2)1056           1057            # Mass flow rates through 'ddd' chamber and 's1' chambers are evaluated using their respective method.1058            CV_s1.mass_flow_su(property_matrix, geo, rho_in, T_in, x_l_in, theta_step)1059            CV_ddd.mass_flow_ex(property_matrix, geo, P_out, theta_step, Reed_valve) 1060            1061            # The 'next_values' method is called for the 's1' chamber. The method updates property attributes based1062            # on mass and temperature derivatives computed through mass and energy balances.1063            CV_s1.next_values(property_matrix, omega, theta_step)1064            alpha = 11065            # The 'next_values' method is called for the compression chambers (in number of N_c_max). The method updates property attributes based1066            # on mass and temperature derivatives computed through mass and energy balances.1067            while alpha <= geo.N_c_max-1:1068                CV_c1[alpha].next_values(property_matrix, omega, theta_step)1069                alpha += 11070            # The 'next_values' method is called for the 'ddd' chamber. The method updates property attributes based1071            # on mass and temperature derivatives computed through mass and energy balances.1072            CV_ddd.next_values(property_matrix, omega, theta_step)1073            1074            # Update of the crank angle1075            theta += theta_step1076        1077        # Mass flow rate, isentropic efficiency, volumetric efficiency and absorbed power are computed after a1078        # complete rotation has been simulated. 1079        (m_dot, eta_is, eta_v, W_dot_loss, W_dot) = performance_computation(geo, property_matrix, CV_ddd, CV_s1, P_out, omega, tau_loss)1080        1081        # The shift angle used for the computation of the 'ddd' volume is set to its initial value of 2*pi.1082        CV_ddd.theta_shift = 2*pi1083        # The initial values of properties in the 'ddd' chambers for the next rotations are set equal to the value1084        # of such properties at theta = 2*pi1085        CV_ddd.initialize(CV_ddd)1086        1087        # The initial values of properties in the innermost compression chambers for the next rotations are set equal to the value1088        # of such properties at theta = 2*pi1089        alpha = geo.N_c_max1090        while alpha > 1:1091            CV_c1[alpha].initialize(CV_c1[alpha-1])1092            alpha -= 11093        # The initial values of properties in the innermost compression chambers for the next rotations are set equal to1094        # the value of these properties in the suction chamber 's1' at theta = 2*pi.1095        CV_c1[1].initialize(CV_s1)1096        # The initial values of properties in the suction chambers for the next rotations are set equal to1097        # the value of these properties in the suction region 'sa' at theta = 2*pi. 1098        CV_s1.initialize(CV_sa)1099    1100        1101    # A final rotation is done after the convergence criteria are respected, the simulation could have been stopped1102    # here but this rotation is done in order to save all the values of the properties for every chamber as a function1103    # of the orbiting angle. The 'record_values' method is called after properties are computed in every chamber.1104    1105    theta = 01106    theta_step = theta_step_base1107    while theta < geo.theta_d - theta_step:1108    1109        leakage_computation(property_matrix, geo, 'ddd', CV_s1, CV_c1, CV_d1, CV_ddd, leakage)1110        CV_s1.mass_flow_su(property_matrix, geo, rho_in, T_in, x_l_in, theta_step)1111        CV_ddd.mass_flow_ex(property_matrix, geo, P_out, theta_step, Reed_valve)         1112        CV_s1.record_values()1113        1114        alpha = 11115        while alpha <= geo.N_c_max:1116            CV_c1[alpha].record_values()1117            alpha += 11118        1119        CV_ddd.record_values()1120        CV_d1.record_values()1121        CV_dd.record_values()1122        1123        CV_s1.next_values(property_matrix, omega, theta_step)11241125        alpha = 11126        while alpha <= geo.N_c_max:1127            CV_c1[alpha].next_values(property_matrix, omega, theta_step)1128            alpha += 11129        CV_ddd.next_values(property_matrix, omega, theta_step)1130              1131        theta += theta_step1132        1133    CV_d1.discharge(CV_c1[geo.N_c_max])1134    CV_dd.discharge(CV_ddd)1135    1136    while abs((CV_d1.P - CV_dd.P)/CV_d1.P) > tol_pressure and theta < geo.theta_d + 2*pi/3:1137    1138        leakage_computation(property_matrix, geo, 'd', CV_s1, CV_c1, CV_d1, CV_ddd, leakage)1139    1140        CV_s1.mass_flow_su(property_matrix, geo, rho_in, T_in, x_l_in, theta_step)1141        CV_dd.mass_flow_ex(property_matrix, geo, P_out, theta_step, Reed_valve)1142        CV_dd.mass_flow_dis_dd(property_matrix, geo, CV_d1)1143              1144        CV_s1.record_values()1145        alpha = 11146        while alpha <= geo.N_c_max:1147            CV_c1[alpha].record_values()1148            alpha += 11149        CV_ddd.record_values()1150        CV_d1.record_values()1151        CV_dd.record_values()1152        1153        CV_s1.next_values(property_matrix, omega, theta_step)    1154        alpha = 11155        while alpha <= geo.N_c_max-1:1156            CV_c1[alpha].next_values(property_matrix, omega, theta_step)  1157            alpha += 11158        CV_dd.next_values(property_matrix, omega, theta_step)1159        CV_d1.next_values(property_matrix, omega, theta_step)1160        1161        theta += theta_step1162    1163    geo.theta_m = theta1164    merging_process(property_matrix, CV_ddd, CV_dd, CV_d1, CV_d1)1165    CV_ddd.theta_shift = 01166            1167    while theta < 2*pi:1168        1169        leakage_computation(property_matrix, geo, 'ddd', CV_s1, CV_c1, CV_d1, CV_ddd, leakage)1170    1171        CV_s1.mass_flow_su(property_matrix, geo, rho_in, T_in, x_l_in, theta_step)1172        CV_ddd.mass_flow_ex(property_matrix, geo, P_out, theta_step, Reed_valve)  1173 1174        CV_s1.record_values()1175        alpha = 11176        while alpha <= geo.N_c_max:1177            CV_c1[alpha].record_values()1178            alpha += 11179        CV_ddd.record_values()1180        CV_d1.record_values()1181        CV_dd.record_values()1182        1183        CV_s1.next_values(property_matrix, omega, theta_step)     1184        alpha = 1  1185        while alpha <= geo.N_c_max-1:1186            CV_c1[alpha].next_values(property_matrix, omega, theta_step)1187            alpha += 11188        CV_ddd.next_values(property_matrix, omega, theta_step)1189        1190        theta += theta_step11911192    (m_dot_final, eta_is_final, eta_v_final, W_dot_loss_final, W_dot_final) =  \1193    performance_computation(geo, property_matrix, CV_ddd, CV_s1, P_out, omega, tau_loss)1194    print('----------------------------------')1195    print('Simulation finished')1196    1197    # 'CV' is a dictionary where al the recorded value 1198    CV = {'s1' : CV_s1, 'c1' : CV_c1, 'd1' : CV_d1, 'dd' : CV_dd, 'ddd' : CV_ddd }1199    1200    return (m_dot_final, eta_is_final, eta_v_final, W_dot_final, CV)12011202120312041205def save_results(name, results):1206    # 'list_file' specifies the precise file directory1207    list_file = ['C:\\Users\\Utente1\\Documents\\Tifeo\\Python\\Compressore\\Results' , name ,'.pkl']1208    1209    # The complete file name is built up1210    filename = "".join(list_file)       1211    1212    # 'file' points at the file in which the pickled object will be written1213    with open(filename, "wb") as file:1214        # The dump() method of the pickle module in Python, converts a Python object hierarchy into a byte stream. 1215        # This process is also called as serilaization.1216        pickle.dump(results, file)1217  121812191220def open_results(name):1221    # 'list_file' specifies the precise file directory1222    list_file = ['C:\\Users\\Utente1\\Documents\\Tifeo\\Python\\Compressore\\Results' , name ,'.pkl']1223    # The complete file name is built up1224    filename = "".join(list_file)       1225    1226    # 'file' points at the file to be opened1227    with open(filename, "rb") as file:1228        # The load() method of Python pickle module reads the pickled byte stream of one or more python objects 1229        # from a file object1230        output = pickle.load(file)1231     1232    return output1233    12341235def file_name(filename):1236    list_path = ['C:\\Users\\Utente1\\Documents\\Tifeo\\Python\\Compressore\\Results\\' , filename ,'.pkl']1237    1238    path = "".join(list_path) 1239    i = 21240    while os.path.isfile(path):1241        1242        filename_list = [ filename , '_' ,str(i)]1243        filename = "".join(filename_list) 1244        1245        list_path = ['C:\\Users\\Utente1\\Documents\\Tifeo\\Python\\Compressore\\Results\\' , filename ,'.pkl']1246        1247        path = "".join(list_path) 1248        i += 11249        1250    return filename12511252
...

processFunctions.py

Source:processFunctions.py

1# -*- coding: utf-8 -*-2"""3Created on Wed Mar 24 11:23:04 202145@author: nicol6"""78import numpy as np9import propertyFunctions as propF1011# This function resolves the energy balance equation given the rotational speed and the dictionaries in which12# mass flow in and out of a given control volume are specified.13def dmCV_dTheta(omega, m_dot_dict):14    """15    Parameters16    ----------17    omega : float.18    m_dot : dict.1920    Returns21    -------22    dmCVdTheta.23    """24    # Initialization of the net mass flow rate25    m_dot_sum =  0 26    27    # Computation of the net mass flow rate28    for i, j in m_dot_dict.items():29        m_dot_sum += m_dot_dict[i] 30    31    # The rate of change of mass inside a control volume with respect to the angle theta is obtained as the ratio32    # between the net mass flow rate and the rotation speed of the machine.33    dmCVdTheta = 1/omega*m_dot_sum34    35    return dmCVdTheta36    37# This function resolves the oil mass balance to compute the rate change of oil mass fraction withing a contro volume,38# given rotational speed, mass inside the control volume, the dictionaries in which mass flow in and out of a given39# control volume are specified and oil mass fraction inside a control volume.40def dxoil_dTheta(omega, m_cv, dmCVdTheta, m_dot_dict, x_l_dict, x_l):41    42    # Initialization of the oil net mass flow rate43    m_dot_x_l_sum = 0 4445    # Computation of the oil net mass flow rate46    for i, j in m_dot_dict.items():47        m_dot_x_l_sum += m_dot_dict[i]* x_l_dict[i]     48   49    # The rate of change of oil mass fraction is obtained as the difference between rate of change of oil mass due to50    # oil mass flow rate and overall mass flow rate, divided by the mass inside the control volume.51    dxoildTheta = 1/m_cv*(1/omega*m_dot_x_l_sum - x_l*dmCVdTheta)5253    return dxoildTheta5455# Temperature derivatives with respect to the crank angle is computed by solving the energy balance inside56# the control volume. 57def dT_dTheta(property_matrix, omega, m_cv, T, rho, dVdTheta, dmCVdTheta, dxoildTheta, m_dot_dict, h_dict, x_l, Q_dot = 0):58  59    # Initialization of the rate of change of temperature60    dTdTheta = 061  62    # The energy balance equation can be decomposed into 5 terms accounting for 5 different mechanisms through63    # which energy is transferred from and to the control volume64  65    # Term 1 - Work exerted by the change of volume66    dPdT = propF.propfunction_mixture(property_matrix, x_l, 'dPdT', T = T, rho = rho)67    term1 = -T*dPdT*(dVdTheta - 1/rho*dmCVdTheta)68    69    # Term 2 - Energy contained in the oil mass fraction70    u_l = propF.propfunction_mixture(property_matrix, 1, 'u', T = T, rho = rho)71    u_g = propF.propfunction_mixture(property_matrix, 0, 'u', T = T, rho = rho)72    term2 = -m_cv*(u_l - u_g)*dxoildTheta73    74    # Term 3 - Energy trasported by the fluid coming into the control volume and going out of it75    h =  propF.propfunction_mixture(property_matrix, x_l, 'h', T = T, rho = rho)   76    term3 = -h*dmCVdTheta77    78    # Term 4 - Heat exchange in the control volume79    term4 = Q_dot/omega80    81    # Term 5 - Rate of change of enthalpy82    m_dot_h_sum = 0 83    for i, j in m_dot_dict.items():84        m_dot_h_sum += m_dot_dict[i]* h_dict[i]     85    term5 = 1/omega* m_dot_h_sum86    87    # Rate of change of internal energy88    dudT = propF.propfunction_mixture(property_matrix, x_l, 'dudT', T = T, rho = rho)89    90    # All the terms are taken into account in a single equation. Temperature derivative is the only unknwon term91    # in the resulting equation92    dTdTheta = (term1 + term2 + term3 + term4 + term5)/dudT/m_cv93    94    return  dTdTheta9596# ***???***97def supply_HT(HT_bool, T_flow_in, m_dot = 0, T_wall = 0):98    99    if HT_bool != True:100        T_flow_out = T_flow_in101        102    return T_flow_out103104# The 'discharge_is_guess' function is meant to compute a first-attempt value of temperature and density in the discharge105# chambers as they are initialized. Temperature and density at machine inlet and pressure at outlet are known.106# The isentropic efficiency used is a first-attempt value that is specified as input (isentropic efficiency for this107# kind of machine is the rage 0.6-0.7)108def discharge_is_guess(property_matrix, rho_in, T_in, P_out, x_l, eta_is_guess):109    110    # Enthalpy and entropy at inlet are computed based on temperature and density at inlet111    h_in = propF.propfunction_mixture(property_matrix, x_l, 'h', T = T_in, rho = rho_in)112    s_in = propF.propfunction_mixture(property_matrix, x_l, 's', T = T_in, rho = rho_in)113    114    # Enthalpy at outlet is computed based on the isentropic assumption115    h_is_guess = propF.propfunction_mixture(property_matrix, x_l, 'h', s = s_in, P = P_out)   116    117    # Enthalpy at outlet is computed using the first-attempt value of the isentropic efficiency provided as input118    h_out_guess = (h_is_guess - h_in)/eta_is_guess + h_in119    120    # Temperature and density at outlet are computed based on the outlet enthalpy guessed121    T_out_guess = propF.propfunction_mixture(property_matrix, x_l, 'T', h = h_out_guess, P = P_out)122    rho_out_guess = propF.propfunction_mixture(property_matrix, x_l, 'rho', h = h_out_guess, P = P_out)123    124    # Solver to implement when adding the oil125    return T_out_guess, rho_out_guess126127128# The 'discharge_is_guess' function is meant to compute a first-attempt value of temperature and density in the compression129# chambers as they are initialized. Temperature and density at machine inlet and pressure at outlet are known.130# The isentropic efficiency used is a first-attempt value that is specified as input (isentropic efficiency for this131# kind of machine is the rage 0.6-0.7) 132def compression_is_guess(property_matrix, dict_V_c, alpha, rho_in, T_in,  x_l, eta_is_guess):133    134    # Entropy at inlet are computed based on temperature and density at inlet135    s_in = propF.propfunction_mixture(property_matrix, x_l, 's', T = T_in, rho = rho_in)136    137    # Volume of the compression chamber as the process begins (alpha = 0 refers to the outermost compression chamber)138    V_init = np.polyval(dict_V_c[0], 0) 139    # Volume of the compression chamber at the time the current function is called (alpha > 1)140    V_final = np.polyval(dict_V_c[alpha - 1], 0)141    142    # To compute the first-attempt density of the compression chamber, the assumption made is that mass inside the 143    # compression chamber is constant throughtout the whole process144    rho_comp_guess = rho_in*(V_init/V_final)145    # First-attempt temperature is computes based on entropy and density.146    T_comp_guess = propF.propfunction_mixture(property_matrix, x_l, 'T', s = s_in, rho = rho_comp_guess)147    148    # Solver to implement when adding the oil149    return T_comp_guess, rho_comp_guess150151# The 'isentropic_nozzle' function computes mass flow rate, enthalpy and Reynolds number of a fluid flowing through 152# an orifice that could be asumed to behave as an isentropic nozzle. Flows though the suction port and the discarge153# port are modeled as isentropic nozzles in the current machine model.154# It is not specified at this stage what is the upstream state and downstream state.155# This function is called within a class method. Properties of the chamber that is the object of the method are156# specified using number 1, properties of the external chamber using number 2.157# i.e.: if we are computing the mass flow rate of the leakage that takes place from compression chamber 'c1' to158# suction chamber 's1', T1 represents temperature of chamber 'c1', T2 represents temperature of chamber 's1'.159def isentropic_nozzle(geo, property_matrix, area, rho1, T1, x_l1, rho2 = None, T2 = None,  x_l2 = None ,h2 = None, P2 = None, leakage_type = None):160    161    # 'flow_function' is defined within the 'isentropic_nozzle' environment and as such they share the same variables.162    # This function is used to compute the mass flow rate of two-phase flow flowing through an isentropic nozzle.163    def flow_function(area, s_up, P_up, P_down , v_up, v_down):   164        #     Parameters for two phase flow165        # area correction factor166        Xd = 1;167        # Morris correction coefficient168        Cd = 0.77;169        # homogenous flow coefficient170        psi = 1;171        # area ratio between throat and suction172        sigma = 0;173        # Number of points for integral calculation174        N_p = 10;175        176        # Discretization of the pressure variation through the nozzle177        D_p = (P_up - P_down)/(N_p-1)178        P = np.linspace(P_up, P_down, N_p, endpoint =True)   179        180        # Initialization of the specific volume181        v = np.ones(len(P))*v_up182 183        # Accordingly with the eqution presented in Ian Bell's thesis, to compute the mass flow rate it184        # is necessary to perform the integral of the specific volume over the whle range of pressure.185        integral = 0186        # The initial value of the specific volume is computed in correspondence of the upstream properties187        v[0] = 1/propF.propfunction_mixture(property_matrix, x_l_up, 'rho', s = s_up , P = P[0])188        for i in range(len(P)-1):189            # The i-th value of the specifiv volume is computed based on the dicretized pressure and assuming190            # an isentropic process191            v[i+1] = 1/propF.propfunction_mixture(property_matrix, x_l_up, 'rho', s = s_up , P = P[i+1]);192            # The value of the integral is updated after each dp considered193            integral = integral + (v[i+1] + v[i])/2*D_p194        195        # The final equation presented in Ian Bell's thesis is solved196        m_dot = Cd*Xd*area*(2*integral/(v_down**2 - sigma*v_up**2))**0.5197    198        return m_dot199   200    # area = 6e-11201    # rho1 = 23202    # rho2 = 22203    # T1 = 302204    # T2 = 303205    # x_l1 = 0206    # x_l2 = 0207    208    # If the flow ares is closed the mass flow rate is zero209    if area == 0: 210        m_dot = 0   211        h = 0212    # If the area is non-zero the mass flow rate computation begins.213    else :214        # Pressure 1 is computed based on temperature and density.215        # Pressure 2 could be computed using either enthalpy and pressure or temperature and density.216        if P2 == None and h2 == None :217            P1 = propF.propfunction_mixture(property_matrix, x_l1, 'P', T = T1 , rho = rho1)218            P2 = propF.propfunction_mixture(property_matrix, x_l2, 'P', T = T2 , rho = rho2)219        elif rho2 == None and T2 == None: 220            P1 = propF.propfunction_mixture(property_matrix, x_l1, 'P', T = T1 , rho = rho1)221            rho2 = propF.propfunction_mixture(property_matrix, x_l2, 'rho', h = h2 , P = P2)222            T2 = propF.propfunction_mixture(property_matrix, x_l2, 'T', h = h2 , P = P2)223        # If input are not coherent the function raises an error.224        else:225            raise Exception('Invalid inputs')226            227        # This is the minimum pressure ratio that is required for existence of flow228        DP_floor = 1; #Pa229    230        # If delta_P is higher than the minimum required the computation takes place.231        if abs(P1 - P2) > DP_floor:232            # Upstream and downstream properties are defined based on the value of pressure in the two chambers233            if P1 > P2:234                x_l_up = x_l1235                P_up = P1236                T_up = T1237                rho_up = rho1 238                h_up = propF.propfunction_mixture(property_matrix, x_l1, 'h', T = T_up , rho = rho_up) 239                visc_up = propF.propfunction_mixture(property_matrix, x_l1, 'mvisc', h = h_up, P = P_up)240                v_up = 1/rho_up241                s_up = propF.propfunction_mixture(property_matrix, x_l1, 's', T = T_up , rho = rho_up)242                x_down = x_l2243                P_down = P2244                T_down = T2245                rho_down = propF.propfunction_mixture(property_matrix, x_l1, 'rho', s = s_up , P = P_down)246                v_down = 1/rho_down247            else:248                x_l_up = x_l2249                P_up = P2250                T_up = T2251                rho_up = rho2252                h_up = propF.propfunction_mixture(property_matrix, x_l1, 'h', T = T_up , rho = rho_up)253                visc_up = propF.propfunction_mixture(property_matrix, x_l1, 'mvisc', h = h_up, P = P_up)254                v_up = 1/rho_up255                s_up = propF.propfunction_mixture(property_matrix, x_l1, 's', T = T_up , rho = rho_up)256                x_down = x_l1257                P_down = P1258                T_down = T1259                rho_down = rho_down = propF.propfunction_mixture(property_matrix, x_l1, 'rho', s = s_up , P = P_down)260                v_down = 1/rho_down261           262            # Computation of the critical pressure ratio. If pressure ratio exceeds this value sonic condition263            # in the nozzle are attained264            k_star = propF.propfunction_mixture(property_matrix, x_l2, 'k_star', T = T_up , rho = rho_up)265            ratio_crit = (1 + (k_star - 1)/2)**(k_star/(1 - k_star))266            # Pressure ratio between downstream chamber and upstream chamber. This ratio is always < 1. 267            ratio = P_down/P_up268            269            # If pressure ratio is lower than the critical pressure ratio (upstream pressure >> downstream pressure)270            # the flow is choked and flow properties do not depend from downstream state anymore.271            if ratio <= ratio_crit:272               P_down = P_up*ratio_crit273               v_down = 1/propF.propfunction_mixture(property_matrix, x_l2, 'rho', s = s_up , P = P_down)   274            # Once all the properties are defined, the mass flow rate is computed by means of the 'flow_function'275            m_dot = flow_function(area, s_up, P_up, P_down , v_up, v_down)       276            # The mass flux is considered to be positive if it enters the chamber on which the method is applied277            if P1 > P2:278                m_dot = -m_dot279                h = h_up            280            if P1 < P2:281                m_dot = m_dot282                h = h_up283        284        # If delta_P is lower than the minimum required the flow through the nozzle is zero.285        else:286            m_dot = 0   287            h = 0288        289    # Computation of Reynolds number of the flow.290    # Re = 0 if mass flow rate is zero291    if leakage_type == None or area == 0 or abs(P1 - P2) <= DP_floor:292        Re = 0 293    # Radial leakage Reynolds number294    elif leakage_type == 'r':295        Re = abs(m_dot)/area*2*geo.delta_r/visc_up296    # Flank leakage Reynolds number297    elif leakage_type == 'fl':298        Re = abs(m_dot)/area*2*geo.delta_f/visc_up299    # Raise error if input are not coherent300    else: 301        raise Exception('Invalid inputs')302        303        304    return (m_dot, h, Re)  305306    307308309    310    311    312    313    314    315    316    317    318    319    320    321    
...

load_db.py

Source:load_db.py

1from phrase_structures import *2from utils import *3import json4import sqlite35def build_spider_dataset(num):6	valid = True  # if there exists empty columns in the database, then the database is assumed invalid, reported7	# and discarded8	with open(TABLE_METADATA_PATH, 'r') as fp:9		json_file = json.load(fp)10	meta = json_file[num]11	if SETTING in ['spider', 'chisp']:12		database_name = DATABASE_PATH + '/%s/%s.sqlite' % (meta['db_id'], meta['db_id'])13	elif SETTING == 'dusql':14		database_name = DATABASE_PATH + '/%s.db' % meta['db_id']15	else:16		raise AssertionError17	print(meta['db_id'])18	propertynps = []19	conn = sqlite3.connect(database_name)20	# conn.create_function('FIXENCODING', 1, lambda s: str(s).encode('latin_1'))21	crsr = conn.cursor()22	typenps = []23	print("")24	print("All table names: ")25	for idx, tab_name in enumerate(meta['table_names']):26		c_english = '@' + tab_name + '@'27		c_chinese = '@' + tab_name + '@'28		z = meta['table_names_original'][idx]29		print(z)30		if ' ' in z or z[0] in ['0', '1', '2', '3', '4', '5', '6', '7', '8', '9'] or '(' in z:31			z = '[' + z + ']'32		properties = []33		for prop_idx, prop_ori in enumerate(meta['column_names']):34			if prop_ori[0] == idx:35				# '*' disgarded at generation time36				properties.append(prop_idx - 1)37		typenps.append(TYPENP(c_english=c_english, c_chinese=c_chinese, z=z, overall_idx=idx, properties=properties))38	# disgard the '*'39	for idx, prop_name in enumerate(meta['column_names']):40		table_idx = prop_name[0]41		if table_idx < 0:42			assert idx == 043			# print('star: ', idx, prop_name)44			continue45		c_english = '{0}$' + prop_name[1] + '$'46		c_chinese = '{0}$' + prop_name[1] + '$'47		col_z = meta['column_names_original'][idx][1]48		# if the name contains ' ', wrap it with brackets49		if ' ' in col_z or col_z[0] in ['0', '1', '2', '3', '4', '5', '6', '7', '8', '9'] or '(' in col_z:50			col_z = '[' + col_z + ']'51		z = '{0}%s' % col_z52		dtype = meta['column_types'][idx]53		if idx in meta['primary_keys']:54			is_primary_key = True55		else:56			is_primary_key = False57		cur_prop = PROPERTYNP(c_english=c_english, c_chinese=c_chinese, z=z, dtype=dtype, table_id=table_idx,58					   overall_idx=idx - 1, meta_idx=idx - 1, is_pk=is_primary_key)59		if idx in meta['primary_keys']:60			cur_prop.is_primary = True61		'''62		query_all = 'select %s from %s' % (cur_prop.z.format(''), typenps[cur_prop.table_id].z)63		try:64			all_values = crsr.execute(query_all).fetchall()65		except Exception as e:66			print(e)67			print(query_all)68			raise69		temp_all_values = []70		for item in all_values:71			try:72				item_str = str(item)73				temp_all_values.append(item_str)74			except Exception as e:75				continue76		if is_unique(temp_all_values):77			cur_prop.is_unique = True78		else:79			cur_prop.is_unique = False80		'''81		propertynps.append(cur_prop)82	fks = meta['foreign_keys']  # a.k.a the prop_rels83	#typenp_unique_fk_num = [[] for _ in range(len(typenps))]84	for pair in fks:85		assert pair[0]-1 >= 0 and pair[1]-1 >= 086		propertynps[pair[0]-1].is_fk_left = True87		propertynps[pair[1]-1].is_fk_right = True88		'''89		prop1 = propertynps[pair[0]-1]90		prop2 = propertynps[pair[1]-1]91		tid1 = prop1.table_id92		tid2 = prop2.table_id93		# if the other side of this foreign key relation is a primary key, it means the column in this side of this94		# foreign key relation is a reference, namely an edge column, two edge columns make up an edge table95		if prop1.is_primary and (pair[1]-1) not in typenp_unique_fk_num[tid2]:96			typenp_unique_fk_num[tid2].append(pair[1]-1)97		if prop2.is_primary and (pair[0]-1) not in typenp_unique_fk_num[tid1]:98			typenp_unique_fk_num[tid1].append(pair[0]-1)99		'''100	#for tid, item in enumerate(typenp_unique_fk_num):101	#	if len(item) == 2:102	#		typenps[tid].is_edge = True103	#		print("edge table: %s" % typenps[tid].c_english)104	#print("")105	for pair in fks:106		if meta['column_names'][pair[0]][0] == meta['column_names'][pair[1]][0]:107			new_typenp = copy.deepcopy(typenps[meta['column_names'][pair[0]][0]])108			new_typenp.overall_idx = len(typenps)109			new_typenp.properties = []110			for pid in typenps[meta['column_names'][pair[0]][0]].properties:111				new_prop = copy.deepcopy(propertynps[pid])112				new_prop.table_id = new_typenp.overall_idx113				new_prop.overall_idx = len(propertynps)114				new_typenp.properties.append(new_prop.overall_idx)115				propertynps.append(new_prop)116				propertynps[pid].clones.append(new_prop.overall_idx)117			typenps.append(new_typenp)118			typenps[meta['column_names'][pair[0]][0]].clones.append(new_typenp.overall_idx)119	# here they're both directed120	fk_type_matrix = [[0] * len(typenps) for i in range(len(typenps))]121	fk_property_matrix = [[0] * len(propertynps) for i in range(len(propertynps))]122	for pair in fks:123		if meta['column_names'][pair[0]][0] == meta['column_names'][pair[1]][0]:124			tid = meta['column_names'][pair[0]][0]125			for i, clone_id in enumerate(typenps[tid].clones):126				fk_type_matrix[clone_id][tid] += 1127				fk_property_matrix[pair[0] - 1][propertynps[pair[1] - 1].clones[i]] += 1128				fk_property_matrix[propertynps[pair[0] - 1].clones[i]][pair[1] - 1] += 1129		else:130			fk_property_matrix[pair[0] - 1][pair[1] - 1] += 1131			fk_type_matrix[meta['column_names'][pair[1]][0]][meta['column_names'][pair[0]][0]] += 1132	# PROBABLY BRING THE NAME MATCHING INTO ACCOUNT?133	nm_type_matrix = [[0] * len(typenps) for i in range(len(typenps))]134	nm_property_matrix = [[0] * len(propertynps) for i in range(len(propertynps))]135	for i in range(len(propertynps)):136		for j in range(len(propertynps)):137			if propertynps[i].table_id in typenps[propertynps[j].table_id].clones or propertynps[j].table_id in typenps[138				propertynps[i].table_id].clones:139				continue140			if propertynps[i].z.format('') == propertynps[j].z.format(''):141				if propertynps[i].z.format('').lower() == 'id' and propertynps[i].overall_idx != propertynps[142					j].overall_idx \143						and propertynps[i].overall_idx not in propertynps[j].clones \144						and propertynps[j].overall_idx not in propertynps[i].clones:145					continue146				nm_property_matrix[i][j] += 1 * NAME_PAIR_WEIGHT147				nm_property_matrix[j][i] += 1 * NAME_PAIR_WEIGHT148				nm_type_matrix[propertynps[j].table_id][propertynps[i].table_id] += 1 * NAME_PAIR_WEIGHT149				nm_type_matrix[propertynps[i].table_id][propertynps[j].table_id] += 1 * NAME_PAIR_WEIGHT150	for prop in propertynps:151		# if prop.dtype == 'varchar(1000)':152		#	conn.execute('UPDATE "{0}" SET "{1}" = FIXENCODING("{1}")'.format(typenps[prop.table_id].z, prop.z.format('')))153		query_all = 'select %s from %s' % (prop.z.format(''), typenps[prop.table_id].z)154		try:155			all_values = crsr.execute(query_all).fetchall()156		except Exception as e:157			print(e)158			print(query_all)159			raise160		# if a property is empty before discarding invalid values, it means whole tables have no entries, assume as161		# invalid database162		if len(all_values) == 0:163			print("Empty property: ", prop.z.format(typenps[prop.table_id].z))164			valid = False165		if len(all_values) > MAX_VALUE_ENTRIES:166			all_values = random.choices(all_values, k=MAX_VALUE_ENTRIES)167		temp_all_values = []168		for idx, item in enumerate(all_values):169			if isinstance(item, tuple) or isinstance(item, list):170				assert (len(item) == 1)171				item = item[0]172			# if isinstance(item, bytes):173			#	item = item.decode('latin_1')174			try:175				item_str = str(item)176			except Exception as e:177				print("!")178				item_str = item179			if prop.dtype in ['datetime', 'bool'] and item_str[0] not in ["'", '"']:180				item_str = "'" + item_str + "'"181			elif ('varchar' in prop.dtype or prop.dtype in ['bit', 'id']) and len(item_str) > 0 and "'" != item_str[0] and '"' != item_str[0] and "'" != item_str[-1] and '"' != item_str[-1]:182				item_str = '"' + item_str + '"'183			elif SETTING == 'dusql':184				item_str = '"' + item_str.strip('\"\'') + '"'185			# do not save those None or '' values to avoid confusion186			if item is None or item == '':187				continue188			temp_all_values.append(189				VALUENP(c_english=item_str, c_chinese=item_str, z=item_str, dtype=prop.dtype))190		# if prop.dtype == 'bool':191		#	print(all_values)192		all_values = temp_all_values193		if is_unique(all_values):194			prop.is_unique = True195		prop.set_values(all_values)196	for idx in range(len(typenps)):197		type_valid = False198		for pid in typenps[idx].properties:199			if propertynps[pid].valid is True:200				type_valid = True201		if not type_valid:202			typenps[idx].valid = False203	type_matrix = [[0] * len(typenps) for i in range(len(typenps))]204	for i in range(len(type_matrix)):205		for j in range(len(type_matrix[0])):206			type_matrix[i][j] = fk_type_matrix[i][j] + nm_type_matrix[i][j]207	property_matrix = [[0] * len(propertynps) for i in range(len(propertynps))]208	for i in range(len(property_matrix)):209		for j in range(len(property_matrix)):210			property_matrix[i][j] = fk_property_matrix[i][j] + nm_property_matrix[i][j]211	property_matrix_np = numpy.matrix(property_matrix)212	# don't tune property edge weight with pagerank, do it through type pagerank.213	'''214	property_scores = pagerank_scores(property_matrix_np)215	# print(property_matrix)216	# print("")217	# print(property_scores)218	for i in range(len(property_matrix)):219		coefficient = property_scores[i] / float(sum(property_matrix[i]))220		for j in range(len(property_matrix)):221			property_matrix[i][j] = coefficient * property_matrix[i][j]222	'''223	# turn the property_matrix un-directioned224	#for i in range(len(property_matrix)):225	#	for j in range(len(property_matrix)):226	#		property_matrix[i][j] += property_matrix[j][i]227	fk_prop_rels = []228	for pair in fks:229		if propertynps[pair[0] - 1].table_id == propertynps[pair[1] - 1].table_id:230			tid = propertynps[pair[0] - 1].table_id231			for i, clone_id in enumerate(typenps[tid].clones):232				fk_prop_rels.append(233					PROP_REL(pair[0] - 1, propertynps[pair[1] - 1].clones[i], tid,234							 clone_id,235							 property_matrix[pair[0] - 1][propertynps[pair[1] - 1].clones[i]]))236				fk_prop_rels.append(237					PROP_REL(propertynps[pair[0] - 1].clones[i], pair[1] - 1, clone_id,238							 tid,239							 property_matrix[propertynps[pair[0] - 1].clones[i]][pair[1] - 1]))240		fk_prop_rels.append(241			PROP_REL(pair[0] - 1, pair[1] - 1, propertynps[pair[0] - 1].table_id, propertynps[pair[1] - 1].table_id,242					 property_matrix[pair[0] - 1][pair[1] - 1]))243	nm_prop_rels = []244	for prop1 in propertynps:245		for prop2 in propertynps:246			if prop1.table_id in typenps[prop2.table_id].clones or prop2.table_id in typenps[247				prop1.table_id].clones or prop1.overall_idx == prop2.overall_idx:248				continue249			if prop1.z == prop2.z:250				# exclude those 'id' queries which are not in the same query251				if prop1.z.format('').lower() == 'id' and prop1.overall_idx != prop2.overall_idx \252						and prop1.overall_idx not in prop2.clones and prop2.overall_idx not in prop1.clones:253					continue254				nm_prop_rels.append(PROP_REL(prop1.overall_idx, prop2.overall_idx, prop1.table_id, prop2.table_id,255											 property_matrix[prop1.overall_idx][prop2.overall_idx] * NAME_PAIR_WEIGHT))256	prop_rels = fk_prop_rels + nm_prop_rels257	type_matrix = numpy.matrix(type_matrix)258	property_matrix = numpy.matrix(property_matrix)259	for p in propertynps:260		if p.values is None:261			print(p.z)262			raise AssertionError263		if len(p.values) == 0:264			p.valid = False265	print("Finished resolving database $ %s $" % meta['db_id'])266	num_queries = calc_num_queries_via_stats(len(fk_prop_rels))267	return database_name, meta[...

Automation Testing Tutorials

Learn to execute automation testing from scratch with LambdaTest Learning Hub. Right from setting up the prerequisites to run your first automation test, to following best practices and diving deeper into advanced test scenarios. LambdaTest Learning Hubs compile a list of step-by-step guides to help you be proficient with different test automation frameworks i.e. Selenium, Cypress, TestNG etc.