Status of triangle functions

Hi there!

First of all thanks for making this very useful piece of software.

I was looking into the FMM method to speed up my electron optics simulation package [Traceon](https://github.com/leon-vv/Traceon). For this purpose I'm very interested in using the triangle methods provided by FMMLIB3D.

After some tweaking I was able to get it to work. First I had to uncomment the `"tria"` line in the `.mako` file. Furthermore I had to change `dipvec` to `dipvec.T` on line 226 of the `__init__.py` file. I have attached an example script below for you to test out.

My question is why the triangle methods are not supported out of the box? It seems like all the infrastructure is in place to enable them. If this is purely a manpower issue, would you accept my pull requests to enable them?

Thanks again!

```Python
import pyfmmlib
import numpy as np
from scipy.integrate import dblquad

N = 5

triangles = np.random.rand(N, 3, 3)
centroids = np.mean(triangles, axis=1)

normals = np.cross(triangles[:, 2]-triangles[:, 0], triangles[:, 1]-triangles[:, 0])
normals /= np.linalg.norm(normals, axis=0)

charges = np.random.rand(N)

class Mesh:
def __init__(self, triangles, normals, centroids):
self.triangles = triangles
self.normals = normals
self.centroids = centroids

def __len__(self):
return len(self.triangles)

# Return the exact value for triangle at index i
def exact(i):

target = centroids[i]

def to_integrate(a, b, j):
v1, v2, v3 = triangles[j]

location = v1 + (v3-v1)*a + (v2-v1)*b
area = np.linalg.norm(np.cross(v3-v1, v2-v1))

return area/np.linalg.norm(location-target)

potential = 0.

for j in range(N):
potential += charges[j]*dblquad(to_integrate, 0, 1, 0, lambda y: 1-y, args=(j,))[0]

return potential

m = Mesh(triangles, normals, centroids)
pot1 = pyfmmlib.fmm_tria("p", 0, pyfmmlib.LaplaceKernel(), m, slp_density=charges).real

e = exact(0)/(4*np.pi)
print(f'Accuracy: {e/pot1[0] - 1:.1e}')
```