# Exploratory Study of a Green Function Based Solver for Nonlinear Partial Differential Equations
## 3 Numerical Implementation
### The Numerical Transposition of the Integral Scheme
#### 3.2.4 Boundary Conditions
]\) operator corresponds to the de-centered finite difference operator, defined over \(2q+1\) points and de-centered a number of points equal to \(i-M\). These coefficients can be computed using the well-known Fornberg algorithm for the computation of finite differences [73], recalling that their grid step must be divided by the homogeneous \(\Delta x\) before implementing them into the numerical scheme. The reflection, on the other hand, only affects the boundary point, \(x_{b}\), as the rest of the discrete reflections will remain outside of the domain of integration, i.e., \(\Omega_{M}\). This effect can be understood looking back at Figure 2, where we can identify the colored Gaussian as the explicit expression of the Green function and the Delta in its center can be seen as the Taylor series expansion of the Green function. Hence, while the reflection activates since the beginning of the boundary domain, the delta term only appears at the final value of the boundary point, imposing or codifying the Green function homogeneous condition over it. As a result, for the Dirichlet boundary condition, all the weights will be canceled at \(x_{b}\) and in the Neumann boundary the weights will be doubled. Furthermore, this effect will be also transmitted at the non-homogeneous terms through the Green function, causing their propagation to also satisfy the homogeneous boundary condition. We present two homogeneous examples in Figures 12 and 13, with constant and variable diffusive coefficients, respectively. To recover these solutions, it is enough to follow the standard methodology we have introduced through this work, adding at the initial and last points the contribution of the reflection over the Green function, i.e., using the weights from (65) doubled at the Neumann Boundary point and of a null value for the Dirichlet one. Even in the situation of the time evolving and spatially variable diffusion coefficient, the scheme will include the weights for each of the propagation towards the boundary. While we will require to introduce the correction terms for the weights, specially to take into account for the temporal variation correction, the algorithm will still hold its simple shape: propagating the information from the stencil points towards the future value, regardless of whether we reach the end of the domain or not. Finally, the last contribution of the boundary at the numerical scheme is its boundary term, that can be extracted from Equation (9). For the Dirichlet problem, its expression is rather straightforward as it is the only term that can affect the value of the unknown function at the boundary, \(u(x_{b},t_{0}+\tau)\) as the rest are null at \(x_{b}\). Therefore, it must be equal to \(u(x_{b},t+\tau)-u(x_{b},t)\), which is what we also recover analytically by developing the integral, i.e., deriving the Dirichlet Green function and applying it to the discretization of the unknown function at that point. This is why most of the numerical algorithms impose directly the Dirichlet boundary condition successfully without any extra regard on the affectation of the boundary over the rest of the domain. The Neumann boundary condition is also straightforward: all the terms of the integral disappear, leaving only the following,  
\[u_{b,Neu}(x,t_{0}+\tau)=\int_{0}^{\tau}\mathcal{G}_{2,\text{sOD}}\big{(}x,t_ {0}+\tau|x_{b},t_{0}+\tau^{\prime}\big{)}f_{u}\big{(}t_{0}+\tau^{\prime}\big{)} d\tau^{\prime} \tag{66}\]  
where \(f_{u}(t)=a_{2}du/dx\) is the usual flux or derivative imposed at the boundary for the unknown function. Hence, this term can be seen as the propagation of the flux that is injected to or extracted from the domain, which can be abstractly thought also as a discrete nonlinear term production/destruction term of the problem. Figure 12: Temporal evolution of a constant coefficient diffusive problem (**left**) and its first-order spatial derivative (**right**) for an implementation of two homogeneous boundaries at the boundaries of the domain: Neumann at \(x=-1\) and Dirichlet at \(x=1\).