I would like to model laminar flow of water from roots to the stem of a plant. At the very end of the roots, the tubes vary from millimeter to centimeter scale in diameter and length. As we get closer to the stem the roots get larger in length and diameter. I want to create random 3D domains that represents the network of roots with varying diameters and lengths. What would be the best way to create this geometry.
Chances are, you don't want something truly random; you want something that has the same abstract 3D structure as a plant root system, but beyond a certain level of abstraction, you don't care what the root system looks like. I'm guessing you want some way to generate 3D fractal domains of the kind mentioned in this paper describing the calculation of fractal dimensions of root systems.
After pulling up this paper on fractal analysis of the efficiency of soil exploration by root systems, I found the SimRoot package that looks like it generates 3D root system geometries that may interest you. Unfortunately, they neglect to provide a way to download their package. However, their web site links to other packages that model root systems, such as PlantGL out of INRIA, which is open source.
Of course, once you create the geometry, you're going to have to figure out how to extract the relevant data in a compatible format and use it in PDE simulations. I leave that part up to you.
You could also borrow code from the folks who do NeuroML for dendrites. I can upload code that generates tubes from NeuroML if you so desire.
I think the answer from Geoff Oxberry is very good. It provides out-of-the-box solutions.
If you want to go on your own:
The mentioned L-systems may generate root-like structures if you provide the right rules. there is this book about "The algorithmic beauty of plants", but it doesn't cover root systems.
Diffusion limited aggregation processes may also generate root-like structures. If you aggregate spheres and once aggregated you produce the boolean union of the structure you will get the volume to mesh directly (smoothing needed, almost for sure).
As I said, not many solutions but maybe the ideas help you. If you implement something do not forget to release it with a free license! :D
Something that might be helpful for you is the following paper:
Olga Wildeotter: "An adaptive numerical method for the Richards equation with root growth", Plant and Soil, 2003
They only treat a 2D model and use a cellular automaton to simulate growth. It does not directly relate to your question, however.