H-S

[JungHulk]

# -*- coding: utf-8 -*-
"""
Created on Wed Jan 25 16:31:38 2023

@author: seungj.jung
"""
## rev3. air temperature = ambient

import os
import pandas as pd
import pickle
import math
import numpy as np

import Unitconv as u
import Fan_curve as f

df_air = pd.read_pickle('Air_prop.p')
unitconv = u.unitconv
fan = f.Fan
fan_press = f.Plot.Curve()[1]

class Chip:
    def __init__(self, power, N, L, W):
        self.power = power
        self.N = N
        self.L = L
        self.W = W
        self.tot_power = self.power * self.N
        # return tot_power
        
class Heatsink: # option 단추로 Heatsink dimension 입력값 선택
    def __init__ (self, H, L, b, t, N, S, k_fin): # H, L, b, t, N, S, k_fin
        if H <= b:
            print("Heatsink dimension (b) should be checked")
            print("H : ", H, " < or =", "b : ",b)
        elif 0 in [H, L, b, t, N, S, k_fin]:
            print("Heatsink dimension should not be zero")
        else:
            self.H = H
            self.L = L
            self.b = b
            self.L_fin = self.H - self.b
            self.t = t
            self.N = N
            self.S = S
            self.W = (self.N-1) * self.S + self.N*self.t
            self.k_fin = k_fin
        
    # def __init__ (self, H, L, b, t, N, W, k_fin): # H, L, b, t, N, W, k_fin
    #     if H <= b:
    #         print("Heatsink dimension (b) should be checked")
    #         print("H : ", H, " < or =", "b : ",b)
    #     elif 0 in [H, L, b, t, N, W, k_fin]:
    #         print("Heatsink dimension should not be zero")
    #     else:
    #         self.H = H
    #         self.L = L
    #         self.b = b
    #         self.L_fin = self.H - self.b
    #         self.t = t
    #         self.N = N
    #         self.W = W
    #         self.S = (self.W - (self.N*self.t)) / (self.N-1) 
    #         self.k_fin = k_fin    

# Emissivity 와 Thermal conductivity of heat sink 테이블로 만들기

class Calculation():
    def setdata(self, T_amb, T_dif, Chip, Heatsink): # first = T_ambient., second = Temperature diff, Chip 정보, Heatsink 정보
        self.T_amb = T_amb
        self.T_dif = T_dif
        
        self.N_s = Chip.N
        self.L_s = Chip.L
        self.W_s = Chip.W
        self.tot_power = Chip.tot_power
        self.H = Heatsink.H 
        self.L = Heatsink.L 
        self.b = Heatsink.b
        self.L_fin = Heatsink.L_fin
        self.t = Heatsink.t 
        self.N = Heatsink.N
        self.S = Heatsink.S
        self.W = Heatsink.W
        self.k_fin = Heatsink.k_fin

    def Air_prop(self, Table, Prop, T): # df_air, input temperature, Properties name(" ")
        
        if T > 500 or T <-80:     #Temperature range 설정 (-40 ~ 80)
            print("Out of range")
            value = []
            prop_list = []
        else:
            if T == math.floor(T*10)/10:
                value = Table.loc[T][Prop] # 해당하는 하나의 물성치만 가져옴
                prop_list = Table.loc[T]    # 온도에 해당되는 모든 물성치 가져옴
            else:
                T_up = math.ceil(T*10)/10
                T_down = math.floor(T*10)/10
                val_list = []
                
                list_prop = Table.columns
                for i in list_prop:
                    val_up   = Table.loc[T_up][i]
                    val_down = Table.loc[T_down][i]
                    if type(val_up) == float or type(val_down) == float:
                        if T_up == T_down:
                            val = T_up
                        else:
                            val = val_down + (T-T_down) * (val_up - val_down)/(T_up-T_down)
                    else:
                        val = Table.loc[T_up][i]
                    
                    val_list.append(val)
                    
                df_val_list = pd.Series(val_list, index = list_prop)
                value = df_val_list[Prop] # Prop 이 몇번째 index 인지 확인 후 prop_list 에서 호출
                prop_list = df_val_list# 온도에 해당되는 모든 물성치 가져옴
        return value, prop_list
        # aa = Air_prop(df_air, 32.6, "Pressure")
        # aa[0]   aa[1]["Cp"] 

    def cal_temp(self, V_flow): # 열전달 후 온도 계산
        # self.V_flow = V_flow          # Fan 유량 받기 # ft3/min
        # T_dif = range(1,500)            # differential temperature between surface and air
        T_amb = self.T_amb  # 외기 값
        # T_dif = self.T_dif # 초기 값
        T_max = self.T_amb + self.T_dif     # Ideal max temperature of heat sink surface
        T_ave = (T_amb+T_max)/2   # air 물성치는 평균값을 사용
        
        
        air = self.Air_prop(df_air, "Pressure", T_ave) # T_max 해당되는 air properties 호출
        
        T    = T_ave + 273.15
        D    = air[1][1]
        Cp   = air[1][4]          # Specific heat, Kj/kg-k
        k    = air[1][6] * 10**-3 # thermal conductivity, W/m-k
        Vis_dyn  = air[1][7] * 10**-6 # dynamic viscosity, Pa-s (=kg/m-s)
        Vis_kin  = air[1][8] * 10**-4 # kinetic viscosity, m^2/s
        Pr   = air[1][11]
        alpha    = air[1][12] * 10**-4 # thermal diffusivity, m^2/s 
        beta     = 1/T   # expension coefficient, 1/K

        sigma = 5.67*10**-8 # Stefan-Boltzmann constant, W/m^2k^4
        epsil = 0 # Emissivity of heat sink 

        # Duct size 가정하기.(Fan size 로 입력 받기?? 현재는 Heat sink size) Heatsink 외부와의 열전달 계수 구하기 위해 필요한 속도계산, 해당 면적 어떻게 받을지 고민하기.
        A_duct = unitconv("length", "mm", "m", self.W).cal() * unitconv("length", "mm", "m", self.H).cal() # m2
        V_flow_meter = unitconv("volume", "ft3", "m3", V_flow).cal() # m3/min
        Vel_inf = V_flow_meter / A_duct /60 # 외기 속도 m/s
        # print(Vel_inf)
        Re = Vel_inf * self.L * 10**-3 / (Vis_kin)
        Nu = 0.664 * Pr**(1/3) * Re**(1/2)
        h1 = Nu * k / (self.L*10**-3)

        A1 = self.H * 10**-3 * self.L * 10**-3 + self.t * 10**-3 * (2 * self.H * 10**-3 + self.L * 10**-3 )
        Q_c1_max = 2 * h1 * A1 * (T_max - T_amb) # Heat loss of external surface
        Q_c1_base_max = 2 * h1 * self.b *10**-3 * (self.L + 2*self.t)*10**-3 *(T_max - T_amb)
        Q_c1_fin_max = Q_c1_max - Q_c1_base_max

        # fin 효율 위한 계산
        L_c = (self.L_fin + self.t /2)*10**-3 # 수정된 fin 길이 (fin 끝단에서도 대류 고려)
        P = 2 * (self.L + self.t)*10**-3 # fin 끝단 둘레 길이
        A_c = self.L*10**-3 * self.t*10**-3 # fin 단면적 ,m2
        m_c1 = math.sqrt(h1 * P / (self.k_fin * A_c))
        eff_c1_fin = math.tanh(m_c1*L_c) / (m_c1*L_c)

        # fin 에서의 convective heat transfer coefficient 계산 위한 과정
        #### 층류 난류 판별식 필요. 이후 계산, fin 사이 Rey 수가 층류-> Nu Eq'n laminar flow
        
        A_s = (self.N - 1)* self.S * (self.L_fin) # mm2, 핀사이 유로의 단면적
        A_s_feet = unitconv("area", "mm2", "ft2", A_s).cal() # ft2
        Vel_mean_feet = V_flow / A_s_feet # ft/min
        Vel_mean = unitconv("length", "ft", "m", Vel_mean_feet).cal() /60 # m/s
        D_hyd = 2 * (self.L_fin*10**-3 * self.S*10**-3) / (self.L_fin*10**-3 + self.S*10**-3)  # Hydraulic diameter
        Re_D  = Vel_mean * D_hyd / Vis_kin # Hydraulic diameter Reynolds Number
        Gz    = Re_D * Pr / (self.L * 10**-3 / D_hyd) # Graetz number : Characteristic Parameter describing flow in a Duct
        
        if Re_D < 2300: ## 아직 검증 안함
            # Nu_m  = 0.664 * Gz**(1/2) / (Pr**(1/2) * Pr**(-1/3)) * (1+7.3 * (Pr/Gz)**(1/2))**(1/2) # Mean Nusselt number for Channel Flow in excel
            Nu_1 = 7.541
            Nu_2 = 1.841 * (Re_D * Pr * D_hyd/self.L)**1/3
            Nu_3 = 0
            Nu_m = (Nu_1**3 + Nu_2**3 + Nu_3**3)**(1/3)
            
        else:
            Nu_m = 0.023 * Re_D**(0.8) * Pr**(0.4) # 잘 맞음 Dittus boelter equation
  
        
        ## Analytical forced convection modeling of plate fin heat sinks
        # Re_D  = Vel_mean * self.S / Vis_kin # Hydraulic diameter Reynolds Number
        # Re_s = Re_D * self.S/self.L
        # Nu_i = ((Re_s*Pr)**(-3)/2 + (0.664*math.sqrt(Re_s)*Pr**(1/3) * math.sqrt(1+(3.65/math.sqrt(Re_s))))**(-3) )**(-1/3)
        # eqt = 2*Nu_i *self.k_fin/k *self.H/self.S *self.H/self.t *(self.t/self.L + 1) 
        # Nu_m = math.tanh(math.sqrt(eqt)) / math.sqrt(eqt) * Nu_i
        
        
        ## Jouhara and Axcell(2009), Hwang and Fan(1963)
        # X = Re_D * Pr
        # Nu_m = 7.55 + (0.024 * X**1.14)/(1+0.0358* Pr**0.17 *X*0.64)
        
        # Nu_m = 0.332 * Re_D**(1/2) * Pr**(1/3) # 
        # Nu_m = 0.9760 * Re_D**(0.3337) * Pr**(1/3) # finned heat exchanger
        # Nu_m = 0.4217 * Re_D**(0.4700) * Pr**(1/3) # finned heat exchanger
        
        # Re_s = Re_D * self.S / self.L # 논문 Nu correlation for tise..
        # Nu_m = ( (Re_s * Pr / 2)**(-3) + 0.664 * math.sqrt(Re_s) * Pr**(1/3) * math.sqrt(1+(3.65/math.sqrt(Re_s)))**(-3))**(-1/3)
        
        
        h2 = Nu_m * k / D_hyd # Heat transfer coefficient of external surface
        A2 = self.L * 10**-3 * (2*(self.L_fin * 10**-3) + self.S*10**-3) + 2*(self.t*10**-3 * self.H*10**-3 + self.S*10**-3 * self.b*10**-3) + self.t*10**-3 * self.L*10**-3 # Total area of external surface
        Q_c2_max = h2 * A2 * (T_max - T_amb)  # Heat loss of external surface
        Q_c2_base_max = h2 * (2 * self.b*10**-3 * self.S*10**-3 + self.L*10**-3*self.S*10**-3 + 2*self.b*10**-3 *self.t*10**-3) * (T_max - T_amb)
        Q_c2_fin_max  = Q_c2_max - Q_c2_base_max

        m_c2 = math.sqrt(h2 * P / (self.k_fin * A_c))
        eff_c2_fin = math.tanh(m_c2*L_c) / (m_c2*L_c)

        Q_r1 = 2 * epsil * sigma * A1 * ((T_max+273.15)**4 - (T_amb+273.15)**4) # Heat loss of side surface
        A_r2 = self.L *10**-3  *(self.t + self.S)*10**-3 + 2* (self.t*10**-3 * self.H*10**-3 + self.S*10**-3 * self.b*10**-3) # Total area of remain surface
        Q_r2 = epsil * sigma * A_r2 * ((T_max+273.15)**4 - (T_amb+273.15)**4)# Heat loss of external surface

        Q_tot = Q_c1_base_max + Q_c1_fin_max * eff_c1_fin + Q_r1 + (self.N-1) * (Q_c2_base_max + Q_c2_fin_max * eff_c2_fin + Q_r2)  # W

        values = np.array([T_max, Q_tot, h1, Re, Nu, A1, Q_c1_max, Q_c1_base_max, Q_c1_fin_max, eff_c1_fin, m_c1, L_c, T, D, Cp, Vis_dyn, Vis_kin, k, alpha, beta, Pr, D_hyd, Re_D, Gz, Nu_m, h2, A2, Q_c2_max, Q_c2_base_max, Q_c2_fin_max, eff_c2_fin, sigma, Q_r1, A_r2, Q_r2])
    
        return Q_tot, T_max, values


    def cal_rsp(self, V_flow): # resistamce spreading
        T_amb = self.T_amb + 273.15
        T_max = T_amb + self.T_dif
        # self.V_flow = V_flow
        # L_s1 = Input.L_s1
        # W_s1 = Input.W_s1
        
        # L = Heat_sink.L
        # W = Heat_sink.W
        # tp = Heat_sink.b
        # t = Heat_sink.t
        # S = Heat_sink.S
        # k_fin = Heat_sink.k_fin
        
        r1 = math.sqrt(self.L_s *10**-3 * self.W_s*10**-3 / math.pi)
        r2 = math.sqrt(self.L *10**-3 * self.W*10**-3 / math.pi)

        r_ratio =  r1/r2
        tau  = self.b*10**-3 / r2
        
        h_ave = self.cal_temp(V_flow)[2][27] / (self.T_dif * self.L *10**-3 * (self.t + self.S)*10**-3)
        epsil = 0 # Emissivity of heat sink 
        sigma = 5.67*10**-8 # Stefan-Boltzmann constant, W/m^2k^4
        h_rad = epsil * sigma * (T_max + T_amb)*(T_max**2 + T_amb**2)
        h_eff = h_ave + h_rad
        
        Bi = h_eff * r2 / self.k_fin
        
        lamda = math.pi + 1/ (r_ratio * math.sqrt(math.pi))
        phi = (math.tanh(lamda * tau) + (lamda/Bi)) / (1 + (lamda/Bi)*math.tanh(lamda * tau))
        psi = r_ratio * tau / math.sqrt(math.pi) + (1-r_ratio) * phi / math.sqrt(math.pi) # dimensionless constriction resistance
        R_sp = psi / (self.k_fin * r1 * math.sqrt(math.pi))
        
        # print(r1, r2, r_ratio, tau, h_ave, epsil, sigma, h_rad, h_eff, Bi, lamda, phi, psi, R_sp)
        
        return R_sp
    
    def cal_dp(self, V_flow): # dP 계산시에는 m로 환산해놨음

        L = unitconv("length", "mm", "m", self.L).cal() # m
        H = unitconv("length", "mm", "m", self.H).cal()
        b = unitconv("length", "mm", "m", self.b).cal()
        L_fin = unitconv("length", "mm", "m", self.L_fin).cal()
        t = unitconv("length", "mm", "m", self.t).cal() # m
        S = unitconv("length", "mm", "m", self.S).cal() # m
        N = self.N
        W = unitconv("length", "mm", "m", self.W).cal() # m
        # T_avg = self.T_amb + 273.15
        
        air = self.Air_prop(df_air, "Density", self.T_amb) # air properties 호출
        
        D_hyd = 2 * S * L_fin / (S + L_fin) # m
        a_ratio = S/L_fin # aspect ratio
        D = air[0] # kg/m3
        Vis_dyn = air[1][7] * 10**-6 # dynamic viscosity, Pa-s (=kg/m-s)
        A_meter = (N-1)*S*L_fin # m2, fin 사이의 속도.
        # A_meter = unitconv("length", "mm", "m", self.W).cal() * unitconv("length", "mm", "m", self.H).cal() # m2, Duct 속도
        A_feet  = unitconv("area", "m2","ft2", A_meter).cal() # ft2
        
        Step = 10
        LastStep = 0.0001
        direction = 1

        # V_flow  = 25.16         # volumetric air flow rate (ft3/min)   Fan 에서 불러오기
        Vel_feet   = V_flow / A_feet # ft/min
        Vel_meter  = unitconv("length", "ft", "m", Vel_feet).cal() / 60 # m/s
        
        Re = D * Vel_meter * D_hyd / Vis_dyn
        f  = (24 - 32.527*a_ratio + 46.721*a_ratio**2 - 40.829*a_ratio**3 + 22.954*a_ratio**4 - 6.089*a_ratio**5) /Re # friction factor for fully developed laminar flow
        L_char = (L/D_hyd) / Re # Characteristic Length
        f_app = ( (3.44/math.sqrt(L_char))**2 + (f * Re)**2)**(1/2) / Re # apparent friction factor
        sigma = 1 - (N*t)/W # ratio of the area of the flow channels to that of the flow approaching the heat sink
        Kc = 0.42 * (1 - sigma**2) # Pressure losses due to sudden constraction of the flow entering the heatsink
        Ke = (1 - sigma**2)**2 # Pressure losses due to sudden expansion of the leaving the heatsink
        
        dp_pa  = (Kc + 4*f_app*L/D_hyd + Ke) * D * Vel_meter**2/2 # Pressure drop across the heatsink, Pa(N/m2)
        dp_mmWTR = unitconv("pressure", "Pa", "mmWTR", dp_pa).cal()
        
        # print(dp_mmWTR, Vel_meter)
        # print("D_hyd :", D_hyd, "a_ratio :", a_ratio, "Vel_feet", Vel_feet, "Re :", Re, "f : ", f, "L_char : ", L_char, "f_app :", f_app, "sigma :", sigma, "Kc :" , Kc, "Ke :", Ke, "dp_pa :", dp_pa)
        return dp_mmWTR, Vel_meter, dp_pa
    
    def T_disch(self, V_flow): # Heatsink discharge temperature
        
        air = self.Air_prop(df_air, "Pressure",self.T_amb) # T_max 해당되는 air properties 호출
        V_flow_meter = unitconv("volume", "ft3", "m3", V_flow).cal() * 60 # m3/h
        
        m_air = V_flow_meter * air[1][1] # kg/h
        
        
        ## find proper Cp values
        T_disch = 0
        Step = 10
        LastStep = 0.001
        direction = 1
        
        for i in range(0, 1000, 1):
            T_disch = T_disch + Step * direction
            T = (self.T_amb + T_disch)/2 # in/out ave. temperature
            
            heat_air = m_air * self.Air_prop(df_air, "Cp", T)[0] * (T_disch - self.T_amb) / 3.6
            # print("Cal :", heat_air, " vs ", "Power : ", self.tot_power)
            # print("T_disch :", T_disch)
            if self.tot_power < heat_air :
                if direction == 1:
                    if Step > LastStep:
                        Step = Step / 10
                        direction = direction * -1
                    else:
                        break
            else:
                if direction == -1:
                    if Step > LastStep:
                        Step = Step / 10
                        direction = direction * -1
                    else:
                        break
        # T_disch = self.T_amb + self.tot_power*3.6 / (m_air * self.Air_prop(df_air, "Cp", T)[0])

        return T_disch

# Target heat sink temperature 
def FindValue(mode, solver, chip, heatsink, v_flow=0): # 0: Tmax 1: dP
##  mode 0 : Tmax, 1:dP
    
    if mode == 0: ## Tmax
        Step = 10
        LastStep = 0.00001
        direction = 1
        T_amb = Cal1.T_amb
        T_dif = Cal1.T_dif
        power = chip.tot_power
    
        for i in range(0, 10000, 1):
            T_dif = T_dif + Step * direction
            solver.setdata(T_amb, T_dif, chip, heatsink)
            solver.cal_temp(v_flow)
            
            cal_value = solver.cal_temp(v_flow)[0]
            # print("Step :", Step, "direction :", direction)
            # print("T_target : ", T_dif)
            # print("Cal :", solver.cal_temp(v_flow)[0], " vs ", "Power : ", power)
    
            if power < cal_value:
                if direction == 1:
                    if Step > LastStep:
                        Step = Step / 10
                        direction = direction * -1
                    else:
                        break
            else:
                if direction == -1:
                    if Step > LastStep:
                        Step = Step / 10
                        direction = direction * -1
                    else:
                        break
    
        sol = T_dif + T_amb
        # print("Temperature(C) :", sol)
        return sol
    
    elif mode == 1: # pressure drop
        Step = 10
        LastStep = 0.00001
        direction = 1
        
        press = 1 # 초기 압력, mmWTR
        
        for i in range(0, 1000, 1):

            # Fan curve 에 해당되는 flow/pressure 와 dp 계산했을 때의 dp가 일치하는 flow 찾기
            press = press + Step * direction
            # 최대나 최소 압력 벗어날 때
            if press >= max(fan_press):
                press = max(fan_press)
            elif press <= min(fan_press):
                press = min(fan_press)
                
            flow = fan(press).Flow() # Fan_curve 에서 press 에 해당하는 flow 찾기
            # print(flow)
            # dp = solver.cal_dp(flow) # Fan_curve 에서 찾은 flow를 적용했을 때의 pressure loss
                        
            cal_value = solver.cal_dp(flow)[0]
            # print("Step :", Step, "direction :", direction)
            # print("dp : ", press)
            # print("Cal :", solver.cal_dp(flow)[0], " vs ", "press : ", press)
    
            if press > cal_value:
                if direction == 1:
                    if Step > LastStep:
                        Step = Step / 10
                        direction = direction * -1
                    else:
                        break
            else:
                if direction == -1:
                    if Step > LastStep:
                        Step = Step / 10
                        direction = direction * -1
                    else:
                        break
        err = press - cal_value
        
        if err > 5:
            print("Please check the Pressure drop and Volumetric flow ")
            print("Please check the Pressure drop and Volumetric flow ")
            print("Please check the Pressure drop and Volumetric flow ")
            print("Please check the Pressure drop and Volumetric flow ")
            print("Please check the Pressure drop and Volumetric flow ")
            
        sol = press
        # print("Pressure(mmH2O) :", press,"Flowrate(ft3/min) :", flow)
        return sol    

 
#####################################################################
######################## Calculation sample #########################
#####################################################################

# ### Chip 정보 입력 (power, N, L, W)
# chip1 = Chip(200, 1, 36, 51)
# ## chip2 = Chip(100, 1, 36, 51)

# ## Heatsink 정보 입력 #H, L, b t, N, S, k_fin # H, L, b, t, N, W, k_fin
# heatsink1 = Heatsink(23, 95, 6, 1, 37, 1.5, 167) # H, L, b, t, N, S, k_fin

# T_ambient = 50 # Input 값

# Cal1 = Calculation()
# Cal1.setdata(T_ambient, 10, chip1, heatsink1) # 초기 temp diff. 10도

# Sol_dp = FindValue(1, Cal1, chip1, heatsink1)     # dP와 Fan flow 일치하는 값 찾기. 결과는 pressure
# V_flow = fan(Sol_dp).Flow()            # 해당 시스템에 적용될 Fan volumetric flow 값
# Sol_temp = FindValue(0, Cal1, chip1, heatsink1, V_flow) # 평균 온도 계산
# HS_therm_resis = (Sol_temp - T_ambient) / chip1.tot_power
# # print(HS_therm_resis)
# Rsp = Cal1.cal_rsp(V_flow)
# T_spot = Sol_temp + (Rsp*chip1.tot_power)# 최대 온도 계산

# ###### Data summary
# ColName = ["H", "L", "b", "L_fin", "t", "N", "S", "W", "Heatsink Temp.","Hot spot temp.", "Rsp", "Therm. Resis.", "Pressure Loss(Pa)", "T_disch(C)", "Vol_flow(m3/s)", "Nu_m", "Pr", "Re_D"]
# init  = np.zeros([1,len(ColName)])
# list_heatsink = pd.DataFrame(init, columns=ColName)

# dp = unitconv("pressure", "mmWTR", "Pa", Sol_dp).cal()
# Vol_f = unitconv("volume", "ft3", "m3", V_flow).cal()/60  # m3/s
# Rsp = Cal1.cal_rsp(V_flow)
# T_spot = Sol_temp + (Rsp*chip1.tot_power)
# Nu_m = Cal1.cal_temp(V_flow)[2][24]
# Pr = Cal1.cal_temp(V_flow)[2][20]
# Re_D = Cal1.cal_temp(V_flow)[2][22]
# T_disch = Cal1.T_disch(V_flow)
# Check_resis = (Sol_temp - T_ambient) / chip1.tot_power

# list_heatsink.loc[1]=[heatsink1.H, heatsink1.L, heatsink1.b, heatsink1.L_fin, heatsink1.t, heatsink1.N, heatsink1.S, heatsink1.W, Sol_temp, T_spot, Rsp, Check_resis, dp, T_disch, Vol_f, Nu_m, Pr, Re_D]
# list_heatsink.drop([0], inplace = True)  # 초기 조건 삭제


#####################################################################
###################### Heatsink 사양 찾기 ############################
#####################################################################
#### Heatsink fin 길이(H), fin 개수만 변경해서 요구 값 찾기

chip1 = Chip(200, 1, 36, 51)
T_ambient = 50 # Input 값

target_temp  = 100 # heatsink 타겟 온도, ℃
s_factor     = 0   # margin, %
target_resis = (target_temp- T_ambient) / chip1.tot_power * (100-s_factor)/100  

Cal1 = Calculation()

ColName = ["H", "L", "b", "L_fin", "t", "N", "S", "W", "Heatsink Temp.","Hot spot temp.", "Rsp", "Therm. Resis.", "Pressure Loss(Pa)", "T_disch(C)", "Vol_flow(m3/s)", "Nu_m", "Pr", "Re_D"]
init  = np.zeros([1,len(ColName)])
list_heatsink = pd.DataFrame(init, columns=ColName)
n = 0

for i in range(10, 17, 3): # (10, 17, 3): (75, 116, 10):
    i = i/10 ## thickness 할때만 int 로 받아야 하므로 i/10  
    
    for j in range(10, 32, 5): # (10, 32, 5): (4, 9, 1):
        # heatsink1 = Heatsink(23, i, j, 1, 37, 91, 167) # H, L, b, t, N, W, k_fin
        
        heatsink1 = Heatsink(23, 95, 6, i, j, 91, 167) # H, L, b, t, N, W, k_fin
        Cal1.setdata(T_ambient, 10, chip1, heatsink1)
        Sol_dp = FindValue(1, Cal1, chip1, heatsink1)
        V_flow = fan(Sol_dp).Flow() #ft3/min
        Sol_temp = FindValue(0, Cal1, chip1, heatsink1, V_flow) # heat sink 평균온도 계산
        Rsp = Cal1.cal_rsp(V_flow)
        T_spot = Sol_temp + (Rsp*chip1.tot_power)   ## 최대 온도로 thermal resistance 계산
        Check_resis = (T_spot - T_ambient) / chip1.tot_power 
        print(str(i).center(5), str(j).center(5), "Check_resis : ", "{:.6f}".format(Check_resis) ,"  ", "Target : ", "{:.6f}".format(target_resis))
        
        if Check_resis < target_resis:
            n = n+1
            dp = unitconv("pressure", "mmWTR", "Pa", Sol_dp).cal()
            Vol_f = unitconv("volume", "ft3", "m3", V_flow).cal()/60  # m3/s
            Nu_m = Cal1.cal_temp(V_flow)[2][24]
            Pr = Cal1.cal_temp(V_flow)[2][20]
            Re_D = Cal1.cal_temp(V_flow)[2][22]
            T_disch = Cal1.T_disch(V_flow)
            list_heatsink.loc[n]=[heatsink1.H, heatsink1.L, heatsink1.b, heatsink1.L_fin, heatsink1.t, heatsink1.N, heatsink1.S, heatsink1.W, Sol_temp, T_spot, Rsp, Check_resis, dp, T_disch, Vol_f, Nu_m, Pr, Re_D]

list_heatsink.drop([0], inplace = True)  # 초기 조건 삭제

#####################################################################
#####################################################################
#####################################################################








# Velocity = Cal1.cal_dp(V_flow)[1]  # 0 :dp, 1:vel

# Sol_temp = FindValue(0, Cal1, chip1)

# dif_press = Cal1.cal_dp(17.62)
# print(dif_press)

### 호출 예제 ###
# aaa = Cal1.Air_prop(df_air, "Pressure") # Properties 
# aaa = Cal1.cal() # Q_tot, T_max
# print(Cal1.cal()[0])


# Solution = FindValue(0, Cal1)

##################################        
## Factor comparison
# Para_colName =  np.array(['T_max', 'Q_tot', 'h1', 'Re', 'Nu', 'A1', 'Q_c1_max', 'Q_c1_base_max',' Q_c1_fin_max',' eff_c1_fin','m_c1', 'L_c', 'T',' D',' Cp',' Vis_dyn',' Vis_kin',' k',' alpha',' beta',' Pr',' D_hyd',' Re_D',' Gz', 'Nu_m', 'h2', 'A2', 'Q_c2_max', 'Q_c2_base_max', 'Q_c2_fin_max', 'eff_c2_fin', 'sigma', 'Q_r1', 'A_r2', 'Q_r2'])
# Para = Cal1.cal_temp()[2]
# df_para = pd.DataFrame(Para)
# df_para = df_para.transpose()
# df_para.columns = Para_colName
# ##################################

[JungHulk]

Convective heat transfer coefficient

https://etd.ohiolink.edu/apexprod/rws_etd/send_file/send?accession=ohiou1292521397&disposition=inline

[JungHulk]

Simple method to estimate air flow pass

https://www.electronics-cooling.com/1997/05/a-simple-method-to-estimate-heat-sink-air-flow-bypass/

[JungHulk]

"air prop table"

# -*- coding: utf-8 -*-
"""
Created on Thu Jan 19 15:39:48 2023

@author: seungj.jung
"""

#### Excel Data base 만들기
import os
import pandas as pd
# from pandas import Series
import pickle


# excel = Dispatch('excel.Application')  # win32com.client로 excel call
openfiledir = os.getcwd() + '\\' # 현재 폴더 위치 
filename = 'Air_properties.xlsx'  # 현재 폴더 위치의 해당 파일명의 .xls 파일 open

df = pd.read_excel(openfiledir + filename)
IndName = df['Temperature']
IndName[0] = "(Units)"
df.index = IndName    # Index 지정

df.drop([1], inplace = True)
df.drop(['Unnamed: 0','Temperature'], axis = 1, inplace = True)
df.drop(df.index[1], inplace = True)
ColName = df.columns  # Column 지정

with open('Air_prop.p','wb') as file:
    pickle.dump(df, file)

[JungHulk]

"Unitconv_pickle_maker"

# -*- coding: utf-8 -*-
"""
Created on Thu Jan 26 08:55:44 2023

@author: seungj.jung
"""

# -*- coding: utf-8 -*-
"""
Created on Thu Jan 19 15:39:48 2023

@author: seungj.jung
"""

#### Excel Data base 만들기
import os
import pandas as pd
import pickle
import numpy as np


# excel = Dispatch('excel.Application')  # win32com.client로 excel call
openfiledir = os.getcwd() + '\\' # 현재 폴더 위치 
filename = 'Unitconv.xlsx'  # 현재 폴더 위치의 해당 파일명의 .xls 파일 open

unit = np.array(["length", "area", "volume", "pressure"])
for i, n in enumerate(unit):
    globals()['df_{}'.format(n)] = pd.read_excel(openfiledir + filename, sheet_name = n)   # 각 sheet name 으로 dataframe 만들기
    globals()['Indname_{}'.format(n)] = globals()['df_{}'.format(n)][format(n)]            # Dataframe 별 index 지정
    globals()['df_{}'.format(n)].index = globals()['Indname_{}'.format(n)]                 # Dataframe 에 index 설정
    globals()['df_{}'.format(n)].drop(format(n), axis = 1, inplace = True)                 # Dataframe 첫 열 삭제
    
    # for 문으로 저장하는건??
    
# df_length = pd.read_excel(openfiledir + filename, sheet_name = "length")
# Indname_length = df_length["length"]
# df_length.index = Indname_length

# df_length.drop('length', axis =1, inplace = True)


with open('Unitconv_length.p','wb') as file:
    pickle.dump(df_length, file)

with open('Unitconv_area.p','wb') as file:
    pickle.dump(df_area, file)    

with open('Unitconv_volume.p','wb') as file:
    pickle.dump(df_volume, file)
    
with open('Unitconv_pressure.p','wb') as file:
    pickle.dump(df_pressure, file)

[JungHulk]

Unitconv.py

# -*- coding: utf-8 -*-
"""
Created on Tue Jan 31 18:05:11 2023

@author: seungj.jung
"""

import os
import pandas as pd
import pickle

unitconv_length = pd.read_pickle('Unitconv_length.p')
unitconv_area = pd.read_pickle('Unitconv_area.p')
unitconv_volume = pd.read_pickle('Unitconv_volume.p')
unitconv_pressure = pd.read_pickle('Unitconv_pressure.p')

class unitconv:

    def __init__(self,selection, fromunit, tounit, value):
        self.type  = selection # length, area, volume
        self.funit = fromunit  # 기존 단위
        self.tunit = tounit    # 변경하고자 하는 단위
        self.value = value     # 값
        
    def cal(self):
        if self.type == "length":
            factor = unitconv_length.loc[self.funit, self.tunit]
            sol = factor * self.value
            
        elif self.type =="area":
            factor = unitconv_area.loc[self.funit, self.tunit]
            sol = factor * self.value
            
        elif self.type =="volume":
            factor = unitconv_volume.loc[self.funit, self.tunit]
            sol = factor * self.value
            
        elif self.type =="pressure":
            factor = unitconv_pressure.loc[self.funit, self.tunit]
            sol = factor * self.value
        
        return sol

      
## a = unitconv("length","mm", "ft", 50)
## b = a.cal()

## a = unitconv("length","mm", "ft", 50).cal()

[JungHulk]

Fan_curve.py

# -*- coding: utf-8 -*-
"""
Created on Wed Feb  1 10:49:05 2023

@author: seungj.jung
"""

## P-Q curve for electric Fan


import os
from openpyxl import load_workbook
import numpy as np
import pandas as pd
import pickle
import matplotlib.pyplot as plt

openfiledir = os.getcwd() + '\\' # 현재 폴더 위치 
filename = 'Fan_curve.xlsx'      # 현재 폴더 위치의 해당 파일명의 .xls 파일 open

df = pd.read_excel(openfiledir + filename)
# 1행, 1열 삭제
df.drop([1], inplace = True)
df.drop(['Unnamed: 0'], axis = 1, inplace = True)


ColName = df.loc[0]
df.columns = ColName
df.drop([0], inplace = True)

IndName = df["Pressure"]
df.drop(["Pressure"], axis = 1, inplace = True)
df.index = IndName

# with open('Fan_curve.p','wb') as file:
#     pickle.dump(df, file)

# # a = df.loc[55.04][0]

class Fan:
    def __init__(self, pressure):
        self.pressure = pressure
        
    def Flow(self):
        
        if self.pressure > max(df.index) or self.pressure <0:     #Temperature range 설정 (-40 ~ 80)
            print("Out of range")
            value = []
            
        else:
            for i in range(len(df.index)):
                if self.pressure == df.index[i]:
                    value = df.loc[self.pressure][0]
                elif self.pressure < df.index[i] and self.pressure > df.index[i+1]:
                    value = df.loc[df.index[i]][0] + (self.pressure - df.index[i]) * (df.loc[df.index[i+1]][0] -df.loc[df.index[i]][0]) / (df.index[i+1] -df.index[i])
                  
        return value

## 호출
# Flow = Fan(pressure).Flow()

class Plot:          
    def Curve():
       
        step = 0.01
        flow = []
        pressure = []

        if df.index[0] > df.index[-1]:
            
            for i in range(int(df.index[-1]*100), int(df.index[0]*100), 1):
                y = i *step
                x = Fan(y).Flow()
 
                flow.append(x)
                pressure.append(y)
        else:
            for i in range(int(df.index[0]*100), int(df.index[-1]*100), 1):
                y = i *step
                x = Fan(y).Flow()
                
                flow.append(x)
                pressure.append(y)
        
        fig = plt.figure(figsize=(12, 8))
        ax = fig.add_subplot(1, 1, 1)
        ax.plot(flow, pressure)
        
        ## plot setting
        plt.rc('xtick', labelsize = 18)
        
        ax.set_xlabel('Volumetric flow, Q (ft3/min)', fontsize=16)
        ax.set_ylabel('Pressure, P (mmH2O)', fontsize=16)
        ax.set_title('Fan P-Q Performance curve', fontsize=20, loc='center') # 'center', 'right'
        
        plt.show()
        
        return flow, pressure
        
## 호출
## Fan_curve = Plot.Curve()
        
                

[JungHulk]

example

start = time.time() # 시작시간
for i in range(10, 26, 5): # H
    for j in range(2, 9, 2): # b
        for k in range(10, 16, 1): # t
            k = k/10 ## thickness 할때만 int 로 받아야 하므로 i/10  
            for l in range(20,41,5): # N
                heatsink1 = Heatsink(i, 100, j, k, l, 90, 167) # H, L, b, t, N, W, k_fin
                Cal1.setdata(T_ambient, 10, chip1, heatsink1)
                Sol_dp = FindValue(1, Cal1, chip1, heatsink1)
                V_flow = fan(Sol_dp).Flow() #ft3/min
                Sol_temp = FindValue(0, Cal1, chip1, heatsink1, V_flow) # heat sink 평균온도 계산
                Rsp = Cal1.cal_rsp(V_flow)
                T_spot = Sol_temp + (Rsp*chip1.tot_power)
                Check_resis = (T_spot - T_ambient) / chip1.tot_power  
                print(str(i).center(5), str(j).center(5), str(k).center(5), str(l).center(5), "Check_resis : ", "{:.6f}".format(Check_resis) ,"  ", "Target : ", "{:.6f}".format(target_resis))

                if Check_resis < target_resis:
                    n = n+1
                    dp = unitconv("pressure", "mmWTR", "Pa", Sol_dp).cal()
                    Vol_f = unitconv("volume", "ft3", "m3", V_flow).cal()/60  # m3/s
                    Nu_m = Cal1.cal_temp(V_flow)[2][24]
                    Pr = Cal1.cal_temp(V_flow)[2][20]
                    Re_D = Cal1.cal_temp(V_flow)[2][22]
                    T_disch = Cal1.T_disch(V_flow)
                    list_heatsink.loc[n]=[heatsink1.H, heatsink1.L, heatsink1.b, heatsink1.L_fin, heatsink1.t, heatsink1.N, heatsink1.S, heatsink1.W, Sol_temp, T_spot, Rsp, Check_resis, dp, T_disch, Vol_f, Nu_m, Pr, Re_D]

end = time.time() # 계산 종료 시간
print(f"{end - start : .5f} sec") # 계산 시간
list_heatsink.drop([0], inplace = True)  # 초기 조건 삭제