CFD Windshield Deicing Simulations for Commercial Vehicle Applications.
Government safety regulations have continued to grow in response to making automotive and commercial vehicles safer for drivers and passengers alike. With an average of 5,760,000 vehicle crashes occurring per year, and 1,259,000 of them being weather-related per the Federal Highway Administration, it is all the more imperative to mitigate any additional hindrances to the tractor driver's ability to navigate in adverse weather, particularly icy conditions . One area that has seen tremendous growth is leveraging computational fluid dynamics (CFD) to assist in the design and validation process of the HVAC system, particularly defrost duct performance. In the past twenty years, the "science" of windshield deicing simulation has been formulated, tested, and proven as a valuable means of not only reducing part design time, but also making the duct design more robust . With CFD, a robust and complete design can be developed prior to ever being prototyped and tested in an experimental facility. Coupled with the exponential growth in computational hardware speed , simulation of this problem has gradually incorporated more geometric details, greater mesh refinement, and more robust physics. Within this paper, the field of deicing simulations will be introduced, followed by the integration of multiple modeling tools that is then demonstrated for a two duct design. Results will show how simulation prediction of experimental performance will meet the design engineer's requirements, specifically the percentage of a given glass surface ("zone") that is free of ice.
SAE Deicing Criteria
The problem being considered in this investigation is extending the CFD work conducted for an automotive vehicle's windshield deicing performance and applying it to a commercial vehicle. As in automotive vehicles, commercial vehicle windshield deicing is also critical to the driver's safety and other motorists' safety. One could argue that it is even more important to address in the commercial vehicle space as the risk to other motorists is much higher due to sheer size and payload of the truck or tractor . Performance of deicing is measured based on SAE J381 test procedures, which dictate what percentage of the windshield and/or side windows need to be clear of ice after a certain time . It should be noted that SAE's test procedure is more stringent than in the real world, where the driver would attempt to scrape the ice off of the glass surfaces prior to driving . The windshield is divided into three "zones" shown schematically in Figure 1 to which the SAE deicing criteria is applied. The SAE J381 test procedure dictates the size and shape of each of these three zones depending on the vehicle type . Zones A and C are specifically called out in the SAE J381 test procedure, while Zone B is an intermediate area generated per the HVAC design engineer's request.
In order of criticality, Section C is the area relating to the driver's visibility. After 30 minutes, this section of the windshield needs to be clear of ice (99% clear). Section B is the next larger area that now takes into account the passenger's view out of the windshield. There are no SAE criteria tied to deicing this section of the windshield. The final section called out in the SAE deicing criteria is Section A, the widest section which will allow the driver the greatest safety when operating the vehicle. In 30 minutes, 80% of this area needs to be clear of ice. These windshield zones were "developed to be compatible with vision requirements necessary to operate trucks, buses, and multipurpose vehicles [and] is based on SAE J941 with certain modifications to accommodate the wide variety of conditions encountered in these vehicles" . Of course, HVAC design engineers may require or demand a higher level of performance than those outlined by SAE's specifications. For instance, the HVAC design engineer who set deicing requirements for the current project required Zone A be 85% free of ice in 30 minutes, and, after 45 minutes, Zone A needs to be 95% clear of ice. Similar criteria exist for the side windows, though not as intensive as outlined above for the windshield. After 45 minutes, the driver side window needs to be 65% clear of ice, and the passenger side window needs to be 70% clear of ice. With the regulatory constraints outlined, the next portion of the problem to consider is the design and packaging of the HVAC system (primarily the ductwork).
HVAC System Packaging
Other than physically packaging the HVAC module, optimizing the pressure drop across the ducting from the HVAC module to the instrument panel (IP) of the vehicle and proper airflow rate into the cabin are the biggest obstacles for an engineer to consider. A basic space claim for the ductwork already exists, as space constrained by the IP is at a premium. The engine and other peripheral components already limit the most direct duct routing options, so pressure drop can only be addressed with subtle changes, such as smoothing bends, adjusting the diameter of the ducting, adjusting outlet count, and other design implementations. The next item requiring the engineer's attention is managing the airflow properly to the cabin from the HVAC unit. Here, the SAE deicing criteria must be taken into consideration as well. If the airflow is not biased properly, the windshield may experience adverse deicing behavior, such as key areas not receiving enough heated air which would lead to a "patchy" deicing pattern rather than clear zones like those shown in Figure 1. As the SAE requirements are biased towards driver visibility, the flow distribution should be biased towards the driver as well; for example, 50-55% driver-side airflow, 45-50% passenger-side airflow. This distribution can be addressed through the physical instrument panel duct design, as the vehicle operator cannot vary the airflow distribution for deicing manually. The outlet duct quantity plays a large role in achieving a reasonable airflow split as a "first pass" measure. The outlet count (located at the base of the windshield), vane angle, opening area, and HVAC module architecture are four key parameters which affect airflow delivery to the windshield. Once an initial deicing pattern is simulated or tested, adjusting the outlet vanes can fine tune the flow to eliminate any remaining frozen zones on the windshield that must be addressed to meet the SAE targets. The aforementioned problems set the stage for leveraging CFD (computational fluid dynamics) simulation as a means to not only frontload the design process (and therefore expedite it), but also allows for a robust design to be utilized in the field.
Time Benefits of CFD Simulation
By addressing this problem in the digital space, engineers are able to determine beforehand how their designs will perform prior to any need for physical parts or testing. Each portion of the problem can be addressed sequentially and holistically. For instance, if the deicing duct does not meet pressure drop requirements, the design can be iterated on until satisfactory performance is achieved. Another benefit of addressing the problem computationally is giving the engineer the ability to gradually add complexity to the simulation as required and "shorten the development cycle" . For instance, if only pressure drop and air distribution ("cold flow" steady state simulation) are required, no transient deicing simulations need to be run. Depending on the computational hardware available, results for this type of simulation can be achieved within hours. This allows the design engineer time to address necessary design changes the same day after reviewing CFD results. Physical testing limits the amount of tests per day as a result of other tests being booked for a given day or time slot. Also, if an engineer desires to see the duct's deicing performance, a transient deicing simulation can be run utilizing the cold flow simulation results as a base model, thus leveraging the knowledge gained during that steady state simulation. In contrast, a physical test offers an "all or nothing" approach. Probes can collect data for HVAC duct pressure drop tests and deicing pattern contours can be collected during a transient deicing test, but a retest would be required to gather any additional data the design engineer may need later on (not to mention rebooking the test chamber for a date potentially much farther into the future).
Visual Benefits of CFD Simulation
CFD also provides a visual benefit that physical testing cannot. A design engineer is able to see results, such as pressure contours and velocity flow fields within the ductwork, allowing for an informed design decision. Physical testing can show a design engineer only that which can be inferred from quantifiable data. A design engineer would not be able to see a flow path of air going through the deicing duct or easily determine a breakthrough point (where ice first begins to melt) on a windshield. The transient results from the simulations can only be achieved by physical testing with intermittent interruption for data collection, which is labor intensive. A major benefit of simulating the transient deicing process is the ability to see the boundary of ice melting on the windshield and truly understand the effectiveness of flow fields and temperature gradients. Video loops of the process during the deicing simulation can be generated and allow a design engineer to see and quantify the ice melted at any time interval as opposed to relying on still images from a physical test that may require manual capturing by a test engineer. CFD puts all these tools at the design engineer's fingertips, and allows them to validate their computational model in the physical test chamber with greater confidence than going to the test chamber directly from a concept. An investigation of the procedural methods developed in simulating the deicing problem will be discussed below.
History of Deicing Simulations
Over the past fifteen years, the deicing problem has been studied in great detail as a result of stringent government regulations with respect to driver/occupant safety. Entities such as SAE and regulations such as Federal Motor Vehicle Safety Standards (FMVSS) have generated deicing criteria for not only passenger vehicles, but also commercial vehicles, requiring design engineers to focus on developing a robust HVAC system . Turning to CFD simulation as a means to reduce cost of physical testing as well as optimize designs at the beginning of the design cycle became a necessity as the years progressed. As stated by Wolfahrt et al., "they [simulation methods] have reached a mature stage where they can contribute successfully to a significant extent to meet the vehicle development targets" . In the early 2000s, hardware limitations forced engineers to be more judicious in selecting not only a proper domain, but a corresponding mesh size and type that would strike a balance between computation time, accuracy, and geometric fidelity. Due to the steady increase in computing power, engineers are no longer as limited by mesh size or number of simulation iterations, as high speed computing and cluster parallelization have enabled simulations to be conducted even faster. Additionally, the ever decreasing cost of CPUs allows engineers even greater access to sheer volume of computing power, which has doubled every two years for the past fifty years according to Moore's Law . In the sections below, the deicing simulation process will be discussed under the lens of historical SAE papers on the topic. These sections will map the process from CAD generation and manipulation, to meshing schemes, steady state (cold flow) simulation setups, and transient deicing simulation setups.
CAD Optimization Simulations
To begin the historical study, a closer look at CAD generation was investigated. For instance, Pak and Simon studied the benefits of utilizing a parametric cabin modeler as opposed to creating raw CAD . Within the cabin modeler, 200 parameters defined the cabin geometry which in theory could limit the fidelity of the design and accuracy of the simulation when compared to actual cabin CAD. "The cabin modeler can capture overall geometry and some detailed cabin configurations; however, very complicated minor details such as knobs, bumps, and levers in the cabin were excluded in the model" . Omissions of intricate features could be particularly dangerous from an accuracy perspective if certain IP topography was not accurately captured (as this would directly impact airflow and therefore deicing performance). The benefits to parameterizing the cabin, however, included its ability to model various automobile type cabins, making it more of a template library that design engineers could draw from instead of starting with a blank canvas. Additionally, Pak and Simon embedded a meshing functionality to the modeler utilizing tetrahedral elements and prism elements to further speed up the meshing process.
Mesh Element Optimization
Prior research compared the deicing results from the cabin modeler with those generated when utilizing true CAD to thermal imaging results from a physical test. It was determined that the cabin modeler method succeeded in achieving an acceptable level of accuracy . The observation was made that "the velocity magnitude of the cabin modeler model was slightly off from the value of the detailed model results. However, the case that the detailed model prediction deviated from the test results cannot be excluded." The paper compares data generated by the cabin modeler-generated CAD, detailed raw CAD, and the test data. The "raw CAD" was generated from scratch, whereas the cabin modeler generates the appropriate CAD based on measurement parameters input into the Cabin Modeler GUI. In other words, the raw CAD will capture all details the engineer creates in the CAD package utilized, while the Cabin Modeler will automate the geometry creation process and capture only key details (such as windshield heights, seat heights, foot well lengths, etc.) to make a reasonable "clone" of the cabin. As to be expected, the simulations did not exactly match experimental test data, but provided an acceptable level of matching trends regardless of each model's shortcomings. One important detail should be noted; the raw CAD model had a total of1.8 million hexahedral elements while the cabin modeler model had 3 million tetrahedral elements . When comparing mesh counts, it makes sense that the cabin modeler would experience no major loss to accuracy as the mesh count was nearly double in size compared to the CAD model. Additionally, tetrahedral elements would return more accurate results for a flow simulation by default, as they discretize the domain more thoroughly than hexahedral elements. In industry today, hexahedral elements are utilized as a means to reduce mesh count in model areas that either experience noncomplex flow behavior or are not value-add areas to computational results. Though not discussed in the paper, it can be inferred that the time saved utilizing the parametric cabin modeler was cancelled out by the computation time required to solve the heavier model. When comparing meshing style utilized by the other sources discussed, using a tetrahedral element volume meshing scheme was more common, similar to one shown in Figure 2. Only  utilized a hexahedral volume meshing scheme.
There are multiple ways to incorporate the flow field into the cabin model. For simplicity, the inlet flow to the cabin from the ductwork can be integrated as a mathematical boundary condition to solve for flow between more complex cabin features . Conversely, one could include the discretized CFD solution ductwork into the simplified cabin model to achieve similar results. However, Farag and Huang stated that either scheme can be utilized for the cabin flow domain, suggesting that the choice of scheme may be more dependent on computational limitations rather than accuracy of results.
Fluid Domain Mesh Count Optimization
Another common method utilized in controlling the mesh count proposed by AbdulNour is cropping the main fluid domain (in this case, the cabin) at a location where important fluid effects no longer play a relevant role in the flow solution. "Usually, it is adequate to use a computational domain that extends from the inlet of the defroster nozzle, where the velocity is assumed to be uniform, to some sufficient downstream location in the vehicle interior, where the velocity becomes negligibly small and the flow is assumed to have constant pressure" . By "cropping" the domain, early deicing simulation development was able to be conducted on hardware that was not as robust as that which is available to CFD engineers today. Another issue encountered by engineers in the early years of CFD modelling was handling design iterations of the defrost duct.
Assuming one fluid domain was utilized for the entire simulation, the volume mesh would need to be regenerated due to the movement of solid ductwork components. To address this problem, splitting the fluid domain into two distinct zones was introduced by . "This approach enables independent control of mesh densities. It also makes the mesh structure very modular without having to mesh the entire passenger compartment every time there is a design change in the defrost duct" . This method requires the use of a non-conformal interface to "link" the two zones together during the simulation. Ideally, using non-conformal interfaces is undesirable from an accuracy perspective due to the required interpolation between the two zones. However, Hill and Sringari have suggested in their article that there is no major impact on accuracy of the results . The CAD cleanup and mesh for the duct can be conducted in a much lighter file, and then migrated into the flow solver during the final stages of simulation setup. Additionally, as Hill and Sringari pointed out, splitting the fluid domain into sub-domains allows for greater control of the meshing scheme, adding mesh density where it is needed (such as the defrost duct) without unduly complicating the mesh in areas where the added level of discretization would not add any meaningful value (such as the rear of the cabin where no complex flow is occurring). Using one fluid domain for the entire problem forces the engineer to make a decision on speed versus accuracy. A lower mesh count would expedite the simulation run time, but may return invalid or unusable results. Too high a mesh count would waste elements on areas not needing ample discretization but would yield a fairly accurate result. While determining what geometry is incorporated in the model is the "what" question of the model, meshing is the "how" question that takes greater time and understanding of the physics to be simulated.
Duct Mesh Preparation for Deicing Simulations
As previously discussed, the major flow domain meshing scheme is a simpler decision when compared to capturing localized flow effects in the cabin. As a majority of the flow physics are occurring in the ducts and near the windshield, discretizing the effects in these regions are extremely important, as the simulations may yield inaccurate values or deicing behavior as a result. Within the scope of deicing simulations, a strong concentration of mesh elements near high pressure gradient locations was required . Therefore, a higher mesh concentration within the defrost duct and near the defrost duct outlet and windshield was required to resolve the fluid flow behavior properly. Two methodologies for defrost duct meshing exist. Patidar recommends using prism layers inside the defrost duct to capture "the air flow boundary layer effect" . Patidar gives a closer view of the defrost duct volume mesh and the associated prism layers for added discretization in the boundary layer along the walls . Farag and Haung, however, make no mention of treating the defrost duct with prism layers inside the duct. This seems to be once again a matter of computational resources available, as Patidar further explains that prism layers were neglected in the defrost nozzle locations due to an element size of 1.5 mm. "[1.5 mm elements are] enough in size of volume mesh to capture the boundary layer effect in that area" . One can infer that, if the mesh density is large enough, once could potentially avoid utilizing prism layers in the duct as well, which would bolster the reasoning why Farag and Huang did not mention this technique in their paper. Not only does meshing the defrost duct require special attention, but meshing the windshield and ice requires additional attention as well.
Windshield Mesh Preparation for Deicing Simulations
Patidar proposes two methods to mesh the windshield, a shell modeling method and a solid modeling method. In the shell modeling method, the thickness of the glass is handled in a 1-D fashion in which only shell elements are created for these surfaces. They are "defined [as] conductive wall[s] with virtual thickness as input conditions" . Only Wolfahrt appears to have utilized this technique for the windshield . All other sources utilize Patidar's solid modeling method (shown in Figure 2) where all glass surfaces which will be covered in ice are modeled as prism layers . As opposed to the shell modeling method, the thickness of the ice and glass is physically represented in the thickness dimension of a prism layer. Farag and Huang recommend utilizing the solid modeling method as "the thickness of [the] ice layer and windshield thickness is small compared to the overall passenger compartment and defroster nozzle" . Though they do not discuss the ideal number of layers to utilize for the windshield, Wolfahrt and Patidar recommend three to five layers and six layers respectively, suggesting a range of three to six layers of prism elements is adequate to represent the windshield properly. Patidar also recommended using ten prism layers "with 0.05 mm thickness" for ice layers on all relevant glass surfaces . It should be noted that Patidar most likely preferred the shell modeling method as it took "40% less CPU time than [the] Solid Modeling Method for computing the same de-icing simulation" . The selection of whether to use the solid modeling method or the shell modeling method largely depends on the availability of computational resources. If resources are at a premium, the shell modeling method would likely be preferred. If computational resources are more than adequate, the computational "cost" associated with the solid modeling method may be worth the additional accuracy. After preparing the model, the last major step was setting up the simulation's boundary conditions and solver parameters.
Steady State (Cold Flow) Simulation Case Setup
Based on research, Ansys Fluent was the most popular flow solver to use for deicing simulations (only Wolfahrt uses another flow solver called SWIFT). In these simulations, the flow is considered to be "steady, incompressible, viscous, Newtonian and isotropic" . To solve for turbulence, the RNG K-epsilon turbulence model was the most popular. Farag and Huang commented that this turbulence model (along with non-equilibrium wall functions) allow for the most robust results during the validation process at that time . Aroussi et al. agree with this selection of turbulence model due to its "robustness and economy of computation", which would reduce computation overhead, a common issue for these simulations during their time of development and validation . As for solver settings, a "SIMPLE algorithm is considered for the pressure-velocity coupling. Green- Gauss Cell Based steady state solver" was proposed by Patidar . While Patidar utilized a first order upwind convection differencing scheme, AbdulNour utilized a second order discretization scheme to solve the flow equations. Farag and Haung discuss that 3D deicing simulations are split into two parts: the steady state and transient simulations. The steady state simulation would optimize the cold flow directional aspect of the airflow in the vehicle cabin while the transient simulation would model the ice melting process on the windshield . For the steady state simulation, only turbulence and continuity equations are solved. The goal of the simulation is to determine the steady state isothermal flow of air going through the defrost duct into the cabin . Once the steady state simulation is completed, the heat transfer solutions can be applied to model the transient deicing process.
Transient Deicing Simulation Case Setup (Method 1)
One proposed method for determining the deicing pattern involves using a Reynolds Analogy. Hill and Sringari state that similarity between the momentum transfer and heat transfer in the windshield's boundary layer (inside the cabin, akin to a flat plate model) can be used to draw a conclusion about the temperature profile along the windshield surface . In this manner, phase-change and energy equations do not require calculation as no accuracy decrease can be perceived in the deicing pattern. "It has been shown that the ice melting patterns closely follow the heat transfer coefficients on the inside of the windshield" . This was conducted as a means to expedite simulation times to increase design iteration assessment throughput.
Transient Deicing Simulation Case Setup (Method 2)
The second method for simulating the deicing performance involves calculating the phase change behavior of the ice on the windshield through the solid model heat transfer. The critical assumption for this transient deicing simulation is that the airflow in the cabin is steady. During this simulation, only the energy equation is solved as a thermal aspect is introduced into the solver utilizing an engine coolant warm-up curve applied to the converged steady-state solution of the flow field. It should be noted that generally, the engine coolant warm-up is the most significant aspect of a vehicle's de-icing performance, but is beyond the scope of this paper (typically outside of the HVAC design engineer's and CFD analyst's control). Because a converged flow solution is the starting point for the transient simulation, the turbulence equation is not required to be solved again (as this would be a redundant computation). This is where the transience in the simulation comes from. "Since energy equation is linear, the analysis can be started as steady state to solve for the flow and turbulence equations and then solve the energy equation as a function of time on the converged flow field" . A common Fluent flow solver setup for the phase change physics is outlined by Farag and Huang. An enthalpy-porosity formulation as opposed to "tracking the liquid-solid front explicitly" is used . "The liquid- solid mushy zone is treated as a porous zone with porosity equal to the liquid fraction, and appropriate momentum sink terms are added to the momentum equations to account for the pressure drop caused by the presence of solid material" . Figure 3 shows a liquid fraction contour plot on the windshield after 25 minutes of transient deicing.
Because deicing testing and simulation is more qualitative in nature, using a liquid fraction to anchor the results mathematically makes the most sense. It should be noted, however, that what liquid fraction is considered to be free of ice will vary based on the individual engineer. With a multitude of principles for developing deicing simulations presented, the current deicing project for a commercial vehicle can be used as a case study. The application of the aforementioned ideas into an integrated design process should help future design analysis.
Extending Historical Deicing Simulation Practices to Present Day Simulations
As previously discussed, many authors developed unique procedures or optimization techniques to address certain hardware limitations or flow physics. However, the overall process, taking CAD models, meshing the internal volumes, running a steady state (cold flow) simulation, and transient deicing simulation is the same. These steps will act as the procedure for the project conducted and discussed in the sections below. However, the opportunity was taken to incorporate techniques addressed above that would provide greater simulation robustness (such as including the full HVAC module CAD geometry and HVAC module velocity profile). At the same time, perceived irrelevant techniques were removed from the simulation, such as limiting volume element size or scheme (i.e., using hexagonal mesh elements over tetrahedral elements due to computational limitations at the time). Due to the significant computational resources and robust computational software currently available, mesh refinement could be pursued without limitations commonly encountered by the authors referenced in the previous section. Following the same cadence as the historical deicing simulation discussion outlined above, the current project discussion will begin with preparing the geometry for the simulations.
Project CAD Preparation
The domain to be simulated is a semi-tractor cabin. To clarify terminology, the defrost duct is the HVAC duct that will direct airflow from the HVAC module to the front windshield, while the demist duct is the HVAC duct that directs airflow to the side windows. In this study, two defrost/demist duct pairs will be assessed. The preprocessing software ANSA was utilized to prepare the raw CAD for meshing and apply the required batch meshing scenarios to the model. The parasolid CAD file of the cabin was watertight and did not require any sealing from a surface perspective. Only a day cab variant cabin was to be simulated, so a parasolid of the curtain used to separate the forward section of the cabin from the sleeping area was merged and stitched into the cabin, creating two separate fluid domains. The sleeper segment of the cabin was removed to lighten the overall mesh count of the model similar to the technique proposed by others .
Per request by the design engineer, the updated IP was merged into the cabin model as well, which included updated demist duct outlet panels. Time was taken to properly smooth and crop areas of the cabin that would require too high a mesh density to resolve, yet would not contribute to the accuracy of the results (such as corners where the windshield and IP met). CAD models that were also merged into the model included the steering wheel, pedals, and seats in an effort to make the model perform as realistically as possible and not artificially disturb the airflow. No mannequins of occupants were utilized in the model to replicate conditions of the physical test. The remaining CAD to be merged into the model were the demist ducts, defrost ducts, and the HVAC module. Only the HVAC module from the heater core outlet to the HVAC duct mating point was utilized. As done by Pak and Simon, a velocity profile was utilized, so no upstream HVAC module CAD was required . In other words, the air flow domain starts at the outlet of the HVAC module heater core. This will be discussed in greater detail during the model setup section of this report. Data planes were added at defrost and demist duct inlets as well as outlets so as to capture relevant air flow data for post-processing. A cabin "outlet" plane was created at the base of the rear wall of the cabin to allow for a pressure outlet. This is accurate, as some air would be escaping the cabin through a vent.
Project Mesh Preparation
Since computational resources were not of major concern for this simulation, a higher mesh density was utilized throughout the model generating approximately 8 million volume elements. Because of ANSA's robust batch meshing tools, a variety of surface mesh densities were applied throughout the model. For instance, a 5-8 millimeter mesh density was utilized for the primary flow domain through the HVAC module and ductwork. Though Farag and Huang recommended utilizing prism layers in the ductwork, no prism layers were utilized in this simulation. The lack of prism layers was compensated with a finer mesh density. This would also provide a means to study if prism layers were required when compared to physical test results. Within the cabin, several surface meshing "scenarios" (as termed in ANSA) were used to vary the refinement as required. Near the windshield, a mesh density of 3-5 millimeters was utilized. ANSA's built-in mesh refinement tool further discretized the mesh in areas that would be susceptible to mesh proximity errors in the volume mesh. Figure 4 displays an example of the mesh refinement between the defrost duct outlet and the IP grille.
While refining this mesh area, one will notice that ANSA refined not only the grille surface elements, but also the defrost duct outlet planes immediately beneath them. The mesh density gradually coarsened to the original mesh criteria that was applied to these surfaces.
Meshing Schemes and Mesh Quality Criteria
Using what was learned in , the mesh on the cabin was coarsened in the foot wells as well as further away from the defrost and demist duct outlets. This would assist in limiting mesh size and thus reduce computation time for areas in the cabin that were non-value adders. Figure 5 displays an image of the surface mesh utilized for the current project.
Special caution was taken to coarsen the mesh gradually and to limit large mesh growth transitions so as to not destabilize the simulation. Once general surface meshing was created, an element skewness check was performed utilizing ANSA's built-in mesh checks. A maximum skewness allowance of 0.85 was imposed as well as a maximum aspect ratio of 15 (Fluent calculations). These values are in-line with the recommendations provided by Ansys (volume mesh skewness not greater than 0.98) and help eliminate common problems involved with simulation stability created by poor quality meshes . ANSA's automatic quality fixes and element restructuring addressed areas that were in violation of the stipulated quality criteria. These same criteria were imposed for the solid (volume) mesh as well. For the solid mesh, ANSA's Tetra Rapid meshing scheme was utilized for the air flow domain. The same fixes and checks were conducted on the volume mesh, yielding a high-quality volume mesh. Figure 6 shows a cutaway view of the volume mesh through the defrost duct and cabin air domain.
As discussed previously, the surface mesh discretization tool assisted in preparing for the 1 mm gap between the IP grille and the defrost duct outlet. Quality criteria was maintained, even with such a small gap between geometry features. Using the Solid Modeling Method discussed in Patidar , prism layers were created to represent the windshield, side windows, and ice layers. Material thicknesses and properties were provided for the glass sections by the design engineer. Prism layers of a constant height were generated for the windshield and side windows to capture each of the different material layers present in the glass. After these prism layers were created, ice layers were applied to each of the aforementioned surfaces. Multiple prism layers were generated to adequately model the thickness of the ice layer created during the physical testing procedure on each of the glass surfaces. It should be noted that the prism layers generated violated the aspect ratio mesh criteria, which was acceptable and did not destabilize the simulation. The overall volume mesh count of the entire model was approximately 8 million elements, nearly 4 times the size of models discussed earlier in the historical review of deicing simulations. Once the mesh was prepared and all quality checks were completed, a mesh file was produced by ANSA to be utilized in Fluent.
Project Steady State (Cold Flow) Simulation
At the time of the simulation, Fluent version 14.5 was used for all deicing simulations. The first case setup to be created was for the steady-state "cold flow" simulation. This would establish the steady- state airflow in the cabin that the engine warm up curve would be applied to in the transient simulation. As there has been significant improvements to CFD turbulence models since the time Farag and Huang conducted their simulations, a realizable K-epsilon turbulence model was utilized. It is a fairly robust solver that covers a wide variety of flows typically encountered within the industrial setting. Additionally, standard wall functions were utilized as opposed to non-equilibrium wall functions. Following suit with the methods proposed by Farag and Huang, the SIMPLE algorithm was used for pressure-velocity coupling and the Green-Gauss Node Based solver were utilized. These settings were also recommended in ANSYS Fluent's tutorial guide for deicing simulations. A new material for the ice was created (called "liquid-solid" from the ANSYS Fluent's tutorial) so as to be able to track the phase change of the ice on the windshield in the transient simulation . Additionally, ANSYS recommended utilizing first order discretization schemes. The only two major boundary conditions utilized for the cold flow simulation were a mass flow inlet (flow from the HVAC module's heat core outlet into the HVAC ductwork) and a pressure outlet to purge the air from the cabin. At the mass flow inlet, a velocity profile was applied to accurately represent the airflow into the cabin as opposed to a uniform mass flow rate condition similar to the method used by  (Figure 7).
It should be noted that this velocity profile was scaled to the design engineer's desired flow rate. Upon completion of the setup, the simulation was run on a Linux cluster using 48 cores (8 Xeon X5690 processors, 96 GB RAM), delivering results within approximately three hours.
Project Steady State (Cold Flow) Results
At the completion of the steady state simulations, pressure drop across the HVAC system to the cabin and flow distribution data was compared. It should be noted that for these simulations, the defrost duct outlet was divided into 8 sub-sections. A schematic is provided below for easy reference (Figure 8).
Outlets 1-4 are on the driver side while Outlets 5-8 are on the passenger side (numerically ordered from the driver's seat, left to right). The divisions that were created were strictly for data assessment purposes. The first set of data was taken for the defrost duct pressure drop, demist duct pressure drop, and then the whole duct (defrost and demist) pressure drop. The defrost duct total pressure drop data can be seen below in Figure 9.
Based on the steady-state CFD results, Iteration 2's design performed marginally better in the way of pressure drop performance. Updated defrost duct design (including vane angle adjustments) and the updated HVAC vane geometry found in Iteration 2 allowed for the reduction in pressure drop over Iteration 1's defrost duct design by approximately 3%. As the pressure drop is smaller, the velocity component of total pressure theoretically should be larger to satisfy the Bernoulli Equation. To assess this hypothesis, airflow out of the defrost duct was compared between Iteration 1 and 2 in Figure 10.
When comparing the total airflow across each iteration's defrost duct, it was determined that Iteration 1 flowed approximately 3.9% more airflow than Iteration 2. This is counterintuitive based on the Bernoulli Equation as discussed above, but this may stem from the uncertainty in the simulation (such as mesh size in the individual duct outlets). Additionally, the revised HVAC module and demist duct design provided by the design engineer generated a larger total pressure drop across Iteration 2's demist duct. Figure 11 shows that a 30% total pressure drop increase in the demist duct exists due to these changes.
Upon inspecting the HVAC module and demist duct more closely, it was determined that the HVAC modules exhibit areas of recirculation. The recirculation zone found in Iteration 1's driver side demist duct plenum was mitigated with Iteration 2's revised HVAC module design, allowing for less restriction of the airflow through the demist duct. Figure 12 shows the change in recirculation zone by reduction in captured volume.
In Iteration 1, one can see two distinct recirculation zones. As the flow enters the demist duct plenum it is already encountering the large recirculation zone. Additionally, most of the flow is pointed directly at the recirculation flow, causing significant flow impingement. With the revisions made in Iteration 2, the air flow enters parallel to the recirculation zone, which allows the flow to experience a "slingshot" effect as it is accelerated around the recirculation zone (vortex), as noted by the increased air velocity directly above the recirculation zone. One can see that the flow is indeed accelerated around this recirculation zone, and assists in driving down the recirculation penalty. With the revised air flow direction, one can see that the passenger-side recirculation zone is no longer an issue. By capitalizing on the vortex effect experienced in the passenger-side demist duct plenum, the flow headed to the driver- side duct plenum is "pulled" upwards to make it more aligned with the flow passage through the driver-side demist duct plenum. However, these benefits are not significant enough to overcome the source of the demist duct total pressure drop, shown in Figure 13.
The cross sectional area change found immediately upon exiting the demist duct plenum is the probable cause of the 30% increase in demist duct total pressure drop. The flow is constricted as the cross-section gradually returns to a circular shape further in the demist duct. As an additional penalty, because of the odd shape of Iteration 2's cross-sectional area (red area in Figure 13 above), the flow is not as easily "guided" through the duct, and takes time to reattach itself to the walls of the duct and pick up lost velocity. Though the defrost duct geometry is fairly well-optimized, clear benefits can be seen from utilizing Iteration 1's demist duct design over Iteration 2's design purely from a pressure drop perspective. Additional views of the airflow through the cabin as a whole were requested by the design engineer for added reference material. An example of this data is shown in Figure 14 below.
With these views at his or her disposal, a design engineer is better able to determine how vane angles in the defrost duct outlets impact steady state airflow as a whole. If a portion of the windshield or side windows are not receiving enough airflow across their surfaces, the design engineer can update his or her design prior to proceeding to transient deicing simulations. With the steady-state post-processing complete, the transient simulation setup becomes the next topic of discussion.
Project Transient Deicing Case Setup
The converged solution of flow internal to the cabin was utilized as the "time zero" case for the transient simulation. Farag and Huang's case setup was utilized. An engine coolant warm-up curve from a 15 liter diesel engine during a physical test was used to apply the thermal behavior of the air flowing through the cabin at the maximum volumetric flow rate dictated by the HVAC module hardware. As mentioned by Patidar, only the energy equation is solved in this simulation as the cold flow field is taken to be steady state . This assists in expediting the transient simulation as a large amount of computation resources are no longer required to manage both the turbulence and energy equation calculations. As done by Farag and Huang, the enthalpy-porosity formulation (utilization of the mushy zone formulation under the solidification and melting model in Fluent) is used for the transient simulation, thus using the porosity of the ice layer (and therefore the liquid fraction) as a measure of how much the windshield and side windows are free from ice. Which liquid fraction value to consider ice melted is up to the individual analyst. In these studies, a liquid fraction of 0.675 is considered free of ice (i.e., melted). This value was based on prior analysis that was calibrated to a physical test. Time steps of 4 seconds were used so as to capture deicing data at all SAE-specified time intervals (for a total of 750 time steps). Each time step utilized 50 iterations for convergence, where convergence was defined as a reduction in energy residual from the start of each time step to the end by at least two orders of magnitude. Because of the refined time step duration, each time step converged within 10 iterations or so. In other words, each "steady state step" in the transient simulation did not have a drastic change to the initial condition value, so the simulation would converge easily as opposed to having time steps that were significantly larger. During the simulation, executable commands were utilized as per the Fluent user manual to output liquid fraction contours of all the necessary glass surfaces at the end of each time step, allowing for an animation of the deicing process to be generated . The two transient simulations were conducted on a remote machine utilizing 12 cores (2 Xeon X5675 processors at 3.06 Ghz), 128 GB RAM, and took approximately 9 hours to complete each transient simulation.
Project Transient Deicing Results
An example of the windshield deiced pattern is shown in Figure 15 (Iteration 1) with a data table showing percent free of ice in particular target areas. Figure 15 is split into two main pieces of information, the quantitative time vs percent of windshield cleared data, and the qualitative deicing pattern images. The tabular data includes the SAE zone target percent cleared values ("Target A" for example) and its associated simulation result (Column "A"). This allows for a quick comparison between design iterations to see immediately what role duct geometry updates have on the windshield deicing performance. One will notice that the Iteration 1 geometry meets the required SAE targets, exceeding target at 30 minutes for both Target A and Target C. In the lower right half of the figure, the simulated windshield is shown. It should be noted that the "cracks" (black lines) in the windshield are only due to the simulation being parallelized on a Linux cluster (i.e., parts of the model are segmented and computed on separate nodes), and should not be considered as part of the deicing pattern data. Each of the red lines border the area that is free from ice (in this case, having a liquid fraction greater than 0.675). The time interval (minutes) is displayed in the blue bubbles and point to the associated deicing contour outline. The design engineer also requested that the breakthrough points (the first location(s) where the windshield is free of ice) be marked on the windshield image, which was estimated based on the deicing animations generated after the transient simulation was completed. Similar images were generated for the side windows as well. It should be noted that an additional sub-criteria was generated for the side windows not called out by the SAE J381 Test Documentation at the request of the design engineer (these areas are called the short zone, tall zone, shared zone, and mirror reflect zones in the images below).
As done in Figure 15 above, the design engineer's additional mirror zones were overlaid onto the window surface as shown in Figure 16. For side windows, only one SAE-mandated criteria exists, the total percentage of the window free of ice at 45 minutes. For the driver side window, this is 65% and for the passenger side window, the value is 70%. As can be seen in the tabular data and window image above, no deicing that met the required liquid fraction criteria occurred during the first 20 minutes. When comparing the contours, there appears to be no design that is far superior to the other. However, when a deep-dive study is conducted into the tabular data, greater insight into the benefits of each HVAC system's design was determined. Zone C (driver sightline) will be studied as it is the most critical area of deicing and has the most stringent SAE criteria. A plot of deicing time (minutes) vs Zone C percent cleared (%) was generated along with a curve for both Iteration 1 and Iteration 2 (Figure 17).
Minutes A Target A B C Target C Whole 15 0% 0% 0% 0% 20 0% 0% 0% 0.3% 25 36 1% 51 5% 81.8% 30.3% 30 93 9% 85% 100% 100% 99% 62.3% 35 100% 100% 100% 72.8% 40 100% 100% 100% 75.7% 45 100% 95% 100% 100% 77 5% 50 100% 100% 100% 78.8%
Though Iteration 1 has a higher pressure drop, the defrost duct (total flow) still flows 3.9% more air than Iteration 2 does. This idea is further supported by the defrost duct individual outlet flow rate comparison shown in Figure 18.
At an individual defrost duct outlet level, the Iteration 1 driver side outlets (Sections 1-4) flow an additional 8% air when compared to Iteration 2's flow. Therefore, Iteration 1's defrost duct nets a gain in deicing speed that is 1.6 times greater than Iteration 2's when compared from 20 to 25 minutes. However, Iteration 1's deicing speed gain is minimized from 25 to 30 minutes, which may suggest that Iteration 2's reduced duct pressure drop (i.e., system flow restriction) begins to surpass the flow rate benefit held by Iteration 1. From 25-30 minutes, Iteration 1's deicing speed slope begins to plateau and meets Iteration 2's deicing speed slope at the 30 minute mark (a slope which remained steady from 20 minutes to 30 minutes).
Minutes Main Tall Short Shared Target Total Target Window Zone Zone Zone Mirror Window Total Reflect Window Zones 15 0% 0% 0% 0% 0% 20 0% 0% 0% 0% 0% 25 25.3% 75.5% 97.8% 100% 30.3% 30 64.2% 100% 100% 100% 70% 67.1% 35 95.7% 100% 100% 100% 96% 40 100% 100% 100% 100% 100% 45 100% 100% 100% 100% 99% 100% 65% 50 100% 100% 100% 100% 100%
Summary and Conclusions
In conclusion, CFD offers a unique solution to the common problem of ice accumulation on automotive and commercial vehicle windshields. Through years of study, the simulating of deicing behavior on windshields and other glass surfaces (such as side windows) has been honed to a repeatable process which can be utilized in part or in whole while designing the HVAC system. While conducting their individual simulations, several innovators introduced their own "style" to the process, such as reducing CAD generation time, optimizing volume and surface mesh schemes, and utilizing the most desirable flow physics in cold flow and transient modelling. Frontloading simulation in the design process gives design engineers the opportunity to benefit from not only large amounts of quantitative data, but also qualitative data generated in the computational domain. Assuming a design engineer does not have ready access to a physical test chamber, or requires several design iterations prior to conducting physical testing, modeling the HVAC system and cabin allow relevant design criteria, such as pressure drop, flow rates, and ultimately transient deicing behavior to be gathered within days at minimal cost. This time benefit has been further compounded with the increased availability of high speed computing and parallelization found in commercial CFD codes available to the industry. By reflecting on the approaches used in deicing simulation over the years, an in-depth review of each simulation process was conducted on a current project which combined the best simulation techniques. By leveraging CFD deicing simulations, a design engineer has more time to iterate on a design, leading to an even more robust product. When starting deicing simulations from scratch, the largest initial time cost is preparing the CAD and mesh for the cabin and duct work (flow domain).
Depending on the experience level of the analyst and complexity of the geometry, the process can take anywhere from several days to three weeks. Going forward, however, the time penalty is small for iterative designs and updates. For the simulations outlined above, the steady state duct flow simulations took roughly three hours per run (48 cores, 8 Xeon X5690 processor at 3.47 Ghz, 96 GB RAM), and the transient deicing simulations took approximately nine hours per run (12 cores, 2 Xeon X5675 processors at 3.06 Ghz, 128 GB RAM). Depending on the design engineer's data requirement, he or she can have steady state and transient CFD simulation data within a few hours or two business days respectively. Iterative updates associated with the cabin geometry, ductwork, or ice thickness would require meshing updates netting approximately a day or so (dependent upon geometry complexity) followed by the simulation time outlined above. Boundary condition updates, such as a new velocity profile, would only require an update to the case file boundary conditions followed by the simulation time noted above.
Because the project was not an exhaustive investigation of all the historical simulation techniques, significant opportunity exists to further investigate simulation procedure and model setup. For instance, prism layers beyond those used on the glass/ice surfaces were not used at the time this project was conducted. Having access to physical test data would allow for follow-up work to determine if there is a clear benefit to using prism layers or eliminating them entirely from the HVAC ductwork mesh.
Furthermore, investigation into using a Cartesian or even hexagonal volume mesh as an additional means of model size control could prove to be beneficial, especially when a cluster or parallelization is unavailable. Additionally, a deeper investigation into "grey ice" (haze on the windshield that may or may not be "ice" based on the passenger's perspective) and its relationship to what liquid fraction should be considered as being "melted" should be further investigated as test data is made available. Lastly, a deeper dive into the heat transfer impact compared to duct pressure drop impact on deicing speed is also worth further investigation. As the need for appropriate windshield deicing measures continues to move from the physical to the computational domain, CFD simulation will continue to play a pivotal role in automotive and commercial vehicle HVAC system design.
Thanks to Navistar, Inc., Ashraf Farag, and Scott Paciero for advocating the importance of this work and its publication.
[1.] U.S. Department of Transportation: Federal Highway Administration, "How Do Weather Events Impact Roads?," http://www.ops.fhwa.dot.gov/weather/q1 roadimpact.htm, accessed May 15, 2015.
[2.] Aroussi, A., Hassan, A., Clayton, B., AbdulNour, B. et al., "Improving Vehicle Windshield Defrosting and Demisting," SAE Technical Paper 2000-01-1278, 2000, doi:10.4271/2000-01-1278.
[3.] Sneed, A., "Moore's Law Keeps Going, Defying Expectations," Scientific American, May 19, 2015, http://www.scientificamerican.com/article/moore-s-law-keeps-going-defying-expectations/.
[4.] SAE International Surface Vehicle Recommended Practice, "Windshield Defrosting Systems Test Procedure and Performance Requirements-Trucks, Buses, and Multipurpose Vehicles," SAE Standard J381, Rev. Jan 2009.
[5.] Farag, A. and Huang, L., "CFD Analysis and Validation of Automotive Windshield De-Icing Simulation," SAE Technical Paper 2003-01-1079, 2003, doi:10.4271/2003-01-1079.
[6.] Wolfahrt, J., Baier, W., Wiesler, B., and Moshammer, T., "Numerical Studies for De-Icing Validation," SAE Technical Paper 2005-01-1883, 2005, doi:10.4271/2005-01-1883.
[7.] Pak, J. and Simon, L., "Cabin Flow and the Windshield De-Icing Simulation with the Parametric Cabin Modeler," SAE Technical Paper 2003-01-1080, 2003, doi:10.4271/2003-01-1080.
[8.] AbdulNour, B., "CFD Prediction of Automotive Windshield Defrost Pattern," SAE Technical Paper 1999-01-1203, 1999, doi:10.4271/1999-01-1203.
[9.] Hill, W. and Sringari, S., "A Parametric Approach for Rapid Design and Analysis of Automotive HVAC Defrost Systems," SAE Technical Paper 2001-01-0584, 2001, doi:10.4271/2001-01-0584.
[10.] Patidar, A., "Comparison of CFD Analysis Methods for Simulating De-Icing Pattern over Automotive Windshield and Windows," SAE Technical Paper 2010-01-0555, 2010, doi:10.4271/2010-01-0555.
[11.] Patidar, A., Natarajan, S., and Pande, M., "CFD Analysis and Validation of an Automotive HVAC System," SAE Technical Paper 2009-01-0535, 2009, doi:10.4271/2009-01-0535.
[12.] Ansys, Inc., "Introduction to ANSYS Fluent Meshing-Module 7: Volume Fill Methods (Version 16.1 Release)," Sept. 22, 2015.
[13.] Ansys, Inc., "Tutorial: Windshield De-Icing Analysis," Dec. 3, 2010.
Krystian J. Link, Navistar, Inc.
Nicholas A. Pohlman, Northern Illinois University
Received: 25 Jan 2018
e-Available: 06 Apr 2018
|Printer friendly Cite/link Email Feedback|
|Title Annotation:||computational fluid dynamic|
|Author:||Link, Krystian J.; Pohlman, Nicholas A.|
|Publication:||SAE International Journal of Commercial Vehicles|
|Date:||Mar 1, 2018|
|Previous Article:||Investigation on Underhood Thermal Analysis of Truck Platooning.|
|Next Article:||Aerodynamic Analysis of Cooling Airflow for Different Front-End Designs of a Heavy-Duty Cab-Over-Engine Truck.|