In this section we shall describe our algorithm to compute the integrals in the discretised variational equation in (3.100) efficiently. The method presented here will be used to solve the discretised integral equation for the crack opening displacement in conjunction with certain iterative solvers for matrix equations. Before we start, we have to define some words precisely. Boundary elements refer to small disjoint patches that cover the crack S. Each boundary element has local nodes on its boundary or in its interior. Local shape functions are defined in the interior of each element, and each of these local shape functions is associated with a local node where it takes the value of 1, while at other local nodes it is equal to 0. Some local nodes of different boundary elements may share a global position called a global node. The union of all the local shape functions that takes the value of 1 at a global node is called a basis function, or a global shape function associated with that global node. The combination of such global shape function and the union of the relevant boundary elements are called the global boundary element associated with the global node. We take Gauss points on each global boundary element to compute the outer integral on the right-hand side of (3.100).
We now move on to the description of the algorithm. We assume that
is known on the crack, and is discretised
Now consider a non-leaf cell C of level l. We compute the multipole moments associated with C by adding all the multipole moments of C's children after shifting the origin from the centroids of C's children to that of C via (3.57) and (3.58). We repeat this procedure tracing the tree structure of cells obtained in step 2 upward (decreasing l) until we reach level 2 cells.
We first compute the coefficients of the local expansion of cells of level 2 according to the definition using (3.60) and (3.61). After that, the coefficients of the local expansion of a level l cell C are computed recursively in 2 steps. Namely, we first add together the contributions of the forms in (3.60) and (3.61) from cells of the level l which are well-separated from C. This is followed by the addition of the coefficients of the local expansion of the parent of C with the origin shifted from the centroid of the parent to that of C via (3.64) or (3.65).
In the actual computation the 4th and 5th steps are carried out at the same time in a single do loop with respect to l. It is seen that the computational cost for the above algorithm is of the order of N if one truncates the infinite sums in (3.60), etc. at a fixed number P. Notice, however, that the size of the boundary elements has to be fine enough so that the conditions of validity of various expansions such as the one stated after (3.105) are not violated.