pyCM: An open-source computational framework for residual stress analysis employing the Contour Method

The contour method for assessing residual stress is a widely accessible method in terms of the equipment required and experimental protocol. However, its application is hampered by the absence of a common computational framework to reconstruct the residual stress field and visualise the final result. In order to help regularise the analysis an open-source reconstruction package is presented: the Python Contour Method (pyCM). This package comprises a graphical user interface framework for the purposes for analysing and resolving stresses from experimental contour method data. The code requires no proprietary software and is expandable and fully transparent, enabling users to evaluate it and extend it to develop best practice. The framework is presented and compared against results stemming from a widely published edge-welded beam dataset produced as part of the EC 5th Framework ENPOWER project to highlight applicability.


Motivation and significance
The contour method has been widely applied to assess and quantify residual stress distributions in predominantly, but not * Correspondence to: B49 Pariser Building, Department of Mechanical, exclusively, weldments since the method's inception by Prime [1,2]. More recently, the technique has been reviewed in depth by Prime and DeWald [3] and a best-practice guide has been developed by Hosseinzadeh et al. [4]. The technique is mature in terms of being widely applied, with recent process development activity centred on assessing uncertainty e.g. with sensitivity analysis [5], as well as protocols for residual stresses with high magnitudes, e.g. necessitating self-restraint features and dealing with plasticity [6][7][8][9] and more complicated geometries e.g. thick-walled pipes and tubular cross-sections [10,11].
Briefly, the method assumes the (ideally) stress-free cutting of a component, and then employing high resolution surface   1. A schematic overview of the processing steps associated with the contour method on a component containing residual stress in the z direction that varies between tension (T) and compression (C). The experimental steps (a) involve stress-free cutting (cut) and then measurement of surface variation (measure). The five analysis and reconstruction steps (b) comprise: registration to identify which coordinate system will be employed (reference vs. floating), alignment and averaging of the mirrored floating surface to the reference for subsequent fitting. Finite element analysis is then used to set up the problem (FEA-Pre) including meshing the object and applying the boundary conditions, including the surface profile. The residual stress is then recovered and post processed (FEA-Post).
profiling metrology to measure the deflection of the cut faces caused by stress relaxation. The load distribution required to generate the measured deflection is then computationally reconstructed ( Fig. 1). In this destructive method, both the stress to be measured as well as the component itself, are permanently affected by the cutting process. For the reconstruction process, a linear elastic finite element analysis (FEA) is employed to apply surface deflections to a numerical fit of the surface profile to return it to being perfectly flat. Here, the stress state required to arrive at the surface deflection is the (i) stresses normal to the surface; and (ii) the shear and in-plane stresses relieved by the cut to arrive at mechanical equilibrium. This approach is based on elastic superposition as applied by Bueckner [12,13], with the technique having been validated many times against other residual stress measurement techniques e.g. X-ray/neutron diffraction [2,14]. Furthermore, on the basis of superposition, a complete stress tensor at a point can be ascertained by using multiple cuts [15] or by a combination of diffraction and contour method [16]. Prime and DeWald [3] have extensively described the mechanics and practical considerations involving the method, as described in Fig. 1.
The experimental aspects of the contour method are very accessible, requiring a wire electro-discharge machining (EDM) facility and a coordinate measurement machine (CMM) or noncontact profilometer, such as confocal, differential intensity or interferometry based instruments. By contrast, access to robust, open-source software for the reconstruction of the residual stress field from acquired experimental data is not currently available. While implementations employing MATLAB 1 and Abaqus 2 [17,18] have been made available, these require dataset-specific changes to the underlying computational approaches to extend their application. As a result, there are often differences between implementations proposed by different practitioners. At the present time, different groups using different analysis software may get different results and it is unclear whether this is due to different practical implementations or different analysis software. In order to make intercomparisons easier and to enable different groups to start to converge in terms of best practice, open-source, fully transparent code is needed so that the whole community can compare their methodologies and develop the most robust and repeatable approaches and integrate them within a publicly available framework.
The present work introduces a new open-source tool for contour method analysis. It encompasses the full work flow of analysis steps from preparing the initial cut surface profile records (registration), to alignment and averaging of these profiles, to numerical fitting, to finite element analysis preprocessing and finally, post-processing and interpretation. This is a tool which uses a single data container to store all aspects of an analysis, capture all qualitative choices regarding the inclusion of raw metrology data and embodies current analysis best practices. Unlike all existing offerings, the full analysis procedure can be conducted without the reliance on commercial software. Available via Ref. [19], the development of the software is motivated by the need to (i) provide a common platform to assess and playback analyses for all practitioners; (ii) provide a teaching tool for new practitioners; and (iii) initiate the process of standardisation of the analysis technique.

Software description
The pyCM analysis approach is presented in Fig. 2. The software has been coded entirely in Python 3. Each step in the analysis is a separate module which can be called independently as a standalone GUI, or coupled together into a comprehensive tabbed workflow. Furthermore, all steps/records are held in a singular MATLAB v7.2 data container ( Fig. 3) accessed by the SciPy library. The choice of the data container is to reduce barriers to adoption for users who are currently employing MATLAB and to encourage comparison between this and existing analysis frameworks.
Here, the focus is only on data analysis and visualisation, not the acquisition of contour method data. For this, the reader is directed to the comprehensive best-practice guide by Hosseinzadeh et al. [4] which covers both practical aspects and qualitative data interpretation steps. The five constituent analysis steps (shown in Fig. 2) are described in turn below.

Registration -point cloud editor
Analysis commences by reading in the outlines and the unregistered point clouds from the two experimentally measured surface profiles in the form of ASCII files. The user is provided with the option of designating a set of data (outline and surface profile) as the reference with the other being the floating. The convex hull of the point cloud can be found to determine an outline if one was not available due to the metrology approach taken. The reference dataset is that to which the floating dataset will be subsequently aligned to, and hence the final analysis will have the reference dataset coordinate system. Several options are provided for discounting points from subsequent steps. These points are not deleted, but are tracked with a boolean array identifying a 'masked' status i.e. when the points are loaded in subsequent steps, if they have been masked, they will be highlighted. Reasons to discount points include: points that have large measurement errors because of their proximity to free edges or an excess of data density. The ability to audit this process is an important aspect because it impacts the repeatability of CM analyses between practitioners.

Alignment and averaging of surface profile data
The objective in alignment and averaging of the reference and floating datasets is to provide a single surface for fitting. The benefit of averaging the surfaces is that this will eliminate antisymmetric cutting artefacts [3]. Once the masked reference and floating surfaces have been identified and saved, then the floating dataset requires mirroring using the appropriate reflection plane (dependent on the metrology arrangement) and alignment with the reference data. This is followed by averaging of the two surfaces to a common grid. Here, both manual alignment and automated iterative closest point (ICP) algorithms have been implemented to operate on the outlines of each dataset. Whilst the contour method assumes that an ideal zero-thickness cut has been taken, discrepancies between the outlines can often prevent accurate alignment via ICP. Once aligned, the z values are linearly interpolated on a common x, y grid between the two surfaces, having the same coordinate system as the reference surface. Values at locations outside the data but inside the outline are extrapolated on the same basis. The user is able to specify the grid spacing, although the software will suggest a spacing based on the mean Euclidean distance between points of the reference surface.

Surface fitting
The aim here is to numerically fit the averaged surface to a function that can smooth the data and further be used in imposing FEA boundary conditions. As averaging can highlight further errors which were not caught in the initial point cloud editing, the averaged surface can be masked in the same manner as that employed above. Afterwards, a bivariate spline of order ranging from 2-6 can be fit to the data. The algorithm will smooth the surface according to points at which discrete polynomials are conjoined (knots), with any combination of spacing of these knots. The FITPACK functionality available in SciPy has been employed (see [20] for a comprehensive overview of the algorithms employed).

Finite element analysis pre-processing
In order to use a FEA to calculate the stresses required for a flat surface to be in equilibrium, pre-processing is required. The geometry of the cut component needs to be meshed, and material properties and displacement boundary conditions need to be imposed prior to a linear elastic analysis. The displacement boundary conditions are those dictated by the fitted surface as well as those needed to prevent rigid body movement. These are required for force and moment equilibrium.
A mesh is first created by Abaqus CAE (proprietary) or Gmsh (open-source) through scripts generated by pyCM. Rigid body boundary conditions, as well as those informed by the surface fitting step are then applied along with singular elastic material properties. This can then be automatically submitted to either the Abaqus (proprietary) or CalculiX (open-source) solver directly from the pyCM interface. Standardised meshing has been adopted to ensure that the mesh conforms to best practices. This includes ensuring that the outline can be respaced such that nodes are regularly placed along the outline while preserving corners. Meshes created with Abaqus are not necessarily required to be processed by the Abaqus solver; meshes comprising fully integrated 8noded quadrilaterals or 4-noded tetrahedra returned from either Abaqus CAE or Gmsh are converted or natively employed (respectively) in a legacy VTK format (unstructured grid) that are then rendered. From this point, boundary conditions and material properties are applied for interpretation by either Abaqus or CalculiX, and the input is submitted with a request to write all stress components at integration points to a data file in an ASCII format.

Finite element analysis post-processing
The original residual stress distribution across the component can be recovered employing the native post-processors of the FE package employed. However, to ensure that identical results are achieved with commercial and open-source solvers, a postprocessor was developed. In this step, the data file generated by the FEA solver is read and stresses at the nodes are calculated using shape functions and Gaussian quadrature [21]. This is required in order for VTK to render a contour plot over the mesh corresponding to the calculated residual stress. A generalised nodal value calculation was adopted for 8-noded, fully integrated elements, with both node and integration numbering identical between Abaqus and CalculiX. Subsequently, the user has the ability to extract stress values at regular intervals based on shape function interpolation, which can be plotted and exported for direct comparison to the results from other residual stress measurement techniques, as demonstrated for the illustrative example below. Furthermore, results are saved in a VTK XML format, which can be loaded and manipulated with ParaView, an open-source FEA post-processor [22].  [27]. A strain gauge was placed on the opposite side of the beam to where the weld was performed for a combined slitting/contour analysis EDM cut.

Illustrative example
In order to illustrate the workflow enabled by pyCM, contour cut metrology data was obtained from Hosseinzadeh et al. [23] for an AISI 316H autogenously edge-welded beam. This surface profile and outline dataset [24] was selected because the residual stress in the artefact has also been measured by a number of methods, facilitating validation. Measurements include synchrotron and neutron diffraction measurements [25]. Furthermore, this particular specimen has been used as a benchmark for experimental and modelling efforts into the development of residual stress due to welding for nuclear applications [9,26]. A depiction of the sample configuration is given in Fig. 4. In the following sections, pyCM is used to process this data and hence to reconstruct the longitudinal residual stress.
For direct comparison purposes, registration, alignment and averaging as well as surface fitting were conducted according to the contour method analysis parameters documented by Hosseinzadeh et al. [23]. This included employing a cubic bivariate spline with a knot spacing of 3 mm across the beam width (transverse) and 7 mm over the height (normal). The fitted surface captures the averaged data with the same efficacy as described by Hosseinzadeh et al., once the spurious data related to the EDM cut completion was masked by manually windowing the data using the tool afforded by the GUI. The basis for omitting data for fitting is that it did not coincide with the overall trend representing the 'obvious outliers' described by Hosseinzadeh et al.; this coincided with approximately the last 4 mm of the cut. Purposefully omitting data can be inadvisable unless the experimenter has good reasons why this should be conducted. In the current analysis, the cut was performed during a slitting operation. The combination of this technique in conjunction with cutting-induced plasticity led to surface displacements measurements that were greater than elsewhere across the slice by over an order of magnitude. Excluding this last 4 mm of cut meant that a knot spacing any smaller than 7 mm could not be used in the longitudinal direction in order for successful fitting. The details relating to this process are available as an example dataset within the pyCM repository [19].
It is found that when viewing results with the pyCM postprocessor, there is no discernible difference in the result (that is the results are identical) irrespective of whether the mesh was made by Gmsh or Abaqus CAE and solved by either CalculiX or Abaqus. This finding is important and necessitated the inclusion of a custom post processor, as each FEA package-specific postprocessor will interpret results slightly differently in terms of colour maps for contour plots, as well as the accuracy of the shape function interpretation. One of the benefits of the pyCM post processor is that it can provide 'line scans' of residual stresses, similar to those obtained by X-ray or neutron diffraction experiments, calculated directly from Gaussian quadrature. This permits direct comparison to results obtained previously by X-ray or neutron diffraction for this component; previously stresses obtained by CM for specific points have been inferred either from integration points or at nodes, which often do not correlate precisely with neutron/X-ray line scans.
The results obtained by pyCM are compared to those published by Hosseinzadeh et al. [23] in Fig. 5 are evidently nearly identical. Differences in peak stress values between the contour method approaches can be attributed to slightly different choices in averaged data, mesh, spline and interpolation of the FEA results.

Impact
The advantages of pyCM over existing publicly available code include portability of results, sharing of analyses and providing practitioners the ability to play out different analysis strategies when analysing a particular dataset more effectively. This functionality was demonstrated recently to validate a materials processing model by Medina-Juárez et al. [28]. Factors which necessitate these decisions are inherent to the technique, and the present contribution lends itself to making the results obtained by the contour method more transparent. This is timely as the contour method has been identified as being a key technique in addressing issues identified by the emergent additive manufacturing field [29,30].
Current advances in the practical execution of the technique currently being introduced to obviate some of the limitations of the contour method remain to be incorporated into the pyCM framework. Predominantly, developments have taken place both on experimental aspects of the technique (cutting strategies), surface fitting approaches and FEA considerations as described in Section 1. The present work affords natural extension of this research, and therefore the current and future user community is large.

Conclusions
An open-source platform for performing contour method analyses to determine residual stress has been developed called pyCM. This tool has been developed to perform initial point cloud processing, surface fitting as well as management of the FEA pre and post processing. The platform has been confirmed to provide the same results whether using a free and open-source finite element package (CalculiX) or a commercial one, Abaqus. An exemplar analysis workflow has been presented examining the stresses in an autogenously welded plate. The results have been compared against those from a published benchmark [23], with the platform providing near-identical results to those produced using closed-source analysis tools.
The presented approach has been developed in such a manner to easily compare the results obtained by different practitioners. The open-source framework which has been employed enables a transparent approach to contour method analysis, allowing convergence between different practitioners for the same incoming dataset, as well as enabling extensions to account for advancements in the contour method.

Declaration of competing interest
The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.