A 3D simulation design method and system for digital molds used in lightweight glass bottle production

Abstract

This invention relates to the field of glass bottle manufacturing technology, and more particularly to a three-dimensional simulation design method and system for a digital mold used in lightweight glass bottle production. The method includes the following steps: constructing a digital mold geometric model for the glass bottle mold; identifying preferred bubble escape paths based on the digital mold geometric model and constructing a digital simulation model of the glass bottle; performing a coupled simulation of molten glass flow-heat transfer-bubble evolution based on the digital simulation model to extract potential bubble nucleation and aggregation regions; calculating bubble defect tendency values ​​for the potential bubble nucleation and aggregation regions and performing spatial mapping of glass bottle production defects to generate glass bottle production defect risk areas; and performing differentiated defect risk optimization on the glass bottle production defect risk areas to achieve the design requirements for lightweight glass bottle production. This invention achieves precise positioning and differentiated optimization of glass bottle bubble defect risks through bubble escape path identification and three-dimensional coupled simulation.

Description

Technical Field

[0001] This invention relates to the field of glass bottle manufacturing technology, and in particular to a three-dimensional simulation design method and system for a digital mold for lightweight glass bottle production. Background Technology

[0002] One of the most serious and common quality problems in glass bottle production is bubble defects. The main reason for these defects is the failure to effectively control the gas dissolution and escape mechanisms during the high-temperature molten glass forming process. In traditional glass bottle forming processes, when molten glass flows at high temperatures, the dissolved gases (mainly oxygen, carbon dioxide, and water vapor) need to escape through specific paths. However, due to unreasonable mold design or inaccurate temperature control, these escape paths are often blocked, causing the gas to be trapped inside the glass, forming bubbles that are difficult to eliminate. This is especially true in the production of thin-walled, lightweight glass bottles, where the reduced material thickness further increases the difficulty of gas escape, making the bubble problem even more prominent. Traditional mold design methods are insufficient to meet the production requirements of high precision and low defect rates. Currently, glass bottle mold design mainly relies on accumulated experience and simplified two-dimensional flow analysis, lacking a precise description of the complex three-dimensional flow state of molten glass within the mold cavity. Existing technologies typically use single-phase flow simulation methods, treating the molten glass as a continuous medium, ignoring the dissolution, precipitation, and migration behavior of gases in the glass melt, resulting in mold designs that cannot provide reasonable gas escape channels. 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, in order to solve at least one of the above-mentioned technical problems.

[0004] To achieve the above objectives, a three-dimensional simulation design method for a digital mold for lightweight glass bottle production includes the following steps:

[0005] Step S1: Perform digital geometric scanning on the glass bottle mold to construct a digital mold geometric model; perform multi-scale adaptive mesh generation of the cavity based on the digital mold geometric model to obtain qualified cavity mesh data; preset the bubble nucleation point positions on the qualified cavity mesh data to generate bubble nucleation point simulation data;

[0006] Step S2: Set the bubble tracking mesh for the qualified cavity mesh data using the bubble nucleation point simulation data to generate bubble tracking mesh data; identify the preferred escape path of bubbles based on the bubble tracking mesh data, and construct a digital simulation model of the glass bottle based on the digital mold geometry model;

[0007] Step S3: Perform coupled simulation of molten glass flow-heat transfer-bubble evolution based on the digital simulation model of the glass bottle to extract potential bubble nucleation and aggregation regions; calculate the bubble defect tendency value for the potential bubble nucleation and aggregation regions, and perform spatial mapping of glass bottle production defects to generate glass bottle production defect risk areas;

[0008] Step S4: Conduct risk attribution analysis on the defect risk areas in glass bottle production and optimize differentiated defect risks to achieve the design requirements for lightweight glass bottle production.

[0009] Preferably, the present invention also provides a three-dimensional simulation design system for a digital mold for lightweight glass bottle production, which executes the three-dimensional simulation design method for a digital mold for lightweight glass bottle production as described above. The three-dimensional simulation design system for a digital mold for lightweight glass bottle production includes:

[0010] The digital mesh construction module is used to perform digital geometric scanning on glass bottle molds to construct a digital mold geometric model; based on the digital mold geometric model, multi-scale adaptive meshing of the cavity is performed to obtain qualified cavity mesh data; the qualified cavity mesh data is used to preset the positions of bubble nucleation points to generate bubble nucleation point simulation data.

[0011] The bubble path tracing module is used to set the bubble tracing mesh on qualified cavity mesh data using bubble nucleation point simulation data, and generate bubble tracing mesh data; it identifies the preferred bubble escape path based on the bubble tracing mesh data, and constructs a digital simulation model of the glass bottle based on the digital mold geometry model;

[0012] The defect risk mapping module is used to perform coupled simulation of molten glass flow-heat transfer-bubble evolution based on the digital simulation model of glass bottles, so as to extract potential bubble nucleation and aggregation regions; calculate the bubble defect tendency value for potential bubble nucleation and aggregation regions, and perform spatial mapping of glass bottle production defects to generate glass bottle production defect risk regions.

[0013] The lightweight optimization module is used to perform risk attribution analysis on defect risk areas in glass bottle production and to optimize differentiated defect risks in order to meet the design requirements for lightweight glass bottle production.

[0014] This invention, by introducing a digital mold geometry model and multi-scale adaptive mesh technology, achieves precise capture of microscopic bubble behavior during glass bottle molding, fundamentally solving the problem that traditional mold design methods cannot effectively control bubble defects. By combining the preset location of bubble nucleation points with bubble tracking mesh technology, a complete prediction mechanism for bubble formation and migration is established, enabling designers to intuitively identify the preferred escape paths of bubbles and effectively avoid the problem of gas being trapped within the glass. Through coupled simulation of molten glass flow-heat transfer-bubble evolution, the limitations of traditional single-phase flow simulation are overcome, achieving full-process simulation of gas dissolution, precipitation, and migration behavior in high-temperature molten glass, significantly improving the accuracy and scientific rigor of mold design. In particular, the calculation and spatial mapping analysis of defect tendency values ​​in potential bubble nucleation and aggregation regions 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 for the bubble risk characteristics of different regions, effectively solving the problem of increased gas escape difficulty under thin-walled 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 regulated, shortening new product development cycles, and improving 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 technological path for the refined and intelligent production of the glass products industry. Therefore, the three-dimensional simulation design method for a lightweight glass bottle production digital mold of the present invention achieves precise digitization of the mold through structured light scanning and feature enhancement technology, improves the calculation accuracy by adopting a multi-scale adaptive mesh generation strategy, and innovatively introduces bubble nucleation point tracking and escape path identification to predict the bubble movement trajectory. It constructs a coupled simulation model of molten glass flow-heat transfer-bubble evolution to analyze the glass bottle forming process, accurately identifies potential bubble nucleation and aggregation regions through spatial superposition method, and quantifies the bubble defect tendency value by combining flow velocity vector processing and solidification front interaction analysis. Based on risk attribution analysis, it implements differentiated defect risk optimization, adds micro venting channels for poor venting areas, optimizes the transition zone geometry curve for areas with excessive flow shear, and adjusts the cooling system parameters for areas with excessively fast solidification, thereby achieving accurate prediction and systematic control of bubble defects in lightweight glass bottles. Attached Figure Description

[0015] Figure 1 This is a flowchart illustrating the steps of the three-dimensional simulation design method for the production of lightweight glass bottles according to the present invention.

[0016] The realization of the objective, functional features and advantages of the present invention will be further explained in conjunction with the embodiments and with reference to the accompanying drawings. Detailed Implementation

[0017] The technical method of the present invention will now be clearly and completely described with reference to the accompanying drawings. Obviously, the described embodiments are only some, not all, of the embodiments of the present invention. All other embodiments obtained by those skilled in the art based on the embodiments of the present invention without inventive effort are within the scope of protection of the present invention.

[0018] Furthermore, the accompanying drawings are merely illustrative of the invention and are not necessarily drawn to scale. The same reference numerals in the drawings denote the same or similar parts, and therefore repeated descriptions of them will be omitted. Some block diagrams shown in the drawings are functional entities and do not necessarily correspond to physically or logically independent entities. These functional entities can be implemented in software, in one or more hardware modules or integrated circuits, or in different network and / or processor methods and / or microcontroller methods.

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

[0020] To achieve the above objectives, please refer to Figure 1 This invention provides a three-dimensional simulation design method for a digital mold for lightweight glass bottle production, comprising the following steps:

[0021] Step S1: Perform digital geometric scanning on the glass bottle mold to construct a digital mold geometric model; perform multi-scale adaptive mesh generation of the cavity based on the digital mold geometric model to obtain qualified cavity mesh data; preset the bubble nucleation point positions on the qualified cavity mesh data to generate bubble nucleation point simulation data;

[0022] Step S2: Set the bubble tracking mesh for the qualified cavity mesh data using the bubble nucleation point simulation data to generate bubble tracking mesh data; identify the preferred escape path of bubbles based on the bubble tracking mesh data, and construct a digital simulation model of the glass bottle based on the digital mold geometry model;

[0023] Step S3: Perform coupled simulation of molten glass flow-heat transfer-bubble evolution based on the digital simulation model of the glass bottle to extract potential bubble nucleation and aggregation regions; calculate the bubble defect tendency value for the potential bubble nucleation and aggregation regions, and perform spatial mapping of glass bottle production defects to generate glass bottle production defect risk areas;

[0024] Step S4: Conduct risk attribution analysis on the defect risk areas in glass bottle production and optimize differentiated defect risks to achieve the design requirements for lightweight glass bottle production.

[0025] In this embodiment of the invention, the three-dimensional simulation design method for the digital mold of lightweight glass bottle production includes the following steps:

[0026] Step S1: Perform digital geometric scanning on the glass bottle mold to construct a digital mold geometric model; perform multi-scale adaptive mesh generation of the cavity based on the digital mold geometric model to obtain qualified cavity mesh data; preset the bubble nucleation point positions on the qualified cavity mesh data to generate bubble nucleation point simulation data;

[0027] In this embodiment of the invention, when performing digital geometric scanning of the glass bottle mold, a blue light structured light scanner with a resolution of 0.05 mm and a wavelength of 465 nm is used to scan the mold assembly fixed on a rotating platform from multiple angles. The platform rotates 15 degrees each time, collecting a total of 24 angle data points. The three-dimensional coordinates are calculated using a phase unwrapping algorithm to generate point cloud data. Noise reduction processing is then performed, with a filter radius of 0.5 mm, deleting isolated points and flying points more than 5 mm from the nearest point group to obtain clean point cloud data. Based on this data, three-dimensional surface reconstruction is performed using a Poisson surface reconstruction algorithm with an octree depth of 10, enhancing features such as bottle neck threads, mold parting lines, and venting grooves. Mesh generation is performed according to the digital model. First, cavity size parameters are extracted, and thickness transition regions are marked based on bottle volume, surface area, and wall thickness change rate. Multi-scale mesh generation is performed using the octree decomposition principle, with a basic mesh of 2 mm, fine regions of 0.5 mm, and thickness transition regions of 0.8 mm, setting 5 layers of structured boundary layer meshes. Mesh quality assessment and optimization are performed to ensure that the Jacobian determinant value is not less than 0.3 and the mesh skewness is less than 0.8, thus obtaining qualified cavity mesh data.

[0028] Step S2: Set the bubble tracking mesh for the qualified cavity mesh data using the bubble nucleation point simulation data to generate bubble tracking mesh data; identify the preferred escape path of bubbles based on the bubble tracking mesh data, and construct a digital simulation model of the glass bottle based on the digital mold geometry model;

[0029] In this embodiment of the invention, when setting the bubble tracking mesh for qualified cavity mesh 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 gravimetric analysis. Nitrogen is The diffusion coefficient was determined by the capillary method. A Lagrange-Eulerian hybrid grid method was used, with 5000 bubble nuclei randomly distributed within the cavity. The particle diameter ranged from 0.05 mm to 0.5 mm, and the flow velocity gradient was greater than [missing value]. The regional mesh was refined to 0.3 mm. Subsequently, the characteristics of the venting channels were analyzed. A feature recognition algorithm was used to extract the position, size, and orientation of 12 venting channels, each with a width of 0.3 mm to 0.5 mm and a depth of 0.2 mm to 0.3 mm. A fast travel method was used to calculate the distance field, establishing a three-dimensional mesh coordinate system. The surface of the venting channels was marked as the zero-distance source point. The Euclidean distance from each mesh point within the cavity to the nearest venting channel was calculated, and a flow resistance weighting factor was introduced for correction. A thermodynamic gradient descent algorithm was used to identify the preferred bubble escape path. Starting from each bubble nucleation point, the bubble trajectory was formed along the steepest direction of potential energy decrease. Three venting channels located in the bottle shoulder region handled 68% of the bubble venting.

[0030] Step S3: Perform coupled simulation of molten glass flow-heat transfer-bubble evolution based on the digital simulation model of the glass bottle to extract potential bubble nucleation and aggregation regions; calculate the bubble defect tendency value for the potential bubble nucleation and aggregation regions, and perform spatial mapping of glass bottle production defects to generate glass bottle production defect risk areas;

[0031] In this embodiment of the invention, when performing a coupled simulation of molten glass flow-heat transfer-bubble evolution based on a digital simulation model of a glass bottle, the multiphysics finite volume method is used for numerical solution. A set of conservation equations for mass, momentum, and energy, and a 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. Temperature, velocity, pressure, and gas phase volume fraction data are collected at six key time windows, and weighted time averaging is used to generate flow field distribution data. The shear strain rate within the cavity is calculated, and values ​​below [a certain value] are identified. Regions are marked as potential flow stagnation areas; temperature gradients less than [a certain value] are identified. Regions were marked as potential slow solidification zones; regions with vapor pressures 10% below saturated vapor pressure were marked as potential negative pressure precipitation zones. These three types of regions were superimposed using spatial Boolean operations to identify four main bubble nucleation and aggregation regions. Bubble quantification characteristics were analyzed in these regions, calculating bubble size and gas phase volume fraction; the bubble escape efficiency index was calculated based on flow field data; and the solidification trapping risk was calculated through solidification front interaction analysis. A multi-factor weighted method was used to calculate the bubble defect tendency value, with weights of: bubble size 0.3, gas phase volume fraction 0.2, bubble escape efficiency index -0.3, and solidification trapping risk 0.2. The results were divided into four risk levels and spatially mapped onto a digital model.

[0032] Step S4: Conduct risk attribution analysis on the defect risk areas in glass bottle production and optimize differentiated defect risks to achieve the design requirements for lightweight glass bottle production.

[0033] In this embodiment of the invention, the dominant factors for each risk area are determined through comparative analysis. Among the eight risk areas, areas 1, 4, 6, and 8 are dominated by poor venting; areas 2 and 5 are dominated by excessive flow shear; and areas 3 and 7 are dominated by excessively rapid solidification leading to bubble retention. For the areas dominated by poor venting, the parting surface offset analysis method is used to adjust the parting surface position, reducing the distance from the risk area to the parting surface to within 10 mm. Simultaneously, micro-venting channels with a diameter of 0.1-0.3 mm are arranged at the gas accumulation location, forming a venting optimization sub-strategy. For the areas dominated by excessive flow shear, the geometry of the cavity neck-shoulder transition area is adjusted, changing the simple circular arc to a Bezier curve to reduce the runner shrinkage rate. The material droplet feeding speed is also pulsed, forming a flow optimization sub-strategy. For the areas dominated by excessively rapid solidification, the mold cooling system parameters are adjusted, reducing the cooling water flow rate and increasing the water temperature. High-frequency induction heating elements are embedded near the risk areas, and differential temperature control technology is implemented, forming a solidification control sub-strategy. The weights of each strategy were determined using the analytic hierarchy process (AHP), and strategy conflict detection and collaborative optimization were performed, ultimately resulting in a differentiated defect risk optimization strategy that includes adjustments to 27 specific parameters. After optimization, the bubble defect tendency value was reduced by an average of 42.5%, and the weight of the glass bottle was reduced by 5.8%, achieving the design requirements for lightweight glass bottle production.

[0034] Preferably, the digital geometric scanning of the glass bottle mold in step S1 includes:

[0035] A structured light scanner was used to perform digital geometric scanning on the cavity, mandrel, and neck mold of the glass bottle mold to obtain point cloud data of the mold assembly surface.

[0036] The point cloud data of the mold assembly surface is denoised by setting the filter radius to 0.5mm and deleting isolated points and flying points to obtain clean mold point cloud data.

[0037] Based on the point cloud data of the clean mold, a three-dimensional surface is reconstructed and the mold structural features are enhanced to obtain a digital mold geometric model; among which, the mold structural features include bottle neck threads, mold parting lines, and venting grooves.

[0038] In this embodiment of the invention, when using a structured light scanner to perform digital geometric scanning on the cavity, mandrel, and neck mold of a glass bottle mold, the mold assembly is first fixed on a rotating platform to ensure mold stability during scanning. A blue light structured light scanner with a resolution of 0.05 mm is used, with a scanning light source wavelength of 465 nm and a projected fringe density of 100 lines / mm. During the scanning process, the rotating platform rotates 15 degrees each time, acquiring mold surface data from 24 angles 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, with the scanning distance maintained within the range of 300 mm to 500 mm. Five different phase fringe images are acquired at each angle, and the three-dimensional coordinates are calculated using a phase unwrapping algorithm to finally generate point cloud data of the mold assembly surface, with a point cloud density reaching [value missing]. When denoising the point cloud data of the mold assembly surface, the average distance from each point to its nearest neighbor is first calculated based on statistical principles. A filtering radius of 0.5 mm is set, and within this radius, the variance of the distance between the target point and its surrounding points is calculated. When the variance exceeds twice the global average variance, the point is marked as a noise point. Isolated point groups that significantly deviate from the main structure are defined as point sets with a distance greater than 5 mm from the nearest point group and are deleted. For flying points, density analysis is used; points with fewer than 5 points within a sphere with a radius of 0.5 mm are identified as flying points and deleted. Through the above processing, clean mold point cloud data is obtained with a point cloud accuracy better than 0.02 mm. When reconstructing a 3D surface based on the clean mold point cloud data, the Poisson surface reconstruction algorithm is used, with an octree depth of 10, a weight factor of 4, and a solution accuracy of 0.01 mm. For geometrically complex regions, the local point cloud density is increased, and the octree depth is increased to 12. To enhance the mold's structural features, a feature recognition algorithm was applied to the bottle neck thread area to extract the thread contour and accurately reconstruct a standard bottle neck thread with a pitch of 2.5mm and a thread height of 1.2mm. When processing the mold parting line, the edge points on both sides were extracted and fitted with a precise plane to ensure the parting line flatness error was controlled within 0.01mm. For venting groove feature enhancement, tiny grooves with a width of 0.5mm and a depth of 0.3mm were identified, and a morphological processing algorithm was used to accurately reconstruct the venting groove geometry, ensuring consistency between the digital model and the actual mold geometry.

[0039] Preferably, step S1, which involves multi-scale adaptive mesh generation of the cavity based on the digital mold geometry model, includes the following steps:

[0040] The three-dimensional dimensional parameters of the cavity are extracted from the geometric model of the digital mold; among them, the three-dimensional dimensional parameters of the cavity include the inner diameter of the bottle mouth, the bottom contour, the radius of curvature of the neck-shoulder transition area, the three-dimensional tilt angle of the shoulder, and the thickness of the bottle wall.

[0041] The bottle volume, bottle surface area, and bottle wall thickness variation rate are evaluated based on the three-dimensional dimensional parameters of the cavity. When the bottle wall thickness variation rate exceeds 15%, it is marked as a thickness transition area to generate bottle geometric feature data.

[0042] The initial cavity mesh data is generated by performing multi-scale adaptive mesh generation on the digital mold geometric model using bottle geometric feature data.

[0043] The initial cavity mesh data is quality assessed, and then the cavity mesh is optimized to ensure that the Jacobian determinant of the smallest mesh unit is not less than 0.3 and the mesh skewness is less than 0.8, thus obtaining qualified cavity mesh data.

[0044] In this embodiment of the 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 a circular outline is extracted and its diameter is calculated. The bottom outline is obtained by creating a horizontal cross-section at the bottom of the model and using an edge extraction algorithm to obtain the closed curve shape. For the neck-shoulder transition area, 10 equally spaced longitudinal cross-sections are created along the longitudinal central axis using a curvature analysis tool. The curvature change is measured on each cross-section to determine the location of the minimum curvature radius point. The three-dimensional tilt angle of the shoulder is determined by creating triangular mesh patches in the shoulder region, calculating the angle between the normal vector of each patch and the vertical direction, and taking the average value. The bottle wall thickness is determined by creating 100 uniformly distributed measurement points on the model surface, measuring the distance between the inner and outer surfaces from each point along the normal vector direction to form a wall thickness distribution map. The three-dimensional closed surface is transformed into a volume integral using Gauss's theorem, dividing the internal space of the bottle into 10,000 tiny tetrahedral units, and the total volume is obtained by summing the volumes of all units. The bottle 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 rate of change of bottle wall thickness is determined by calculating the ratio of the difference in wall thickness between adjacent measurement points to the distance between the two points. The specific calculation formula is as follows: ,in The first point represents the wall thickness (in millimeters). The second point represents the wall thickness value (in millimeters). The distance between two points is expressed in millimeters. When the wall thickness change rate between any two adjacent points exceeds 15%, the area where the line connecting these two points lies is marked as a thickness transition region and highlighted in red on the 3D model. A hierarchical octree partitioning strategy is adopted to divide the entire model space into initial mesh cells, with the side length set to 1 / 10 of the maximum size of the bottle. Subsequently, mesh refinement is performed based on the bottle's geometric feature data: in the threaded area of ​​the bottle neck, the mesh cells are subdivided to 0.5 mm; in the thickness transition region (areas with a wall thickness change rate exceeding 15%), the mesh cells are subdivided to 0.8 mm; in the neck-shoulder transition region (areas with a curvature radius less than 5 mm), the mesh cells are subdivided to 0.6 mm; in areas where the shoulder's 3D tilt angle changes by more than 30 degrees, the mesh cells are subdivided to 1.0 mm; in thin-walled areas with a wall thickness of less than 2 mm, the mesh cells are subdivided to 1.2 mm; the remaining areas maintain the original mesh size. Smooth transition zones are set between mesh cells, with the size ratio of adjacent cells not exceeding 1.5 to ensure a smooth mesh transition. After mesh generation, an initial cavity mesh containing approximately 1.5 million elements is generated. The following criteria are used: the Jacobian determinant value of each mesh element is calculated; the mesh element skewness is calculated using the formula... ,in The maximum angle (in degrees) within the unit. The minimum angle (in degrees) within the element is defined. Through global scanning, elements with a Jacobian determinant value below 0.3 (approximately 8% of the total) and a skewness greater than 0.8 (approximately 5% of the total) were identified. For these substandard elements, local mesh optimization was performed: for elements with low Jacobian values, topology optimization was used to adjust node connections; for elements with high skewness, a spring smoothing method was used to redistribute node positions, with the spring stiffness inversely proportional to the element edge length; for elements exhibiting both problems, node redistribution was performed first, followed by topology adjustment. After optimization, the mesh quality was re-evaluated to ensure that all elements had a Jacobian value no lower than 0.3 and a skewness less than 0.8, ultimately yielding qualified cavity mesh data that met simulation accuracy requirements.

[0045] Preferably, step S2 includes the following steps:

[0046] Step S21: Obtain the glass type; set the solubility coefficient and diffusion coefficient of the molten gas in the glass melt according to the glass type, and generate the characteristic coefficient of the molten gas;

[0047] Step S22: Based on the characteristic coefficients of molten gas and the simulation data of bubble nucleation points, set the bubble tracking mesh for the qualified cavity mesh data to generate bubble tracking mesh data;

[0048] Step S23: Analyze the position, size and orientation of the venting grooves based on the digital mold geometric model to obtain the mold venting structure data;

[0049] Step S24: Analyze the venting channel position of the bubble tracking grid data using the mold venting structure data, calculate the shortest distance from each bubble tracking grid data to the venting groove, and generate distance field distribution data;

[0050] Step S25: Identify the preferred bubble escape path based on the distance field distribution data; use the preferred bubble escape path and bubble tracking grid data as the simulation configuration data for glass melt-gas two-phase flow;

[0051] Step S26: Obtain glass bottle production process parameters;

[0052] 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 the boundary conditions and initial conditions, and construct the digital simulation model of the glass bottle based on the digital mold geometric model.

[0053] In this embodiment of the invention, when obtaining the glass type, an X-ray fluorescence spectrometer was used to determine the composition of the glass sample, accurately detecting the content of major components such as silicon dioxide, sodium oxide, and calcium oxide. The analysis results showed that the glass used in this embodiment was soda-lime glass, with a silicon dioxide content of 72.5%, a sodium oxide content of 14.2%, a calcium oxide content of 10.1%, and the remainder being trace elements such as magnesium oxide and aluminum oxide. Based on this glass type, the solubility coefficient of oxygen in molten glass was set as follows: This coefficient was obtained through equilibrium solubility testing at 1450°C; the diffusion coefficient of oxygen in this glass melt was set. This coefficient was determined by measuring the bubble movement rate in the glass melt at 1450°C using the rotating disk method. Simultaneously, the solubility and diffusion coefficients of other molten gases such as carbon dioxide and sulfur dioxide were set. Based on the qualified cavity mesh data from the preceding steps, the temperature gradient during the glass bottle forming process was greater than... Regions with pressure gradients greater than 0.2 MPa / cm were marked as potential bubble nucleation areas. Within these regions, the mesh cells were further refined to 0.3 mm to improve the computational accuracy during the initial bubble formation stage. Subsequently, based on the gas diffusion equation and temperature field distribution, the main bubble movement paths were predicted, and mesh density gradients were set along these paths, with a mesh size of 0.3 mm near the nucleation point, gradually increasing to 0.8 mm along the bubble movement direction. Furthermore, in predicted bubble aggregation areas, such as the upper corners of the bottle and the transition area of ​​the bottle shoulder, the mesh size was uniformly set to 0.5 mm. To track the merging behavior of microbubbles, connecting mesh bands with a cell size of 0.2 mm were created in areas where the distance between adjacent bubbles was less than 1 mm. The mold parting surface contour line was extracted, which represents the boundary where the two halves of the mold meet. A detection point was set every 5 mm along the parting surface contour line, for a total of 86 detection points. At each detection point, a depth analysis was performed perpendicular to the parting surface inside the mold to detect the groove structure. Linear structures with a groove width less than 1 mm, a depth between 0.5 and 1.2 mm, and a length greater than 5 mm were identified as venting grooves. The centerlines of all venting grooves were extracted, and the three-dimensional coordinate sequence of each centerline was recorded. The width variation of each venting groove was measured, and the narrowest and widest dimensions were precisely recorded as 0.6 mm and 0.9 mm, respectively. The depth distribution of the venting grooves was measured, and the shallowest and deepest values ​​were recorded as 0.5 mm and 1.1 mm, respectively. The orientation angle of each venting groove, i.e., the angle between the venting groove centerline and the horizontal plane, was determined, ranging from 15 degrees to 75 degrees. All measurement data were integrated to form structured mold venting structure data, including the spatial distribution, geometric dimensions, and orientation information of the venting grooves. The Euclidean distance from each grid node to the nearest venting groove was calculated using the fast traversal method. The specific implementation method is as follows: First, the surface of the exhaust channel 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 nearest exhaust channel reference point to that node, and the Euclidean distance D between the two points is calculated. The calculation formula is as follows: ,in The coordinates of the grid nodes are in millimeters. The coordinates (in millimeters) of the reference point for the exhaust groove are given. For cases with multiple exhaust grooves, the minimum distance value is taken as the exhaust distance for that node. The exhaust distance data of all nodes constitute a scalar field in three-dimensional space. Trilinear interpolation is used to interpolate the distance values ​​at the nodes to the center of the grid cells, forming a continuous distance field distribution. The distance field data is normalized to generate dimensionless distance field distribution data with values ​​ranging from 0 to 1, where 0 represents being located on the surface of the exhaust groove and 1 represents the farthest position from the exhaust groove. At each location marked as a bubble nucleation point, the negative gradient vector of the distance field is calculated, expressed as: ,in Let G be the distance field function (dimensionless), and G be the vector pointing in the direction of the fastest decrease in distance (1 / mm). The bubble's path is traced along the negative gradient direction, with a step size of 0.2 mm, until the distance field value is less than 0.05 (i.e., very close to the exhaust channel) or the outer surface of the glass bottle is reached. The drag coefficient is calculated for each tracing path, considering glass viscosity, temperature distribution, and pressure gradient. The drag coefficient is expressed as: ,in The viscosity of glass (Pa·s) The path length is in millimeters. The equivalent diameter of the pipeline path (in millimeters) The pressure correction factor is 1 / Pa. The pressure gradient is represented in Pascals per millimeter (Pa). The path with the lowest drag coefficient is marked as the preferred bubble escape path. All preferred paths are integrated with the bubble tracking mesh data to generate glass melt-gas two-phase flow simulation configuration data, including geometric mesh, bubble nucleation locations, and escape paths. The glass melting temperature is recorded as 1450°C, obtained by measuring a high-temperature thermocouple installed at the furnace outlet. The glass droplet temperature is measured as 1200°C, obtained by measuring the droplet at the moment of droplet cut using an infrared thermometer. The glass droplet weight is determined to be 355 grams, determined by weighing 30 consecutive droplets using a precision balance and taking the average value. The mold preheating temperature is recorded as 550°C, obtained by measuring a thermocouple embedded in the mold wall. The compressed air pressure is measured as 0.6 MPa, obtained by measuring the pressure on the compressed air pipeline using a pressure sensor. The blowing time parameters are recorded as follows: initial blowing time 0.8 seconds, reblowing time 1.2 seconds, and interval between two blowing times 0.5 seconds, recorded by a timer in the production line control system. Mold opening and closing times were measured: closing time 0.3 seconds, opening time 0.4 seconds, using displacement sensors to measure mold action time. Cooling parameters were recorded: mold cooling water temperature 35℃, flow rate 15 liters / minute, measured using temperature sensors and flow meters. All parameters were recorded through the production line data acquisition system, forming a complete dataset of glass bottle production process parameters. Boundary conditions were set: a mass inflow boundary was set at the glass inlet, with a flow rate of 355 g / 1.5 seconds; a temperature boundary was set on the mold wall, with a temperature of 550℃ and a thermal conductivity of 25 W / (m·K); a pressure outlet boundary was set at the venting channel, with a pressure of 0.1 MPa. Initial conditions were set: initial glass temperature 1200℃, initial pressure 0.6 MPa, initial bubble volume fraction 0.001. Based on the digital mold geometry model, a flow field and temperature field computational mesh was constructed, coupling the bubble tracking mesh with the flow field mesh. Material properties are defined as follows: glass density is 2500 kg / m³, specific heat capacity is 1200 J / (kg·K), thermal conductivity is 1.5 W / (m·K), and viscosity-temperature relationship follows the VFT equation. Where μ is the glass viscosity (Pa·s), Temperature (Kelvin), The material constants are -2.5, 4500, and 250, respectively. Gas properties and interaction parameters are set, including surface tension, contact angle, and bubble coalescence and rupture criteria. Through these settings, a digital simulation model of the glass bottle is constructed.

[0054] Preferably, identifying the preferred bubble escape path based on distance field distribution data includes:

[0055] Based on the bubble tracking grid data, bubble grid node analysis is performed, and the gradient vector of each bubble grid node is calculated based on the distance field distribution data to generate flow direction gradient field data;

[0056] Based on the digital mold geometry model, the streamline path from each point in the mold cavity to the venting groove is traced using flow gradient field data to obtain complete streamline path data; the streamline integration step size is set to 0.1 mm.

[0057] The path length and turning angle of each streamline are calculated based on the complete streamline path data to obtain streamline characteristic analysis data;

[0058] Based on the streamline characteristic analysis data, the top 20% of streamlines with the shortest path length and turning angle of less than 45 degrees were selected as the main escape channels to obtain the main channel streamline data.

[0059] The flow channel cross-sectional area change is analyzed based on the main channel streamline data, and the flow channel contraction rate in the neck-shoulder transition zone is calculated. When the contraction rate exceeds 40%, it is marked as a flow resistance area.

[0060] By analyzing the flow resistance region, the main channel streamline data is used to select bubble escape paths, generating priority escape paths for bubbles.

[0061] In this embodiment of the invention, the coordinates of all nodes in the grid data, totaling approximately 350,000 nodes, are extracted. For each node, the gradient vector of the distance field distribution data is calculated using the central difference method, with the following formula: Where D is the distance field function value (dimensionless), In spatial coordinates (millimeters), in the specific implementation, six neighboring points around the node are taken, located at... The gradient components of the central node are calculated using the distance field values ​​from the positive and negative sides of the three directions, with a spacing of 0.2 mm. , and The same method is used for calculation. For boundary nodes, one-sided difference is used to calculate the gradient. The calculated gradient vector is normalized to obtain the unit direction vector, representing the local preferred direction of bubble movement. 5000 starting points are evenly distributed within the mold cavity, with a spacing of 1.5 mm between them. Starting from each starting point, the streamline integration step size is set to 0.1 mm, and the streamline is traced along the gradient direction. The integration calculation formula for each step is as follows: ,in Let k be the current position coordinate vector (in millimeters), and k be the step size (0.1 millimeters). This is the normalized gradient vector (dimensionless) at that point. The fourth-order Runge-Kutta method improves integration accuracy by calculating four intermediate points: , , , , The integration process continues until the streamline reaches the surface of the exhaust channel (distance field value less than 0.05) or reaches the preset maximum integration step count of 500 steps. The coordinates of all integration points for each streamline are recorded sequentially to form complete streamline path data, including the three-dimensional coordinates of the starting point, ending point, and all intermediate tracking points. When calculating 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 distances between adjacent points on the streamline, using the following formula: ,in Let be the coordinates (in millimeters) of the i-th point on the streamline, and n be the number of points on the streamline. The turning angle is obtained by calculating and summing the angles formed by three adjacent points on the streamline. The formula is: ,in Let i be the direction vector of the i-th line segment. , The total turning angle (in degrees) of the streamline is given. For more detailed analysis, each streamline is divided into five segments, and the local turning angle and length of each segment are calculated. The average curvature of the streamline is also calculated using the following formula: The unit is degrees per millimeter. All calculation results are summarized to form streamline characteristic analysis data, including the total length, total turning angle, segment length, segment turning angle, and average curvature value of each streamline. The 5000 streamlines are sorted in ascending order of path length. The average length (35.8 mm) and standard deviation (12.3 mm) of all streamline lengths are calculated, and abnormal streamlines exceeding the average plus twice the standard deviation (60.4 mm) are removed. For the remaining streamlines, streamlines with turning angles less than 45 degrees are further screened, resulting in 3280 qualified streamlines. These streamlines are sorted in ascending order of path length, and the top 20% (656 streamlines) are selected as candidates for main escape channels. To ensure uniform spatial distribution, the mold cavity space is divided into 8 quadrants, and the same proportion of streamlines are selected in each quadrant to avoid excessive concentration of streamlines in local areas. Spatial clustering analysis was performed on the selected streamlines using the density clustering algorithm DBSCAN, with a cluster radius of 3 mm and a minimum number of points of 5. Streamlines with similar spatial locations were grouped into the same escape channel, resulting in 22 main escape channels containing a total of 656 streamlines, forming the main channel streamline data. Twenty cross-sections were uniformly sampled on each main channel streamline, with the cross-sectional plane perpendicular to the streamline tangent. On each cross-section, the equivalent channel diameter was calculated using fluid dynamics cross-sectional analysis methods. Specifically, 16 directions were radiated outwards from the streamline point, advancing along the direction with the steepest distance gradient until a wall was encountered or the distance gradient changed significantly (rate of change exceeding 200%). The lengths of the 16 radiating directions were measured, and the average value was taken 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: The analysis focuses on the cross-sectional area changes in the neck-shoulder transition zone (20%-40% of the streamline length), and the contraction rate is calculated. When the contraction rate exceeds 40%, the region is marked as a flow resistance region and highlighted in yellow in the 3D model. Analysis shows that among the 22 main escape channels, 8 channels have significant flow resistance regions, with contraction rates of 42.3%, 45.7%, 49.2%, 51.8%, 53.5%, 56.1%, 58.4%, and 62.9%, respectively. A weighted score is assigned to each main channel streamline, calculated using the following formula: Where L is the normalized streamline length (dimensionless, 0-1 range), θ is the normalized turning angle (dimensionless, 0-1 range), and CR is the maximum contraction rate (percentage). , , The weighting coefficients are set to 20, 30, and 50, respectively. For streamlines without significant flow resistance regions (contraction rate less than 40%), The value for this item is 0. All main channel streamlines are sorted in descending order of weight score, and the top 30% (197 paths) are selected as the preferred bubble escape paths. These paths are mainly distributed in the upper region of the bottle, with an average length of 22.6 mm, an average turning angle of 28.3 degrees, and an average contraction rate of 22.7%. To improve bubble removal efficiency, these preferred paths are clustered into 5 main bubble escape channels based on spatial location, and the inlet position of each channel is marked to form a bubble escape priority area distribution map, guiding subsequent mold structure optimization. The final generated bubble escape priority path data will be used for glass melt-gas two-phase flow simulation configuration.

[0062] Preferably, step S3, which involves performing a coupled simulation of molten glass flow-heat transfer-bubble evolution based on a digital simulation model of the glass bottle, includes:

[0063] Based on the digital simulation model of the glass bottle, a coupled simulation of molten glass flow-heat transfer-bubble evolution was performed. Temperature field, velocity field, pressure field and gas phase volume fraction field were collected according to the preset periodic time window to obtain mold flow field-bubble simulation test data.

[0064] Time averaging is performed on the mold flow field-bubble simulation test data to generate flow field distribution state data;

[0065] The shear strain rate inside the cavity is calculated based on the velocity field in the flow field distribution data. Then, low-speed flow regions with shear strain rates below 5 s⁻¹ are identified and marked as potential flow stagnation regions.

[0066] Based on the temperature field in the flow field distribution data, identify slow cooling regions with a temperature gradient of less than 2 ℃ / mm and mark them as potential slow solidification regions.

[0067] Based on the pressure field in the flow field distribution data, regions below 10% of the saturated vapor pressure of the glass melt are identified and marked as potential negative pressure precipitation zones.

[0068] By spatially overlaying data from potential flow stagnation zones, potential slow solidification zones, and potential negative pressure precipitation zones, potential bubble nucleation and aggregation regions can be obtained.

[0069] In this embodiment of the invention, the finite volume method is used to discretize and solve the governing equations. The flow equations employ the mass conservation equation and the momentum conservation equation: and ,in The density of glass (kg / m³) The velocity vector (m / s) Pressure (Pa), For viscous stress tensor (Pa), Let be the acceleration due to gravity (m / s²). The heat transfer equation uses the energy conservation equation: C p Specific heat capacity (J / (kg·K)), Let λ be the temperature (K) and λ be the thermal conductivity (W / (m·K)). This represents the viscous dissipation term (W / m³). Bubble evolution is governed by the population balance equation: ,in Bubble number density distribution function , Let the bubble radius be (m). The bubble growth rate is (m / s). For source item The time step was set to 0.001 seconds, with a total simulation time of 10 seconds. Data was collected in 0.1-second intervals. Temperature, velocity, pressure, and gas volume fraction were recorded at 25 evenly distributed monitoring points, generating mold flow field-bubble simulation test data containing 4000 sets of spatiotemporal data points. The simulation process was divided into 100 equal time intervals, each 0.1 seconds. For each grid node, the weighted average of each physical quantity at each time point was calculated, using an exponential decay function as the weighting function. ,in The current time (in seconds). To evaluate the time point (seconds), The feature time scale is 0.5 seconds. The specific calculation formula is as follows: ,in These are physical quantities (temperature, velocity, pressure, or gas volume fraction). The time average is used. For locally rapidly changing regions, an adaptive time window is used, with the window width decreasing as the physical quantity gradient increases to ensure the capture of transient phenomena. For the entire flow field, the time averages of the velocity field, temperature field, pressure field, and gas volume fraction field are calculated, where x, y, and z are spatial coordinates (millimeters). All time-averaged data are integrated to form the flow field distribution data. The velocity gradient tensor is calculated using the central difference method. The velocity gradient tensor is defined as: ,in For speed in i Directional component (m / s), for j The coordinates of the direction (in meters). 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. The formula for calculating the scalar value of shear strain rate is: Where ":" represents the tensor double dot product, , The unit is In actual calculations, for the center point of each cell in the 3D mesh, the velocity gradient in each direction is calculated using the velocity values ​​of the six adjacent points around that point, and then the shear strain rate is calculated. After calculating the shear strain rate distribution for the entire flow field, regions with shear strain rate values ​​below 5 s⁻¹ are identified. These regions exhibit slow glass melt flow and are prone to forming stagnation zones. Adjacent low shear strain rate points are connected to form regions. When the volume of a continuous region exceeds 0.5 cubic centimeters, it is marked as a potential flow stagnation zone, generating potential flow stagnation zone data, including region location, shape, and volume information. The temperature gradient vector is calculated. ,in For temperature ( ), Here are the spatial coordinates (in millimeters). The formula for calculating the temperature gradient amplitude is: The unit is / mm. A fourth-order precision finite difference scheme is used in the calculation: Where h is the grid spacing (millimeters), and A similar formula was used for calculation. The temperature gradient distribution was calculated for the entire temperature field, identifying regions with a temperature gradient amplitude less than 2℃ / mm. These regions exhibit slow cooling rates and prolonged glass solidification times. Adjacent low temperature gradient points were grouped into continuous regions; when the volume of a continuous region exceeded 0.8 cubic centimeters, it was marked as a potential slow solidification region. Feature extraction was performed on all marked regions, recording the region center location, shape parameters (length-width-height ratio), and volume, generating potential slow solidification region data containing 15 main slow-cooling areas. The saturated vapor pressure of the glass melt at the local temperature was determined. Based on experimental measurements, the saturated vapor pressure P of the soda-lime glass used in the operating temperature range (800℃-1200℃) was determined. sat The relationship between Pascal (Pa) and temperature T (Kelvin) is as follows: ,in =10.2, =8500. For each grid point in the model, the corresponding saturated vapor pressure is calculated based on the temperature at that point, and then the actual pressure is... With saturated vapor pressure The ratio is calculated as follows: Identification In regions where the vapor pressure is below 90%, or 10% below the saturated vapor pressure, gas is highly susceptible to precipitation and bubble formation. =90% threshold was used to extract isosurfaces, and the regions inside the isosurfaces were marked as potential negative pressure precipitation zones. The volume, surface area, and location of the lowest pressure point of each negative pressure precipitation zone were calculated. Regions with a volume greater than 0.3 cubic centimeters were identified as primary negative pressure precipitation zones, generating potential negative pressure precipitation zone data that includes regional geometric features and pressure distribution. The mold cavity space was discretized into cubic voxels with a side length of 0.2 mm, totaling approximately Individual voxels. For each voxel, check whether it simultaneously belongs to three potential problem regions and assign different weight coefficients: Stagnation region weight. =0.4, weight of slow solidification zone =0.3, weight of negative pressure precipitate zone =0.3. Calculate the overall risk coefficient for each voxel: ,in This is an indicator function; it takes a value of 1 when the voxel belongs to the corresponding region, and 0 otherwise. The risk coefficient R ranges from 0 to 1. When the volumetric density is ≥0.7, the voxel is marked as a potential bubble nucleation and aggregation region. A three-dimensional morphological closing operation (using a sphere with a radius of 1 mm as the structuring element) is used to connect the marked regions, eliminating isolated points and small holes, making the regions more continuous. The boundary surfaces of the continuous regions are extracted, and the volume, surface area, centroid location, and shape factor of each region are calculated. Regions with a volume greater than 0.5 cubic centimeters are identified as the main bubble nucleation and aggregation regions. A total of 8 main regions were identified, accounting for 4.2% of the total volume of the mold cavity, mainly distributed in the bottle shoulder transition area and the bottle bottom edge area.

[0070] Preferably, step S3, which involves calculating the bubble defect tendency value for the potential bubble nucleation and aggregation region and performing spatial mapping of glass bottle production defects, includes:

[0071] Based on the gas phase volume fraction field in the flow field distribution state data, the regional bubble characteristics of potential bubble nucleation and aggregation regions are quantified to generate regional bubble quantification feature data; among which, the regional bubble quantification feature data includes bubble size and gas phase volume fraction.

[0072] The velocity field in the flow field distribution data is used to perform local velocity vector processing on the potential bubble nucleation and aggregation region, and the bubble escape efficiency index is calculated based on the preferred bubble escape path.

[0073] Interactive analysis of solidification fronts was conducted on potential bubble nucleation and aggregation regions to obtain solidification capture risk data;

[0074] Based on bubble size, gas volume fraction, bubble escape efficiency index, and solidification capture risk data, a weighted defect tendency value is calculated for potential bubble nucleation and aggregation regions to generate a bubble defect tendency value. Specifically, bubble size is weighted at 0.3, gas volume fraction at 0.2, bubble escape efficiency index at -0.3, and solidification capture risk data at 0.2.

[0075] Based on the bubble defect tendency value, the bubble risk concern level is classified, and the spatial mapping of glass bottle production defects is performed according to the corresponding potential bubble nucleation and aggregation areas to generate glass bottle production defect risk areas.

[0076] In this embodiment of the invention, each potential bubble nucleation and aggregation region is independently numbered, marking a total of 8 main regions. For each region, the three-dimensional distribution data of the gas phase volume fraction α is extracted. The density clustering algorithm DBSCAN is used to analyze the data. Clustering is performed on grid points with values ​​greater than 0.001 to identify individual bubbles. The cluster radius is set to 0.5 mm, and the minimum number of points is 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 range: Statistically determine the number density n (bubbles / cubic centimeter) and average size of bubbles in each region. (mm), standard deviation of size distribution (mm), and maximum bubble size (mm). Calculate the region-average gas volume fraction: ,in For the first Gas phase volume fraction per grid cell Let this be the volume of the unit. Generate a bubble size probability density function for each region, and record the bubble size values ​​at the 10th, 50th, and 90th quantiles. Extract the velocity field data within each potential bubble nucleation and aggregation region. Calculate the average velocity vector over the region: ,in Let be the velocity vector (m / s) of the i-th grid cell. Let the volume of this unit be (cubic meters). Calculate the flow velocity direction consistency index: The value ranges from 0 to 1, with larger values ​​indicating more consistent flow directions. Subsequently, the velocity vector within each region is multiplied by the direction of the preferred bubble escape path to calculate the degree of directional consistency. ,in This is the unit direction vector of the preferred escape path of the bubble at this point. Calculate the average directional consistency of the region: Based on the magnitude of the flow velocity, consistency of its direction, and degree of consistency with the escape path, the bubble escape efficiency index is calculated: ,in The distance (in millimeters) from the area to the nearest exhaust vent. For reference distance (10 mm). The value range is 0 to A higher value indicates a stronger bubble escape capability. The bubble escape efficiency index was calculated for eight potential bubble nucleation and aggregation regions, with results of 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℃ to 680℃. Using the 700℃ isotherm as the solidification front, the advancement of the solidification front over time was tracked, and the positions of the solidification front at 100 time points were extracted, with a time interval of 0.1 seconds. The advancement velocity of the solidification front was calculated. ,in The distance (in millimeters) of the solidification front displacement at adjacent time points. The time interval is 0.1 seconds. The bubble's rising velocity is calculated using Stokes' law: Where g is the acceleration due to gravity (9.8 m / s²), Where is the bubble radius (meters) The density of liquid glass is 2500 kg / m³. The density of the gas inside the bubble is 1.2 kg / m³. The viscosity of the glass melt (Pa·s) follows the temperature-dependent property. equation: For each potential bubble nucleation and aggregation region, calculate the solidification trap risk: ,in For size The summation condition is the number of bubbles, where N is the total number of bubbles in the region, and the summation condition is that the solidification velocity is greater than the bubble rise velocity. The F value ranges from 0 to 1, and the larger the value, the higher the risk of the bubble being captured by the solidification front. All parameters are normalized to ensure their values ​​are uniformly between 0 and 1. Bubble size normalization formula: ,in The average bubble diameter (mm) for the region. The maximum bubble diameter (1.8 mm) is given across all regions. The normalized formula for gas phase volume fraction is: ,in This represents the regional average gas volume fraction. The maximum gas phase volume fraction (0.025) is given across all regions. Normalized formula for bubble escape efficiency index: Where E is the regional bubble escape efficiency index (m / s), The maximum efficiency index (0.067 m / s) is used across all regions. Since the solidification capture risk is already in the range of 0 to 1, normalization is not required; the F-value is used directly. Based on the set weighting coefficients, the bubble defect tendency value is calculated: , where 0.2 is added at the end as the basic risk value 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 defect risk regions of the glass bottle. 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 bottle shoulder transition area (accounting for 42% of the total high risk regions), the bottle bottom corner (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 structures that need to be optimized are determined, and a dataset of the defect risk regions of glass bottle production containing detailed spatial coordinates and risk levels is generated.

[0077] Especially importantly, local flow velocity vector processing is performed on the potential bubble nucleation and aggregation regions through the velocity field in the flow field distribution state data, and the bubble escape efficiency index is calculated according to the preferred bubble escape path, specifically:

[0078] Grid cell positioning is performed on the velocity field in the flow field distribution state data in each potential bubble nucleation and aggregation region, and the three-dimensional flow velocity vector components in the corresponding region are extracted to obtain regionalized local flow velocity vector field data;

[0079] Based on the preferred bubble escape path, path segment identification is performed on each potential bubble nucleation and aggregation region, and the average tangent direction of the path segment in the region is calculated to obtain the regional escape path reference direction vector;

[0080] Vector dot product operation is performed on the regionalized local flow velocity vector field data and the corresponding regional escape path reference direction vector, and the projection length in the escape path reference direction is extracted to obtain the driving flow velocity component of the bubble along the preferred path;

[0081] Statistical averaging is performed on the driving flow velocity components of the bubbles along the preferred path in each potential bubble nucleation and aggregation region to obtain the regional average escape driving flow velocity;

[0082] The average escape driving velocity in the region is compared with a preset critical escape rate threshold and normalized to generate a bubble escape efficiency index for each potential bubble nucleation and aggregation region.

[0083] In this embodiment of the invention, an octree spatial 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 regions obtained in the preceding steps, region bounding boxes are established. For each bounding box, point-polyhedral inclusion relationship detection is performed to identify the grid cells completely located within the region, generating a grid cell index list for 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 cells are extracted from the global velocity field database. Trilinear interpolation was used during extraction to ensure coordinate transformation accuracy. Considering the region boundary effect, a weighted average method was used to calculate the velocity value for mesh cells less than 0.5 mm from the boundary, with a weighting factor... ,in The shortest distance to the boundary is calculated. Statistical analysis is performed on the extracted velocity vectors to calculate the average velocity magnitude, standard deviation, maximum velocity, and minimum velocity within each region. The extracted complete three-dimensional velocity vector data is categorized and stored according to region numbering, forming regionalized local velocity vector field data. When identifying path segments for each potential bubble nucleation and aggregation region based on the bubble escape priority path, the spatial intersection relationship between the region and the bubble escape priority path is first determined. The entry points of each streamline in the bubble escape priority path through the potential bubble nucleation and aggregation region are calculated. and exit point The sets of entry and exit points are denoted as follows: and For each region, calculate the center point. ,in Let be the coordinates of all grid nodes within the region, and n be the total number of nodes. Based on the center point. right Sort the points in the region and select the m entry points closest to the center point (m is the square root of the region's volume divided by 2.5, in mm³) to form a representative set of entry points. For each entry point Find its associated exit point. This forms a local path segment traversing the region. The geometry of each path segment is reconstructed using cubic spline interpolation. Twenty equidistant points are sampled on each path segment, and the tangent vector between each adjacent pair of points is calculated. The tangent vectors of each path segment are weighted and averaged. The weighting coefficients are proportional to the flow capacity of the preferred bubble escape path, as shown in the formula. ,in This represents the flow capacity value corresponding to the streamline. The combined result is the region's escape path reference direction vector. This serves as an optimal direction indicator for bubble escape from the region. When performing a vector dot product operation between the regionalized local velocity vector field data and the corresponding regional escape path reference direction vector, a projection decomposition method is employed. This is used for the velocity vector of each grid cell j within the potential bubble nucleation and aggregation region. Calculate its escape path reference direction vector in the region. The projection on, i.e. ,in Let t_ref be the angle between the two vectors. Considering t_ref is a unit vector, the projection value simplifies to... Positive projection values ​​indicate that the flow velocity direction facilitates bubble escape along the preferred path, while negative values ​​indicate that escape is hindered. Calculate the driving component. , This forms a dataset of the driving velocity components for all grid cells within the region. To assess the consistency of the local flow field, the directional consistency coefficient of the driving velocity components is calculated. The value ranges from [-1, 1]. A value closer to 1 indicates that the local flow field is more favorable for bubbles to escape along the preferred path. For grid cells with excessively low flow velocities... The driving velocity component is multiplied by the attenuation factor. This reflects the limiting effect of low flow velocity on bubble migration. For the case where the flow velocity direction is perpendicular to the reference direction... Setting the driving velocity component to zero indicates that it neither promotes nor hinders bubble escape. The driving velocity component data for each region are then processed. Sort by size and divide into positive value groups. and negative value group Calculate the average of the positive value group. average of negative value groups ,in These represent the number of positive and negative values, respectively. Considering the complexity of the flow field structure within the region, a mixing factor is introduced. This represents the proportion of the favorable flow field in the region. The effective average driving velocity is calculated. The combined effects of promoting and hindering factors were considered. A size correction coefficient was introduced to address the sensitivity of bubble size to the flow field response. Where d̄ is the average bubble size in the region, The reference bubble size is 0.5 mm. The correction factor reflects the physical law that smaller bubbles are more easily driven by the flow field. For highly anisotropic flow fields (standard deviation of driving velocity component)... Introducing a turbulence correction factor This reflects the interference of turbulent fluctuations on bubble migration. The final region-averaged escape driving velocity... The calculation results show that the average escape driving velocity is 3.2 mm / s in the transition zone of the bottle shoulder (region 1), 1.8 mm / s in the center of the bottle bottom (region 2), 0.9 mm / s at the inner corner of the bottle neck (region 3), and 2.7 mm / s on the side wall of the bottle body (region 4). When comparing and normalizing the average escape driving velocity of each region with a preset critical escape rate threshold, the critical threshold is determined based on the theory of bubble kinematics. The critical escape rate threshold v_crit is calculated according to Stokes' law and the formula for bubble buoyancy. Where g is the acceleration due to gravity, 9.8 m / s². The average bubble diameter, The density of the gas inside the bubble is approximately 0. The density of the glass melt is 2500 kg / m³. Let be the viscosity of the glass melt (approximately 100 Pa·s at the operating temperature). The critical escape rate threshold is calculated for an average bubble diameter of 0.5 mm. The formula for calculating the bubble escape efficiency index E is as follows: ,in This is a path tortuosity correction factor. , The average turning angle of streamlines within the region; This is the pressure field correction factor. , The average pressure in the region, The critical pressure value is taken as 95% of the saturated vapor pressure of the glass melt. Normalization ensures that the index E is limited to the range of [0, 1]. The closer the value is to 1, the higher the bubble escape efficiency in that region. The calculation results show that the bubble escape efficiency index is 0.72 in the transition zone of the bottle shoulder (region 1), 0.48 in the center of the bottle bottom (region 2), 0.31 in the inner corner of the bottle neck (region 3), and 0.65 in the side wall of the bottle body (region 4).

[0084] Preferably, the solidification front interaction analysis of potential bubble nucleation and aggregation regions includes:

[0085] Based on the temperature field in the flow field distribution data, the equivalent solidification temperature of the glass melt in each potential bubble nucleation and aggregation region and the spatial location of the isotherm are identified. Then, the normal moving velocity is calculated to generate the local solidification front moving velocity.

[0086] Based on the velocity field and bubble size in the flow field distribution data, the relative migration velocity of bubbles under the action of buoyancy and drag force in the potential bubble nucleation and aggregation region is estimated, and the equivalent migration rate data of regional bubbles is generated.

[0087] The bubble migration-solidification rate ratio is generated by processing the local solidification front movement velocity and the regional bubble equivalent migration rate data.

[0088] When the bubble migration-solidification rate ratio is less than 1.1 and the shear strain rate is less than When the bubble is captured by the slowly moving solidification front in the corresponding potential bubble nucleation and aggregation region, the capture correction factor is set to 1.2; otherwise, it is 1.0.

[0089] The reciprocal of the bubble migration-solidification rate ratio is used as the basic capture risk, and multiplied by the capture correction factor to obtain the preliminary quantitative capture index;

[0090] The preliminary quantitative capture index is normalized to generate solidification capture risk data.

[0091] In this embodiment of the invention, differential scanning calorimetry was used to perform thermal analysis on the soda-lime glass used, determining its glass transition temperature to be 520℃ and its crystallization initiation temperature to be 680℃. In the actual forming process, the equivalent solidification temperature of the glass is defined as the temperature at which the glass viscosity reaches 10^7 Pa·s, determined through the viscosity-temperature relationship... The calculation shows that, among which Let be the temperature (Kelvin), and μ be the viscosity (Pa·s). The calculated equivalent solidification temperature is 700℃. Temperature values ​​of all grid nodes are extracted from the flow field distribution data to construct a continuous temperature field function. Extracting from the temperature field The isothermal surface, i.e., the solidification front, is determined. The locations of the solidification fronts at 10 different time intervals (0.1 seconds apart) are extracted around each potential bubble nucleation and aggregation region. For each region, the normal direction of the front is calculated. ,in This represents the temperature gradient vector. Calculate the displacement of the solidification front along the normal direction at adjacent time points. ,in Let be the position vector of a point on the front at time t. Extract the bubble size distribution of each region from the regional bubble quantization feature data, and calculate the 25th, 50th, and 75th percentile values ​​of the bubble diameter, denoted as [missing information]. For each region, the local flow field average temperature is extracted. Calculate the local glass melt viscosity based on temperature. Pa·s. Calculate the buoyancy-induced upward velocity of a bubble in a stationary melt. Where g is the acceleration due to gravity (9.8 m / s²), and d is the diameter of the bubble (meters). The density of the glass melt (2500 kg / m³) Let the gas density inside the bubble be 1.2 kg / m³. Considering the non-Newtonian fluid properties of the glass melt, a shear thinning correction factor is introduced. ,in The local shear strain rate (1 / s), α=0.05, n=0.3, is the corrected buoyancy rise rate. Extract the average velocity vector of the local flow field. Calculate the overall migration velocity of the bubble relative to the melt. Where g / |g| is the unit direction vector of gravity. Based on the local solidification front movement velocity... When processing the bubble migration-solidification rate ratio using the equivalent bubble migration rate v_m_eff, the spatial relationship between the two velocity vectors is first determined. For each potential bubble nucleation and aggregation region, the normal vector of the solidification front is extracted. and bubble migration direction vector Calculate the angle between two vectors. When θ < 90°, the bubble migration direction is basically the same as the solidification front's advancing direction, and the solidification front chases the bubble; when θ > 90°, the bubble migration direction is basically opposite to the solidification front's advancing direction, and the bubble escapes the solidification front. For the case where θ < 90°, the effective migration velocity is calculated. , representing the projected velocity of the bubble in the direction normal to the solidification front. For the case θ>90°, set... This indicates that the bubble has completely escaped the solidification front. The bubble migration-solidification rate ratio is calculated as follows: R>1 indicates that the bubble's migration speed exceeds the solidification front's speed, giving the bubble a chance to escape; R<1 indicates that the solidification front's speed exceeds the bubble's speed, and the bubble will be trapped. When the bubble migration-solidification rate ratio is less than 1.1 and the shear strain rate is less than... In determining whether bubbles are captured by a solidification front, the average shear strain rate of each potential bubble nucleation and aggregation region is first extracted from the flow field distribution data. The formula for calculating shear strain rate is: = Where D is the strain rate tensor. This represents the tensor double dot product. For each region, examine the bubble migration-solidification rate ratio R and the mean shear strain rate. Does the capture condition meet? and When the conditions are met, it is determined that the bubbles in that region will be captured by the slowly moving solidification front, and a capture correction factor C = 1.2 is set; otherwise, C = 1.0 is set. The physical meaning of the correction factor is that when the shear strain rate is low, the melt flow is slow, which is insufficient to enhance the migration ability of the bubbles, thus increasing the risk of bubble capture. It was determined that 6 out of the 8 potential regions meet the capture conditions, namely region... ,area ,area ,area ,area Areas that do not meet the conditions are designated as regions. ,area ,area Calculate the basic capture risk for each region. ,in This represents the ratio of bubble migration rate to solidification rate. The physical meaning of is the ratio of the solidification front's moving velocity to the bubble's migration velocity. A larger value indicates that the solidification front moves faster relative to the bubble, and the higher the risk of the bubble being trapped. When When the bubble's movement speed exceeds the solidification front speed, the risk of capture is low; when When the velocity of the solidification front exceeds the velocity of the bubble, the capture risk is relatively high. Multiplying the basic capture risk by the capture correction factor C yields the preliminary quantitative capture index. For regions meeting the capture criteria, a correction factor C = 1.2 is applied, indicating that under low shear rate conditions, the bubble escape ability is further reduced, increasing the capture risk by 20%. For regions not meeting the criteria, a correction factor C = 1.0 is applied, without additionally increasing the capture risk. Statistics are compiled for all regions. Value, determine the maximum value (Region 1) and minimum value (Region 3). The transformation is performed using a linear normalization formula: ,in This represents the normalized solidification trapping risk level, ranging from 0 to 1. A higher value indicates a higher risk of bubbles being trapped by the solidification front. In extreme cases, when... At that time, set =0, indicating that the bubble is almost impossible to capture; when At that time, set This indicates that the bubble is almost certain to be captured. Calculations are performed on 8 potential regions. 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), respectively. The final solidification capture risk data will be used for subsequent calculation of bubble defect tendency values.

[0092] Preferably, step S4 includes the following steps:

[0093] Step S41: Perform risk attribution analysis on the risk areas of glass bottle production defects and generate data on the dominant factors of the risk areas; among which, the data on the dominant factors of the risk areas include the dominant factors of poor venting, excessive flow shear, and excessively fast solidification leading to bubble retention.

[0094] Step S42: Develop differentiated defect risk optimization strategies based on the dominant factor data of the risk area;

[0095] Step S43: Optimize design parameters based on differentiated defect risk optimization strategy to obtain defect optimization design parameters;

[0096] Step S44: Adjust the three-dimensional simulation model based on the defect optimization design parameters, and then re-execute the coupled simulation of molten glass flow-heat transfer-bubble evolution to achieve the lightweight glass bottle production design requirements.

[0097] In this embodiment of the invention, a multidimensional feature vector is constructed. Characterizing the risk characteristics of each region, among which For exhaust coefficient, For the flow shear coefficient, This is the solidification rate coefficient. The exhaust coefficient is... The calculation formula is: ,in The bubble escape efficiency index. The maximum efficiency index value (0.78) is given across all regions; the flow shear coefficient is... The calculation formula is: ,in The region-average shear strain rate (1 / s) is given. The optimal shear strain rate is 7.5 / s; the solidification rate coefficient is... The calculation formula is: This refers to the risk level of solidification capture. For each risk region, an eigenvector is calculated, and the dominant factor index is determined. ,when When the value is greater than 0.8, factor i is the dominant factor. Analysis shows that in the eight risk areas, poor exhaust is the dominant factor in areas 1, 4, 6, and 8. Regions 2 and 5 are dominated by excessive flow shear. In regions 3 and 7, the primary factor was excessively rapid solidification leading to bubble retention. The generated risk area dominant factor data includes 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 venting (areas 1, 4, 6, and 8), an venting channel optimization strategy is adopted: venting channels are added at locations where the solidification capture risk exceeds 0.8, with a width of 0.8 mm and a depth of 1.0 mm; the width of existing venting channels is increased by 20%, and the depth by 15%; the orientation of the venting channels is adjusted so that the angle between them and the preferred bubble escape path 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 transition curve length is extended by 20%. For the regions dominated by rapid solidification (regions 3 and 7), a temperature field optimization strategy was adopted: the mold preheating temperature was increased from 550℃ to 580℃; two heating points were added to the rapid solidification region with a power density of 2.5 W / cm²; and the thickness of the chrome plating layer on the mold cavity surface was increased from 0.05 mm to 0.08 mm to improve heat reflectivity. Through these three differentiated strategies, a targeted defect optimization scheme was formed. The key parameter set requiring optimization was determined: the venting 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 preheating temperature Heating point power and coating thickness Define a multi-objective optimization function: ,in This is the bubble defect index. For the weight increment of the glass bottle, To increase the complexity of mold manufacturing, The weighting coefficients are set to 0.5, 0.3, and 0.2 respectively. The constraint condition is set as follows: exhaust groove width. Millimeters, depth 0.8≤d_v≤1.2 millimeters, orientation angle Shoulder angle Preheating temperature A genetic algorithm was used for 100 generations of optimization iterations, with a population size of 50, a crossover probability of 0.85, and a mutation probability of 0.1, to obtain the final defect optimization design parameters. The geometric model of the digital mold was precisely modified. For the exhaust system, eight new exhaust channels were created on the STEP format model, with a width of 0.85 mm and a depth of 1.05 mm; the width of existing exhaust channels was increased from 0.75 mm to 0.9 mm, and the depth from 0.8 mm to 0.92 mm; the exhaust channel orientation was adjusted to be parallel to the preferred path of local bubble escape. For the flow control region, the transition curve of the bottle shoulder was changed from a circular arc to an elliptical curve, with a major-to-minor axis ratio of 1.4 and an inclination angle reduced 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 region, 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℃. The adjusted model was re-meshed, with the number of mesh cells increased from 2.5 million to 2.8 million to ensure mesh quality in critical areas. Based on the modified model, a coupled simulation of molten glass flow, heat transfer, and bubble evolution was re-executed. The total simulation time was set to 12 seconds, with a time step of 0.0008 seconds, and data was collected every 0.08 seconds. The optimized model resulted in an average reduction of 42.5% in bubble defect tendency and a 5.8% reduction in glass bottle weight, achieving the design requirements for lightweight glass bottle production.

[0098] Of particular importance is the development of differentiated defect risk optimization strategies based on data on the dominant factors in risk areas.

[0099] When the dominant factor classification data determines that the risk area is dominated by poor venting, the mold parting surface is optimized and the addition of micro venting channels is planned for the gas accumulation location corresponding to the glass bottle production defect risk area; the diameter of the venting channel is set to 0.1-0.3mm, forming a sub-strategy for venting optimization in this area;

[0100] When the risk area is determined to be dominated by excessive flow shear in the dominant factor classification data, the flow optimization sub-strategy for the high shear position corresponding to the glass bottle production defect risk area is planned to adjust the flow channel shrinkage rate of the mold cavity neck-shoulder transition area and the material droplet feeding speed, forming a flow optimization sub-strategy for this area.

[0101] When the dominant factor classification data determines that the risk area is dominated by excessively rapid solidification leading to bubble retention, the mold cooling system parameter adjustment plan is carried out to form a solidification control sub-strategy for this area.

[0102] Based on the exhaust optimization sub-strategy, flow optimization sub-strategy and solidification control sub-strategy, the strategy priority is sorted and coordinated to obtain a differentiated defect risk optimization strategy.

[0103] In this embodiment of the invention, the parting line position of an existing mold is evaluated. The geometric center coordinates of regions 1, 4, 6, and 8, which are the four main areas with poor venting, are extracted. Measure the distance d from each center point to the nearest parting surface. When When the parting line is in millimeters, a local redesign of the parting line is performed, shifting the parting line to the risk area, so that the shifted distance is... In specific operations, a local offset of 12.4 mm was applied to the parting surface of region 6, and a local offset of 8.3 mm was applied to region 8. Next, the streamline tracing method was used to determine the gas accumulation location. Five to eight micro-exhaust channels with diameters ranging from 0.1 to 0.3 mm were arranged in each risk region. Specifically, six exhaust channels with a diameter of 0.22 mm were arranged in region 1, eight with a diameter of 0.18 mm in region 4, five with a diameter of 0.25 mm in region 6, and seven with a diameter of 0.15 mm in region 8. The exhaust channels were arranged in a ring layout with a spacing of 3.5 mm and a depth of 15 mm to ensure connection to the parting surface. The resulting exhaust optimization sub-strategy included specific parting surface offset parameters and exhaust channel specification parameters. A detailed analysis of the shear strain rate distribution was performed on regions 2 and 5. The shear strain rate gradient was calculated using a high-order flow field interpolation method. The location of the most drastic shear rate change was determined. A quantitative analysis of the geometry of the cavity neck-shoulder transition zone was performed by measuring the rate of change of the flow channel cross-sectional area. ,in The cross-sectional area of ​​the neck (square millimeters) The area is the shoulder cross-sectional area (square millimeters). For region 2, the original shrinkage rate was 58%. The transition curve was changed from a simple circular arc (radius 12 mm) to a Bezier curve (control point coordinate offsets of 3.5, 4.2, and 2.8 mm), reducing the shrinkage rate to 42%. For region 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, reducing the shrinkage rate to 38%. Simultaneously, the feed rate was adjusted using pulsed control, changing from a constant speed... mm / s changed to time-varying speed Where t is time (seconds) and T is the pulse period (0.4 seconds), the peak local shear strain rate is reduced. The resulting flow optimization sub-strategy clarifies the geometric adjustment parameters and feed control curve. Three-dimensional heat flow analysis of regions 3 and 7 is performed using thermal field analysis software. The local cooling rate is 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. The mold cooling system parameters were adjusted, reducing the cooling water flow rate from 15 liters / minute to 12 liters / minute; increasing the cooling water temperature from 35°C to 42°C; and reducing the diameter of the cooling water channel near the risk area from 8 mm to 6 mm to decrease cooling intensity. Simultaneously, two high-frequency induction heating elements with a power density of 2.2 W / cm² were embedded in the mold structure near the risk area. The heating element dimensions were 15 mm × 10 mm × 3 mm, and the depth was located 10 mm from the mold surface. Differential temperature control technology was implemented, setting the target cooling rate for area 3 to 20°C / second and the target cooling rate for area 7 to 18°C / second. To prevent overheating from causing new problems in other areas, four temperature monitoring points were set 25 mm away from the edge of the risk area. Auxiliary cooling was activated when the temperature exceeded 600°C. The resulting solidification control sub-strategy includes cooling parameter adjustment data and local heating element configuration schemes. The weights of each strategy are determined using the analytic hierarchy process (AHP). A judgment matrix A is established, and the weights are obtained by calculating eigenvectors based on scores for bubble defect sensitivity, implementation complexity, and cost factors: [Weight of exhaust optimization strategy]. Flow optimization strategy weights Solidification control strategy weights Subsequently, strategy conflict detection was conducted, identifying three potential conflict points: a conflict between flow optimization and solidification control requirements at the boundary between Region 2 and Region 3; spatial overlap in mold structure between Region 1 and Region 5; and temperature field interference between venting requirements in Region 8 and solidification control in Region 7. For these conflict points, a constraint reconciliation method was used for collaborative optimization: at the boundary between Region 2 and Region 3, flow optimization was prioritized, shifting the solidification control heating element position 7 mm backward; in the spatially overlapping area between Region 1 and Region 5, the venting channel design was retained, reducing the channel adjustment range by 25%; in the area where Region 8 and Region 7 interfered with each other, a heat insulation ring (2 mm wide and 1.5 mm thick) was added around the venting channel to reduce thermal field interference. Finally, a collaborative and integrated differentiated defect risk optimization strategy was formed, including 27 specific parameter adjustments, executed in the following order: first, improving the venting system; then, optimizing the channel structure; and finally, adjusting the cooling system.

[0104] Preferably, the present invention also provides a three-dimensional simulation design system for a digital mold for lightweight glass bottle production, which executes the three-dimensional simulation design method for a digital mold for lightweight glass bottle production as described above. The three-dimensional simulation design system for a digital mold for lightweight glass bottle production includes:

[0105] The digital mesh construction module is used to perform digital geometric scanning on glass bottle molds to construct a digital mold geometric model; based on the digital mold geometric model, multi-scale adaptive meshing of the cavity is performed to obtain qualified cavity mesh data; the qualified cavity mesh data is used to preset the positions of bubble nucleation points to generate bubble nucleation point simulation data.

[0106] The bubble path tracing module is used to set the bubble tracing mesh on qualified cavity mesh data using bubble nucleation point simulation data, and generate bubble tracing mesh data; it identifies the preferred bubble escape path based on the bubble tracing mesh data, and constructs a digital simulation model of the glass bottle based on the digital mold geometry model;

[0107] The defect risk mapping module is used to perform coupled simulation of molten glass flow-heat transfer-bubble evolution based on the digital simulation model of glass bottles, so as to extract potential bubble nucleation and aggregation regions; calculate the bubble defect tendency value for potential bubble nucleation and aggregation regions, and perform spatial mapping of glass bottle production defects to generate glass bottle production defect risk regions.

[0108] The lightweight optimization module is used to perform risk attribution analysis on defect risk areas in glass bottle production and to optimize differentiated defect risks in order to meet the design requirements for lightweight glass bottle production.

[0109] Therefore, the embodiments should be considered as exemplary and non-limiting in all respects, and the scope of the invention is defined by the appended claims rather than the foregoing description. Thus, all variations falling within the meaning and scope of the equivalents of the application are intended to be included within the invention.

[0110] The above description is merely a specific embodiment of the present invention, enabling those skilled in the art to understand or implement the 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 invention. Therefore, the present invention is not to be limited to the embodiments shown herein, but is to be accorded the widest scope consistent with the principles and novel features of the invention herein.

Claims

1. A three-dimensional simulation design method for a digital mold for lightweight glass bottle production, characterized in that, Includes the following steps: Step S1: Perform digital geometric scanning on the glass bottle mold to construct a digital mold geometric model; Based on the digital mold geometric model, multi-scale adaptive mesh generation of the cavity is performed to obtain qualified cavity mesh data; Pre-set the positions of bubble nucleation points in qualified cavity mesh data to generate bubble nucleation point simulation data; Step S2: Set up the bubble tracking mesh for the qualified cavity mesh data using the bubble nucleation point simulation data to generate bubble tracking mesh data; Based on the bubble tracking grid data, the preferred escape path of the bubbles is identified, and a digital simulation model of the glass bottle is constructed based on the geometric model of the digital mold. Step S3: Perform coupled simulation of molten glass flow-heat transfer-bubble evolution based on the digital simulation model of the glass bottle to extract potential bubble nucleation and aggregation regions; Calculate the bubble defect tendency value for potential bubble nucleation and aggregation regions, and perform spatial mapping of glass bottle production defects to generate glass bottle production defect risk areas; Step S4: Conduct risk attribution analysis on the defect risk areas in glass bottle production and optimize differentiated defect risks to achieve the design requirements for lightweight glass bottle production.

2. The three-dimensional simulation design method for the digital mold of lightweight glass bottle production according to claim 1, characterized in that, Step S1, which involves digitally scanning the glass bottle mold, includes: A structured light scanner was used to perform digital geometric scanning on the cavity, mandrel, and neck mold of the glass bottle mold to obtain point cloud data of the mold assembly surface. The point cloud data of the mold assembly surface was denoised by setting the filter radius to 0.5mm and deleting isolated points and flying points to obtain clean mold point cloud data. Based on the point cloud data of the clean mold, a three-dimensional surface is reconstructed and the mold structural features are enhanced to obtain a digital mold geometric model; among which, the mold structural features include bottle neck threads, mold parting lines, and venting grooves.

3. The three-dimensional simulation design method for the digital mold of lightweight glass bottle production according to claim 1, characterized in that, Step S1, which involves multi-scale adaptive mesh generation of the cavity based on the digital mold geometry model, includes the following steps: The three-dimensional dimensional parameters of the cavity are extracted from the geometric model of the digital mold; among them, the three-dimensional dimensional parameters of the cavity include the inner diameter of the bottle mouth, the bottom contour, the radius of curvature of the neck-shoulder transition area, the three-dimensional tilt angle of the shoulder, and the thickness of the bottle wall. The bottle volume, bottle surface area, and bottle wall thickness variation rate are evaluated based on the three-dimensional dimensional parameters of the cavity. When the bottle wall thickness variation rate exceeds 15%, it is marked as a thickness transition area to generate bottle geometric feature data. The initial cavity mesh data is generated by performing multi-scale adaptive mesh generation on the digital mold geometric model using bottle geometric feature data. The initial cavity mesh data is quality assessed, and then the cavity mesh is optimized to ensure that the Jacobian determinant of the smallest mesh unit is not less than 0.3 and the mesh skewness is less than 0.8, thus obtaining qualified cavity mesh data.

4. The three-dimensional simulation design method for the digital mold of lightweight glass bottle production according to claim 1, characterized in that, Step S2 includes the following steps: Step S21: Obtain the glass type; set the solubility coefficient and diffusion coefficient of the molten gas in the glass melt according to the glass type, and generate the characteristic coefficient of the molten gas; Step S22: Based on the characteristic coefficients of molten gas and the simulation data of bubble nucleation points, set the bubble tracking mesh for the qualified cavity mesh data to generate bubble tracking mesh data; Step S23: Analyze the position, size and orientation of the venting grooves based on the digital mold geometric model to obtain the mold venting structure data; Step S24: Analyze the venting channel position of the bubble tracking grid data using the mold venting structure data, calculate the shortest distance from each bubble tracking grid data to the venting groove, and generate distance field distribution data; Step S25: Identify the preferred bubble escape path based on the distance field distribution data; use the preferred bubble escape path and bubble tracking grid data as the simulation configuration data for glass melt-gas two-phase flow; Step S26: Obtain 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 the boundary conditions and initial conditions, and construct the digital simulation model of the glass bottle based on the digital mold geometric model.

5. The three-dimensional simulation design method for the digital mold of lightweight glass bottle production according to claim 4, characterized in that, Based on the distance field distribution data, the preferred escape paths for bubbles include: Based on the bubble tracking grid data, bubble grid node analysis is performed, and the gradient vector of each bubble grid node is calculated based on the distance field distribution data to generate flow direction gradient field data; Based on the digital mold geometry model, the streamline path from each point in the mold cavity to the venting groove is traced using flow gradient field data to obtain complete streamline path data; the streamline integration step size is set to 0.1 mm. The path length and turning angle of each streamline are calculated based on the complete streamline path data to obtain streamline characteristic analysis data; Based on the streamline characteristic analysis data, the top 20% of streamlines with the shortest path length and turning angle of less than 45 degrees were selected as the main escape channels to obtain the main channel streamline data. The flow channel cross-sectional area change is analyzed based on the main channel streamline data, and the flow channel contraction rate in the neck-shoulder transition zone is calculated. When the contraction rate exceeds 40%, it is marked as a flow resistance area. By analyzing the flow resistance region, the main channel streamline data is used to select bubble escape paths, generating priority escape paths for bubbles.

6. The three-dimensional simulation design method for the digital mold of lightweight glass bottle production according to claim 1, characterized in that, Step S3, which involves performing a coupled simulation of molten glass flow, heat transfer, and bubble evolution based on a 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-bubble evolution was performed. Temperature field, velocity field, pressure field and gas phase volume fraction field were collected according to the preset periodic time window to obtain mold flow field-bubble simulation test data. Time averaging is performed on the mold flow field-bubble simulation test data to generate flow field distribution state data; The shear strain rate inside the cavity is calculated based on the velocity field in the flow field distribution data. Then, low-speed flow regions with shear strain rates below 5 s⁻¹ are identified and marked as potential flow stagnation regions. Based on the temperature field in the flow field distribution data, identify slow cooling regions with a temperature gradient of less than 2 ℃ / mm and mark them as potential slow solidification regions. Based on the pressure field in the flow field distribution data, regions below 10% of the saturated vapor pressure of the glass melt are identified and marked as potential negative pressure precipitation zones. By spatially overlaying data from potential flow stagnation zones, potential slow solidification zones, and potential negative pressure precipitation zones, potential bubble nucleation and aggregation regions can be obtained.

7. The three-dimensional simulation design method for the digital mold of lightweight glass bottle production according to claim 6, characterized in that, Step S3 involves calculating the bubble defect tendency value for potential bubble nucleation and aggregation regions and performing spatial mapping of defects in glass bottle production, including: Based on the gas phase volume fraction field in the flow field distribution state data, the regional bubble characteristics of potential bubble nucleation and aggregation regions are quantified to generate regional bubble quantification feature data; among which, the regional bubble quantification feature data includes bubble size and gas phase volume fraction. The velocity field in the flow field distribution data is used to perform local velocity vector processing on the potential bubble nucleation and aggregation region, and the bubble escape efficiency index is calculated based on the preferred bubble escape path. Interactive analysis of solidification fronts was conducted on potential bubble nucleation and aggregation regions to obtain solidification capture risk data; Based on bubble size, gas volume fraction, bubble escape efficiency index, and solidification capture risk data, a weighted defect tendency value is calculated for potential bubble nucleation and aggregation regions to generate a bubble defect tendency value. Specifically, bubble size is weighted at 0.3, gas volume fraction at 0.2, bubble escape efficiency index at -0.3, and solidification capture risk data at 0.

2. Based on the bubble defect tendency value, the bubble risk concern level is classified, and the spatial mapping of glass bottle production defects is performed according to the corresponding potential bubble nucleation and aggregation areas to generate glass bottle production defect risk areas.

8. The three-dimensional simulation design method for the digital mold of lightweight glass bottle production according to claim 7, characterized in that, Interaction analysis of solidification fronts in potential bubble nucleation and aggregation regions includes: Based on the temperature field in the flow field distribution data, the equivalent solidification temperature of the glass melt in each potential bubble nucleation and aggregation region and the spatial location of the isotherm are identified. Then, the normal moving velocity is calculated to generate the local solidification front moving velocity. Based on the velocity field and bubble size in the flow field distribution data, the relative migration velocity of bubbles under the action of buoyancy and drag force in the potential bubble nucleation and aggregation region is estimated, and the equivalent migration rate data of regional bubbles is generated. The bubble migration-solidification rate ratio is generated by processing the local solidification front movement velocity and the regional bubble equivalent migration rate data. When the bubble migration-solidification rate ratio is less than 1.1 and the shear strain rate is less than 3 s⁻¹, it is determined that the bubbles in the corresponding potential bubble nucleation and aggregation region are captured by the slowly moving solidification front, and the capture correction factor is set to 1.2; otherwise, it is 1.

0. The reciprocal of the bubble migration-solidification rate ratio is used as the basic capture risk, and multiplied by the capture correction factor to obtain the 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 the digital mold of lightweight glass bottle production according to claim 1, characterized in that, Step S4 includes the following steps: Step S41: Perform risk attribution analysis on the risk areas of glass bottle production defects and generate data on the dominant factors of the risk areas; among which, the data on the dominant factors of the risk areas include the dominant factors of poor venting, excessive flow shear, and excessively fast solidification leading to bubble retention. Step S42: Develop differentiated defect risk optimization strategies based on the dominant factor data of the risk area; Step S43: Optimize design parameters based on differentiated defect risk optimization strategy to obtain defect optimization design parameters; Step S44: Adjust the three-dimensional simulation model based on the defect optimization design parameters, and then re-execute the coupled simulation of molten glass flow-heat transfer-bubble evolution to achieve the lightweight glass bottle production design requirements.

10. A three-dimensional simulation design system for a digital mold for lightweight glass bottle production, characterized in that, The three-dimensional simulation design system for executing the digital mold for lightweight glass bottle production as described in claim 1 includes: The digital mesh construction module is used to perform digital geometric scanning on glass bottle molds to construct a digital mold geometric model; based on the digital mold geometric model, multi-scale adaptive meshing of the cavity is performed to obtain qualified cavity mesh data; the qualified cavity mesh data is used to preset the positions of bubble nucleation points to generate bubble nucleation point simulation data. The bubble path tracing module is used to set the bubble tracing mesh on qualified cavity mesh data using bubble nucleation point simulation data, and generate bubble tracing mesh data; it identifies the preferred bubble escape path based on the bubble tracing mesh data, and constructs 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 simulation of molten glass flow-heat transfer-bubble evolution based on the digital simulation model of glass bottles, so as to extract potential bubble nucleation and aggregation regions; calculate the bubble defect tendency value for potential bubble nucleation and aggregation regions, and perform spatial mapping of glass bottle production defects to generate glass bottle production defect risk regions. The lightweight optimization module is used to perform risk attribution analysis on defect risk areas in glass bottle production and to optimize differentiated defect risks in order to meet the design requirements for lightweight glass bottle production.

Images