A contact problem simulation method based on the virtual element method

Abstract

The application provides a contact problem simulation method based on a virtual unit method, which comprises: (1) mesh generation, (1.1) data extraction and mesh construction; (1.2) boundary information extraction; (2) material definition: for different materials, material attribute parameters are given, a constitutive model of the material is defined, and material attribute information of a contact body is provided for subsequent iterative solving; (3) contact search, (3.1) global contact detection; (3.2) local contact detection and non-matching mesh processing; (3.3) contact constraint construction and friction parameter setting; (4) boundary condition application, (4.1) displacement constraint setting; (4.2) stress boundary setting; (5) iterative solving, (5.1) stiffness matrix assembly and residual calculation; (5.2) nonlinear solving strategy; (6) post-processing, (6.1) data extraction and result export; (6.2) visualization processing. The application has reasonable concept, high mesh adaptability, low algorithm complexity, and can realize unified solving of multiple models.

Description

Technical Field

[0001] This invention relates to the field of contact problem simulation technology, and more specifically to a contact problem simulation method based on the virtual element method. Background Technology

[0002] Contact problems are prevalent in engineering fields (such as machining and aerospace structural analysis) and natural phenomena (such as interactions between geological materials). Among traditional numerical methods, the finite element method (FEM) can accurately calculate stress and strain distributions, but its requirement for mesh matching limits its efficiency in handling irregular geometries and complex contact interfaces. While the discrete element method (DEM) is suitable for multi-body contacts, it cannot capture internal deformations of continuums; meshless methods simplify contact detection, but their accuracy is insufficient. Key technologies for solving contact problems in existing solutions include:

[0003] Contact detection technology consists of a multi-stage search strategy and a dynamic contact update mechanism. The multi-stage search strategy first performs a global search, using algorithms such as hierarchical bounding volumes and spatial segmentation to quickly filter potential contact pairs. Secondly, it performs a local precise search, determining the contact point location through nearest-neighbor mapping combined with a point-to-surface algorithm. In large deformation problems, an adaptive small-sliding logic is used to reconfigure the connection relationship between the contact detection points and the target elements, balancing accuracy and efficiency.

[0004] Contact discretization methods are divided into two categories: matched meshes and unmatched meshes. For matched meshes, node-to-node discretization schemes directly and accurately transmit forces, but this is only applicable to matched meshes. For unmatched meshes, node-to-face or face-to-face schemes can effectively handle unmatched mesh problems, but their accuracy depends on the choice of master and slave faces, and force transmission is not accurate enough.

[0005] Constraint application strategies typically include the penalty function method and the Lagrange multiplier method. The former is computationally efficient, but the penetration amount depends on the stiffness parameter and requires empirical adjustment; the latter strictly satisfies the no-penetration condition, but introduces additional degrees of freedom, leading to increased computational complexity.

[0006] Inadequacies of existing technology and problems to be solved:

[0007] Mesh dependency: Traditional finite element method requires fine meshing of the contact interface, which makes it difficult to accurately handle mismatched meshes;

[0008] Algorithm complexity: Non-matching mesh contact requires the use of face-to-face or Mortar methods, which have high computational complexity;

[0009] Model scalability: Existing frameworks require complex stiffness matrix derivations to handle nonlinear material model problems, making them difficult to extend easily.

[0010] In conclusion, it is necessary to further innovate existing technologies. Summary of the Invention

[0011] To address the technical problems existing in the background art mentioned above, this invention proposes a contact problem simulation method based on the virtual element method. It has a reasonable concept, strong mesh adaptability, low algorithm complexity, and can realize unified solution of multiple models, making it suitable for promotion and application.

[0012] To address the aforementioned technical problems, this invention provides a contact problem simulation method based on the virtual element method, which specifically includes the following steps:

[0013] (1) Mesh generation

[0014] (1.1) Data extraction and grid construction

[0015] Extract mesh data from a standard VTK format mesh file in the specified path and construct a virtual element method mesh. This mesh supports the generation of regular polygonal elements or custom arbitrary-shaped polygonal elements and calculates the geometric features of each mesh element.

[0016] (1.2) Boundary Information Extraction

[0017] Extract the boundary nodes of the contact body from the virtual element mesh constructed in the previous step, and construct a set of contact body boundary nodes;

[0018] (2) Material definition

[0019] For different materials, given material property parameters, define the constitutive model of the material;

[0020] (3) Contact Search

[0021] (3.1) Global Contact Detection

[0022] Rapidly filter potential contact pairs and optimize search efficiency using hierarchical spatial partitioning techniques from computer graphics;

[0023] (3.2) Local contact detection and mismatched mesh processing

[0024] Calculate the potential contact nodes of a potential contact body pair, traverse all boundary nodes of the potential contact body pair, mark the boundary nodes with negative or zero signed distance as potential contact nodes, construct the contact node pair between the master surface and the subordinate surface; process the non-matching mesh interface, insert the non-matching contact nodes into the original contact to construct point-to-point contact constraints, and calculate the normal distance between the contact nodes.

[0025] (3.3) Contact constraint construction and friction parameter setting

[0026] Based on the coordinates of the contact node pair constructed in the previous step, calculate the relative displacement and tangential vector of the contact node pair, specify the Coulomb friction coefficient, and generate the normal and tangential contact force transmission rules.

[0027] (4) Apply boundary conditions

[0028] (4.1) Displacement constraint settings

[0029] By marking the boundary surfaces of the contact body that are subject to displacement constraints, displacement boundary conditions are applied to the boundary surfaces of the contact body.

[0030] (4.2) Stress boundary setting

[0031] By marking the boundary surfaces of the contact body subjected to stress loads, stress boundary conditions are applied to the boundary surfaces of the contact body, supporting a step-by-step loading strategy to improve the convergence of nonlinear problems.

[0032] (5) Iterative solution

[0033] (5.1) Stiffness matrix assembly and residual calculation

[0034] The stiffness matrix and residual of each element are calculated based on the virtual element method, and then assembled into the global stiffness matrix and residual.

[0035] (5.2) Nonlinear solution strategy

[0036] The nonlinear solution process involves first setting the residual tolerance and the maximum number of iterations. For large deformation problems, the total load is decomposed into multiple incremental steps, and then the decomposed load values ​​are applied step by step to avoid divergence. Finally, the convergence path is optimized through an adaptive step size adjustment strategy.

[0037] (6) Post-processing

[0038] (6.1) Data extraction and result export

[0039] After iteratively solving the previous step (5), the displacement degree of freedom solution at each grid node is obtained. The strain and stress field are then obtained through geometric equations and constitutive equations. The above results are converted into a standardized VTK file data structure that is compatible with mainstream visualization software.

[0040] (6.2) Visualization Processing

[0041] Based on the VTK file output in the previous step, a cloud map is further generated using post-processing tools to extract the time history curves of physical quantities at specific nodes for analysis and reference.

[0042] The contact problem simulation method based on the virtual element method, wherein step (3.1) is based on the axis-aligned bounding box tree algorithm to quickly screen potential contact body pairs.

[0043] The contact problem simulation method based on the virtual element method, wherein: in step (3.2), the potential contact nodes of a potential contact body pair are calculated using the symbolic distance field technique.

[0044] The contact problem simulation method based on the virtual element method, wherein: in step (3.2), the non-matching mesh interface is processed by the nearest projection algorithm.

[0045] The contact problem simulation method based on the virtual element method, wherein the specific process of step (5.1) is as follows:

[0046] First, define the overall energy U of the contact system:

[0047]

[0048] In the above formula, U p For the strain potential energy of all contacting bodies, U c The energy contributed by the contact constraint to all contact bodies; u is the displacement degree of freedom of the contact node, h is the plastic internal variable of the contact node, Ψ is the constitutive model expression of each contact body, f is the body force on each contact body, and t is the surface force on each contact body.

[0049] Then, volume integrals are performed on Ψ and f·u, and area integrals are performed on t·u. Finally, the potential energies of the n contact bodies are summed. Since the energy of each body is composed of the energy of each element in the mesh, the stiffness matrix and residuals are calculated and assembled in the order from element to whole. The calculation process of the stiffness matrix and residuals is as follows:

[0050] First, calculate the potential energy of each element based on the coordinates and displacement of each element's nodes. Then, the unified construction from energy to stiffness matrix is ​​achieved through automatic differentiation technology, that is, calculating the element potential energy with respect to the element displacement u. e The partial derivatives are used to obtain the residual R. e ;Calculate the residual with respect to the element displacement u e The partial derivatives are used to obtain the stiffness matrix K. e :

[0051]

[0052] The stiffness matrices and residuals of all elements are assembled to obtain the overall stiffness matrix and residuals. For elastoplastic problems, the plastic internal variables are updated first, and then the above calculation process of stiffness matrix and residuals is performed.

[0053] The contact problem simulation method based on the virtual element method, wherein step (5.2) is a nonlinear solution process using the Newton-Raphson iterative method.

[0054] By adopting the above technical solution, the present invention has the following beneficial effects:

[0055] The contact problem simulation method based on the virtual element method of this invention is reasonably conceived and has the following advantages:

[0056] (1) Enhanced mesh adaptability: It supports polygonal units of arbitrary shape and accurately transmits forces between contact interfaces under complex mismatched meshes, thus expanding the mesh adaptability boundary of the original technical solution.

[0057] (2) Reduced algorithm complexity: By discretizing polygonal units, adaptive point-to-point contact format discretization is achieved. Only the unit where the contact interface is located is modified, without the need to add a virtual interface or re-divide the mesh.

[0058] (3) Unified solution of multiple models: Under the automatic differential framework, the calculation can be completed by only modifying the energy function and stress expression when adding a new material model, without the need to manually derive the stiffness matrix.

[0059] The present invention also has the following features and advantages:

[0060] ① A node insertion algorithm based on the virtual element method, which achieves an adaptive point-to-point contact format through polygon discretization and nearest point projection;

[0061] ② The SDF-AABB hybrid contact search algorithm is adopted, which combines the efficiency of AABB trees with the accuracy of the symbolic distance field (SDF) to achieve global search and local matching of multi-body contact;

[0062] ③ A unified solution framework for contact bodies using multiple material models is adopted, namely, iterative solution is achieved through automatic differentiation-driven stiffness matrix construction, which supports coupled analysis of material nonlinearity (Neo-Hookean, Mooney-Rivlin, J2 plasticity) with geometric nonlinearity and contact nonlinearity, eliminates manual derivation errors, and supports arbitrary constitutive model extension.

[0063] In defining materials, this invention automatically associates yield criteria and hardening models for elastoplastic materials, generating corresponding potential energy functions and stress update algorithms; for elastoplastic stress updates, it implicitly iterates to solve for plastic internal variables, ensuring the convergence of the material response. Attached Figure Description

[0064] To more clearly illustrate the specific embodiments of the present invention or the technical solutions in the prior art, the drawings used in the description of the specific embodiments or the prior art will be briefly introduced below. Obviously, the drawings described below are some embodiments of the present invention. For those skilled in the art, other drawings can be obtained from these drawings without creative effort.

[0065] Figure 1 This is a schematic diagram of the polygonal mesh and its symbolic distance field involved in the contact problem simulation method based on the virtual element method of the present invention;

[0066] Figure 2 This is a schematic diagram of the contact detection and contact format involved in the contact problem simulation method based on the virtual element method of this invention. Figure 2 In the figure, (a) is a schematic diagram of contact detection based on AABB tree, and (b) is a schematic diagram of adaptive node-to-node contact format;

[0067] Figure 3 This is a schematic diagram illustrating the derivation of the energy-stiffness matrix under the automatic differentiation framework involved in the contact problem simulation method based on the virtual element method of this invention.

[0068] Figure 4 This is a flowchart of the nested Newton-Raphson iterative algorithm involved in the contact problem simulation method based on the virtual element method of this invention;

[0069] Figure 5 This is an example diagram illustrating the stress and contact force during the post-processing stage in the contact problem simulation method based on the virtual element method of this invention. Figure 5 In the figure, (a) is a stress cloud diagram and (b) is an example diagram of contact force distribution. Detailed Implementation

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

[0071] The present invention will be further explained below with reference to specific embodiments.

[0072] This embodiment provides a contact problem simulation method based on the virtual element method, which includes the following steps:

[0073] (1) Mesh generation:

[0074] (1.1) Data extraction and grid construction

[0075] Extract mesh data from a standard VTK format mesh file located at the specified path and construct a virtual element method mesh (i.e., directly read mesh data and element topology from a standard VTK format file, see...). Figure 1 This mesh supports the generation of regular polygonal cells or custom arbitrary-shaped polygonal cells, and calculates the geometric features (area, centroid coordinates, etc.) of each mesh cell, providing basic data for subsequent contact detection and mechanical analysis.

[0076] (1.2) Boundary Information Extraction

[0077] Extract the boundary nodes of the contact body from the virtual element mesh constructed in the previous step, and construct a set of contact body boundary nodes (that is, store the mesh nodes located on the boundary into a set) to provide interface information for subsequent contact search.

[0078] (2) Material definition

[0079] It supports defining various constitutive models such as linear elastic, hyperelastic, and elastoplastic. Among them, for different materials, given material property parameters such as Young's modulus and Poisson's ratio, the constitutive model of the material is defined (that is, the given material property parameters are substituted into the generally accepted constitutive model formula), which provides contact body material property information for subsequent iterative solutions.

[0080] For elastoplastic materials, the yield criterion (such as J2 plasticity) and hardening model are automatically associated to generate the corresponding potential energy function and stress update algorithm; for elastoplastic stress update, the plastic internal variables are solved implicitly to ensure the convergence of the material response.

[0081] (3) Contact Search

[0082] (3.1) Global Contact Detection

[0083] A fast screening of potential contact body pairs based on the axis-aligned bounding box (AABB) tree algorithm, see [link / reference]. Figure 2 Figure (a) in the figure shows a significant reduction in computational cost; hierarchical spatial partitioning techniques from computer graphics (such as spatial hashing or octrees) are used to optimize search efficiency.

[0084] (3.2) Local contact detection and mismatched mesh processing

[0085] The potential contact nodes of a given potential contact body pair are calculated using the Signed Distance Field (SDF) technique. All boundary nodes of the potential contact body pair are traversed, and boundary nodes with negative or zero signed distances are marked as potential contact nodes. This constructs contact node pairs between the controlling and subordinate surfaces. The nearest projection algorithm is used to process mismatched mesh interfaces, inserting the mismatched contact nodes into the existing contacts to construct point-to-point contact constraints. (See...) Figure 2 In Figure (b), calculate the normal distance between the contact nodes; Figure 2 Figure (b) shows the local element and node insertion of two contact bodies, where Ω1 and Ω2 are two potential contact bodies. and These are the potential contact surfaces of two contacting bodies, where i-1, i, and i+1 are the node numbers of Ω1, and k-1, k, and k+1 are the node numbers of Ω2. and For the newly inserted node of Ω1, and For the newly inserted node of Ω2, n 1Let ξ be the outward normal vector of the contact point, and ξ be the local coordinate system.

[0086] (3.3) Contact constraint construction and friction parameter setting

[0087] Based on the coordinates of the contact node pair constructed in the previous step, calculate the relative displacement and tangential vector of the contact node pair, specify the Coulomb friction coefficient, and generate the normal and tangential contact force transmission rules.

[0088] (4) Apply boundary conditions

[0089] (4.1) Displacement constraint settings

[0090] By marking the boundary surfaces of the contact body that are subject to displacement constraints, displacement boundary conditions are applied to the boundary surfaces of the contact body.

[0091] (4.2) Stress boundary setting

[0092] By marking the boundary surfaces of the contact body subjected to stress loads, stress boundary conditions are applied to the boundary surfaces of the contact body; a step-by-step loading strategy is supported to improve the convergence of nonlinear problems.

[0093] (5) Iterative solution (see...) Figure 3 )

[0094] (5.1) Stiffness matrix assembly and residual calculation

[0095] The stiffness matrix and residuals of each element are calculated using the virtual element method, and then assembled into the global stiffness matrix and residuals. The specific process is as follows:

[0096] First, define the overall energy U of the contact system:

[0097]

[0098] Here, U p For the strain potential energy of all contacting bodies, U c The energy contributed by the contact constraints to all contact bodies; u is the displacement degree of freedom of the contact node, h is the plastic internal variable of the contact node, Ψ is the constitutive model expression of each contact body, f is the body force on each contact body, and t is the surface force on each contact body; then, volume integrals are performed on Ψ and f·u, and surface integrals are performed on t·u, and finally the potential energies of the n contact bodies are summed; since the energy of each body is composed of the energy of each element in the mesh, the stiffness matrix and residuals are calculated and assembled in the order from element to whole, and the calculation process of stiffness matrix and residuals is as follows:

[0099] First, calculate the potential energy of each element based on the coordinates and displacement of each element's nodes. Then, the unified construction from energy to stiffness matrix is ​​achieved through automatic differentiation technology, that is, calculating the element potential energy with respect to the element displacement u.e The partial derivatives are used to obtain the residual R. e ;Calculate the residual with respect to the element displacement u e The partial derivatives are used to obtain the stiffness matrix K. e :

[0100]

[0101] The stiffness matrices and residuals of all elements are assembled to obtain the overall stiffness matrix and residuals; for elastoplastic problems, the plastic internal variables are updated first, and then the above calculation process is performed.

[0102] (5.2) Nonlinear solution strategy

[0103] The Newton-Raphson iterative method is used to handle the nonlinear solution process (see...). Figure 4 Set the residual tolerance and maximum number of iterations; for large deformation problems, decompose the total load into multiple incremental steps and apply the decomposed load values ​​step by step to avoid divergence; optimize the convergence path through an adaptive step size adjustment strategy.

[0104] (6) Post-processing

[0105] (6.1) Data extraction and result export

[0106] Following the iterative solution in the previous step, the displacement degrees of freedom solution at each mesh node is obtained. The strain and stress fields are then derived using geometric and constitutive equations. These results are then converted into a standardized VTK file data structure, compatible with mainstream visualization software (see [link to documentation]). Figure 5 );

[0107] (6.2) Visualization Processing

[0108] Based on the VTK file output in the previous step, a cloud map is further generated using a post-processing tool (such as Paraview software), and the time history curves of physical quantities of specific nodes are extracted for analysis and reference.

[0109] This invention has the following features and advantages:

[0110] ① Adaptive processing of non-matching meshes: Discretized using arbitrary polygonal elements (see...) Figure 1 It supports adaptive point-to-point contact format and eliminates the limitation of non-matching meshes; among them, the virtual element method achieves efficient discretization of complex geometry (such as custom polygon meshes) by expanding the shape function space.

[0111] ② Efficient Contact Search Algorithm: Combining Axis-Aligned Bounding Box (AABB) and Signed Distance Function (SDF), a two-stage contact detection is achieved (see...). Figure 2 Global search for potential contact objects: Fast filtering of potential contact object pairs based on the AABB tree of the CGAL library; Precise matching of local contact points: Calculation of contact point distance and penetration through SDF to construct a set of contact node pairs.

[0112] ③ Unified solution framework: Automatic differentiation technology is used to automatically derive the energy to stiffness matrix, supporting constitutive models such as linear elastic, hyperelastic, and elastoplastic (see...). Figure 3 ), and solves for the plastic internal variables and displacement degrees of freedom through nested Newton-Raphson iterations (see Figure 4 ).

[0113] This invention has a reasonable concept, strong grid adaptability, low algorithm complexity, and can achieve unified solution for multiple models, making it suitable for promotion and application.

[0114] Finally, it should be noted that the above embodiments are only used to illustrate the technical solutions of the present invention, and not to limit them; although the present invention has been described in detail with reference to the foregoing embodiments, those skilled in the art should understand that modifications can still be made to the technical solutions described in the foregoing embodiments, or equivalent substitutions can be made to some or all of the technical features; and these modifications or substitutions do not cause the essence of the corresponding technical solutions to deviate from the scope of the technical solutions of the embodiments of the present invention.

Claims

1. A simulation method for contact problems based on the virtual element method, characterized in that, Specifically, the steps include: (1) Mesh generation (1.1) Data extraction and grid construction Extract mesh data from a standard VTK format mesh file in a specified path and construct a virtual element method mesh. This mesh supports the generation of regular polygonal elements or custom arbitrary-shaped polygonal elements and calculates the geometric features of each mesh element. (1.2) Boundary Information Extraction Extract the boundary nodes of the contact body from the virtual element mesh constructed in the previous step, and construct a set of contact body boundary nodes; (2) Material definition For different materials, given material property parameters, define the constitutive model of the material; (3) Contact Search (3.1) Global Contact Detection Rapidly filter potential contact pairs and optimize search efficiency using hierarchical spatial partitioning techniques from computer graphics; (3.2) Local contact detection and mismatched mesh processing Calculate the potential contact nodes of a potential contact body pair, traverse all boundary nodes of the potential contact body pair, mark the boundary nodes with negative or zero signed distance as potential contact nodes, construct the contact node pair between the master surface and the subordinate surface; process the non-matching mesh interface, insert the non-matching contact nodes into the original contact to construct point-to-point contact constraints, and calculate the normal distance between the contact nodes. (3.3) Contact constraint construction and friction parameter setting Based on the coordinates of the contact node pair constructed in the previous step, calculate the relative displacement and tangential vector of the contact node pair, specify the Coulomb friction coefficient, and generate the normal and tangential contact force transmission rules. (4) Applying boundary conditions (4.1) Displacement constraint settings By marking the boundary surfaces of the contact body that are subject to displacement constraints, displacement boundary conditions are applied to the boundary surfaces of the contact body. (4.2) Stress boundary setting By marking the boundary surfaces of the contact body subjected to stress loads, stress boundary conditions are applied to the boundary surfaces of the contact body, supporting a step-by-step loading strategy to improve the convergence of nonlinear problems. (5) Iterative solution (5.1) Stiffness matrix assembly and residual calculation The stiffness matrix and residual of each element are calculated based on the virtual element method, and then assembled into the global stiffness matrix and residual. (5.2) Nonlinear solution strategy The nonlinear solution process involves first setting the residual tolerance and the maximum number of iterations. For large deformation problems, the total load is decomposed into multiple incremental steps, and then the decomposed load values ​​are applied step by step to avoid divergence. Finally, the convergence path is optimized through an adaptive step size adjustment strategy. (6) Post-processing (6.1) Data extraction and result export After iteratively solving the previous step (5), the displacement degree of freedom solution at each grid node is obtained. The strain and stress field are then obtained through geometric equations and constitutive equations. The above results are then converted into a standardized VTK file data structure that is compatible with mainstream visualization software. (6.2) Visualization Processing Based on the VTK file output in the previous step, a cloud map is further generated using post-processing tools to extract the time history curves of physical quantities at specific nodes for analysis and reference.

2. The contact problem simulation method based on the virtual element method as described in claim 1, characterized in that: In step (3.1), potential contact pairs are quickly screened based on the axis-aligned bounding box tree algorithm.

3. The contact problem simulation method based on the virtual element method as described in claim 1, characterized in that: In step (3.2), the potential contact nodes of a potential contact pair are calculated using the symbolic distance field technique.

4. The contact problem simulation method based on the virtual element method as described in claim 1, characterized in that: In step (3.2), the non-matching mesh interface is processed using the nearest projection algorithm.

5. The contact problem simulation method based on the virtual element method as described in claim 1, characterized in that, The specific process of step (5.1) is as follows: First, define the overall energy U of the contact system: In the above formula, U p For the strain potential energy of all contacting bodies, U c The energy contributed by the contact constraint to all contact bodies; u is the displacement degree of freedom of the contact node, h is the plastic internal variable of the contact node, Ψ is the constitutive model expression of each contact body, f is the body force on each contact body, and t is the surface force on each contact body. Then, volume integrals are performed on Ψ and f·u, and area integrals are performed on t·u. Finally, the potential energies of the n contact bodies are summed. Since the energy of each body is composed of the energy of each element in the mesh, the stiffness matrix and residuals are calculated and assembled in the order from element to whole. The calculation process of the stiffness matrix and residuals is as follows: First, calculate the potential energy of each element based on the coordinates and displacement of each element's nodes. Then, the unified construction from energy to stiffness matrix is ​​achieved through automatic differentiation technology, that is, calculating the element potential energy with respect to the element displacement u. e The partial derivatives are used to obtain the residual R. e ;Calculate the residual with respect to the element displacement u e The partial derivatives are used to obtain the stiffness matrix K. e : The stiffness matrices and residuals of all elements are assembled to obtain the overall stiffness matrix and residuals. For elastoplastic problems, the plastic internal variables are updated first, and then the above calculation process of stiffness matrix and residuals is performed.

6. The contact problem simulation method based on the virtual element method as described in claim 1, characterized in that: Step (5.2) is a nonlinear solution process using the Newton-Raphson iterative method.

Images