Organizing a CFD program written in python

Question

I've spent some time writing a spectral element solver aimed at solving fluid flow problems concerning the motion of self-propelled bodies. I've gotten to the point where the code has grown long enough that I'm having trouble organizing things in a way that is simple and understandable.

Presently, I have written it in python (using numpy/scipy etc.) and I'm trying to take advantage of the use of modules and OOP concepts to break things down better. Unfortunately, I'm still a bit new with OOP concepts as my background is mainly in engineering. I have things organized follows:

1. A "preprocess" module that contains two classes: "MasterElement", "Element" and "Mesh". Basically, the Mesh class contains an array of nodes and a list of Element objects which are each linked to a MasterElement object. Objects instantiated from MasterElement contain attributes for interpolation, differentiation and integration within each element.

2. A "gmsh" module that parses files generated by the Gmsh automatic meshing program.

Now is where I get a little confused as far as the best way to proceed. Basically, I need 4 things to happen:

1. An "initialization" step where all necessary arrays are initialized into memory (a lot of some_array = np.zeros(...)).

2. A "processing" step where each element in the mesh is looped over. The Jacobian is calculated and operators (weak Laplacian matrix, mass matrix etc.) are subsequently computed for each element. Boundary integrals associated with the "right-hand-side" of the global system are also computed and arranged into a global "contour integral" vector.

3. A solution step where the global system is assembled and solved (directly for now) inside a Newton-Raphson iteration loop.

4. Post processing where the solution is plotted etc.

Currently, I have steps 1, 2, and 3 mashed together in a long python script. Should I just break parts of the script into functions? This is straight forward, but I would be passing a ton of arguments between functions which might be just as confusing. Or is there a better way that I haven't considered?

Even if you could just refer me to some literature on this, that would be great. Most of the CFD books I have focus on the algorithms and not necessarily how to write good, reusable code.

Thanks!

Answer 1

I suggest you to have steps 1, 2, 3 in different files. Since each step is conceptually different, and can be replaced by external routines/software it is better to have a modular architecture.

Splitting the whole code into functions is useful. Regarding the passing of arguments, Python has different ways of doing this: optional parameters, kwargs. You can take a look at SfePy code that is a FEM code written in Python, it could be useful as example.