Three-dimensional simulation design method and system for light-weight glass bottle production digital mold

Abstract

The invention relates to the technical field of glass bottle production, in particular to a three-dimensional simulation design method and system for a light-weight glass bottle production digital mold. The method comprises the following steps: constructing a digital mold geometric model for a glass bottle mold; a bubble escape priority path is identified according to the digital mold geometric model, and a glass bottle digital simulation model is constructed; executing molten glass flow-heat transfer-bubble evolution coupling simulation based on the glass bottle digital simulation model to extract a potential bubble nucleation and aggregation area; calculating a bubble defect tendency value of the potential bubble nucleation and aggregation area, and performing glass bottle production defect space mapping to generate a glass bottle production defect risk area; differentiated defect risk optimization is carried out on the glass bottle production defect risk area, so that the production design requirement of the light-weight glass bottle is met. According to the invention, through bubble escape path identification and three-dimensional coupling simulation, accurate positioning and differential optimization of the bubble defect risk of the glass bottle are realized.

Description

Technical Field

[0001] The present invention relates to the technical field of glass bottle production, and in particular to a three-dimensional simulation design method and system for a digital mold for lightweight glass bottle production. Background Art

[0002] One of the most serious and common quality issues in glass bottle production is bubble defects. The main reason for bubble defects is that the gas dissolution and escape mechanism is not effectively controlled during the high-temperature molten glass molding process. In traditional glass bottle molding processes, as the molten glass flows at high temperatures, the dissolved gases (mainly oxygen, carbon dioxide, and water vapor) within it need to escape through specific paths. However, due to unreasonable mold design or imprecise temperature control, the gas escape path is often blocked, resulting in the gas being trapped in the glass body and forming bubbles that are difficult to eliminate. Especially in the production of thin-walled, lightweight glass bottles, the reduction in material thickness further increases the difficulty of gas escape, making the bubble problem more prominent. Traditional mold design methods cannot meet the production requirements of high precision and low defect rate. Currently, glass bottle mold design mainly relies on accumulated experience and simplified two-dimensional flow analysis, lacking an accurate description of the complex three-dimensional flow state of molten glass within the mold cavity. Existing technologies generally use single-phase flow simulation methods, treating molten glass as a continuous medium and ignoring the dissolution, precipitation, and migration behavior of gases in the glass melt. As a result, the designed mold cannot provide a reasonable gas escape path. Summary of the Invention

[0003] Based on this, the present invention provides a three-dimensional simulation design method and system for digital molds for lightweight glass bottle production to solve at least one of the above technical problems.

[0004] To achieve the above objectives, a 3D simulation design method for a digital mold for lightweight glass bottle production is provided, comprising the following steps: Step S1: Performing a digital geometric scan on the glass bottle mold to construct a digital mold geometric model; performing multi-scale adaptive meshing of the cavity based on the digital mold geometric model to obtain qualified cavity mesh data; and presetting the bubble nucleation point position on the qualified cavity mesh data to generate bubble nucleation point simulation data; Step S2: performing bubble tracking mesh setting on the qualified cavity mesh data using the bubble nucleation point simulation data to generate bubble tracking mesh data; identifying the preferred bubble escape path based on the bubble tracking mesh data, and constructing a digital simulation model of the glass bottle based on the digital mold geometry model; Step S3: Based on the digital simulation model of the glass bottle, a coupled simulation of molten glass flow, heat transfer, and bubble evolution is performed to extract potential bubble nucleation and aggregation areas; bubble defect propensity values ​​are calculated for the potential bubble nucleation and aggregation areas, and a glass bottle production defect space is mapped to generate glass bottle production defect risk areas; Step S4: Conduct risk attribution analysis on the defect risk areas of glass bottle production and perform differentiated defect risk optimization to achieve lightweight glass bottle production design requirements.

[0005] Preferably, the present invention further provides a three-dimensional simulation design system for a digital mold for producing lightweight glass bottles, which executes the three-dimensional simulation design method for a digital mold for producing lightweight glass bottles as described above. The three-dimensional simulation design system for a digital mold for producing lightweight glass bottles comprises: The digital mesh construction module is used to perform digital geometric scanning of the glass bottle mold to construct a digital mold geometric model; perform multi-scale adaptive meshing of the cavity based on the digital mold geometric model to obtain qualified cavity mesh data; and preset the bubble nucleation point position for the qualified cavity mesh data to generate bubble nucleation point simulation data; The bubble path tracking module is used to perform bubble tracking meshing on qualified cavity mesh data using bubble nucleation point simulation data to generate bubble tracking mesh data; identify the preferred bubble escape path based on the bubble tracking mesh data, and construct a digital simulation model of the glass bottle based on the digital mold geometry model; The defect risk mapping module is used to perform coupled simulations of molten glass flow, heat transfer, and bubble evolution based on a digital simulation model of glass bottles to extract potential bubble nucleation and aggregation areas. The module also calculates bubble defect propensity values ​​for these potential bubble nucleation and aggregation areas, and performs spatial mapping of glass bottle production defects to generate glass bottle production defect risk areas. The lightweight optimization module is used to perform risk attribution analysis on defect risk areas in glass bottle production and conduct differentiated defect risk optimization to achieve lightweight glass bottle production design requirements.

[0006] The present invention achieves accurate capture of microscopic bubble behavior during the glass bottle molding process by introducing a digital mold geometry model and multi-scale adaptive grid technology, fundamentally solving the problem that traditional mold design methods cannot effectively control bubble defects. By combining the preset position of the bubble nucleation point with the bubble tracking grid technology, a complete prediction mechanism for bubble formation and migration is established, allowing designers to intuitively identify the priority path for bubble escape, effectively avoiding the problem of gas being trapped in the glass body. Through the coupled simulation of molten glass flow-heat transfer-bubble evolution, the limitations of traditional single-phase flow simulation are broken through, and the full process simulation of gas dissolution, precipitation and migration behavior in high-temperature molten glass is achieved, which greatly improves the accuracy and scientificity of mold design. In particular, the defect tendency value calculation and spatial mapping analysis of potential bubble nucleation and aggregation areas provide designers with an intuitive defect risk distribution map, giving mold optimization a clear direction and basis. In the production of lightweight glass bottles, through differentiated defect risk optimization strategies, corresponding design adjustments are made according to the bubble risk characteristics of different areas, effectively solving the problem of increased difficulty in gas escape under thin-wall conditions. This design method significantly reduces the bubble defect rate in glass bottle production, improves product qualification rate and production efficiency, and reduces material and energy waste. Furthermore, the digital mold design method reduces reliance on process experience, making mold design more standardized and regularized, shortening new product development cycles and increasing companies' flexibility in responding to market changes. This 3D simulation design method not only solves the bubble defect problem in lightweight glass bottle production but also provides a new technical path for refined and intelligent production in the glass products industry. Therefore, the three-dimensional simulation design method of the digital mold for lightweight glass bottle production of the present invention realizes the precise digitization of the mold through structured light scanning and feature enhancement technology, adopts multi-scale adaptive grid division strategy to improve calculation accuracy, innovatively introduces bubble nucleation point tracking and escape path identification to predict the bubble movement trajectory, constructs a molten glass flow-heat transfer-bubble evolution coupling simulation model to analyze the glass bottle forming process, accurately identifies potential bubble nucleation and aggregation areas through spatial superposition method, combines flow velocity vector processing and solidification front interaction analysis to quantify the bubble defect tendency value, and implements differentiated defect risk optimization based on risk attribution analysis, adds tiny exhaust ducts for poor exhaust areas, optimizes the transition zone geometric curve for areas with excessive flow shear, and adjusts the cooling system parameters for areas with too fast solidification, so as to realize accurate prediction and systematic control of bubble defects in lightweight glass bottles. BRIEF DESCRIPTION OF THE DRAWINGS

[0007] Figure 1 A schematic flow chart of the steps of a three-dimensional simulation design method for a digital mold for producing lightweight glass bottles according to the present invention; The purpose, features and advantages of the present invention will be further described with reference to the accompanying drawings and in conjunction with the embodiments. DETAILED DESCRIPTION

[0008] The following is a clear and complete description of the technical method of the present invention in conjunction with the accompanying drawings. It is obvious that the embodiments described are part of the embodiments of the present invention, but not all of the embodiments. Based on the embodiments of the present invention, all other embodiments obtained by those skilled in the art without making any creative efforts are within the scope of protection of the present invention.

[0009] In addition, the accompanying drawings are merely schematic illustrations of the present invention and are not necessarily drawn to scale. Identical reference numerals in the figures denote identical or similar parts, and thus repetitive descriptions thereof will be omitted. Some of the block diagrams shown in the accompanying drawings are functional entities that do not necessarily correspond to physically or logically separate entities. These functional entities may be implemented in software, in one or more hardware modules or integrated circuits, or in different network and / or processor and / or microcontroller approaches.

[0010] It should be understood that although the terms "first," "second," and the like may be used herein to describe various elements, these elements should not be limited by these terms. These terms are used solely to distinguish one element from another. For example, a first element may be referred to as a second element, and similarly, a second element may be referred to as a first element, without departing from the scope of the exemplary embodiments. The term "and / or" as used herein includes any and all combinations of one or more of the listed associated items.

[0011] To achieve this, please refer to Figure 1 The present invention provides a three-dimensional simulation design method for a digital mold for lightweight glass bottle production, comprising the following steps: Step S1: Performing a digital geometric scan on the glass bottle mold to construct a digital mold geometric model; performing multi-scale adaptive meshing of the cavity based on the digital mold geometric model to obtain qualified cavity mesh data; and presetting the bubble nucleation point position on the qualified cavity mesh data to generate bubble nucleation point simulation data; Step S2: performing bubble tracking mesh setting on the qualified cavity mesh data using the bubble nucleation point simulation data to generate bubble tracking mesh data; identifying the preferred bubble escape path based on the bubble tracking mesh data, and constructing a digital simulation model of the glass bottle based on the digital mold geometry model; Step S3: Based on the digital simulation model of the glass bottle, a coupled simulation of molten glass flow, heat transfer, and bubble evolution is performed to extract potential bubble nucleation and aggregation areas; bubble defect propensity values ​​are calculated for the potential bubble nucleation and aggregation areas, and a glass bottle production defect space is mapped to generate glass bottle production defect risk areas; Step S4: Conduct risk attribution analysis on the defect risk areas of glass bottle production and perform differentiated defect risk optimization to achieve lightweight glass bottle production design requirements.

[0012] In an embodiment of the present invention, the three-dimensional simulation design method for a digital mold for producing lightweight glass bottles includes the following steps: Step S1: Performing a digital geometric scan on the glass bottle mold to construct a digital mold geometric model; performing multi-scale adaptive meshing of the cavity based on the digital mold geometric model to obtain qualified cavity mesh data; and presetting the bubble nucleation point position on the qualified cavity mesh data to generate bubble nucleation point simulation data; In an embodiment of the present invention, a digital geometric scan of a glass bottle mold was performed using a blue structured light scanner with a resolution of 0.05mm and a wavelength of 465nm. The mold assembly, mounted on a rotating platform, was scanned from multiple angles. Each 15-degree rotation of the platform collected 24 angles of data. Three-dimensional coordinates were calculated using a phase unwrapping algorithm to generate point cloud data. De-noising was then performed using a 0.5mm filter radius. Isolated points and outliers exceeding 5mm from the nearest point cluster were removed to obtain clean point cloud data. Three-dimensional surface reconstruction was performed based on this data using a Poisson surface reconstruction algorithm with an octree depth of 10. Features such as the bottle thread, mold parting line, and venting grooves were enhanced. Meshing was performed based on the digital model. Cavity dimensions were first extracted. Thickness transition regions were then marked based on bottle volume, surface area, and wall thickness change rate. Multi-scale meshing was performed using an octree decomposition principle, with a base mesh of 2mm, a fine area of ​​0.5mm, and a thickness transition region of 0.8mm. Five structured boundary layer meshes were used. Perform mesh quality evaluation and optimization to ensure that the Jacobian value is not less than 0.3 and the mesh skewness is less than 0.8, so as to obtain qualified cavity mesh data.

[0013] Step S2: performing bubble tracking mesh setting on the qualified cavity mesh data using the bubble nucleation point simulation data to generate bubble tracking mesh data; identifying the preferred bubble escape path based on the bubble tracking mesh data, and constructing a digital simulation model of the glass bottle based on the digital mold geometry model; In the embodiment of the present invention, when bubble tracking grid setting is performed on qualified cavity grid data using bubble nucleation point simulation data, the gas solubility and diffusion coefficient of soda-lime glass are first obtained. The oxygen solubility coefficient is measured by high temperature differential thermogravimetric analysis. , nitrogen is The diffusion coefficient is determined by capillary method. Using Lagrangian-Eulerian hybrid grid method, 5000 bubble cores are randomly distributed in the cavity, with particle diameter ranging from 0.05mm to 0.5mm, and the flow velocity gradient is greater than The regional grid is encrypted to 0.3mm. Then the exhaust groove features are analyzed, and the feature recognition algorithm is used to extract the position, size and direction of 12 exhaust grooves, each with a width of 0.3mm to 0.5mm and a depth of 0.2mm to 0.3mm. The fast marching method is used to calculate the distance field, and a three-dimensional grid coordinate system is established. The exhaust groove surface is marked as a zero-distance source point, and the Euclidean distance from each grid point in the cavity to the nearest exhaust groove is calculated, and the flow resistance weight factor is introduced for correction. The priority path for bubble escape is identified based on the thermodynamic gradient descent algorithm. Starting from each bubble nucleation point, the bubble movement trajectory is formed along the direction of the steepest potential energy drop. The three exhaust grooves located in the bottle shoulder area undertake 68% of the bubble discharge tasks.

[0014] Step S3: Based on the digital simulation model of the glass bottle, a coupled simulation of molten glass flow, heat transfer, and bubble evolution is performed to extract potential bubble nucleation and aggregation areas; bubble defect propensity values ​​are calculated for the potential bubble nucleation and aggregation areas, and a glass bottle production defect space is mapped to generate glass bottle production defect risk areas; In the embodiment of the present invention, when performing the coupled simulation of molten glass flow-heat transfer-bubble evolution based on the digital simulation model of the glass bottle, the multi-physics field finite volume method is used for numerical solution, and the mass, momentum, and energy conservation equations and the gas phase volume fraction tracking equation are established. The calculation time step is 0.005 seconds and the total calculation time is 8.5 seconds. The temperature field, velocity field, pressure field and gas phase volume fraction field data are collected in 6 key time windows, and the flow field distribution state data is generated by weighted time averaging. The shear strain rate in the cavity is calculated to identify the shear strain rate below The areas with temperature gradient less than The areas with a pressure lower than 10% of the saturated vapor pressure were marked as potential slow solidification areas; the areas with a pressure lower than 10% of the saturated vapor pressure were identified as potential negative pressure precipitation areas. The three types of areas were superimposed through spatial Boolean operations to identify four main bubble nucleation and aggregation areas. Bubble quantitative characteristics of these areas were analyzed to calculate the bubble size and gas volume fraction; the bubble escape efficiency index was calculated based on the flow field data; and the solidification capture risk was calculated through solidification front interaction analysis. The bubble defect tendency value was calculated using a multi-factor weighted method with the following weights: bubble size 0.3, gas volume fraction 0.2, bubble escape efficiency index -0.3, solidification capture risk 0.2. The results were divided into four risk levels and spatially mapped on the digital model.

[0015] Step S4: Conduct risk attribution analysis on the defect risk areas of glass bottle production and perform differentiated defect risk optimization to achieve lightweight glass bottle production design requirements.

[0016] In an embodiment of the present invention, the dominant factor of each risk area is determined through comparative analysis. Among the eight risk areas, areas 1, 4, 6, and 8 are dominated by poor exhaust; areas 2 and 5 are dominated by excessive flow shear; and areas 3 and 7 are dominated by bubble retention caused by rapid solidification. For areas dominated by poor exhaust, the parting surface offset analysis method is used to adjust the parting surface position so that the distance from the risk area to the parting surface is reduced to within 10 mm. At the same time, tiny exhaust channels with a diameter of 0.1-0.3 mm are arranged at the gas gathering position to form an exhaust optimization sub-strategy. For areas dominated by excessive flow shear, the geometry of the neck-shoulder transition zone of the cavity is adjusted, the simple arc is changed to a Bezier curve, the flow channel shrinkage rate is reduced, and the droplet feeding speed is pulsed to form a flow optimization sub-strategy. For areas dominated by rapid solidification, the mold cooling system parameters are adjusted, the cooling water flow rate is reduced, the water temperature is increased, and high-frequency induction heating elements are buried near the risk area. Differential temperature control technology is implemented to form a solidification control sub-strategy. Using the Analytic Hierarchy Process (AHP) to determine the weights of each strategy, conduct conflict detection and collaborative optimization, and ultimately develop a differentiated defect risk optimization strategy that includes 27 specific parameter adjustments. After implementation, the bubble defect propensity value decreased by an average of 42.5%, and the weight of the glass bottles was reduced by 5.8%, achieving the design requirements for lightweight glass bottle production.

[0017] Preferably, the digital geometric scanning of the glass bottle mold in step S1 includes: A structured light scanner is used to perform digital geometric scanning of the glass bottle mold cavity, core rod, and neck mold to obtain point cloud data of the mold assembly surface; The point cloud data of the mold assembly surface is denoised, the filter radius is set to 0.5 mm, isolated points and flying points are deleted, and clean mold point cloud data is obtained; Based on the clean mold point cloud data, three-dimensional surface reconstruction is performed and the mold structural features are enhanced to obtain a digital mold geometric model; among them, the mold structural features include bottle mouth threads, mold parting lines, and exhaust grooves.

[0018] In an embodiment of the present invention, when a structured light scanner is used to perform digital geometric scanning on the glass bottle mold cavity, core rod and neck mold, the mold assembly is first fixed on a rotating platform to ensure the stability of the mold during scanning. A blue light structured light scanner with a resolution of 0.05mm is used, the scanning light source wavelength is 465nm, and the projection fringe density is 100 lines / mm. During the scanning process, the rotating platform rotates 15 degrees each time, and a total of 24 angles of mold surface data are collected to ensure that each surface of the mold is fully exposed to the scanning field of view. For deep holes and recessed areas, multi-angle supplementary scanning is used, and the scanning distance is maintained within the range of 300mm to 500mm. Five fringe images with different phases are collected at each angle, and the three-dimensional coordinates are calculated by the phase unwrapping algorithm to finally generate the surface point cloud data of the mold assembly, and the point cloud density reaches . When denoising the point cloud data of the mold assembly surface, first calculate the average distance of each point to the adjacent points based on statistical principles. Set the filter radius to 0.5mm, and within this radius, calculate the distance variance between the target point and the surrounding points. When the distance variance exceeds twice the global average variance, the point is marked as a noise point. For isolated point groups that obviously deviate from the main structure, they are defined as point sets with a distance greater than 5mm from the nearest point group and are deleted. For flying points, the density analysis method is used to calculate the points with less than 5 points within a sphere with a radius of 0.5mm, and determine them as flying points and delete them. Through the above processing, clean mold point cloud data is obtained, and the point cloud accuracy is better than 0.02mm. When reconstructing three-dimensional surfaces based on clean mold point cloud data, the Poisson surface reconstruction algorithm is used, the octree depth is set to 10, the weight factor is 4, and the solution accuracy is 0.01mm. For areas with complex geometric features, the local point cloud density is increased, and the octree depth is increased to 12. To enhance mold structural features, a feature recognition algorithm was applied to the bottle thread area to extract the thread contours and accurately reconstruct a standard bottle thread with a pitch of 2.5mm and a thread height of 1.2mm. When processing the mold parting line, the point cloud edges on both sides were extracted and fitted with precise planes to ensure that the parting line flatness error was controlled within 0.01mm. When enhancing the vent groove feature, small grooves with a width of 0.5mm and a depth of 0.3mm were identified and the morphological processing algorithm was used to accurately reconstruct the vent groove geometry, ensuring consistency between the digital model and the actual mold geometry.

[0019] Preferably, performing multi-scale adaptive meshing of the cavity according to the digital mold geometric model in step S1 includes the following steps: Extracting the three-dimensional dimensional parameters of the cavity based on the digital mold geometry model; the three-dimensional dimensional parameters of the cavity include the inner diameter of the bottle mouth, the bottom contour, the curvature radius of the neck-shoulder transition zone, the three-dimensional inclination angle of the shoulder, and the wall thickness of the bottle body; The bottle volume, surface area, and wall thickness variation rate are evaluated based on the three-dimensional cavity parameters. When the wall thickness variation rate exceeds 15%, it is marked as a thickness transition area to generate bottle geometric feature data. Perform cavity multi-scale adaptive meshing on the digital mold geometry model using the bottle body geometric feature data to generate initial cavity mesh data; The quality of the initial cavity mesh data is evaluated, and then the cavity mesh is optimized to ensure that the Jacobian value of the minimum mesh unit is not less than 0.3 and the mesh skewness is less than 0.8, so as to obtain qualified cavity mesh data.

[0020] In this embodiment of the present invention, the inner diameter of the bottle mouth is measured using a cross-sectional analysis method. A horizontal cross-section is created at the top of the model, and the circular contour is extracted and its diameter is calculated. The bottom contour is calculated by creating a horizontal cross-section at the bottom of the model and using an edge extraction algorithm to obtain a closed curve shape. For the neck-shoulder transition zone, a curvature analysis tool is used to create 10 equally spaced longitudinal cross-sections along the longitudinal midline. The curvature change is measured at each cross-section to determine the location of the minimum curvature radius. The three-dimensional shoulder inclination angle is determined by creating triangular mesh patches in the shoulder area, calculating the angle between each patch's normal vector and the vertical direction, and taking the average value. The bottle wall thickness is determined by creating 100 evenly distributed measurement points on the model surface. From each point, the distance between the inner and outer surfaces is measured along the normal vector direction to generate a wall thickness distribution map. Gauss's theorem is used to convert the three-dimensional closed surface into a volume integral. The bottle's internal space is divided into 10,000 tiny tetrahedral elements, and the volumes of all elements are accumulated to obtain the total volume value. The bottle's surface area is calculated by discretizing the outer surface into 20,000 triangular patches and summing the areas of all patches to obtain the total surface area. The bottle wall thickness change rate is determined by calculating the ratio of the wall thickness difference between adjacent measuring points to the distance between the two points. The specific calculation formula is: ,in is the wall thickness value of the first point (mm), is the wall thickness value of the second point (mm), is the distance between two points (mm). When the wall thickness variation rate between any two adjacent points exceeds 15%, the area where the line connecting these two points is located is marked as the thickness transition area and highlighted in red on the three-dimensional model. A hierarchical octree partitioning strategy is used to divide the entire model space into initial grid units, with the side length set to 1 / 10 of the maximum size of the bottle. The grid is then encrypted based on the geometric feature data of the bottle: in the bottle thread area, the grid units are subdivided to 0.5 mm; in the thickness transition area (area where the wall thickness variation rate exceeds 15%), the grid units are subdivided to 0.8 mm; in the neck-shoulder transition area (area where the curvature radius is less than 5 mm), the grid units are subdivided to 0.6 mm; in the area where the three-dimensional inclination angle of the shoulder changes by more than 30 degrees, the grid units are subdivided to 1.0 mm; in the thin-walled area where the bottle wall thickness is less than 2 mm, the grid units are subdivided to 1.2 mm; the original grid size is maintained in the remaining areas. Smooth transition zones are set between grid units, and the ratio of adjacent unit sizes does not exceed 1.5 to ensure smooth transition of the grid. After the division is completed, the initial cavity mesh data containing about 1.5 million cells is generated. The following criteria are used: the Jacobian value of each mesh cell is calculated; the mesh cell skewness is calculated using the formula ,in is the maximum angle in the unit (degrees), is the minimum angle (degrees) within the unit. Through global scanning, units with Jacobian values ​​lower than 0.3 (about 8% of the total) and units with skewness greater than 0.8 (about 5% of the total) are identified. For these unqualified units, local mesh optimization is performed: for units with low Jacobian values, the topology optimization method is used to adjust the node connection relationship; for units with large skewness, the spring smoothing method is used to redistribute the node position, and the spring stiffness is inversely proportional to the unit side length; for units with both problems, node redistribution is performed first and then the topology structure is adjusted. After optimization, the mesh quality is evaluated again to ensure that the Jacobian value of all units is not lower than 0.3 and the skewness is less than 0.8, and finally qualified cavity mesh data that meets the simulation accuracy requirements is obtained.

[0021] Preferably, step S2 includes the following steps: Step S21: Obtain the type of glass; set the solubility coefficient and diffusion coefficient of the molten gas in the glass melt according to the type of glass to generate a molten gas characteristic coefficient; Step S22: performing bubble tracking grid setting on the qualified cavity grid data based on the melt gas characteristic coefficient and the bubble nucleation point simulation data to generate bubble tracking grid data; Step S23: Analyze the position, size, and direction of the exhaust groove according to the digital mold geometric model to obtain mold exhaust structure data; Step S24: analyzing the exhaust channel position of the bubble tracking grid data using the mold exhaust structure data, and calculating the shortest distance from each bubble tracking grid data to the exhaust slot to generate distance field distribution data; Step S25: identifying a preferred bubble escape path based on the distance field distribution data; using the preferred bubble escape path and the bubble tracking grid data as glass melt-gas two-phase flow simulation configuration data; Step S26: Obtaining glass bottle production process parameters; Step S27: Import the glass melt-gas two-phase flow simulation configuration data, glass bottle production process parameters, and molten gas characteristic coefficients into the multiphase flow simulation software, set boundary conditions and initial conditions, and build a digital simulation model of the glass bottle based on the digital mold geometry model.

[0022] In the embodiment of the present invention, when obtaining the type of glass, an X-ray fluorescence spectrometer is used to measure the composition of the glass sample, accurately detecting the content of main components such as silicon dioxide, sodium oxide, and calcium oxide. The analysis results show that the glass used in this embodiment is soda-lime glass, in which the silicon dioxide content is 72.5%, the sodium oxide content is 14.2%, the calcium oxide content is 10.1%, and the rest are trace elements such as magnesium oxide and aluminum oxide. According to the type of glass, the solubility coefficient of oxygen in the molten glass is set to , which is obtained by conducting equilibrium solubility tests at 1450°C; the diffusion coefficient of oxygen in the glass melt is set to The coefficient is determined by measuring the migration rate of bubbles in the glass melt at 1450℃ using the rotating disk method. At the same time, the solubility coefficient and diffusion coefficient of other molten gases such as carbon dioxide and sulfur dioxide are set. According to the qualified cavity grid data in the previous step, the temperature gradient is greater than Areas with a pressure gradient greater than 0.2 MPa / cm are marked as potential bubble nucleation zones. Within these areas, the grid cells are further refined to 0.3 mm to improve computational accuracy during the initial bubble formation phase. Subsequently, based on the gas diffusion equation and the temperature field distribution, the main bubble movement paths are predicted. A grid density gradient is set along these paths, with a grid size of 0.3 mm near the nucleation point and gradually increasing to 0.8 mm along the direction of bubble movement. Furthermore, a uniform grid size of 0.5 mm is used in predicted bubble aggregation areas, such as the upper corners of the bottle body and the transition zone between the bottle shoulders. To track the merging behavior of tiny bubbles, connecting grid strips with a cell size of 0.2 mm are created in areas where the distance between adjacent bubbles is less than 1 mm. The mold parting surface contour is extracted, marking the boundary where the two halves of the mold meet. Inspection points are set every 5 mm along the parting surface contour, for a total of 86 inspection points. At each inspection point, a depth analysis is performed perpendicular to the parting surface into the mold interior to examine the groove structure. Linear structures with a groove width of less than 1 mm, a depth between 0.5 and 1.2 mm, and a length greater than 5 mm are identified as venting grooves. The center lines of all venting grooves are extracted, and the three-dimensional coordinate sequence of each center line is recorded. The width change of each venting groove is measured, and the narrowest and widest dimensions are accurately recorded, which are 0.6 mm and 0.9 mm, respectively. The depth distribution of the venting grooves is measured, and the shallowest and deepest values ​​are recorded, which are 0.5 mm and 1.1 mm, respectively. The strike angle of each venting groove is determined, that is, the angle between the center line of the venting groove and the horizontal plane, with an angle range of 15 degrees to 75 degrees. All measurement data are integrated to form structured mold venting structure data, which includes the spatial distribution, geometric dimensions and strike information of the venting grooves. The fast marching method is used to calculate the Euclidean distance from each grid node to the nearest venting groove. The specific implementation method is as follows: first, the exhaust slot surface is discretized into 5000 reference points and a kd-tree spatial index structure is constructed; then, each node in the bubble tracking grid is traversed and the nearest neighbor search algorithm is used to find the exhaust slot reference point closest to the node, and the Euclidean distance D between the two points is calculated using the following formula: ,in are the grid node coordinates (mm), is the coordinate of the exhaust groove reference point (mm). In the case of multiple exhaust grooves, the minimum distance value is taken as the exhaust distance of the node. The exhaust distance data of all nodes constitute a scalar field in three-dimensional space. The distance value on the node is interpolated to the center of the grid unit through trilinear interpolation to form a continuous distance field distribution. The distance field data is normalized to generate dimensionless distance field distribution data with a value range of 0 to 1, where 0 represents the location on the exhaust groove surface and 1 represents the farthest position from the exhaust groove. At each position marked as a bubble nucleation point, the negative gradient vector of the distance field is calculated, which is expressed as: ,in is the distance field function (dimensionless), and G is the vector pointing in the direction of the fastest distance decrease (1 / mm). The bubble's path is traced along the negative gradient direction, with each step length of 0.2 mm, until the distance field value reaches less than 0.05 (i.e., very close to the exhaust groove) or reaches the outer surface of the glass bottle. The drag coefficient is calculated for each traced path, taking into account the glass viscosity, temperature distribution, and pressure gradient. The drag coefficient is expressed as: ,in is the glass viscosity (Pa·s), is the path length (mm), is the equivalent diameter of the path pipe (mm), is the pressure correction factor (1 / Pa), is the pressure gradient (Pa / mm). The path with the smallest drag coefficient is marked as the preferred bubble escape path. All preferred paths are integrated with the bubble tracking grid data to generate glass melt-gas two-phase flow simulation configuration data including the geometric grid, bubble nucleation location, and escape path. The glass melting temperature is recorded as 1450°C, which is obtained by measuring with a high-temperature thermocouple installed at the furnace outlet. The glass droplet temperature is measured at 1200°C, which is measured by an infrared thermometer when the droplet is cut off. The glass droplet weight is measured to be 355 grams, which is determined by weighing 30 consecutive droplets using a precision balance and taking the average value. The mold preheat temperature is recorded as 550°C, which is obtained by measuring with a thermocouple embedded in the mold wall. The compressed air pressure is measured to be 0.6 MPa, which is measured on the compressed air pipe by a pressure sensor. The blowing time parameters are recorded: initial blowing time 0.8 seconds, re-blowing time 1.2 seconds, and the interval between two blows 0.5 seconds, which are recorded by the timer of the production line control system. Mold opening and closing times were measured: mold closing time was 0.3 seconds, mold opening time was 0.4 seconds, and mold movement time was measured using a displacement sensor. Cooling parameters were recorded: mold cooling water temperature was 35°C, flow rate was 15 liters / minute, measured using a temperature sensor and flowmeter. All parameters were recorded by the production line data acquisition system, forming a complete data set of glass bottle production process parameters. Boundary conditions were set: a mass inflow boundary was set at the glass feed inlet, with a flow rate of 355 grams / 1.5 seconds; a temperature boundary was set at the mold wall, with a temperature of 550°C and a thermal conductivity of 25 W / (m·K); and a pressure outlet boundary was set at the venting slot, with a pressure of 0.1 MPa. Initial conditions were set: the initial glass temperature was 1200°C, the initial pressure was 0.6 MPa, and the initial bubble volume fraction was 0.001. Based on the digital mold geometry model, a flow field domain and temperature field computational grid were constructed, and the bubble tracking grid was coupled to the flow field grid. Set the material properties: glass density is 2500 kg / m³, specific heat capacity is 1200 J / (kg·K), thermal conductivity is 1.5 W / (m·K), and the viscosity-temperature relationship follows the VFT equation: , where μ is the glass viscosity (Pa·s), is the temperature (Kelvin), are material constants, which are -2.5, 4500, and 250, respectively. Set gas properties and interaction parameters, including surface tension, contact angle, and bubble coalescence and breakup criteria. With these settings, the digital simulation model of the glass bottle is constructed.

[0023] Preferably, identifying the preferred bubble escape path according to the distance field distribution data includes: Bubble grid node analysis is performed based on bubble tracking grid data, and the gradient vector of each bubble grid node is calculated based on the distance field distribution data to generate flow gradient field data; Based on the digital mold geometry model, the streamline path from each point in the mold cavity to the exhaust groove is tracked using the flow gradient field data to obtain the complete streamline path data; wherein, the streamline integration step size is set to 0.1mm; Calculate the path length and turning angle of each streamline based on the complete streamline path data to obtain streamline characteristic analysis data; Based on the streamline feature analysis data, the top 20% streamlines with the shortest path length and turning angle less than 45 degrees are selected as the main escape channels to obtain the main channel streamline data; Based on the main channel streamline data, the flow channel cross-sectional area changes were analyzed, and the flow channel shrinkage rate in the neck-shoulder transition area was calculated. When the shrinkage rate exceeded 40%, it was marked as a flow resistance area. The bubble escape path is selected based on the main channel streamline data through the flow resistance area to generate the bubble escape priority path.

[0024] In this embodiment of the present invention, the coordinates of all nodes in the grid data are extracted, totaling approximately 350,000 nodes. For each node, the gradient vector of the distance field distribution data is calculated using the central difference method. The calculation formula is: , where D is the distance field function value (dimensionless), is the spatial coordinate (mm). In the specific implementation process, six adjacent points around the node are taken, respectively located at The positive and negative sides of the three directions are spaced 0.2 mm apart, and the gradient component of the central node is calculated using the distance field values ​​of these six points: , and The same method is used for calculation. For the boundary nodes, the gradient is calculated using a one-sided difference. The calculated gradient vector is normalized to obtain a unit direction vector, which represents the local preferred direction of bubble movement. 5000 starting points are evenly arranged in the mold cavity, with a starting point spacing of 1.5 mm. Starting from each starting point, the streamline integration step is set to 0.1 mm, and the streamline is traced along the gradient direction. The integral calculation formula for each step is: ,in is the current position coordinate vector (mm), k is the step size (0.1 mm), is the normalized gradient vector (dimensionless) at this point. The fourth-order Runge-Kutta method improves the integration accuracy by calculating four intermediate points: , , , , . The integration process continues until the streamline reaches the exhaust slot surface (the distance field value is less than 0.05) or the preset maximum number of integration steps of 500 is reached. The coordinates of all integration points of each streamline are recorded in sequence to form complete streamline path data, including the starting point position, the end point position and the three-dimensional coordinates of all tracking points in the middle. When calculating the streamline characteristics based on the complete streamline path data, two key indicators are analyzed for each streamline: path length and turning angle. The path length is calculated by accumulating the Euclidean distance between adjacent points on the streamline, and the formula is: ,in is the coordinate of the i-th point on the streamline (in millimeters), and n is the number of points on the streamline. The steering angle is calculated by calculating the angle formed by three adjacent points on the streamline and summing them up. The calculation formula is: ,in is the direction vector of the i-th line segment, , is the total turning angle of the streamline (degrees). To refine the analysis, each streamline is divided into five segments, and the local turning angle and length of each segment are calculated. At the same time, the average curvature of the streamline is calculated using the formula: , in degrees / mm. All calculation results are summarized to form streamline feature analysis data, including the total length, total turning angle, segment length, segment turning angle and average curvature value of each streamline. 5000 streamlines are sorted in ascending order by path length. The average value (35.8 mm) and standard deviation (12.3 mm) of all streamline lengths are calculated, and abnormal streamlines with a length exceeding the average value plus twice the standard deviation (60.4 mm) are eliminated. For the remaining streamlines, streamlines with a turning angle of less than 45 degrees are further screened, and a total of 3280 streamlines that meet the conditions are obtained. These streamlines are arranged in ascending order by path length, and the top 20% (656) are selected as candidates for the main escape channels. To ensure uniformity of spatial distribution, the mold cavity space is divided into 8 quadrants, and streamlines of the same proportion are selected in each quadrant to avoid excessive concentration of streamlines in local areas. A spatial clustering analysis was performed on the selected streamlines. The density clustering algorithm DBSCAN was used. The cluster radius was set to 3 mm and the minimum number of points was 5. The streamlines with close spatial positions were classified as the same escape channel. Finally, 22 main escape channels were obtained, containing a total of 656 streamlines, forming the main channel streamline data. 20 cross sections were evenly sampled on each main channel streamline, and the cross-sectional plane was perpendicular to the tangent direction of the streamline. On each cross section, the equivalent flow channel diameter was calculated using the fluid mechanics cross-sectional analysis method. The specific implementation method is: with the streamline point as the center, radiate in 16 directions around, and advance in each direction along the direction of the steepest distance field gradient until encountering a wall or the distance field gradient changes significantly (the rate of change exceeds 200%), measure the length of the 16 radiation directions and take the average value as the equivalent radius. (mm), flow channel cross-sectional area (square millimeters). Calculate the rate of change of cross-sectional area along the streamline direction: Focus on analyzing the cross-sectional area change of the neck-shoulder transition zone (20%-40% of the streamline length) and calculate the shrinkage rate: When the shrinkage rate exceeds 40%, the area is marked as a flow resistance area and highlighted in yellow in the 3D model. Analysis shows that among the 22 main escape channels, 8 channels have significant flow resistance areas, with shrinkage rates of 42.3%, 45.7%, 49.2%, 51.8%, 53.5%, 56.1%, 58.4% and 62.9% respectively. A weight score is assigned to each main channel streamline, and the calculation formula is: , where L is the normalized streamline length (dimensionless, 0-1 range), θ is the normalized turning angle (dimensionless, 0-1 range), CR is the maximum shrinkage (percentage), 、 、 are weight coefficients, which are set to 20, 30, and 50 respectively. For streamlines with no significant flow resistance (shrinkage rate less than 40%), The value of the item is 0. All main channel streamlines are arranged in descending order according to the weight scores, and the top 30% (197) are taken as the priority paths for bubble escape. These paths are mainly distributed in the upper area of ​​the bottle body, with an average length of 22.6 mm, an average turning angle of 28.3 degrees, and an average shrinkage rate of 22.7%. In order to improve the efficiency of bubble removal, these priority paths are clustered into 5 main bubble escape channels according to their spatial positions, and the entrance position of each channel is marked to form a bubble escape priority area distribution map to guide subsequent mold structure optimization. The final generated bubble escape priority path data will be used for the glass melt-gas two-phase flow simulation configuration.

[0025] Preferably, in step S3, performing coupled simulation of molten glass flow-heat transfer-bubble evolution based on the digital simulation model of the glass bottle includes: Based on the digital simulation model of the glass bottle, a coupled simulation of molten glass flow, heat transfer, and bubble evolution is performed. The temperature field, velocity field, pressure field, and gas volume fraction field are collected according to the preset cycle time window to obtain mold flow field and bubble simulation test data. Perform time averaging on the mold flow field-bubble simulation test data to generate flow field distribution state data; The shear strain rate in the cavity is calculated based on the velocity field in the flow field distribution data. Then, the low-speed flow area with a shear strain rate lower than 5 s⁻¹ is identified and marked as potential flow stagnation area data. According to the temperature field in the flow field distribution data, the slow cooling area with a temperature gradient less than 2 °C / mm is identified and marked as potential slow solidification area data; According to the pressure field in the flow field distribution data, the area below 10% of the saturated vapor pressure of the glass melt is identified and marked as the potential negative pressure precipitation area data; The data of potential flow stagnation zone, potential slow solidification zone and potential negative pressure precipitation zone are spatially superimposed to obtain the potential bubble nucleation and aggregation area.

[0026] In the embodiment of the present invention, the finite volume method is used to discretely solve the control equations. The flow equations use the mass conservation equation and the momentum conservation equation: and ,in is the density of glass (kg / m³), is the velocity vector (m / s), is the pressure (Pa), is the viscous stress tensor (Pa), is the acceleration due to gravity (m / s²). The heat transfer equation uses the energy conservation equation: , where C p is the specific heat capacity (J / (kg·K)), is the temperature (K), λ is the thermal conductivity (W / (m·K)), is the viscous dissipation term (W / m³). The bubble evolution adopts the population balance equation: ,in is the bubble number density distribution function , is the bubble radius (m), is the bubble growth rate (m / s), is the source term . Set the time step to 0.001 seconds, the total simulation time to 10 seconds, and collect data every 0.1 seconds as a periodic time window. Record the temperature, velocity, pressure and gas volume fraction through 25 evenly distributed monitoring points to generate mold flow field-bubble simulation test data containing 4000 sets of time-space data points. Divide the simulation process into 100 equal time intervals, each interval is 0.1 seconds. For each grid node, calculate the weighted average value of each physical quantity at each time point, and the weight function is an exponential decay function ,in is the current time (seconds), is the evaluation time point (seconds), is the characteristic time scale (0.5 seconds). The specific calculation formula is: ,in is a physical quantity (temperature, velocity, pressure or gas volume fraction), is the time average. For local rapid change areas, an adaptive time window is used, and the window width decreases as the gradient of the physical quantity increases to ensure that transient phenomena are captured. For the entire flow field, the time average of the velocity field, the time average of the temperature field, the time average of the pressure field, and the time average of the gas phase volume fraction field are calculated, where x, y, and z are spatial coordinates (mm). All time-averaged data are integrated to form flow field distribution state data. The velocity gradient tensor is calculated using the central difference method. The velocity gradient tensor is defined as: ,in For speed i Directional component (m / s), for j The shear strain rate tensor is defined as the symmetric part of the velocity gradient tensor: ,in for L The transpose of . The shear strain rate scalar value calculation formula is: , where ":" represents the tensor double dot product, , The unit is . In actual calculations, for the center point of each unit in the three-dimensional grid, the velocity gradients in all directions are calculated through the velocity values ​​of the six adjacent points around this point, and then the shear strain rate is calculated. After calculating the shear strain rate distribution of the entire flow field, the areas with shear strain rate values ​​lower than 5 s⁻¹ are identified. In these areas, the glass melt flows slowly and stagnation areas are easily formed. Adjacent low shear strain rate points are connected into regions. When the volume of the continuous region exceeds 0.5 cubic centimeters, it is marked as a potential flow stagnation area, and the potential flow stagnation area data is generated, which includes the area position, shape and volume information. Calculate the temperature gradient vector: ,in is the temperature ( ), is the spatial coordinate (mm). The temperature gradient amplitude calculation formula is: , the unit is / mm. The fourth-order precision finite difference scheme is used in the calculation: , where h is the grid spacing (mm), and A similar formula is used for calculation. The temperature gradient distribution is calculated for the entire temperature field, and areas with a temperature gradient amplitude of less than 2°C / mm are identified. These areas have a slow cooling rate and the glass solidification time is prolonged. Adjacent low temperature gradient points are combined into continuous areas. When the volume of the continuous area exceeds 0.8 cubic centimeters, it is marked as a potential slow solidification area. Feature extraction is performed on all marked areas, and the center position, shape parameters (length-to-width ratio) and volume of the area are recorded to generate potential slow solidification area data containing 15 major slow cooling areas. Determine the saturated vapor pressure of the glass melt at local temperature. According to experimental measurements, the saturated vapor pressure of the soda-lime glass used in the working temperature range (800°C-1200°C) is The relationship between (Pa) and temperature T (Kelvin) is: ,in =10.2, =8500. For each grid point in the model, the corresponding saturated vapor pressure is calculated according to the temperature of the point, and then the actual pressure is converted to and saturated vapor pressure The ratio is calculated as: Identification The area where the value is lower than 90%, that is, the area where the pressure is lower than 10% of the saturated vapor pressure, is prone to gas precipitation and bubble formation. = 90% is used as the threshold to extract the isosurface, and the inner area of ​​the isosurface is marked as the potential negative pressure precipitation area. The volume, surface area and the position of the lowest pressure point of each negative pressure precipitation area are calculated. The area with a volume greater than 0.3 cubic centimeters is regarded as the main negative pressure precipitation area, and the potential negative pressure precipitation area data including the regional geometric characteristics and pressure distribution is generated. The mold cavity space is discretized into cubic voxels with a side length of 0.2 mm, totaling about For each voxel, check whether it belongs to three potential problem areas at the same time and assign different weight coefficients: flow stagnation area weight =0.4, weight of slow solidification zone =0.3, weight of negative pressure precipitation zone = 0.3. Calculate the comprehensive risk factor for each voxel: ,in is an indicator function, which takes the value 1 when the voxel belongs to the corresponding region, otherwise it takes the value 0. The risk coefficient R ranges from 0 to 1. When the value is ≥0.7, the voxel is marked as a potential bubble nucleation and aggregation area. A three-dimensional morphological closing operation (the structural element is a sphere with a radius of 1 mm) is used to connect the marked areas, eliminate isolated points and small holes, and make the area more continuous. The boundary surface of the continuous area is extracted, and the volume, surface area, center of mass position and shape factor of each area are calculated. The area with a volume greater than 0.5 cubic centimeters is identified as the main bubble nucleation and aggregation area. A total of 8 main areas are identified, with a total volume accounting for 4.2% of the mold cavity volume, mainly distributed in the bottle shoulder transition area and the bottle bottom edge area.

[0027] Preferably, in step S3, calculating the bubble defect tendency value for the potential bubble nucleation and aggregation area and performing glass bottle production defect space mapping includes: Based on the gas phase volume fraction field in the flow field distribution state data, regional bubble characteristics of potential bubble nucleation and aggregation areas are quantified to generate regional bubble quantitative characteristic data; wherein the regional bubble quantitative characteristic data includes bubble size and gas phase volume fraction; The velocity field in the flow field distribution state data is used to process the local velocity vector of the potential bubble nucleation and aggregation area, and the bubble escape efficiency index is calculated based on the bubble escape priority path. Conducting solidification front interaction analysis on potential bubble nucleation and accumulation areas to obtain solidification capture risk data; The weighted defect tendency value of the potential bubble nucleation and aggregation area is calculated based on the bubble size, gas phase volume fraction, bubble escape efficiency index and solidification capture risk data to generate the bubble defect tendency value. The bubble size weight is set to 0.3, the gas phase volume fraction weight is set to 0.2, the bubble escape efficiency index weight is set to -0.3 and the solidification capture risk data weight is set to 0.2. The bubble risk attention level is divided according to the bubble defect propensity value, and the glass bottle production defect space is mapped according to the corresponding potential bubble nucleation and aggregation areas to generate the glass bottle production defect risk area.

[0028] In the embodiment of the present invention, each potential bubble nucleation and aggregation area is numbered independently, and a total of 8 main areas are marked. For each area, the three-dimensional distribution data of the gas phase volume fraction α value is extracted. The density clustering algorithm DBSCAN is used to Grid points with values ​​greater than 0.001 are clustered to identify independent bubbles, with a cluster radius of 0.5 mm and a minimum number of points of 10. For each identified bubble, its equivalent diameter is calculated: , where V is the bubble volume (cubic millimeters), obtained by integrating the gas phase volume within the bubble: . Count the number density n (per cubic centimeter) and average size of bubbles in each area (mm), size distribution standard deviation (mm), and maximum bubble size (mm). Calculate the regional average gas volume fraction: ,in For the The gas phase volume fraction of each grid cell, For this unit volume, generate a bubble size probability density function for each region and record the 10%, 50%, and 90% quantile bubble size values. Extract the velocity field data within each potential bubble nucleation and aggregation region. . Calculate the regional average velocity vector: ,in is the velocity vector of the ith grid cell (m / s), is the unit volume (cubic meters). Calculate the flow velocity direction consistency index: , with a value range of 0 to 1. A larger value indicates a more consistent flow direction. Subsequently, the dot product operation is performed on the velocity vector in each region and the direction of the bubble escape priority path to calculate the direction consistency: ,in is the unit direction vector of the bubble escape priority path at this point. Calculate the average direction consistency of the region: The bubble escape efficiency index is calculated based on the flow velocity, directional consistency, and directional consistency with the escape path: ,in is the distance from the area to the nearest exhaust slot (mm), is the reference distance (10 mm). The value range is 0 to , the larger the value, the stronger the bubble escape ability. The bubble escape efficiency index was calculated for 8 potential bubble nucleation and aggregation areas, and the results were 0.023, 0.045, 0.067, 0.031, 0.052, 0.018, 0.039, and 0.026 m / s, respectively. For the soda-lime glass used, the solidification temperature range was 720°C to 680°C, and the 700°C isothermal surface was used as the solidification front. The advancement process of the solidification front with time was tracked, and the solidification front positions at 100 time points were extracted with a time interval of 0.1 second. Calculation of the solidification front advancement speed ,in is the displacement distance of the solidification front at adjacent time points (mm), is the time interval (0.1 seconds). For the bubble rising speed, Stokes' law is applied to calculate: , where g is the acceleration due to gravity (9.8 m / s²), is the bubble radius (m), is the density of liquid glass (2500 kg / m3), is the gas density in the bubble (1.2 kg / m3), is the viscosity of the glass melt (Pa·s), and the viscosity changes with temperature in the following manner equation: For each potential bubble nucleation and accumulation area, the risk of solidification capture is calculated: ,in For size The number of bubbles in the region is N, and the summation condition is that the solidification rate of the bubbles is greater than the rising rate of the bubbles. The F value ranges from 0 to 1, and the larger the value, the higher the risk of the bubbles being captured by the solidification front. Each parameter is normalized so that its value range is unified between 0 and 1. The formula for bubble size normalization is: ,in is the regional average bubble diameter (mm), is the maximum bubble diameter in all regions (1.8 mm). The normalized formula for gas phase volume fraction is: ,in is the regional average gas volume fraction, is the maximum gas phase volume fraction in all regions (0.025). The normalized formula for the bubble escape efficiency index is: , where E is the regional bubble escape efficiency index (m / s), The maximum efficiency index in all areas (0.067 m / s). The solidification capture risk is already in the range of 0 to 1, so no normalization is required and the F value is used directly. Based on the set weight coefficient, the bubble defect tendency value is calculated: , where 0.2 is added as the basic risk value at the end to ensure that the defect tendency value is within a meaningful range. The theoretical range of the D value is from 0 to 1, and the larger the value, the higher the tendency to form bubble defects. When classifying the bubble risk attention levels and mapping the defect space of glass bottle production according to the bubble defect tendency value, a four-level risk classification standard is adopted: D ≤ 0.4 is low risk (green), 0.4 < D ≤ 0.55 is medium risk (yellow), 0.55 < D ≤ 0.7 is high risk (orange), and D > 0.7 is extremely high risk (red). According to this standard, the risk levels of the 8 potential bubble nucleation and aggregation regions are as follows: Region 1 is medium risk (D = 0.51), Region 2 is high risk (D = 0.67), Region 3 is low risk (D = 0.38), Region 4 is extremely high risk (D = 0.72), Region 5 is medium risk (D = 0.55), Region 6 is high risk (D = 0.68), Region 7 is medium risk (D = 0.49), and Region 8 is high risk (D = 0.62). Using three-dimensional rendering technology, the regions with different risk levels are marked with corresponding colors on the three-dimensional model of the glass bottle to generate a visual model of the glass bottle defect risk region. Combining the geometric structure characteristics of the glass bottle, the distribution laws of the extremely high risk and high risk regions are analyzed, and it is found that these regions are mainly concentrated in the shoulder transition area of the bottle (accounting for 42% of the total high risk regions), the bottom corner of the bottle (accounting for 28%), and the uneven wall thickness area of the bottle body (accounting for 30%). According to the spatial distribution of the risk regions, the positions of the mold structure that need to be optimized are determined, and a dataset of the glass bottle production defect risk regions containing detailed spatial coordinates and risk levels is generated.

[0029] Especially important is to perform local flow velocity vector processing on the potential bubble nucleation and aggregation regions through the velocity field in the flow field distribution state data, and calculate the bubble escape efficiency index according to the preferential path of bubble escape, specifically as follows: Locate the grid cells of the velocity field in the flow field distribution state data in each potential bubble nucleation and aggregation region, and extract the three-dimensional flow velocity vector components in the corresponding region to obtain regionalized local flow velocity vector field data; Based on the preferential path of bubble escape, identify the path segments for each potential bubble nucleation and aggregation region, and calculate the average tangent direction of the path segment in the region to obtain the regional escape path reference direction vector; Perform vector dot product operations on the regionalized local flow velocity vector field data and the corresponding regional escape path reference direction vector, and extract the projection length in the escape path reference direction to obtain the driving flow velocity component of the bubble along the preferential path; Perform statistical averaging on the driving flow velocity components of the bubbles along the preferential path in each potential bubble nucleation and aggregation region to obtain the regional average escape driving flow velocity; The regional average escape driving velocity is compared with the preset critical escape rate threshold and normalized to generate the bubble escape efficiency index of each potential bubble nucleation and aggregation area.

[0030] In the embodiment of the present invention, the octree space partitioning algorithm is used to divide the entire computational domain into 2048 sub-regions with a grid resolution of 0.25 mm. Based on the spatial coordinate range of the potential bubble nucleation and aggregation region obtained in the above steps, a regional bounding box is established. For each bounding box, a point-polyhedron inclusion relationship check is performed to identify the grid cells that are completely located in the region and generate an index list of the grid cells in the region. Taking the inner corner of the bottleneck (region 3) as an example, it contains 12635 grid cells. According to the index list, the velocity vector data of the corresponding cell is extracted from the global velocity field database. The trilinear interpolation method is used during extraction to ensure the accuracy of coordinate transformation. Taking into account the regional boundary effect, for grid cells less than 0.5 mm from the boundary, the weighted average method is used to calculate the velocity value, and the weight factor ,in is the shortest distance to the boundary. The extracted velocity vectors are statistically analyzed to calculate the average velocity, standard deviation, maximum velocity, and minimum velocity in each area. The extracted complete three-dimensional velocity vector data are classified and stored according to the area number to form regionalized local velocity vector field data. When identifying the path segment of each potential bubble nucleation and aggregation area based on the bubble escape priority path, the spatial intersection relationship between the area and the bubble escape priority path is first determined. Calculate the entry point of each streamline in the bubble escape priority path through the potential bubble nucleation and aggregation area. and exit points , the entry and exit point sets are respectively recorded as and For each region, calculate the center point ,in is the coordinate of all grid nodes in the area, and n is the total number of nodes. right Sort the points in the image and select the m entry points closest to the center point (m is the square root of the area volume divided by 2.5, in mm³) to form a representative entry point set. For each entry point , find the exit point associated with it , forming a local path segment that passes through the area. The geometric shape of each path segment is reconstructed using cubic spline interpolation. 20 equidistant points are sampled on the path segment, and the tangent vector between each two adjacent points is calculated. The tangent vector of each path segment is weighted averaged, and the weight coefficient is proportional to the flow capacity of the bubble escape priority path. The formula is: ,in is the flow capacity value of the corresponding streamline. The reference direction vector of the regional escape path is obtained by comprehensive analysis. , as the optimal direction for bubbles to escape from the area. The projection decomposition method is used to perform vector dot product calculations on the regionalized local velocity vector field data and the corresponding regional escape path reference direction vector. The velocity vector of each grid cell j in the potential bubble nucleation and aggregation area is calculated. , calculate the reference direction vector of its escape path in the area The projection on ,in is the angle between the two vectors. Considering that t_ref is a unit vector, the projection value is simplified to A positive projection value indicates that the flow direction helps the bubbles escape along the preferred path, while a negative value hinders the escape. , , forming a driving velocity component dataset for all grid cells in the region To evaluate the consistency of the local flow field, the directional consistency coefficient of the driving velocity component is calculated. , the value range is [-1, 1]. The closer the value is to 1, the more favorable the local flow field is for bubbles to escape along the preferred path. , the driving velocity component multiplied by the attenuation factor , to reflect the limiting effect of low flow rate on bubble migration. For the case where the flow direction is perpendicular to the reference direction , set the driving velocity component to zero, indicating that it neither promotes nor hinders the escape of bubbles. Sort by size and divide into positive value groups and negative value groups Calculate the average value of the positive group , the average value of the negative group ,in Considering the complexity of the flow field structure in the region, the mixing factor is introduced. , represents the proportion of favorable flow field in the region. Calculate the effective average driving velocity , taking into account the combined effects of promoting and hindering factors. In view of the sensitivity of bubble size to flow field response, a size correction factor is introduced , where d̄ is the regional average bubble size, The reference bubble size is 0.5 mm. The correction factor reflects the physical law that small bubbles are more susceptible to flow field drive. For a flow field with strong anisotropy (the standard deviation of the driving velocity component , introducing the turbulence correction factor , reflecting the interference of turbulent pulsation on bubble migration. The final regional average escape driving velocity The calculation results show that the average escape driving flow rate in the bottle shoulder transition zone (region 1) is 3.2 mm / s, the center of the bottle bottom (region 2) is 1.8 mm / s, the inner corner of the bottleneck (region 3) is 0.9 mm / s, and the side wall of the bottle (region 4) is 2.7 mm / s. When the regional average escape driving flow rate is compared with the preset critical escape rate threshold and normalized, the critical threshold is determined based on the theory of bubble kinematic mechanics. The critical escape rate threshold v_crit is calculated according to Stokes' law and the bubble floating velocity formula: , where g is the acceleration due to gravity, 9.8 m / s², is the average bubble diameter, is the gas density inside the bubble (approximately 0), The density of glass melt is 2500kg / m³, is the viscosity of the glass melt (about 100 Pa·s at the operating temperature). For an average bubble diameter of 0.5 mm, the critical escape rate threshold is calculated. The calculation formula of bubble escape efficiency index E is: ,in is the path tortuosity correction coefficient, , is the average turning angle of the streamlines in the region; is the pressure field correction coefficient, , is the average pressure in the region, is the critical pressure (95% of the saturated vapor pressure of the glass melt). Normalization ensures that the exponent E is confined to the interval [0, 1]. Values ​​closer to 1 indicate higher bubble escape efficiency in that region. Calculations show that the bubble escape efficiency index is 0.72 in the shoulder transition zone (Region 1), 0.48 in the center of the bottom (Region 2), 0.31 at the inner corner of the neck (Region 3), and 0.65 at the sidewall (Region 4).

[0031] Preferably, performing a solidification front interaction analysis on potential bubble nucleation and accumulation areas comprises: According to the temperature field in the flow field distribution data, the equivalent solidification temperature of the glass melt in each potential bubble nucleation and aggregation area and the spatial position of the isotherm are identified, and then the normal movement speed is calculated to generate the local solidification front movement speed; The relative migration velocity of bubbles under the action of buoyancy and drag forces in the potential bubble nucleation and aggregation areas is estimated based on the velocity field and bubble size in the flow field distribution data, and the regional bubble equivalent migration rate data is generated; The bubble migration-solidification rate ratio is processed according to the local solidification front movement speed and the regional bubble equivalent migration rate data to generate the bubble migration-solidification rate ratio; When the bubble migration-solidification rate ratio is lower than 1.1 and the shear strain rate is less than When , it is determined that the bubbles in the corresponding potential bubble nucleation and aggregation area are captured by the slowly moving solidification front, and the capture correction factor is set to 1.2, otherwise it is 1.0; The inverse of the bubble migration-solidification rate ratio is used as the basic capture risk and multiplied by the capture correction factor to obtain a preliminary quantitative capture index; The preliminary quantitative capture index is normalized to generate solidification capture risk data.

[0032] In the embodiment of the present invention, the soda lime glass used was thermally analyzed by differential scanning calorimetry, and its glass transition temperature was determined to be 520°C and the crystallization onset temperature was 680°C. In the actual molding process, the equivalent solidification temperature of the glass is defined as the temperature at which the glass viscosity reaches 10^7 Pa·s. According to the viscosity-temperature relationship Calculated, where is the temperature (Kelvin), μ is the viscosity (Pa·s), and the calculated equivalent solidification temperature is 700°C. The temperature values ​​of all grid nodes are extracted from the flow field distribution state data to construct a continuous temperature field function Extracted from the temperature field The isothermal surface is the solidification front. The solidification front positions at 10 different times (0.1 second intervals) are extracted around each potential bubble nucleation and aggregation area. For each area, the normal direction of the front is calculated. ,in is the temperature gradient vector. Calculate the displacement distance of the solidification front along the normal direction at adjacent moments ,in is the position vector of a point on the front at time t. The bubble size distribution of each region is extracted from the regional bubble quantitative characteristic data, and the 25%, 50%, and 75% quantile values ​​of the bubble diameter are calculated and recorded as For each region, extract the local flow field average temperature , calculate the local glass melt viscosity based on temperature Pa·s. Calculate the buoyant rise velocity of a bubble in a stationary melt , where g is the acceleration due to gravity (9.8 m / s²), d is the bubble diameter (m), is the density of glass melt (2500 kg / m3), is the gas density in the bubble (1.2 kg / m3). Considering the non-Newtonian fluid characteristics of glass melt, a shear thinning correction factor is introduced ,in is the local shear strain rate (1 / s), α = 0.05, n = 0.3, and the corrected buoyancy rise velocity Extract the local flow field average velocity vector , calculate the integrated migration velocity of bubbles relative to the melt , where g / |g| is the unit direction vector of gravity. According to the moving speed of the local solidification front When processing the bubble migration-solidification rate ratio using the regional bubble equivalent migration rate v_m_eff, the spatial relationship between the two velocity vectors is first determined. For each potential bubble nucleation and aggregation area, the solidification front normal vector is extracted. and the bubble migration direction vector . Calculate the angle between two vectors When θ<90°, the bubble migration direction is basically consistent with the solidification front advancing direction, and the solidification front chases the bubble; when θ>90°, the bubble migration direction is basically opposite to the solidification front advancing direction, and the bubble escapes from the solidification front. For the case of θ<90°, calculate the effective migration speed , represents the projection velocity of the bubble in the normal direction of the solidification front. For the case of θ>90°, set , indicating that the bubbles have completely escaped the solidification front. The bubble migration-solidification rate ratio is calculated as , R>1 means that the bubble migration speed exceeds the solidification front movement speed, and the bubble has a chance to escape; R<1 means that the solidification front movement speed exceeds the bubble movement speed, and the bubble will be captured. When the bubble migration-solidification rate ratio is lower than 1.1 and the shear strain rate is less than When determining whether the bubble is captured by the solidification front, the average shear strain rate of each potential bubble nucleation and aggregation area is first extracted from the flow field distribution state data. The shear strain rate calculation formula is = , where D is the strain rate tensor, For each region, the bubble migration-solidification rate ratio R and the average shear strain rate are examined. Whether the capture conditions are met: and When the conditions are met, it is determined that the bubbles in this area will be captured by the slowly moving solidification front, and the capture correction factor C is set to 1.2; otherwise, C is set to 1.0. The physical meaning of the correction factor is that when the shear strain rate is low, the melt flows slowly, which is not enough to enhance the migration ability of bubbles, increasing the risk of bubbles being captured. It was determined that 6 of the 8 potential areas meet the capture conditions, namely area ,area ,area ,area ,area ; The area that does not meet the conditions is the area ,area ,area Calculate the base capture risk for each region ,in is the bubble migration-solidification rate ratio. The physical meaning of is the ratio of the solidification front moving speed to the bubble migration speed. The larger the value, the faster the solidification front moves relative to the bubble, and the higher the risk of the bubble being captured. When , it means that the bubble moving speed exceeds the solidification front speed and the capture risk is low; when When the solidification front speed exceeds the bubble movement speed, the capture risk is high. Multiply the basic capture risk by the capture correction factor C to obtain the preliminary quantitative capture index. For areas that meet the capture conditions, the correction factor C=1.2, which means that under low shear rate conditions, the bubble escape ability is further reduced and the capture risk is increased by 20%; for areas that do not meet the conditions, the correction factor C=1.0, which does not increase the capture risk. Statistics of all areas value, determine the maximum value (Area 1) and minimum (Region 3). Use the linear normalization formula for conversion: ,in is the normalized solidification capture risk, ranging from 0 to 1. A larger value indicates a higher risk of the bubble being captured by the solidification front. When setting =0, indicating that bubbles are almost impossible to be captured; when When setting , indicating that the bubble is almost certainly captured. Calculated for 8 potential areas The values ​​are 1.00 (region 1), 0.77 (region 2), 0.00 (region 3), 0.89 (region 4), 0.12 (region 5), 0.75 (region 6), 0.11 (region 7), and 0.82 (region 8). The resulting solidification capture risk data will be used to calculate the subsequent bubble defect propensity value.

[0033] Preferably, step S4 includes the following steps: Step S41: Perform risk attribution analysis on the risk areas of glass bottle production defects to generate risk area dominant factor data; wherein the risk area dominant factor data includes the dominant factor of poor exhaust, the dominant factor of excessive flow shear, and the dominant factor of bubble retention caused by rapid solidification; Step S42: formulating a differentiated defect risk optimization strategy based on the risk area dominant factor data; Step S43: Optimizing design parameters based on the differentiated defect risk optimization strategy to obtain defect optimized design parameters; Step S44: Adjust the three-dimensional simulation model based on the defect optimization design parameters, and then re-execute the molten glass flow-heat transfer-bubble evolution coupled simulation to achieve the lightweight glass bottle production design requirements.

[0034] In the embodiment of the present invention, a multidimensional feature vector is constructed Characterize the risk characteristics of each region, where is the exhaust coefficient, is the flow shear coefficient, is the solidification rate coefficient. Exhaust coefficient The calculation formula is: ,in is the bubble escape efficiency index, The maximum efficiency index value (0.78) in all regions; the flow shear coefficient The calculation formula is: ,in is the regional average shear strain rate (1 / s), is the optimal shear strain rate (7.5 / s); solidification rate coefficient The calculation formula is: , i.e. solidification capture risk. For each risk area, calculate the eigenvector and determine the dominant factor index ,when When >0.8, the judgment factor i is the dominant factor. After analysis, among the 8 risk areas, poor exhaust is the dominant factor in area 1, area 4, area 6 and area 8 ( ); in regions 2 and 5, excessive flow shear is the dominant factor ( ); The dominant factor in areas 3 and 7 is the rapid solidification leading to bubble retention ( ). The generated data of dominant factors of risk areas include area number, spatial location, risk level and dominant factor type. A multi-level decision matrix method is used for systematic planning. For the dominant areas of poor exhaust (areas 1, 4, 6, and 8), an exhaust channel optimization strategy is adopted: an exhaust groove is added at the position where the solidification capture risk exceeds 0.8, and the exhaust groove width is set to 0.8 mm and the depth is 1.0 mm; the width of the existing exhaust groove is increased by 20% and the depth is increased by 15%; the direction of the exhaust groove is adjusted so that the angle between it and the priority path of bubble escape does not exceed 15 degrees. For the dominant areas of excessive flow shear (areas 2 and 5), a shear control strategy is adopted: the surface roughness of the mold cavity is reduced from Ra3.2 to Ra1.6; a flow channel buffer zone with a radius of 5 mm is set at the peak of the shear strain rate; the inclination angle of the bottle shoulder transition zone is adjusted from 45 degrees to 35 degrees, and the length of the transition curve is extended by 20%. For the areas dominated by rapid solidification (areas 3 and 7), a temperature field optimization strategy is adopted: the mold preheating temperature is adjusted from 550°C to 580°C; two heating points are added in the rapid solidification area with a power density of 2.5 watts / square centimeter; the thickness of the chrome plating layer on the mold cavity surface is increased from 0.05 mm to 0.08 mm to improve the thermal reflectivity. Through the above three differentiation strategies, a targeted defect optimization plan is formed. Determine the key parameter set that needs to be optimized: the exhaust groove parameter set , including width ,depth , direction angle and quantity ;Flow control parameter set , including surface roughness , buffer radius , shoulder angle and transition length ; Temperature control parameter set , including preheat temperature , heating point power and coating thickness . Define the multi-objective optimization function: ,in is the bubble defect index, is the weight increment of the glass bottle, For the complexity of mold manufacturing, are weight coefficients, set to 0.5, 0.3, and 0.2 respectively. Set constraints: Exhaust slot width mm, depth 0.8≤d_v≤1.2 mm, direction angle Shoulder angle ; Preheat temperature . 100 generations of optimization iterations were performed using a genetic algorithm, with a population size of 50, a crossover probability of 0.85, and a mutation probability of 0.1, ultimately obtaining the defect-optimized design parameters. The digital mold geometry model was precisely modified. For the exhaust system, eight new exhaust grooves were created on the STEP format model, with a width of 0.85 mm and a depth of 1.05 mm; the width of the existing exhaust grooves was increased from 0.75 mm to 0.9 mm, and the depth was increased from 0.8 mm to 0.92 mm; the direction of the exhaust grooves was adjusted to be parallel to the preferred path for local bubble escape. For the flow control area, the shoulder transition zone curve was adjusted from an arc transition to an elliptical curve with a major-minor axis ratio of 1.4 and a reduced inclination angle from 45 degrees to 32 degrees; four buffer zones with a radius of 5.5 mm were added at locations where the shear strain rate exceeded 15 / s. For the temperature control area, the thermal conductivity coefficient in the model material properties was adjusted from 25 W / (m·K) to 22 W / (m·K), and the mold preheating temperature was set to 585°C. The adjusted model was re-meshed, with the number of mesh elements increased from 2.5 million to 2.8 million to ensure mesh quality in critical areas. A coupled simulation of molten glass flow, heat transfer, and bubble evolution was re-executed based on the modified model. The simulation time was set to 12 seconds, with a time step of 0.0008 seconds and data collected every 0.08 seconds. The optimized bubble defect propensity value was reduced by an average of 42.5%, and the weight of the glass bottles was reduced by 5.8%, achieving the design requirements for lightweight glass bottle production.

[0035] It is particularly important to formulate differentiated defect risk optimization strategies based on the dominant factor data of risk areas, specifically: When the dominant factor classification data identifies a risk area dominated by poor exhaust, the mold parting surface is optimized and the addition of tiny exhaust ducts is planned for the gas accumulation location corresponding to the glass bottle production defect risk area. The exhaust duct diameter is set to 0.1-0.3mm, forming an exhaust optimization sub-strategy for this area. When the dominant factor classification data identifies a risk area dominated by excessive flow shear, adjustments are made to the runner shrinkage rate and gob feed rate in the mold cavity neck-shoulder transition zone at the high shear location corresponding to the glass bottle production defect risk area, forming a flow optimization sub-strategy for that area. When the dominant factor classification data identifies a risk area dominated by bubble retention due to rapid solidification, the mold cooling system parameter adjustment plan is carried out to form a solidification control sub-strategy for this area; Based on the exhaust optimization sub-strategy, flow optimization sub-strategy, and solidification control sub-strategy, strategy priority sorting and collaborative integration are carried out to obtain a differentiated defect risk optimization strategy; In the embodiment of the present invention, the position of the existing mold parting surface is evaluated. The geometric center coordinates of the four poor exhaust dominant areas, areas 1, 4, 6, and 8, are extracted. , measure the distance d from each center point to the nearest parting surface. When the parting surface is partially redesigned, the parting surface is offset to the risk area so that the offset distance mm. In the specific operation, a local offset of 12.4 mm was implemented on the parting surface of area 6, and a local offset of 8.3 mm was implemented on area 8. Secondly, the streamline tracking method was used to determine the gas accumulation position, and 5-8 tiny exhaust channels were arranged in each risk area, with a channel diameter range of 0.1-0.3 mm. Among them, area 1 was arranged with 6 exhaust channels with a diameter of 0.22 mm, area 4 was arranged with 8 exhaust channels with a diameter of 0.18 mm, area 6 was arranged with 5 exhaust channels with a diameter of 0.25 mm, and area 8 was arranged with 7 exhaust channels with a diameter of 0.15 mm. The exhaust channels are arranged in a circular layout, with a channel spacing of 3.5 mm and a channel depth designed to be 15 mm to ensure connectivity to the parting surface. The exhaust optimization sub-strategy formed includes specific parting surface offset parameters and exhaust channel specification parameters. A detailed analysis of the shear strain rate distribution was performed on areas 2 and 5. The shear strain rate gradient was calculated using the high-order interpolation method of the flow field. , determine the location where the shear rate changes most dramatically. Quantitatively analyze the geometry of the transition zone between the neck and shoulder of the cavity, and measure the rate of change of the flow channel cross-sectional area. ,in is the cross-sectional area of ​​the neck (square millimeters), is the shoulder cross-sectional area (square millimeters). For area 2, the original shrinkage rate was 58%. The transition curve was adjusted from a simple arc (radius of 12 mm) to a Bezier curve (control point coordinate offsets of 3.5, 4.2, and 2.8 mm), and the shrinkage rate was reduced to 42%. For area 5, the original shrinkage rate was 51%. The transition angle was reduced from 42 degrees to 35 degrees, and the curve length was increased from 18 mm to 24 mm, and the shrinkage rate was reduced to 38%. At the same time, the gob feeding speed was pulsed and the original constant speed was adjusted. mm / s to time-varying speed , where t is time (seconds) and T is the pulse period (0.4 seconds), reducing the local shear strain rate peak. The resulting flow optimization sub-strategy defines the geometric adjustment parameters and feeding control curve. A three-dimensional heat flow analysis of regions 3 and 7 was performed using thermal field analysis software. The local cooling rate was calculated. , where ΔT is the temperature change (°C) and Δt is the time change (seconds). The original cooling rates were 28°C / second and 24°C / second, respectively. A mold cooling system parameter adjustment plan was implemented, reducing the cooling water flow rate from 15 liters / minute to 12 liters / minute; the cooling water temperature was increased from 35°C to 42°C; and the diameter of the cooling water channel near the risk area was reduced from 8 mm to 6 mm to reduce the cooling intensity. Simultaneously, two high-frequency induction heating elements were embedded in the mold structure near the risk area, with a power density of 2.2 watts / square centimeter. The heating elements measured 15 mm × 10 mm × 3 mm and were located 10 mm from the mold surface. Differential temperature control technology was implemented, setting the target cooling rate for zone 3 to 20°C / second and for zone 7 to 18°C / second. To prevent overheating from causing new problems in other areas, four temperature monitoring points were set 25 mm from the edge of the risk area. Auxiliary cooling was initiated when the temperature exceeded 600°C. The resulting solidification control sub-strategy includes cooling parameter adjustment data and local heating element configuration schemes. The weight of each strategy is determined using the hierarchical analysis method. A judgment matrix A is established, and based on the bubble defect sensitivity, implementation complexity, and cost factor scores, the characteristic vector is calculated to obtain the weight: exhaust optimization strategy weight , flow optimization strategy weight , solidification control strategy weight . Subsequently, a strategy conflict detection was performed, and three potential conflict points were identified: there was a conflict between flow optimization and solidification control requirements at the boundary between region 2 and region 3; there was spatial overlap between region 1 and region 5 in the mold structure; there was temperature field interference between the exhaust requirements of region 8 and the solidification control of region 7. For the conflict points, the constraint reconciliation method was used for collaborative optimization: at the junction of region 2 and region 3, flow optimization was performed first, and the position of the solidification control heating element was moved back 7 mm; for the spatial overlap area between region 1 and region 5, the exhaust duct design was retained, and the flow duct adjustment range was reduced by 25%; for the mutually interfering area between region 8 and region 7, an insulation ring (width of 2 mm, thickness of 1.5 mm) was added to the periphery of the exhaust duct to reduce thermal field interference. Finally, a synergistically integrated differentiated defect risk optimization strategy was formed, which included 27 specific parameter adjustments. The execution order was to first implement exhaust system improvements, then optimize the flow duct structure, and finally adjust the cooling system.

[0036] Preferably, the present invention further provides a three-dimensional simulation design system for a digital mold for producing lightweight glass bottles, which executes the three-dimensional simulation design method for a digital mold for producing lightweight glass bottles as described above. The three-dimensional simulation design system for a digital mold for producing lightweight glass bottles comprises: The digital mesh construction module is used to perform digital geometric scanning of the glass bottle mold to construct a digital mold geometric model; perform multi-scale adaptive meshing of the cavity based on the digital mold geometric model to obtain qualified cavity mesh data; and preset the bubble nucleation point position for the qualified cavity mesh data to generate bubble nucleation point simulation data; The bubble path tracking module is used to perform bubble tracking meshing on qualified cavity mesh data using bubble nucleation point simulation data to generate bubble tracking mesh data; identify the preferred bubble escape path based on the bubble tracking mesh data, and construct a digital simulation model of the glass bottle based on the digital mold geometry model; The defect risk mapping module is used to perform coupled simulations of molten glass flow, heat transfer, and bubble evolution based on a digital simulation model of glass bottles to extract potential bubble nucleation and aggregation areas. The module also calculates bubble defect propensity values ​​for these potential bubble nucleation and aggregation areas, and performs spatial mapping of glass bottle production defects to generate glass bottle production defect risk areas. The lightweight optimization module is used to perform risk attribution analysis on defect risk areas in glass bottle production and conduct differentiated defect risk optimization to achieve lightweight glass bottle production design requirements.

[0037] The present invention is therefore intended to be illustrative and non-restrictive in all respects, with the scope of the invention being defined by the appended claims rather than the foregoing description, and all changes that come within the meaning and range of equivalents of the application documents are intended to be embraced therein.

[0038] The foregoing description is intended only to provide specific embodiments of the present invention, which will enable those skilled in the art to understand and implement the present invention. Various modifications to these embodiments will be readily apparent to those skilled in the art, and the general principles defined herein may be implemented in other embodiments without departing from the spirit or scope of the present invention. Therefore, the present invention is not intended to be limited to the embodiments shown herein, but is to be construed in the widest possible manner consistent with the principles and novel features disclosed herein.

Claims

1. A three-dimensional simulation design method for digital molds for lightweight glass bottle production, characterized in that: The following steps are involved: Step S1: performing digital geometric scanning on the glass bottle mold to construct a digital mold geometric model; Perform multi-scale adaptive meshing of the cavity based on the digital mold geometry model to obtain qualified cavity mesh data; Preset the bubble nucleation point position for the qualified cavity mesh data and generate bubble nucleation point simulation data; Step S2: performing bubble tracking grid setting on the qualified cavity grid data using the bubble nucleation point simulation data to generate bubble tracking grid data; Identify the preferred bubble escape path based on bubble tracking grid data, and build a digital simulation model of the glass bottle based on the digital mold geometry model; Step S3: performing a coupled simulation of molten glass flow, heat transfer, and bubble evolution based on the digital simulation model of the glass bottle to extract potential bubble nucleation and aggregation areas; Calculate the bubble defect tendency value for potential bubble nucleation and aggregation areas, and perform glass bottle production defect space mapping to generate glass bottle production defect risk areas; Step S4: Conduct risk attribution analysis on the defect risk areas of glass bottle production and perform differentiated defect risk optimization to achieve lightweight glass bottle production design requirements.

2. The three-dimensional simulation design method for digital molds for lightweight glass bottle production according to claim 1 is characterized in that: The digital geometric scanning of the glass bottle mold in step S1 includes: A structured light scanner is used to digitally scan the glass bottle mold cavity, core rod, and neck mold to obtain point cloud data of the mold assembly surface. De-noising was performed on the point cloud data of the mold assembly surface, the filter radius was set to 0.5 mm, isolated points and flying points were deleted, and clean mold point cloud data was obtained; Based on the clean mold point cloud data, three-dimensional surface reconstruction is performed and the mold structural features are enhanced to obtain a digital mold geometric model; among them, the mold structural features include bottle mouth threads, mold parting lines, and exhaust grooves.

3. The three-dimensional simulation design method for digital molds for lightweight glass bottle production according to claim 1 is characterized in that: In step S1, multi-scale adaptive meshing of the cavity is performed based on the digital mold geometry model, including the following steps: Extracting the three-dimensional dimensional parameters of the cavity based on the digital mold geometry model; the three-dimensional dimensional parameters of the cavity include the inner diameter of the bottle mouth, the bottom contour, the curvature radius of the neck-shoulder transition zone, the three-dimensional inclination angle of the shoulder, and the wall thickness of the bottle body; The bottle volume, surface area, and wall thickness variation rate are evaluated based on the three-dimensional cavity dimensions. When the wall thickness variation rate exceeds 15%, it is marked as a thickness transition area to generate bottle geometric feature data. Perform cavity multi-scale adaptive meshing on the digital mold geometry model using the bottle body geometric feature data to generate initial cavity mesh data; The quality of the initial cavity mesh data is evaluated, and then the cavity mesh is optimized to ensure that the Jacobian value of the minimum mesh unit is not less than 0.3 and the mesh skewness is less than 0.8, so as to obtain qualified cavity mesh data.

4. The three-dimensional simulation design method for digital molds for lightweight glass bottle production according to claim 1 is characterized in that: Step S2 includes the following steps: Step S21: Obtain the type of glass; set the solubility coefficient and diffusion coefficient of the molten gas in the glass melt according to the type of glass to generate a molten gas characteristic coefficient; Step S22: performing bubble tracking grid setting on the qualified cavity grid data based on the melt gas characteristic coefficient and the bubble nucleation point simulation data to generate bubble tracking grid data; Step S23: Analyze the position, size, and direction of the exhaust groove according to the digital mold geometric model to obtain mold exhaust structure data; Step S24: analyzing the exhaust channel position of the bubble tracking grid data using the mold exhaust structure data, and calculating the shortest distance from each bubble tracking grid data to the exhaust slot to generate distance field distribution data; Step S25: identifying a preferred bubble escape path based on the distance field distribution data; using the preferred bubble escape path and the bubble tracking grid data as glass melt-gas two-phase flow simulation configuration data; Step S26: Obtaining glass bottle production process parameters; Step S27: Import the glass melt-gas two-phase flow simulation configuration data, glass bottle production process parameters, and molten gas characteristic coefficients into the multiphase flow simulation software, set boundary conditions and initial conditions, and build a digital simulation model of the glass bottle based on the digital mold geometry model.

5. The three-dimensional simulation design method for digital molds for lightweight glass bottle production according to claim 4 is characterized in that: Identifying the preferred bubble escape paths based on distance field distribution data includes: Bubble grid node analysis is performed based on bubble tracking grid data, and the gradient vector of each bubble grid node is calculated based on the distance field distribution data to generate flow gradient field data; Based on the digital mold geometry model, the streamline path from each point in the mold cavity to the exhaust groove is tracked using the flow gradient field data to obtain the complete streamline path data; wherein, the streamline integration step size is set to 0.1mm; Calculate the path length and turning angle of each streamline based on the complete streamline path data to obtain streamline characteristic analysis data; Based on the streamline feature analysis data, the top 20% streamlines with the shortest path length and turning angle less than 45 degrees are selected as the main escape channels to obtain the main channel streamline data; Based on the main channel streamline data, the flow channel cross-sectional area changes were analyzed, and the flow channel shrinkage rate in the neck-shoulder transition area was calculated. When the shrinkage rate exceeded 40%, it was marked as a flow resistance area. The bubble escape path is selected based on the main channel streamline data through the flow resistance area to generate the bubble escape priority path.

6. The three-dimensional simulation design method for digital molds for lightweight glass bottle production according to claim 1 is characterized in that: In step S3, the coupled simulation of molten glass flow, heat transfer, and bubble evolution is performed based on the digital simulation model of the glass bottle, including: Based on the digital simulation model of the glass bottle, a coupled simulation of molten glass flow, heat transfer, and bubble evolution is performed. The temperature field, velocity field, pressure field, and gas volume fraction field are collected according to the preset cycle time window to obtain mold flow field and bubble simulation test data. Perform time averaging on the mold flow field-bubble simulation test data to generate flow field distribution state data; The shear strain rate in the cavity is calculated based on the velocity field in the flow field distribution data. Then, the low-speed flow area with a shear strain rate lower than 5s⁻¹ is identified and marked as potential flow stagnation area data. According to the temperature field in the flow field distribution data, the slow cooling area with a temperature gradient less than 2 °C / mm is identified and marked as potential slow solidification area data; According to the pressure field in the flow field distribution data, the area below 10% of the saturated vapor pressure of the glass melt is identified and marked as the potential negative pressure precipitation area data; The data of potential flow stagnation zone, potential slow solidification zone and potential negative pressure precipitation zone are spatially superimposed to obtain the potential bubble nucleation and aggregation area.

7. The three-dimensional simulation design method for digital molds for lightweight glass bottle production according to claim 6 is characterized in that: In step S3, the bubble defect tendency value is calculated for the potential bubble nucleation and aggregation area, and the glass bottle production defect space mapping is performed, including: Based on the gas phase volume fraction field in the flow field distribution state data, regional bubble characteristics of potential bubble nucleation and aggregation areas are quantified to generate regional bubble quantitative characteristic data; wherein the regional bubble quantitative characteristic data includes bubble size and gas phase volume fraction; The velocity field in the flow field distribution state data is used to process the local velocity vector of the potential bubble nucleation and aggregation area, and the bubble escape efficiency index is calculated based on the bubble escape priority path. Conducting solidification front interaction analysis on potential bubble nucleation and accumulation areas to obtain solidification capture risk data; The weighted defect tendency value of the potential bubble nucleation and aggregation area is calculated based on the bubble size, gas phase volume fraction, bubble escape efficiency index and solidification capture risk data to generate the bubble defect tendency value. The bubble size weight is set to 0.3, the gas phase volume fraction weight is set to 0.2, the bubble escape efficiency index weight is set to -0.3 and the solidification capture risk data weight is set to 0.

2. The bubble risk attention level is divided according to the bubble defect propensity value, and the glass bottle production defect space is mapped according to the corresponding potential bubble nucleation and aggregation areas to generate the glass bottle production defect risk area.

8. The three-dimensional simulation design method for digital molds for lightweight glass bottle production according to claim 7 is characterized in that: Analysis of solidification front interactions in potential bubble nucleation and accumulation areas includes: According to the temperature field in the flow field distribution data, the equivalent solidification temperature of the glass melt in each potential bubble nucleation and aggregation area and the spatial position of the isotherm are identified, and then the normal movement speed is calculated to generate the local solidification front movement speed; The relative migration velocity of bubbles under the action of buoyancy and drag forces in the potential bubble nucleation and aggregation areas is estimated based on the velocity field and bubble size in the flow field distribution data, and the regional bubble equivalent migration rate data is generated; The bubble migration-solidification rate ratio is processed according to the local solidification front movement speed and the regional bubble equivalent migration rate data to generate the bubble migration-solidification rate ratio; When the bubble migration-solidification rate ratio is lower than 1.1 and the shear strain rate is less than 3 s⁻¹, the bubbles in the corresponding potential bubble nucleation and aggregation area are determined to be captured by the slowly moving solidification front, and the capture correction factor is set to 1.2; otherwise, it is set to 1.0; The inverse of the bubble migration-solidification rate ratio is used as the basic capture risk and multiplied by the capture correction factor to obtain a preliminary quantitative capture index; The preliminary quantitative capture index is normalized to generate solidification capture risk data.

9. The three-dimensional simulation design method for digital molds for lightweight glass bottle production according to claim 1 is characterized in that: Step S4 includes the following steps: Step S41: Perform risk attribution analysis on the risk areas of glass bottle production defects to generate risk area dominant factor data; wherein the risk area dominant factor data includes the dominant factor of poor exhaust, the dominant factor of excessive flow shear, and the dominant factor of bubble retention caused by rapid solidification; Step S42: formulating a differentiated defect risk optimization strategy based on the risk area dominant factor data; Step S43: Optimizing design parameters based on the differentiated defect risk optimization strategy to obtain defect optimized design parameters; Step S44: Adjust the three-dimensional simulation model based on the defect optimization design parameters, and then re-execute the molten glass flow-heat transfer-bubble evolution coupled simulation to achieve the lightweight glass bottle production design requirements.

10. A 3D simulation design system for digital molds for lightweight glass bottle production, characterized in that: For executing the three-dimensional simulation design method of the digital mold for lightweight glass bottle production according to claim 1, the three-dimensional simulation design system for the digital mold for lightweight glass bottle production comprises: The digital mesh construction module is used to perform digital geometric scanning of the glass bottle mold to construct a digital mold geometric model; perform multi-scale adaptive meshing of the cavity based on the digital mold geometric model to obtain qualified cavity mesh data; and preset the bubble nucleation point position for the qualified cavity mesh data to generate bubble nucleation point simulation data; The bubble path tracking module is used to perform bubble tracking meshing on qualified cavity mesh data using bubble nucleation point simulation data to generate bubble tracking mesh data; identify the preferred bubble escape path based on the bubble tracking mesh data, and construct a digital simulation model of the glass bottle based on the digital mold geometry model; The defect risk mapping module is used to perform coupled simulations of molten glass flow, heat transfer, and bubble evolution based on a digital simulation model of glass bottles to extract potential bubble nucleation and aggregation areas. The module also calculates bubble defect propensity values ​​for these potential bubble nucleation and aggregation areas, and performs spatial mapping of glass bottle production defects to generate glass bottle production defect risk areas. The lightweight optimization module is used to perform risk attribution analysis on defect risk areas in glass bottle production and conduct differentiated defect risk optimization to achieve lightweight glass bottle production design requirements.

Images