Optimization of An Air-Based Heat Management System for Dusty Particulate Matter-Covered Lithium-Ion Battery Packs

Summary

Here, we present the adaptive simulated annealing method (ASAM) to optimize an approximate quadratic response surface model (QRSM) corresponding to a dusty particulate matter-covered battery heat management system and fulfill the temperature drops back by adjusting the airflow velocities combination of system inlets.

Abstract

This study aims to solve the problem of the cell temperature rise and performance decline caused by dusty particulate matter covering the surface of the cell through the allocation of airflow velocities at the inlets of the battery cooling box under the goal of low energy consumption. We take the maximum temperature of the battery pack at a specified airflow velocity and dust-free environment as the expected temperature in a dusty environment. The maximum temperature of the battery pack in a dusty environment is solved at different inlet airflow velocities, which are the boundary conditions of the analysis model constructed in the simulation software. The arrays representing the different airflow velocity combinations of inlets are generated randomly through the optimal Latin hypercube algorithm (OLHA), where the lower and upper limits of velocities corresponding to the temperatures above the desired temperature are set in the optimization software. We establish an approximate QRSM between the velocity combination and the maximum temperature using the fitting module of the optimization software. The QRSM is optimized based on the ASAM, and the optimal result is in good agreement with the analysis result obtained by the simulation software. After optimization, the flow rate of the middle inlet is changed from 5.5 m/s to 5 m/s, and the total airflow velocity is decreased by 3%. The protocol here presents an optimization method simultaneously considering energy consumption and thermal performance of the battery management system that has been established, and it can be widely used to improve the life cycle of the battery pack with minimum operating cost.

Introduction

With the rapid development of the automobile industry, traditional fuel vehicles consume a lot of non-renewable resources, resulting in serious environmental pollution and energy shortage. One of the most promising solutions is the development of electric vehicles (EVs)^1,2.

The power batteries used for EVs can store electrochemical energy, which is the key to replacing traditional fuel vehicles. Power batteries used in EVs include lithium-ion battery (LIB), nickel-metal hydride battery (NiMH), and electric double-layer capacitor (EDLC)³. Compared to the other batteries, lithium-ion batteries are currently widely used as energy storage units in EVs owing to their advantages such as high energy density, high efficiency, and long life cycle^4,5,6,7.

However, due to chemical reaction heat and Joule heat, it is easy to accumulate a large amount of heat and increase the battery temperature during rapid charging and high-intensity discharging. The ideal operating temperature of LIB is 20-40 °C^8,9. The maximum temperature difference between the batteries in a battery string should not exceed 5 °C^10,11. Otherwise, it may lead to a series of risks such as temperature imbalance between the batteries, accelerated aging, even overheating, fire, explosion, and so on¹². Therefore, the critical issue to be resolved is designing and optimizing an efficient battery thermal management system (BTMS) that can control the temperature and temperature difference of the battery pack within a narrow.

Typical BTMS include air cooling, water cooling, and phase change material cooling¹³. Among these cooling methods, the air cooling type is widely used because of its low cost and simplicity of the structure¹⁴. Due to the limited specific heat capacity of air, high temperature and large temperature differences are easy to occur between battery cells in air-cooled systems. In order to improve the cooling performance of air-cooled BTMS, it is necessary to design an efficient system^15,16,17. Qian et al.¹⁸ collected the battery pack's maximum temperature and temperature difference to train the corresponding Bayesian neural network model, which is used to optimize cell spacings of the series air-cooled battery pack. Chen et al.¹⁹ reported using the Newton method and the flow resistance network model for optimization of the widths of the inlet divergence plenum and the outlet convergence plenum in the Z-type parallel air-cooled system. The results showed a 45% reduction in the temperature difference of the system. Liu et al.²⁰ sampled five groups of the cooling ducts in the J-BTMS and obtained the best combination of cell spacings by the ensemble surrogate-based optimization algorithm. Baveja et al.²¹ modeled a passively balanced battery module, and the study described the effects of thermal prediction on module-level passive balancing and vice versa. Singh et al.²² investigated a battery thermal management system (BTMS) that used encapsulated phase change material along with forced convective air cooling designed using the coupled electrochemical-thermal modeling. Fan et al.²³ proposed a liquid cooling plate comprising a multi-stage Tesla valve configuration to provide a safer temperature range for a prismatic-type lithium-ion battery with high recognition in microfluidic applications. Feng et al. ²⁴ used the coefficient of variation method to evaluate the schemes with different inlet flow rates and battery clearances. Talele et al.²⁵ introduced wall-enhanced pyro lining thermal insulation to store potential generated heating based on optimal placement of heating films.

When one uses air-cooling BTMS, metal dust particles, mineral dust particles, building materials dust particles, and other particles in the external environment will be brought into the air-cooling BTMS by the blower, which can cause the surface of the batteries to be covered with DPM. If there is no heat dissipation plan, it may cause accidents due to the excessively high battery temperature. After simulation, we take the maximum temperature of the battery pack in a specified airflow velocity and dust-free environment as the expected temperature in a dusty environment. First, C-rate refers to the current value required when the battery releases its rated capacity within the specified time, which is equal to a multiple of the battery's rated capacity in the data value. In this paper, the simulation uses 2C rate discharge. The rated capacity is 10 Ah, and the nominal voltage is 3.2 V. Lithium iron phosphate (LiFePO4) is used as the positive electrode material, and carbon is used as the negative electrode material. The electrolyte has electrolyte lithium salt, a high-purity organic solvent, necessary additives, and other raw materials. The random array representing the different velocity combinations at the inlets was determined through the OLHA, and a 2nd order function between the maximum temperature of the battery pack and the inlet flow velocity combination was set up under the condition of checking the accuracy of the curve fitting. Latin hypercube (LH) designs have been applied in many computer experiments since they were proposed by McKay et al.²⁶. An LH is given by an N x p-matrix L, where each column of L consists of a permutation of the integers 1 to N. In this paper, the optimal Latin hypercube sampling method is used to reduce the computational burden. The method uses stratified sampling to ensure that the sampling points can cover all the sampling internals.

In the following step, the inlet flow velocity combination was optimized to decrease the maximum temperature of the battery pack in a dusty environment based on the ASAM under the condition of considering energy consumption simultaneously. The adaptive simulated annealing algorithm has been extensively developed and widely used in many optimization problems^27,28. This algorithm can avoid getting trapped in a local optimum by accepting the worst solution with a certain probability. The global optimum is achieved by defining the acceptance probability and temperature; the calculation speed can also be adjusted by using these two parameters. Finally, for checking the accuracy of the optimization, the optimal result was compared with the analysis result obtained from the simulation software.

In this paper, an optimization method for the inlet flow rate of the battery box is proposed for the battery pack whose temperature rises due to dust cover. The purpose is to reduce the maximum temperature of the dust-covered battery pack to below the maximum temperature of the non-dust-covered battery pack in the case of low energy consumption.

Protocol

NOTE: The research technology roadmap is shown in Figure 1, where the modeling, simulation, and optimization software are used. The materials required are shown in the Table of Materials.

1. Creating the 3D model

NOTE: We used Solidworks to create the 3D model.

1. Draw a 252 mm x 175 mm rectangle, click Extrude Boss/Base, and enter 73. Create a new plane 4 mm from the outer surface.
2. Draw a rectangle 131 mm x 16 mm and click Linear Sketch Pattern. Enter 22 and 6 in spacing and number of instances, respectively. Select all four sides of the rectangle and click OK. Enter 180 in angle and run it again. This step is for symmetry in the center of the model.
3. Click Extrude Cut, enter 65, and click OK. Click Extrude Boss/Base and enter 65, uncheck Merge result, and click Reverse Direction and OK.
   NOTE: When the merge result is unchecked, the stretched entity becomes a separate entity. There are 23 parts in total, including 11 batteries, 11 dusty particulate matter, and 1 air domain.
4. Draw a rectangle 16 mm x 1 mm. Repeat steps 1.2 and 1.3.
5. Draw a rectangle 63 mm x 15 mm, click the rectangle's top edge and Linear Sketch Pattern. Enter 21, 3, and 270, and click OK. Click Split Line and the face of the cube, click OK.
6. Draw a rectangle 63 mm x 15 mm. Click Split Line and the face of the cube, click OK.
7. Click File and save it as an X_T file.
   NOTE: The specified size: L_box:73 mm; W_box:252 mm; H_box:175 mm; L_b, L_d:65 mm; W_b, W_d:10 mm; H_b:131 mm; H_d:1 mm; L_i:63 mm; W_i:15 mm; d₁, d₂:5 mm, d₃:6 mm is shown in Figure 2.
8. Drag the mesh component by clicking Toolbox > Component Systems > Mesh to the project schematic zone. Import the previously saved X_T file by clicking Geometry.
9. Enter the mesh-design modeler window, and the battery pack model, including 23 parts as independent bodies, is displayed again by clicking Generate.
10. Select all 23 parts of the battery to be a new part named as battery part, all dusty particulate matters of 23 parts as dust part, and air cavity as air part, in the tree outline for the convenience of subsequent hiding and naming objects.
11. First, right-click on BatteryPart and DustPart and select the Hide Part so that the pop-up will only show the air part.
12. Move the mouse to the selection toolbar to select Selection Filter: Bodies, right-click on the air cavity model on the graphics zone to select Named Selection, and rename the air cavity model on the details view zone as air domain.
13. Switch to Selection Filter: Faces, right-click and rename the surface that has been split into three pieces, from bottom to top, as inlet1, inlet2, and inlet3, the separate surface to the right of these three faces is named outlet, the remaining outer surface is named outerBorder, respectively.
14. Switch the Select Mode to the Box Select, click the Y axis to obtain the suitable view of the air cavity model for the convenience of box selecting, rename and number all inner surfaces as cavity surface1 to cavity surface11 using box selection.
15. In order to show only the batteryPart, right-click airPart and select Hide Part. Right-click batteryPart and select Show Part on the pop-up shortcut menu.
16. Move the mouse to the selection toolbar to select Selection Filter: Bodies, switch the Select Mode to Single Select, right-click on Each Battery Model on the Graphics zone to select the Named Selection, rename and number the 11 battery models on the details view zone as batteryDomain1 to batteryDomain11, respectively.
17. Furthermore, each battery model has six sides, then switch to Selection Filter: Faces, right-click on each Side of the Numbered batteryDomains to select Named Selection and rename them according to the orientation of the battery side. For example, rename six sides of the numbered batteryDomain1 as batteryDomain1_Upper, batteryDomain1_Lower, batteryDomain1_Left, batteryDomain1_Right, batteryDomain1_Front and batteryDomain1_Back.
18. In order to show only the dustPart, right-click batteryPart and select Hide Part. Right-click dustPart and select Show Part on the pop-up shortcut menu.
19. Move the mouse to the selection toolbar to select the Selection Filter: Bodies, right-click each dusty particulate matter model on the Graphics zone to select Named Selection, rename and number the 11 dusty particulate matter models on the details view zone as dpmDomain1 to dpmDomain11, respectively.
20. Furthermore, each dusty particulate matter model has six sides; then switch to the Selection Filter: Faces, right-click on Each Side of the Numbered dpmDomains to select Named Selection and rename them according to the orientation of the dusty particulate matter side. For example, rename six sides of numbered dpmDomain1 as dpmDomain1_Upper, dpmDomain1_Lower, dpmDomain1_Left, dpmDomain1_Right, dpmDomain1_Front, and dpmDomain1_Back.
21. Show all bodies and return to the initial window again.

2. Generate the mesh model

NOTE: Finite element meshing is a very important step in finite element numerical simulation analysis, which directly affects the accuracy of subsequent numerical analysis results. The renamed entities are then meshed.

1. In order to mesh the air domain, battery domain, and dpm domain independently, drag two Mesh components again from Toolbox > Component Systems > Mesh to the project schematic zone and rename them as airFEM, batteryFEM, and dpmFEM, respectively. Hold the airFEM > Geometry with the left mouse button and drag it to the batteryFEM > Geometry.
2. Next, hold the batteryFEM > Geometry with the left mouse button and drag it to the dpmFEM > Geometry. Right-click the Lines among the three mesh components and select Delete to disassociate them from each other.
3. Double-click airFEM's Mesh, enter the meshing window, right-click batteryPart and dustPart to select the Suppress Body, and change the physical preference from mechanical to CFD. Generate the FEM air domain model through the face sizing of 2 mm and the body sizing of 4 mm by clicking Update and return to the initial window.
4. Double-click batteryFEM's Mesh, enter the meshing window, right-click airPart and dustPart to select the Suppress Body, and change the physical preference from Mechanical to CFD. Generate the FEM battery domain model through the body sizing 2 mm by clicking Update and return to the initial window.
5. Double-click dpmFEM's Mesh, enter the meshing window, right-click airPart and batteryPart to select Suppress Body, and change the physical preference from Mechanical to CFD. Generate the FEM dpm domain model through body sizing 2 mm by clicking Update, return to the initial window.
   NOTE: Figure 3A shows the grid of the air domain, Figure 3B shows the grid of the battery domain, and Figure 3C shows the grid of the dpm domain.
6. Set the minimum size of the air grid to 4 mm and the minimum size of the battery and dusty particulate matter grid to 2 mm. Ensure that the grid is solution independent, change the minimum cell size of the grid, and perform a grid sensitivity study.
   NOTE: As shown in Figure 4, with the number of grids increasing from 519343 to 1053849, the maximum battery temperature changes are less than 0.6 K. Considering the computation ability and accuracy, the following analysis is based on the grid model with 931189 grids.

3. Simulation analysis

1. Drag Fluid Flow from Toolbox > Analysis Systems > Fluid Flow into the project schematic zone. Hold airFEM > Mesh, then batteryFEM > Mesh and dpmFEM > Mesh with the left mouse button and drag them to Fluid Flow > Setup. Right-click Fluid Flow > Setup and select Update to enter the set window.
2. Verify the validity of the FEM model and check whether the mesh has a negative volume. The software automatically suggests the volume of the model, and a reasonable model value is positive. If there is any problem with the divided grid or model settings, an error message will pop up to tell.
3. Activate the energy equation in heat transfer models. Enter the setting interface of the viscous model and the radiation model and select the K-epsilon Model and the Discrete Ordinates Model.
   NOTE: As shown in Figure 5, comparing four viscous models, the calculation results of the Spalart-Allmaras model are quite different from those of other models. The results of the Standard K-epsilon model are like those of other K-epsilon models. The Standard K-epsilon model with higher stability and economy is widely used; the following analysis is based on the Standard K-epsilon model.
4. Set the new materials with different attributes for air material, battery material, dpm material, and battery box material based on Table 1.
   NOTE: Inside the battery pack, there are three different physical materials: air as a fluid and the rest as a solid. Next, set up the material.
   1. Change the fluid type of the numbered battery domains to the Solid type and change the dpm material to the battery material on the Solid window by double-clicking each Battery Domain. Subsequently, choose the Source Terms item and click the highlighted Source Terms to add an energy source by assigning the number in the number of energy sources and selecting Constant type to input the value of 209993 w/m³.
   2. Change the fluid type of the numbered dpm domains to Solid type.
5. Next, set the interface for simulation calculation of several different domains according to the actual setting flow rate and heat transfer coefficient as described below.
   1. Convert the type of all renamed surfaces, including the inner surfaces of the air domain and all sides of the battery domains, as well as dpm domains from the default wall to the interface. Once the above steps are finished successfully, the mesh interfaces will be generated immediately.
   2. Click the Mesh Interfaces and enter the Create/Edit Mesh Interfaces window. Match the cavity surfaces to all sides except the battery domains' upper sides and the dpm domian's lower sides. Next, name and number them as interface1 to interface11, respectively. So, the 11 mesh interfaces can be created among the air domain and battery domians as well as dpm domains.
   3. Match the battery domains' upper sides and the dpm domains' lower sides. Next, name and number them as interface12 to interface22, respectively. Then, the 11 mesh interfaces are created between the battery domains and dpm domains.
   4. Assign the surface of the outer border as the wall thermal boundary condition by setting the heat transfer coefficient as 5 in the mixed thermal condition and changing its material from default aluminum to the previously self-defined battery box material.
   5. Set airflow velocities of all inlets as 5 m/s in the velocity inlet window and the gauge pressure of the outlet as zero in the pressure outlet window.
6. Next, set the state of the computing domain at the initial moment, such as the initial temperature of 300 K, which will affect the process of computing convergence.
   1. Set the type of solution initialization as the standard initialization before initializing.
   2. Set the Number of Iterations as 2000.
   3. Click Calculate to simulate. Return to the initial window until the simulation is finished.
7. The above part completes the simulation calculation of the temperature and air velocity inside the battery pack and then displays the simulation result in Result. Perform the following steps in the results displayed.
   1. Double-click Fluid Flow > Results to enter the CFD post window, then click the icon of Contour in the toolbox.
   2. Select All Sides of the Batteries in the location selector and change pressure to temperature. Then click Apply to generate the temperature contour of the batteries.
   3. Click File > Export to select the temperature of the selected variable(s). Click the Dropdown button of the locations to pop up the location selector window where all battery domains should be selected. Click OK and Save button to quit.
      NOTE: A spreadsheet whose data corresponds to the temperatures of all batteries' mesh nodes will be saved automatically when the save button is clicked.
   4. Open the spreadsheet to find the maximum value, which indicates the maximum temperature of the batteries in a dusty environment at 5 m/s of all airflow inlets.
   5. Acquire the maximum temperature of the batteries under the free-dust state as the expected temperature and compare it with the maximum temperature under the dusty state; the result shows the entire temperature increasing.
      NOTE: To acquire the maximum temperature of batteries in a dust-free environment, the new battery pack model shown in Figure 6 should be re-established, and all the steps 1.1-3.4.3 should be repeated.
   6. In order to lower the maximum temperature inside the battery pack, set the airflow velocities at the inlets from 5 m/s to 6 m/s, increase by 5%, and calculate the corresponding maximum temperatures of the dusty-covered batteries.
      NOTE: The sensitivity analysis of airflow velocity parameters should be done well in advance before changing the parameter values. As shown in Figure 7 and Table 2, we have kept the same total flow for each of the seven groups of different inlet airflow velocity combinations. There still be an obvious variation in the maximum temperature owing to the difference in airflow velocity allocation. In other words, there is somehow a strong correlation between the airflow velocity combination and the maximum temperature. Therefore, those velocity parameters can be used as design variables.
   7. Plot the temperature-velocity curve as shown in Figure 8, where the red line indicates the temperature characteristic curve decreases with the increase of airflow velocity, and the blue line represents the expected temperature.
   8. Maintain an increase in airflow velocity of 10%. When the velocity increment is more than 10%, the maximum temperature is already lower than the expected temperature, but this does not meet the purpose of low energy consumption. For the remaining air flow rate, reduce the maximum temperature of the battery pack to the expected temperature through optimization, thus achieving the goal of low energy consumption.

4. Optimal Latin hypercube sampling and response surface modeling

NOTE: For the retained flow rates of 5 m/s-5.5 m/s, samples are selected to construct different flow rate combinations within this flow rate range. The velocity combinations are simulated to obtain the maximum temperature. Construct the function of velocity and maximum temperature.

1. Open a new empty spreadsheet to create a table whose rows in the first column are named inlet1, inlet2, and inlet3, and save the file as sampling.xlsx.
2. Run the optimization software and drag the Spreadsheet icon onto the single arrow of Task 1. Next, double-click the Spreadsheet icon to pop up the Component Editor-Excel window.
3. Import the sampling.xlsx by clicking the Browse button and map the inlet1, inlet2, and inlet3 to the A1, A2, and A3 as parameters by clicking the Add this mapping. Click the OK button to return to the initial window.
4. Drag the DOE icon into Task1 and double-click it to pop up the Component Editor-DOE window. Select the OptimOKal Latin Hypercube and set the Number of Points as 15 in the General window.
5. Switch to the Factors window and set 5.5 as the upper limit and 5 as the lower limit for A1, A2, and A3.
6. Switch to the Design Matrix window and click Generate to generate the random sampling points corresponding to the different inlet velocities. Shut down optimization software.
7. Take the velocity combinations arrays of the random sampling points back to calculate and repeat steps 3.5.5-3.7.5 to obtain the corresponding temperature array composed of the maximum temperatures of batteries.
8. Combine the predictor variables x1, x2, and x3 of the velocity combinations arrays and y of the temperature arrays to form a new table of variables, as shown in Table 3, and save it as a sample.txt file. Import the file to fit a response surface model.
9. Rerun the optimization software and drag the Approximation icon onto the single arrow of Task1. Double-click the Task1 icon to pop up the component editor-approximation window to select the Response Surface Model.
10. Switch to the Data File window and import the sample.txt file containing the prediction variables.
11. Switch to the Parameters window and click the Scan to open the parameters in the data file window where the predictor variables of x1, x2, and x3 are defined as input and y as output.
12. Switch to the Technique Options window and select the Quadratic in polynomial order. Switch to the Error Analysis Options window and select the Cross-Validation in the error analysis method.
13. Switch to the View Data window and click Initialize Now to obtain the coefficients of the quadratic linear regression equation.
14. Click the Error Analysis button to pop up the approximation error analysis window to check whether the errors can satisfy with the acceptable standards for each error type. Close the approximation component window. If the arbitrary error cannot satisfy the corresponding acceptable standards, then add more sample points to participate in the model fitting.

5. Adaptive simulated annealing algorithm-based approximate fitting model

NOTE: Next, software and algorithm are used to find the optimal value of the approximate model

1. Drag the Optimization icon into Task1 and double-click it to pop up the component editor-optimization window. Select the Adaptive Simulated Annealing (ASA) in the optimization technique.
2. Switch to the Variables window to set 5.5 as the upper limit and 5 as the lower limit.
3. Switch to the Objectives window and select the Y parameter before closing the component editor-optimization window.
4. Click the Run Optimization button and wait for the optimization result.

Representative Results

Following the protocol, the first three parts are the most important, which include modeling, meshing, and simulation, all in order to get the maximum temperature of the battery pack. Then, the airflow velocity is adjusted by sampling, and finally, the optimal flow rate combination is obtained by optimization.

Figure 9 shows the comparison of battery pack temperature distribution in different environments, and Figure 10 shows the comparison of the second battery temperature distribution in different environments. As shown in Figure 9 and Figure 10, the temperature of the battery under the dusty state is increased to a certain level due to the low thermal conductivity of DPM (dusty particulate matter). 

In order to adjust the battery temperature distribution, set the airflow velocities at the inlets from 5 m/s to 6 m/s, increase by 5% under the dusty model, and obtain the maximum temperatures at each airflow velocity. When the airflow velocity was increased by 15% and 20%, the maximum temperature of the battery pack under the dusty state dropped below the maximum temperature of the battery pack under the free-dust state, as shown in Figure 8. Considering energy consumption, the maximum inlet velocity is set as 5.5 m/s (increased by 10%) to decrease the maximum temperature of the battery pack in the dusty state.

When establishing the quadratic QRSM, the minimum number of samples is calculated by (N + 1) x (N + 2)/2, where N is the number of test variables. There are three design variables in this article, which are the inlet velocities and the minimum number of samples is 10. In order to establish a response surface model with high fitting accuracy, 15 samples were selected using the DOE component of the optimization software platform. The least square method is used to complete the fitting of the response surface between the maximum temperature of the battery pack obtained by the simulation software and three inlet velocities. The approximated response surface model is established as follows:

R2 measures the overall fit of the regression equation and expresses the overall relationship between the dependent variable and all independent variables. R2 is equal to the ratio of the regression sum of squares to the total sum of squares, that is, the percentage of the variability of the dependent variable that the regression equation can explain. The closer the value of R2 is to 1, the better the fit of the regression curve to the observed value.

The error analysis of the calculation results shows that R2 is 0.93127, as shown in Figure 11, which shows that the second-order polynomial response surface approximation model has a good fitting accuracy.

In the end, adaptive simulated annealing (ASA) is used as the optimization method for finding optimal inlet flow velocity combinations. The maximum number of generated designs is 10,000, the number of designs for convergence check is 5, and the convergence epsilon is 1.0 x 10^-8. The relative rate of parameter annealing, cost annealing, parameter quenching, and cost quenching were the same value of 1.

The maximum temperature of the battery pack obtained by optimization was 309.391420 K. The inlets' air flow velocities are 5.5 m/s, 5m/s, and 5.5 m/s. To confirm the accuracy, the optimal case was analyzed by the simulation software. Table 4 shows the comparison between the optimization and simulation verification results. It can be seen that the error of the maximum temperature of the battery pack is within 0.001% under three inlet airflow velocities conditions, which indicates that the optimization method adopted in this work is effective and feasible.

The comparison of the second battery temperature distribution under the different inlet airflow velocities is shown in Figure 12, and the comparison of battery pack temperature distribution before and after optimization is shown in Figure 13. Table 5 shows the specific values of the maximum temperatures and the combinations of airflow velocities. When the airflow velocities of inlets 1-3 are 5.5 m/s, 5.5 m/s, and 5.5 m/s, respectively, the maximum temperature of the battery pack is 309.426208 K. After optimization, the airflow velocity of inlets 1-3 are 5.5 m/s, 5m/s, and 5.5 m/s, and the maximum temperature of the battery pack is 309.392853 K. It should be noted that the sum of airflow velocities of the optimized case shown in Figure 12B is less than the sum of airflow velocities of the case shown in Figure 12A. However, the maximum temperature does not increase with decreasing airflow velocity. Also, the optimized battery pack is compared with the initial battery pack (that is, the airflow velocities of the three inlets are all 5 m/s, and the batteries are covered with DPM). Figure 14 compares the flow line distribution before and after optimization, and it can be seen that the flow line distribution after optimization is wider. Figure 15 compares the effects of each factor on temperature; factor x1 has the greatest influence on temperature. Factors x1 and x3 have similar effects on temperature. In a word, the total airflow velocity decreases by 3%, and the maximum temperature of the battery pack is decreased to the expected temperature (that is, the maximum temperature of the battery pack under a dust-free state).

The optimization method can be widely used to improve the life cycle of the battery pack with low energy consumption.

Figure 1: The technical roadmap. This figure describes the detailed simulation and optimization process according to the research content, including research objects, methods, solutions, modeling, simulation, and optimization software. Please click here to view a larger version of this figure.

Figure 2: A 3D model of lithium-ion battery pack in a dusty environment. The 3D model of the LIB pack, which can be saved as an X_T file and imported into simulation software to simulate, is drawn by modeling software. Please click here to view a larger version of this figure.

Figure 3: Grid diagram. (A) This figure shows the grid of the air domain. (B) This figure shows the grid of the battery domain. (C) This figure shows the grid of the dpm domain. Please click here to view a larger version of this figure.

Figure 4: Grid independence test. The X-axis is the different total number of grids in the mesh model, and the Y-axis is temperature. Please click here to view a larger version of this figure.

Figure 5: Viscous model test. The X-axis is the type of viscous model, the number 1 represents the Standard k-epsilon model, the number 2 represents the RNG k-epsilon model, the number 3 represents the Realizable k-epsilon model, the number 4 represents the Spalart-Allmaras model, the Y-axis is temperature. Please click here to view a larger version of this figure.

Figure 6: 3D model of lithium-ion battery pack in a dust-free environment. The 3D model of the LIB pack, which can be saved as an X_T file and imported into simulation software to simulate, is drawn by modeling software. Please click here to view a larger version of this figure.

Figure 7: Parameter sensitivity analysis. The number on the x-axis represents the nth combination of inlet airflow velocities. For example, the number 5 represents the velocity combination (3,5,7) corresponding to 3 m/s at inlet1, 5 m/s at inlet2, 7 m/s at inlet3. Similarly, number 1,2,3,4,6 represents the different inlets air flow velocity combination of (5,5,5), (4,5,6), (5,6,4), (5,4,6), (3,5,7), (5,3,7), (5,7,3), respectively. The Y-axis is temperature. Please click here to view a larger version of this figure.

Figure 8: Battery pack temperature variation at different inlet airflow velocities. The figure shows the maximum battery pack temperature decreasing with the increase of inlet airflow velocity. The x-axis is the rate of airflow velocity increase at inlets. The Y-axis is temperature. Please click here to view a larger version of this figure.

Figure 9: Comparison of battery pack temperature distribution in different environments. (A) This figure shows the temperature distribution of the battery pack in a dust-free environment. (B) This figure shows the temperature distribution of the battery pack in a dusty environment, from which the temperature is highest in the number 2 battery. Please click here to view a larger version of this figure.

Figure 10: Comparison of number 2 battery temperature distribution in different environments. (A) This figure shows the temperature distribution of the number 2 battery in a dust-free environment. (B) This figure shows the temperature distribution of the number 2 battery in a dusty environment. Please click here to view a larger version of this figure.

Figure 11: Error analysis of approximation response surface model. The figure indicates the quadratic polynomial response surface approximation model has good fitting accuracy. Please click here to view a larger version of this figure.

Figure 12: Comparison of the number 2 battery temperature distribution under different inlet airflow velocities. (A) This figure shows the temperature distribution of the number 2 battery by just increasing the inlet airflow velocity itself. (B) This figure shows the temperature distribution of the number 2 battery after optimization of inlet airflow velocity. Please click here to view a larger version of this figure.

Figure 13: Comparison of battery pack temperature distribution before and after optimization. (A) This figure shows the temperature distribution of the battery pack in a dusty environment without optimization. (B) This figure shows the temperature distribution of the battery pack in a dusty environment after optimization. Please click here to view a larger version of this figure.

Figure 14: Comparison of battery pack streamline distribution before and after optimization. (A) This figure shows the streamlined distribution of the battery pack in a dusty environment without optimization. (B) This figure shows the streamlined distribution of the battery pack in a dusty environment after optimization. Please click here to view a larger version of this figure.

Figure 15: Influence of three factors on temperature. (A) This figure shows the effects of x1 and x2 on temperature. (B) This figure shows the effects of x1 and x3 on temperature. Please click here to view a larger version of this figure.

Name of the medium	ρ/kg·m-3	C/J·(kg·K)-1	K/W (m·K)-1
air Material	1.225	1006.43	0.0242
battery Material	1958.7	733	kx=3.6,ky=kz=10.8
dpm Material	2870	910	1.75
batterybox Material	7930	500	16.3

Table 1: Material properties. The material properties corresponding to the air, battery, dusty particulate matter, and battery box will be used in the parameter Settings of the simulation software.

Number	Inlet1(m/s)	Inlet2(m/s)	Inlet3(m/s)	Maximum temperature of battery pack (K)
1	5	5	5	309.72049
2	4	5	6	309.26413
3	5	6	4	309.703369
4	5	4	6	309.389038
5	3	5	7	311.54599
6	5	3	7	308.858704
7	5	7	3	309.801086

Table 2: Parameter sensitivity analysis. The table shows the seven combinations of inlet airflow velocities and the corresponding maximum temperature of the battery pack. For example, the number 5 represents the velocity combination (3,5,7) corresponding to 3 m/s at inlet1, 5m/s at inlet2, 7 m/s at inlet3, and the corresponding battery pack maximum temperature of 311.54599 K.

Number	Inlet1(m/s)	Inlet2(m/s)	Inlet3(m/s)	Maximum temperature of battery pack (K)
1	5.071	5.429	5.179	309.58725
2	5.286	5.071	5.036	309.59982
3	5.393	5.143	5.429	309.48029
4	5.464	5.25	5.071	309.52237
5	5.179	5.036	5.25	309.59082
6	5.143	5.107	5.5	309.50894
7	5.5	5.357	5.321	309.46039
8	5.107	5.393	5.464	309.52564
9	5.036	5.179	5.107	309.64923
10	5.214	5.321	5	309.59052
11	5.321	5.5	5.393	309.48645
12	5.357	5.464	5.143	309.5264
13	5.429	5	5.214	309.50253
14	5	5.214	5.357	309.58344
15	5.25	5.286	5.286	309.54627

Table 3: Velocity and temperature arrays used for quadratic response surface model. The different airflow velocity combinations at inlets can be randomly generated by the OLHA, and the corresponding maximum temperatures are calculated by the simulation software.

Name	Inlet1(m/s)	Inlet2(m/s)	Inlet3(m/s)	Maximum temperature of battery pack (K)
Optimization result	5.5	5	5.5	309.39142
Simulation verification result	5.5	5	5.5	309.392853

Table 4: Comparison between the optimization and simulation verification results. The suitable airflow velocity combination at inlets and corresponding temperature can be obtained by optimizing, which is also proved to be accurate by the simulation verification.

Name	Inlet1(m/s)	Inlet2(m/s)	Inlet3(m/s)	Maximum temperature of battery pack (K)
A	5	5	5	309.412537
B	5	5	5	309.72049
C	5.5	5.5	5.5	309.426208
D	5.5	5	5.5	309.392853

Table 5: Comparisons of inlets air flow velocity and maximum temperature of the battery pack under different conditions. (A) The battery pack under the normal inlets air flow velocity and free-dust environment. (B) The battery pack under the normal inlets air flow velocity and dusty environment. (C) The battery pack under the inlets air flow velocities increase and dusty environment. (D) The battery pack under the optimized airflow velocities and dusty environment.

Discussion

 The BTMS used in this study was established based on the air-cooling system due to its low cost and simplicity of the structure. Because of the low heat transfer capacity, the performance of the air-cooling system is lower than that of the liquid cooling system and phase change material cooling system. However, the liquid cooling system has the disadvantage of refrigerant leakage, and the phase change material cooling system has high mass and low energy density²⁹. These cooling systems have their advantages and disadvantages. Therefore, the BTMS can be established by combining an air-cooling system with a liquid cooling system or a phase change material cooling system to promote cooling performance.

A CFD solver was implemented to simulate the flow and temperature profile of the model. The governing equations^30, such as continuity (2) and the energy conservation equation (3), were employed to solve the time-dependent thermal problem of the airflow.

Where p, k, and c are the properties of air employed, which are density, thermal conductivity, and specific heat, respectively; T, and  are the static pressure, temperature, and velocity of the cooling air.

Momentum equations³¹

Where u_i and u_j are Reynolds-averaged velocity components; x_i and x_j are cartesian coordinates; P is Reynolds-averaged pressure; μ is dynamic viscosity; μ_t is turbulent dynamic viscosity. k is turbulent kinetic energy; ε is turbulent kinetic energy dissipation rate.

The Reynolds number based on the inlet flow velocity (v=5 m/s) and the equivalent diameter was estimated to be 0.0242308; the Reynolds number is calculated as 9894, thereby a turbulence model of the standard k-e model was selected.

Reynolds number equation³²

Where P_l is the density, V_max is the maximum flow velocity of the liquid, D is the equivalent diameter of the container, and u_l is the dynamic viscosity of liquid.

Turbulent kinetic energy equation³³

Where k_t and ε is the turbulent kinetic energy and turbulence dissipation rate, respectively; u_j is the jth component of the velocity vector, and μ and u_t are the molecular and turbulent dynamic viscosity, respectively; the G_kt and G_b are the turbulent kinetic energy generation caused by mean velocity and the turbulent kinetic energy generation as a result of buoyancy effects, respectively; Y_M represents the influence of the fluctuating dilation incompressible turbulent to the sum of dissipation rates; S_kt is source term of k_t; α_kt is the inverse effective Prandtl number for k_t.

Turbulent kinetic energy dissipation equation³³

Where S_ε is source term of ε; α_t is the inverse effective Prandtl number for ε; C_1ε , C_2ε and C_3ε are empirical constants.

For the battery cells, the energy conservation equation³⁴

Where Q, k_b, c_b; and P_b represent the generated heat, thermal conductivity, specific heat capacity, and density of the battery, respectively.

Heat convection formula³⁵

Where h_f represents convection heat transfer coefficient; T_s represents the surface temperature of LIBs; T_B represents the temperature of ambient air; and q* represents convection heat transfer rate.

The inlet of the BTMS was set to a velocity-inlet boundary condition of 5 m/s and temperature of 300 K while the system outlet was conditioned to pressure-outlet with the surrounding pressure set to atmospheric pressure. The walls around the system are set for natural convection.

This paper began the research under the condition that the structure of battery pack model was determined, dust covering the surface of the battery will cause the temperature of the battery to rise. Then we present the ASAM to optimize an approximate QRSM and fulfill the temperature drops back through the optimal air flow velocities combination of system inlets for solving the problem of DPM effect. It should be mentioned that the positions of the air inlet and outlet of the battery pack also have a great influence on the temperature of the BTMS¹⁴.

There are some critical steps in protocol. When creating the 3D model of the battery pack, give each body and surface in the model a recognizable name for subsequent material addition material, creation of mesh interface and setting of boundary conditions. When operating the simulation software, it is necessary to set each parameter accurately, especially the unit of the parameter.

In terms of fitting model, error analysis is significant in response surface modeling, if the arbitrary error could not satisfy the corresponding acceptable standards, then more sample points should be added to participate in the model fitting until the error reaches the acceptable standards. After the simulation software imports the grid model, troubleshoot the mesh model, click Check to check whether the mesh has a negative volume. If there is any problem with the divided grid or model settings, an error message will pop up.

The main limitation of this study is that the geometric model used in the simulation is derived by simplifying the realistic battery pack model, it's almost impossible to fully reflect reality. Then, the boundary conditions imposed are unlikely to be consistent with the actual situation. The calculation results are also different according to different calculation theories. To facilitate the simulation, we simplified the heat generation model of the battery, the average heat generation rate of the battery is 20.993 kW/m³ as the internal heat source^36,37.

The significance concerning existing methods and any future applications of the technique:

This protocol will help establish an optimization method while simultaneously considering energy consumption and thermal performance of the battery management system, and it can be widely used to improve the life cycle of the battery pack with minimum operating cost. This technique can also be used in mechanical design, architectural design and other fields.

Materials

Name	Company	Catalog Number	Comments
Ansys-Workbench	ANSYS	N/A	Multi-purpose finite element method computer design program software.https://www.ansys.com
Isight	Engineous Sogtware	N/A	Comprehensive computer-aided engineering software.https://www.3ds.com
NVIDIA GPU	NVIDIA	N/A	An NVIDIA GPU is needed as some of the software frameworks below will not work otherwise. https://www.nvidia.com
Software
SOLIDWORKS	Dassault Systemes	N/A	SolidWorks provides different design solutions, reduces errors in the design process, and improves product quality www.solidworks.com