Herbert W. Franke – historic Research paper from 1999 – published in:
“Die höhle, zeitschrift für karst- und höhlenkunDe”, 51. Jg., Vol. 2/2000 – Download of article
Digital Simulation and Visualization of Dripstone Caves
The formative processes involved in the formation of cave carbonates (also known as speleothems) are now known in principle [Dreybrodt 1980, Buhmann and Dreybrodt 1985, Dreybrodt and Franke 1987]. One of the most important insights is that flowstones and speleothems forming on the cave ceiling have a fundamentally different set of forms. The layered structure is particularly pronounced in flowstone [Franke 1956 and 1961]. What appears on the outside is nothing other than the surface of the uppermost layer. However, the shape of the underlying layers, which normally remains hidden from view, is also of scientific interest, as it reveals climate-dependent sedimentation processes, a kind of climate calendar. The position of the layers is also particularly important when it comes to finding the best places to take samples for chemical investigations or for dating purposes.
For various reasons, it is difficult or impossible to calculate the layer pattern, but computer-aided simulation can provide a good insight into the internal structure. The first experiments were carried out by Dreybrodt and Lamprecht [Dreybrodt and Lamprecht, approximately 1980], who used the special simulation program Simula for this purpose; as a result, they obtained cross-sections through stalagmites. The following describes how such a simulation of the layer pattern of stalagmites can also be carried out with the aid of conventional PCs. The software Mathematica used for this purpose also makes it possible to switch from the cross-sectional images to spatial perspective views. This assumes a rotationally symmetrical structure.
Following Dreybrodt’s example, it is assumed that crystal growth is always perpendicular to the base surface. This results in a system of orthogonal trajectories that can be represented in any approximation depending on the selected resolution. If rotational symmetry is assumed, the problem can be restricted to two dimensions. To capture the geometry, a decay curve is required, which describes the layer thickness deposited per time unit. Due to the many locally changing influencing factors, it is not possible to derive this; instead, it is necessary to look for the one that provides the best description. Since the settling of the carbonate is proportional to the imbalance between the carbon dioxide contained in the solution or in the air and this is equalized as the solution evaporates, the deposited layers are greatest at the axis and become thinner and thinner towards the periphery. This case corresponds to many other processes known in physics and chemistry and is normally described with the aid of a so-called negative exponential function with base e. While Dreybrodt opted for a linear function, a quadratic one is used here – since the solution spreads over a surface, a quadratic decay law is likely. The program described below makes it possible to use various decay functions on a trial basis and to check the results. Incidentally, it turns out that – as already recognized earlier [Franke 1956] – different functions do not lead to any significant differences; a stationary situation soon arises in all of them during growth – in the sense that cap-shaped layers of uniform form are placed vertically on top of each other.
The example images show that the simulation can also be adapted to changing situations, particularly a change in the solution supply rate, which is interesting because it reflects the humidity of the relevant climatic phase: The diameter of a stalagmite is proportional to the amount of solution supplied per unit time – regardless of its concentration. Decreasing diameters during growth, which lead to cone shapes, indicate a climatic phase of increasing dryness. Increasing diameters that build up cones at the top indicate an increase in humidity; however, this results in overhanging boundary surfaces of the stalagmite, for which the growth rules of the floor sinter are no longer decisive, but those of the ceiling formations. The consequence of this is that such formations do occur in nature but cannot be observed as they are covered by curtains.
Simplifications in Simulation
Simplifications are indispensable in computer simulations of natural processes; it then depends on the extent to which they describe the essential properties of the phenomenon under consideration. One simplification in the modeling of stalagmites is the assumed rotational symmetry. If only individual speleothem formations are considered and superpositions are disregarded, then this measure only has an effect in the initial phase of growth; the subsequent stationary form proves, as mentioned, to be a clear function of the parameters considered, regardless of the shape of the subsurface and the lowest layers. The model stalagmite can therefore be allowed to grow on a flat base surface and still correctly represent the shapes that are visible in natural caves. However, the program also allows the assumption of an uneven base as long as the rotational symmetry is maintained. However, care must be taken to ensure that the height coordinates of the support points fall steadily; otherwise, the case would occur in which the solution collects in depressions, so that the basic prerequisite for carbonate precipitation in flowstones – those from thin films – would no longer be given. In this way, it is also possible to convince oneself through simulation tests that the stationary forms are independent of the support form. It would also be possible to extend the program to the asymmetric case, but with a considerable increase in computational effort. Such programs could then also be used to represent ensembles of stalagmites with superimpositions, alternating source points, etc.
Another simplification concerns the accuracy of the representation, which in this case depends on the number of support points and the thickness of the superimposed layers. In principle, it is possible to calculate with arbitrarily fine orthogonal meshes and thereby approach the real state as closely as desired – however, this is soon limited by increasing computing time and memory requirements.
Description of the Program
The selectable initial conditions are expressed by the following parameters:
dic ……. thickness increase of the layer per time unit
pz …….. number of support points
rep ……. number of layers
xfolge … x-coordinates of the interpolation points
yfolge … y-coordinates of the interpolation points
In the iterative section (do loop), vertical lines are drawn from the interpolation points to the next generation of interpolation points. (editorial not: the complete code is available in the published version of the historic article from 2000). For this purpose, the thickness increase per time unit dicke[ ] is calculated, where j is the number of the current interpolation point. In the decay function dic * Exp[-(bogenlänge[ j ] / r)^2], r is equal to the radius of the stationary sections of the stalagmites. The part of the program introduced with if allows the adaptation to different formation conditions, e.g. with constant, continuously variable or graduated water supply velocity. In the final section, the interpolation points for the next run of the loop are brought into a suitable form. As the trajectories close to the axes tilt more and more outwards during the iteration, their distances increase so much that the lines describing the layer surface become too coarse. In a similar way to Dreybrodt, new trajectories are therefore interposed; in order not to increase the number of support points, the outermost of the series is then removed.
The result is cross-sections through stalagmites in which the layer structure and the direction of crystallization can be seen. Of course, each layer surface can be calculated and output individually, especially that of the outermost layer, which represents the surface of the stalagmite.
Photo-realistic Images of Speleothem-filled Caverns
The procedure described, the details of which can be seen in the program, can also be applied to other mathematical program systems, most of which also offer corresponding graphical visualization options.
Another aspect has recently been added to the scientific aspect of speleothems, namely the realistic representation of cave chambers. For educational purposes, for example, the possibility of depicting the growth of stalactites in fast motion would be extremely interesting, but there is also interest in computer-animated journeys through caves for films or computer games. This raises the fundamental question of the extent to which computer-generated reproductions of natural objects must correspond to the given set of forms. In the case of mountain landscapes, for example, one is usually content with a fractally structured relief – a simplification that even some laymen find inauthentic. Clifford Pickover, who created very impressive computer-generated simulations of cave spaces [Pickover 1998], also limits his representation to the mathematically and computer-graphically much easier to grasp ceiling shapes and depicts the floor shapes as mirror images of the ceiling shapes.
In terms of programming, it is easy to import the stalagmites created with the Mathematica software into computer-generated images of cave rooms and thus also realistically depict flowstones. The software Bryce, which can be used to create realistic terrain forms – e.g. eroded ground sections – offers good prerequisites for this. It also offers a wide range of colors and textures that can be assigned to individual objects or object groups. Finally, light sources can be positioned in the three-dimensional, perspective rooms in order to illuminate the scene well – a task that is strikingly similar to that of a cave photographer looking for the most favorable positions for the flash bulbs.
Editorial Notes
The article with the complete Mathematica code can be downloaded here:
Herbert W. Franke: “Computersimulationen zum Schichtenbild der Stalagmiten” , Die Höhle, Zeitschrift für Karst- und Höhlenkunde, Vol. 2/2000.
More informations about coding dynamic processes with Mathematica at the border of art and science:
Herbert W. Franke: “Animation mit Mathematica”, Springer, Berlin, New York, Tokio 2002
Franke has been working on the question of determining the age of sediments in caves since the early 1950s. His first research paper on the subject Altersbestimmung von Kalzitkonkretionen mit radioaktivem Kohlenstoff”, Naturwissenschaftliche Rundschau, Springer-Verlag Berlin, Göttingen, Heidelberg, Vol. 22 (1951)