Hi, I am currently a little stuck in attempting to implement the 2D boundary element method simulation for electrostatics, as outlined in Chapter 3 of Banerjee "The boundary element methods in engineering". I am implementing the solution in "scilab" (scilab.org), and have constructed a simple 4 panel interior problem, with unit potentials on each panel, and the boundary as a 1x1 square. The implementation is attached) I am a little unsure if I am constructing the matrix entries for each of the two matrices Fu=Gq, correctly as I am getting odd results, where the panel fluxes (q_i) are nonzero. In the Banerjee formulation, each panel is placed in its own reference frame, with the normal rotated to align with the X axis, and the AB vector pointing along the positive Y axis For the F matrix, I use -0.5 down the main diagonal, as this corresponds to the "workaround" for the singularity when the field point approaches the panel. For the off-diagonal terms, I am using the [(theta/2*pi)]_A^B solution, theta being the angle formed from the dot product of the (field point to panel) normal and the x axis. For the G matrix I use the -b/(pi*k) ( ln (b) - 1) solution given in Banerjee (3.30, if its useful), on the off-diagonal, and 1/(2*pi*k)*(r*sin(theta)*(log(r)-1) + theta*h) solution. For my problem, k=1. I think I am either misinterpreting the equations, or the coordinate transformations somewhere, but I am unsure. any advice that anyone might have would be very useful! I've tried reversing out a few other sample codes, such as the one by Kirkup, but these each seem to use different formulations of the problem (eg kirkup using numerical integration without coordinate transformations, and the implementation is by no means transparent), which is confusing when trying to match up the matrices. I am of course operating from the assumption that an equi potential problem should yield zero flux - which maybe for some reason in BEM i am wrong?