Transport of Geothermal Fluids along Dikes and Fault Zones

: Field observations of active and fossil natural geothermal fields indicate that geothermal fluids are primarily transported along dikes and fault zones. Fluid transport along dikes (com-monly through fractures at their margins) is controlled by the cubic law where the volumetric flow rate depends on the aperture of the fracture in the 3rd power. Dikes (and inclined sheets) also act as heat sources for geothermal fields. In high-temperature fields in volcanoes in Iceland dikes and inclined sheets constitute 80–100% of the rock at crustal depths of 1.5–2 km. Holocene feeder-dikes are known to have increased the activity of associated geothermal fields. Fault zones transport geothermal fluids along their two main hydromechanical units, the core and the damage zone. The core is comparatively thin and primarily composed of breccia, gouge, and clay and related low-permeability porous materials. By contrast, the fault damage zone is characterised by fractures whose frequency is normally highest at the contact between the core and the damage zone. Fluid transport in the damage zone, and in the core following fault slip, is controlled by the cubic law. During non-slip periods fluid transport in the core is primarily controlled by Darcy’s law. Secondary mineralisation (forming mineral veins and amygdales) tends to reduce the fault-zone permeability. Repeated earthquake activity is thus needed to maintain the permeability of fault zones in active natural geothermal fields.


Introduction
Geothermal fluids, like other crustal fluids, are transported in two basic ways. One is along rock fractures, the other along interconnected pore networks. The main crustal fluids are oil and gas, magma, groundwater, and geothermal fluids. Here the focus is on geothermal fluids, which include hot water (geothermal water), steam and, at the surface, mud. Pores occur in all rocks and are of different origin. Some are primary (such as vesicles in volcanic rocks), while others are secondary (such as those formed by exsolution). Similarly, some fractures are primary (such as cooling or columnar joints in igneous rocks) while others are secondary (such as faults generated in later tectonic stress fields).
In volcanic areas-and formerly volcanically active areas-the main channels for the transport of water in the solid crust are dikes and fault zones; this applies to the groundwater of volcanic islands such as Tenerife, where dikes are important collectors of water [1], and also to continental desert areas such as Sudan [2]. The effects of faults in collecting and transporting groundwater are well known and apply everywhere [3][4][5][6][7][8][9][10][11][12].
In Iceland, it was established long time ago that dikes conduct geothermal water [13]; this channelling is particularly noticeable in eroded palaeo-geothermal fields where mineral veins, primarily located at the margins of the dikes but also inside the dikes, indicate the former paths of the geothermal fluids ( Figure 1). Faults are also well-known channels for geothermal fluids. The channelling of water along active faults is abundantly clear in paleo-geothermal fields, such as in sedimentary basins ( Figure 2). Additionally, geothermal fields in volcanic areas remain active for a long time primarily if seismogenic faulting occurs from time to time. Faulting maintains the permeability and without repeated fault slip the activity of a high-temperature geothermal field gradually decreases. For example, the erupting geysers in Iceland increase or renew their activity following earthquakes [14][15][16]. Low-temperature fields outside active volcanic areas, however, may be active for long periods (and often associated with dikes) without renewed seismogenic faulting. Figure 1. Basaltic dike, providing paths for geothermal fluids, in the damage zone of the Husavik-Flatey Fault, a transform fault partly exposed on land in North Iceland [10,14]. The mineral veins dissect the central part of the dike but occur mainly at the dike contacts with the host rock (one vein is indicated). There stress concentrations encourage fracture formation [17,18], including the propagation of hydrofractures (seen here as mineral veins). View east, the dike is about 0.4 m thick (the height of the red-black notebook at the right margin of the dike is about 15 cm).
Dikes collect and channel geothermal water primarily because of their commonly low dike-transverse and high dike-parallel permeability [1,2,13,17]; they act as barriers to dike-transverse flow while acting as conduits for dike-parallel flow, primarily along their contacts with the host rock ( Figure 1). Many fault zones act as barriers to transverse flow while allowing fault-parallel flow. However, an entire active fault zone may conduct fluids during certain periods, particularly following a large fault slip in the core. The fluid transport is rarely along the boundary between the fault zone and the host rock, primarily because there is normally no clear mechanical contact at that boundary [9,17]. A fault zone consists of two main hydromechanical units, namely a fault core in the centre and fault damage zones on either side of the core (Figure 2). Fluid transport in the damage zone is primarily through fractures, whereas transport in the core is commonly through porous media (the core material is mostly breccia, gouge, and other particulate/porous materials). Fluid transport in geothermal fields is thus partly through porous media, and partly through factures. Fluid transport in porous media is controlled by Darcy's law, whereas fluid transport along fractures is controlled by the cubic law [17].
The aim of this paper is to provide an overview of the effects of dikes and faults on fluid transport in geothermal fields. The paper first summarises the relevant extensive field studies on the structures of dikes and faults in relation to fluid transport, particularly the transport of geothermal fluids. Then the mechanisms of fluid transport along and within dikes and fault zones are discussed, with particular reference to Darcy's law and the cubic law. The focus is on analytical models of fluid transport, but the use of numerical models is reviewed at the end of the paper. Examples of palaeo-geothermal fields are also presented where the effects of dikes and faults on geothermal fluid transport can be observed. Additionally, there is a discussion of the effects of active seismogenic faulting on the maintenance of permeability in geothermal fields. The last sections of the paper focus on application of the field and theoretical results to explain observed active geothermal fields and the implications of the results in a wider context. The damage zone has a thickness of a few metres and is much thicker than the core.

Dike Structure
The term dike referred originally to a solidified or 'frozen' magma-filled fracture that cuts through the host-rock layers ( Figure 3). Because the strike-dimension (length) and the dip-dimension (height) of the dike are many times greater than its thickness (Figure 3), the dike is also referred to as a tabular intrusion, or a sheet-like intrusion. And because the dike dissects the layers, it is additionally described as a discordant intrusion, in contrast with sills, which are concordant tabular or sheet-like intrusions; this means that sills are commonly emplaced at the contacts between layers/strata ( Figure 4). The third type of tabular intrusions is the inclined sheet, a structure mostly found close to the shallow chambers of polygenetic (central) volcanoes ( Figure 5). All these intrusive structures act as heat sources and (when completely solidified and at much lower temperatures than initially) may conduct geothermal fluids; they are described in detail elsewhere [15,17,18]. Here, however, the focus is on the description of dikes, while inclined sheets are also mentioned, in relation to their function as conduits for geothermal fluids.

Figure 3.
Part of a basaltic dike exposed in the walls of the caldera of Santorini, Greece. Dike is a discordant sheet-like intrusion, dissecting the layers at close to a right angle. Indicated are the directions of the strike-dimension (length) and dip-dimension (height) of the dike. Following its emplacement, a dike may act as a heat source for geothermal fluids in addition to acting as a conduit (when the dike has cooled down) for the fluids. . Basaltic sill in East Iceland exposed at the depth of about 1200 m below the original surface of the lava pile. The sill is about 30 m thick and has a well-developed system of (mostly vertical) columnar joints. Sills may develop into shallow magma chambers [18] and act as heat sources in geothermal fields.
Most dikes, sills, and inclined sheets are extension fractures; this means that their thickness (originally, roughly the aperture of the magma-filled fracture) is entirely the result of opening displacement in the direction of the minimum compressive (maximum tensile) principal stress, σ3. Furthermore, they propagate in a direction parallel with the maximum compressive principal stress σ1, the trajectories (ticks) of which can therefore be used to forecast the likely propagation paths of dikes (and inclined sheets and sills); this has been done, using numerical models, for numerous dike paths [18].
All dikes are segmented ( Figure 6). In some dikes, the offsets (lateral separations) are large in comparison with the dike thickness (Figure 7), while in others the offset is so small that the segments are in physical contact ( Figure 8). While in the small offsets there may be thin igneous veins or dikelets connecting the segments, the connection between the segments need not occur in the studied outcrop but rather at a different location along the segments. If, in addition, there is large underlap between the segments ( Figure  6), then the dike cannot act as a continuous barrier to the flow of geothermal fluids, because the fluid will find paths, escape, through the 'gaps' in the dike wall, that is, between segments where the underlap and offset are comparatively large.

Figure 5.
Inclined basaltic sheets in a fossil high-temperature geothermal field (the alteration is indicated) in the mountain Esja in Southwest Iceland. View east, the sheets are mostly about 1 m thick and dipping towards their source, a shallow magma chamber of the fossil Videy Volcano, which was active at about 2.5 Ma [15]. Above the swarm of sheets is a pile of subhorizontal basaltic lava flows. The sheets provided heat and acted as channels for geothermal fluids.
If, however, the offsets are comparatively small and showing overlap ( Figure 8) rather than underlap, the dikes-particularly thick dikes-can act as barriers to the flow of geothermal fluids ( Figure 9). The dikes then tend to collect the fluids at their margins, especially at the margin facing the upstream part of the fluid flow [2,17,18]. How effective dikes are in collecting geothermal fluids depends partly on their thickness and (commonly) low transverse permeability, but also on their attitude, that is, strike and dip, in relation to the hydraulic gradient controlling the flow of the geothermal fluids. Figure 6. Schematic illustration of two dike segments (one whose dip-dimension is denoted by a) with an offset b and underlap c. While the present illustration is for dike segments in a vertical section (a cliff section), it could also be viewed as being in a lateral section, in which case a would be the strike-dimension of one of the segments. When either b or c are significant in a dike, it does not function well as a collector of geothermal fluids because the fluids tend to flow through the 'gaps' in the dike wall, generated by b and/or c.

Figure 7.
Offset dike segments of a basaltic dike in a vertical section as seen in the caldera walls of Santorini, Greece. The maximum thickness of the upper (left) segment is about 1 m. The segments have a large lateral offset in comparison with their thicknesses while the underlap is close to zero. The large offset (and the thinness of the segments) means that these segments would not have acted as a major collector of geothermal fluids-the fluids would normally find a path through the 'gap' made by the offset.
While dikes can have a variety of strikes in a given active or recently active volcanic area-in such areas high-temperature geothermal fluids are most common-the strike of regional dikes-dikes injected from a deep-seated reservoir rather than a shallow magma chamber [18]-often follows (crudely) a normal (Gaussian) distribution (Figure 10a). For a given regional hydraulic gradient that partly controls the flow of geothermal fluids-as is common in low-temperature fields such as in Iceland-the dikes striking about perpendicular to the hydraulic gradient are, other things being equal, most likely to collect such fluids. Dike dip has also effects on the collection of fluids. Regional dikes tend to be subvertical and close to perpendicular to the hosting rock layers (Figures 9 and 10b). Those dikes that dissect the layers at close to right angles and are also striking perpen-dicular to the general dip of the layers, and thus the hydraulic gradient, are most likely to collect, and act as channels for, geothermal water. Figure 8. Basaltic dike segments A and B (in Tenerife, Canary Islands); they are offset and with some overlap but in physical contact. The maximum thickness of segment A is 1.3 m, and the host rock is pyroclastic. Since the segments are in physical contact and overlapping, they form a continuous 'wall' that acts as a barrier to the transverse flow of water. The segmented dike could thus have acted as a collector of and conduit for geothermal fluids.
Dike thickness is also important. It normally follows a heavy-tailed size distribution and, often, a power-law distribution ( Figure 10c); this means that in a given dike area or swarm most of the dikes are comparatively thin, while a few ones are very thick. The thick dikes are the ones most likely to collect and conduct geothermal water, that is, to act as sources for hot or warm springs; this is partly because thick dikes provide a more effective barrier to transverse flow (flow across the dike thickness), and partly because thick dikes tend to be of larger lateral dimensions. For example, aspect ratios, that is, the ratios of strike-dimension (length) to the thickness of regional dikes in Iceland vary from 300 to 1500, with an average of about 803 [19]. In Sudan, the same aspect ratios of basaltic dikes range from 243 to 2331, with an average of 862 [2]. The most detailed data on the aspect ratios of regional dikes is from swarms in Sinai (Egypt). For three regional swarms of basaltic dikes [20], the strike-dimension/thickness aspect ratios range from 594 to 4312, the average aspect ratios in the swarms being 1523, 1540, and 2487. For felsic dikes in the same swarms, the aspect ratios range from 330 to 2487 with an average aspect ratio of 609 [20]. The aspect ratios of these regional dikes, while showing a variation of about an order of magnitude, have averages that mostly range from about 800 to about 2500, so by a factor of about 3.
All the above regional dikes are exposed at significant depths below the original surface of the associated volcanic areas-from depths of many hundred metres to some kilometres. The aspect ratios tend to be much smaller at shallower crustal depths. For example, Becerril et al. [21] analysed many feeder-dike in El Hierro (Canary Islands) and Kusumoto et al. [22] studied feeder-dikes and non-feeders at very shallow depths below the surface in the volcano Miyakejima in Japan. For the dikes in El Hierro the aspect ratios range from 141 to 227, with an average of 171. In Miyakejima the aspect ratio is not for the strike-dimension/thickness ratio, but rather for the dip-dimension (as seen in the caldera walls of the volcano) to the dike thickness; these aspect ratios range from 61 to 246, with an average value of 136. Thus, the aspect-ratio results from El Hierro and Miyakejima are very similar and are much lower than for the regional dikes above, mainly because the dikes are exposed at a very shallow depth-at the surface in the case of El Hierro-where the (mostly young volcanic) layers hosting the dike are very compliant, that is, with a low Young's modulus (stiffness). Generally, however, in a particular area or swarm, the comparatively thick dikes tend to be of larger lateral and/or vertical dimensions than thin dikes and therefore have a greater chance of collecting and conducting geothermal fluids. Figure 9. Part of a swarm of regional basaltic dikes (10 of which are indicated by numbers) dissecting a 400 m thick pile of (mostly) basaltic lava flows in Southeast Iceland. Many of the dikes are thick (up to 13 m) and made of segments in physical contact with each other, thereby offering a barrier to transverse flow of water-which is down the dip of the lava flows (the direction of the hydraulic gradient), as indicated by the arrows. Many of the dikes may have collected and acted as conduits for geothermal fluids of initially high-temperature and later low-temperature fields. View northeast, the valley floor is at about 2 km beneath the original top of the lava pile. Many dikes are formed through repeated magma injections into the same fracture. When the magma composition in the injections is the same or very similar, say basalt, the dike is referred to as multiple. When, however, the injections are of widely different composition, such as rhyolite and basalt, the dike is referred to as composite [18]. Multi-ple/composite dikes are very common and can have great effects on the mode of transport of geothermal fluids.
The main difference as regards geothermal fluid transport along single dikes (Figures 3 and 8) and multiple/composite dikes ( Figure 11) has to do with the potential flow channels associated with the dikes. For a single dike, while some geothermal fluids may be transported through the centre-the last part of the dike to solidify-the main along-dike channels are normally at the margins of the dike. The dike functions as an elastic inclusion (that is, with different mechanical properties from those of the host rock) and concentrates stresses at its contacts with the host rock, at the dike margins. By contrast, a multiple/composite dike has many contacts inside the dike body ( Figure 11); these are contacts between individual sheet intrusions or, as they are also named, columnar rows. The contacts are common paths for the along-dike migration of geothermal fluids. The migration of the fluids into the dike is then mainly through the columnar joints, which are generally perpendicular to the margins of the individual injections, the individual columnar rows. Figure 11. Basaltic dike in North Iceland composed of seven rows of columns (numbered), that is, a multiple dike. View north, the total exposed thickness of the dike is 20 m (a person provides a scale). The contacts between the columnar rows may provide conduits for geothermal fluids.
The columnar or cooling joints contribute significantly to dike permeability and storage ( Figure 12). The cooling joints form networks that generate fracture porosity and, thereby, storage for geothermal fluids. Additionally, the joints are connected to the contacts between columnar rows (injections) in multiple dikes. That connection, however, depends on the time between successive injections. If the time between injections is long in relation to the thickness of the injections, then a chilled selvage (glassy surface) may form between the columnar rows [18]. The glassy surface tends to reduce the flow of geothermal water from one columnar row to another and into the contacts between the rows. If, however, the rows form in a rapid succession, so that the earlier rows are still hot when the next one is injected, then there will be no chilled selvage at the contacts and flow between the rows and into their contacts will be easier.
The structure of inclined sheets ( Figure 5) is basically the same as that of dikes. Inclined sheets are very common in polygenetic (central) volcanoes and act as heat sources for geothermal fields. In addition, they, like dikes, may collect geothermal fluids and channel them to the surface. There are two principal differences between local swarms of inclined sheets and those of regional dikes, however. One is that, as the name implies, the sheets are commonly neither vertical nor horizontal but rather inclined, mostly with dips from 20° to 50° [18]. The local swarms, however, are not composed only of inclined sheets, but also of radial dikes and sills ( Figure 13). The radial dikes are mostly subvertical but, like the inclined sheets, they are generally much thinner than the regional dikes-commonly about a metre or less. It follows that individual sheets and radial dikes have strike-dimensions and dip-dimensions that are only a fraction of similar dimensions for regional dikes. Thus, the potential of the sheets/radial dikes to collect geothermal fluids from a large region, in comparison with that of regional dikes, is much more limited. However, the sheet swarms are always connected with shallow magma chambers ( Figure 13) that act as major sources of heat for high-temperature geothermal reservoirs and fields ( Figure 5). Thus, the availability of high-temperature geothermal fluids is normally much greater inside the sheet swarms than inside regional dike swarms. The sources of major high-temperature geothermal fields can be studied through drilling into active volcanoes and through field studies of deeply eroded fossil polygenetic/central volcanoes. Field studies of fossil central volcanoes show that the main heat sources of geothermal fields associated with such volcanoes are shallow magma chambers and local sheet swarms ( Figure 14)-as is illustrated in Figure 13. The principal heat source is the long-lived shallow chamber, exposed as a pluton in a fossil volcano, but the inclined sheets and radial dikes act as heat sources at shallower levels [18]. Shallow magma chambers worldwide are commonly with roofs (top contacts between magma and rock) at depths from 0.5 km to 6 km [18]. The tops of observed fossil shallow chambers exposed in Iceland range from 0.5 km to about 2 km in depth, but 2 km is also the deepest erosional level in Iceland.
Drilling into high-temperature geothermal fields in active volcanoes, such as the Hengill Volcano in Southwest Iceland, shows that inclined sheets and radial dikes increase in frequency with depth and reach 80-100% of the rock at depths of about 1.5-2 km [24,25]. The Hengill Volcano and the surrounding areas, primarily Hellisheidi, provide energy for the largest geothermal power plants in Iceland. Studies of the geothermal fields in Hengill and Hellisheidi indicate that the three most recent feeder-dikes (dikes supplying magma to a volcanic eruption) in the area are major contributors to the currently active geothermal fields. The three feeder-dikes were formed 2000, 5000, and 9000 years ago. In particular, there is clear connection between the feeder-dikes formed 2000 and 5000 years ago and increased or renewed geothermal activity at the location of the dikes [25]; this renewed activity is presumably related to the feeder-dikes both as recently added heat sources for the geothermal fields as well as the feeders acting as new conduits of geothermal fluids. Figure 13. The main sheet-like intrusions are dikes, inclined sheets, and sills. Common sources are deep-seated reservoirs, located in the lower crust or crust-mantle boundary, and shallow chambers, located in the upper crust and, in turn, supplied with magma from (one or more) deep-seated reservoirs [18]. Some dikes are injected directly from deep-seated reservoirs and tend to be thicker and more primitive in composition than the radial dikes and inclined sheets injected from shallow chambers. The inclined sheets are primarily characterised by their dips, while sills tend to be parallel with the layering of the host rock. The time for solidification of a sheet-like intrusion depends on its thickness in the 2nd power, so that thicker intrusions act as heat sources for geothermal fields much longer than thin ones. Figure 14. Swarm of inclined sheets and radial dikes close to the contact (indicated) with a fossil magma chamber (a gabbro pluton) in the extinct Geitafell Volcano, Southeast Iceland. At the contact about 80% of the rock is composed of inclined sheets and radial dikes, but the intensity of the swarm falls off rapidly with distance from the contact. The top of the gabbro pluton/fossil chamber is at about 2 km beneath the initial top of the Geitafell Volcano when it was active. The chamber and the sheet-intrusions acted as heat sources for high-temperature geothermal fields, as indicated by numerous mineral veins that dissect the rocks.

Fault-Zone Structure
The internal structure of a fault zone is characterised by two hydromechanical units: a fault core and a fault damage zone [26]. The internal structure is the same for all fault types, namely normal, reverse, and strike-slip [9,10,17,26]. The main difference is that in dip-slip faults (normal and reverse) the thicknesses (widths) of the zones need not be symmetric about the fault plane. More specifically, the thicknesses of the core and the damage zone in the hanging wall of a dip-slip fault zone may be greater than the thicknesses of the same in the footwall ( Figure 15); this is because the hanging wall is subject to greater loading (stress) than the footwall. By contrast, in a strike-slip fault zone the loading is normally the same on either side of the fault plane, so that the deformation, hence the thicknesses of the core and the damage zones, is the same on either side of (is symmetric about) the fault plane ( Figure 16).
Field observations show that fault zones normally consist of a fault core and a fault damage zone (Figures 15 and 16). The core accommodates the bulk of the fault displacement and it is also named the fault-slip zone [26]. There are many fractures in the core, but it is characterised by fault gouge, breccia, and other particulate materials. Thus, the core rock is crushed and altered into porous and fragmented materials ( Figure 17) that, in active faults, may behave as ductile (and creep) except at very high strain rates such as during seismogenic faulting. In the core, there are commonly numerous veins filled with secondary minerals spaced at centimetres or millimetres, that form dense networks ( Figure 2). These networks, when actively transporting geothermal fluids, resulting in the core having a granular-media structure at the millimetre or centimetre scale, thereby supporting flow of geothermal fluids in the core being modelled using Darcy's law [17].
In large fault zones, the core thickness is commonly from several metres to a few tens of metres ( Figure 17; [27][28][29][30][31]). Transform faults and other very large fault zones sometimes develop more than one fault core and many damage zones [32,33]. The core permeability is usually low during most interseismic (non-slip) periods, in which case the core commonly acts as a barrier to, particularly transverse, geothermal fluid transport. However, during periods of high strain rates, such as are associated with fault slip, the core may fracture and very greatly increase its permeability [9].
The damage zone is also known at the transition zone [26]. The damage zone consists partly of lenses of breccias but primarily of sets of extension fractures and, to a lesser degree, shear fractures [14,34,35]. Many, and presumably most, of the extension fractures are hydrofractures that eventually become mineral-filled veins (Figures 1 and 2). The extension fractures and faults generally make the damage zone much more permeable than the fault core. Laboratory measurements indicate that hydraulic conductivity in the damage zone may be several orders of magnitude greater than that in either the fault core or in the host rock [36,37]-results that have great effects on the potential of the core versus that of the damage zone for transporting geothermal fluids.
The boundary between the damage zone and the host rock is usually diffuse. The main difference is that the fracture frequency in the host rock-also known as the protolith [26,37]-is normally less than that in the damage zone ( Figure 18a; [38]). In the outermost parts or subzones of the damage zone, namely on approaching the host rock, this difference in frequency may, however, be small ( Figure 18b). The boundaries between the fault core and the damage zone are generally much sharper than those between the damage zone and the host rock. All these boundaries vary along the length of the fault and change, in time and space, with the evolution of the fault zone ( Figure 19). The core is mostly composed of (commonly low-permeability) breccia, gouge, and clay, whereas the damage zone is characterised by fractures whose frequency increases towards the contact with the core. In dip-slip faults the damage zone in the hanging wall is normally thicker (subject to greater strain) than that in the footwall [17,26,27].
The field descriptions presented here of the fault core and damage zone focus on outcrop-scale fault zones (Figures 15-19). It should be made clear, however, that that cores and damage zones are seen in much smaller fault zones. For example, laboratory experiments on small rock samples may produce very similar units, that is, a thin core and a thicker damage zone; in the latter, the frequency of fractures changes irregularly, but generally decreases, with increasing distance from the core [38].

Figure 16.
Internal hydromechanical structure of a strike-slip (here a dextral) fault. As for a dip-slip fault, the damage zone is normally much thicker than the core. However, in a strike-slip fault the thickness of the damage zone on either side of the core is usually the same. During non-slip (no earthquake) periods, the transport of geothermal fluids is normally easiest in the part of the damage zone that is close to its contacts with the core, namely where the fracture frequency is highest. Figure 17. Fault core and parts of the damage zone of the Husavik-Flatey Fault, an oceanic transform fault exposed partly on land in North Iceland [33]. View west, the core is 10 m thick and mainly composed of breccias of crushed basaltic lava flows. The damage zone, by contrast, is characterised by tilted basaltic lava flows (initially with zero dip but now dipping around 40°), as indicated by white arrows, and dense sets of mineral veins. The person provides a scale. Figure 18. Fracture frequency in a damage zone as a function of distance from the contact between the core and the damage zone. (a) Fracture frequency decreases (irregularly) with distance from the core in a small fault in West Norway. The host rock (gneiss) takes over from the damage zone at a distance of about 10 m from the core. (b) Fracture frequency at specific distances from the core of a large normal fault in West Norway. Closest to the core/damage zone contact there are 24 fractures/metre which falls to 14 fractures/metre at a distance of 20 m, to 4 fractures/metre at a distance of 80 m, and then stays the same to a distance of 110 m from the core/damage zone contact. The fracture frequency at 80 m (and thus 110 m) is regarded as that of the host rock [39]. Neither fault shows a sharp boundary between the damage zone and the host rock. Figure 19. The fault core and the damage zone both become thicker as the fault develops with time (a-d). As indicated in d, Young's modulus (stiffness) is inversely proportional to the fracture frequency in the damage zone and reaches its minimum value in the damage zone at its contact with the core. Variations in fracture frequency result in different Young's moduli within fault zones, thereby in their having local stresses (controlling fault slip and permeability) that may be widely different from the associated regional stress fields [17,27].
As regards the potential for fault zones to transport geothermal fluids, the field observations summarised above indicate as follows. The damage zone is generally the main conduit for transport of geothermal fluids. Measurements on laboratory specimens indicate a permeability of the damage zone that is up to 10,000-times higher than that of the core and the host rock [36] While the laboratory specimens from fault cores yield hydraulic conductivity values that may be similar (but somewhat lower than) the bulk in-situ values, the laboratory values of hydraulic conductivity for the damage zones are likely to be considerably lower than the corresponding bulk in-situ values. This is so because the larger, highly conductive fractures that are common in the damage zone of major fault zones (Figures 2 and 15-18) are not represented in the small, relatively non-fractured laboratory specimens. This conclusion finds support in reported in-situ measurements giving fault-zone permeabilities up to 1000-times greater than the maximum laboratory values of the same fault rocks [36]. It follows that the in-situ permeability difference between the core and the damage zone may be greater than the cited laboratory values would indicate.

Mineral Fillings and Permeability
Circulation of geothermal fluids along dikes, fault zones, and through crustal layers results in the precipitation of secondary minerals; these accumulate in cavities and holes (pores) in the rock, as well as in fractures. Mineral accumulations in the pores of igneous rocks, such as in vesicles, are named amygdales. Mineral accumulations in fractures are named mineral veins. Amygdales and mineral veins are both are indication of palaeo-geothermal fluid transport and gradual reduction in permeability ( Figure 20).
Amygdales form in a lava pile everywhere where geothermal fluids circulate. Common secondary minerals that form amygdales are calcite, quartz, and zeolites. As the amygdales accumulate in the vesicles and other pores in the host rocks of the geo-thermal reservoir, the porosity and, therefore, the storage of the reservoir gradually diminishes. The porosity can, however, be renewed partly through new sheet injections (dikes, inclined sheet, or sills) and earthquakes.
Sheet intrusions (Figures 1, 4, 5 and 14) not only provide new, temporary heat sources at shallow depths but also numerous new vesicles that in due course become filled with secondary minerals to form amygdales. Additionally, sheet intrusions develop cooling or columnar joints (Figures 11 and 12) and contacts ( Figure 1) that increase permeability and storage for geothermal fluids [18]. Earthquakes reactivate or generate new fractures and also shake the ground, with the overall effect that the porosity and permeability of the rocks increases. In fact, many geothermal fields, such as in South Iceland, maintain their permeability primarily through repeated earthquake activity [15]. When secondary minerals precipitate in fractures the resulting structures are referred to as mineral veins (Figures 1, 2 and 20). The veins commonly form dense networks that are, to a large degree, simultaneously active ( Figure 21). This means many or most of the veins may be transporting geothermal fluids at the same time. This is in contrast with networks or swarms of dikes and inclined sheets, as are common in central volcanoes (Figures 5,9,13 and 14), which are, as a rule, not active (fluid) at the same time (normally only one or a few sheets/dikes are fluid at the same time). Thus, the magma in a given dike/inclined sheet normally solidifies before a new dike/sheet is injected-even in multiple and composite dikes ( Figure 11). Thus, dikes and sheets do not normally form networks of simultaneously fluid and interacting intrusions.
But mineral veins do. While active as hydrofractures, the veins commonly form large networks that transport geothermal fluids from the heat sources (such as magma chambers or sheet intrusions) towards or to the surface ( Figure 21). Cross-cutting relationships ( Figure 22) and the orientation of fibres show that the great majority of mineral veins are extension fractures; this applies even to mineral veins in fault zones, where detailed field studies show that the great majority of veins are extension fractures [10,14,34,35].

Figure 21.
Dense network of mineral veins, representing a fossil high-temperature geothermal system, in the damage zone of the Husavik-Flatey Fault in North Iceland. The vein networks were measured in profiles (scan lines) such as the one along the white measuring tape. The heat sources of the high-temperature systems in this fault zone were primarily thick regional dikes [40]. The person provides a scale.
As an example of the size distributions of mineral veins, we may consider the aperture/thickness and length of veins in a fault zone in North Iceland, namely the Husavik-Flatey Fault, an oceanic transform fault initiated some 8 Ma ago and partly exposed on land [10,14,33]. The veins are exposed at coast of North Iceland at the depth of about 1500 m below the original top of the 10-12-Ma-old lava pile-the erosion is partly due to the ice sheets of the Pleistocene and partly due to the action of the sea. A total of 1717 mineral veins were measured, most of which are composed of quartz, chalcedony, and zeolites (Figures 1, 20 and 21). Of 829 cross-cutting relationships among the veins, 79% indicate pure extension, while 21% show evidence of strike-slip movement, that is, are partly shear fractures-in agreement with the general observations above. The common vein intensity in several measured profiles within the networks ( Figure 21) is 5-10 veins per metre. Most of the veins strike either subparallel or sub-perpendicular to the strike of the Husavik-Flatey Fault. The thickness-size distribution of 1333 mineral-vein thicknesses follows a heavy-tailed distribution (Figure 23), more specifically a power-law. The veins range in thickness from 0.1 mm to 85 mm within the networks (Figure 21), the thickness of the great majority of the veins being less than 2.5 mm. Outside the networks, however, some veins reach thicknesses as large as 100 mm (Figure 24). Many veins show evidence of repeated injections of geothermal fluids, that is, they were reactivated several times, as is commonly observed [41,42]. For the veins that are free to grow, the growth is mostly self-similar. This means that the aspect ratio (between length or strike-dimension and thickness) is maintained even if the final dimensions of the vein are the result of many fluid injections.
Most of the mineral veins have restricted lengths (strike-dimensions). This means that they are not free to grow laterally because they meet other veins or discontinuities within the rock body that may arrest or deflect them, thereby restricting their expansion ( Figure 21). Yet, in the dataset of 1717 veins, 384 veins were both unrestricted as regards length, so that both lateral ends were exposed and seen as tapering away (thinning out) without meeting other veins or discontinuities, and pure extension fractures (an excellent example of such non-restricted extension veins, from a different area, is seen in Figure  25). The measured veins range in length from 2.5 to 400 cm and have thicknesses, that is, palaeo-apertures of 0.1 to 85 mm.
There is a linear correlation between the lengths and aperture/thickness of the 384 veins ( Figure 26). The coefficient of determination 66 . 0 2 = R , which means that 66% of the variation in aperture can be explained through variation in length or strike-dimension. As discussed below, this linear correlation between aperture and length for extension fractures modelled as mode I crack is a direct consequence of the appropriate fracture-mechanics equations. That the coefficient of determination is not higher than observed is because some of the factors that enter those equations themselves vary, in particular the Young's moduli of the host rocks.     (Figures 22a and 25) is also heavy-tailed where a power-law or a log-normal function fits the distribution [35]. Similarly, the size distribution of lengths or strike-dimensions of these veins also follows log-normal or power-law functions; these veins are similar in size to those in Iceland. Thus, the calcite veins reach a maximum length of 10.3 m and a maximum thickness of 28 mm. The coefficient of determination for the linear correlation between length and thickness of the veins in the Bristol Channel is 74 .
, so that 74% of the variation in aperture or thickness can be explained in terms of variation in vein length [35].

Transport of Geothermal Fluids along Dikes
Fluid transport along dikes is partly inside the dikes-particularly along the contacts between columnar rows on sheets in multiple dikes ( Figure 11)-but primarily along the contacts between the dike and the host rock ( Figure 1). The flow of geothermal fluids along these contacts, as well as along columnar joints inside a dike (Figure 12), is thus modelled as transport along fractures.
The two basic models for the flow of fluids along rock fractures are (1) the elliptical-crack-opening model, and (2) the parallel-plate model. Both models are useful-in most non-restricted veins the aperture/thickness distribution is as in a flat ellipse ( Figure  27)-and are reviewed here. Both derive from the Navier-Stokes equations [43,44].

Elliptical Crack/Conduit Model
Elliptical models for the flow of geothermal fluids may be used for fractures and contacts, but also for conduits that supply water to geysers and hot springs in general. Consider first the real shape of a mineral vein. The veins in Figure 27 show typical variation in thickness/palaeo-aperture seen in many mineral veins and other extension fractures [45]. The thicknesses of the veins vary as in a very flat ellipse. In these cases, and many other similar thickness/aperture variations of extension fractures, we can thus either use the elliptical crack model, or the parallel-plate model. Here we consider first the elliptical crack/conduit model. An elliptical cross-section of a hydrofracture ( Figure 28) or a conduit is given by: The circular conduit model applies also approximately to the uppermost conduits of some erupting geysers [16,17].
Equations (2)-(5) apply to rocks that behave as rigid, that is, do not deform however much the fluid pressure in, or stress on, the fracture network and its source changes. While widely used in the modelling of fluid transport in the Earth's crust [17], the assumption of a rigid behaviour of crustal rocks is generally not a realistic one. In particular, that assumption means that the important contribution of buoyancy to the driving pressure (overpressure) of fluid transport up (vertically or inclined) through hydrofractures is ignored.
When modelling the behaviour of the host rock as elastic, as the brittle (seismogenic) crust normally behaves to a first approximation, the difference in density between that of the fluid and the rock is reflected in buoyancy contributing to the driving pressure (overpressure). For fluid transport in a vertical hydrofracture hosted by crustal rocks that behave as elastic, the overpressure o p is given by [18]: where e p is the excess pressure in the hydrofracture source at the time of rupture and start of fluid transport up towards the surface, r ρ is the average host-rock density, f ρ is the average geothermal fluid density (unit: kg m −3 ), g is acceleration due to gravity (unit: m s −2 ), and h is the dip-dimension (height) of the hydrofracture as measured from its point of contact with its source. The symbol d σ here denotes the differential stress ( may be taken as zero and is ignored in the derivation below [18]. Using Equation (6), and thus assuming the rocks to behave as elastic during the fluid transport, Equation (2) can be rewritten on the form: where all the symbols have been defined above. Equation (7) can then be recast as: Because most mineral veins are much longer than they are thick, it implies that the strike-dimension of a vein/hydrofracture is normally hundreds of times the size of its aperture ( Figure 26). Moreover, theoretical considerations and measurements of mineral veins show that they are commonly approximately flat ellipses where the aperture/thickness varies little along the greater part of the strike dimension or surface length ( Figure 27). It follows that a fluid transport up through a mineral vein can be modelled using the parallel-plate model, to which we turn now.

Parallel-Plate Model
The parallel plate model is the basis of the so-called cubic law in hydrogeology, volcanology, and related fields. The name 'cubic law' stems from the volumetric flow rate (effusion rate in volcanology) being a function of the aperture or opening of the hydrofracture in the 3rd power, that is, the cube of the aperture. This is the standard model in hydrogeology [17], and is also the most common model for dikes, inclined sheets, and sills in volcanotectonics and hydrofractures/mineral veins in geothermal fields [18].
Because the cubic law implies that the volumetric flow rate depends on the aperture in the third power any slight variations in the aperture or opening of a fluid-conducting fracture (Figures 27 and 28) can give rise to great variations in the volumetric flow rate. It follows that there is channelling (focusing) of the flow into those the parts of the hydrofracture where the aperture is larger than in the neighbouring parts.
Consider first a vertical hydrofracture/mineral vein in a rigid host rock (Figure 1). From the Navier-Stokes equations [43,44] it follows that the volumetric flow Q up through the fracture is given by: where u ∆ is the opening or aperture of the hydrofracture, and W its width perpendicular to the flow direction, that is, the strike-dimension of the fracture ( Figure 29). Thus, the cross-sectional area perpendicular to the vertical fluid transport direction is where the strike dimension is many times (commonly 10 2−3 times) the opening, so that W >> u ∆ . The other symbols in Equation (9) are the geothermal fluid dynamic viscosity f µ and density f ρ , the acceleration due to gravity g, and the excess-pressure gradient / e p z ∂ ∂ in the direction of the fluid flow.
As indicated above, a rigid model, while commonly used, is not very realistic for fluid transport in the Earth's crust: an elastic model is more realistic. Then the walls of the hydrofracture and its source are free to deform elastically due to fluid-pressure changes as the fluid is transported (here vertically) out of the source and (here) to the surface (Figure 29b). Since the host-rock density r ρ normally differs from the fluid density f ρ , a buoyancy term becomes added to the excess-pressure gradient. Thus, for a vertical hydrofracture in a host rock that behaves as elastic, and by analogy with Equation (8), then Equation (9) may be rewritten in the form: where all the symbols are as defined above.
Equations (9) and (10) apply to geothermal fluid transport up through a vertical fracture. In a network of hydrofractures/mineral veins (Figures 1, 2, 21 and 22) there are many non-vertical, that is, inclined and horizontal, hydrofractures (Figures 1, 2, 21 and 22b). Let the dip-dimension (the path of the geothermal fluid from its source) be denoted by L and its dip by α . Then, by analogy with Equation (9), for a rigid host rock the volumetric flow rate Q up the inclined hydrofracture is given by: Similarly, by analogy with Equation (10), for an inclined hydrofracture in a host rock that behaves as elastic the volumetric flow rate Q becomes: Equations (11) and (12) are more general than Equations (9) and (10). More specifically, Equations (9) and (10) are special cases of Equations (11) and (12) where α = 90°, that is, the hydrofracture is vertical so that sin α = 1. Similarly, for a horizontal hydrofracture (a sill-like fracture) we have α = 0° and sin α = 0, so that the first term in the brackets of Equations (11) and (12) drops out. For a hydrofracture in the horizontal xy-plane with width W measured along the y-axis and the length L measured along the x-axis, we may substitute x for L and obtain, from Equations (11) and (12), with sin α = 0, the volumetric flow rate through the horizontal hydrofracture as: Here the only 'force' driving the flow of geothermal fluid along the horizontal fracture is the excess pressure gradient / e p x ∂ ∂ . When a vertical hydrofracture supplies geothermal fluids to a horizontal hydrofracture, a water sill, then the excess pressure in the sill is equal to the overpressure (Equation (6)) that develops in the vertical hydrofracture on its path towards the contact with the water sill. Figure 28. Ideal mode I crack geometry used to model extensional rock fractures such as hydrofractures. Only half of the crack is shown. The letters x, y, and z denote the coordinate axes. Here the crack is in the y-z plane and opens up parallel with the direction of the x-axis. The total opening displacement is 2u = Δu. The crack is opened by a stress σ3 due to the fluid pressure p, which varies as a function of y, that is, σ3 = −p(y).

Figure 29.
For geothermal fluid transport up through rock facture we can model the host rock (and thus the fracture walls) either as rigid or elastic. Only haf the fracture is seen here. The vertical fracture direction parallel with Q is the dip-dimenion, while the fracture direction at the surface is (here) half the strike-dimension. (a) In a rigid model the vertical hydrofracture and its sill-like source maintain their size and shape when the fluid pressure changes. Q is the volumetric flow rate and pe is the excess pressure in the sill-like reservoir before rupture and hydrofracture formation (Equation (9)). (b) In an elastic model both the hydrofracture and its sill-like source change their shape and size when the fluid overpressure in the fracture and the fluid excess pressure in the source change (Equations (10)-(12)).

Transport of Geothermal Fluids along Fault Zones
Transport of geothermal fluids in the damage zone is primarily along fractures (Figures 15 and 16). It follows that the elliptical crack/conduit model (Section 5.1) and the parallel-plate model (Section 5.2) are generally applicable to fluid transport in the dam-age zone; this conclusion is supported by studies of palaeo-geothermal reservoirs in the damage zones of fault zones. Such zones show dense networks of mineral veins, indicating that much of the fluid transport when the geothermal reservoirs within the fault zone were active was through hydrofracture networks (Figures 2, 21 and 22a).
The way geothermal fluids are transported within the fault core is more complex (Figures 2, 15 and 16). When there is a fault slip in the core, there may be, for a while, flow of geothermal fluids in the core along the slip surface, that is, the fault. The temporary permeability can then increase by a factor of 10 6−9 because all the pores and fractures that the slip-surface meets become interconnected following the slip [17]. More specifically, the percolation threshold is reached along the slip surface, thereby giving rise to this enormous increase in temporary permeability.
For most of the time, even in an active fault zone, the fluid transport within the fault core is through porous media flow; this means that the flow of geothermal fluids in the core is governed by Darcy's law, which can be stated as follows. The velocity of fluid transport is proportional to the hydraulic gradient given by explained below. The law may be presented by the following equation [17] Here q is the specific discharge (or discharge velocity) in m s −1 , Q is the volumetric flow rate (flow volume per unit time) in m 3 s −1 , A is the cross-sectional area normal to the flow in m 2 , K is the coefficient of permeability in m s −1 , Δh is the difference in total head along the length of the fluid transport ΔL. Δh and ΔL have both the unit of length, so that the hydraulic gradient / h L ∆ ∆ (or / h L ∂ ∂ ) has no units. The minus sign is to indicate that q is in the direction of decreasing total head, but this sign is commonly omitted.
Although it was initially derived from experiments on water flowing thorough sand, Darcy's law is valid for most granular media so long as the flow is laminar, which is also known as viscous or streamline flow. In particular, Darcy's law is generally valid for sediments and for porous, solid rocks, such as scoria layers between lava flows and fault-zone breccia. It can also be used as a crude approximation for some densely fractured rocks, particularly if the fractures have various attitudes and form a well-interconnected fracture-pore space. Thus, geothermal fluid transport in parts of the highly fractured damage zone in some fault zones ( Figure 21) could possibly be modelled using Darcy's law, but, as indicated, fluid transport in the damage zone is generally more appropriately modelled using flow in fractured rocks whose basis is the cubic law (Section 5).
The physical meaning of Darcy's law is that the flow is driven towards regions of lower potential energy; only in the special case of horizontal flow is the transport necessarily from higher to lower pressure. Darcy's law can be presented in various forms. The letter i is commonly used for the hydraulic gradient / h L ∆ ∆ , in which case alternative expressions for Darcy's law include the following: q Ki = (15) and Q KiA = (16) where K is the hydraulic conductivity, i is the hydraulic gradient, and A is the cross-sectional area perpendicular to the flow, as defined above. A more general version of Darcy's law, for isotropic hydraulic conductivity K, may be given as: where q is the specific discharge (a vector), the minus sign is omitted, and grad h is the gradient of the scalar field h given by: where ∇ is del or nabla, and , and are unit vectors in the directions of the coordinate axes x, y, and z, respectively. The gradient ∇h represents the rate of change of h. The direction of ∇h coincides with that in which h changes fastest; namely for isotropic K, the direction of flow of the geothermal fluid.

Application
Here I present some applications of the theoretical framework established above as to fluid transport in geothermal fields. I first consider active geothermal fields, and then apply the general theory and field observations to palaeo-geothermal fields.

Fluid Transport in Active Geothermal Fields
There are many active high-temperature geothermal fields (and source reservoirs for these fields) in Iceland. Most of these are supplied with geothermal fluids through fractures. Many of the fractures are related to earthquake activity [15,16], and some to dikes ( Figure 1). If the geothermal water first collects in a reservoir where it builds up fluid pressure, then at the time of rupture of the reservoir and injection of a hydrofracture the excess pressure pe in the reservoir is normally the same as the tensile strength of the host rock, 0 T , at the rupture site. The general condition for rupture of the source is then [17]: Here, l p is the overburden pressure or lithostatic stress at the rupture site in the roof of the reservoir/source, l t e p p p − = is the excess pressure, that is, the difference between the total fluid pressure t p in the reservoir and the lithostatic stress at the time of reservoir rupture, and 3 σ is the minimum compressive (maximum tensile) principal stress.
If a hydrofracture becomes injected into the roof of the reservoir and begins to propagate up into the crustal layers above, its overpressure po becomes the one shown by Equation (6).
Consider now a hot spring supplied with water through a fracture. If there is accumulation of water in a reservoir which then ruptures (for example, during stress changes related to earthquakes in the area, Figure 21), then the excess pressure should be equal to the tensile strength. The in-situ tensile strength of rocks is generally from 0.5 MPa to 9 MPa, primarily based on hydraulic fracture measurements in drill holes, with the most common values being 2-4 MPa [18]. Thus, we can use the average value as = = 0 T p e 3 MPa. The surface temperature of many hot springs such as in the Geysir field in Iceland is from about 70 °C to 100 °C [47]. For Geysir itself the water temperature rises quickly with depth, and reaches about 120 °C at the depth of 20 m in the main pipewhich is still not enough to cause boiling, for which overheating is needed [15]. More specifically, in the pipe of Geysir the temperature is 90-95 °C from the surface to a depth of 10 m in the pipe, where it suddenly rises first to 118 °C at 12 m depth, and then to about 124 °C at 15 m depth and remains at that temperature to the depth of 22 m, which is the bottom of the pipe or 'conduit' [16]. Similar temperature profile is obtained for the other (and much more frequently) erupting geyser in the field, name Strokkur. The temperature of geothermal fluids increases with further depth; in fact, in so-called high-temperature fields in Iceland the temperature at 1 km depth is above 200 °C. Additionally, using geochemical measures, the reservoir temperatures of Geysir itself are estimated at about 240 °C [47]. In the low-temperature fields close to the Geysir geothermal field, however, the water temperature is mostly between 45 °C and 70 °C [47].
The sources can be at various depths. For example, studies of mineral veins in fault zones in the Bristol Channel Basin in the UK as well as in fault zones in Iceland both indicate that the depth of origin of the geothermal fluids in exposed mineral veins is of the order of several hundred metres to a maximum of about 1200 m [10,35]. In both cases the figures refer to the depth below the present exposures of the veins. For the mineral veins of the Husavik-Flatey Fault in North Iceland, the veins are currently at about 1.5 km below the original surface, suggesting that the depths of their sources were mostly between 2 and 3 km below the original surface at the time when the geothermal fields generating the veins were active. The exact depth of the exposures of the veins in the UK is not known, but appears similar (1-2 km) to that of the veins in Iceland [35,48].
For the present purpose, consider geothermal water at temperature of 100 °C. The dynamic viscosity may be taken as 2.8 × 10 −4 Pa s and the density as 958 kg m −3 [49]. Both values would depend on the composition of the water, but we may use this as approximate values for water at this temperature. Let us apply this to mineral veins in Iceland. If most of the veins originate at depths of 2 km or less, then the average density ρr of the crust in Iceland to that depth is about 2540 kg m −3 [50]. Additionally, the tensile strength of the source, thereby the excess pressure pe, may be taken as 3 MPa. Since most of the veins and the strike-slip fault hosting them are vertical, we use Equation (10) where factor B is about 5 × 10 6 m −1 s −1 . Here we have used the acceleration of gravity g as 9.8 m s −2 . Moreover, since the excess pressure gradient itself is negative (excess pressure decreases from 3 MPa in the source to 0.1 MPa at the Earth's surface, which we take as 0 MPa) the pressure gradient becomes added to the buoyancy term to drive the geothermal water towards the surface. While the aspect ratios ( where the factor C = 2 × 10 9 m −1 s −1 . We can now apply these results to yields in active geothermal fields. For example, a well at Helludalur just north of the Geysir geothermal field in South Iceland met with a reservoir at a depth of 386 m providing a yield of 5 L s −1 [47]. Using Equation (21) we can estimate what size of a hydrofracture would be needed to provide that flow into the well (and, through the well, up to the surface). From Equation (21) a yield of 5 L s −1 = 5 × 10 −3 m 3 s −1 gives Δu = 0.0012 m, that is, 1.2 mm. Here we have used the average aspect ratio, namely W/Δu = 400, in which case the length W corresponding to the aperture of 1.2 mm would be 480 mm or close to 0.5 m. There is a considerable range in the aspect ratios, from less than 100 to over 600 (excluding extreme outliers), so by a factor of at least 6 ( Figure 26). Thus, for the maximum length, the hydrofracture with the aperture of 1.2 mm could be as long as 0.7 m. This is similar to the common length of mineral veins in the palaeo-geothermal network from the Husavik-Flatey Fault in Iceland (Figure 26).
The main result here, however, is that a single small hydrofracture can theoretically supply the volume of geothermal water provided by a well. While it may be a single fracture that provides much or, or all, the geothermal water through the walls of the well at a given time, the transport of the water from the source, the reservoir or aquifer, is normally through a network of fractures ( Figure 21); this network transports geothermal fluids within the reservoir as well as from deeper levels in the crust.
The entire Geysir geothermal field, which contains many hot springs (only two of which are erupting, Strokkur, every 5-10 min, and the Great Geysir, very rarely nowadays), has an area of about 3 km 2 [16] and a total yield of about 12 × 10 −3 m 3 s −1 during the past decades [47]. Based on Equation (21), the entire field could theoretically be supplied with geothermal water through a single hydrofracture with an aperture of about 0.0016 m, or about 1.6 mm. Using the same aspect ratio as before, 400, the corresponding length or strike dimension of the fracture could be 640 mm or about 0.6 m. It may come as a surprise that so small increase in aperture and length of a hydrofracture is all that is needed to provide a yield that is greater by a factor of 2.5 than that of the well above. However, this follows from the cubic law for fractures in rocks behaving as elastic (Equations (10), (12) and (20)), which shows that a comparatively slight increase in aperture can greatly increase the volumetric flow rate Q through a hydrofracture.
Many active geothermal fields maintain their permeability through earthquake activity [6,10,14]. The same applies to the Geysir geothermal field, and the Great Geysir itself in particular [15]. In South Iceland, at the location of the Geysir geothermal field, there are primarily three sets of faults: NNE-trending, NE-trending, and ESE-trending [15]. All these fault sets are seen in the vicinity of the Geysir geothermal field, such as in the canyon of the river Hvita at the location of the waterfall Gullfoss (the Golden Waterfall). The NNE-strikes are primarily associated with dextral faults, the NE-strikes with normal faults, and the ENE-strikes with sinistral faults.
Following the two M6.6 earthquakes in the South Iceland Seismic Zone in the year 2000, there was an increased eruption activity in the Great Geysir [15]. The main earthquake faults did not reach so far north that they ruptured the Geysir geothermal field but the earthquakes triggered microseismicity primarily on NNE-trending faults that dissect the geothermal field and are also found south of the field [16,51]. Many of the hot springs in the Geysir geothermal field and in its surroundings align with lineaments of similar strikes [16]. As indicated, most of the NNE-striking lineaments are dextral strike-slip faults. In fact, detailed studies of the pipe supplying hot water to the Great Geysir and to Strokkur show that the pipes gradually change from essentially circular at the surface to elliptical fractures with NNE-strike at depths exceeding about 20 m [16]. This change is as expected and in perfect harmony with the hot springs being supplied with water through networks of hydrofractures associated with fault zones. Figure 21 shows part of such a network in a fault zone which, however, is orders of magnitude larger than the fault zones associated with the Geysir geothermal field.

Fluid Transport in Palaeo-Geothermal Fields
When the aspect (W/Δu) ratio is known, as well as the elastic properties of the rock at the time of hydrofracture emplacement, the fluid overpressure or driving pressure po can be estimated from standard fracture-mechanics equations [52]. For application to hydrofractures, the standard equation for fluid overpressure po at the time of hydrofracture formation is [17]: where Δu is the mode I opening displacement (palaeo-aperture, vein thickness), E is Young's modulus, ν is Poisson's ratio, and W is the strike-dimension or length of the fracture/vein in a direction perpendicular to the fluid transport (only half of W is supposed to be shown in Figure 29). This fracture model assumes that the fracture is a through-crack, that is, extends from one free surface to another one. The model also assumes that the strike-dimension W of the fracture is smaller than its dip-dimension. These assumptions may not always hold true, but are commonly good first approximations for hydrofractures in a dense network. This follows because many such fractures would meet contacts between flow units/lava flows in volcanic areas, and contacts between limestone (sandstone) and shale (mudrock) layers in sedimentary basins, where the contacts would be filled with fluids and thus act as free surfaces. Moreover, many of the hydrofractures meet other active hydrofractures, which then act as free surfaces (Figures 21 and 22a). The thicknesses of lava flows are generally greater than the strike-dimensions of mineral veins ( Figure 6) and the availability of cooling joints that may be partly used by the hydrofractures both would favour the dip-dimension being larger than the strike-dimension, in agreement with the assumptions in Equation (22).
The mineral vein in Figure 30 is located at about 1.5 km depth below the initial surface of the lava pile and belongs to the vein network ( Figure 21) of the damage zone ( Figure 16) of the Husavik-Flatey Fault. We can now use Equation (22) to estimate the fluid overpressure or driving pressure at the time of vein formation. The strike-dimension (outcrop length) of the vein is W = 80 cm or 0.8 m, and a maximum thickness is 2 mm or 0.002 m. At the depth of 1.5 km in the Icelandic lava pile, a typical Young's modulus is about 15 GPa and Poisson's ratio is 0.25 [50]. There is no evidence of multiple injections, so we assume that the entire vein was formed in a single fluid injection. For veins that grow in a self-similar manner, however, even if they form in several injections, the approach used here would still be valid.
This result is very similar to other overpressure estimates from sets of mineral veins. For example, as indicated above, the aspect ratio (W/Δu) for the veins in the Husavik-Flatey Fault is, on average, about 400 ( Figure 27). That is the same as the aspect ratio in Equation (23). Thus, assuming that the elastic constants E and ν did not change much during the formation of the vein network (Figure 21), the average overpressure associated with the formation of that network would be the same as above, namely 20 MPa.
In a very different environment, namely the Bristol Channel Basin in Southwest England, the aspect ratios of mineral (mostly calcite) veins resulted in an average geothermal fluid overpressure at the time of vein formation of 18 MPa [35]. While these re-sults from Iceland and England cannot be generalised, they indicate that overpressures of this order may be common in vein networks that supply fluids to geothermal fields while active.

Discussions
Understanding the transport of geothermal fluids is of increasing importance for several reasons. The mechanics of transport of geothermal fluids in the crust adds to our knowledge of the flow of crustal fluids in general; these fluids include groundwater, oil and gas, and magma, in addition to geothermal fluids. Although the interest in understanding the mechanics of flow of oil and gas is likely to decline in the coming years and decades in connection with the energy transition to renewables, the interest in understanding the flow and accumulation of groundwater, geothermal water, and magma are all going to increase. While the theoretical aspects of the flow of fluids in porous rocks are now well understood [53][54][55], those of flow in fractured rocks are still far from fully understood and the subject of intensive research efforts [17,56,57].
Fluid transport in natural geothermal fields is primarily along dikes and faults. Dikes are confined to areas of current or earlier volcanic activity whereas faults occur in geothermal fields both in volcanic areas as well as in non-volcanic sedimentary basins. In the paper it is demonstrated that flow of geothermal fluids (and groundwater as well) is very common along dikes and inclined sheets. The evidence from palaeo-geothermal fields is very clear on this matter: there are numerous mineral veins along the margins, and inside, dikes and inclined sheets (Figure 1). The central parts of dikes are the last ones to cool down and may encourage fluid transport, particularly if the central parts contain numerous vesicles that become interconnected, as is also common in lava flows (cf. Figure 21). The mineral veins commonly use cooling (columnar) joints as parts of their paths (Figure 12).
Perhaps the main paths for geothermal fluids along dikes and sheets are at their margins, that is, at the contacts between the dike rock and the host rock; these paths are commonly seen in the field (Figure 1). The reason for the contacts being common flow paths is primarily the stress concentration that occurs at the contacts and results in associated fracturing. Dikes and other intrusions act as elastic inclusions, namely solid bodies whose elastic properties differ from those of the host rock [58]. Elastic inclusions generally modify the stresses in their vicinity and, when the inclusion has a higher Young's modulus (stiffness) than its surroundings (as is normally the case for a basaltic intrusion), the stresses concentrate (are raised) at the inclusion. Concentration of stresses at contacts between the dike and the host rock favours fracture development, hence paths for geothermal fluids.
Dikes act as channels for geothermal fluids not only in high-temperature fields, but also in low-temperature fields. In Iceland the distinction between these is the temperature at 1 km depth. If the temperature is above 200 °C, it is classified as a high-temperature field, but if the temperature at this depth is below 150 °C it is classified as low-temperature fields; those with temperatures in the range of 150-200 °C are rare and classified as intermediate-temperature fields [59]. Most low-temperature fields in Iceland are associated with dikes. The dikes commonly strike at right angles to the dipping lava flows ( Figure 9) and tend to collect the water and channel it towards the surface along the fractures at the contact between the dike and the host rock.
The function of faults as geothermal fluid conduits is well established. Not only is this function clear from field observations of faults in palaeo-geothermal fields (Figures 2,  20 and 21), but is also confirmed by observations of active geothermal fields. The activity of individual hot springs is commonly enhanced following earthquakes. In particular, dormant geysers may resume their erupting activities during and following earthquakes [15,16]. In fact, many natural geothermal fields are associated with faults, both normal faults [6,7] and strike-slip faults [10,14].
Fault zones are generally composed of two main hydromechanical units, namely a core and a damage zone (Figures 15-17 and 19). In the absence of recent fault slip, geothermal fluid transport in the core is controlled by Darcy's law (Equations (14)- (17)). This law applies to flow in all particulate or porous materials. These are the materials that characterise the fault core; in particular, they are fault breccia, gouge, and clay. Thus so long as there is no renewed fault slip along the core, the flow in the core is porous-media flow and occurs through interconnected pores of various types.
Since large parts of the fault core may not reach the percolation threshold in the absence of renewed fault slip, the permeability along the core may be very low [36]. In particular, the transverse permeability may be so low as to allow the formation of 'compartments' in the geothermal reservoir. This means that the pore-fluid pressure on one side of a fault may differ significantly from that on the other side. Fault-related compartments are well known in hydrocarbon reservoirs and may occur in some geothermal reservoirs as well. The difference between these types of reservoirs as regards fault activity, however, is that most natural high-temperature geothermal reservoirs occur in tectonically active areas-they are mainly confined to the boundaries of the tectonic plate. This implies that many of the faults associated with geothermal reservoirs are active (but rarely so in hydrocarbon reservoirs), so that fault slip occurs from time to time.
Following fault slip the temporary permeability increase of a fault core may be of the order of 10 6 to 10 9 [17]. The primary reason for this increase is that the percolation threshold is reached along the entire fault-slip plane ( Figure 31). For a while, the geothermal fluid transport in the core is then no longer primarily controlled by Darcy's law but rather by the cubic law (Equations (9)-(13)). It follows that there will be an enormous increase in permeability along the slip plane. Gradually, however, during periods without further fault slips, healing and sealing (secondary mineralisation) of the fractures opened during the fault slip, as well as of the slip plane itself, will reduce the permeability and bring it back to porous-media dominated flow (Darcy's law). During subsequent fault slips, if they occur, the 'cycle' will repeat itself-porous media controlled flow, following by fracture-controlled flow, followed again by porous media controlled flow. Figure 31. Schematic illustration of the internal hydromechanical structure of a dextral strike-slip fault zone. When fault slip occurs, all the numerous small cracks and pores in the core that meet with the slip plane become interconnected along that plane, so that the percolation threshold is reached along that plane. This is the main reason why the permeability of the core may temporarily increase by many orders of magnitude following slip.
As for the damage zone, the fluid flow is normally fracture-controlled, as expressed by Equations (9)- (13). Since the fracture frequency in the damage zone increases towards the contact between the damage zone and the core [17], the most highly conductive parts of the damage zone are commonly close to the core or at the boundary between the core and the damage zone. In large fault zones, such as the Husavik-Flatey Fault in Iceland, very extensive networks of geothermal-fluid-conducting fractures may develop in the damage zone ( Figure 21). These can supply fluids to large and long-lasting geothermal fields. Understanding the way hydrofracture networks form and evolve over time, particularly in the damage zone of fault zones in general, is of very great importance for improving and refining our knowledge of natural geothermal systems.
In this paper the focus has been on analytical models in relation to fluid transport along dikes and faults. Analytical models are very useful; not only as regards actual calculations but also in bringing into focus the physical principles that control fluid transport. Numerical methods, however, have been much used to model fluid transport in the crust. For example, the numerical software Modflow has been available for several decades to model the transport of groundwater flow [60]. This software is mostly for porous media, or fracture networks that can be modelled as analogous to porous media, and therefore, as regards flow along dikes and fault zone, primarily applicable to fluid transport along the fault core. There have been many numerical studies on fluid transport and pressure in fault zones [14,[61][62][63][64][65][66][67]. Some of these studies focus on the interaction between fluid pressure and transport, in the fault core and damage zone, and fault slip [63]. Others are more concerned with fluid-pressure development in fault zones in relation to the earthquake cycle [62]. Still others model the effects of fault intersections in channelling geothermal fluids to the surface to form hot springs [64,66].
There are also numerous numerical models on the propagating paths of hydrofractures driven by various types of fluids, such as magma, water, and geothermal fluids [14,17,18]. Many hydrofractures form their paths through rupturing layered and faulted rocks, while others follow fault planes and dikes. Studies of these kind include how hydrofractures become arrested [68][69][70], the effects of fractures in general [71] and faults in particular [72,73] on hydrofracture paths, as well as general energy considerations for the way fluid-driven fractures select their paths in layered and faulted rocks [18,71]. Numerical models for understanding fluid transport in the crust, and in particular that of geothermal fluids, are very useful and complementary to the basic analytical models which are emphasised the present paper.

Summary and Conclusions
The focus of this paper is on the transport of geothermal fluids in geothermal systems. The emphasis is on fluid transport in natural geothermal systems, but many of the results apply as well to Enhanced Geothermal Systems (EGS), particularly those that use fault zones for fluid transport. The main conclusions of the paper may be summarised as follows:

•
Geothermal fluids in natural systems flow primarily through two types of large-scale structures: dikes and faults. Faults are the universal type of conduits of fluids and occur in all geothermal systems, including sedimentary basins and volcanic fields. Dikes (and inclined sheets) occur in active (and inactive) volcanic areas and in addition to acting as conduits for geothermal fluids also function as their heat sources. • Some geothermal fluids are conducted in the central parts of dikes, but the main conduits are normally at the dike margins, that is, at the dike contacts with the host rock. This is due to stress concentration at the margins that gives rise to fracture formation and provides potential paths for fluids. The stress concentration is because dikes act as elastic inclusions that modify and mainly raise the stresses at their contacts with the host rock.

•
Fluid transport along dikes (and inclined sheets) at their margins is controlled by the cubic law. This law states that the volumetric flow rate depends on the aperture of the fracture in the power of 3. This means that slight changes in the aperture may result in very large changes in the volumetric flow rate. Fluid transport to the surface of many geothermal fields is largely controlled by flow along dikes. In particular, in low-temperature fields many dikes strike at right angles to the dip direction of the lava pile (and the flow direction of the associated geothermal fluids) and tend to collect geothermal fluids and transport them to the surface as warm springs. • Dikes (and inclined sheets) also act as heat sources for hundreds or thousands of years following their emplacement. The rate of cooling depends on the thickness of the dike in the power of 2 so that the thick dikes act as heat sources for a much longer time than the thin dikes. Drilling into high-temperature fields in active vol-canoes, such as in Iceland, shows that the frequency of dikes and inclined sheets increases rapidly with depth and the dike/sheet intensity reaches 80-100% of the rock at depths of 1.5-2 km. Moreover, feeder-dikes have been shown to have increased the geothermal activity in fields such as in the active Hengill Volcano in Iceland.

•
Fault zones are the other main type of structures conducting geothermal fluids. The connection between fault zones and the transport of geothermal fluids is well established, both in palaeo-geothermal fields as well as in active geothermal fields. Fault zones are composed of two main hydromechanical units, namely, the fault core and the fault damage zone. The core is comparatively thin and made of particulate materials, primarily breccia, gouge, and clay and related low-permeability materials. By contrast, the fault damage zone is characterised by fractures whose frequency normally increases (irregularly) from the contact between the damage zone and the host rock to the contact between the damage zone and the core.

•
Geothermal fluid transport in the core during non-slipping periods is controlled by Darcy's law for fluid flow in porous media. During fault slip in the core, and for a while following slip, the transport of geothermal fluids may by mostly controlled by the cubic law for flow in fractures. Furthermore, fault slip results in an enormous number of pores and fractures becoming interconnected along the slip plane, meaning that the percolation threshold is reached along the plane. Gradually, however, healing and sealing of the slip plane and the associated fractures occurs so as to reduce the permeability. Thus, in the absence of renewed fault slip, the permeability diminishes and the flow in the core again becomes primarily controlled by Darcy's law.

•
Geothermal fluid transport in the damage zone is normally controlled by the cubic law; this follows because the damage zone is characterised by fractures. However, the fracture frequency generally varies much within the damage zone, resulting in the local Young's modulus (stiffness) and thus the stress field controlling fracture formation and slip also varying inside the zone. It follows that the local stress field inside the damage zone (and the core) may be very different from the regional stress field outside the fault zone. It is the local stress field that primarily controls fracture formation and reactivation (including fault slip) and, thereby, fluid transport in the damage zone.

•
Healing and sealing of fractures and pores in the fault core and the fault damage zone are mainly due to secondary mineralisation associated with the transport of geothermal fluids thorough these zones. The mineralisation gradually decreases the permeability and porosity through the formation of mineral veins and amygdales. Thus, if no further fracturing or fault slip takes place the permeability, and thus the activity of the hot springs, mud pots, and/or steam vents of the associated geothermal fields, decreases. It follows that earthquake activity is needed to maintain the permeability in active natural geothermal fields.
Funding: This research received no external funding.

Data Availability Statement:
The available data is mostly from the cited publications.