Topological optimization method for heat flow problem of cooling runner of cooling element of high-speed precision machine tool

Abstract

The invention discloses a topological optimization method for a heat flow problem of a cooling channel of a cooling element of a high-speed precision machine tool based on self-adaptive mesh refinement and parallel calculation. The method aims at maximizing the heat transfer capacity in the cooling channel and ensuring that the minimum flow requirement is met at the same time. In the process, the pressure drop is used as a constraint function to facilitate adjustment of flow characteristics, and it is ensured that fluid dynamics cannot excessively hinder cooling efficiency; the effectiveness of the design channel is quantitatively evaluated through the average temperature in the design domain; the lower average temperature represents the higher heat dissipation capability; this approach balances the dual requirements of thermal control and hydrodynamics to design efficient and effective cooling channels. According to the method, the finite volume method is combined with adaptive mesh refinement and parallel calculation, fluid-solid boundary definition is sharpened, the calculation process is accelerated, the calculation efficiency and precision are improved, and the cost related to thermal fluid topological optimization is reduced, so that an effective solution is provided for design of cooling elements.

Description

Technical Field

[0001] The present invention belongs to the technical field of heat flow problems, and specifically relates to a topology optimization method for the heat flow problem of the cooling flow path of a cooling element of a high-speed precision machine tool. Background Art

[0002] There are many heat sources with significant heat load intensity in precision boring machine tools. Due to the increase in temperature, obvious thermal errors are caused, as Figure 1 shown. The thermal error of a precision boring machine tool is the main cause of machining errors. Therefore, when using a boring machine for high-precision machining, effective temperature control measures must be taken. The cooling system and ambient temperature control are used to reduce the influence of thermal errors on machining accuracy, and traditional cooling elements with serpentine and spiral cooling channels are designed. The coolant is driven by an industrial cooler, which effectively dissipates the internal heat and reduces the thermal error. Therefore, designing a cooling element with strong heat transfer ability and low pressure drop loss is crucial for reducing temperature rise and thermal error and improving the machining accuracy of a precision boring machine. In addition, cooling capacity and efficiency are the key criteria for evaluating the performance of a cooling element, and the balance between these two criteria is crucial for efficient heat dissipation and thermal error control. Therefore, it is crucial to determine the optimal balance between heat dissipation efficiency and pressure drop when designing a cooling element.

[0003] Design strategies such as shape optimization, optimization based on fractal theory, and topology optimization (TO) are often used to achieve this balance, aiming to optimize the heat transfer performance and fluid dynamics in the cooling channel. Shape optimization and optimization based on fractal theory rely heavily on the designer's empirical knowledge, which often makes it complicated to achieve the optimal solution in the design of cooling elements. For example, common designs such as spiral and serpentine cooling channels, although widely used, have disadvantages such as uneven cooling temperature and large pressure drop loss, making them suboptimal components for efficient cooling. These designs often strongly rely on experience or bionic strategies.

[0004] With the progress of computing power, TO has become a revolutionary method in radiator design. TO is particularly effective in applications such as electronic chips and battery energy management. The heat transfer ability and efficiency of TO channels are far stronger than those of traditional serpentine and spiral cooling channels. Due to its wide design freedom, TO also greatly shortens the design cycle. In addition, as additive manufacturing develops towards more cost-effective and scalable production, the manufacturing of complex structures obtained through TO has become increasingly feasible, paving the way for innovative cooling channel designs. Existing research has emphasized the concentrated efforts to solve the laminar flow problem in thermal fluid TO. However, high heat flux radiators often encounter turbulence in the cooling medium, which is a major challenge. TO designs optimized for laminar flow conditions cannot be seamlessly transferred to applications involving turbulence in cooling channel designs, indicating the need for tailored methods to handle these more complex flow situations.

[0005] Turbulent TO is still in the initial stage of development, and the complexity of the turbulent problem leads to huge computational requirements due to its high nonlinearity. The Darcy flow model has been proven effective in reducing computational time and achieving stable numerical convergence in the challenges of thermal-fluid TO. However, most existing thermal-fluid TO methods are proposed based on the finite element method (FEM), which may not be as robust and stable as the finite volume method (FVM) when dealing with computational fluid dynamics (CFD) problems.

[0006] Density-based methods have been widely used in TO, but they face challenges such as unclear boundaries between solid and fluid domains. Coarse grids often produce blurred boundaries, and only fine grids can ensure the smoothness of the boundaries, thus generating TO results with enhanced manufacturability. However, grid refinement greatly increases the computational load. Accurately depicting the fluid-solid boundary is crucial for effectively capturing the heat transfer effect in the design of thermal-fluid systems. Therefore, adopting an appropriate grid resolution is crucial for ensuring an accurate numerical solution of heat transfer.

[0007] As TO develops towards solving large-scale problems and benefits from improvements in computer hardware and software performance, parallel computing technology has become increasingly important in this field. The MPI function in the Pstream library of OpenFOAM is used to implement parallel computing, thus improving the computational efficiency of TO applications. This strategic integration of parallel computing frameworks promotes a more effective and efficient optimization process, which is crucial for solving the increasingly complex and large-scale problems of TO challenges. Summary of the Invention

[0008] In view of this, the purpose of the present invention is to provide a topology optimization method for the heat flow problem of the cooling channel of a cooling element of a high-speed precision machine tool, so that the designed cooling element has stronger heat transfer ability and lower pressure drop loss.

[0009] To achieve the above purpose, the present invention provides the following technical solutions:

[0010] A topology optimization method for the heat flow problem of the cooling channel of a cooling element of a high-speed precision machine tool, comprising the following steps:

[0011] Step 1: Construct the objective function and optimization conditions for the topology optimization of the heat flow problem of the cooling channel to maximize the heat transfer capacity in the cooling channel while ensuring that the minimum flow rate requirement is met:

[0012] Find{γ1,γ2,...,γ i} T

[0013]

[0014] Where: Ψ is the objective function; J1 is the volume constraint of the fluid; J2 is the pressure drop constraint of the cooling channel; V is the volume of the fluid region in the design domain; Ω is the design domain; R is the control equation; P in is the inlet pressure; ∫ inlet *dΓ is the integral of the computational domain on the inlet boundary; u is the velocity vector; p is the pressure; T is the temperature; γ is the design density; γ i is the i-th objective to be optimized; d is the design domain variable;

[0015] Step 2: Solve the objective function

[0016] 21) Initialization: Decompose the computational domain for parallel computing and initialize the design domain Ω;

[0017] 22) Determine whether the absolute value of the difference between the objective function Ψ n at the current iteration step and the objective function Ψ n―1 at the previous iteration step is greater than the preset threshold, or, determine whether the current iteration number n is less than or equal to the preset maximum iteration number N: If so, execute step 23); If not, the design optimization is completed;

[0018] 23) Calculate the state variables by solving the original problem;

[0019] 24) Calculate the adjoint variables by solving the adjoint problem;

[0020] 25) Calculate the objective function and sensitivity; If the current projection slope β is less than the preset maximum projection slope β Max , then update the projection slope β;

[0021] 26) Update the design variables by the MMA algorithm; If the following conditions are met:

[0022]

[0023] Where: ε Ψ is the objective function error threshold; ε V is the volume constraint error threshold; ε p is the pressure drop error threshold; ∫ Ω *dΩ is the integral in the design domain; V frac is the flow channel ratio; p n is the pressure at the current iteration step; p n―1 is the pressure value at the previous iteration step;

[0024] Then, adopt the adaptive mesh refinement method to update the mesh;

[0025] 27) If the permeability penalty parameter p κ is less than or equal to the preset maximum permeability penalty parameter p κMax , then update the permeability penalty parameter p κ;

[0026] 28) Calculate the material permeability κ and thermal diffusivity D T ;

[0027] 29) Let n = n + 1, and loop to execute step 22).

[0028] Furthermore, in the said step 23), the state variables include the velocity vector u, pressure p, and temperature T; the method for calculating the state variables by solving the original problem is as follows:

[0029] Construct the continuity equation for a steady incompressible flow with mass conservation:

[0030]

[0031] where: u is the velocity vector; represents the gradient operator;

[0032] Adopt the Darcy model for turbulence modeling:

[0033]

[0034] where: p is the pressure, and κ and μ represent the fluid permeability and dynamic viscosity respectively;

[0035] For an incompressible flow, when considering the transport phenomena involving temperature, the description can be extended to include the transfer and diffusion of temperature; the equations for controlling these processes in an incompressible fluid are:

[0036]

[0037] where: T is the temperature; D t is the thermal diffusivity;

[0038] The heat source q(X) depends on the design domain x:

[0039]

[0040] where: Ω is the design domain; θ is the heat source domain; x ∈ Ω\θ is the design domain excluding other regions except the heat source.

[0041] Furthermore, in the said step 24), the adjoint variables include the adjoint velocity u A 、the adjoint pressure p A and the adjoint temperature T A ; the method for calculating the adjoint variables by solving the adjoint problem is as follows:

[0042] Combine the objective function and the constraint conditions into a Lagrangian function:

[0043] L(u,u A ,p,p A,T,T A ,λ1,λ2,Ω)

[0044] =Ψ+∫ Ω u A ·R u dΩ+∫ Ω p A R p dΩ+∫ Ω T A R T dΩ+λ1J1+λ2J2

[0045] where: u A , p A , T A , λ1 and λ2 represent the adjoint velocity, adjoint pressure, adjoint temperature, the Lagrange multiplier related to the volume constraint, and the Lagrange multiplier related to the pressure drop constraint, respectively; u A , p A , T A , λ1 and λ2 are all Lagrange multipliers, introduced as dummy variables;

[0046] The adjoint model is obtained as:

[0047]

[0048] Since the governing equations and the constraint conditions of J1 and J2 are theoretically constant, so L = Ψ; then the sensitivity of the objective function is The following relationship is obtained:

[0049]

[0050] Find the appropriate adjoint variables u A , p A , T A , λ1 and λ2 such that Then the final Lagrangian equation is constructed as follows:

[0051]

[0052] where: is the governing equation;

[0053] To obtain the optimal topology optimization result, the following conditions should be satisfied:

[0054]

[0055] where: δ ζ Ψ, δ ζ J and are the objective function Ψ, the constraint conditions, respectivelyJ and the differential of the control equation i with respect to the variable ξ; ξ is a variable and ξ ∈ (u, p, T);

[0056] The differential form is obtained as follows:

[0057]

[0058] where: ∫ Ω *dΩ and ∫ Γ *dΓ represent the integral in the design domain and the integral on the boundary, respectively;

[0059] For any δu, δp, and δT, the calculation requirements within the design domain are satisfied, and the adjoint equation is obtained as follows:

[0060]

[0061] Consider the adjoint boundary conditions:

[0062]

[0063] For the inlet and the wall, the velocity and temperature are fixed values, so δu = 0 and δT = 0, and the adjoint velocity and adjoint temperature at the inlet and the wall are:

[0064]

[0065] T A = 0

[0066] For the outlet, the pressure is 0 and the temperature has a zero gradient, so δp = 0 and The adjoint velocity and adjoint temperature at the outlet are:

[0067] p A = 0

[0068]

[0069] where: is the adjoint normal velocity.

[0070] Furthermore, in step 25), the sensitivity calculation method for the thermal-fluid topology optimization problem is:

[0071]

[0072] where: δ γ is the volume factor; k s and κ f represent the solid permeability and the fluid permeability, respectively; p κ is the penalty parameter for permeability; is the penalty parameter for the thermal diffusivity; and respectively represent the solid thermal diffusivity and the fluid thermal diffusivity.

[0073] Further, in the step 28), the calculation method of the permeability κ is as follows:

[0074]

[0075] where: κ is the permeability; κ s and κ f respectively represent the solid permeability and the fluid permeability; p κ is the penalty parameter of the permeability;

[0076] The calculation method of the thermal diffusivity D T is as follows:

[0077]

[0078] where: is the penalty parameter of the thermal diffusivity; and respectively represent the solid thermal diffusivity and the fluid thermal diffusivity.

[0079] The beneficial effects of the present invention are as follows:

[0080] The best cooling channel design aims to achieve a large amount of heat transfer while minimizing the flow resistance to the greatest extent. These two goals often conflict with each other. Therefore, the goal of the topology optimization method for the heat flow problem of the cooling channel of the cooling element of the high-speed precision machine tool of the present invention is to maximize the heat transfer capacity in the cooling channel while ensuring that the minimum flow requirement is met. In this process, using the pressure drop as the constraint function helps to adjust the flow characteristics and ensure that the fluid dynamics does not overly impede the cooling efficiency; the effectiveness of the designed channel is quantitatively evaluated by the average temperature in the design domain; a lower average temperature indicates a higher heat dissipation capacity; this method balances the dual requirements of thermal control and fluid dynamics to design an efficient and effective cooling channel. The present invention combines the finite volume method with adaptive grid refinement and parallel computing, sharpens the definition of the fluid-solid boundary and accelerates the calculation process, improves the calculation efficiency and accuracy, and reduces the cost associated with the thermal-fluid topology optimization, thereby providing an effective solution for the design of the cooling element. Description of the Drawings

[0081] In order to make the objectives, technical solutions, and beneficial effects of the present invention clearer, the present invention provides the following drawings for description:

[0082] Figure 1 is the thermal field of the precision boring machine tool; Figure 2 is for different cooling applications; Figure 3 is the adaptive grid refinement of 1 unit; Figure 4 is the distribution diagram of the unit processors;

[0083] Figure 5 Results comparison between the finite volume method solver and the finite element method solver; Figure 6 Cooling plate design domain diagram;

[0084] Figure 7 Iteration process of the cooling plate objective function, constraints, and number of meshes; Figure 8 Topology optimization results for different adaptive mesh refinement levels;

[0085] Figure 9 Numerical comparison of different adaptive mesh refinement levels on the cutting line; Figure 10 Temperature gradient distribution on the cross-section for different mesh refinement levels;

[0086] Figure 11 3D design results for different temperature and pressure drop combinations; Figure 12 Heat dissipation comparison between the rectangular channel and the topology optimized channel;

[0087] Figure 13 Pressure drop and average temperature of different topology optimized channels; Figure 14 Water jacket design domain; Figure 15 Application results of adaptive mesh refinement;

[0088] Figure 16 Iteration process of the cylinder objective function switching, constraint switching, volume, and number of meshes;

[0089] Figure 17 TO water jacket design results; Figure 18 Material distribution of different cross-sectional areas of TO cooling channel B-6;

[0090] Figure 19 Cooling capacity and efficiency comparison of different TO cooling channels; Figure 20 Topology optimization results calculated using two different solvers;

[0091] Figure 21 Pressure field and temperature field distributions calculated using two different solvers;

[0092] Figure 22 Freely installed cooling plate; Figure 23 Thermal behavior measurement; Figure 24 Finite element model of the boring machine;

[0093] Figure 25 Temperature, temperature gradient distribution, and thermal deformation of the boring machine without cooling;

[0094] Figure 26 Cooling plate installation; Figure 27 Boring machine with cooling plate; Figure 28is the temperature of the skateboard and the spindle system;

[0095] Figure 29 is the dynamic temperature change of key components; Figure 30 is the thermal deformation of the boring machine; Figure 31 is the thermal error of the spindle system and the center point of the workbench;

[0096] Figure 32 is the thermal error of the boring machine; Figure 33 is the topological optimization of the water jacket; Figure 34 is the mesh division result; Figure 35 is the temperature field of the spindle system of the boring machine; Figure 36 is the average temperature of the spindle of the boring machine; Figure 37 are the highest temperature and the stator temperature; Figure 38 is the coolant pressure field; Figure 39 is the pressure drop under different inlet Reynolds numbers; Figure 40 is the thermal deformation;

[0097] Figure 41 is the thermal deformation under different inlet Reynolds numbers.

[0098]

[0099] Specific implementation manners

[0100] The present invention will be further described below in conjunction with the accompanying drawings and specific embodiments, so that those skilled in the art can better understand the present invention and be able to implement it, but the embodiments given are not intended to limit the present invention.

[0101] 1. Topological optimization process of the heat flow problem of the cooling channel

[0102] In this embodiment, based on the thermal fluid TO model and numerical methods, a thermal fluid TO solver was developed using the open-source platform OpenFOAM. The main program includes four basic components: (1) dividing the computational domain for multi-processor parallel computing; (2) calculating the original state variable equation and the adjoint state variable equation to obtain the objective function and sensitivity analysis; (3) applying the moving asymptote (MMA) algorithm to refine the design variables; (4) facilitating communication and implementing the AMR strategy between each processor. A high projection slope β may cause the design variable γ to tend to a binary distribution; in order to maintain the stability of the iterative process, the projection slope β is adjusted to enhance the optimization step. The high value of p κ in the interpolation model exacerbates the nonlinearity of the TO model; in order to avoid premature convergence to a local optimal value in the initial stage, p κIt is initially set to 1 and then gradually increased to 5. The module for solving the original and adjoint state variable equations is encapsulated in a library. In addition, the MMA algorithm is used as an optimization tool and compiled into a library for integration with the main program. The termination of the iterative optimization process is triggered either by the relative change of the objective function dropping below 0.00001 or by reaching a predetermined number of iterations.

[0103] Specifically, the topology optimization method for the heat flow problem of the cooling element cooling flow path of the high-speed precision machine tool in this embodiment includes the following steps:

[0104] Step 1: To maximize the heat transfer capacity in the cooling channel while ensuring that the minimum flow rate requirement is met, construct the objective function and optimization conditions for the topology optimization of the heat flow problem of the cooling channel:

[0105] Find{γ1,γ2,...,γ i} T

[0106]

[0107] where: Ψ is the objective function; J1 is the volume constraint of the fluid; J2 is the pressure drop constraint of the cooling channel; V is the volume of the fluid region in the design domain; Ω is the design domain; R is the control equation; P in is the inlet pressure; ∫ inlet *dΓ is the integral of the computational domain on the inlet boundary; u is the velocity vector; p is the pressure; T is the temperature; γ is the design density; γ i is the i-th optimization target to be optimized; d is the design domain variable.

[0108] Step 2: Solve the objective function

[0109] 21) Initialization: Decompose the computational domain for parallel computing and initialize the design domain Ω.

[0110] 22) Determine whether the absolute value of the difference between the objective function Ψ n at the current iteration step and the objective function Ψ n―1 at the previous iteration step is greater than a preset threshold, or, determine whether the current iteration number n is less than or equal to the preset maximum iteration number N: If so, execute step 23); if not, the design optimization is completed.

[0111] 23) Calculate the state variables by solving the original problem. Specifically, the state variables include the velocity vector u, the pressure p, and the temperature T.

[0112] 24) Calculate the adjoint variables by solving the adjoint problem. Specifically, the adjoint variables include the adjoint velocity u A 、the adjoint pressure p A and the adjoint temperature T A .

[0113] 25) Calculate the objective function and sensitivity; if the current projection slope β is less than the preset maximum projection slope β Max , then update the projection slope β. In this embodiment, the method for updating the projection slope β is: control the projection slope β to gradually increase as the optimization process progresses, and the initial value of β is 1 and the maximum value is 200.

[0114] 26) Update the design variables through the MMA algorithm; if the following conditions are met:

[0115]

[0116] where: ε Ψ is the objective function error threshold; ε v is the volume constraint error threshold; ε p is the pressure drop error threshold; ∫ Ω *dΩ is the integral in the design domain; V frac is the flow channel ratio; p n is the pressure at the current iteration step; p n―1 is the pressure value at the previous iteration step;

[0117] then, adopt the adaptive mesh refinement method to update the mesh.

[0118] 27) If the permeability penalty parameter p κ is less than or equal to the preset maximum permeability penalty parameter p κMax , then update the permeability penalty parameter p κ . In this embodiment, the update method of the permeability penalty parameter p κ is: p κ is initially set to 1 and then gradually increased to 5.

[0119] 28) Calculate the material permeability κ and thermal diffusivity D T .

[0120] 29) Let n = n + 1 and loop to execute step 22).

[0121]

[0122] 2. Topological optimization of the thermo-fluid problem

[0123] Figure 2(a) shows the application of a freely mounted cooling plate in a boring machine, aiming to improve thermal stability and reduce thermal errors. The freely mounted cooling plate is connected to the boring machine slide using a magnetic base, which can effectively dissipate the heat generated during the operation of the boring machine. An industrial chiller provides a stable coolant supply to the freely mounted cooling plate, helping to remove the internal heat of the boring machine. Data acquisition sensors are strategically placed to continuously monitor the temperatures of various components of the boring machine, including the machine bed and column. Then, the real-time temperature data collected by these sensors is transmitted to a computer system. The collected temperature data is used to regulate the temperature and flow rate of the coolant supplied by the industrial chiller, controlling the temperature field and thermal errors. This integrated approach involves several key principles of temperature and thermal error control: (1) The freely mounted cooling plate effectively dissipates the heat of the boring machine, preventing overheating and ensuring stable operation; (2) The industrial chiller provides a stable coolant flow, keeping the machine components at the optimal working temperature; (3) The sensors continuously collect the temperature data of various components of the boring machine, providing crucial real-time monitoring; (4) The computer system dynamically adjusts the temperature and flow rate of the coolant based on the real-time data, ensuring precise temperature control. By combining these elements, the cooling system effectively controls heat dissipation and minimizes thermal errors, thereby reducing thermal errors, improving machining accuracy, and the operating stability of the boring machine.

[0124] Figure 2 (b) shows the application of a water jacket in the boring machine spindle system, and details the principles of temperature rise and thermal error control. The temperature control system connected to the cooler regulates and maintains the coolant temperature. The cooler provides the coolant, ensuring it remains within a preset range for effective cooling. The water jacket is installed on the high-speed spindle system of the boring machine and absorbs and dissipates the heat generated during the operation of the high-speed spindle system through circulating coolant. Temperature sensors are located at the inlet and outlet of the cooling channels to monitor the temperature changes of the cooling water when it enters and exits. Then, the temperature data is monitored and recorded in real time through a computer and monitor, enabling the operator to view the temperature changes and adjust the temperature control system settings as needed. In addition, the traditional water jacket with serpentine and spiral channels is replaced by a TO water jacket. The principles of temperature and thermal error control are based on several key mechanisms. Temperature regulation is achieved by regulating the temperature and flow rate of the cooling water through the temperature control system, ensuring that the cooling water effectively removes the heat generated during the operation of the boring machine spindle. The temperature sensors at the inlet and outlet contribute to real-time monitoring, ensuring stable and continuous cooling performance. Feedback regulation allows the operator to modify the cooling water temperature and flow rate based on the real-time data, keeping the boring machine spindle temperature within an appropriate range, thereby reducing thermal errors caused by temperature fluctuations. By precisely controlling and monitoring the temperature, the temperature and thermal errors of the spindle system are effectively reduced, improving the machining accuracy and stability of the boring machine.

[0125] 2.1. Control equations

[0126] The continuity equation for a steady incompressible flow that ensures mass conservation is mathematically expressed as:

[0127]

[0128] where: u is the velocity vector; denotes the gradient operator.

[0129] In turbulence modeling, the Reynolds-averaged Navier-Stokes (RANS) model has extensive applications in simulating fluid flow. Due to the non-linear nature of the solution process, its complexity is notable. However, applying the RANS model in the TO of thermo-fluid coupling problems may require a large amount of computation and may lead to non-convergence issues. Therefore, to address these challenges, a simpler and more computationally efficient Darcy model is usually adopted to replace the fluid flow in this case, thus greatly reducing the computational requirements and improving the convergence reliability, making it a practical choice for the TO of complex thermo-fluid systems.

[0130]

[0131] where: p is the pressure, and κ and μ represent the fluid permeability and dynamic viscosity, respectively.

[0132] For incompressible flow, especially when considering transport phenomena involving temperature, the description can be extended to include the transfer and diffusion of temperature. The fundamental equations governing these processes in an incompressible fluid are:

[0133]

[0134] where: T is the temperature; D T is the thermal diffusivity.

[0135] The heat source Q(x) depends on the design domain x:

[0136]

[0137] where: Ω is the design domain; Θ is the heat source domain; x ∈ Ω\Θ is the design domain excluding other regions except the heat source.

[0138] 2.2. Topology Optimization

[0139] 2.2.1. Material Interpolation

[0140] Interpolate the permeability in the flow equation. κ s and κ f represent the solid permeability and fluid permeability, respectively. When γ = 0 and κ = κ f , it is a fluid. When γ = 1 and κ = κ s , it is a solid. Here, γ represents the design density.

[0141]

[0142] Among them: It is necessary to calculate κ f to satisfy and represent the pressure drops of the Darcy and RANS models respectively; p κ is the penalty parameter of the permeability.

[0143] represents fluid flow. Assign a small value to κ s while u≈0, indicating that there is no flow at this position. The thermal diffusivity is interpolated into the energy equation:

[0144]

[0145] Among them: is the penalty parameter of the thermal diffusivity; D Ts and D Tf represent the solid thermal diffusivity and the fluid thermal diffusivity respectively.

[0146] The thermal conductivity can be calculated by where K and rhoCp represent the thermal conductivity and the product of density and heat capacity respectively. In addition, when γ = 0, D T = D Tf , the following can be obtained: represents fluid heat transfer. When γ = 1, D T = D Ts , the following can be obtained: represents solid heat transfer. In the TO model, due to the use of the penalty strategy for material properties, it may lead to mesh dependence or checkerboard elements in the optimization results, thus affecting the convergence performance of the TO process. Therefore, a Helmholtz-type partial differential equation is adopted as a regularization filter:

[0147]

[0148] Among them: γ f represents the pseudo-density before filtering; γ a represents the pseudo-density after filtering; R min is the filtering radius.

[0149] The larger the filtering radius, the easier it is to erase the shape features of γ a . In addition, using the filter usually generates gray elements, which are fuzzy elements between fluid and solid material elements (0 < γ a < 1). Therefore, a projection function is used to achieve a clear fluid-solid boundary. In this embodiment, the hyperbolic tangent function is used for projection:

[0150]

[0151] Among them: η and β represent the projection point and the projection slope respectively. The projection function parameter η takes the value of V to ensure the volume ratio of the material. The projection slope β is controlled to gradually increase as the optimization process progresses. Specifically, the initial value of β is 1 and the maximum value is 200. The optimized variables after projection almost converge to a set close to {0, 1}, and then a clear fluid-solid boundary is obtained.

[0152] 3. Sensitivity calculation and implementation of topology optimization program

[0153] 3.1 Adjoint model

[0154] The continuous adjoint method is used to calculate the sensitivity of the thermal-fluid problem, which is crucial for understanding how changes in design parameters affect the objective function. In this method, the objective function and the constraints are combined into a Lagrangian function. This integration allows the simultaneous consideration of the objective and the constraints, thus promoting efficient optimization. The adjoint method effectively calculates the gradient of the objective function with respect to the design variables, even for complex systems, thus providing a powerful tool for TO in engineering applications that are crucial for both performance and compliance with constraints.

[0155] L(u, u A , p, p A , T, T A , λ1, θ2, Ω)

[0156] = Ψ + ∫ Ω u A ·R u dΩ + ∫ Ω p A R p dΩ + ∫ Ω T A R T dΩ + λ1J1 + λ2J2

[0157] Among them: u A , p A , T A , λ1 and λ2 represent the adjoint velocity, the adjoint pressure, the adjoint temperature, the Lagrange multiplier related to the volume constraint, and the Lagrange multiplier related to the pressure drop constraint respectively. u A , p A , T A , λ1 and λ2 are all Lagrange multipliers, which are introduced as dummy variables.

[0158] The adjoint model obtained is:

[0159]

[0160] Since the governing equations And the constraints of J1 and J2 are theoretically constant, so L = Ψ; then the sensitivity of the objective function is The following relational expressions are obtained:

[0161]

[0162] Find the appropriate adjoint variables u A , p A , T A , λ1 and λ2, such that Then the final Lagrangian equation is constructed as follows:

[0163]

[0164] Where: is the control equation.

[0165] 3.2. Sensitivity calculation

[0166] To obtain the optimal TO result, the following conditions must be satisfied, which can be proven mathematically.

[0167]

[0168] Where: δ ζ Ψ, δ ζ J and are the differentials of the objective function Ψ, the constraint condition J and the control equation i with respect to the change of the variable ξ; ξ is a variable, and ξ ∈ (u, p, T).

[0169] The differential form is obtained:

[0170]

[0171] Where: ∫ Ω *dΩ and ∫ Γ *dΓ represent the integral in the design domain and the integral on the boundary respectively;

[0172] For any δu, δp and δT, the calculation requirements in the design domain are satisfied, and the adjoint equation is obtained as follows:

[0173]

[0174] Considering the adjoint boundary conditions:

[0175]

[0176] For the inlet and the wall, the velocity and temperature are fixed values, so δu = 0 and δT = 0, and the adjoint velocity and adjoint temperature of the inlet and the wall are:

[0177]

[0178] T A = 0

[0179] For the outlet, the pressure is 0 and the temperature is a zero gradient, so δp = 0 and the adjoint velocity and adjoint temperature at the outlet are:[[]]

[0180] p A = 0

[0181]

[0182] where:[[]] is the adjoint normal velocity.[[]]

[0183] Finally, the sensitivity calculation method for the thermal-fluid topology optimization problem is:[[]]

[0184]

[0185] where: δ γ is the volume factor; κ s and κ f represent the solid permeability and fluid permeability respectively; p k is the penalty parameter for permeability; is the penalty parameter for thermal diffusivity; and represent the solid thermal diffusivity and fluid thermal diffusivity respectively.[[]]

[0186] 3.3. Adaptive Mesh Refinement

[0187] AMR (Adaptive Mesh Refinement) is used to promote precise mesh refinement in key areas of interest, especially at the fluid-solid boundary. This targeted approach reduces the use of unnecessary mesh elements, thus saving mesh resources and improving computational efficiency. In addition, refining the mesh at the fluid-solid interface ensures a more accurate simulation of the physical phenomena occurring near these critical boundaries. AMR employs the h-type refinement method, which involves dynamically adjusting the mesh by adding or deleting nodes. This process allows the elements to be divided or aggregated according to the specific needs of the simulation region. For example, as Figure 3 shown, introducing an octree structure can generate sub-grids from the original hexahedral mesh. This operation divides the parent mesh into eight smaller sub-grids, creating 36 faces, including 12 internal faces that are part of the original parent mesh. This method significantly enhances the mesh's ability to adapt to complex geometries and flow dynamics, improving the overall fidelity and accuracy of the thermal-fluid TO simulation.[[]]

[0188] Density Gradient Used to represent the rate of change of design variables. This gradient is a key metric for determining where adaptive meshing should be applied. Specifically, when the density gradient falls within a predefined range, the adaptive meshing process in these regions is triggered, ensuring that the mesh is finely tuned to capture subtle changes in the material distribution. Conversely, when the density change exceeds this specified range, indicating a more uniform material distribution, the mesh refinement is reduced and the mesh reverts to its initial coarse state. This adaptive approach allows for an increase in mesh density where the design variables change rapidly, providing higher resolution and accuracy in critical areas while conserving computational resources in areas with less change. This dynamic adjustment of mesh density based on the local requirements of the TO model optimizes computational efficiency and simulation accuracy.

[0189]

[0190] Where: AMR l and AMR h represent the lower and upper bounds of the adaptive physical field, respectively.

[0191] When the results of the latest calculated objective function and constraint conditions are below a predefined error threshold, the AMR strategy is activated. This threshold is established to determine when the resolution of the current mesh is insufficient to accurately capture the dynamics, thus requiring improvement.

[0192]

[0193] Where: ε Ψ is the objective function error threshold; ε V is the volume constraint error threshold; ε p is the pressure drop error threshold; ∫ Ω *dΩ is the integral in the design domain; V frac is the flow channel ratio; p n is the pressure at the current iteration step; p n―1 is the pressure value at the previous iteration step.

[0194] 3.4. Parallel Computing

[0195] The domain decomposition method is used in OpenFOAM to facilitate parallel computing, improving the efficiency of computational tasks. In this method, as Figure 4As shown, the structured grid of the computational domain is divided into multiple subdomains, each subdomain being distinguished by a different color. This segmentation allows the computational tasks within each subdomain to be executed simultaneously on various computational processors, significantly accelerating the simulation process. The steps of parallel computing are as follows: (1) Parallel computing mechanism. (2) Communication and processor management. (3) Efficient data processing. By using these tools and methods, OpenFOAM can effectively handle the complexity of parallel processing in CFD, achieving faster and more accurate simulations, suitable for complex engineering analyses encountered in thermal fluid TO.

[0196] The parallel solution processor in OpenFOAM includes establishing a parallel execution environment, decomposing the grid, running a parallel solver, and performing post-processing and analysis. Generally, the configuration of parallel parameters is described in the decomposeParDict dictionary file located in the system subdirectory of the case folder. The comprehensive settings are shown in Table 1.

[0197] Table 1 Parallel computing settings

[0198]

[0199] where n x 、n y and n z represent the number of decompositions in the x, y, and z directions respectively, and np = n x ×n y ×n z .

[0200] 4. Numerical verification of thermal fluid topology optimization

[0201] 4.1 Verification of the finite volume method solver

[0202] Before performing thermal fluid TO using the Darcy model, the accuracy of FVM programming must be verified. For this purpose, the thermal fluid simulation results using the Darcy and RANS models in COMSOL Multiphysics 6.1 were evaluated. The simulation was benchmarked and compared using the computational domain model shown in Figure 5 (a), which depicts a cuboid of 3m × 1m × 1m. Within this domain, the Darcy model simulated the flow of a blue cuboid cross-section (dimensions 3m × 0.2m × 0.2 m) composed of a material with a permeability of 1×10 - 9 m 2 、rhoCp of 5000. The rest of the region was composed of a material with a permeability of 1×10 -14 m 2, composed of a material with rhoCp of 20000. In contrast, in the RANS model, the central blue cuboid region represents fluid flow, and the remaining space is considered solid, with a kinematic viscosity and thermal conductivity of 1×10 -5 Pa·s and 0.01 W / (m·K), respectively. The inlet temperature is set at 273.15 K, and the heat flux from the heat sources at the upper and lower boundaries is 1 W / m 2 , while the remaining walls are insulated and maintained with a non-slip boundary condition. The mesh parameters in COMSOL and OpenFOAM are kept consistent. Figure 5 (b) shows the cross-sectional temperature field distribution at y = 0.5 m. Clearly, the calculation results of COMSOL-Dacy (Darcy model in COMSOL), COMSOL-RANS (RANS model in COMSOL-RANS), and the proposed OpenFOAM model are in very good agreement, showing similar overall temperature distributions. To further evaluate the numerical differences, the pressure and temperature along the specified red line for different solutions were analyzed, as shown in Figure 5 (c). Inside the fluid domain, the temperature remains constant at 273.15 K. Inside the solid domain, the temperature results of the COMSOL-Dacy and COMSOL-RANS models are exactly the same, while the temperature calculated by the proposed OpenFOAM model is slightly higher compared to the temperatures of the other two models. Since the COMSOL-RANS model does not provide the pressure distribution inside the solid domain, the pressure distributions between COMSOL-Dacy and OpenFOAM calculations were compared. Although slight deviations were observed at the start and end segments of the pressure distribution line, the overall similarity is still well maintained.

[0203] Table 2 provides a comparative analysis of the results obtained from the above three numerical models, with the results produced by the COMSOL-RANS model as the baseline. The differences in the average temperature, maximum temperature, and pressure drop calculated by the COMSOL-Dacy model are 0.03%, 0.05%, and 0.28%, respectively. These minimal errors emphasize the feasibility of using the Darcy model to simulate fluid flow, confirming the findings discussed previously in [5]

[18] . The changes in the average and maximum temperatures calculated by OpenFOAM are slightly larger than those observed in the COMSOL-Dacy model, with a maximum deviation reaching 0.86%. Nevertheless, this difference is completely within the acceptable error range. The pressure drop error matches the error observed in the COMSOL-Dacy model. Therefore, it is certain that the numerical results obtained from COMSOL and OpenFOAM have a high degree of consistency, thus verifying the accuracy and effectiveness of the FVM programming used in these calculations.

[0204] Table 2 Comparison of Results of Three Numerical Models

[0205]

[0206] 4.2, Design of the Cooling Plate Freely Installed on the Boring Machine Slide Block

[0207] According to Figure 2 (a), the design of the cooling plate freely installed on the drilling machine slide block proves the effectiveness of the thermal fluid TO solver. Figure 6 The schematic diagram of the 3D TO model of the freely installed cooling plate is shown. The main design domain is described as a blue cuboid with dimensions of 35 mm × 7.5 mm × 1 mm. The fluid inlet and outlet are represented by two gray cuboids, each with dimensions of 5 mm × 2.5 mm × 1 mm. The dashed line in the middle marks a symmetry surface. The fluid is vertically guided into the inlet from the x direction, while the outlet boundary is set to p in = 0. For hydrodynamics, all walls are subject to the no-slip boundary condition, ensuring that the velocity on the wall is zero, i.e., u out = 0 and wall All walls are treated with the adiabatic boundary condition which means that no heat transfer occurs through these boundaries. A volume constraint of 0.5 is imposed on the design domain, and the refinement strategy of AMR includes specific parameters, which are ε = 10 Ψ = 10 ―4 , ε p = 10 ―4 and ε V = 10 ―3 . The material properties related to the simulation are shown in Table 3 in detail.

[0208] Table 3 Material Properties

[0209]

[0210] In the TO study, the changes in the objective function, constraints, and the number of grids during the iterative process are as shown in Figure 7 . Initially, both the average temperature and the pressure drop showed significant fluctuations. However, as the iteration progresses, these variables begin to stabilize, indicating convergence to the optimal solution. When the conditions of the AMR strategy are met, the grid will undergo its initial refinement, rapidly upgrading from 22528 elements to 88139 elements. The significant increase in grid density leads to significant changes in the simulation results, especially a sharp rise in the average temperature and a small but noticeable change in the pressure drop. In the 370th iteration of the optimization process, another round of grid refinement was carried out. This adjustment expanded the number of elements from 88139 to 433974 elements. This further improvement shows that the average temperature rises again while the pressure drop remains relatively stable, indicating the influence of grid density on the heat transfer performance of the cooling plate without a corresponding increase in fluid resistance.

[0211] To further investigate the influence of AMR on the TO results, the density gradients and velocity distributions across different cross-sections of the freely mounted cooling plate were analyzed at 200, 360, and 400 iterations. The TO results showed that the trends for different levels of mesh refinement were consistent. Figure 8 (a) shows the TO results without any refinement, where the constraints imposed by the original mesh size led to a jagged and uneven surface at the fluid-solid boundary. This irregularity highlights the limitations of the unrefined mesh in accurately defining complex boundaries. Figure 8 (b) and Figure 8 (c) show the TO results after one and two levels of refinement, respectively. With each additional level of refinement, the boundary between the solid and fluid domains becomes increasingly clear and smooth. This improvement significantly enhances the manufacturability of the design, ensuring that the refined geometry is more suitable for machining and manufacturing processes. Importantly, all refinement efforts were strategically focused on the focal domain, maintaining the original mesh size in the area. This targeted approach ensures that the computational cost remains controllable, avoiding excessive increases while still significantly improving the accuracy and quality of the TO results.

[0212] According to Figure 7 , AMR has a significant impact on the average temperature and pressure drop throughout the TO process. To investigate this influence in more detail, the TO result data for different levels of refinement along the specified red line in Figure 8 (c) were analyzed, as shown in Figure 9 . Initially, in the regions along the line where there is no fluid-solid boundary domain, all meshes remained unrefined, resulting in three overlapping data lines. As the line extends, the middle section intersects a total of seven fluid-solid boundaries, causing the temperature curve to exhibit seven distinct inflection points. At the highest refinement level (level 2), the temperature reaches its peak, followed by level 1, with level 0 representing the lowest temperature. Notably, the high-temperature regions within the refined mesh are mainly located in the solid region, indicating the presence of local heating effects due to restricted fluid flow. Regarding the fluid pressure, as the pipeline length increases, the trend generally decreases, with the smallest fluctuations observed. Slightly higher pressures correspond to increased refinement levels, although the overall differences are small, confirming the minor impact of mesh refinement on the overall pressure drop, as shown in Figure 9 (b). Figure 9 (c) shows the velocity distribution along the line, indicating that due to the extremely low permeability of the material, the velocity is highest at the inlet and approaches zero in the region. This supports the applicability of using the Darcy model to represent fluid flow within the TO framework, as the decrease in velocity is consistent with the expected results based on the model assumptions. Similar to the pressure distribution, changes in the AMR level have a negligible impact on the fluid velocity distribution, highlighting the robustness of the flow characteristics to changes in mesh granularity.

[0213] The temperature gradient values on the cross-section are as Figure 10 shown, illustrating the significant variations due to the mesh refinement levels in the solid domain (grey) and the fluid domain (blue). It is noteworthy that without mesh refinement, the temperature gradients in these domains are nearly identical. However, in the seven fluid-solid boundary domains (dark grey), significant changes in the temperature gradient occur due to the influence of AMR. When using AMR, the proximity of the mesh elements decreases, which greatly amplifies the temperature gradient. Specifically, in the first fluid-solid boundary domain, the maximum temperature gradients for refinement levels 0, 1, and 2 are 3560, 6732, and 13870 respectively. The gradient at level 2 is approximately four times that at level 0 and twice that at level 1, indicating an inverse relationship with the element size. This significant increase explains why Figure 7 (a) shows significant temperature variations in the level 2 TO channels at the boundary domain. Additionally, the temperature difference in the solid region is more pronounced due to the significantly higher thermal diffusivity of the solid compared to the fluid. Therefore, the implementation of AMR results in a finer boundary layer that can capture temperature gradient fluctuations more precisely. This improvement significantly enhances the accuracy of the numerical simulation in terms of heat conduction processes, which is crucial for accurately predicting heat transfer at the fluid-solid boundary.

[0214] Figure 11 shows the TO results obtained under different pressure drop constraints. In Figure 11 (a), the streamline intuitively represents the pressure distribution of different TO channel designs, labeled a-1 to a-6. As the complexity of the TO cooling channel structure and the number of branches increase, the pressure drop also increases accordingly. When the fluid flows through the cooling channel, the pressure in all TO channels shows a downward trend, starting from the maximum pressure at the inlet and reaching 0 Pa at the outlet. Figure 11 (b) depicts the temperature field of the cooling channel. The temperature in the TO cooling channel gradually decreases and remains uniformly distributed from the top to the bottom of the design domain. In particular, for channels A-1 and A-2, the temperature rises sharply in the region without fluid flow. As the cooling channels are more evenly distributed throughout the design domain, both the heat dissipation efficiency and the uniformity of the cooling plate temperature distribution are significantly improved, thereby enhancing the overall cooling performance of the entire cooling plate.

[0215] To evaluate the heat dissipation performance of the TO channels, a traditional rectangular channel with the same fluid volume ratio was designed for comparison. Figure 12Shows the temperature and flow rate distributions within these cooling channels. Notably, the temperature in the traditional rectangular channels increases significantly at the end of the side cooling channels, while the temperature within the TO channels remains notably more uniform. Specifically, the maximum temperature in the traditional rectangular cooling channels is 19.88 °C, and the average temperature is 7.65 °C. In contrast, the maximum temperature in the TO cooling channels is 12.88 °C, and the average temperature is 5.91 °C. Compared with the traditional design, the average temperature and the maximum temperature of the TO cooling channels are reduced by 35.21% and 22.75% respectively. This improvement not only enhances the overall heat dissipation efficiency of the cooling channels but also significantly improves the temperature uniformity.

[0216] Streamlines show that the rectangular and TO cooling channels maintain similar flow velocities in the main channel. However, the TO cooling channel is characterized by a dispersed design throughout the area, with large and small branches in the main channel, increasing the fluid-solid heat transfer area. This configuration allows the TO cooling channel to dissipate heat more effectively and maintain a lower overall temperature.

[0217] Regarding the pressure drop, the recorded pressure drop in the rectangular cooling channel is 786.89 Pa, while the pressure drop in the TO cooling channel is slightly higher, at 789.68 Pa, indicating a minimum difference of 3.5%. Although the TO cooling channel has a larger heat transfer area, which may result in greater fluid flow losses, its streamlined shape minimizes the flow resistance. In contrast, the rectangular cooling channel and its sharp corners generate secondary flows, which enhance fluid perturbation and result in a slightly higher pressure drop.

[0218] Ultimately, the comparative analysis confirms that the TO cooling channel has significant advantages in terms of heat dissipation and flow performance. This case study effectively validates the effectiveness of the TO cooling method in the design of cooling channels, demonstrating its superiority over the traditional rectangular cooling channel design in terms of thermal control and fluid dynamics.

[0219] The pressure drop and average temperature of various TO cooling channels are as Figure 13 shown, highlighting the inherent trade-off between these two conflicting physical parameters. As the cooling channel contains numerous branches and is more widely dispersed, the average temperature decreases. This indicates that the strong cooling performance of the cooling channel corresponds to a high pressure drop, suggesting a compromise where the enhanced energy dissipation of fluid flow requires increased pressure support. As the pressure drop decreases, the design of the TO cooling channel focuses on the center of the cooling channel at the inlet and outlet. This configuration uses the shortest flow distance to minimize flow losses, which in turn leads to a significant increase in the temperature of the solid region without the cooling plate. The average temperature of the TO channel rises from 5.91 °C to 42.07 °C, while the pressure drop decreases from 960.57 Pa to 664.41 Pa.

[0220] By fitting Figure 13The curves on it yielded the Pareto front for the thermal-fluid TO problem related to the cooling channels. Supported by an R-squared value of 0.9353, the curve fitting showed high precision. This Pareto front allows designers to select the most suitable cooling channel structure for specific practical engineering applications to balance heat dissipation and hydraulic performance. This method simplifies the optimization process and reduces the need to continuously adjust design parameters during the TO process, thus simplifying and enhancing the design workflow.

[0221] 4.3. Spindle System Water Jacket Design

[0222] According to Figure 2 (b), as Figure 14 shown, the geometric model of the water jacket design domain was obtained. The design domain was constructed as a semi-cylinder with an inner diameter of 83 mm, a height of 86 mm, and a thickness of 4.5 mm. To optimize the calculation efficiency, the middle surface of the water jacket was considered as the symmetry surface. The design domain also had a uniform heat source with a strength of 10000000 W / m3. The fluid dynamics was initiated by the fluid entering perpendicular to the red surface at a speed of 0.5 m / s and leaving perpendicular to the green surface. The inlet temperature Tin was set at 0 °C, and the outlet pressure pout was set at 0 Pa. All cylinder walls were designed to be anti-slip and adiabatic to prevent heat transfer through these surfaces. A volume constraint was applied, set as V = 0.6. The refinement strategy of AMR included specific parameters that were crucial for improving the accuracy of the simulation by refining the mesh in critical areas while effectively managing computational resources. That is, ε Ψ = 5×10 ―5 , ε p = 5×10 ―5 and ε V = 10 ―3 . Table 3 lists the material properties crucial for the simulation, providing the necessary data for the accurate modeling and analysis of the thermal-fluid dynamics in this specific design domain.

[0223] Figure 15 shows the TO results focusing on the pressure drop Lagrange multiplier. The upper part shows the distribution of the density gradient, while the lower part shows the velocity distribution of the entire design domain. In Figure 15 (a), before applying AMR, the flow channels were evenly distributed within the design domain, with the highest velocities at the inlet and outlet. Additionally, the visible jagged edges of the mesh indicated the initial unrefined state. As Figure 15 (b) shows, the TO results after implementing AMR. The overall distribution structure of the flow channels remained basically unchanged. However, in the areas where AMR was applied, especially in those places, there were significant improvements. Figure 15(b) The magnified details on the right side are highlighted. The mesh at the fluid-solid boundary is refined in two layers, effectively reducing the element length of the level-2 mesh to one-fourth of the level-0 mesh. The use of AMR enables a clearer boundary definition, thus significantly improving the accuracy of the manufacturing process. The fine mesh not only enhances the visual clarity and accuracy of the thermal-fluid TO model but also ensures a more faithful representation of the design in practical applications, thereby promoting more accurate and efficient manufacturing results.

[0224] Figure 16 shows the iterative history of the TO results corresponding to the Figure 15 scenario shown. The convergence of temperature and pressure drop during the entire optimization process is shown in Figure 16 (a) and (b) respectively. When the predefined parameters of AMR are met, AMR is initiated at the 305th iteration step. This results in an increase in the number of meshes from 13,750 to 69,071. Subsequently, the mesh count is further refined and upgraded to 404,742 at the 450th iteration. After an additional five iterations, the criteria for AMR are met again, and the number of meshes slightly decreases to 367,005. The convergence tolerances for the average temperature, pressure drop, and volume ratio change significantly at the beginning of the iterative process, reflecting the initial adjustment and stabilizing as the iteration progresses. Notably, each time AMR is applied, the convergence tolerance for the average temperature changes significantly, indicating the sensitivity of temperature calculation to mesh refinement. In contrast, the convergence tolerances for pressure drop and volume ratio do not exhibit substantial fluctuations, indicating that these parameters are relatively stable despite the change in mesh density. This iterative history emphasizes the dynamic nature of the TO results under AMR conditions, illustrating how mesh adjustment significantly affects the convergence behavior of different optimization metrics.

[0225] The pressure and temperature distributions in the TO results are shown in Figure 17 as shown. The pressure and temperature distributions vary with different pressure drop constraints. As the constraints intensify, the TO channels labeled B-1 to B-6 develop increasingly complex structures. This complexity is manifested as a proliferation of flow channel branches and the simultaneous narrowing of these channels. The design is adjusted by dispersing the flow from the central axis to the periphery, thereby achieving a more uniform temperature distribution throughout the design domain.

[0226] The material distribution within the cross-section of the TO cooling channel B-6 is shown in Figure 18As shown. At specific axial planes -z = 10 mm, z = 30 mm, and z = 40 mm, the cross-section shows a pattern where the fluid is distributed along the outer periphery while the solid material is concentrated on the inner periphery. This arrangement results in a non-uniform material distribution around the circumference, which is in stark contrast to the more consistent configurations typically shown in cooling channels. These observations highlight how variable pressure drop constraints can significantly affect the structural and functional characteristics of the TO channels, thereby optimizing the design to better distribute heat and fluid pressure throughout the cooling channels. This configuration is crucial for achieving the target performance goals in a cooling system, where maintaining a balanced temperature and pressure distribution is essential for efficient operation.

[0227] A comparative analysis of the cooling capacity and efficiency of various TO cooling channels is as Figure 19 shown. The cooling capacity and efficiency are consistent with the Figure 17 results shown. By adjusting the magnitude of the pressure drop Lagrange multiplier, the pressure drop across these channels can be effectively controlled, denoted as. As the pressure drop Lagrange multiplier increases, the pressure drop decreases, while simplifying the complexity of the cooling channels. For example, when the pressure drop Lagrange multiplier is set to 0.06, the pressure drop is 40.41 Pa, the average temperature is 3.22 °C, and the maximum temperature is 8.01 °C. Conversely, increasing the pressure drop Lagrange multiplier to 0.6 results in a pressure drop reduction of 27.98 Pa. Correspondingly, the average temperature and the maximum temperature increase to 9.66 °C and 25.84 °C, respectively. This trend is consistent with what is typically observed in the field of cooling plate design, where enhanced heat dissipation performance corresponds to higher pumping power. Therefore, the selection of pumping power can be controlled according to specific operating conditions, thereby optimizing the energy efficiency and cooling effect in the designed TO cooling channels. This approach helps to strike a balance between minimizing energy consumption and maximizing heat dissipation performance, which is crucial for achieving optimal performance in practical engineering applications.

[0228] 4.4. Comparison of Different Solvers

[0229] This example evaluates the computational efficiency and resource utilization of two solvers: one using the finite element method in the commercial software COMSOL Multiphysics 6.1, and the other using the FVM through the open-source platform OpenFOAM. The computing hardware includes an Intel(R) Core(TM) i7-1700 CPU 2.90GHz, 64GB RAM (63.8GB available). The operating system is a 64-bit version running on an X64 processor architecture, which can efficiently handle complex computational tasks. COMSOL Multiphysics 6.1 runs on a Windows 10 64-bit operating system, providing a stable and familiar environment for simulations requiring detailed finite element analysis. In contrast, OpenFOAM 10 runs in a Ubuntu 20.04 virtual machine hosted on VMware Workstation Pro, which is allocated 16.0GB of memory and configured to use an 8-core processor setup. This arrangement ensures that OpenFOAM has sufficient resources to effectively perform the thermal-fluid TO simulation, although within the scope of the virtualized environment, this may introduce some performance overhead.

[0230] To ensure the comparability of the simulations conducted by the two solvers, the same geometric model described in Section 4.3 is used. Table 4 details the mesh statistics and solver settings, and all other parameter settings are kept consistent across the two platforms. This setup allows for a direct comparison of the two solvers in terms of computational time, efficiency, and resource consumption, providing valuable insights into the applicability of each method to specific types of thermal-fluid TO simulations.

[0231] Table 4 Mesh details and solver settings

[0232]

[0233] The TO results obtained from the above two different solvers are as Figure 20Shown: one is implemented through the OpenFOAM platform and the other through COMSOL Multiphysics. The results of OpenFOAM are post-processed using Paraview 5.6, enhancing the visualization capabilities, while the results of COMSOL are analyzed using its integrated post-processing module, which provides native tools tailored for such evaluations. The cooling channel distributions calculated by the two software platforms show a high degree of overall similarity, indicating that both solvers effectively applied the TO principle under the given pressure and volume constraints. Generally, these results have a wide main cooling channel extending between the inlet and outlet, with two narrower branches on both sides, optimizing fluid flow and heat dissipation within the constraints. Despite these similarities, some differences in the cooling channel configurations can be observed. In the TO results obtained by COMSOL, the main cooling channel is almost parallel to the x-axis, with a symmetric distribution on both sides. This indicates a very structured approach to cooling channel optimization, influenced by the specific handling of the optimization algorithm by the solver. On the other hand, the TO results obtained by OpenFOAM show that the main cooling channel protrudes upward, and the branched cooling channels are slightly backward relative to the main cooling path. This variation may reflect differences in each software's handling of mesh refinement, optimization criteria, or solver dynamics.

[0234] These differences were further explored by presenting the fluid pressure field distributions and temperature fields calculated by the two software platforms, as Figure 21 shown. Despite the minor differences in the cooling channel layouts, the calculated physical fields exhibit similar trends, further confirming that both solvers are capable of obtaining comparable and reliable results in thermal-fluid TO. Cross-validation between the two different computational tools not only increases confidence in the simulation results but also demonstrates the flexibility and robustness of advanced numerical methods in addressing complex engineering challenges.

[0235] A detailed comparison of the computational results obtained by the above two different solvers is shown in Table 5. The average temperature and the maximum temperature obtained by OpenFOAM are 3.27% and 0.55% lower, respectively, than those obtained by COMSOL. Conversely, the pressure drop obtained by OpenFOAM is 0.14% higher than that calculated by COMSOL. These findings emphasize that, despite the minor differences in the physical field values calculated by the finite element method used in COMSOL and the FVM used in OpenFOAM, the differences between the two different solvers are not significant. This indicates that both of these computational methods are robust and reliable for engineering applications, providing reliable results within reasonable tolerances. The minor differences in the computational results stem from the inherent characteristics of each numerical method, including how they handle boundary conditions, mesh discretization, and solver algorithms.

[0236] Table 5 Comparison of Computational Results Using Two Different Solvers

[0237]

[0238] Table 6 conducts a comparative analysis of the computational efficiency and resource utilization of OpenFOAM and COMSOL in terms of thermal-fluid TO. OpenFOAM demonstrates excellent computational efficiency, which is mainly attributed to its ability to utilize physical and logical processing units, as well as its strong support for parallel computing. This makes the processing speed of OpenFOAM significantly faster than that of COMSOL, which only relies on physical cores. In addition, COMSOL uses automatic differentiation, especially in the "black box" method, which results in a large amount of memory requirements due to the need to store intermediate variables during the backpropagation of derivatives

[33] . This significantly affects the memory efficiency, especially for those involving complex derivative calculations. In terms of mesh processing, the AMR implementation of OpenFOAM manages approximately 200,000 meshes, and the mesh changes very little under the influence of different partition boundaries. This AMR method improves the computational performance by optimizing the mesh density according to the needs of thermal-fluid TO simulations, thereby reducing unnecessary computational loads.

[0239] The application of parallel computing in OpenFOAM significantly reduces the computational time of thermal-fluid TO models with a larger number of degrees of freedom, by 93.89%, highlighting its significant efficiency advantage over COMSOL. However, it is worth noting that the efficiency improvement of parallel computing is not linearly inversely proportional to the number of processors used. The data listed in Table 6 shows that while increasing the number of processors generally leads to faster calculations, the rate of increase in computational speed is not always directly related to the increase in the number of processors. This is partly due to the overhead associated with communicating between processors using MPI, as well as potential inefficiencies such as high central processing unit (CPU) utilization and memory saturation. For example, when using 2, 4, and 8-core processors, the speedup ratios are 2.53, 3.55, and 3.27 respectively. This indicates that as the number of core processors increases, the returns gradually decrease, which illustrates the importance of optimizing the number of core processors to achieve optimal performance without consuming unnecessary resources. Therefore, before fully deploying a parallel computing strategy in a practical scenario, it is necessary to carefully evaluate the parallel speedup ratio of the thermal-fluid TO model. This ensures that the resource allocation is optimal and the computational execution efficiency is maximized without compromising system stability or incurring excessive operating costs.

[0240] Table 6 Execution Time and Resource Occupancy

[0241]

[0242] 5. Practical Applications of Thermal-Fluid Topology Optimization

[0243] 5.1 Application of a Freely Mounted Cooling Plate in a Boring Machine

[0244] Manufacture the freely mounted cooling plate designed in Section 4.2, as Figure 22 shown. Then mount the freely mounted cooling plate on the slide of the boring machine to study its influence on the thermal behavior of the boring machine.

[0245] Circulating cooling water is supplied to the cooling plate, and the cooling plate is in contact with the slide. Then, when the feed rates of the X, Y, and Z axes are 10 m / min, the thermal behavior of the boring machine is measured. An infrared thermal imager is used to measure the temperature field distribution of the boring machine, as Figure 23 (a) shown. Temperature sensors are used to measure the temperatures of key points, including the left bearing of the Y axis, the upper bearing of the X axis, and the left screw of the Y axis. Five displacement sensors are used to measure the thermal error between the spindle system and the center point of the worktable, as Figure 23 (b) shown.

[0246] The thermal behavior of the boring machine is analyzed without installing the cooling plate, as Figure 24 shown, with a total of 896731 elements and 235411 nodes. The boring machine bed is fixed to the ground. The ball screw feed drive systems of the X, Y, and Z axes are of double drive structure, with a total of 6 ball screws and 6 motors. In this paper, the temperature and thermal deformation of the boring machine at a feed rate of 10 m / min for the X, Y, and Z axes are analyzed. According to

[34]

[35] , the thermal boundary conditions of the boring machine are calculated, as shown in Table 7.

[0247] Table 7 Thermal boundary conditions at a feed rate of 10 m / s

[0248]

[0249] The temperature field of the precision boring machine is as Figure 25 (a) shown. The temperatures of the motors of the X, Y, and Z axes are higher than those of other components, and the highest temperature of the upper motor of the X axis is 52.6 °C. The temperature rises of heating components such as nut pairs, rolling guide pairs, and bearing pairs are obvious. However, it is not easy to install the cooling plate on these components. Therefore, it is considered to install the cooling plate on the large basic components of the boring machine, including the slide, the bed, and the column. The temperature rise of the X axis is higher than that of the feed drive systems of the Y and Z axes because the external load applied to the X axis is greater than the loads applied to the slide and the Y axis. The temperature gradient distribution of the boring machine is as Figure 25 (b) shown. The regions with larger temperature gradients are concentrated near the heating components. In addition, the average temperature gradient on both sides of the Y-axis slide is 4.36 K / m. Installing cooling plates on both sides of the slide is expected to reduce thermal deformation. The thermal deformation field of the boring machine is as Figure 25(c). Since the bottom surface of the machine tool bed is fixed, the deformation of the bed is minimal, the column tilts backward away from the machine tool bed, the slide tilts forward towards the bed, and the spindle system exhibits obvious thermal deformation. The bed and the worktable slightly expand upward. For the entire boring machine, the maximum deformation is approximately 252 μm and occurs at the motor.

[0250] As Figure 26 shown, in order to adapt to the dimensions of the slider structure, the cooling plate designed in Section 4.2 was magnified ten times and then arranged on the slider of the boring machine. The heat transfer coefficient of the cooling plate was calculated to be 1000 W / (m 2 ·K). The thermal behaviors of the precision boring machine with and without cooling were compared and analyzed.

[0251] 5.1.1. Temperature Field

[0252] The temperature of the boring machine with cooling is as Figure 27 shown. The overall temperature distribution pattern is similar to Figure 25 (a). The highest temperature is also 52.6 °C because the highest temperature appears at the X-axis motor and the cooling plate does not directly cool the X-axis motor. Therefore, the highest temperature of the boring machine remains unchanged. The difference is that the temperature at the positions where four cooling plates are arranged on the slider is significantly reduced.

[0253] As Figure 28 shown, the temperatures of the slides with and without cooling were obtained. The highest temperature still remains at the motor. The lowest temperature of the slider with cooling is 20 °C, which is the temperature of the coolant. The lowest temperature of the slide without cooling is 20.4 °C. The reduction in the temperature rise of the slider is beneficial for reducing thermal errors.

[0254] The dynamic changes in the temperatures of the key components of the boring machine are as Figure 29 shown. Under the two cooling conditions, the temperatures of the left nut of the Y-axis and the front upper bearing of the X-axis are almost the same. The temperature of the left front bearing of the Y-axis is 24.22 °C without cooling and 23.88 °C with cooling. Since the cooling plate is installed on the slider equipped with the Y-axis, close to the motor and the bearing, the heat generated by the Y-axis can be removed. Therefore, under the action of the cooling plate, the temperature of the left front bearing of the Y-axis is significantly reduced. The average temperature of the boring machine is 20.94 °C (without cooling) and 20.90 °C (with cooling), and the thermal equilibrium times of the boring machine with and without cooling are 4.1 hours and 4.3 hours respectively. The installation of the cooling plate shortens the thermal equilibrium time of the boring machine.

[0255] 5.1.2. Thermal Deformation

[0256] The thermal deformation of the precision boring machine after cooling is as Figure 30As shown. For the entire boring machine, the thermal deformation is asymmetric because the motor of the X-axis is installed on the left side of the boring machine, resulting in more obvious thermal deformation on the left side than on the right side. The upper motor of the X-axis experiences the largest thermal deformation. The maximum thermal deformations of the precision boring machine are 252 μm without cooling and 249 μm with cooling, respectively. Compared with the uncooled boring machine, the thermal deformations of the slide block, Y-axis, and spindle system of the cooled boring machine are significantly reduced. Then, the effectiveness of the cooling strategy and the efficient TO method for the thermal fluid problem in the cooling elements of the precision boring machine is verified. The deformation on the left side of the slide block is greater than that on the right side. The connection between the slide block and the X-axis nut is on the left side of the slide block, and then the heat entering the left side of the slide block is greater than that entering the right side of the slide block.

[0257] The average thermal deformations of the key components of the precision boring machine are shown in Table 8. The thermal deformation of the workbench is small, and the reduction rate is 0.14%. The reduction rates of the thermal deformations of the front and rear bearings on the X-axis are 3.52% and 6.06%, respectively. The reduction rates of the thermal deformations of the nut on the X-axis, the left front bearing of the Y-axis, and the spindle system are 15.02%, 15.07%, and 12.93%, respectively. The reduction rates of the thermal deformations of the nut on the X-axis, the left front bearing of the Y-axis, and the spindle system are more significant than those of the workbench and the front and rear bearings on the X-axis because the nut on the X-axis, the left front bearing on the Y-axis, and the spindle system are closer to the cooling plate than the workbench and the front and rear bearings on the X-axis. The average thermal deformation of the entire boring machine is 17.45 μm without cooling and 16.35 μm after cooling. When the designed cooling plate is used, the average thermal deformation of the entire boring machine is reduced by 6.29%.

[0258] Table 8 Average Thermal Deformations of Key Components

[0259]

[0260] To verify the effectiveness of the cooling plate, the thermal errors at the center points of the spindle system and the workbench are extracted. The total thermal errors and thermal error components at the center points of the spindle system and the workbench are obtained, as Figure 31 (a) shown. The thermal error component in the Z direction is the largest, followed by that in the X direction, and the smallest in the Y direction. By using the cooling plate, the thermal error components in the X, Y, and Z directions are reduced to a certain extent. Without cooling, the thermal errors of the spindle center point in the X, Y, and Z directions are 15.23 μm, 13.41 μm, and 64.68 μm, respectively. After cooling, these thermal errors are reduced to 14.27 μm, 11.53 μm, and 56.90 μm, respectively. Compared with when not in use, the thermal errors of the spindle center point in the X, Y, and Z directions are reduced by 6.30%, 14.02%, and 12.03%, respectively, when using the cooling plate. For the center point of the workbench, since the cooling plate is far from the workbench, the thermal error is not significantly reduced, as Figure 31As shown in (b). Therefore, the cooling plate has little influence on the thermal error of the center point of the workbench. With and without the cooling plate, the total thermal error and thermal error components at the center point of the workbench are consistent.

[0261] As Figure 32 shown, the thermal error changes of the precision boring machine are obtained under cooling and non-cooling conditions. With cooling, the thermal error is much smaller compared to the non-cooling case. Specifically, the total thermal error during cooling remains within a lower range, while the overall thermal error without cooling is significantly higher. Under cooling and non-cooling conditions, the experimental data is consistent with the simulation results. This indicates that the simulation model can accurately predict the thermal error of the precision boring machine under different conditions. The thermal error gradually increases with the running time and tends to be stable. Under cooling conditions, the thermal error rapidly increases to a certain level in the initial stage and then remains stable. Under non-cooling conditions, the thermal error continues to increase and reaches a higher stable value. In summary, implementing cooling measures can effectively reduce the thermal error of the precision boring machine. The experimental data is highly consistent with the simulation results, verifying the accuracy of the thermal-fluid-solid simulation model. Compared with the boring machine slider without a cooling plate, the total thermal error of the boring machine is reduced by 16.67%.

[0262] 5.2 Application of the water jacket in the high-speed spindle system

[0263] Taking the TO water jacket obtained in Section 4.3 as an example, verify the effectiveness of the efficient TO method for the thermal-fluid problems in the cooling components of the high-speed precision boring machine, as Figure 33 shown.

[0264] Due to the extremely complex three-dimensional structure of the TO water jacket, use 3D printing technology to manufacture an optimized TO water jacket according to the dimensions and structure of the original spiral water jacket in the boring machine spindle, as Figure 33 (c) shown. The 3D printing preparation includes slicing the CAD model, setting printing parameters, and calibrating the printer. The actual printing process requires depositing materials layer by layer to build complex structures, including complex cooling channels. Quality inspection is carried out to verify the dimensional accuracy and structural integrity. Then the TO water jacket is embedded in the boring machine spindle system to replace the original spiral water jacket. Then, thermal-fluid-solid behavior analysis and experimental research are carried out. Figure 34 Shows the mesh division results of the boring machine spindle system, all using tetrahedral meshes. The model has a total of 1,989,096 elements and 484,440 nodes. The boundary conditions obtained in [5] are used for the thermal-fluid-solid behavior model. The spindle speed is 10,000 r / min, and the maximum cooling water inlet flow rate is 11 L / min. The room temperature and the coolant inlet temperature are set to 20 °C. The heat generation of the stator, rotor, front bearing, and rear bearing are 830 W, 415 W, 62 W, and 58 W respectively. The forced convection coefficient between the rotating surface and the air is 60 W / (m 2·K), the natural convection coefficient between the housing and the air is 9.7 W / (m 2 ·K). Finally, in [5], the cooling capacity of the TO water jacket obtained by the proposed thermal-fluid TO method was compared with the refrigeration capacity of the traditional spiral water jacket.

[0265] 5.2.1 Heat dissipation performance

[0266] Figure 35 Shows the temperature distributions of the boring machine spindle system with embedded spiral and TO water jackets at Reynolds numbers of 5000 and 15000. Generally speaking, under various cooling conditions, the highest temperature is observed at the stator of the boring machine spindle system. Although the heating power of the rotor is higher than that of the stator, the direct contact between the water jacket and the stator enables the coolant to effectively remove the heat generated by the stator. Since the rotor is in an enclosed space, its internal heat can only be dissipated through heat transfer across the gap between the stator and the rotor and forced convection heat transfer at the rotating end face, resulting in the highest temperature at the stator of the boring machine spindle system. As the Reynolds number increases to 15000, the cooling capacity of the water jacket increases, leading to a decrease in the temperature of the boring machine spindle system. However, the temperature drop near the water jacket of the boring machine spindle system with a TO water jacket is more significant.

[0267] Figure 36 Shows the average temperature of the boring machine spindle system at different inlet Reynolds numbers. As the inlet Reynolds number increases, the heat dissipation performance of the water jacket is improved, and the average temperature of the boring machine spindle system decreases. In addition, the average temperature difference between the spindle systems with spiral water jackets and TO water jackets also increases as the Reynolds number increases because the coolant flow in the TO cooling channels is more complex than that in the spiral cooling channels, making its advantage more obvious under turbulent conditions. Therefore, when the flow rate is high, the cooling performance of the spindle system with a TO water jacket is much better than that of the boring machine spindle system with a spiral water jacket.

[0268] The highest and lowest temperatures of the stator at different Re are as Figure 37 shown. Both the highest and lowest temperatures decrease as Re decreases. For the same inlet temperature, there is no significant difference in the lowest temperature of the stator when using a TO water jacket and when using a spiral water jacket. In addition, due to the more uniform distribution of the TO channels on the outer surface of the stator, the highest temperature of the stator cooled by the TO water jacket is significantly lower than the temperature of the stator cooled by the spiral water jacket.

[0269] 5.2.2 Pressure drop

[0270] Figure 38It shows the pressure distribution of the coolant in the water jacket when the inlet Reynolds number is 15000. Due to the single flow path of the coolant in the spiral cooling channel, there is a significant pumping pressure requirement. The maximum pressure of the coolant at the inlet of the spiral channel is 11872 Pa, and the fluid pressure decreases as the flow distance increases. The coolant pressure in the TO cooling channel is much lower than that in the spiral cooling channel, and the pressure distribution of the coolant in the TO cooling channel is uniform. The maximum pressure at the inlet is 3479 Pa, and the fluid pressure in the main channel is about 1500 Pa.

[0271] The pressure drops at different inlet Reynolds numbers are as Figure 39 shown. The TO cooling channel has a double-inlet and double-outlet structure, so the pressure drop is the sum of the pressures at the two inlets. The results show that as Re increases, the pressure drop of the coolant in the TO cooling channel gradually increases, and the pressure drop difference between the spiral and TO cooling channels also increases. At a Reynolds number of 25000, the pressure drops of the coolant in the spiral and TO cooling channels are 28850 Pa and 9298 Pa respectively. The pressure drop of the coolant in the TO cooling channel is reduced by 67.77%, which not only reduces the pumping pressure requirement of the boring machine spindle cooling system but also reduces the risk of coolant leakage.

[0272] 5.2.3. Thermal deformation

[0273] The thermal deformation of the shaft core in the boring machine spindle system is as Figure 40 shown. The results show that when the spiral water jacket and the TO water jacket are used as cooling elements, the maximum thermal deformations of the shaft core in the boring machine spindle system are 115.19 μm and 97.83 μm respectively. The TO water jacket has much stronger control ability for thermal deformation than the spiral water jacket. Then, the effectiveness of the proposed efficient TO method for the thermal fluid problems in the cooling elements of the precision boring machine spindle is verified.

[0274] The thermal deformations of the spindle at different Re are as Figure 41 shown. The cooling performance of the water jacket improves as the Reynolds number increases, effectively controlling the temperature rise of the boring machine spindle, thereby reducing the thermal deformation of the shaft core. When Re is 5000, when the spiral and TO water jackets are used as the cooling elements of the boring machine spindle, the thermal deformations of the shaft core are 149.94 μm and 126.41 μm respectively, and the thermal deformation is reduced by 15.69%. When Re is 25000, when the spiral and TO cooling water jackets are used as the cooling elements of the boring machine spindle, the thermal deformations of the shaft core are 110.63 μm and 85.33 μm respectively, and the thermal deformation is reduced by 22.87%. However, as Re further increases, the heat dissipation capacity of the cooling water jacket will reach the limit, resulting in the thermal deformation gradually stabilizing. In addition, this shows that the thermal deformation of the boring machine spindle cannot be continuously reduced only by increasing the coolant flow rate. Therefore, it is more critical to improve the heat dissipation and thermal deformation control ability by redesigning the cooling elements.

[0275] 6. Conclusions

[0276] On the OpenFOAM platform, an innovative and efficient topology optimization method is proposed using the finite volume method for the thermal-fluid problems in the cooling elements of high-speed precision boring machines. The linear Darcy model is used to simulate fluid flow, and an adaptive mesh refinement strategy is adopted to enhance the mesh distribution near the fluid-solid boundary. This strategy ensures high accuracy in the physical field calculation while maintaining controllable computational costs. In addition, a parallel computing framework significantly reduces the computational time required for the topology optimization process to achieve efficient topology optimization. Then, a thermal-fluid topology optimization program developed based on the efficient topology optimization method for the thermal-fluid problems of high-speed precision boring machine cooling elements is embedded in the OpenFOAM platform. Then, the thermal-fluid topology optimization program is verified through various case studies, demonstrating the feasibility and effectiveness of the proposed efficient topology optimization method. By comparing the results calculated by OpenFOAM and the finite element method-based COMSOL software, the effectiveness of the finite volume method in thermal-fluid topology optimization is confirmed. Finally, the proposed efficient thermal-fluid topology optimization method is applied to the design of the cooling elements of high-speed precision boring machines. Then, the heat dissipation capacity and efficiency of the topology-optimized cooling elements are compared with those of traditional cooling elements. The following conclusions are drawn:

[0277] (1) The proposed thermal-fluid topology optimization method can be used to design the cooling elements of precision boring machines. By defining the objective function to minimize the temperature and constraining the pressure drop, various configurations of the topology-optimized cooling channels are derived. Adjusting the Lagrange multiplier of the pressure drop helps generate the Pareto front for the cooling channel design, which is helpful for balancing the pumping power and heat dissipation requirements in practical engineering applications. The feasibility of this method is confirmed through the design of the boring machine slider cooling plate and the spindle system cooling water jacket.

[0278] (2) The proposed adaptive mesh refinement strategy effectively improves the manufacturability. Implementing the adaptive mesh refinement strategy can obtain clearer and more manufacturable topology optimization boundaries while maintaining reasonable computational costs. Compared with traditional finite element methods, the computational efficiency of the proposed efficient topology optimization method is significantly improved. Compared with the results obtained by COMSOL Multiphysics based on the finite element method, using the OpenFOAM platform with four processors reduces the computational time by 93.89%, the CPU utilization rate by 41.73%, and the memory usage rate by 82.6%.

[0279] (3) The topology-optimized cooling elements designed using the proposed thermo-fluid topology optimization method can dissipate internal heat and reduce the thermal deformation of the boring machine. For the application of the designed cooling plates in the large base components of the boring machine, the thermal equilibrium time of the boring machine with freely installed cooling plates designed using the proposed thermo-fluid topology optimization method is shorter. Compared with the boring machine slide without a cooling plate, the total thermal error of the boring machine is reduced by 16.67%. For the application of the designed cooling water jacket in the boring machine spindle system, compared with the traditional spiral cooling water jacket, the average temperature of the boring machine spindle with an embedded topology-optimized water jacket is significantly reduced. The pressure drop of the coolant in the topology-optimized water jacket is significantly lower than that in the spiral water jacket. Compared with the shaft core in the spiral water jacket boring machine spindle, the thermal deformation of the shaft core in the boring machine spindle with an embedded topology-optimized water jacket is reduced by 22.87%.

[0280] These findings highlight the effectiveness and efficiency of the proposed thermo-fluid topology optimization method in designing cooling elements for precision machine tools, providing a new method for improving the thermal behavior and machining accuracy of precision boring machines.

[0281] The above-described embodiments are only preferred embodiments given to fully illustrate the present invention, and the protection scope of the present invention is not limited thereto. Equivalent substitutions or transformations made by those skilled in the art on the basis of the present invention are within the protection scope of the present invention. The protection scope of the present invention is subject to the claims.

Claims

1. A topological optimization method for the heat flow problem of the cooling flow path of a cooling element in a high-speed precision machine tool, characterized in that: It includes the following steps: Step 1: To maximize the heat transfer capacity in the cooling channels while ensuring the minimum flow rate requirement is met, construct the objective function and optimization conditions for the topological optimization of the heat flow problem in the cooling channels: Find{γ1,γ2,...,γ i} T Where: Ψ is the objective function; J1 is the volume constraint of the fluid; J2 is the pressure drop constraint of the cooling channel; V is the volume of the fluid region in the design domain; Ω is the design domain; R is the control equation; P in is the inlet pressure; ∫ inlet *dΓ is the integral of the computational domain on the inlet boundary; u is the velocity vector; p is the pressure; T is the temperature; γ is the design density; γ i is the i-th objective to be optimized; d is the design domain variable. Step 2: Solve the objective function 21) Initialization: Decompose the computational domain for parallel computing and initialize the design domain Ω; 22) Determine the objective function Ψ at the current iteration step n and the objective function Ψ at the previous iteration step n-1 and check if the absolute value of their difference is greater than a preset threshold, or check if the current iteration number n is less than or equal to the preset maximum iteration number N: If so, execute step 23); if not, the design optimization is completed; 23) Calculate the state variables by solving the primal problem; 24) Calculate the adjoint variables by solving the adjoint problem; 25) Calculate the objective function and sensitivity; if the current projection slope β is less than the preset maximum projection slope β Max , then update the projection slope β; 26) Update the design variables by the MMA algorithm; If satisfied: Among them: ε Ψ is the target function error threshold; ε V is the volume constraint error threshold; ε p is the pressure drop error threshold; ∫ Ω *dΩ is the integral in the design domain; V frac is the flow channel ratio; p n is the pressure at the current iteration step; p n-1 is the pressure value at the previous iteration step; Then, adopt the adaptive mesh refinement method to update the mesh; 27) If the permeability penalty parameter p κ is less than or equal to the preset maximum value of the permeability penalty parameter p κMax , then update the permeability penalty parameter p κ ; 28) Calculate the material permeability κ and thermal diffusivity D T ; 29) Let n = n + 1 and loop to execute step 22).

2. The topological optimization method for the heat flow problem of the cooling flow path of the cooling element of a high-speed precision machine tool according to claim 1, wherein: In the said step 23), the state variables include the velocity vector u, pressure p, and temperature T; The method for calculating the state variables by solving the primal problem is: Construct the continuity equation for a steady incompressible flow that conserves mass: where: u is the velocity vector; represents the gradient operator; Adopt the Darcy model for turbulence modeling: Where: p is the pressure, and κ and μ represent the fluid permeability and dynamic viscosity respectively; For an incompressible flow, when considering the transport phenomena involving temperature, the description can be extended to include the transfer and diffusion of temperature; The equations governing these processes in an incompressible fluid are: Where: T is the temperature; D T is the thermal diffusivity; The heat source Q(x) depends on the design domain x: Where: Ω is the design domain; Θ is the heat source domain; x ∈ Ω\Θ is the design domain excluding other regions except the heat source.

3. The topology optimization method for the heat flow problem of the cooling flow path of the cooling element of a high-speed precision machine tool according to claim 2, characterized in that: In the said step 24), the adjoint variables include the adjoint velocity u A , the adjoint pressure p A and the adjoint temperature T A ; The method for calculating the adjoint variables by solving the adjoint problem is as follows: Combine the objective function and the constraints into a Lagrangian function: L(u, u A , p, p A , T, T A , λ1, λ2, Ω) = Ψ + ∫ Ω u A ·R u dΩ + ∫ Ω p A R p dΩ + ∫ Ω T A R T dΩ + λ1J1 + λ2J2 where: u A , p A , T A , λ1 and λ2 represent the adjoint velocity, adjoint pressure, adjoint temperature, the Lagrange multiplier associated with the volume constraint, and the Lagrange multiplier associated with the pressure drop constraint, respectively; u A , p A , T A , λ1 and λ2 are all Lagrange multipliers, introduced as dummy variables; Obtain the adjoint model as: Due to the governing equations and the constraint conditions of J1 and J2 are theoretically constant, so L = Ψ; then the sensitivity of the objective function is The following relationship is obtained: Find the appropriate adjoint variables u A , p A , T A , λ1 and λ2 such that Then the final Lagrangian equation is constructed as follows: Wherein: is the control equation; To obtain the optimal topological optimization result, the following conditions should be met: where: δ ζ Ψ, δ ζ J and are the differentials of the objective function Ψ , constraint conditions J and governing equations i with respect to the change of the variable ξ; ξ is a variable and ξ ∈ (u, p, T); Obtain the differential form: where: ∫ Ω *dΩ and ∫ Γ *dΓ denote the integral in the design domain and the integral on the boundary, respectively; For any δu, δp, and δT, all meet the computational requirements within the design domain, and obtain the adjoint equation as follows: Consider the adjoint boundary conditions: For the inlet and the wall, the velocity and temperature are fixed values, then δu = 0 and δT = 0, and the adjoint velocity and adjoint temperature of the inlet and the wall are: T A ＝0 For the outlet, the pressure is 0 and the temperature has a zero gradient, so δp = 0 and The adjoint velocity and adjoint temperature at the outlet are: p A ＝0 Wherein: is the accompanying normal velocity.

4. The topological optimization method for the heat flow problem of the cooling flow path of the cooling element of a high-speed precision machine tool according to claim 3, characterized in that: In the said step 25), the sensitivity calculation method for the heat fluid topological optimization problem is: where: δ γ is the volume factor; κ s and κ f represent the solid permeability and the fluid permeability, respectively; p κ is the penalty parameter for permeability; is the penalty parameter for thermal diffusivity; and represent the solid thermal diffusivity and the fluid thermal diffusivity, respectively.

5. The topological optimization method for the heat flow problem of the cooling flow path of the cooling element of the high-speed precision machine tool according to claim 1, characterized in that: In the said step 28), the calculation method for the permeability κ is: Where: κ is the permeability; κ s and κ f represent the solid permeability and the fluid permeability, respectively; p κ is the penalty parameter of the permeability; Thermal diffusivity D T is calculated as follows: Wherein: is the penalty parameter of the thermal diffusivity; and respectively represent the solid thermal diffusivity and the fluid thermal diffusivity.

Images