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 PDF

Info

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
Application number
CN202210180332.7A
Other languages
Chinese (zh)
Other versions
CN114626263B (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.)
Beihang University
Original Assignee
Beihang University
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 Beihang University filed Critical Beihang University
Priority to CN202210180332.7A priority Critical patent/CN114626263B/en
Publication of CN114626263A publication Critical patent/CN114626263A/en
Application granted granted Critical
Publication of CN114626263B publication Critical patent/CN114626263B/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
    • G06F30/23Design optimisation, verification or simulation using finite element methods [FEM] or finite difference methods [FDM]
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F2111/00Details relating to CAD techniques
    • G06F2111/10Numerical modelling
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F2119/00Details relating to the type or aim of the analysis or the optimisation
    • G06F2119/14Force 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

High-temperature alloy material short crack propagation numerical simulation method based on crystal plasticity
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:
Figure BDA0003520405080000021
Figure BDA0003520405080000022
wherein the content of the first and second substances,
Figure BDA0003520405080000023
the plastic shear strain rate of the slip system alpha,
Figure BDA0003520405080000024
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,
Figure BDA0003520405080000025
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:
Figure BDA0003520405080000026
the hardening model used for the model was as follows:
Figure BDA0003520405080000027
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
Figure BDA0003520405080000031
φ、
Figure BDA0003520405080000032
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
Figure BDA0003520405080000033
φ、
Figure BDA0003520405080000034
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:
Figure BDA0003520405080000035
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:
Figure BDA0003520405080000041
wherein t is time, γαCumulative shear strain of slip system alpha,
Figure BDA0003520405080000042
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:
Figure BDA0003520405080000043
wherein, γαIs the accumulated shear strain of the slip system alpha,
Figure BDA0003520405080000044
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:
Figure BDA0003520405080000045
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,
Figure BDA0003520405080000046
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:
Figure BDA0003520405080000047
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:
Figure BDA0003520405080000051
Figure BDA0003520405080000052
wherein the content of the first and second substances,
Figure BDA0003520405080000053
the slip strength and the accumulated shear strain of the slip system calculated for the current cycle j +1,
Figure BDA0003520405080000054
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:
Figure BDA0003520405080000061
Figure BDA0003520405080000062
wherein, the first and the second end of the pipe are connected with each other,
Figure BDA0003520405080000063
the plastic shear strain rate of the slip system alpha,
Figure BDA0003520405080000064
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,
Figure BDA0003520405080000065
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:
Figure BDA0003520405080000066
the hardening model used for the model was as follows:
Figure BDA0003520405080000067
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
Figure BDA0003520405080000071
Wherein, C11、C12、C44Respectively, the elastic constant parameters of the orthotropic material FGH96,
Figure BDA0003520405080000077
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
Figure BDA0003520405080000072
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
Figure BDA0003520405080000073
φ、
Figure BDA0003520405080000074
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
Figure BDA0003520405080000075
φ、
Figure BDA0003520405080000076
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:
Figure BDA0003520405080000081
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:
Figure BDA0003520405080000082
wherein t is time, γαIs the accumulated shear strain of the slip system alpha,
Figure BDA0003520405080000083
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:
Figure BDA0003520405080000084
wherein, γαIs the cumulative shear strain of the slip system alpha,
Figure BDA0003520405080000085
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:
Figure BDA0003520405080000091
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,
Figure BDA0003520405080000092
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:
Figure BDA0003520405080000093
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:
Figure BDA0003520405080000094
Figure BDA0003520405080000095
wherein, the first and the second end of the pipe are connected with each other,
Figure BDA0003520405080000096
the slip strength and the accumulated shear strain of the slip system calculated for the current cycle j +1,
Figure BDA0003520405080000097
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:
Figure FDA0003520405070000021
Figure FDA0003520405070000022
wherein, the first and the second end of the pipe are connected with each other,
Figure FDA0003520405070000023
the plastic shear strain rate of the slip system alpha,
Figure FDA0003520405070000024
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,
Figure FDA0003520405070000027
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:
Figure FDA0003520405070000025
the hardening model used for the model was as follows:
Figure FDA0003520405070000026
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
Figure FDA0003520405070000034
φ、
Figure FDA0003520405070000035
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
Figure FDA0003520405070000036
φ、
Figure FDA0003520405070000037
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:
Figure FDA0003520405070000031
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:
Figure FDA0003520405070000032
wherein t is time, γαIs the accumulated shear strain of the slip system alpha,
Figure FDA0003520405070000038
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:
Figure FDA0003520405070000033
wherein, γαIs the accumulated shear strain of the slip system alpha,
Figure FDA0003520405070000048
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.
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:
Figure FDA0003520405070000041
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,
Figure FDA0003520405070000042
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:
Figure FDA0003520405070000043
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:
Figure FDA0003520405070000044
Figure FDA0003520405070000045
wherein the content of the first and second substances,
Figure FDA0003520405070000046
the slip strength and the accumulated shear strain of the slip system calculated for the current cycle j +1,
Figure FDA0003520405070000047
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.
CN202210180332.7A 2022-02-25 2022-02-25 High-temperature alloy material short crack propagation numerical simulation method based on crystal plasticity Active CN114626263B (en)

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)

* Cited by examiner, † Cited by third party
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

Patent Citations (6)

* Cited by examiner, † Cited by third party
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