Abaqus XFEM - Stress intensity factor calculation

Greetings all,

I've been working a lot on simulating stationary 3D cracks in Abaqus XFEM. Trying to benchmark the method against well-known stress intensity factor(SIF) cases. But the convergence of the SIF contour integrals has been the biggest problem so far. Even if the average of SIFs agree reasonably good with benchmark cases, the difference between the contour integrals (contour number n) is bigger. More or less, I get oscillating results for the contours included (using up to 5-15 depending on the model). Theoretically the integral should converge as the region is increased. So my main question in this context is, how are the contour integrals constructed/calculated in Abaqus for XFEM?

At the present time I decided that when too many contours are included, and the contours are far away from the crack tip, might not give such a realistic result since it is too far away from the crack tip behavior, and other parts of the model starts to effect it. Usually the SIFs diverge down from the previous values when too many are included (~12-15). And as well, the first 1-2 (3) have to be disregarded due to large variations. A very fine mesh at the crack-tip does usually not give the diverging SIFs for large contours.

But then is the question, how do I use the remaining contour integral values (3-12)? Averaging? Excluding more? At the moment I'm trying to find a "rule", deciding in which interval the contour integral can be accept with certain meshes where the result converges compared to benchmark cases. To be able to perform further analysis using XFEM.

So now to the real questions to help my analysis of the contour integrals:

How are the elements included in the contour integrals around the crack tip, depending on the number of contours chosen?

According to Abaqus manual:

"To evaluate these integrals, Abaqus defines the domain in terms of rings of elements surrounding the crack tip. Different "contours" (domains) are created. The first contour consists of those elements directly connected to crack-tip nodes. The next contour consists of the ring of elements that share nodes with the elements in the first contour as well as the elements in the first contour. Each subsequent contour is defined by adding the next ring of elements that share nodes with the elements in the previous contour. q is chosen to have a magnitude of zero at the nodes on the outside of the contour and to be one (in the crack direction) at all nodes inside the contour except for the midside nodes (if they exist) in the outer ring of elements. These midside nodes are assigned a value between zero and one according to the position of the node on the side of the element."

So how I have interpreted the manual, closest rings of elements are used for each contour integral. Starting with an easy geometry (see attached picture), I've marked how I expect the rings to be constructed around the tip. Every second ring is marked. This is a uniform case/mesh in the z-direction, and just a partition of the model. Also an explanation of a general 3D geometry would be appreciated of possible (but I think that would be much harder...).

The reason I want to know is that I want to keep the radius fixed for the total number of contour integrals constant, then refine the mesh until convergence (where I simultaneously increases the number of contours to match the radius).

Another problem related to my interpretation of the construction of the element rings, is that I dont get the same contour integral values when I change the number of contours included? For example, comparing the first 4 values when 4 and 5 contours are included."

Im thankful for all answers and as well other tips or material regarding xfem and contour integrals.

[img_assist|nid=12037|title=Contour integrals|desc=|link=none|align=left|width=1600|height=1200]