CN109190260B - Laser-arc hybrid welding three-dimensional transient numerical simulation method - Google Patents

Laser-arc hybrid welding three-dimensional transient numerical simulation method Download PDF

Info

Publication number
CN109190260B
CN109190260B CN201811044792.7A CN201811044792A CN109190260B CN 109190260 B CN109190260 B CN 109190260B CN 201811044792 A CN201811044792 A CN 201811044792A CN 109190260 B CN109190260 B CN 109190260B
Authority
CN
China
Prior art keywords
area
gas
workpiece
calculating
calculation
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
CN201811044792.7A
Other languages
Chinese (zh)
Other versions
CN109190260A (en
Inventor
庞盛永
陈鑫
母中彦
黄安国
胡仁志
罗曼乐兰
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Huazhong University of Science and Technology
Original Assignee
Huazhong University of Science and Technology
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Huazhong University of Science and Technology filed Critical Huazhong University of Science and Technology
Priority to CN201811044792.7A priority Critical patent/CN109190260B/en
Publication of CN109190260A publication Critical patent/CN109190260A/en
Application granted granted Critical
Publication of CN109190260B publication Critical patent/CN109190260B/en
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00Computer-aided design [CAD]
    • G06F30/20Design optimisation, verification or simulation
    • BPERFORMING OPERATIONS; TRANSPORTING
    • B23MACHINE TOOLS; METAL-WORKING NOT OTHERWISE PROVIDED FOR
    • B23KSOLDERING OR UNSOLDERING; WELDING; CLADDING OR PLATING BY SOLDERING OR WELDING; CUTTING BY APPLYING HEAT LOCALLY, e.g. FLAME CUTTING; WORKING BY LASER BEAM
    • B23K26/00Working by laser beam, e.g. welding, cutting or boring
    • B23K26/346Working by laser beam, e.g. welding, cutting or boring in combination with welding or cutting covered by groups B23K5/00 - B23K25/00, e.g. in combination with resistance welding
    • B23K26/348Working by laser beam, e.g. welding, cutting or boring in combination with welding or cutting covered by groups B23K5/00 - B23K25/00, e.g. in combination with resistance welding in combination with arc heating, e.g. TIG [tungsten inert gas], MIG [metal inert gas] or plasma welding
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F2111/00Details relating to CAD techniques
    • G06F2111/10Numerical modelling

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • Optics & Photonics (AREA)
  • Theoretical Computer Science (AREA)
  • Plasma & Fusion (AREA)
  • Mechanical Engineering (AREA)
  • Computer Hardware Design (AREA)
  • Evolutionary Computation (AREA)
  • Geometry (AREA)
  • General Engineering & Computer Science (AREA)
  • General Physics & Mathematics (AREA)
  • Arc Welding In General (AREA)

Abstract

The invention relates to a three-dimensional transient numerical simulation method for laser-arc hybrid welding, which solves the problems that the prior laser arc welding simulation method does not realize the overall calculation and can not truly reproduce the hybrid welding process. The invention realizes the overall accurate numerical solution of the evolution law of the pores, the molten pool, the plasma and the metal steam of the laser-TIG electric arc hybrid welding, and can be used for the theoretical research and the process optimization of the laser-TIG electric arc hybrid welding.

Description

Laser-arc hybrid welding three-dimensional transient numerical simulation method
Technical Field
The invention belongs to a welding numerical simulation method, and particularly relates to a laser-arc hybrid welding three-dimensional transient numerical simulation method.
Technical Field
The laser-electric arc hybrid welding can combine the advantages of two heat sources, namely laser and electric arc, and realize advantage complementation; compared with electric arc welding, the laser-electric arc hybrid welding has the advantages of large fusion depth, high speed, small residual stress and deformation; compared with laser welding, the laser-electric arc hybrid welding has the advantages of enhanced bridging capacity, easy adjustment of welding seam components, reduced cold brittleness phase tendency and improved pinhole stability. The laser-arc hybrid welding has unique advantages, can realize high-quality and high-efficiency welding manufacture, becomes a welding technology with wide application prospect, and has important application in the fields of ocean, nuclear power and the like.
As shown in fig. 1, the laser and the electric arc are directly used as heat sources to act on the material, so that the material is subjected to multi-state transformation such as melting, gasification and ionization, a molten pool, a small hole, photoinduced metal vapor/plasma, electric arc plasma and the like coexist, severe physical changes such as extremely high temperature, high pressure, high speed, instantaneous and the like exist in a millimeter-scale area, and the complex physical processes and characteristics make it difficult to deeply research the composite welding mechanism by adopting an experimental method. Meanwhile, the laser-arc hybrid welding combines two heat sources, the technological parameters of the laser-arc hybrid welding are further increased, and more difficulties are brought to obtaining the optimal technological parameters through an experimental method. The numerical simulation method can be used for carrying out visual quantitative analysis on all physical quantities and evolution rules thereof in the laser-arc hybrid welding process, and can reveal the evolution rules of pinholes, a molten pool, plasma and metal steam along with time in the welding process, so that the numerical simulation has important significance in the application of laser-arc hybrid welding research.
The laser-arc hybrid welding process is extremely complex, physical modeling is difficult, the density, temperature, potential, speed and the like on two sides of an interface are suddenly changed in an extremely small space scale, so that the characteristics of large gradient and high nonlinearity are caused near the interface, a unified solution dual-phase flow model is difficult to solve, and accurate solution of the interface position is difficult to obtain. In view of the above difficulties, those skilled in the art either neglect to calculate the arc and metal vapor only for the molten pool, see Z Gao, et al, analysis of the world porous ceramic stationary Laser-MIG hybrid welding, the International journal Advanced Manufacturing Technology,2009,44(9-10):870, or neglect to calculate the arc only for the molten pool, see Y.T.Cho, et al, numerical analysis of the hybrid planar Manufacturing by Nd: YAGlaser and gas porous arc & Laser Technology,2011,43(3):711 and 720, and there is no current report on the overall numerical calculation implementation of Laser-arc hybrid welding. However, the behavior of the plasma during hybrid welding directly determines the thermal effect on the weld pool, and the metal vapor generated by the local violent evaporation of the weld pool will also change the composition and properties of the plasma. Therefore, the method of only calculating the molten pool or only calculating the electric arc can not truly reproduce the composite welding process, and is not beneficial to the research of the interaction mechanism of the laser-electric arc composite welding.
Disclosure of Invention
The invention provides a laser-arc hybrid welding three-dimensional transient numerical simulation method, which solves the problems that the conventional laser arc welding simulation method does not realize integral calculation and can not truly reproduce the hybrid welding process, can realize the integrated simulation calculation of melting, solidification and evaporation behaviors, heat transfer flow of arc plasma, laser plume gas and metal molten pool and electromagnetic field evolution in the laser-TIG arc hybrid welding process, can reproduce the physical process more comprehensively and truly, is convenient to carry out deeper theoretical research, and can further promote the engineering application of laser-arc hybrid welding.
The invention provides a laser-arc hybrid welding three-dimensional transient numerical simulation method, which comprises the steps of establishing a geometric model, setting an initial value, updating a time step, updating physical parameters, solving an electromagnetic field, decomposing a space discontinuity, calculating a gas region state, calculating a workpiece region state, updating a molten pool free interface and judging, and is characterized in that:
firstly, establishing a geometric model:
establishing a geometric model by adopting three-dimensional modeling software, wherein the geometric model is a cuboid, the geometric model comprises a workpiece area, an electrode area and a gas area, and the workpiece area is the cuboid and has the same length and width as the geometric model; the electrode area is a cylinder or a cylinder with one end being a cone, when the electrode area is a cylinder, the included angle between the central axis of the electrode area and the vertical direction is 0-60 degrees, and when the electrode area is a cylinder with one end being a cone, the diameter of the bottom surface of the cone is the same as that of the cylinder; the vertical distance between the end part of the electrode area and the workpiece area is 1 mm-10 mm; the gas area is the rest area of the geometric model except the workpiece area and the electrode area;
the region of the geometric model except the electrode region is a calculation region;
carrying out grid division on the calculation region by adopting grid division software, wherein the grid is a cuboid or a cube, and the side length of the cuboid or the cube is 10-500 mu m;
secondly, setting an initial value:
given the boundary conditions of the computation region: the environment temperature is 300-1000K, the electrode temperature is 3000-3500K, the environment pressure is 0.1-10 atm, and the bottom potential of the workpiece is 0V;
setting the initial value of the physical quantity of each grid center point in the whole calculation area: potential 0V, current density 0A/m2The magnetic field intensity is 0T, the speed is 0m/s, the temperature is 300-1000K, the pressure is 0.1-10 atm, and the mass percentage content of the metal steam is 0;
setting an accumulation time variable tsumCalculating the cutoff time t as 0end=0.1~10s;
Setting the time step Δ tiIs an initial time step Δ t0,Δt0=1×10-7~1×10-6s, performing the fourth step;
step three, updating the time step:
calculating the maximum time step for gas region stability calculation
Figure BDA0001793010580000041
Then calculating the maximum time step length of stable calculation of the workpiece region
Figure BDA0001793010580000042
Take (Δ t)0,Δtg,Δtm) The minimum value of the three is used as the updated time step delta tiCarrying out the step four;
where, i is 1,2 …, which indicates the i-th step, the maximum condition number Cmax=0.5~1.0,ug_x,ug_y,ug_zGas velocity at the center point of the gas area grid
Figure BDA0001793010580000043
Component in three directions x, y, z, um_x,um_y,um_zThe velocity of the molten pool being the central point of the grid of the region of the workpiece
Figure BDA0001793010580000044
Components in the x, y, z directions; the sizes of the grids in the x direction, the y direction and the z direction are respectively delta x, delta y and delta z;
fourthly, updating physical property parameters:
setting physical property parameters in the whole calculation region:
the gas area is filled with mixed gas of protective gas and metal steam, the protective gas is one or more of argon, helium and carbon dioxide, and the metal steam is the metal steam of the workpiece material; calculating the mass percentage content and the temperature of the metal steam in the mixed gas as initial values for the first time, and acquiring the mass percentage content and the temperature of the metal steam in the mixed gas in the next time step from the substeps (7.2) and (7.3);
conductivity σ of the gas mixtureeDensity rhogViscosity [ mu ] ofgSpecific heat capacity CpgCoefficient of thermal conductivity kgAnd net gas emissivityNThe composition and temperature of the mixed gas at the center of each grid are directly determined;
the workpiece material of the workpiece area is iron-based alloy or non-ferrous metal alloy; corresponding conductivity σeDensity rhomViscosity [ mu ] ofmSpecific heat capacity CpmCoefficient of thermal conductivity kmThermal expansion coefficient β, surface tension coefficient gamma, and thermal capillary coefficient
Figure BDA0001793010580000051
Work potential
Figure BDA0001793010580000052
Emissivity of radiationrDirectly determined by the workpiece material and temperature at each grid center;
fifthly, solving the electromagnetic field:
solving for the electromagnetic field throughout the calculation area:
welding current is introduced in the welding process, an electromagnetic field is generated in the whole calculation area, welding current of 30-300A is given through the electrodes, and potential of the center points of all grids in the whole calculation area is calculated and obtained
Figure BDA0001793010580000053
Current density
Figure BDA0001793010580000054
And magnetic induction intensity
Figure BDA0001793010580000055
Distribution, the calculation equation is as follows:
Figure BDA0001793010580000056
Figure BDA0001793010580000057
Figure BDA0001793010580000058
Figure BDA0001793010580000059
wherein the vacuum dielectric constant is0Magnetic loss of potential at center point of grid
Figure BDA00017930105800000510
Gradient operator
Figure BDA00017930105800000511
Divergence operator
Figure BDA00017930105800000512
Rotation operator
Figure BDA00017930105800000513
Laplacian operator
Figure BDA00017930105800000514
Electrical conductivity sigmaeRespectively adopting the mixed gas in the gas area and the corresponding conductivity of the workpiece material in the workpiece area during calculation;
sixthly, spatial discontinuous decomposition:
performing space discontinuous decomposition by taking the upper surface of a workpiece area as an interface in the first calculation and taking an updated free interface of a molten pool as an interface in the subsequent calculation, and dividing the whole calculation area into a gas area and a workpiece area, wherein the gas area is an arc-metal vapor area and the workpiece area is a workpiece solid-liquid area; defining interface boundary points as intersection points of grid central point connecting lines at two sides of a free interface of the molten pool and the free interface of the molten pool;
seventhly, calculating the gas area state:
the method comprises the following substeps of calculating the flow field, the temperature field and the metal vapor distribution of the gas region in sequence:
(7.1) calculating the gas area flow field: calculating the gas velocity of the center points of all grids in the gas area
Figure BDA0001793010580000061
Gas pressure PgThe calculation formula is as follows:
Figure BDA0001793010580000062
Figure BDA0001793010580000063
wherein
Figure BDA0001793010580000064
Is the acceleration of gravity;
setting the gas velocity of the boundary point of the gas zone and the workpiece zone at the free interface of the molten bath
Figure BDA0001793010580000065
And gas pressure Pg_sRespectively as follows:
Figure BDA0001793010580000066
wherein the evaporation recoil pressure PrBy using pressure based on the environmentObtaining a recoil pressure model;
(7.2) calculating the gas zone temperature field: calculating the gas temperature T of the central point of all grids in the gas areagThe calculation formula is as follows:
Figure BDA0001793010580000071
wherein k isbBoltzmann constant, e is elementary charge, t is time,
Figure BDA0001793010580000075
a differential operator;
calculating the heat flow q of the boundary point of the interface into the gas region at the free interface of the molten pool of the gas region and the workpiece regiongas
Figure BDA0001793010580000072
Wherein T isgs、TmsAcquiring the temperatures of a gas area and a workpiece area which are boundary points of an interface respectively by adopting a virtual flow method; k is a radical ofeffIs the effective thermal conductivity, k, of the boundary point of the interfaceeff=0.5×(kg(Tgs)+kg(Tms)),kg(Tgs)、kg(Tms) Respectively is the temperature Tgs、TmsThe corresponding gas thermal conductivity; the effective thickness of the boundary layer of the free interface of the molten pool is 0.1 mm-0.4 mm;
(7.3) calculating the gas zone metal vapor diffusion: calculating the mass percentage content C of the metal vapor at the central points of all grids in the gas areavaporThe calculation formula is as follows:
Figure BDA0001793010580000073
wherein D is the mass diffusion coefficient of the metal vapor, and the value can be obtained by a second viscosity approximation method;
setting the metal vapour quality of boundary point of each interface at the free interface of molten pool in gas region and workpiece regionAmount percentage of Cvapor
Figure BDA0001793010580000074
Eighthly, calculating the state of the workpiece area:
sequentially calculating a flow field and a temperature field of a workpiece region, and comprising the following substeps:
(8.1) calculating the flow field of the workpiece area: calculating the molten pool speed of the central points of all grids in the workpiece area
Figure BDA0001793010580000081
Pressure P of the molten bathmThe calculation formula is as follows:
Figure BDA0001793010580000082
Figure BDA0001793010580000083
wherein the temperature T of the central point of the grid of the workpiece areamTaking the ambient temperature during the first calculation; the next time step is calculated by adopting the calculated value of the last time step (7.2);
reference temperature TrefTaking the melting point temperature of a workpiece material; K. c is the seepage coefficient, relaxation coefficient and liquidus fraction flClosely related:
Figure BDA0001793010580000084
wherein d is a constant closely related to the size of dendrite spacing, and d is 0.05mm to 0.5 mm; fraction of liquid phase flLinear with temperature:
Figure BDA0001793010580000085
wherein, Tl、TsThe liquidus temperature and solidus temperature of the workpiece material, respectively;
calculating the normal pressure P of the molten pool at the boundary point of the interface at the free interface of the molten pool in the gas area and the workpiece areafAnd the tangential force of the molten pool
Figure BDA0001793010580000086
Figure BDA0001793010580000087
Wherein
Figure BDA0001793010580000088
Normal vector and tangent vector of the free interface of the molten pool, kappa is the curvature of the free interface of the molten pool, PeffFor effective pressure of gas region to workpiece region, take Peff=max(Pr,Pg);
(8.2) calculating the temperature field of the workpiece area: calculating the temperature T of the central point of all grids in the workpiece areamThe calculation formula is as follows:
Figure BDA0001793010580000091
wherein the workpiece absorbs the energy of the laser beam qlaserAcquiring by adopting a ray tracing method; the laser is fiber laser or Nd, namely YAG laser, the power is 100W-10000W, and the spot radius is 0.1 mm-0.5 mm;
calculating the heat flow q of the boundary point of the interface into the workpiece region at the free interface of the molten pool of the gas region and the workpiece regionmatel
Figure BDA0001793010580000092
Wherein σ is Stefan-Boltzmann constant, TIs the ambient temperature, qevapHeat loss for evaporation;
ninthly, updating a free interface of the molten pool:
molten pool velocity through all grid center points of workpiece area
Figure BDA0001793010580000093
Calculating the distance phi from the central point of each grid to the interface in the area, wherein the obtained equivalent surface with phi equal to 0 is the updated free interface of the molten pool, and the calculation formula is as follows:
Figure BDA0001793010580000094
the normal vector and curvature of the free interface of the molten pool are solved as follows:
Figure BDA0001793010580000095
in the fifth step to the ninth step, all the involved calculation equations are dispersed by adopting a finite difference method or a finite volume method, and the time step delta t isiSolving is carried out internally;
tenthly, judging:
will tsum+ΔtiIs given to tsumJudging whether t is presentsum<tendIf yes, turning to the step three, updating the time step, otherwise, ending the calculation.
In the first step, the three-dimensional modeling Software is AutoCAD of American Autodesk Software company, UG of Germany Siemens PLM Software company or Pro/E of American parameter technology company; the cuboid of the geometric model is 5 mm-500 mm long, 5 mm-500 mm wide and 5 mm-200 mm high; the cuboid of the workpiece area is 1 mm-100 mm high; when the electrode area is a cylinder, the diameter of the cylinder is 1 mm-4 mm, the height is 1 mm-150 mm, when the electrode area is a cylinder with one end being a cone, the cone angle is 30-120 degrees, the diameter of the cylinder is 1 mm-4 mm, and the height is 1 mm-150 mm; the gridding software is Hypermesh of Altair corporation, Gridgen of Pointwise corporation, or Patran of MSC.
6. In the fourth step, the electric conductivity σ of the mixed gas is determined by the composition and temperature of the mixed gaseDensity rhogViscosity [ mu ] ofgSpecific heat capacity CpgCoefficient of thermal conductivity kgAnd gas purifying agentCoefficient of radiationNReference is made to the substances from Mougenot J, et al, plasma-assisted interaction in structured input-gas configuration, journal of Physics D: Applied Physics,2013,46(13):135206, or Murphy A B, et al, modeling of thermal plasma for arc welding, the roll of the shifting gas properties and of metallic force. journal of Physics D: Applied Physics,2009,42(19):194006, or Murphy AB, the effects of metallic in arc welding, journal of Physics D: Applied Physics,2010,43(43): 434001;
the corresponding electrical conductivity σ of the workpiece material being determined directly from the workpiece material and the temperatureeDensity rhomViscosity [ mu ] ofmSpecific heat capacity CpmCoefficient of thermal conductivity kmThermal expansion coefficient β, surface tension coefficient gamma, and thermal capillary coefficient
Figure BDA0001793010580000111
Work potential
Figure BDA0001793010580000112
Emissivity of radiationrObtained by calculation using the Material Performance calculation Software JMatPro, Inc. of Sente Software, UK, either cited by Wei HL, et al, original of grain orientation reduction solution of an aluminium alloy, acta materials, 2016,115:123-131, or by Wu C S.computer correlation of this-dimensional correlation in transforming MIG well locations, engineering compositions, 1992,9(5): 529-537.
In the sub-step (7.1), the recoil pressure model based on the environmental pressure is as follows: s. Pang, et. Explosition of duration reduction Laser welding indirect representations. journal of Laser Applications 2015,27(2): 022007; in the sub-step (7.2), the virtual stream method is as follows: fedkiw R P, et al, A non-oscilloscopic Eulerapaach to interfaces in multimaterial flows Journalo f computational physics,1999,152(2): 457-; in the sub-step (7.3), the second viscosity approximation method is as follows: A.B.Murphy.A composition of transactions of differentiation in thermoplasmas. journal of Physics D: Applied Physics,1996,29(7): 1922.
In the sub-step (8.2), the ray tracing method is as follows: computing the energy density of the laser deep fusion welding small holes with arbitrary shapes, the laser technology 2010,34(05) 614-618; in the step (8.2), the heat loss q is evaporatedevapThe calculation method is as follows: wu D, et al.infraduction of scatter formation in fiber laser welding of5083aluminum alloy&Mass Transfer,2017,113:730-740。
The invention can be realized by high-level programming languages such as C + + computer, and can simulate the three-dimensional dynamic behaviors of electric arc, molten pool, small hole and metal vapor in the laser-electric arc welding process by the computer.
The invention solves the whole calculation area uniformly when solving the electromagnetic field, divides the integral calculation area with multiple phases (liquid phase, solid phase, electric arc and metal vapor) into a gas area (electric arc-metal vapor area) and a workpiece area (workpiece solid-liquid area) through the instantaneous free interface of the workpiece molten pool when solving the concentration of the flow field, the temperature field and the metal vapor, respectively calculates two areas with physical and mathematical discontinuity in a minimum time step, and carries out bidirectional sequential coupling at the interface position in a boundary loading mode through respective calculation data. The method provided by the invention can well solve the problem of nonlinear solution near the interface, realize the heat transfer flow of a molten pool and arc plasma, the diffusion behavior of metal vapor, free interface evolution and electromagnetic field in the laser-arc hybrid welding process, realize integrated simulation calculation, more comprehensively and truly reproduce the physical process, facilitate deeper theoretical research and further promote the engineering application of the laser-arc hybrid welding.
Drawings
FIG. 1 is a schematic view of a laser-arc hybrid welding process;
FIG. 2 is a block flow diagram of the present invention;
FIG. 3 is a schematic view of an overall spatial intermittent decomposition with dashed lines indicating the free interface of the molten bath;
FIG. 4 shows the simulation results of the temperature field in example 1 of the present invention, wherein the unit of the numerical quantity is K;
FIG. 5 shows the velocity field simulation results of example 1 of the present invention, in which the numerical quantity is m/s;
FIG. 6 shows the results of the mass percentage simulation of the metal vapor in example 1 of the present invention;
FIG. 7 shows the simulation results of current density distribution in example 1 of the present invention, in which the unit of digital quantity is A/m2
Detailed Description
The invention is further illustrated below with reference to specific embodiments and the accompanying drawings.
Embodiment 1, the flowchart of which is shown in fig. 2, includes a geometric model building step, an initial value setting step, a time step updating step, a physical property parameter updating step, an electromagnetic field solving step, a spatial intermittent decomposition step, a gas region state calculating step, a workpiece region state calculating step, a molten pool free interface updating step, and a judgment step:
firstly, establishing a geometric model:
establishing a geometric model by using three-dimensional modeling software, wherein the geometric model is a cuboid, and is 14mm in length, 10mm in width and 12mm in height; the geometric model comprises a workpiece area, an electrode area and a gas area, wherein the workpiece area is a cuboid, the length and the width of the workpiece area are the same as those of the geometric model, and the height of the workpiece area is 4 mm; the electrode area is a cylinder, the diameter of the cylinder is 2mm, the height of the cylinder is 3mm, the included angle between the central axis of the electrode area and the vertical direction is 30 degrees, and the vertical distance between the end part of the electrode area and the workpiece area is 5 mm; the gas area is the rest area of the geometric model except the workpiece area and the electrode area; the three-dimensional modeling Software is UG of Siemens PLM Software company in Germany;
the region of the geometric model except the electrode region is a calculation region;
carrying out grid division on the calculation region by adopting grid division software, wherein the grid is a cube, and the side length of the cube is 100 micrometers; the grid division software is Hypermesh of Altair company in America;
secondly, setting an initial value:
given the boundary conditions of the computation region: the environment temperature is 300K, the electrode temperature is 3000K, the environment pressure is 1atm, and the bottom potential of the workpiece is 0V;
setting the initial value of the physical quantity of each grid center point in the whole calculation area: potential 0V, current density 0A/m2The magnetic field intensity is 0T, the speed is 0m/s, the temperature is 300K, the pressure is 1atm, and the mass percentage of the metal vapor is 0;
setting an accumulation time variable tsumCalculating the cutoff time t as 0end=1.5s;
Setting the time step Δ tiIs an initial time step Δ t0,Δt0=1×10-7s, performing the fourth step;
step three, updating the time step:
calculating the maximum time step for gas region stability calculation
Figure BDA0001793010580000141
Then calculating the maximum time step length of stable calculation of the workpiece region
Figure BDA0001793010580000142
Take (Δ t)0,Δtg,Δtm) The minimum value of the three is used as the updated time step delta tiCarrying out the step four;
where, i is 1,2 …, which indicates the i-th step, the maximum condition number Cmax=0.8,ug_x,ug_y,ug_zGas velocity at the center point of the gas area grid
Figure BDA0001793010580000143
Component in three directions x, y, z, um_x,um_y,um_zThe velocity of the molten pool being the central point of the grid of the region of the workpiece
Figure BDA0001793010580000144
Components in the x, y, z directions; the sizes of the grids in the x direction, the y direction and the z direction are respectively delta x, delta y and delta z;
fourthly, updating physical property parameters:
setting physical property parameters in the whole calculation region:
the gas area is filled with mixed gas of protective gas and metal steam, the protective gas is argon, and the metal steam is metal steam of a workpiece material; calculating the mass percentage content and the temperature of the metal steam in the mixed gas as initial values for the first time, and acquiring the mass percentage content and the temperature of the metal steam in the mixed gas in the next time step from the substeps (7.2) and (7.3);
conductivity σ of the gas mixtureeDensity rhogViscosity [ mu ] ofgSpecific heat capacity CpgCoefficient of thermal conductivity kgAnd net gas emissivityNThe composition and temperature of the mixed gas at the center of each grid are directly determined; the electrical conductivity σ of the gas mixture being determined by the composition and temperature of the gas mixtureeDensity rhogViscosity [ mu ] ofgSpecific heat capacity CpgCoefficient of thermal conductivity kgAnd net gas emissivityNCited from MougenotJ, et al, plasma-field pore interaction in tubular input-gassconfiguration, journal of Physics D: Applied Physics,2013,46(13): 135206;
the workpiece material of the workpiece area is 316L stainless steel; corresponding conductivity σeDensity rhomViscosity [ mu ] ofmSpecific heat capacity CpmCoefficient of thermal conductivity kmThermal expansion coefficient β, surface tension coefficient gamma, and thermal capillary coefficient
Figure BDA0001793010580000151
Work potential
Figure BDA0001793010580000158
Emissivity of radiationrDirectly determined by the workpiece material and temperature at each grid center; the corresponding electrical conductivity σ of the workpiece material being determined directly from the workpiece material and the temperatureeDensity rhomViscosity [ mu ] ofmSpecific heat capacity CpmCoefficient of thermal conductivity kmThermal expansion coefficient β, surface tension coefficient gamma, and thermal capillary coefficient
Figure BDA0001793010580000152
Work byPotential energy
Figure BDA0001793010580000153
Emissivity of radiationrThe material property calculation Software JMatPro of Sente Software company in England is adopted for calculation and acquisition;
fifthly, solving the electromagnetic field:
solving for the electromagnetic field throughout the calculation area:
welding current is introduced in the welding process, an electromagnetic field is generated in the whole calculation area, the welding current of 150A is given through the electrodes, and the potential of the center points of all grids in the whole calculation area is obtained through calculation
Figure BDA0001793010580000154
Current density
Figure BDA0001793010580000155
And magnetic induction intensity
Figure BDA0001793010580000156
Distribution, the calculation equation is as follows:
Figure BDA0001793010580000157
Figure BDA0001793010580000161
Figure BDA0001793010580000162
Figure BDA0001793010580000163
wherein the vacuum dielectric constant is0Magnetic loss of potential at center point of grid
Figure BDA0001793010580000164
Gradient operator
Figure BDA0001793010580000165
Divergence operator
Figure BDA0001793010580000166
Rotation operator
Figure BDA0001793010580000167
Laplacian operator
Figure BDA0001793010580000168
Electrical conductivity sigmaeRespectively adopting the mixed gas in the gas area and the corresponding conductivity of the workpiece material in the workpiece area during calculation;
sixthly, spatial discontinuous decomposition:
performing spatial intermittent decomposition by taking the upper surface of a workpiece area as an interface in the first calculation and taking an updated free interface of a molten pool as an interface in the subsequent calculation, and dividing the whole calculation area into a gas area and a workpiece area as shown in fig. 3, wherein the gas area is an arc-metal vapor area, and the workpiece area is a workpiece solid-liquid area; defining interface boundary points as intersection points of grid central point connecting lines at two sides of a free interface of the molten pool and the free interface of the molten pool;
seventhly, calculating the gas area state:
the method comprises the following substeps of calculating the flow field, the temperature field and the metal vapor distribution of the gas region in sequence:
(7.1) calculating the gas area flow field: calculating the gas velocity of the center points of all grids in the gas area
Figure BDA0001793010580000169
Gas pressure PgThe calculation formula is as follows:
Figure BDA00017930105800001610
Figure BDA00017930105800001611
wherein
Figure BDA00017930105800001612
Acceleration of gravity;
setting the gas velocity of the boundary point of the gas zone and the workpiece zone at the free interface of the molten bath
Figure BDA0001793010580000171
And gas pressure Pg_sRespectively as follows:
Figure BDA0001793010580000172
wherein the evaporation recoil pressure PrObtaining by adopting a recoil pressure model based on the environmental pressure; the recoil pressure model based on the environmental pressure is as follows: s. Pang, et al. Explosition of duration depth variant laser welding under variable temporal representation. journal of laser applications 2015,27(2): 022007;
(7.2) calculating the gas zone temperature field: calculating the gas temperature T of the central point of all grids in the gas areagThe calculation formula is as follows:
Figure BDA0001793010580000173
wherein k isbBoltzmann constant, e is elementary charge, t is time,
Figure BDA0001793010580000175
a differential operator;
calculating the heat flow q of the boundary point of the interface into the gas region at the free interface of the molten pool of the gas region and the workpiece regiongas
Figure BDA0001793010580000174
Wherein T isgs、TmsAcquiring the temperatures of a gas area and a workpiece area which are boundary points of an interface respectively by adopting a virtual flow method; k is a radical ofeffIs the effective thermal conductivity, k, of the boundary point of the interfaceeff=0.5×(kg(Tgs)+kg(Tms)),kg(Tgs)、kg(Tms) Respectively is the temperature Tgs、TmsThe corresponding gas thermal conductivity; the effective thickness of the boundary layer of the free interface of the molten pool is 0.1 mm; the virtual stream method is as follows: fedkiw R P, et al, A non-oscillotory Eulerian approach to interfaces in multimaterial flows [ J ] method].Journal ofcomputational physics,1999,152(2):457-492;
(7.3) calculating the gas zone metal vapor diffusion: calculating the mass percentage content C of the metal vapor at the central points of all grids in the gas areavaporThe calculation formula is as follows:
Figure BDA0001793010580000181
wherein D is the mass diffusion coefficient of the metal vapor, and the value can be obtained by a second viscosity approximation method; the second viscosity approximation method is as follows: A.B.Murphy.A composition of transactions of differentiation in thermoplasmas. journal of Physics D: Applied Physics,1996,29(7): 1922;
setting the mass percentage content C of metal steam at each interface boundary point at the free interface of a molten pool of a gas area and a workpiece areavapor
Figure BDA0001793010580000182
Eighthly, calculating the state of the workpiece area:
sequentially calculating a flow field and a temperature field of a workpiece region, and comprising the following substeps:
(8.1) calculating the flow field of the workpiece area: calculating the molten pool speed of the central points of all grids in the workpiece area
Figure BDA0001793010580000183
Pressure P of the molten bathmThe calculation formula is as follows:
Figure BDA0001793010580000184
Figure BDA0001793010580000185
wherein the temperature T of the central point of the grid of the workpiece areamTaking the ambient temperature during the first calculation; the next time step is calculated by adopting the calculated value of the last time step (7.2);
reference temperature TrefTaking the melting point temperature of a workpiece material; K. c is the seepage coefficient, relaxation coefficient and liquidus fraction flClosely related:
Figure BDA0001793010580000191
wherein d is a constant closely related to the size of the dendrite spacing, and d is 0.1 mm; fraction of liquid phase flLinear with temperature:
Figure BDA0001793010580000192
wherein, Tl、TsThe liquidus temperature and solidus temperature of the workpiece material, respectively;
calculating the normal pressure P of the molten pool at the boundary point of the interface at the free interface of the molten pool in the gas area and the workpiece areafAnd the tangential force of the molten pool
Figure BDA0001793010580000193
Figure BDA0001793010580000194
Wherein
Figure BDA0001793010580000195
Normal vector and tangent vector of the free interface of the molten pool, kappa is the curvature of the free interface of the molten pool, PeffFor effective pressure of gas region to workpiece region, take Peff=max(Pr,Pg);
(8.2) calculating the temperature field of the workpiece area: calculating the temperature T of the central point of all grids in the workpiece areamThe calculation formula is as follows:
Figure BDA0001793010580000196
wherein the workpiece absorbs the energy of the laser beam qlaserThe laser is obtained by adopting a ray tracing method, the laser is an optical fiber laser, the power is 2000W, and the radius of a light spot is 0.3 mm; the ray tracing method is as follows: computing the energy density of the laser deep fusion welding small holes with arbitrary shapes, the laser technology 2010,34(05) 614-618;
calculating the heat flow q of the boundary point of the interface into the workpiece region at the free interface of the molten pool of the gas region and the workpiece regionmatel
Figure BDA0001793010580000201
Wherein σ is Stefan-Boltzmann constant, TIs the ambient temperature, qevapHeat loss for evaporation; heat loss q by evaporationevapThe calculation method is as follows: wu D, et al, understanding of spray formation in fiber laser welding of5083aluminum alloy&MassTransfer,2017,113:730-740;
Ninthly, updating a free interface of the molten pool:
molten pool velocity through all grid center points of workpiece area
Figure BDA0001793010580000202
Calculating the distance phi from the central point of each grid to the interface in the area, wherein the obtained equivalent surface with phi equal to 0 is the updated free interface of the molten pool, and the calculation formula is as follows:
Figure BDA0001793010580000203
the normal vector and curvature of the free interface of the molten pool are solved as follows:
Figure BDA0001793010580000204
in the fifth step to the ninth step, all the related calculation equations are dispersed by adopting a finite difference method, and the time step delta t isiSolving is carried out internally;
tenthly, judging:
will tsum+ΔtiIs given to tsumJudging whether t is presentsum<tendIf yes, turning to the step three, updating the time step, otherwise, ending the calculation.
After the calculation is finished, an open source visualization program Paraview is adopted to perform visualization processing on simulation data obtained in the simulation process, and fig. 4 to 7 are transient results of the laser-arc hybrid welding temperature field, the speed field, the metal vapor distribution and the current density distribution simulated in the embodiment.
Embodiment 2, the flowchart of which is shown in fig. 1, includes a geometric model building step, an initial value setting step, a time step updating step, a physical property parameter updating step, an electromagnetic field solving step, a spatial intermittent decomposition step, a gas region state calculating step, a workpiece region state calculating step, a molten pool free interface updating step, and a judgment step:
firstly, establishing a geometric model:
same as example 1;
secondly, setting an initial value:
given the boundary conditions of the computation region: the environment temperature is 1000K, the electrode temperature is 3500, the environment pressure is 10atm, and the bottom potential of the workpiece is 0V;
setting the initial value of the physical quantity of each grid center point in the whole calculation area: potential 0V, current density 0A/m2The magnetic field intensity is 0T, the speed is 0m/s, the temperature is 1000K, the pressure is 10atm, and the mass percentage of the metal vapor is 0;
setting an accumulation time variable tsumCalculating the cutoff time t as 0end=10s;
Setting the time step Δ tiIs an initialStep of time Δ t0,Δt0=1×10-6s, performing the fourth step;
step three, updating the time step:
same as example 1;
fourthly, updating physical property parameters:
setting physical property parameters in the whole calculation region:
the gas area is filled with mixed gas of protective gas and metal vapor, the protective gas is helium, and the metal vapor is metal vapor of a workpiece material; calculating the mass percentage content and the temperature of the metal steam in the mixed gas as initial values for the first time, and acquiring the mass percentage content and the temperature of the metal steam in the mixed gas in the next time step from the substeps (7.2) and (7.3);
conductivity σ of the gas mixtureeDensity rhogViscosity [ mu ] ofgSpecific heat capacity CpgCoefficient of thermal conductivity kgAnd net gas emissivityNThe composition and temperature of the mixed gas at the center of each grid are directly determined; the electrical conductivity σ of the gas mixture being determined by the composition and temperature of the gas mixtureeDensity rhogViscosity [ mu ] ofgSpecific heat capacity CpgCoefficient of thermal conductivity kgAnd net gas emissivityNFrom Murphy AB, et al modeling of thermal plastics for arc welding, the roll of the shifting gas properties and of the metal vacuum. journal of Physics D: applied Physics,2009,42(19): 194006;
the workpiece material of the workpiece area is 304 stainless steel; corresponding conductivity σeDensity rhomViscosity [ mu ] ofmSpecific heat capacity CpmCoefficient of thermal conductivity kmThermal expansion coefficient β, surface tension coefficient gamma, and thermal capillary coefficient
Figure BDA0001793010580000221
Work potential
Figure BDA0001793010580000222
Emissivity of radiationrDirectly determined by the workpiece material and temperature at each grid center; from workpiece materialThe electrical conductivity σ corresponding to the material of the workpiece determined directly by the temperatureeDensity rhomViscosity [ mu ] ofmSpecific heat capacity CpmCoefficient of thermal conductivity kmThermal expansion coefficient β, surface tension coefficient gamma, and thermal capillary coefficient
Figure BDA0001793010580000223
Work potential
Figure BDA0001793010580000224
Emissivity of radiationrThe citation is from: wu C S.computer correlation of thread-dimensional correlation in transforming MIG wells engineering computations,1992,9(5): 529-;
fifthly, solving the electromagnetic field:
same as example 1;
sixthly, spatial discontinuous decomposition:
same as example 1;
seventhly, calculating the gas area state:
(7.1) same as in example 1;
(7.2) the same procedure as in example 1 was repeated except that the thickness was changed to 0.4 mm;
(7.3) same as in example 1;
eighthly, calculating the state of the workpiece area:
same as example 1;
ninthly, updating a free interface of the molten pool:
same as example 1;
tenthly, judging:
same as example 1;
after the calculation is finished, the simulation data obtained in the simulation process is visualized by adopting visualization processing software Tecplot 360 of the American Tecplot company, and the transient results of the laser-arc hybrid welding temperature field, the speed field, the metal vapor distribution and the current density distribution can be obtained.
Embodiment 3, the flowchart of which is shown in fig. 1, includes a geometric model building step, an initial value setting step, a time step updating step, a physical property parameter updating step, an electromagnetic field solving step, a spatial intermittent decomposition step, a gas region state calculating step, a workpiece region state calculating step, a molten pool free interface updating step, and a judgment step:
firstly, establishing a geometric model:
establishing a geometric model by using three-dimensional modeling software, wherein the geometric model is a cuboid, and is 500mm long, 400mm wide and 200mm high; the geometric model comprises a workpiece area, an electrode area and a gas area, wherein the workpiece area is a cuboid, the length and the width of the workpiece area are the same as those of the geometric model, and the height of the workpiece area is 40 mm; the electrode area is a cone, the cone angle is 120 degrees, the diameter of the bottom surface of the cone is the same as that of the cylinder, the diameter of the cylinder is 4mm, the height of the cylinder is 150mm, the included angle between the axis in the electrode area and the vertical direction is 60 degrees, and the vertical distance between the end part of the electrode area and the workpiece area is 10 mm; the region of the geometric model except the electrode region is a calculation region; the three-dimensional modeling software is Pro/E of American parameter technology company;
the gas area is the rest area of the geometric model except the workpiece area and the electrode area;
carrying out grid division on the calculation region by adopting grid division software, wherein the grid is a cuboid, and the side length of the cuboid is 500 microns, the width of the cuboid is 200 microns, and the height of the cuboid is 200 microns; the grid division software is Gridgen of the American Pointwise company;
the second to fourth steps are the same as in example 1;
fifthly, solving the electromagnetic field:
welding current is introduced in the welding process, an electromagnetic field is generated in the whole calculation area, and the welding current of 300A is given through the electrode; otherwise, the same as example 1;
the sixth step to the seventh step are the same as those in the embodiment 1;
eighthly, calculating the state of the workpiece area:
(8.1) the same procedure as in example 1 was repeated except that d was 0.5 mm;
(8.2) the laser is Nd, namely YAG laser, the power is 10000W, and the spot radius is 0.5 mm; otherwise, the same as example 1;
nine to ten steps are the same as in example 1;
embodiment 4, the flowchart of which is shown in fig. 1, includes a geometric model building step, an initial value setting step, a time step updating step, a physical property parameter updating step, an electromagnetic field solving step, a spatial intermittent decomposition step, a gas region state calculating step, a workpiece region state calculating step, a molten pool free interface updating step, and a judgment step:
firstly, establishing a geometric model:
the geometric model is a cuboid, the length is 5mm, the width is 5mm, and the height is 4 mm; the height of the workpiece area is 1 mm; the electrode area is a cone, the cone angle is 30 degrees, the diameter of the bottom surface of the cone is the same as that of the cylinder, the diameter of the cylinder is 1mm, the height of the cylinder is 1mm, the included angle between the central axis of the electrode area and the vertical direction is 60 degrees, and the vertical distance between the end part of the electrode area and the workpiece area is 2 mm;
the grid is a cube, and the side length of the cube is 20 micrometers;
otherwise, the same as example 1;
secondly, setting an initial value:
calculating the cut-off time tend0.1 s; otherwise, the same as example 1;
the third to fourth steps are the same as in example 1;
fifthly, solving the electromagnetic field:
welding current is introduced in the welding process, an electromagnetic field is generated in the whole calculation area, and 30A welding current is given through the electrodes; otherwise, the same as example 1;
the sixth step to the seventh step are the same as those in the embodiment 1;
eighthly, calculating the state of the workpiece area:
(8.1) taking d as 0.05mm, and the rest is the same as the example 1;
(8.2) the laser is a fiber laser, the power is 100W, and the radius of a light spot is 0.1 mm; otherwise, the same as example 1;
nine to ten steps are the same as in example 1;
the method realizes the numerical simulation calculation of the laser-TIG hybrid welding under the condition of the process parameters of the embodiment, and can simultaneously solve the evolution law of solid phase, liquid phase, metal vapor and plasma in the hybrid welding process. By simply modifying the embodiment and the implementation mode in the implementation scheme, the numerical simulation of the laser-TIG hybrid welding of different shielding gas types, different workpiece materials and different process parameters can be realized by adopting the method disclosed by the invention, the expansibility is good, and the method has application potential in the aspects of mechanism research of the hybrid welding and engineering process parameter formulation.
The above description is only a preferred embodiment of the present invention, but the scope of the present invention is not limited thereto, and any changes or substitutions that can be easily conceived by those skilled in the art within the technical scope of the present invention are included in the scope of the present invention. Therefore, the protection scope of the present invention shall be subject to the protection scope of the claims.

Claims (5)

1. A three-dimensional transient numerical simulation method for laser-arc hybrid welding comprises the steps of establishing a geometric model, setting an initial value, updating a time step, updating physical property parameters, solving an electromagnetic field, performing spatial discontinuous decomposition, calculating a gas region state, calculating a workpiece region state, updating a molten pool free interface and judging, and is characterized in that:
firstly, establishing a geometric model:
establishing a geometric model by adopting three-dimensional modeling software, wherein the geometric model is a cuboid, the geometric model comprises a workpiece area, an electrode area and a gas area, and the workpiece area is the cuboid and has the same length and width as the geometric model; the electrode area is a cylinder or a cylinder with one end being a cone, when the electrode area is a cylinder, the included angle between the central axis of the electrode area and the vertical direction is 0-60 degrees, and when the electrode area is a cylinder with one end being a cone, the diameter of the bottom surface of the cone is the same as that of the cylinder; the vertical distance between the end part of the electrode area and the workpiece area is 1 mm-10 mm; the gas area is the rest area of the geometric model except the workpiece area and the electrode area;
the region of the geometric model except the electrode region is a calculation region;
carrying out grid division on the calculation region by adopting grid division software, wherein the grid is a cuboid or a cube, and the side length of the cuboid or the cube is 10-500 mu m;
secondly, setting an initial value:
given the boundary conditions of the computation region: the environment temperature is 300-1000K, the electrode temperature is 3000-3500K, the environment pressure is 0.1-10 atm, and the bottom potential of the workpiece is 0V;
setting the initial value of the physical quantity of each grid center point in the whole calculation area: potential 0V, current density 0A/m2The magnetic field intensity is 0T, the speed is 0m/s, the temperature is 300-1000K, the pressure is 0.1-10 atm, and the mass percentage content of the metal steam is 0;
setting an accumulation time variable tsumCalculating the cutoff time t as 0end=0.1~10s;
Setting the time step Δ tiIs an initial time step Δ t0,Δt0=1×10-7~1×10-6s, performing the fourth step;
step three, updating the time step:
calculating the maximum time step for gas region stability calculation
Figure FDA0001793010570000021
Then calculating the maximum time step length of stable calculation of the workpiece region
Figure FDA0001793010570000022
Take (Δ t)0,Δtg,Δtm) The minimum value of the three is used as the updated time step delta tiCarrying out the step four;
where, i is 1,2 …, which indicates the i-th step, the maximum condition number Cmax=0.5~1.0,ug_x,ug_y,ug_zGas velocity at the center point of the gas area grid
Figure FDA0001793010570000023
At three of x, y and zComponent of direction, um_x,um_y,um_zThe velocity of the molten pool being the central point of the grid of the region of the workpiece
Figure FDA0001793010570000024
Components in the x, y, z directions; the sizes of the grids in the x direction, the y direction and the z direction are respectively delta x, delta y and delta z;
fourthly, updating physical property parameters:
setting physical property parameters in the whole calculation region:
the gas area is filled with mixed gas of protective gas and metal steam, the protective gas is one or more of argon, helium and carbon dioxide, and the metal steam is the metal steam of the workpiece material; calculating the mass percentage content and the temperature of the metal steam in the mixed gas as initial values for the first time, and acquiring the mass percentage content and the temperature of the metal steam in the mixed gas in the next time step from the substeps (7.2) and (7.3);
conductivity σ of the gas mixtureeDensity rhogViscosity [ mu ] ofgSpecific heat capacity CpgCoefficient of thermal conductivity kgAnd net gas emissivityNThe composition and temperature of the mixed gas at the center of each grid are directly determined;
the workpiece material of the workpiece area is iron-based alloy or non-ferrous metal alloy; corresponding conductivity σeDensity rhomViscosity [ mu ] ofmSpecific heat capacity CpmCoefficient of thermal conductivity kmThermal expansion coefficient β, surface tension coefficient gamma, and thermal capillary coefficient
Figure FDA0001793010570000031
Work potential
Figure FDA0001793010570000032
Emissivity of radiationrDirectly determined by the workpiece material and temperature at each grid center;
fifthly, solving the electromagnetic field:
solving for the electromagnetic field throughout the calculation area:
in the welding processWelding current is input, an electromagnetic field is generated in the whole calculation area, welding current of 30-300A is given through the electrodes, and the potential of the center points of all grids in the whole calculation area is obtained through calculation
Figure FDA0001793010570000033
Current density
Figure FDA0001793010570000034
And magnetic induction intensity
Figure FDA0001793010570000035
Distribution, the calculation equation is as follows:
Figure FDA0001793010570000036
Figure FDA0001793010570000037
Figure FDA0001793010570000038
Figure FDA0001793010570000039
wherein the vacuum dielectric constant is0Magnetic loss of potential at center point of grid
Figure FDA00017930105700000310
Gradient operator
Figure FDA00017930105700000311
Divergence operator
Figure FDA00017930105700000312
Rotation operator
Figure FDA00017930105700000313
Laplacian operator
Figure FDA00017930105700000314
Electrical conductivity sigmaeRespectively adopting the mixed gas in the gas area and the corresponding conductivity of the workpiece material in the workpiece area during calculation;
sixthly, spatial discontinuous decomposition:
performing space discontinuous decomposition by taking the upper surface of a workpiece area as an interface in the first calculation and taking an updated free interface of a molten pool as an interface in the subsequent calculation, and dividing the whole calculation area into a gas area and a workpiece area, wherein the gas area is an arc-metal vapor area and the workpiece area is a workpiece solid-liquid area; defining interface boundary points as intersection points of grid central point connecting lines at two sides of a free interface of the molten pool and the free interface of the molten pool;
seventhly, calculating the gas area state:
the method comprises the following substeps of calculating the flow field, the temperature field and the metal vapor distribution of the gas region in sequence:
(7.1) calculating the gas area flow field: calculating the gas velocity of the center points of all grids in the gas area
Figure FDA0001793010570000041
Gas pressure PgThe calculation formula is as follows:
Figure FDA0001793010570000042
Figure FDA0001793010570000043
wherein
Figure FDA0001793010570000044
Is the acceleration of gravity;
setting the gas velocity of the boundary point of the gas zone and the workpiece zone at the free interface of the molten bath
Figure FDA0001793010570000048
And gas pressure Pg_sRespectively as follows:
Figure FDA0001793010570000045
wherein the evaporation recoil pressure PrObtaining by adopting a recoil pressure model based on the environmental pressure;
(7.2) calculating the gas zone temperature field: calculating the gas temperature T of the central point of all grids in the gas areagThe calculation formula is as follows:
Figure FDA0001793010570000046
wherein k isbBoltzmann constant, e is elementary charge, t is time,
Figure FDA0001793010570000047
a differential operator;
calculating the heat flow q of the boundary point of the interface into the gas region at the free interface of the molten pool of the gas region and the workpiece regiongas
Figure FDA0001793010570000051
Wherein T isgs、TmsAcquiring the temperatures of a gas area and a workpiece area which are boundary points of an interface respectively by adopting a virtual flow method; k is a radical ofeffIs the effective thermal conductivity, k, of the boundary point of the interfaceeff=0.5×(kg(Tgs)+kg(Tms)),kg(Tgs)、kg(Tms) Respectively is the temperature Tgs、TmsThe corresponding gas thermal conductivity; the effective thickness of the boundary layer of the free interface of the molten pool is 0.1 mm-0.4 mm;
(7.3) calculating the gas zone metal vapor diffusion: calculating the mass percentage content C of the metal vapor at the central points of all grids in the gas areavaporMeter for measuringThe calculation formula is as follows:
Figure FDA0001793010570000052
wherein D is the mass diffusion coefficient of the metal vapor, and the value can be obtained by a second viscosity approximation method;
setting the mass percentage content C of metal steam at each interface boundary point at the free interface of a molten pool of a gas area and a workpiece areavapor
Figure FDA0001793010570000053
Eighthly, calculating the state of the workpiece area:
sequentially calculating a flow field and a temperature field of a workpiece region, and comprising the following substeps:
(8.1) calculating the flow field of the workpiece area: calculating the molten pool speed of the central points of all grids in the workpiece area
Figure FDA0001793010570000054
Pressure P of the molten bathmThe calculation formula is as follows:
Figure FDA0001793010570000055
Figure FDA0001793010570000061
wherein the temperature T of the central point of the grid of the workpiece areamTaking the ambient temperature during the first calculation; the next time step is calculated by adopting the calculated value of the last time step (7.2);
reference temperature TrefTaking the melting point temperature of a workpiece material; K. c is the seepage coefficient, relaxation coefficient and liquidus fraction flClosely related:
Figure FDA0001793010570000067
wherein d is a constant closely related to the size of dendrite spacing, and d is 0.05mm to 0.5 mm; fraction of liquid phase flLinear with temperature:
Figure FDA0001793010570000062
wherein, Tl、TsThe liquidus temperature and solidus temperature of the workpiece material, respectively;
calculating the normal pressure P of the molten pool at the boundary point of the interface at the free interface of the molten pool in the gas area and the workpiece areafAnd the tangential force of the molten pool
Figure FDA0001793010570000063
Figure FDA0001793010570000064
Wherein
Figure FDA0001793010570000065
Normal vector and tangent vector of the free interface of the molten pool, kappa is the curvature of the free interface of the molten pool, PeffFor effective pressure of gas region to workpiece region, take Peff=max(Pr,Pg);
(8.2) calculating the temperature field of the workpiece area: calculating the temperature T of the central point of all grids in the workpiece areamThe calculation formula is as follows:
Figure FDA0001793010570000066
wherein the workpiece absorbs the energy of the laser beam qlaserAcquiring by adopting a ray tracing method; the laser is fiber laser or Nd, namely YAG laser, the power is 100W-10000W, and the spot radius is 0.1 mm-0.5 mm;
calculating the boundary point of the interface into the workpiece zone at the free interface of the molten pool between the gas zone and the workpiece zoneOf heat flow qmatel
Figure FDA0001793010570000071
Wherein σ is Stefan-Boltzmann constant, TIs the ambient temperature, qevapHeat loss for evaporation;
ninthly, updating a free interface of the molten pool:
molten pool velocity through all grid center points of workpiece area
Figure FDA0001793010570000072
Calculating the distance phi from the central point of each grid to the interface in the area, wherein the obtained equivalent surface with phi equal to 0 is the updated free interface of the molten pool, and the calculation formula is as follows:
Figure FDA0001793010570000073
the normal vector and curvature of the free interface of the molten pool are solved as follows:
Figure FDA0001793010570000074
in the fifth step to the ninth step, all the involved calculation equations are dispersed by adopting a finite difference method or a finite volume method, and the time step delta t isiSolving is carried out internally;
tenthly, judging:
will tsum+ΔtiIs given to tsumJudging whether t is presentsum<tendIf yes, turning to the step three, updating the time step, otherwise, ending the calculation.
2. The laser-arc hybrid welding three-dimensional transient numerical simulation method of claim 1, characterized in that: in the first step, the three-dimensional modeling software is AutoCAD of American Autodesk software company, UG of Germany Siemens PLMSofware company or Pro/E of American parameter technology company; the cuboid of the geometric model is 5 mm-500 mm long, 5 mm-500 mm wide and 5 mm-200 mm high; the cuboid of the workpiece area is 1 mm-100 mm high; when the electrode area is a cylinder, the diameter of the cylinder is 1 mm-4 mm, the height is 1 mm-150 mm, when the electrode area is a cylinder with one end being a cone, the cone angle is 30-120 degrees, the diameter of the cylinder is 1 mm-4 mm, and the height is 1 mm-150 mm;
the gridding software is Hypermesh of Altair corporation, Gridgen of Pointwise corporation, or Patran of MSC.
3. The laser-arc hybrid welding three-dimensional transient numerical simulation method of claim 1, characterized in that: in the fourth step, the electric conductivity σ of the mixed gas is determined by the composition and temperature of the mixed gaseDensity rhogViscosity [ mu ] ofgSpecific heat capacity CpgCoefficient of thermal conductivity kgAnd net gas emissivityNReference is made to the substances from Mougenot J, et al, plasma-field poolvination in structured input-gas configuration, journal of Physics D: Applied Physics,2013,46(13):135206, or Murphy A B, et al, modeling of thermal plasma for arc welding: the roll of the shifting gas properties and of metallic publication. journal of Physics D: Applied Physics,2009,42(19):194006, or Murphy AB.the effects of metallic publication in arc welding. journal of Physics D: Applied Physics,2010,43(43): 434001;
the corresponding electrical conductivity σ of the workpiece material being determined directly from the workpiece material and the temperatureeDensity rhomViscosity [ mu ] ofmSpecific heat capacity CpmCoefficient of thermal conductivity kmThermal expansion coefficient β, surface tension coefficient gamma, and thermal capillary coefficient
Figure FDA0001793010570000091
Work potential
Figure FDA0001793010570000092
Emissivity of radiationrThe material property calculation Software JMatPro of Sente Software company in England is adopted for calculationObtained from either Wei H L, the equivalent, the original, the orientation reducing, the identification of an aluminum alloy, acta materials, 2016,115: 123-.
4. The laser-arc hybrid welding three-dimensional transient numerical simulation method of claim 1, characterized in that: in the sub-step (7.1), the recoil pressure model based on the environmental pressure is as follows: s. Pang, et al. exhibition of duration reduction Laser welding under variable electromagnetic Applications. journal of Laser Applications 2015,27(2): 022007; in the sub-step (7.2), the virtual stream method is as follows: fedkiw R P, et al, A non-oscillotory Eulerian approach in multimaterial flows, Journal of computational physics,1999,152(2) 457-; in the sub-step (7.3), the second viscosity approximation method is as follows: A.B.Murphy.A composition of transactions of differentiation in thermoplasmas. journal of Physics D: Applied Physics,1996,29(7): 1922.
5. The laser-arc hybrid welding three-dimensional transient numerical simulation method of claim 1, characterized in that: in the sub-step (8.2), the ray tracing method is as follows: computing the energy density of the laser deep fusion welding small holes with arbitrary shapes, the laser technology 2010,34(05) 614-618; in the step (8.2), the heat loss q is evaporatedevapThe calculation method is as follows: wu D, et al.infraduction of scatter formation in fiber laser welding of5083aluminum alloy&Mass Transfer,2017,113:730-740。
CN201811044792.7A 2018-09-07 2018-09-07 Laser-arc hybrid welding three-dimensional transient numerical simulation method Active CN109190260B (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201811044792.7A CN109190260B (en) 2018-09-07 2018-09-07 Laser-arc hybrid welding three-dimensional transient numerical simulation method

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201811044792.7A CN109190260B (en) 2018-09-07 2018-09-07 Laser-arc hybrid welding three-dimensional transient numerical simulation method

Publications (2)

Publication Number Publication Date
CN109190260A CN109190260A (en) 2019-01-11
CN109190260B true CN109190260B (en) 2020-08-18

Family

ID=64915480

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201811044792.7A Active CN109190260B (en) 2018-09-07 2018-09-07 Laser-arc hybrid welding three-dimensional transient numerical simulation method

Country Status (1)

Country Link
CN (1) CN109190260B (en)

Families Citing this family (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN110625229A (en) * 2019-08-05 2019-12-31 江苏大学 Method for calculating minimum consumption of welding shielding gas
CN112214863B (en) * 2020-07-31 2024-01-02 中国科学院金属研究所 Method for acquiring characteristic thermal cycle curve of coarse-grain region of austenitic stainless steel welding heat affected zone
CN113221412B (en) * 2021-05-08 2023-03-21 桂林理工大学 Method, system, terminal and medium for calculating charging potential data of lossy medium
CN114523205B (en) * 2022-02-11 2023-11-10 东北电力大学 Dynamic tracking method for pasty area of magnesium alloy Nd-YAG pulse laser spot welding molten pool
CN116713560B (en) * 2023-08-10 2023-10-20 中建安装集团西安建设投资有限公司 Nondestructive testing control welding method based on digital model
CN117195663B (en) * 2023-11-03 2024-02-20 山东理工大学 Simulation method for removing electric spark machining materials in liquid based on three-phase flow interface tracking

Family Cites Families (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
DE4321713C2 (en) * 1992-07-07 1994-08-25 Ford Werke Ag Composite disc brake rotor and method for its manufacture
CN100411799C (en) * 2005-10-20 2008-08-20 武汉理工大学 Magnetic control melting electrode welding method, and its developed application, and its universal equipment
CN100551603C (en) * 2007-04-28 2009-10-21 重庆大学 A kind of consumable electrode surfacing method of electromagnetic complex field, equipment and expansion thereof are used
CN101462196B (en) * 2008-12-31 2011-02-09 重庆大学 Method and equipment for electromagnetic composite double-face submerged arc welding of diphase stainless steel thick plate
CN103737176B (en) * 2013-12-30 2016-08-31 华中科技大学 A kind of laser and electromagnetic pulse complex welding method and equipment
CN104708177B (en) * 2015-01-28 2017-03-01 辽宁工程技术大学 A kind of measurement apparatus of tungsten electrode argon-arc welding electric arc power and measuring method
CN106709176A (en) * 2016-11-29 2017-05-24 中国航空工业集团公司沈阳飞机设计研究所 Dynamic numerical simulation technology for laser melting deposition formed molten pool

Also Published As

Publication number Publication date
CN109190260A (en) 2019-01-11

Similar Documents

Publication Publication Date Title
CN109190260B (en) Laser-arc hybrid welding three-dimensional transient numerical simulation method
Fan et al. Numerical modeling of the additive manufacturing (AM) processes of titanium alloy
Koepf et al. 3D multi-layer grain structure simulation of powder bed fusion additive manufacturing
Zhang et al. 3-Dimensional heat transfer modeling for laser powder-bed fusion additive manufacturing with volumetric heat sources based on varied thermal conductivity and absorptivity
Lopez-Botello et al. Two-dimensional simulation of grain structure growth within selective laser melted AA-2024
Michaleris Modeling metal deposition in heat transfer analyses of additive manufacturing processes
Zinoviev et al. Evolution of grain structure during laser additive manufacturing. Simulation by a cellular automata method
Lee et al. Modeling of heat transfer, fluid flow and solidification microstructure of nickel-base superalloy fabricated by laser powder bed fusion
Raghavan et al. Numerical modeling of heat-transfer and the influence of process parameters on tailoring the grain morphology of IN718 in electron beam additive manufacturing
Gürtler et al. Influence of powder distribution on process stability in laser beam melting: Analysis of melt pool dynamics by numerical simulations
Liu et al. Micro scale 3D FEM simulation on thermal evolution within the porous structure in selective laser sintering
Fallah et al. Phase-field simulation of solidification morphology in laser powder deposition of Ti–Nb alloys
Shahabad et al. Heat source model calibration for thermal analysis of laser powder-bed fusion
Kumar et al. Three-dimensional finite element modeling of gas metal-arc welding
Veyhl et al. On the thermal conductivity of sintered metallic fibre structures
Geng et al. Modelling and experimental observation of the deposition geometry and microstructure evolution of aluminum alloy fabricated by wire-arc additive manufacturing
Ding et al. The well-distributed volumetric heat source model for numerical simulation of wire arc additive manufacturing process
Tong et al. A dynamic welding heat source model in pulsed current gas tungsten arc welding
Wu et al. The virtual node polygonal element method for nonlinear thermal analysis with application to hybrid laser welding
Vrána et al. Contour laser strategy and its benefits for lattice structure manufacturing by selective laser melting technology
Long et al. Numerical simulation of thermal behavior during laser metal deposition shaping
Philo et al. A pragmatic continuum level model for the prediction of the onset of keyholing in laser powder bed fusion
Unni et al. Numerical simulation of the influence of oxygen content on the weld pool depth during activated TIG welding
CN115238605A (en) Numerical simulation method for predicting surface quality of SLM (Selective laser melting) single melting channel
Darabi et al. Thermal study of a cladding layer of Inconel 625 in Directed Energy Deposition (DED) process using a phase-field model

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant