CN114626263A - High-temperature alloy material short crack propagation numerical simulation method based on crystal plasticity - Google Patents
High-temperature alloy material short crack propagation numerical simulation method based on crystal plasticity Download PDFInfo
- Publication number
- CN114626263A CN114626263A CN202210180332.7A CN202210180332A CN114626263A CN 114626263 A CN114626263 A CN 114626263A CN 202210180332 A CN202210180332 A CN 202210180332A CN 114626263 A CN114626263 A CN 114626263A
- Authority
- CN
- China
- Prior art keywords
- crystal
- crack
- model
- crack propagation
- shear strain
- 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.)
- Granted
Links
- 239000013078 crystal Substances 0.000 title claims abstract description 110
- 238000000034 method Methods 0.000 title claims abstract description 60
- 239000000956 alloy Substances 0.000 title claims abstract description 28
- 238000004088 simulation Methods 0.000 title claims abstract description 24
- 239000000463 material Substances 0.000 claims abstract description 94
- 229920003023 plastic Polymers 0.000 claims abstract description 56
- 238000012360 testing method Methods 0.000 claims abstract description 28
- 230000008569 process Effects 0.000 claims abstract description 13
- 238000011065 in-situ storage Methods 0.000 claims abstract description 5
- 238000000638 solvent extraction Methods 0.000 claims abstract description 4
- 230000009191 jumping Effects 0.000 claims abstract 3
- 238000005336 cracking Methods 0.000 claims description 32
- 230000001186 cumulative effect Effects 0.000 claims description 26
- 230000010354 integration Effects 0.000 claims description 17
- 238000001514 detection method Methods 0.000 claims description 13
- 229920006395 saturated elastomer Polymers 0.000 claims description 8
- 230000000694 effects Effects 0.000 claims description 7
- 230000007246 mechanism Effects 0.000 claims description 7
- 229920001887 crystalline plastic Polymers 0.000 claims description 6
- 230000000977 initiatory effect Effects 0.000 claims description 5
- 230000000737 periodic effect Effects 0.000 claims description 5
- 238000012545 processing Methods 0.000 claims description 4
- 230000035945 sensitivity Effects 0.000 claims description 4
- 238000012800 visualization Methods 0.000 claims description 4
- 238000012512 characterization method Methods 0.000 claims description 3
- 238000006243 chemical reaction Methods 0.000 claims description 3
- 238000000354 decomposition reaction Methods 0.000 claims description 3
- 230000006870 function Effects 0.000 claims description 3
- 238000005728 strengthening Methods 0.000 claims description 3
- 239000000126 substance Substances 0.000 claims description 3
- 238000009825 accumulation Methods 0.000 claims 1
- 238000012805 post-processing Methods 0.000 claims 1
- 230000035882 stress Effects 0.000 description 11
- 229910000601 superalloy Inorganic materials 0.000 description 6
- 230000006399 behavior Effects 0.000 description 5
- 238000004364 calculation method Methods 0.000 description 5
- 102100021503 ATP-binding cassette sub-family B member 6 Human genes 0.000 description 2
- 101100000375 Homo sapiens ABCB6 gene Proteins 0.000 description 2
- 125000004122 cyclic group Chemical group 0.000 description 2
- 238000010894 electron beam technology Methods 0.000 description 2
- 238000009661 fatigue test Methods 0.000 description 2
- 238000011160 research Methods 0.000 description 2
- 230000009471 action Effects 0.000 description 1
- 238000004458 analytical method Methods 0.000 description 1
- 239000003086 colorant Substances 0.000 description 1
- 238000005520 cutting process Methods 0.000 description 1
- 230000007547 defect Effects 0.000 description 1
- 230000007812 deficiency Effects 0.000 description 1
- 238000013461 design Methods 0.000 description 1
- 238000011161 development Methods 0.000 description 1
- 238000010586 diagram Methods 0.000 description 1
- 238000001887 electron backscatter diffraction Methods 0.000 description 1
- 238000011156 evaluation Methods 0.000 description 1
- 230000006872 improvement Effects 0.000 description 1
- 238000011835 investigation Methods 0.000 description 1
- 238000012986 modification Methods 0.000 description 1
- 230000004048 modification Effects 0.000 description 1
- 238000006467 substitution reaction Methods 0.000 description 1
- 230000008646 thermal stress Effects 0.000 description 1
- 238000012546 transfer Methods 0.000 description 1
Images
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F30/00—Computer-aided design [CAD]
- G06F30/20—Design optimisation, verification or simulation
- G06F30/23—Design optimisation, verification or simulation using finite element methods [FEM] or finite difference methods [FDM]
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F2111/00—Details relating to CAD techniques
- G06F2111/10—Numerical modelling
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F2119/00—Details relating to the type or aim of the analysis or the optimisation
- G06F2119/14—Force analysis or force optimisation, e.g. static or dynamic forces
Landscapes
- Engineering & Computer Science (AREA)
- Physics & Mathematics (AREA)
- Theoretical Computer Science (AREA)
- Computer Hardware Design (AREA)
- Evolutionary Computation (AREA)
- Geometry (AREA)
- General Engineering & Computer Science (AREA)
- General Physics & Mathematics (AREA)
- Investigating Strength Of Materials By Application Of Mechanical Stress (AREA)
Abstract
The invention relates to a crystal plasticity-based high-temperature alloy material short crack propagation numerical simulation method, which comprises the following steps of: (1) establishing a short crack propagation finite element model based on ABAQUS software and based on a short crack propagation in-situ test of a high-temperature alloy material, and partitioning the model; (2) based on the uniaxial tensile stress-strain curve of the material, the stress-strain curve predicted by the model is accurately matched with a uniaxial tensile result; (3) deriving the grid unit and node information of the model, and completing the microstructure modeling of the model crystal plastic material area; (4) introducing a crystal plasticity theory into an expansion finite element method, customizing a short crack expansion criterion, and simulating a short crack expansion behavior; (5) and simulating the short crack propagation process by adopting a quasi-periodic jumping method. The method accurately characterizes the microstructure of the high-temperature alloy material, introduces the microstructure into a short crack propagation model, focuses on the short crack propagation characteristics, and improves the accuracy of the numerical simulation method of the short crack propagation of the high-temperature alloy material.
Description
Technical Field
The invention belongs to the technical field of aerospace engines, and particularly relates to a crystal plasticity-based numerical simulation method for short crack propagation of a high-temperature alloy material.
Background
The turbine disc is used as a key important part of an aircraft engine, works in a complex environment with high temperature and high rotating speed for a long time, and bears large centrifugal load and thermal stress. With the improvement of the performance of the engine, the working conditions of the turbine disc are more severe, which becomes a great problem restricting the development of the aircraft engine. High temperature alloy materials are the material of choice for advanced engine turbine disks and fatigue fracture is one of the major failure modes.
In the structural fatigue failure process, the crack propagation can be divided into two stages, namely a short crack propagation stage and a long crack propagation stage, wherein the length of the short crack is in the range of 10-20 crystal grain sizes, and the propagation life of the short crack accounts for more than 70% of the fatigue life, so that the accurate evaluation of the short crack propagation stage is very important for the analysis of the structural fatigue life.
The crack length and the grain size in the short crack propagation stage are in the same order, and the propagation process is greatly influenced by the microstructure of the crack tip. And in a longer crack propagation stage, a short crack propagation path is tortuous, the crack can also propagate when the stress intensity factor is lower than a macroscopic threshold value, the propagation rate has volatility, and the propagation rate is obviously reduced under the action of obstruction when the crack passes through a grain boundary.
The traditional crack propagation simulation method is mostly based on linear elastic fracture mechanics, is not suitable for short crack propagation behaviors aiming at the conventional crack propagation stage expansion, and cannot explain the fluctuation characteristics of the deflection and the propagation rate of a short crack propagation path, so that the research of the high-temperature alloy material short crack propagation numerical simulation method based on the crystal plasticity theory is innovative from the physical mechanism level by focusing on the short crack propagation characteristics and considering the crystal level influence.
Disclosure of Invention
In order to overcome the defects of the prior art, the invention provides a numerical simulation method for the short crack propagation of the high-temperature alloy material, which is based on the crystal plasticity theory, can consider the influence of the microstructure of the material on the short crack propagation, realizes the simulation of the short crack propagation behavior of the high-temperature alloy material, and lays a foundation for the subsequent fatigue life prediction and damage tolerance design.
The technical scheme of the invention is as follows: a numerical simulation method for the short crack propagation of a high-temperature alloy material based on crystal plasticity is characterized by writing an MATLAB script program based on an Electron Back Scattering Diffraction (EBSD) detection result of the high-temperature alloy material, representing the shape and crystal orientation of material grains, simulating the short crack propagation by adopting an expansion finite element method, customizing a short crack propagation criterion from a physical mechanism level based on a crystal plasticity theory, and obtaining the short crack propagation rate and the propagation direction.
The method comprises the following implementation steps:
firstly, establishing a short crack propagation finite element model of the high-temperature alloy material based on ABAQUS software. The finite element model is based on the size of a test piece in a high-temperature alloy material in-situ fatigue test, and only an assessment section part with the minimum middle width of the test piece is reserved by the finite element model in combination with the Saint-Venn principle. And partitioning the model, dividing the vicinity of a crack initiation region into a crystalline plastic material region, and dividing the rest part into an isotropic material region. The model grid unit type adopts 8-node hexahedron linear units (C3D8) to subdivide grid units in the crystal plastic material area and establish local fine grids.
And secondly, fitting the parameters of the crystal plastic material by adopting a representative volume unit model of the polycrystalline material. The strengthening effect of the deformation rate is considered in the crystal plastic material, and the rule for controlling the plastic shear strain rate of a single slip system in the model is as follows:
wherein the content of the first and second substances,the plastic shear strain rate of the slip system alpha,to reference the plastic shear strain rate, ταThe breaking shear stress of slip system alpha, gαAs the current strength of the slip system alpha, hαβIn order to harden the modulus,plastic shear strain rate of slip system beta, NtotalFor the total number of slip trains, n is the rate sensitivity parameter for the slip shear rate, sgn (τ)α) As a function of the sign:
the hardening model used for the model was as follows:
hαβ=qh(γ)(α≠β)
wherein h isααDescribed is the effect of self-hardening, in which the slip system α ═ β, hαβDescribed is the influence of latent hardening, where the slip system α ≠ β, q characterizes the relationship between latent hardening and self-hardening behaviour of the material, h0、hs、τsRespectively is initial hardening modulus, saturated hardening modulus and saturated shear stress, and gamma is the accumulated shear strain of all slip systems;
introducing a representative volume unit model which accords with the statistical information description of the microstructure of the material, and fitting an orthotropic elastic constant and a slip system hardening parameter in a crystal plastic constitutive model by a uniaxial tensile stress-strain curve of the material to ensure that the stress-strain curve predicted by simulation is accurately matched with a uniaxial tensile result and reflect the mechanical property of the polycrystalline material;
and endowing the fitted crystal plastic material attribute to the region near the crack initiation of the model, and endowing the isotropic material attribute to the rest regions.
And thirdly, finishing the microstructure modeling of the crystal plastic material area of the short crack propagation finite element model. And deriving the grid unit number, the node number in the unit and the global coordinate of each node of the crystal plastic material region of the finite element model in the first step, and writing an MATLAB script program to represent the microstructure of the test piece based on the electron backscatter diffraction (EBSD) detection result of the test piece. The EBSD detection result of the test piece can obtain the position coordinates and Euler angle orientation information of a scanning point, wherein the scanning point refers to a corresponding position point of a sample obtained when an incident electron beam continuously moves to scan a specified area of the sample in the detection process;
the MATLAB script program processing steps for characterizing the microstructure of the test piece include:
comparing the position coordinates of the scanning points in the EBSD detection result with the information of the grid cells in the model crystal plastic material area, and grouping the scanning points according to whether each scanning point belongs to the same cell to obtain the scanning points contained in each cell and the corresponding Euler angle orientation thereofφ、
Combining the similar Euler angles in the grouped units, determining the orientation with the maximum number of scanning points in the unit as the unit dominant orientation, and determining the unit dominant orientationφ、Converting the crystal orientation g into a crystal orientation g in an hkl crystal coordinate system, wherein the hkl crystal coordinate system is a face-centered cubic crystal triaxial coordinate system, h, k and l are three crystal axes of the crystal, and the crystal orientation conversion expression is as follows:
the processed sizes of the test pieces have dispersity, the geometric shapes of the test pieces are not completely consistent with the model, when scanning points with similar Euler angles are combined, non-oriented 'empty' units inevitably appear, the positions of the 'empty' units in the crystalline plastic material area are judged, and the crystal orientation closest to the geometric positions of the 'empty' units is endowed to the 'empty' units;
comparing the crystal orientations of different units, combining the units with the orientation difference smaller than 10 degrees in the adjacent units to form the same crystal grain, obtaining the appearance and the corresponding crystal orientation of the crystal grain, and storing the related information in a python language form which can be identified by ABAQUS commercial finite element software, thereby finishing the microstructure characterization of the high-temperature alloy material;
and (3) importing the crystal grain morphology and crystal orientation information obtained by running an MATLAB script program into an INP file of the short crack propagation model to complete the microstructure modeling of the model crystal plastic material region.
And fourthly, introducing a crystal plasticity theory into an expansion finite element method, and customizing a crack expansion criterion reflecting a short crack expansion physical mechanism. The criterion for controlling the short crack propagation is self-defined by adopting a first subprogram of ABAQUS commercial finite element software, the crack propagation of the unit is controlled by the accumulated shear strain of a single slip system of the enrichment unit in front of the crack tip, and the propagation direction of the short crack is vertical to the normal direction of a slip plane where the maximum accumulated shear strain of the single slip system is located;
the cumulative shear strain of the single slip system of the unit is calculated by a second subprogram of ABAQUS commercial finite element software and is stored in a state variable statev, and the cumulative shear strain of the single slip system is:
wherein t is time, γαCumulative shear strain of slip system alpha,The plastic shear strain rate of the slip system α.
The first subprogram of the self-defining crack propagation criterion only acts on the pre-crack-tip enrichment unit, and the accumulated shear strain of the single slip system calculated and stored by the second subprogram can be called in the first subprogram, so as to judge the cracking time and the crack propagation direction of the crack in the unit. The conditions for crack propagation of the enrichment unit are as follows:
wherein, γαIs the accumulated shear strain of the slip system alpha,plastic shear strain rate, N, of slip system alphatotalIs the total number of slip trains, t is time, γcrCritical cumulative shear strain values for individual slip systems for unit cracking.
And fifthly, simulating the short crack propagation process of the material by adopting a high-efficiency calculation method of short crack propagation, namely a similar periodic jump method. Considering the relatively large number of cycles (about 10) actually involved in the short crack propagation process4Magnitude), in order to save the calculation cost, a similar periodic jump method is adopted to carry out short crack propagation simulation, namely, the cracking of the unit in front of the crack tip is realized only by calculating the accumulated shear strain value of a slip system under a single cycle in a second subprogram, and the method specifically comprises the following steps:
the second subprogram calculates the decomposition shear stress, the slip strength and the slip system accumulated shear strain value of all units in the crystal plastic material region under a single cycle, the calculated single slip system accumulated shear strain under the current increment step and the slip system shear strain value accumulated in the single cycle are transmitted to the first subprogram, the first subprogram judges the maximum accumulated shear strain of the single slip system under the current increment step, and the first subprogram calculates the integration points of the enrichment units, and the cycle number required by the integration points to reach the cracking critical value can be obtained by taking the Gaussian integration points of the C3D8 unit as objects:
wherein N isiThe number of cycles, γ, required for the ith integration point of the cell to reach the cracking thresholdcrCritical cumulative shear strain value of single slip system for unit cracking,the maximum cumulative shear strain value of each slip system is obtained at the current increment step, and the Delta gamma is the cumulative shear strain value of the slip system in a single cycle.
The number of cycles N required to bring the unit integration point to the cracking thresholdiStored in an external file (. dat), in conjunction with a third subroutine for processing the external file, the number of cycles required for the enrichment unit to crack before the crack tip is calculated:
wherein N is the number of cycles required for the cracking of the enrichment unit before the crack tip, NiThe number of cycles required for the unit ith integration point to reach the cracking threshold.
Then, the number N of cycles needed for the enrichment unit to crack before the crack tip is transmitted to a second sub-program, and the slip strength and slip system accumulated shear strain related physical quantity related to the accumulated quantity in the second sub-program are updated:
wherein the content of the first and second substances,the slip strength and the accumulated shear strain of the slip system calculated for the current cycle j +1,the calculated slip strength and slip system cumulative shear strain for the previous cycle j, N is the number of cycles required for unit cracking before the crack tip, Δ gj、ΔγjThe cumulative slip strength, cumulative slip shear strain in the previous cycle.
By the method, the accumulated shear strain value of the slip system under a single cycle is calculated, so that the enrichment unit before the crack tip cracks. The crack will propagate forward along the propagation direction calculated by the first subroutine in the fourth step;
then, as the crack expands forwards to generate a new crack tip, the model repeats the fourth step and the fifth step corresponding to a new crack tip enrichment unit to perform the next crack expansion judgment;
and after the crack propagation simulation is finished, obtaining a short crack propagation path in a visualization module, and further calculating the short crack propagation rate by adopting a secant method.
Further, the first sub-routine is a UDMGINI sub-routine, the second sub-routine is a UMAT sub-routine, and the third sub-routine is a UEXTERNALDB sub-routine.
Compared with the prior art, the invention has the advantages that:
the method accurately characterizes the microstructure of the high-temperature alloy material, introduces the microstructure into a short crack propagation model, focuses on the short crack propagation characteristics, determines the short crack propagation criterion from the physical mechanism level based on the crystal plasticity theory, improves the accuracy of the high-temperature alloy material short crack propagation numerical simulation method, and makes up for the related research deficiencies.
Drawings
FIG. 1 is a flow chart of the method for simulating the short crack propagation of the superalloy material according to the present invention;
FIG. 2 is a schematic view of a short crack propagation model of an exemplary superalloy material;
FIG. 3(a) is a simulation result of a short crack propagation path of an exemplary superalloy material;
FIG. 3(b) is a graph showing the simulation result of the short crack growth rate of the superalloy material of the example.
Detailed Description
In order to make the objects, technical solutions and advantages of the present invention more apparent, the present invention is described in further detail below with reference to the accompanying drawings and embodiments. It should be understood that the specific embodiments described herein are merely illustrative of the invention and are not intended to limit the invention. In addition, the technical features involved in the respective embodiments of the present invention described below may be combined with each other as long as they do not conflict with each other.
The technical scheme of the numerical simulation method for the short crack propagation of the high-temperature alloy material is further explained below with reference to the attached drawings. The material of investigation of this example was powdered superalloy material FGH 96.
As shown in fig. 1, the specific implementation process of the present invention is as follows:
firstly, establishing a FGH96 material short crack propagation finite element model based on ABAQUS software. Based on a notch test piece in an FGH96 material in-situ fatigue test, in combination with the Saint-Venn principle, the model only keeps the examination section part with the minimum middle width of the test piece, and in order to improve the calculation efficiency, only the surface crack condition is considered, the thickness direction is compressed, and the size of the short crack propagation model is 8mm multiplied by 3mm multiplied by 0.005 mm. In order to ensure that cracks can be initiated in the fixed area and can be expanded in a single crack mode, the U-shaped notch with the processing depth of 0.2mm enables stress to be concentrated on the edge of the notch. The model was partitioned, and the region with dimensions of 3mm × 1.5mm near the notch was a crystalline plastic material region, and the remaining region was an isotropic material region. The model grid unit type adopts 8-node hexahedron linear units (C3D8) to generate 12607 grid units, wherein the grid unit subdivision is carried out on the crystalline plastic material area near the notch, the unit size is set to be 0.003mm, 6384 units are generated, and local fine grids are established.
And secondly, fitting the parameters of the crystal plastic material by adopting a representative volume unit model of the polycrystalline material. The strengthening effect of the deformation rate is considered in the crystal plastic material, and the rule for controlling the plastic shear strain rate of a single slip system in the model is as follows:
wherein, the first and the second end of the pipe are connected with each other,the plastic shear strain rate of the slip system alpha,to reference the plastic shear strain rate, ταThe breaking shear stress of slip system alpha, gαAs the current strength of the slip system alpha, hαβIn order to harden the modulus,plastic shear strain rate of slip system beta, NtotalFor the total number of slip trains, n is the rate sensitivity parameter for the slip shear rate, sgn (τ)α) As a function of the sign:
the hardening model used for the model was as follows:
hαβ=qh(γ)(α≠β)
wherein h isααDescribed is the effect of self-hardening, in which the slip system α ═ β, hαβDescribed is the effect of latent hardening, in whichThe transfer system alpha is not equal to beta, q represents the relationship between latent hardening and self hardening behavior of the material, h0、hs、τsThe initial hardening modulus, the saturated hardening modulus and the saturated shear stress are respectively, and gamma is the accumulated shear strain of all slip systems.
Introducing a representative volume unit model which accords with the statistical information description of the microstructure of the FGH96 material, and fitting an orthotropic elastic constant and a slip system hardening parameter in a crystal plastic constitutive model by a uniaxial tensile stress-strain curve of the material to ensure that the stress-strain curve predicted by simulation is accurately matched with a uniaxial tensile result and reflect the mechanical property of the polycrystalline material;
the fitted crystalline plastic material properties were assigned to the vicinity of the model crack initiation, as shown in table 1; the remaining regions were assigned isotropic material properties as obtained from a query in the handbook of Chinese superalloy, as shown in Table 2.
TABLE 1 FGH96 Material Crystal Plastic Material parameters
Wherein, C11、C12、C44Respectively, the elastic constant parameters of the orthotropic material FGH96,for reference to the plastic shear strain rate, n is the rate sensitivity parameter for the slip shear strain rate, h0、hs、τsThe initial hardening modulus, the saturated hardening modulus and the saturated shear stress are respectively shown, and q represents the relationship between latent hardening and self-hardening behaviors of the material.
TABLE 2 isotropic Material Properties of FGH96 Material
And thirdly, finishing the microstructure modeling of a crystal plastic material area of 3mm multiplied by 1.5mm near the model gap. And deriving the unit number of the region, the number of 8 nodes of each unit and the global coordinate of each node, and writing an MATLAB script program to represent the microstructure of the notch test piece based on the EBSD detection result of the FGH96 material notch test piece. The EBSD detection result of the test piece can obtain the position coordinates and Euler angle orientation information of a scanning point, wherein the scanning point refers to a corresponding position point of the sample obtained when an incident electron beam continuously moves to scan a specified area of the sample in the detection process;
the MATLAB script program processing steps for characterizing the microstructure of the notched test piece included:
comparing the position coordinates of the scanning points in the EBSD detection result of the FGH96 material with the information of the grid cells of the model crystal plastic material region, and finishing grouping the scanning points according to whether the scanning points belong to the same cell to obtain the scanning points contained in each cell and the corresponding Euler angle orientation thereofφ、
Combining the similar Euler angles in the grouped units, determining the orientation with the maximum number of scanning points in the unit as the unit dominant orientation, and determining the unit dominant orientationφ、Converting the crystal orientation into a crystal orientation g in an hkl crystal coordinate system, wherein the hkl crystal coordinate system is a face-centered cubic crystal triaxial coordinate system, h, k and l are three crystal axes of the crystal, and the crystal orientation conversion expression is as follows:
as the size of the U-shaped notch of the FGH96 material test piece subjected to wire cutting has dispersity, the geometric shape of the U-shaped notch is not completely consistent with that of the model, when scanning points with similar Euler angles are combined, a non-oriented 'empty' unit is inevitably generated, the position of the 'empty' unit in the crystal plastic material area is judged, and the crystal orientation which is closest to the geometric position of the 'empty' unit is endowed to the 'empty' unit;
comparing the crystal orientations of different units, combining the units with the orientation difference smaller than 10 degrees in adjacent units to form the same crystal grain, obtaining the crystal grain morphology and the corresponding crystal orientation, and storing the related information in a python language form recognizable by ABAQUS commercial finite element software, thereby finishing the microstructure characterization of the FGH96 material;
and (3) importing the grain morphology and the crystal orientation information obtained by running an MATLAB script program into an INP file of the FGH96 material short crack propagation notch model to complete the microstructure modeling of the model crystal plastic material region. The model was created as shown in figure 2. The upper part of the graph shows the whole of the short crack propagation gap model; the left side of the lower part of the figure shows a crystal plastic material area after crystal grains are introduced near the middle gap in an enlarged manner, and unit areas with different colors and gray levels represent different crystal grains; the lower right side of the figure shows a polar diagram of the crystal orientation of this region of the model.
And fourthly, introducing a crystal plasticity theory into an expansion finite element method, and customizing a crack expansion criterion reflecting a short crack expansion physical mechanism. The criterion for controlling the short crack propagation is self-defined by adopting a first subprogram of ABAQUS commercial finite element software, the crack propagation of the unit is controlled by the accumulated shear strain of a single slip system of the enrichment unit in front of the crack tip, and the propagation direction of the short crack is vertical to the normal direction of a slip plane where the maximum accumulated shear strain of the single slip system is located;
the cumulative shear strain of the single slip system of the unit is calculated by the second subprogram of ABAQUS commercial finite element software and stored in the state variable statev, and the cumulative shear strain of the single slip system is:
wherein t is time, γαIs the accumulated shear strain of the slip system alpha,the plastic shear strain rate of the slip system α.
The first subprogram of the self-defining crack propagation criterion only acts on the pre-crack-tip enrichment unit, and the accumulated shear strain of the single slip system calculated and stored by the second subprogram can be called in the first subprogram, so as to judge the cracking time and the crack propagation direction of the crack in the unit. The conditions for crack propagation of the enrichment unit are as follows:
wherein, γαIs the cumulative shear strain of the slip system alpha,plastic shear strain rate, N, of slip system alphatotalIs the total number of slip trains, t is time, γcrCritical cumulative shear strain values for individual slip systems for unit cracking.
And fifthly, simulating the FGH96 material short crack propagation process by adopting a short crack propagation efficient calculation method, namely a similar periodic jump method. Considering the relatively large number of cycles (about 10) actually involved in the short crack propagation process4Magnitude), in order to save the calculation cost, a similar periodic jump method is adopted to carry out short crack propagation simulation, namely, the cracking of the unit in front of the crack tip is realized only by calculating the accumulated shear strain value of a slip system under a single cycle in a second subprogram, and the method specifically comprises the following steps:
the second subprogram calculates the decomposition shear stress, the slip strength and the slip system accumulated shear strain value of all units in the crystal plastic material region under a single cycle, the calculated single slip system accumulated shear strain under the current increment step and the slip system shear strain value accumulated in the single cycle are transmitted to the first subprogram, the first subprogram judges the maximum accumulated shear strain of the single slip system under the current increment step, and the first subprogram calculates the integration points of the enrichment units, and the cycle number required by the integration points to reach the cracking critical value can be obtained by taking the Gaussian integration points of the C3D8 unit as objects:
wherein, NiThe number of cycles, γ, required for the ith integration point of the cell to reach the cracking thresholdcrCritical cumulative shear strain value of single slip system for unit cracking,the maximum cumulative shear strain value of each slip system is obtained in the current increment step, and the delta gamma is the cumulative shear strain value of the slip system in a single cycle.
The number of cycles N required to bring the unit integration point to the cracking thresholdiStored in an external file (. dat), in conjunction with a third subroutine for processing the external file, the number of cycles required for the enrichment unit to crack before the crack tip is calculated:
wherein N is the number of cycles required for the cracking of the enrichment unit before the crack tip, NiThe number of cycles required for the unit ith integration point to reach the cracking threshold.
Then, the number N of cycles needed for the enrichment unit to crack before the crack tip is transmitted to a second sub-program, and the slip strength and slip system accumulated shear strain related physical quantity related to the accumulated quantity in the second sub-program are updated:
wherein, the first and the second end of the pipe are connected with each other,the slip strength and the accumulated shear strain of the slip system calculated for the current cycle j +1,the calculated slip strength and slip system cumulative shear strain for the previous cycle j, N is the number of cycles required for unit cracking before the crack tip, Δ gj、ΔγjThe cumulative slip strength, cumulative slip shear strain in the previous cycle.
By the method, the accumulated shear strain value of the slip system under a single cycle is calculated, so that the enrichment unit before the crack tip cracks. The crack will propagate forward along the propagation direction calculated by the first subroutine in the fourth step;
then, as the crack expands forwards to generate a new crack tip, the model repeats the fourth step and the fifth step corresponding to a new crack tip enrichment unit to perform the next crack expansion judgment;
setting boundary and load conditions according to the actual conditions of the FGH96 material short crack propagation in-situ test, wherein the boundary conditions of the model are left-side hinge joint, right-side cyclic loading and cyclic load maximum stress sigmamax800MPa, stress ratio R is 0.1, and loading frequency is 1 Hz;
and after the crack propagation simulation is finished, obtaining a short crack propagation path in a visualization module, and further calculating the short crack propagation rate by adopting a secant method. The simulation results are shown in fig. 3(a) and 3 (b).
Further, the first sub-routine is a UDMGINI sub-routine, the second sub-routine is a UMAT sub-routine, and the third sub-routine is a UEXTERNALDB sub-routine.
The above examples are provided only for the purpose of describing the present invention, and are not intended to limit the scope of the present invention. The scope of the invention is defined by the appended claims. Various equivalent substitutions and modifications can be made without departing from the spirit and principles of the invention, and are intended to be within the scope of the invention.
Claims (5)
1. A high-temperature alloy material short crack propagation numerical simulation method based on crystal plasticity is characterized by comprising the following steps:
step (1): establishing a short crack propagation finite element model of a high-temperature alloy material, and partitioning the model; the high-temperature alloy material short crack propagation finite element model is established based on the size and test conditions of a test piece in a short crack propagation in-situ test of the high-temperature alloy material; the step of partitioning the model refers to the step of dividing the model into an isotropic material area and a crystal plastic material area, wherein the crystal plastic material attribute is given to the vicinity of the crack initiation area, and the isotropic material attribute is given to the rest areas;
step (2): fitting crystal plastic material parameters by adopting a representative volume unit model of the polycrystalline material; the representative volume unit model of the polycrystalline material is a model which accords with the statistical information description of the microstructure of the material and can accurately reflect the macroscopic mechanical properties of the polycrystalline material; the fitted crystal plastic material parameters refer to orthogonal anisotropic elastic constants and slip system hardening parameters in a crystal plastic constitutive model fitted by a uniaxial tensile stress-strain curve of the material, so that the stress-strain curve predicted by simulation is accurately matched with a uniaxial tensile result;
and (3): characterizing the microstructure of the test piece to complete the microstructure modeling of the crystal plastic material area of the model; the representation of the microstructure of the test piece refers to deriving grid unit numbers, node numbers in the unit and global coordinates of all nodes of a crystal plastic material area of the model, endowing the units of the crystal plastic material area with corresponding crystal orientations based on an electron back scattering diffraction detection result of the test piece, combining the units with smaller crystal orientation angle difference, and dividing the units into the same crystal grains to realize the representation of the microstructure of the material; the step of finishing the microstructure modeling of the crystal plastic material area of the model refers to introducing the obtained material grain morphology and crystal orientation information into the model so as to realize the microstructure modeling of the crystal plastic material area;
and (4): introducing a crystal plasticity theory into an expansion finite element method to simulate a short crack expansion behavior, and customizing a crack expansion criterion reflecting a short crack expansion physical mechanism; the introduction of the crystal plasticity theory into the expansion finite element method is to calculate the relevant physical quantities of model stress and slip system shear strain based on the crystal plasticity theory and simulate the short crack expansion behavior of the material by adopting the expansion finite element method without grid repartition; the self-defined crack propagation criterion reflecting the physical mechanism of short crack propagation refers to a self-defined crack propagation criterion, and the criterion and the crack propagation direction for controlling the cracking of a front unit at the tip of the crack are determined from a crystal layer;
and (5): simulating the short crack propagation process of the material by adopting a quasi-periodic jumping method, and acquiring a short crack propagation path and a propagation rate; the periodic jump-like method is that the crack is expanded in a unit by only simulating one cycle to replace the cycle number required by the unit cracking before the crack tip; the step of obtaining the short crack propagation path and the propagation rate refers to directly observing the short crack propagation path in a visualization module and carrying out post-processing on the short crack propagation path to obtain the corresponding crack length, and then calculating the crack propagation rate by adopting a secant method.
2. The method for simulating the short crack propagation numerical value of the high-temperature alloy material according to claim 1, characterized in that:
the crystal plastic material in the step (2) considers the strengthening effect of the deformation rate, and the rule for controlling the plastic shear strain rate of a single slip system in the model is as follows:
wherein, the first and the second end of the pipe are connected with each other,the plastic shear strain rate of the slip system alpha,to reference the plastic shear strain rate, ταIs the resolved shear stress of slip system alpha, gαAs the current strength of the slip system alpha, hαβIn order to harden the modulus,plastic shear strain rate of slip system beta, NtotalFor the total number of slip trains, n is the rate sensitivity parameter for the slip shear rate, sgn (τ)α) As a function of the sign:
the hardening model used for the model was as follows:
hαβ=qh(γ)(α≠β)
wherein h isααDescribed is the effect of self-hardening, in which the slip system α ═ β, hαβDescribed is the influence of latent hardening, where the slip system α ≠ β, q characterizes the relationship between latent hardening and self-hardening behaviour of the material, h0、hs、τsRespectively is initial hardening modulus, saturated hardening modulus and saturated shear stress, and gamma is the accumulated shear strain of all slip systems;
and endowing the fitted attribute of the crystal plastic material to the area near the crack initiation of the model, and endowing the attribute of the isotropic material to the rest area.
3. The method for simulating the short crack propagation numerical value of the high-temperature alloy material according to claim 1, wherein the method comprises the following steps:
the microstructure modeling of the crystal plastic material area of the model in the step (3) is to compile a script program to reproduce the grain morphology and the crystal orientation of the test piece based on the electron back scattering diffraction detection result of the test piece, so as to realize the microstructure modeling of the model, wherein the electron back scattering diffraction detection result of the test piece obtains the position coordinate of a scanning point and the orientation information of an Euler angle;
the script program processing step of the micro-structure modeling of the model comprises the following steps:
comparing the position coordinates of the scanning points in the electron back scattering diffraction detection result with the information of the grid cells of the crystal plastic material area of the model, and grouping the scanning points according to whether the scanning points belong to the same cell to obtain the scanning points contained in each cell and the corresponding Euler angle orientation of the scanning pointsφ、
Combining the similar Euler angles in the grouped units, determining the orientation with the maximum number of scanning points in the unit as the unit dominant orientation, and determining the unit dominant orientationφ、Converting the crystal orientation g into a crystal orientation g in an hkl crystal coordinate system, wherein the hkl crystal coordinate system is a face-centered cubic crystal triaxial coordinate system, h, k and l are three crystal axes of the crystal, and the crystal orientation conversion expression is as follows:
the processed sizes of the test pieces have dispersity, the geometric shapes of the test pieces are not completely consistent with the model, when scanning points with similar Euler angles are combined, non-oriented 'empty' units inevitably appear, the positions of the 'empty' units in the crystalline plastic material area are judged, and the crystal orientation of the unit closest to the geometric positions of the 'empty' units is endowed to the 'empty' units;
comparing the crystal orientations of different units, combining the units with the orientation difference smaller than 10 degrees in the adjacent units to form the same crystal grain, obtaining the appearance and the corresponding crystal orientation of the crystal grain, and storing the related information in a language form which can be identified by finite element software, so as to finish the microstructure characterization of the high-temperature alloy material;
and (4) importing the crystal grain morphology and crystal orientation information obtained by running the script program into an INP file of the short crack propagation model to complete the microstructure modeling of the model crystal plastic material region.
4. The method for simulating the short crack propagation numerical value of the high-temperature alloy material according to claim 1, wherein the method comprises the following steps:
introducing a crystal plasticity theory into an expansion finite element method to simulate the material short crack expansion behavior in the step (4), wherein a criterion for controlling the short crack expansion is customized by a first subprogram of finite element software, the crack expansion of the unit is controlled by the accumulated shear strain of a single slip system of an enrichment unit in front of the crack tip, and the expansion direction of the short crack is vertical to the normal direction of the slip plane where the maximum accumulated shear strain of the single slip system is located;
the cumulative shear strain of the individual slip trains of the unit is calculated by the second subroutine of the finite element software and stored in the state variables, which are:
wherein t is time, γαIs the accumulated shear strain of the slip system alpha,is the plastic shear strain rate of the slip system alpha;
the first subprogram of the self-defined crack propagation criterion only acts on the front enrichment unit of the crack tip, and the accumulated shear strain of the single slip system calculated and stored by the second subprogram can be called in the first subprogram, so that the cracking time and the crack propagation direction of the crack in the unit are judged; the conditions for crack propagation of the enrichment unit are as follows:
5. The method for simulating the short crack propagation numerical value of the high-temperature alloy material according to claim 1, wherein the method comprises the following steps:
in the step (5), a material short crack propagation process is simulated by adopting a quasi-periodic jumping method, and considering that the number of actually involved cycles in the short crack propagation process is large, the cracking of the unit before the crack tip is realized by calculating the accumulated shear strain value of the slip system under a single cycle in the second subprogram, specifically:
the second subprogram calculates the decomposition shear stress, the slip strength and the slip system accumulated shear strain value of all units in the crystal plastic material region under a single cycle, transmits the calculated single slip system accumulated shear strain under the current increment step and the slip system shear strain value accumulated in the single cycle to the first subprogram, judges the maximum accumulated shear strain of the single slip system under the current increment step through the first subprogram, and obtains the cycle number required by reaching the cracking critical value by taking each Gaussian integration point of the C3D8 unit as an object because the first subprogram calculates the integration point of the enrichment unit:
wherein N isiThe number of cycles, γ, required for the unit ith integration point to reach the cracking thresholdcrThe critical cumulative shear strain value of the individual slip system for unit cracking,the maximum cumulative shear strain value of each slip system is obtained in the current increment step, and the delta gamma is the cumulative shear strain value of the slip system in a single cycle;
the number of cycles N required to bring the unit integration point to the cracking thresholdiStored in an external file, and in conjunction with a third subroutine for processing the external file, the number of cycles required for the enrichment unit to crack before the crack tip is calculated:
wherein N is the number of cycles required for the cracking of the enrichment unit before the crack tip, NiThe number of cycles required for the ith integration point of the unit to reach the cracking critical value;
and then, transmitting the cycle number N required by the enrichment unit before the crack tip to a second sub-program, and updating the slip strength related to the accumulation amount and the slip system accumulated shear strain related physical quantity in the second sub-program:
wherein the content of the first and second substances,the slip strength and the accumulated shear strain of the slip system calculated for the current cycle j +1,the calculated slip strength and the accumulated shear strain of the slip system for the previous cycle j, N is the number of cycles required for unit cracking before the crack tip, and deltagj、ΔγjThe cumulative slip strength and the cumulative slip shear strain in the previous cycle;
therefore, the accumulated shear strain value of the slip system under a single cycle is calculated, so that the enrichment unit in front of the tip of the crack cracks, and the crack forwards expands along the expansion direction calculated by the first subprogram in the fourth step;
then, as the crack is expanded forwards to generate a new crack tip, corresponding to a new crack tip enrichment unit, the model repeats the fourth step and the fifth step to perform the next crack expansion judgment;
and after the crack propagation simulation is finished, obtaining a short crack propagation path in a visualization module, and calculating the short crack propagation rate by adopting a secant method.
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202210180332.7A CN114626263B (en) | 2022-02-25 | 2022-02-25 | High-temperature alloy material short crack propagation numerical simulation method based on crystal plasticity |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202210180332.7A CN114626263B (en) | 2022-02-25 | 2022-02-25 | High-temperature alloy material short crack propagation numerical simulation method based on crystal plasticity |
Publications (2)
Publication Number | Publication Date |
---|---|
CN114626263A true CN114626263A (en) | 2022-06-14 |
CN114626263B CN114626263B (en) | 2024-06-11 |
Family
ID=81901036
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN202210180332.7A Active CN114626263B (en) | 2022-02-25 | 2022-02-25 | High-temperature alloy material short crack propagation numerical simulation method based on crystal plasticity |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN114626263B (en) |
Citations (6)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
US20090240646A1 (en) * | 2000-10-26 | 2009-09-24 | Vextec Corporation | Method and Apparatus for Predicting the Failure of a Component |
CN109725123A (en) * | 2019-02-28 | 2019-05-07 | 北京航空航天大学 | It is a kind of consider shot peening strengthening surface layer grain refinement crack propagation life determine method |
CN111721787A (en) * | 2020-06-24 | 2020-09-29 | 四川大学 | Damage life evaluation method for fatigue crack initiation and propagation based on crystal plasticity |
CN112100702A (en) * | 2020-09-09 | 2020-12-18 | 北京航空航天大学 | Additive material small crack propagation numerical simulation method considering microstructure |
CN113094946A (en) * | 2021-03-23 | 2021-07-09 | 武汉大学 | Phase field model localization self-adaptive algorithm for simulating material cracking |
CN114047210A (en) * | 2021-10-28 | 2022-02-15 | 北京理工大学 | Fatigue crack initiation prediction method considering surface integrity |
-
2022
- 2022-02-25 CN CN202210180332.7A patent/CN114626263B/en active Active
Patent Citations (6)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
US20090240646A1 (en) * | 2000-10-26 | 2009-09-24 | Vextec Corporation | Method and Apparatus for Predicting the Failure of a Component |
CN109725123A (en) * | 2019-02-28 | 2019-05-07 | 北京航空航天大学 | It is a kind of consider shot peening strengthening surface layer grain refinement crack propagation life determine method |
CN111721787A (en) * | 2020-06-24 | 2020-09-29 | 四川大学 | Damage life evaluation method for fatigue crack initiation and propagation based on crystal plasticity |
CN112100702A (en) * | 2020-09-09 | 2020-12-18 | 北京航空航天大学 | Additive material small crack propagation numerical simulation method considering microstructure |
CN113094946A (en) * | 2021-03-23 | 2021-07-09 | 武汉大学 | Phase field model localization self-adaptive algorithm for simulating material cracking |
CN114047210A (en) * | 2021-10-28 | 2022-02-15 | 北京理工大学 | Fatigue crack initiation prediction method considering surface integrity |
Also Published As
Publication number | Publication date |
---|---|
CN114626263B (en) | 2024-06-11 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
Park et al. | Monitoring impact events using a system-identification method | |
Guida et al. | SPH–Lagrangian study of bird impact on leading edge wing | |
Schijve | Fatigue crack growth under variable-amplitude loading | |
CN103970965A (en) | Test run method for accelerated life test of gas turbine engine | |
CN112580228B (en) | Fan blade structural design optimization method for mixed structure of turbofan engine | |
Liu et al. | A numerical approach to simulate 3D crack propagation in turbine blades | |
Zhu et al. | Probabilistic fatigue assessment of notched components under size effect using generalized weakest-link model | |
Chen et al. | Numerical prediction based on XFEM for mixed-mode crack growth path and fatigue life under cyclic overload | |
CN115270239A (en) | Bridge reliability prediction method based on dynamic characteristics and intelligent algorithm response surface method | |
Bunnell et al. | Rapid visualization of compressor blade finite element models using surrogate modeling | |
CN115292849A (en) | Mechanical structure residual life prediction method based on phase field method and BP neural network | |
Machniewicz | Fatigue crack growth prediction models for metallic materials Part II: Strip yield model–choices and decisions | |
Heng et al. | Machine Learning-Assisted probabilistic fatigue evaluation of Rib-to-Deck joints in orthotropic steel decks | |
CN109583037A (en) | A kind of parameter control method of blade of aviation engine shot-peening machining deformation | |
Liu et al. | Modeling fatigue crack growth for a through thickness crack: An out-of-plane constraint-based approach considering thickness effect | |
Pandey et al. | A stress triaxiality based modified Liu–Murakami creep damage model for creep crack growth life prediction in different specimens | |
Tianshuai et al. | Detail fatigue rating method based on bimodal weibull distribution for DED Ti-6.5 Al-2Zr-1Mo-1V titanium alloy | |
Rakowitz et al. | Structured and unstructured computations on the DLR-F4 wing-body configuration | |
CN114626263A (en) | High-temperature alloy material short crack propagation numerical simulation method based on crystal plasticity | |
CN110826182B (en) | Aeroelastic design method of aircraft structure based on vertex method and sequential optimization strategy | |
CN116542177A (en) | Water turbine service life assessment method and system based on start-up and shutdown condition analysis and judgment | |
Hallmann et al. | A method for analyzing the influence of process and design parameters on the build time of additively manufactured components | |
Tvergaard | Discrete modelling of ductile crack growth by void growth to coalescence | |
Kovalev et al. | Multi-criteria analysis of aircraft structures fracture | |
Collins | A multi-fidelity framework for physics based rotor blade simulation and optimization |
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 |