The Local Formulation for the Modified Green’s Function Method

The Modified Global Green's Function Method (MGGFM) is an integral technique that is characterized by good accuracy in the evaluation of boundary fluxes. This method uses only projections of the Green's Function for the solution of the discrete problem and this is the origin of the term 'Modified' of its name. In this paper the local strategy for calculating the projections of Green's function using de Finite Element Method (FEM) are detailed. The numerical examples show some aspects of the method that had not yet been observed and good results for the flux in all nodes of the mesh.


INTRODUCTION
The Modified Green's Function Method (MGFM) is an integral method that is characterized mainly by the following factors: 1 -Contrary to the BEM (Boundary Element Method) and the Trefftz Method, this method does not use an analytical fundamental solution and/or a Green's function; 2 -Only projections of the Green's Function (denoted by GFp from now on) are required in the calculations.These projections can be calculated with appropriate numerical techniques.The Finite Element Method (FEM) is an interesting technique to obtain the discrete values for the GFp because they do not use fundamental solutions and/or Green's Functions.However, it is necessary to use an additional operator, N + , prescribed on the boundary by the user in order to avoid the singularity of the final system of equations; 3 -Precise values for the boundary generalized fluxes; and Latin American Journal of Solids and Structures 12 (2015) 883-904 4-When using the FEM for the evaluation of the GFp, the MGFM results for the potential (displacements) are identical to the FEM results.Therefore, all measures of the errors based on the displacements of the FEM also apply to the MGFM.
Using the FEM as a numerical method to obtain the discrete values for the GFp at least two approaches are possible: the global and the local formulations.
The global approach (MGGFM) is the formulation that has been most used by several authors in the last two decades and its implementation has been extensively studied.Therefore, this will not be shown in this paper.In this global formulation the auxiliary problems used to calculate the GFp are solved after the assembly of all elements of the mesh and one additional differential operator N + is preferably prescribed at the boundary nodes with homogeneous Dirichlet boundary conditions.The time spent to calculate the GFp is higher than the conventional finite element calculations and the good flux convergence is observed only at the boundary nodes.The following works are examples of this approach: membranes, Barcellos and Silva (1987); Mindlin´s plate, Barbieri andBarcellos (1991, 1993); singular potential, Barcellos and Barbieri (1991); h and p convergence, Filippin et al. (1992); semi-thick shell, Barbieri et al. (1993), non-homogenous potentials problems, Barbieri and Barcellos (1993); subregions, Barbieri et al. (1994), laminated plates, Machado et al. (1994Machado et al. ( , 2008Machado et al. ( , 2012););3D elasticity, Meira Jr. et al. (1995) and mathematical formulation, Barbieri et al. (1998,b).
Using the local formulation (MLGFM) the GFp are obtained in each element of the finite element mesh.The additional operator N + is prescribed by the user at the element boundary nodes.The advantages of this formulation are the accurate flux at all nodes of the mesh and the reduction of processing time of the analysis compared with the global formulation.
This paper was organized in a way of providing the reader with just a review of the MGFM basics mathematical concepts and showing the details of the MLGFM implementation.
The examples solved with the local approach (MLGFM) included in the text show accurate results for the flux in all nodes of the mesh.The solution of the one-dimensional Burgers equation (nonlinear problem) using the Hopf-Cole transformation is another application of this local technique to approximate the GFp due to the flux accurate evaluation at all nodes of the mesh.
The examples solved with the global approach (MGGFM) show new results obtained with high-order (p=4) elements for two-dimensional singular potential problem and new observations for the h and p convergences for regular potential problems.

MGFM: BASIC CONCEPTS AND ITS NUMERICAL IMPLEMENTATION
The mathematical aspects shown in this section and all procedures for the numerical implementation of the MGFM have been previously described by Barcellos and Silva (1987), Silva (1988) and Barbieri et al. (1998).Small reviews of the main concepts of this method are included in the text to show some details that will be discussed throughout this paper.
The u(Q) solution for the linear equation system Latin American Journal of Solids and Structures 12 (2015) 883-904 With boundary conditions properly prescribed can be obtained using the following integral equation: where P and Q are  domain points, (P,Q); p and q are  boundary points, (p,q); N and N * are the Neumann operators associated with the differential operator A and with its adjunct operator A * , respectively; b(P) is the generalized force vector and G(,) is a fundamental solution.
Adding and subtracting the quantity in Eq. ( 2) we obtain: where N + is an additional operator prescribed by the user.The choice and specification of this operator was detailed by Barbieri et al. (1998a).
If it is prescribed the boundary condition (N*+N + )G(p,Q)=0, the function G(,) becomes a Green Function and the solution for u(Q) can be rewritten as follows: where F(p)=(N+N + )u(p).
After the conventional approximations of the variables u(Q), b(P) and F(p) with the finite elements or the boundary elements interpolations, the process of residues orthogonalization with the Galerkin's method of Eq.(4) results: being the set of interpolation functions of the domain (base of the finite elements subspace) and [] is the set of interpolation functions of the boundary (base of the boundary elements subspace).The vector components {u D }, {b} and {f} represent nodal values for the generalized displacement, for the generalized body forces and for the reactions on the boundary, respectively.
Eq.( 6) can be rewritten as: where G d () represents the GFp on the subspace generated by the finite element basis, Latin American Journal of Solids and Structures 12 (2015) 883-904 Repeating the same procedure for Eq.( 5), we have: where u(q)= []{u b }, [] is the basis of the boundary elements subspace (trace of []) and {u b } is the displacements vector on the boundary.According to the convenience, the forces on the boundary F(p) can be interpolated using Again, it is convenient to rewrite Eq.( 8) in the following way: where G b () represents the GFp on the subspace generated by the boundary elements basis, . The numerical calculation of the Green's Function projections G d () and G b () is carried out by solving the following associated problems, Barcellos and Silva (1987): Finally, after calculating the nodal values of G d () and G b (), the nodal values for u(Q) and F(p) are obtained solving Eqs.( 7) and ( 9).

THE MLGFM FORMULATION
The formulation shown in sequence is an extension of the concepts used by Barbieri et al. (1994) and Barbieri and Muñoz (1998) for subregions.A one dimensional problem modeled with linear finite elements is taken as an example to better understand the concepts of the local formulation of the method.Thus, when N + =0 the finite element equations for the problems ( 10) and ( 11) can be written as: Latin American Journal of Solids and Structures 12 (2015) 883-904 To obtain discrete values of the GFp for each node of the element it is necessary to prescribe the additional operator N + 0 because the system of Eq ( 12) is singular.Figure 1 shows the scheme used in this work to prescribe the additional operator N + at the element boundary nodes.The scheme illustrated in Figure 2 shows the boundary conditions and the loading used for the discrete calculation of G d () and G b () for each element.The GFp is calculated solving the following system of equations: Barcellos and Silva (1987) have shown that Eq. ( 7) can be rewritten in the following matricial form:

THE FINAL SYSTEM OF EQUATIONS
where [G c ] and [G d ] are matrices of discrete values of the GFp obtained by solving equations ( 10) and ( 11), respectively.
The assembled final system of equation is obtained considering the potencial (or displacement) continuity and the flux balance for the internal nodes of the mesh.The final system of equations is solved after imposing the boundary conditions of the problem.At the end of this step all nodal values for {u} and {f} are known.At this point, we have to remember that fictitious fluxes were added when prescribing the additional operator N + u=k 0 u at all internal nodes of the mesh during the evaluation of the GFp.Therefore, the correct flux is obtained through subtraction of the additional part k 0 u j for any internal node j of the mesh.

RESULTS
If not mentioned, the flux results of the FEM are obtained using the global smoothing technique described by Hinton and Campbell (1974) and using Lagrangian elements.

One-dimensional Wave Equation
Problem: Find the solution u(x) for: with boundary conditions: u(0)=0 and du(1)/dx=0.The analytical solution for this problem is u(x)= sin(x)/cos(1) -x.Homogeneous finite element meshes (N elements with degree p=1 to p=4) were used to obtain the GFp and the additional operator 'N + u=1.u' was prescribed according to the scheme illustrated in Figure 1.The objective of this analysis is to evaluate the flux convergence at all nodes of the mesh and compare these results with the errors obtained with the FEM.
In Figure 3 are shown the results obtained for the h convergence of the flux at x=0 and in Figure 4 are shown the curves for the relative error, E*=E% FEM /E% MLGFM .All curves of these figures can be approximated with functions of the type, E*C(p)N m , where C(p) and m are two constants that depend on the degree of the element.The adjusted values for these two parameters are shown in Table 1.
Latin American Journal of Solids and Structures 12 (2015) 883-904 Contrary to expectations, the flux results obtained with the MLGFM using the linear element are close to the FEM results and with similar rates of convergence.For higher order elements (p2) the results obtained with the MLGFM are always better than those obtained with the FEM and with higher convergence rates.Whereas the results MLGFM have convergence rate with a value close to m-2p, the results of the FEM present convergence rates near to m-2 for linear and quadratic elements and m-4 for cubic and fourth order elements.More detailed analysis on h and p convergence of finite element solutions for the problem can be found in the literature, for example Ihlenburg (1998), Bouillard (1999), Ihlenburg and Babuska (1997) and Bouillard and Ihlenburg (1999).
In Tables 2 and 3 are shown the results for the flux at x=0 obtained with a homogeneous mesh with 5 elements.The first observation to be made is that the flux errors obtained with MLGFM are approximately constant for all nodes of the mesh and the Neumann condition is automatically satisfied when the boundary condition is imposed for the solution of the final system of equations.Moreover, the error in the flux obtained with the FEM varies from node to node and the Neumann boundary condition at x=1 is only achieved when h0.   1.05310 -9 1.93610 -12 0,8 5.40610 -3
To conclude this example, other two results are shown in Table 4 and Figure 5.The results in Table 4 show the values of flux errors for all nodes of the ends of the elements (N=10, p=2) and the results of Figure 5 illustrate the errors in the potential and flux obtained with the MLGFM (N=20, p=2).It is noteworthy the fact that the errors for the flux are lower than the errors for the potential in Figure 5.
x Figure 5: Errors for potential and flux (N=20 and p=2).

One-dimensional Burgers Equation
The one-dimensional Burgers equation can be written as follows: where u(x,t) denotes the variable of interest,  is a positive constant and its value and physical interpretation depends on the problem being studied and f(x,t) is a source term.There are no analytical solutions for the Burgers equation, but for some combinations of boundary conditions and initial conditions it is possible to obtain the 'exact' solution usually written as infinite series.
For >0 values, it is assumed that there is a function (x,t) such that u(x,t)=(x,t)/x and its derivatives are limited for all (x,t).Substituting this transformation of variables in Eq.( 16) and employing the Hopf-Cole transformation (x,t)=-2ln((x,t)), Cole (1951), then u(x,t) is obtained with the following expression and the non-linear equation ( 16) becomes linear for (x,t), i.e: Employing the Crank-Nicolson rule to approximate the time derivative and writing the Eq. ( 18) for the time (n+1)t we have: which is a one-dimensional linear differential equation with constant coefficients and time-varying excitation.
However, although the Hopf-Cole transformation probably makes the solution of the Burgers equation less laborious, one of the inconveniencies of this method is that the calculation of the partial derivative of the auxiliary function (x,t) with respect to x must be precise to obtain the correct value of u(x,t).Obtaining an accurate value for the spatial derivatives has always been the weak point of many numerical methods used for solving this type of problem.
To verify the efficiency of MLGFM in solving this type of problem, the procedure described above was used in the solution of the Burgers equation with initial conditions given by u(x,0)=sin (x), x[0,1] and boundary conditions u(0,t)=u(1,t)=0, t[0,T].Although an analytical solution exists for this problem it is written in the form of a series with infinite terms and according to Kutluay at al. (1999) 'It is known that the exact solutions for <0.01 fail because of the slow convergence of the infinite series'.
The results shown in the sequence were obtained from the MLGFM using homogeneous meshes with 800 quadratic elements (to compare with others results) for the calculation of the GFp and =0.01.The analytical solution for some instants of time is shown in Figure 6 and we note that some curves have a strong gradient around x=1.The values shown in Table 5 illustrate the influence of the t value in the accuracy of the method.In Figure 7 are carried out comparative studies with literature data.As expected, the errors obtained with the MLGFM are better than those obtained by other techniques due the fact that the partial derivatives of (x, t) with respect to x are better evaluated by this method.Table 5: Influence of the t value in the MLGFM solutions (=0.01,p=2 and N=800).
Error % and presents singularity at r=0.
Fourth order finite elements (p=4) are used to discretize the domain according to the pattern shown in Figure 8a.Since this is a singular problem, special finite elements (obtained from the analytical behavior of the solution) were used around the singular point.These special elements are perfectly compatible with the Lagrangian elements and approach u(r,) with behavior r  for 01.The formulation of the special elements was detailed by Barbieri (2005) and its geometry is shown schematically in Figure 8b.The node 1 is always placed at the singular point and the other nodes are numbered counterclockwise and placed on the semicircle (or circle).The radius R of such semicircle is small enough to ensure that the singularity is quite well represented.The results of this work was obtained using a mesh with 12 subdivisions in the  direction and arranged in pg (geometric progression) subdivisions with ratio equal to 0.5 in the radial direction from r=0.5 to r=R=1.52587910-5 (shaded region in the Figure 8a).
(a) Finite Element Mesh Contrary to observed results in other applications, Barcellos and Barbieri (1991) and Maldaner (1993), the oscillations of the flux results around the analytical solution are very small, even close to the singular point.For (r,)   3 the analytical expression for the flux is u(r,0)/n=0.5r -0.5 and the flux obtained with the MGGFM was u(r,0)/n0.500005r-0.499999 , see Figure 9. Observing the results for the flux of Figure 10 we note that the best values around the singular point were obtained at the inter-element nodes (boundary element mesh).Barcellos and Barbieri (1991) also observed this behavior with elements of lower order.Figure 10: Error % for the flux at the nodes in  3 ={(r,): 0 r1;=}.

Two-dimensional Steady State Heat Transfer
Problem: Find the solution for: The analytical solution of this problem is: with R e =10 and R i =1.The radial flux u(r,)/r, calculated at r=R i =1 is the parameter used in this example for comparisons with the numerical results.
The global formulation of the MGGFM is used to solve this problem and, due to the characteristic of the analytical solution, only a circular sector (05/180) is modeled using only one element in the  direction and N homogeneous elements in the radial direction with variable p (p=2,…,p=8).
The curves plotted in Figure 11 show the h convergence for the radial flux u(1,)/r obtained with different elements (p=2, p=3 and p=4).Besides showing that the values obtained with the MGGFM are much better than the values of the FEM, these results show significant reduction of errors when increasing the order of the element (for both methods).The curves of Figure 12 show the relative errors, E*=Error% FEM /Error% MGGFM , obtained with the elements of order 2, 3 and 4. A mathematical expression of type E*=C(p)N m can be used to represent these errors and the fitted curves valid for N>4 are also shown in this figure.Finally, to complete this example, the problem being studied was again solved with homogeneous mesh with N=5 and p varying from 2 to 8. The results in Figure 13 show the p convergence for the radial flux u(1,)/r.The mathematical expression for the relative error, E*, obtained with these data was E*=3.44942e 1.46107p with values varying from 64.08 for p=2 to 410941 for p=8.

Two-dimensional Potential Problem
Find the solution u(x) for the equation: with n being an integer (n=5 in this example) and with homogeneous Dirichlet boundary conditions.This problem was solved by Barbieri et al. (1998) and he was once again solved with MGGFM to show new conclusions based on the results of the flux.Its analytical solution is u(x)=(1-x 2n )(1-y 2n ) and the analytical flux at x=(1,0) is equal to -10.
Due to the problem´s symmetry, only a quarter of the domain is discretized with uniform finite element meshes and their respective boundary meshes are employed for the evaluation of the GFp using the MGGFM, Figure 14.10.From these results we note that while the rates of p convergence for the FEM grow with approximately constant rate, this same behavior does not happen for the results of the MGGFM where m(2)=3.861;m(3)=3.939;m(4)=5.879;m(5)=5.996and m(6)=7.871.
This behavior is best visualized from the results of the p convergence to the flux in (1,0) illustrated in Figure 17.The results of the FEM are always improved with increasing the order of the elements and this same behavior does not apply to the MGGFM results.However, the values for the relative error, E*, indicated in Table 10 show that for the same element, the results of MGGFM are always better than the results of the FEM.

CONCLUSIONS
In this paper the local formulation for the MGFM was shown.With this formulation it is possible to obtain good results for the flux at all nodes of the mesh.The formulation shown was meant for one-dimensional problems.The extension of these concepts to two-dimensional and threedimensional cases is still a topic under study.
In Example 5.1 the solution of the wave equation is displayed using the local formulation of the method and the results showed: -good convergence for the in all nodes of the mesh for p> 1; -Contrary to our expectations, for p=1 the flux results of the MLGFM are better than the results of the FEM, but with small difference; -the MLGFM errors for the flux are approximately constant for all nodes of the mesh, and -E * values indicate that the change of the quadratic element (p=2) for the cubic element (p=3) does not improve its rate of convergence.Example 5.2 was directed to the solution of the Burgers equation.This is a nonlinear equation and the Hopf-Cole transformation and the Crank-Nicolson rule was used to obtain a discrete time linear problem.Local MGFM approach was used to solve this equation and the results show the influence of t on the accuracy of the method and the good convergence of the method for any instant of time.In Example 5.3 the solution of a singular potential problem using the global formulation of MGFM is displayed.The meshes used to solve this example were prepared with elements of fourth order (p=4) and using specific elements (singular) around the singular point.The results showed that: -the quality of approximations of GFp influences the flux results; -the approximation of F(p) with appropriate functions (singular) produce accurate results for the flux; -no flux fluctuation around the singular point was observed.
Problems 5.4 and 5.5 were solved with MGGFM with the aim of showing some new observations on the flux convergence.More research is still needed on this topic to study the small reduction in error of the flux when odd elements are used to obtain the GFp.
In this paper the FEM was used to obtain GFp and for this reason the MGFM results were compared with the FEM results.Comparisons with flux results obtained with other techniques are still needed to extend the conclusions on the quality of the flux results shown in this paper.

Figure 1 :
Figure 1: Boundary Conditions and the additional operator N + for the local formulation.

Figure 2 :
Figure 2: Loads and boundary conditions for the local formulation.

Figure 8 :
Figure 8: Finite Element Mesh for evaluation of the GFp.
FEM results for the Error % in du(x)/dx (N=5).

Table 10 :
Adjusted curves for the error% in u(1,0)/n and for the Relative Error, E*.