The other day I was talking with the grad student I work for and we were discussing flux boundary conditions. It might help to know more about what it is we’ve been doing lately. As I mentioned in my last post, the grad student has been working on an FEM package. What it does is takes, as input, a file describing a kind of problem to solve. Among other things the file has sections that need to describe the boundary conditions for the problem.

For example, if you have a problem where you want to know the temperature distribution throughout a solid you need to impose certain conditions along the boundary. The two main types of BCs are called Dirichlet and Neumann boundary conditions. A Dirichlet BC is where you impose the actual value of the solution at a particular point (e.g. if one part of the material had its temperature fixed at 100 degrees). A Neumann BC (or flux BC) is where you impose the value of a derivative of the solution (again, for temperature, this would be imposing the amount of heat flux coming into or out of part of the boundary).

What we were talking about was a way to write these flux boundary conditions into the file in the most flexible way possible. What if your flux is a function of position along an edge (which it almost certainly would be)? The way the FEM works it would have to take that value to be constant over the edge (that constant being the average value of the flux through that edge) and it would be kind of a pain in the ass to integrate and find all these average values yourself, especially for a very fine mesh with a lot of elements along the edge.

We decided that it would be nice if the user could simply specify a mathematical expression and have the program integrate that on its own (numerically, of course). Fortunately, as part of the utilities he had written in the past, he had a class capable of evaluating any mathematical expression at run time. It wasn’t capable of accepting parameters to these expressions, though. It only understood expressions that would have a constant value. I figured if we wanted to use this class, we could invent a syntax for denoting an argument and perform a string substitution to get the argument in there. This was a decent solution, but I was a bit wary of the performance of this.

Each time we wanted to evaluate this expression with different parameters, we’d have to turn the parameters into strings (from doubles), substitute them into the string (which is an insertion operation), and then parse and evaluate the expression. A few days prior I had been explaining byte code to him since I had being showing off how easy it is to use Python (and also how you don’t have to wait 10 minutes to compile code each time you change something) and noted how it would execute slower than C/C++ code. Anyway, I brought up the byte code solution for our new problem. We could, instead, parse the expression only once and turn it into a list of tokens representing an RPN expression.

Now, I have written (and rewritten several times) code that is capable of doing this. I told him about it as he left for lunch and when he came back he sounded pretty excited about implementing it. Since I had done this so many times, I said I could write the code pretty quickly and decided to write a new Calculator class for this project. Now, the reason I didn’t use the old code I’ve written before is because that’s for a much more general (and unfinished) expression language. I just needed a simpler subset of that. I also wanted to redesign how it worked.

Now after some tinkering we have a Calculator that can even take a matrix argument (though it can’t do matrix math, it just expects the entire string to be a single matrix). This allows us to input spatially varying material properties and other such things into problems.

Now, I’m not sure I can post the code here since, strictly speaking, it belongs to the university (since I got paid for the hours I worked on it). If people here are interested enough I might actually go through the trouble of finding out whether or not I can post it here. Otherwise I won’t bother.