Method for predicting fatigue crack propagation path of outer container of marine cryogenic tank

Abstract

This invention relates to the field of fatigue crack prediction technology, specifically a method for predicting the fatigue crack propagation path of external containers in marine cryogenic tanks. The method includes the following steps: extracting acceleration and liquid level to construct a kinetic energy scalar; combining heat flux and cryogenic gradient to generate a cryogenic impact tensor; using transient temperature and principal stress to obtain distortion values ​​and generate dynamic propagation increments; establishing perturbation displacement vectors using strain energy and stiffness inverse matrix; fusing reference displacements to obtain step-crack displacements at each stage and establishing a trajectory coordinate set. In this invention, swaying kinetic energy is established by calculating the inner product of hull acceleration and liquid level; a cryogenic impact tensor is constructed by combining heat flux and cryogenic gradient; propagation increments are generated based on the distortion characteristics of temperature and principal stress fluctuations exceeding limits; perturbation displacements are extracted using energy release rate and stiffness inverse matrix; step-crack displacements are obtained by fusing with reference displacements and then spliced ​​and adjusted to avoid the accumulation of errors from mesh re-division, thus achieving accurate tracking of the propagation trajectory.

Description

Technical Field

[0001] This invention relates to the field of fatigue crack prediction technology, and in particular to a method for predicting the fatigue crack propagation path of external containers for marine cryogenic tanks. Background Technology

[0002] The field of fatigue crack prediction technology mainly involves the study of damage evolution laws of materials and load-bearing structures under alternating loads. Its core aspects include stress-strain field analysis, fracture mechanics parameter calculation, and crack tip driving force assessment. This overall technical field quantifies the expansion behavior of initial micro-defects inside the structure when subjected to repeated tensile compression, bending, or thermal stress cycles through theoretical formulas, numerical simulations, or physical experiments, thereby providing a basic mechanical basis for the structural physical state assessment of large-scale engineering equipment. Among them, the traditional method for predicting the fatigue crack propagation path of marine cryogenic tank containers refers to the technical evaluation scheme for the fatigue damage evolution process of the outer shell structure of cryogenic storage tanks carrying liquefied natural gas or hydrogen under the action of alternating stress caused by drastic temperature gradient changes, alternating internal and external pressure differences, and ship swaying. Such methods usually use finite element programs to construct the shell geometry model and mesh nodes of the outer container. After setting the initial surface crack or penetrating crack size, wave load, fluid sloshing load, and thermal stress are applied as boundary conditions to the model nodes. Then, the J integral rule or displacement extrapolation method is used to solve the type I, type II, and type III stress intensity factors at the crack tip. Subsequently, the crack propagation length increment under a single stress cycle is calculated based on the Paris formula. The propagation trajectory of the crack is depicted step by step by modifying the coordinates of the crack tip node and re-grinding the local mesh topology of the crack front end.

[0003] Existing technologies apply various load boundaries using finite element programs and employ conventional methods of calculating stress intensity factors through integration to assess length increments. However, these methods fail to accurately reflect the transient impact effects of alternating stress and extreme temperature gradients on structures under severe sea conditions. Calculation models based on unidirectional derivation using ideal theoretical formulas cannot precisely capture the nonlinear distortion characteristics of the strain field at the crack tip under multi-source disturbances. Relying on repeated modifications to node coordinates and re-division of the mesh topology leads to error accumulation and consumes computational resources, failing to meet the need for accurate tracking of structural damage evolution trajectories under multidimensional and complex working conditions. Summary of the Invention

[0004] To address the technical problems existing in the prior art, embodiments of the present invention provide a method for predicting the fatigue crack propagation path of the outer container of a marine cryogenic tank.

[0005] To achieve the above objectives, the present invention adopts the following technical solution: a method for predicting the fatigue crack propagation path of the outer container of a marine cryogenic tank, comprising the following steps:

[0006] S1: Obtain the values ​​of the ship's three-axis acceleration and the vertical height of the tank liquid surface. Perform square root calculation on the ship's three-axis acceleration to extract the total acceleration scalar. Perform product calculation on the total acceleration scalar and the vertical height of the tank liquid surface to establish the swaying impact kinetic energy scalar.

[0007] S2: Extract the nodal thermal conductivity and cryogenic gradient values, introduce a temperature variable, perform a dot product on the swaying impact kinetic energy scalar and nodal thermal conductivity to extract the heat flux value, and obtain the cryogenic impact tensor.

[0008] S3: Call the cryogenic impact tensor to extract the transient temperature and principal stress values ​​of the two phases at the crack tip, and perform absolute deviation extraction of temperature time-series fluctuation values ​​and stress fluctuation values ​​respectively. Take the inner product of the temperature gradient and stress fluctuation values ​​to obtain the dimensionless thermodynamic distortion value, and perform over-limit judgment with the preset material yield benchmark value to generate dynamic expansion increment.

[0009] S4: Based on the dynamic expansion increment, obtain the strain energy value of the grid element and divide it by the dynamic expansion increment to extract the energy release rate. Extract the local stiffness inverse matrix and perform tensor mapping on the energy release rate and the local stiffness inverse matrix to establish the perturbation displacement vector.

[0010] S5: Call the perturbation displacement vector to extract the node reference displacement value, perform vector fusion with the node reference displacement and the perturbation displacement vector to extract the step cracking displacement, and perform array splicing adjustment of the step cracking displacement of each period according to the time sequence to obtain the dynamic trajectory coordinate set of crack propagation.

[0011] As a further aspect of the present invention, the swaying impact kinetic energy scalar includes fluid impact dissipation energy, bulkhead hysteresis work, and surface wave surge potential energy; the cryogenic impact tensor includes the cold contraction principal stress component, low-temperature embrittlement normal stress, and phase transformation-induced crack shear stress; the dynamic propagation increment includes the plastic passivation leading edge step size, microscopic cleavage extension, and grain boundary slip distance; the perturbation displacement vector includes the structural weakening offset, fracture work normal offset, and interface instability slip distance; and the crack propagation dynamic trajectory coordinate set includes the crack leader abrupt change point, the main crack trajectory positioning marker, and the secondary micro-fracture plane vertex.

[0012] As a further aspect of the present invention, the specific steps of S1 are as follows:

[0013] S101: Obtain the values ​​of the ship's three-axis acceleration and the vertical height of the liquid surface in the tank, perform exponentiation on the three-dimensional axial components of the ship's three-axis acceleration, perform summation on the exponentiation values, perform root extraction on the summation values, extract the amplitude of the spatial motion vector, strip the directional attributes, extract the independent measure of the motion magnitude, and generate the total acceleration scalar.

[0014] S102: Based on the total acceleration scalar and the vertical height of the liquid surface in the tank, perform an alignment operation on the total acceleration scalar and the vertical height of the liquid surface, extract the time index elements of the same sequence and perform multiplication calculation, perform traversal and accumulation processing on the obtained multiplication item sequence, aggregate the related data element items, and obtain the liquid surface sloshing impact interaction item.

[0015] S103: Perform energy dimension mapping calculation on the liquid surface sloshing impact interaction term, extract impact energy feature parameters from the mapped values, perform modulus constraint processing on the impact energy feature parameters, strip the residual directional attributes, extract independent absolute numerical parameters and perform scalar encapsulation processing to establish a sloshing impact kinetic energy scalar.

[0016] As a further aspect of the present invention, the specific steps of S2 are as follows:

[0017] S201: Obtain the node thermal conductivity and cryogenic gradient values, perform joint solution on the node thermal conductivity and cryogenic gradient values ​​based on the heat transfer equation to obtain the basic heat conduction flux, and perform energy dissipation conversion processing on the swaying impact kinetic energy scalar to extract the equivalent heat source term of the corresponding node. Align and couple the basic heat conduction flux and the equivalent heat source term to extract the spatial direction mapping parameter and aggregate the heat conduction state feature term to obtain the node heat flux distribution parameter.

[0018] S202: Based on the node heat flow distribution parameters and the cryogenic gradient values, perform coordinate mapping transformation on the cryogenic gradient values, extract the partial derivative parameters of the orthogonal coordinate axes, perform Cartesian product operation on the partial derivative parameters and the node heat flow distribution parameters, construct a two-dimensional mesh feature arrangement sequence, perform topological connection configuration on the feature node sequence, and generate a cryogenic heat flow coupling matrix.

[0019] S203: For the cryogenic thermal flux coupling matrix, extract the internal row list feature elements, perform reshaping processing on the feature elements according to the spatial recombination configuration term to obtain a one-dimensional continuous feature sequence, introduce a time evolution dimension term to the feature sequence for order expansion reconstruction, perform rank transformation constraint operation on the data sequence, and establish a cryogenic shock tensor.

[0020] As a further aspect of the present invention, the specific steps of S3 are as follows:

[0021] S301: Call the cryogenic impact tensor to extract the transient temperature and principal stress values ​​of the two phases at the crack tip. Perform time-series subtraction and absolute value calculation on the transient temperature of the two phases at the crack tip to extract the temperature time-series fluctuation value. Perform subtraction calculation on the principal stress values ​​according to the same source axis to extract the stress fluctuation value. Combine the temperature gradient and stress fluctuation values ​​and perform sequence dimension encapsulation operation to obtain multi-dimensional gradient fluctuation parameters.

[0022] S302: Perform sequence unpacking and splitting on the multidimensional gradient fluctuation parameter to obtain the temperature gradient sequence and stress fluctuation sequence. Perform dimensionless division transformation on the temperature gradient sequence and stress fluctuation sequence respectively. Perform full element traversal and accumulation operation on the obtained product term sequence. Extract the energy distortion characterization parameter and perform scalar mapping processing to generate dimensionless thermodynamic distortion value.

[0023] S303: Extract the preset material yield benchmark value, perform a numerical domain comparison operation between the dimensionless thermal distortion value and the preset material yield benchmark value, determine whether the dimensionless thermal distortion value exceeds the limit of the preset material yield benchmark value, extract the nonlinear evolution factor based on the limit difference term, multiply the nonlinear evolution factor into the benchmark step size term to perform scale amplification calculation, and establish a dynamic expansion increment.

[0024] As a further aspect of the present invention, the process of extracting the nonlinear evolution factor based on the over-limit difference term specifically involves: extracting a preset material hardening index parameter and a reference stress correction term;

[0025] Perform a division operation between the over-limit difference term and the reference stress correction term to obtain the dimensionless difference ratio parameter; then perform a power operation on the dimensionless difference ratio parameter according to the material hardening index parameter.

[0026] The numerical result of the power operation is extracted as a nonlinear evolution factor.

[0027] As a further aspect of the present invention, the specific steps of S4 are as follows:

[0028] S401: Based on the dynamic expansion increment, obtain the strain energy value of the grid cell, perform element alignment operation between the strain energy value of the grid cell and the dynamic expansion increment, perform quotient division operation between the grid strain energy node and the dynamic expansion increment node, extract the energy decay ratio element, perform sequence cascade encapsulation processing, and obtain the energy release rate value.

[0029] S402: Extract the local stiffness inverse matrix, perform dimensional expansion on the local stiffness inverse matrix and the energy release rate value, perform Kronecker product calculation on the expanded matrix and the energy release rate value to construct a multidimensional data space, extract spatial feature parameters, perform coordinate transformation mapping on the feature parameters to obtain the stiffness release mapping tensor;

[0030] S403: Call the stiffness release mapping tensor, perform tensor expansion operation to convert it into a two-dimensional mapping matrix, perform singular value decomposition operation on the two-dimensional mapping matrix, extract non-zero singular value principal component direction elements, perform weighted summation calculation on the principal component direction elements, remove redundant component data, aggregate principal component elements and perform row vector reshaping processing to establish perturbation displacement vector.

[0031] As a further aspect of the present invention, the process of performing a quotient division operation on the mesh strain energy nodes and the dynamically expanded incremental nodes to extract the energy decay ratio element is specifically as follows:

[0032] The values ​​corresponding to the grid strain energy nodes are set as the dividend elements, and the values ​​corresponding to the dynamically expanded incremental nodes are set as the divisors. A division operation is performed on the dividend and divisors, and the quotient data of the division operation is extracted as the energy decay ratio element.

[0033] As a further aspect of the present invention, the specific steps of S5 are as follows:

[0034] S501: Extract the node reference displacement value, perform spatial coordinate alignment operation on the node reference displacement value and the perturbation displacement vector, match the same dimension vector elements, perform element-by-element addition calculation on the coordinate aligned node reference displacement and perturbation displacement vector, aggregate the vector superposition components and perform serialized vector encapsulation to obtain the step cracking displacement vector.

[0035] S502: Based on the step cracking displacement vector, perform timestamp tag parsing, extract real-time periodic displacement parameters and previous periodic evolution parameters, sort the displacement parameters of each period in ascending order according to the time dimension, construct a one-dimensional time-series evolution sequence, perform multi-dimensional matrix expansion and reconstruction assembly on the sorted sequence, and generate a time-series cracking displacement array.

[0036] S503: For the time-series cracking displacement array, extract the internal splicing feature elements, perform smooth interpolation calculation on the feature elements and adjacent nodes, remove overlapping displacement coordinate terms in the array, perform spatial topology mapping processing on the interpolation node data, reconstruct the discrete displacement nodes into continuous spatial geometric points and perform set normalization encapsulation, and establish a dynamic trajectory coordinate set for crack propagation.

[0037] As a further aspect of the present invention, the process of performing multidimensional matrix expansion and reconstruction assembly on the sorted sequence to generate a temporal crack displacement array specifically involves: extracting a preset coordinate axis dimension quantity parameter; performing a division calculation between the total number of elements in the sorted sequence and the coordinate axis dimension quantity parameter to obtain the row segmentation boundary parameter;

[0038] A two-dimensional empty matrix is ​​established based on the row partitioning boundary parameter and the number of coordinate axis dimensions; the data contained in the sorted sequence are filled into the corresponding arrangement node positions of the two-dimensional empty matrix, and the matrix assembly data is extracted as a temporal crack displacement array.

[0039] Compared with the prior art, the advantages and positive effects of the present invention are as follows:

[0040] In this invention, the swaying impact kinetic energy is established by performing an inner product calculation on the ship's acceleration and the liquid level height. The cryogenic impact tensor is constructed by combining the heat flow and cryogenic gradient values. The multidimensional environmental load is transformed into a dynamic boundary reflecting the microscopic forced state. Based on the inner product distortion characteristics of the transient temperature at the crack tip and the principal stress fluctuation, the limit judgment is performed to generate a dynamic propagation increment. The perturbation displacement is extracted by mapping the energy release rate of the mesh element and the local stiffness inverse matrix. The reference displacement and the perturbation displacement are fused to obtain the step cracking displacement and spliced ​​and adjusted according to the time sequence. This avoids the error accumulation and computational consumption caused by mesh re-division, and realizes the accurate reconstruction of the fatigue crack propagation evolution trajectory under various working conditions. Attached Figure Description

[0041] To more clearly illustrate the technical solutions in the embodiments of the present invention, the accompanying drawings used in the description of the embodiments will be briefly introduced below. Obviously, the accompanying drawings described below are only some embodiments of the present invention. For those skilled in the art, other drawings can be obtained based on these drawings without creative effort.

[0042] Figure 1 This is a schematic diagram of the steps of the present invention;

[0043] Figure 2 This is a detailed schematic diagram of S1 of the present invention;

[0044] Figure 3 This is a detailed schematic diagram of S2 of the present invention;

[0045] Figure 4 This is a detailed schematic diagram of S3 of the present invention;

[0046] Figure 5 This is a detailed schematic diagram of S4 of the present invention;

[0047] Figure 6 This is a detailed schematic diagram of S5 of the present invention. Detailed Implementation

[0048] The technical solution of the present invention will now be described with reference to the accompanying drawings.

[0049] To make the technical problems, technical solutions and advantages of the present invention clearer, a detailed description will be given below in conjunction with the accompanying drawings and specific embodiments.

[0050] Please see Figure 1 This invention provides a method for predicting the fatigue crack propagation path of the outer container of a marine cryogenic tank, comprising the following steps:

[0051] S1: Obtain the values ​​of the ship's three-axis acceleration and the vertical height of the tank liquid surface. Perform square root calculation on the ship's three-axis acceleration to extract the total acceleration scalar. Perform product calculation on the total acceleration scalar and the vertical height of the tank liquid surface to establish the swaying impact kinetic energy scalar.

[0052] S2: Extract the nodal thermal conductivity and cryogenic gradient values, introduce the temperature variable, perform a dot product on the swaying impact kinetic energy scalar and nodal thermal conductivity to extract the heat flux value, and obtain the cryogenic impact tensor.

[0053] S3: Call the cryogenic impact tensor to extract the transient temperature and principal stress values ​​of the two phases at the crack tip. Perform absolute deviation extraction of temperature time-series fluctuation values ​​and stress fluctuation values ​​respectively. Take the inner product of the temperature gradient and stress fluctuation values ​​to obtain the dimensionless thermodynamic distortion value. Perform over-limit judgment on the dimensionless thermodynamic distortion value and the preset material yield benchmark value to generate dynamic expansion increment.

[0054] S4: Based on the dynamic expansion increment, obtain the strain energy value of the mesh element and divide it by the dynamic expansion increment to extract the energy release rate. Extract the local stiffness inverse matrix, perform tensor mapping on the energy release rate and the local stiffness inverse matrix, and establish the perturbation displacement vector.

[0055] S5: Call the perturbation displacement vector, extract the node reference displacement value, perform vector fusion with the node reference displacement and the perturbation displacement vector to extract the step cracking displacement, and perform array splicing adjustment of the step cracking displacement of each period according to the time sequence to obtain the dynamic trajectory coordinate set of crack propagation.

[0056] The swaying impact kinetic energy scalar includes fluid impact dissipation energy, bulkhead hysteresis work, and surface wave surge potential energy; the cryogenic impact tensor includes the cold contraction principal stress component, low-temperature embrittlement normal stress, and phase transformation-induced crack shear stress; the dynamic propagation increment includes the plastic passivation leading edge step size, microscopic cleavage extension, and grain boundary slip distance; the perturbation displacement vector includes the structural weakening offset, fracture work normal offset, and interface instability slip distance; and the crack propagation dynamic trajectory coordinate set includes the crack leader abrupt change point, the main crack trajectory location marker, and the secondary micro-fracture plane vertex.

[0057] Please see Figure 2 The specific steps of S1 are as follows:

[0058] S101: Obtain the values ​​of the ship's three-axis acceleration and the vertical height of the liquid surface in the tank, perform exponentiation on the three-dimensional axial components of the ship's three-axis acceleration, perform summation on the exponentiation values, perform root extraction on the summation values, extract the amplitude of the spatial motion vector, strip the directional attributes, extract the independent measure of the motion magnitude, and generate the total acceleration scalar.

[0059] Continuous time-series data was acquired from a piezoelectric 3-axis accelerometer deployed within the target load-bearing structure at a sampling frequency of 100 Hz. Simultaneously, the vertical height of the tank liquid level was acquired from an ultrasonic level gauge. Data cleaning was performed on the time-series data, removing extreme values ​​and null values ​​that deviated more than 20% from the nearest mean. A sliding time window mean filtering operation with five time steps was used to smooth the data. The filtered lateral, longitudinal, and vertical acceleration components were extracted, and each component was exponentiated and summed. The sum was then square-rooted to extract the spatial motion vector amplitude, and the direction attribute was removed to generate a total acceleration scalar. For example, if the lateral, longitudinal, and vertical acceleration values ​​were 3, 4, and 12 m / s², respectively, the exponentiation yielded 9, 16, and 144, the sum was 169, and the square root extracted the vector amplitude to be 13 m / s². After removing the direction attribute, a total acceleration scalar of 13 m / s² was generated.

[0060] S102: Based on the total acceleration scalar and the vertical height of the liquid surface in the tank, perform an alignment operation on the total acceleration scalar and the vertical height of the liquid surface, extract the time index elements of the same sequence and perform multiplication calculation, perform traversal and accumulation processing on the obtained multiplication item sequence, aggregate the related data element items, and obtain the liquid surface sloshing impact interaction item.

[0061] Based on a timestamp tagging mechanism, the continuous time-series total acceleration scalar sequence is aligned with the vertical height of the liquid surface numerical sequence. Elements with timestamp differences within a preset tolerance range of 0.01 seconds are identified as synchronization nodes. The aligned time index elements of the same sequence are extracted and multiplied to generate discrete multiplication terms. A full-sequence traversal and accumulation process is performed on the multiplication term sequences of each time node throughout the entire analysis period, aggregating the time-domain correlated data to generate a liquid surface sloshing impact interaction term. For example, three consecutive aligned nodes are extracted: node 1 has a scalar of 13 m / s² and a height of 5 m, with a product of 65; node 2 has a scalar of 10 m / s² and a height of 6 m, with a product of 60; node 3 has a scalar of 15 m / s² and a height of 4 m, with a product of 60. A traversal and accumulation process is performed on this sequence, i.e., summing 65, 60, and 60, ultimately establishing a liquid surface sloshing impact interaction term with a value of 185.

[0062] S103: Perform energy dimension mapping calculation on the liquid surface sloshing impact interaction term, extract impact energy feature parameters from the mapped values, perform modulus constraint processing on the impact energy feature parameters, strip the residual directional attributes, extract independent absolute numerical parameters and perform scalar encapsulation processing to establish a sloshing impact kinetic energy scalar.

[0063] The system retrieves the corresponding medium density parameter and the equivalent bearing area parameter of the container bottom from the preset physical property database. It then multiplies the liquid surface sloshing impact interaction term value with the medium density and equivalent bearing area parameters to obtain a Joule energy value. This value is extracted as an impact energy feature parameter. If it is negative, the absolute value is extracted and subjected to modulus constraint processing to remove residual directional attributes and obtain an independent absolute numerical parameter. A one-dimensional scalar data label is added to the absolute numerical parameter for encapsulation, establishing a sloshing impact kinetic energy scalar. Specifically, taking the previously calculated interaction term 185 as an example, the system retrieves a medium density of 400 kg / m³ and an equivalent bearing area of ​​2 m². Multiplying 185 by 400 and 2 yields 148000. This value is positive; the original value is preserved, and the directional attribute is removed. A scalar encapsulation is then used to establish a sloshing impact kinetic energy scalar with a value of 148000 Joules.

[0064] Please see Figure 3 The specific steps of S2 are as follows:

[0065] S201: Obtain the node thermal conductivity and cryogenic gradient values, perform joint solution of the node thermal conductivity and cryogenic gradient values ​​based on the heat transfer equation to obtain the basic heat conduction flux, and perform energy dissipation conversion processing on the swaying impact kinetic energy scalar to extract the equivalent heat source term of the corresponding node. Align and couple the basic heat conduction flux and the equivalent heat source term, extract the spatial direction mapping parameter and aggregate the heat conduction state feature term to obtain the node heat flux distribution parameter;

[0066] The extracted nodal thermal conductivity is multiplied by the cryogenic gradient value to derive the basic heat conduction flux. For example, when the extracted nodal thermal conductivity is 0.15 W / m Kelvin and the calculated cryogenic gradient value is 200 Kelvin / m, the two are directly multiplied to obtain a basic heat conduction flux of 30 W / m². The execution node of this calculation logic is to calculate the heat flux result based on the thermal conductivity and temperature gradient parameters, thereby triggering subsequent processing operations on the impact energy. Simultaneously, the step retrieves the swaying impact kinetic energy scalar recorded by the liquid tank pressure sensor, multiplies this swaying impact kinetic energy scalar by a preset fluid viscosity dissipation coefficient to obtain an energy dissipation conversion value, and directly uses this conversion value as the equivalent heat source term for the corresponding node. Based on this, the step extracts the three-dimensional coordinate information of the node, and based on the same node spatial coordinate matching relationship, performs a summation operation on the calculated basic heat conduction flux and the corresponding equivalent heat source term to complete the alignment coupling and obtain the coupled heat flux data. Subsequently, the spatial orientation mapping parameters of each node, namely the node normal vector, are extracted. The coupled heat flux data are multiplied by the normal vector values ​​of the corresponding dimensions to obtain the heat flux components of each dimension. The heat flux component values ​​of each dimension, the node spatial coordinate values, and the current sampling timestamp are combined and concatenated into a data matrix to generate the final node heat flux distribution parameters.

[0067] S202: Based on the node heat flux distribution parameters and the cryogenic gradient values, coordinate mapping transformation is performed on the cryogenic gradient values ​​to extract the partial derivative parameters of the orthogonal coordinate axes. The partial derivative parameters and the node heat flux distribution parameters are then subjected to Cartesian product operation to construct a two-dimensional mesh feature arrangement sequence. Topological connection configuration is performed on the feature node sequence to generate a cryogenic heat flux coupling matrix.

[0068] The cryogenic gradient values ​​are converted into standard rectangular coordinate data using the Jacobian matrix, and the partial derivative parameters along the first and second orthogonal coordinate axes are extracted. The set of partial derivative parameters is then exhaustively multiplied by the Cartesian product operation with the set of nodal heat flux distribution parameters to construct a 2D mesh feature arrangement sequence.

[0069] Table 1: Nodal Heat Flow and Cryogenic Gradient Mapping Characteristics

[0070]

[0071] As shown in Table 1, the heat flow distribution parameter 222000 and the first orthogonal partial derivative parameter 1.5 are extracted from grid node 1. These are multiplied by a Cartesian product to obtain the characteristic sequence combination value 333000. This value is filled into a 2D sequence and participates in network topology connectivity, constructing the core elements of the cryogenic heat flow coupling matrix.

[0072] S203: For the cryogenic thermal flux coupling matrix, extract the internal row list feature elements, perform reshaping processing on the feature elements according to the spatial recombination configuration term to obtain a one-dimensional continuous feature sequence, introduce the time evolution dimension term to the feature sequence for order expansion reconstruction, perform rank transformation constraint operation on the data sequence, and establish the cryogenic shock tensor.

[0073] Extract the eigenvalues ​​of the main diagonal and adjacent diagonal rows of the cryogenic thermal-fluid coupling matrix. Perform reshaping by connecting the first and last elements sequentially from left to right and top to bottom to obtain a 1D continuous feature sequence. Introduce time-step slice copies according to a preset evolution period for expanded reconstruction, forming a 2D data layer with temporal and spatial dimensions, and add material principal axis deformation stress properties to the elements. Taking the feature value 333000 of node 1 generated in Table 1 as an example, reshape it to the beginning of the 1D sequence, expand it by 5 time-step slices to become a 2D square matrix, and establish the cryogenic shock tensor.

[0074] Please see Figure 4 The specific steps of S3 are as follows:

[0075] S301: Call the cryogenic shock tensor to extract the transient temperature and principal stress values ​​of the two phases at the crack tip. Perform time-series subtraction on the transient temperature of the two phases at the crack tip to extract the temperature time-series fluctuation value. Perform subtraction calculation on the principal stress values ​​according to the same source axis to extract the stress fluctuation value. Combine the temperature gradient and stress fluctuation values ​​and perform sequence dimension encapsulation operation to obtain multi-dimensional gradient fluctuation parameters.

[0076] The transient temperature values ​​and principal stress values ​​of the material at the crack tip are extracted for the first and second consecutive phases. The transient temperature value of the second phase is subtracted from that of the first phase, and the absolute value of the difference is taken to extract the temperature time-series fluctuation value. The principal stress value of the second phase is subtracted from that of the first phase by performing a subtraction along the same axis to extract the stress fluctuation value. These two values ​​are bundled in index order, and a sequence dimension encapsulation operation is performed to generate a multidimensional gradient fluctuation parameter. Based on this, combined with representation learning methods from the field of deep learning, the multidimensional gradient fluctuation parameter is mapped to a continuous low-dimensional latent space, thereby stripping away surface physical noise and effectively extracting the deep nonlinear latent feature components most sensitive to the crack initiation stage of the structure. For example, the transient temperature at the crack tip is extracted as -160 degrees Celsius in the first phase and -145 degrees Celsius in the second phase, and the absolute value of the difference is taken to obtain a temperature gradient of 15 degrees Celsius. Simultaneously, the principal stresses of the first phase (120 MPa) and the second phase (180 MPa) are extracted, and the difference is taken to obtain a stress fluctuation value of 60 MPa. Aggregate 15 and 60, and perform a wrapper operation to generate multidimensional gradient fluctuation parameters.

[0077] S302: Perform sequence unpacking and splitting on the multidimensional gradient fluctuation parameter to obtain the temperature gradient sequence and stress fluctuation sequence. Perform dimensionless division transformation on the temperature gradient sequence and stress fluctuation sequence respectively. Perform full element traversal and accumulation operation on the obtained product term sequence. Extract the energy distortion characterization parameter and perform scalar mapping processing to generate dimensionless thermodynamic distortion values.

[0078] The encapsulated multidimensional gradient fluctuation parameters are unpacked and separated into temperature gradient and stress fluctuation sequences. Temperature gradient and stress fluctuation elements at the same grid location are directly multiplied to generate a local product term sequence. All product term elements within the target region are summed to extract the energy distortion characterization parameter. A quotient mapping process is performed using a reference distortion variable provided by the cryogenic steel standard database to generate a dimensionless thermal distortion value. Taking the unpacked temperature gradient 15 and stress fluctuation 60 as an example, multiplying them at corresponding positions yields a product term of 900. Assuming five identical product terms are extracted within the region, summing them yields an energy distortion characterization parameter of 4500. The reference distortion variable parameter 10 set in the database is retrieved, and 4500 is divided by 10 to calculate the scalarized dimensionless thermal distortion value of 450.

[0079] S303: Extract the preset material yield benchmark value, perform a numerical domain comparison operation between the dimensionless thermal distortion value and the preset material yield benchmark value, determine the state where the dimensionless thermal distortion value exceeds the limit of the preset material yield benchmark value, extract the nonlinear evolution factor based on the limit difference term, multiply the nonlinear evolution factor into the benchmark step size term to perform scale amplification calculation, and establish a dynamic expansion increment.

[0080] The material property configuration parameter file is retrieved, and the preset material yield benchmark value for the current working condition is extracted. The dimensionless thermal distortion value is compared with this yield benchmark value. If the dimensionless thermal distortion value is greater than the benchmark value, the condition is determined to exceed the limit. The difference between the two is extracted as the excess value term. This difference term is multiplied by a sensitivity coefficient based on historical test fitting to extract the nonlinear evolution factor. This factor is multiplied into the default benchmark step size term to perform amplification calculations and establish a dynamic expansion increment. In the example demonstration, the preset material yield benchmark value is retrieved as 355, and the previously calculated actual thermal distortion is 450. Comparison shows that 450 is greater than 355, confirming yielding. The difference between the two yields an excess value term of 95. Multiplying 95 by the preset sensitivity coefficient 0.02 yields a nonlinear evolution factor of 1.9. Multiplying 1.9 by the set benchmark step size term of 10 mm yields a dynamic expansion increment of 19 mm.

[0081] Please see Figure 5 The specific steps of S4 are as follows:

[0082] S401: Based on the dynamic expansion increment, obtain the strain energy value of the grid cell, perform element alignment operation between the strain energy value of the grid cell and the dynamic expansion increment, perform quotient division operation between the grid strain energy node and the dynamic expansion increment node, extract the energy decay ratio element, perform sequence cascade encapsulation processing, and obtain the energy release rate value.

[0083] Access the material stiffness evolution analysis program to obtain the strain energy values ​​of discrete mesh elements at the crack tip. Align the set of strain energy values ​​with the set of increment distributions generated based on the dynamic propagation increment in spatial coordinates. Using the mesh element strain energy values ​​as the dividend and the aligned dynamic propagation increment values ​​as the divisor, perform a quotient division operation to extract the energy decay ratio element within a unit step size. Concatenate this ratio element according to the spatial arrangement of discrete nodes to obtain the energy release rate value. For example, if the extracted strain energy of a mesh is 7600 joules, and the corresponding dynamic propagation increment is 19 millimeters, performing a quotient operation by dividing 7600 by 19 yields a nodal energy decay ratio element of 400 joules per millimeter. If the crack tip contains four nodes with the same decay ratio, cascade them according to the propagation direction to obtain a comprehensive energy release rate value of 400 for the surface region.

[0084] S402: Extract the local stiffness inverse matrix, perform dimensional expansion on the local stiffness inverse matrix and the energy release rate value, perform Kronecker product calculation on the expanded matrix and the energy release rate value to construct a multidimensional data space, extract spatial feature parameters, perform coordinate transformation mapping on the feature parameters to obtain the stiffness release mapping tensor;

[0085] The local stiffness inverse matrix of the material is extracted, and its dimension is expanded by adding empty rows and columns, along with the set of energy release rate values. The elements of the expanded stiffness inverse matrix and their corresponding energy release rate scalars are then multiplied using the Kronecker product. Spatial feature parameters are extracted from the multidimensional data space constructed by the product calculation, and coordinate transformation is performed on them using the Euler angle rotation matrix to translate them to the global coordinate system and establish a stiffness release mapping tensor.

[0086] Table 2: Derivation of the Coupled Stiffness Inverse Matrix and Release Rate

[0087]

[0088] As shown in Table 2, at row 1, column 1 of the local stiffness inverse matrix, the original compliance component of 0.05 is matched with an energy release rate of 400. Performing the Kronecker product multiplication yields a mapping result of 20.0. After Euler angular coordinate transformation and direction attenuation calculation, a global feature term of 18.5 is generated. The mapping tensor is then established by traversing according to this rule.

[0089] S403: Call the stiffness release mapping tensor, perform tensor expansion operation to convert it into a two-dimensional mapping matrix, perform singular value decomposition operation on the two-dimensional mapping matrix, extract non-zero singular value principal component direction elements, perform weighted summation calculation on the principal component direction elements, remove redundant component data, aggregate principal component elements and perform row vector reshaping processing, and establish perturbation displacement vector.

[0090] A two-dimensional mapping matrix is ​​extracted from the stiffness release mapping tensor set. Singular value decomposition (SVD) is performed on the covariance matrix, decomposing it into a product of an orthogonal matrix and a diagonal matrix. The principal component direction element corresponding to the largest non-zero singular value is extracted from the diagonal matrix. This element is then weighted and summed according to the corresponding singular value, while redundant components corresponding to smaller singular values ​​are removed. The corrected principal component elements are then reshaped into row vectors in a row-column inverted format to establish perturbation displacement vectors. For example, the first largest singular value obtained from SVD of the covariance matrix is ​​36.5, and the second largest singular value, 2.1, is discarded as redundancy. The principal component corresponding to 36.5, with a horizontal component of 0.8 mm and a vertical component of 0.6 mm, is selected and weighted with 1 for weighted retention. The extracted principal component elements are then inverted and reshaped to generate a standardized perturbation displacement vector parameter set.

[0091] Please see Figure 6 The specific steps of S5 are as follows:

[0092] S501: Extract the node reference displacement value, perform spatial coordinate alignment operation on the node reference displacement value and the perturbation displacement vector, match the same dimension vector elements, perform element-wise addition calculation on the coordinate aligned node reference displacement and perturbation displacement vector, aggregate the vector superposition components and perform serialized vector encapsulation to obtain the step cracking displacement vector.

[0093] The initial static node reference displacement values ​​of the nodes to be analyzed are extracted and aligned with the previously output dynamic perturbation displacement vectors according to the X, Y, and Z axes of the three-dimensional coordinate system. Element-by-element addition is performed on displacement vector elements of the same source dimension. The superimposed components of the composite vectors in each direction are annotated with bracketed data identifier arrays, and serialized vector encapsulation is performed to obtain the step cracking displacement vector. For example, the lateral reference displacement value of 2 mm for the key node is extracted and added to the lateral perturbation displacement value of 0.8 mm to obtain a superimposed component of 2.8 mm. The vertical reference displacement of 1.5 mm is extracted and added to the vertical perturbation displacement of 0.6 mm to obtain 2.1 mm. The lateral 2.8 mm and vertical 2.1 mm are aggregated and encapsulated with an outer identifier to generate the step cracking displacement vector.

[0094] S502: Based on the step cracking displacement vector, perform timestamp tag parsing, extract real-time periodic displacement parameters and previous periodic evolution parameters, sort the displacement parameters of each period in ascending order according to the time dimension, construct a one-dimensional time-series evolution sequence, perform multi-dimensional matrix expansion and reconstruction assembly on the sorted sequence, and generate a time-series cracking displacement array.

[0095] The time series parser is invoked to deeply analyze the timestamp labels attached to the step crack displacement vectors. The system period is compared to separate the real-time periodic displacement parameters of the current analysis step from the historical previous periodic evolution parameters. The parameters are sorted in ascending order of time flow to construct a 1D time series evolution sequence. The time step index and 3D displacement components of the sequence are expanded into independent matrix dimensions, and a multi-dimensional matrix expansion and reconstruction assembly is performed to generate a time series crack displacement array.

[0096] S503: For time-series cracking displacement array, extract internal splicing feature elements, perform smooth interpolation calculation on feature elements and adjacent nodes, remove overlapping displacement coordinate terms in the array, perform spatial topology mapping processing on interpolation node data, reconstruct discrete displacement nodes into continuous spatial geometric points and perform set normalization encapsulation, and establish a dynamic trajectory coordinate set for crack propagation.

[0097] Internal splicing feature elements located at grid boundaries in the temporal crack displacement array are extracted, and a cubic spline interpolation algorithm is used to perform smooth interpolation calculations between them and adjacent nodes. Coordinate values ​​are compared, and redundant overlapping terms are removed when the difference is less than 0.001 mm. Spatial topology mapping is performed on the data according to preset surface rules to generate continuous point connections. The data is uniformly scaled to a reference window with dimensions of 0 to 1 for normalization and encapsulation, establishing a dynamic trajectory coordinate set for crack propagation.

[0098] The above description is merely a specific embodiment of the present invention, but the scope of protection of the present invention is not limited thereto. Any variations or substitutions that can be easily conceived by those skilled in the art within the technical scope disclosed in the present invention should be included within the scope of protection of the present invention. Therefore, the scope of protection of the present invention should be determined by the scope of protection of the described technical solutions.

Claims

1. A method for predicting the fatigue crack propagation path of the outer container of a marine cryogenic tank, characterized in that, Includes the following steps: S1: Obtain the values ​​of the ship's three-axis acceleration and the vertical height of the tank liquid surface. Perform square root calculation on the ship's three-axis acceleration to extract the total acceleration scalar. Perform product calculation on the total acceleration scalar and the vertical height of the tank liquid surface to establish the swaying impact kinetic energy scalar. S2: Extract the nodal thermal conductivity and cryogenic gradient values, introduce a temperature variable, perform a dot product on the swaying impact kinetic energy scalar and nodal thermal conductivity to extract the heat flux value, and obtain the cryogenic impact tensor. S3: Call the cryogenic impact tensor to extract the transient temperature and principal stress values ​​of the two phases at the crack tip, and perform absolute deviation extraction of temperature time-series fluctuation values ​​and stress fluctuation values ​​respectively. Take the inner product of the temperature gradient and stress fluctuation values ​​to obtain the dimensionless thermodynamic distortion value, and perform over-limit judgment with the preset material yield benchmark value to generate dynamic expansion increment. S4: Based on the dynamic expansion increment, obtain the strain energy value of the grid element and divide it by the dynamic expansion increment to extract the energy release rate. Extract the local stiffness inverse matrix and perform tensor mapping on the energy release rate and the local stiffness inverse matrix to establish the perturbation displacement vector. S5: Call the perturbation displacement vector to extract the node reference displacement value, perform vector fusion with the node reference displacement and the perturbation displacement vector to extract the step cracking displacement, and perform array splicing adjustment of the step cracking displacement of each period according to the time sequence to obtain the dynamic trajectory coordinate set of crack propagation.

2. The fatigue crack growth path prediction method for a marine cryogenic tank box outer container according to claim 1, characterized by, The swaying impact kinetic energy scalar includes fluid impact dissipation energy, bulkhead hysteresis work, and surface wave surge potential energy; the cryogenic impact tensor includes the cold contraction principal stress component, low-temperature embrittlement normal stress, and phase transformation-induced crack shear stress; the dynamic propagation increment includes the plastic passivation leading edge step size, microscopic cleavage extension, and grain boundary slip distance; the perturbation displacement vector includes the structural weakening offset, fracture work normal offset, and interface instability slip distance; and the crack propagation dynamic trajectory coordinate set includes the crack leader abrupt change point, the main crack trajectory positioning marker, and the apex of the secondary micro-fracture plane.

3. The method of claim 1, wherein the method is characterized by: The specific steps of S1 are as follows: S101: Obtain the values ​​of the ship's three-axis acceleration and the vertical height of the liquid surface in the tank, perform exponentiation on the three-dimensional axial components of the ship's three-axis acceleration, perform summation on the exponentiation values, perform root extraction on the summation values, extract the amplitude of the spatial motion vector, strip the directional attributes, extract the independent measure of the motion magnitude, and generate the total acceleration scalar. S102: Based on the total acceleration scalar and the vertical height of the liquid surface in the tank, perform an alignment operation on the total acceleration scalar and the vertical height of the liquid surface, extract the time index elements of the same sequence and perform multiplication calculation, perform traversal and accumulation processing on the obtained multiplication item sequence, aggregate the related data element items, and obtain the liquid surface sloshing impact interaction item. S103: Perform energy dimension mapping calculation on the liquid surface sloshing impact interaction term, extract impact energy feature parameters from the mapped values, perform modulus constraint processing on the impact energy feature parameters, strip the residual directional attributes, extract independent absolute numerical parameters and perform scalar encapsulation processing to establish a sloshing impact kinetic energy scalar.

4. The fatigue crack growth path prediction method for a marine cryogenic tank box outer container according to claim 3, characterized by, The specific steps of S2 are as follows: S201: Obtain the node thermal conductivity and cryogenic gradient values, perform joint solution on the node thermal conductivity and cryogenic gradient values ​​based on the heat transfer equation to obtain the basic heat conduction flux, and perform energy dissipation conversion processing on the swaying impact kinetic energy scalar to extract the equivalent heat source term of the corresponding node. Align and couple the basic heat conduction flux and the equivalent heat source term to extract the spatial direction mapping parameter and aggregate the heat conduction state feature term to obtain the node heat flux distribution parameter. S202: Based on the node heat flow distribution parameters and the cryogenic gradient values, perform coordinate mapping transformation on the cryogenic gradient values, extract the partial derivative parameters of the orthogonal coordinate axes, perform Cartesian product operation on the partial derivative parameters and the node heat flow distribution parameters, construct a two-dimensional mesh feature arrangement sequence, perform topological connection configuration on the feature node sequence, and generate a cryogenic heat flow coupling matrix. S203: For the cryogenic thermal flux coupling matrix, extract the internal row list feature elements, perform reshaping processing on the feature elements according to the spatial recombination configuration term to obtain a one-dimensional continuous feature sequence, introduce a time evolution dimension term to the feature sequence for order expansion reconstruction, perform rank transformation constraint operation on the data sequence, and establish a cryogenic shock tensor.

5. The method of claim 4, wherein the method is characterized by: The specific steps for S3 are as follows: S301: Call the cryogenic impact tensor to extract the transient temperature and principal stress values ​​of the two phases at the crack tip. Perform time-series subtraction and absolute value calculation on the transient temperature of the two phases at the crack tip to extract the temperature time-series fluctuation value. Perform subtraction calculation on the principal stress values ​​according to the same source axis to extract the stress fluctuation value. Combine the temperature gradient and stress fluctuation values ​​and perform sequence dimension encapsulation operation to obtain multi-dimensional gradient fluctuation parameters. S302: Perform sequence unpacking and splitting on the multidimensional gradient fluctuation parameter to obtain the temperature gradient sequence and stress fluctuation sequence. Perform dimensionless division transformation on the temperature gradient sequence and stress fluctuation sequence respectively. Perform full element traversal and accumulation operation on the obtained product term sequence. Extract the energy distortion characterization parameter and perform scalar mapping processing to generate dimensionless thermodynamic distortion value. S303: Extract the preset material yield benchmark value, perform a numerical domain comparison operation between the dimensionless thermal distortion value and the preset material yield benchmark value, determine whether the dimensionless thermal distortion value exceeds the limit of the preset material yield benchmark value, extract the nonlinear evolution factor based on the limit difference term, multiply the nonlinear evolution factor into the benchmark step size term to perform scale amplification calculation, and establish a dynamic expansion increment.

6. The method of claim 5, wherein the method is characterized by: The process of extracting the nonlinear evolution factor based on the over-limit difference term specifically involves: extracting the preset material hardening index parameter and the reference stress correction term; Perform a division operation between the over-limit difference term and the reference stress correction term to obtain the dimensionless difference ratio parameter; Perform a power operation on the dimensionless difference ratio parameter according to the material hardening index parameter; The numerical result of the power operation is extracted as a nonlinear evolution factor.

7. The fatigue crack propagation path prediction method for the outer container of a marine cryogenic tank as described in claim 5, characterized in that, The specific steps of S4 are as follows: S401: Based on the dynamic expansion increment, obtain the strain energy value of the grid cell, perform element alignment operation between the strain energy value of the grid cell and the dynamic expansion increment, perform quotient division operation between the grid strain energy node and the dynamic expansion increment node, extract the energy decay ratio element, perform sequence cascade encapsulation processing, and obtain the energy release rate value. S402: Extract the local stiffness inverse matrix, perform dimensional expansion on the local stiffness inverse matrix and the energy release rate value, perform Kronecker product calculation on the expanded matrix and the energy release rate value to construct a multidimensional data space, extract spatial feature parameters, perform coordinate transformation mapping on the feature parameters to obtain the stiffness release mapping tensor; S403: Call the stiffness release mapping tensor, perform tensor expansion operation to convert it into a two-dimensional mapping matrix, perform singular value decomposition operation on the two-dimensional mapping matrix, extract non-zero singular value principal component direction elements, perform weighted summation calculation on the principal component direction elements, remove redundant component data, aggregate principal component elements and perform row vector reshaping processing to establish perturbation displacement vector.

8. The method of claim 7, wherein the method is characterized by: The process of performing a quotient division operation between the grid strain energy nodes and the dynamically expanded incremental nodes to extract the energy decay ratio element is as follows: The values ​​corresponding to the grid strain energy nodes are set as the dividend elements, and the values ​​corresponding to the dynamically expanded incremental nodes are set as the divisors. A division operation is performed on the dividend and divisors, and the quotient data of the division operation is extracted as the energy decay ratio element.

9. The method of claim 7, wherein the method is characterized by: The specific steps of S5 are as follows: S501: Extract the node reference displacement value, perform spatial coordinate alignment operation on the node reference displacement value and the perturbation displacement vector, match the same dimension vector elements, perform element-by-element addition calculation on the coordinate aligned node reference displacement and perturbation displacement vector, aggregate the vector superposition components and perform serialized vector encapsulation to obtain the step cracking displacement vector. S502: Based on the step cracking displacement vector, perform timestamp tag parsing, extract real-time periodic displacement parameters and previous periodic evolution parameters, sort the displacement parameters of each period in ascending order according to the time dimension, construct a one-dimensional time-series evolution sequence, perform multi-dimensional matrix expansion and reconstruction assembly on the sorted sequence, and generate a time-series cracking displacement array. S503: For the time-series cracking displacement array, extract the internal splicing feature elements, perform smooth interpolation calculation on the feature elements and adjacent nodes, remove overlapping displacement coordinate terms in the array, perform spatial topology mapping processing on the interpolation node data, reconstruct the discrete displacement nodes into continuous spatial geometric points and perform set normalization encapsulation, and establish a dynamic trajectory coordinate set for crack propagation.

10. The method of claim 9, wherein the method is characterized by: The process of performing multidimensional matrix expansion and reconstruction assembly on the sorted sequence to generate a temporal crack displacement array is as follows: extract the preset coordinate axis dimension quantity parameter; perform a division calculation between the total number of elements in the sorted sequence and the coordinate axis dimension quantity parameter to obtain the row segmentation boundary parameter; A two-dimensional empty matrix is ​​established based on the row partitioning boundary parameter and the number of coordinate axis dimensions; the data contained in the sorted sequence are filled into the corresponding arrangement node positions of the two-dimensional empty matrix, and the matrix assembly data is extracted as a temporal crack displacement array.

Images