Alice Namias

Mathieu Sicard

Mylène Weill
email: mylene.weill@umontpellier.fr

Sylvain Charlat
email: sylvain.charlat@univ-lyon1.fr

From Wolbachia genomics to phenotype: molecular models of cytoplasmic incompatibility must account for the multiplicity of compatibility types

Keywords: Wolbachia, toxin-antidote, symbiosis, Culex, Drosophila

Wolbachia endosymbionts commonly induce cytoplasmic incompatibility, making infected males' sperm lethal to the embryos unless these are rescued by the same bacterium, inherited from their mother. Causal genes were recently identified but two families of mechanistic models are still opposed. In the toxin-antidote model, interaction between the toxin and the antidote is required for rescuing the embryos. In host modification models, a host factor is misregulated in sperm and rescue occurs through compensation or withdrawal of this modification. While these models have been thoroughly discussed, the multiplicity of compatibility types, i.e., the existence of many mutually incompatible strains, as seen in Culex mosquitoes, has not received sufficient attention. To explain such a fact, host modification models must posit that the same embryonic defects can be induced and rescued through a large variety of host targets. Conversely, the toxin-antidote model simply accommodates this pattern in a lock-key fashion, through variations in the toxin-antidote interaction sites.

Introduction

Cytoplasmic incompatibility (CI) denotes a phenomenon of conditional sterility induced by various endosymbiotic bacteria, of which Wolbachia is the most renowned. In its simplest form (Figure 1), CI occurs when males carrying Wolbachia mate with uninfected females. In such crosses, fertilization takes place normally but most or all embryos die before hatching [START_REF] Marshall | Some notes regarding the morphological and biological differentiation of Culex pipiens Linnaeus and Culex molestus Forskal (Diptera, Culicidae)[END_REF][START_REF] Yen | New hypothesis of the cause of cytoplasmic incompatibility in Culex pipiens L[END_REF][START_REF] Yen | The etiological agent of cytoplasmic incompatibility in Culex pipiens[END_REF]. In contrast, development proceeds if the female is infected, regardless of the male's infection status. This means that Wolbachia makes infected females more fertile than infected ones on average and on the contrary renders infected males less fertile than uninfected ones. However, the male side does not matter as far as Wolbachia is concerned since only females transmit the infection to their offspring, through the egg cytoplasm. The infected lineage frequency thus tends to increase, and more effectively so if the infection has no negative side effect and is perfectly transmitted across generations [START_REF] Turelli | Dynamics of cytoplasmic incompatibility and mtDNA variation in natural Drosophila simulans populations[END_REF]. In this perspective, CI is understood as a selfish adaptive feature of the symbiont, increasing its chances of invading new species following rare events of horizontal transmission [START_REF] Sanaei | Wolbachia host shifts: routes, mechanisms, constraints and evolutionary consequences[END_REF]. Cytological studies indicate that CI induces similar cellular phenotypes in a wide range of hosts [START_REF] Ryan | Post-fertilization effect of incompatibility factors in Mormoniella[END_REF][START_REF] Lassy | Cytological analysis of fertilization and early embryonic development in incompatible crosses of Drosophila simulans[END_REF][START_REF] Callaini | Mitotic defects associated with cytoplasmic incompatibility in Drosophila simulans[END_REF][START_REF] Tram | Role of delayed nuclear envelope breakdown and mitosis in Wolbachia-induced cytoplasmic incompatibility[END_REF][START_REF] Landmann | Wolbachia-mediated cytoplasmic incompatibility is associated with impaired histone deposition in the male pronucleus[END_REF][START_REF] Bonneau | The cellular phenotype of cytoplasmic incompatibility in Culex pipiens in the light of cidB diversity[END_REF]. During the first embryonic division, paternal chromosomes exhibit defects in condensation and improperly segregate during anaphase. The resulting embryos are aneuploid or haploid and thus inviable, except in haplodiploid species where they may survive as males [START_REF] Bordenstein | Host genotype determines cytoplasmic incompatibility type in the haplodiploid genus Nasonia[END_REF]. While these developmental failures were described a while ago, the underlying genetics have remained elusive until recent years. Because embryonic viability is impeded by the sperm of infected males and rescued by the presence of Wolbachia in the eggs, it has long been formalized that CI is a two-sided phenomenon, implying some kind of "modification" in mature sperm and some kind of "rescue" taking place in the eggs [START_REF] Hurst | The evolution of cytoplasmic incompatibility or when spite can be successful[END_REF][START_REF] Werren | Endless forms: Species and Speciaton[END_REF]. This modification/rescue or mod/resc model is generally acknowledged as a flexible framework that could accommodate any underlying molecular mechanism (box 1).

Two causal genes were recently identified concomitantly in the Wolbachia genomes of major CI model systems: the wPip strain from Culex pipiens mosquitoes and the wMel strain from Drosophila melanogaster flies [15,16] (although other genes may also be involved, e.g. Beckmann, eLife 2019, Scholz Nat Com 2021). In both study systems, these genes, generally referred to as "CI factors" (cif) form an operon-like structure in the prophage WO region [16]. Transgenic expression of the upstream gene (cifA) is necessary and sufficient to prevent CI onset in fly embryos, that is, to recapitulate the resc function. Transgenic expression of cifB is toxic in yeasts and this effect is rescued by its coexpression with cifA, making CifB a good "mod factor" candidate. Yet only when cifA and cifB are coexpressed in male flies is CI induction effectively recapitulated. Intriguingly, transgenic expression of cifB from another strain (wRec, naturally present in Drosophila recens) is sufficient to induce strong CI-like embryonic mortality, but this effect is not rescuable, casting doubt on the hypothesis that it induces the same molecular defects as the dual cifA/B expression [17]. Putting aside this particular case, transgenic expression of the CifA protein thus seems paradoxically required for both CI rescue and induction. This may indicate that CifA displays two distinct and somewhat opposite molecular activities depending on the context, making it a mod co-factor in maturing sperm and a resc factor in embryos, as proposed in the so-called "two-by-one model" [18]. Alternatively CifA may just act as a rescue factor, required not only for the normal development of embryos but also for preventing the inherent toxicity of CifB during spermatogenesis [19]. A recent study in Anopheles gambiae mosquitoes appears to fit this view: transgenic expression in males of a single cifB gene from the wPip was sufficient to induce CI [20]. Besides, rather than being a necessary mod factor in this system, high expression of cifA in males attenuates the CI penetrance. Within the cif genes family, at least two active operons have been identified, that can be specifically referred to as cid and cin, based on the enzymatic activity, deubiquitinase or nuclease, of their downstream protein.

Yet, comparative genomics studies led to the description of additional homologs distributed into five clades [16,[21][22][23]. Thus far, CI induction has been experimentally associated with cid and cin operons, respectively falling into clades I and IV, and also with one operon from group II, but putative functional domains have been identified in all five groups [15,17,18,24]. Now that the key CI effectors have been identified, understanding their molecular mechanism will require deciphering what these proteins actually do. The scope of possibilities remains very large, making it relevant to evaluate the various plausible mechanistic models in light of the current data. Several recent review papers thoroughly discussed these models [19, [START_REF] Shropshire | Models and Nomenclature for Cytoplasmic Incompatibility: Caution over Premature Conclusions -A Response to Beckmann et al[END_REF][START_REF] Beckmann | Caution Does Not Preclude Predictive and Testable Models of Cytoplasmic Incompatibility: A Reply to Shropshire et al[END_REF][START_REF] Chen | The Biochemistry of Cytoplasmic Incompatibility Caused by Endosymbiotic Bacteria[END_REF][START_REF] Shropshire | Symbiont-mediated cytoplasmic incompatibility: What have we learned in 50 years?[END_REF] but perhaps did not put enough emphasis on one specific question that may help guiding future research: based on current knowledge, should we expect the CI rescue to stem from a direct interaction between CifA and CifB in the embryos? Answering this question will essentially distinguish two families of mechanistic models that are still opposed, namely host-modification (HM) versus toxin-antidote (TA) models (Figure 2) [17,19,[START_REF] Shropshire | Symbiont-mediated cytoplasmic incompatibility: What have we learned in 50 years?[END_REF][START_REF] Poinsot | On the mechanism of Wolbachia-induced cytoplasmic incompatibility: Confronting the models with the facts[END_REF].

As recently summarized (e.g. [17,19]), under TA models, the mod and resc factors act like a lock and key: the toxin (mod) affecting the paternal chromosomes is transported through sperm into the embryo where its interaction with the antidote (resc) is needed to prevent the modification and thereby let the first mitosis proceed (Figure 2A). On the contrary, HM models propose that the mod/resc interaction is only indirect and happens through host effectors: the mod factor(s) induces paternal DNA defects but is not transferred into the embryo so that rescue happens by reversing its effects (Figure 2B). In fact, in one particular HM version, called the "mistiming model", the mod and resc sides of CI stem from a single process: slowing down the mitotic dynamics; under this model, embryos are only viable if the paternal and maternal pronuclei are synchronous. A second, perhaps more abstract HM version, called the "goalkeeper model", assumes modification takes place in males in a strain-specific quantity, needing to be precisely remedied to rescue embryos. Finally, the "titration-restitution" model posits that Wolbachia alters the concentration of a crucial host effector in sperm and thus induces CI unless this is counter-balanced in the embryo. inhibited by an antidote produced by the maternal Wolbachia. The paternal antidote may also be required during sperm development but is presumably less stable than the toxin and would thus be degraded prior to fertilization. b) In Host Modification (HM) models, a host factor is modified in the sperm. This modification causes developmental defects unless CifA in the egg inhibits/reverses its effects. In both frameworks, the maternal CifA proteins would constitute the rescue factor, while paternal CifB (with a debated contribution of paternal CifA) would be responsible for CI induction.

Which of these two families of models is most likely correct? Binding assays provide some elements of answer by showing that the CidA and CidB proteins, as well as CinA and CinB, tightly bind in vitro, while CidA only binds to few other elements ( [START_REF] Beckmann | The Wolbachia cytoplasmic incompatibility enzyme CidB targets nuclear import and protaminehistone exchange factors[END_REF], resp. [24]). While not ruling out HM models, this pattern indicates that a central requirement of TA models is likely fulfilled although it still awaits in vivo confirmation. Notably, another central requirement of the TA model, namely, the production of the CidB protein in sperm and its conveyance into the eggs has not been documented yet, neither in Culex nor in Drosophila. However, such absence of evidence should not be taken as evidence of absence, because this protein has not been specifically searched for in these tissues. In any case, we now turn to consider a different line of arguments, building on the multiplicity of Wolbachia compatibility types.

Compatibility types: a simple case study

We have so far described CI in its simplest form, that occurring between infected males and uninfected females, a pattern often referred to as unidirectional CI, because the reverse cross produces viable progeny (Figure 1). Yet it has long been known that CI may also occur between males and females that are both infected, but by distinct Wolbachia strains (Fig. 3). The Drosophila simulans system provides a rather simple situation to apprehend this phenomenon (reviewed in [START_REF] Merçot | Wolbachia infections in Drosophila melanogaster and D. simulans: Polymorphism and levels of cytoplasmic incompatibility[END_REF]). This species naturally hosts several Wolbachia lineages, usually called "strains", easily distinguished on the basis of standard molecular markers such as 16S rRNA. Some of these strains happen to have lost the ability to induce CI [START_REF] Hoffmann | Naturally-occuring Wolbachia infection in Drosophila simulans that does not cause cytoplasmic incompatibility[END_REF][START_REF] Rousset | Evolution of single and double Wolbachia symbioses during speciation in the Drosophila simulans complex[END_REF] but three of them are still capable of both CI induction and rescue: wRi, wHa, and wNo, originally described in lines from Riverside (California), Hawaii and Noumea (New Caledonia), respectively. These three strains display quantitative differences with regard to CI penetrance but most importantly, they are all bidirectionally incompatible. In that sense, they are said to harbor distinct "compatibility types". To what extent can the TA and HM frameworks account for such mutual incompatibilities (Fig. 3)? Within the TA model, the multiplicity of compatibility types is readily explained by the possibility that distinct mod/resc pairs could differ, not in their targets, but in their interaction sites, as the Lock-Key analogy makes obvious (Fig. 3A). Within the HM family, the mistiming model, while elegant in proposing that CI induction and rescue may stem from a single phenomenon, just fails to explain the existence of more than one compatibility type, as previously noted [START_REF] Poinsot | On the mechanism of Wolbachia-induced cytoplasmic incompatibility: Confronting the models with the facts[END_REF]. In contrast, other HM models may account for mutual incompatibilities, but this would imply that the mod factor(s) of incompatible Wolbachia strains achieve the same cytological defects through the activation of different host effectors (Figure 3B). Furthermore, the resc factors should similarly vary in their targets to specifically inhibit the consequences of their mod counterpart (Figure 3B). The above-mentioned existence of more than one functional family within cifB genes (cidB, cinB and possibly other yet uncharacterized enzymes) could form the basis of such a diversity. Transgenic expressions in yeast provide indirect evidence that the cid and cin operons from wPip may be mutually incompatible: CidA and CinA can only rescue the toxicity respectively caused by CidB and CinB [15]. Yet the same approach suggested that compatibility types do not simply come down to enzymatic families: while the toxicity associated with cidB from wHa was fully rescued by its cognate cidA, expression of the cidA gene from wPip only lead to imperfect rescue [START_REF] Beckmann | The Wolbachia cytoplasmic incompatibility enzyme CidB targets nuclear import and protaminehistone exchange factors[END_REF] Coming back to the aforementioned Wolbachia strains from D. simulans, it is noticeable that the wRi genome carries both a type I (cid) and a type II operon, while wHa harbors only a type I operon and wNo only a type III operon [16,21]. Assuming compatibility is preserved within functional families, embryos carrying the wRi strain should be immune to the CI induced by wHa. Yet they are not, which further highlights that incompatibility patterns cannot be readily reduced to enzymatic families.

In brief, the TA model can rather simply explain bidirectional incompatibilities through changes in the CifA-CifB interaction sites. Notably, this model does not rule out the possibility that different Wolbachia strains may induce the same cellular defects through different effectors (as suggested by the fact that both nuclease and deubiquitinase activities produce similar cytological phenotypes) but it does not rely on this enzymatic diversity to explain mutual incompatibilities (Figure 3A). The HM framework also remains compatible with the existence of more than one compatibility type. Mutual incompatibilities would then stem from differences in the host targets of the distinct mod factor(s), appropriately matched by their resc counterparts (Figure 3B). The existence of more than one functional family offers at best a partial explanation for the observed compatibility patterns. While still plausible, the HM framework then faces the challenge of explaining that a variety of distinct pathways would lead to the same cytological phenotype. As we shall now discuss, a further challenge is to explain that compatibility types are substantially more diverse than suggested by the Drosophila case.

Compatibility types: a bigger picture

While the D. simulans case is sufficient to demonstrate that different Wolbachia strains may be mutually incompatible, the scale of this phenomenon is better apprehended by considering data from Culex pipiens. Even before Wolbachia was identified as the causal agent, this species was known for displaying complex patterns of incompatibilities among distinct geographic populations [START_REF] Laven | Eradication of Culex pipiens fatigans through Cytoplasmic Incompatibility[END_REF][START_REF] Duron | Tracking factors modulating cytoplasmic incompatibilities in the mosquito Culex pipiens[END_REF][START_REF] Atyame | Multiple Wolbachia determinants control the evolution of cytoplasmic incompatibilities in Culex pipiens mosquito populations[END_REF]. Following the observation that, in contrast with D. simulans, C. pipiens only hosts one Wolbachia "strain" (named wPip) as far as standard molecular markers such as 16S rRNA can tell, it has been hypothesized that host genetic variation likely contributed to this complexity [START_REF] Rousset | Cytoplasmic incompatibilities in the mosquito Culex pipiens: How to explain a cytotype polymorphism?[END_REF]. Yet, more discriminant molecular markers soon revealed that the wPip clade is in fact composed of several lineages [START_REF] Duron | Transposable element polymorphism of Wolbachia in the mosquito Culex pipiens: evidence of genetic diversity, superinfection and recombination[END_REF] and more importantly, experiments controlling for the host genetic background repeatedly indicated that Wolbachia alone is responsible for the observed incompatibilities [START_REF] Duron | Tracking factors modulating cytoplasmic incompatibilities in the mosquito Culex pipiens[END_REF][START_REF] Atyame | Multiple Wolbachia determinants control the evolution of cytoplasmic incompatibilities in Culex pipiens mosquito populations[END_REF][START_REF] Sicard | Cytoplasmic incompatibility variations in relation with Wolbachia cid genes divergence in Culex pipiens[END_REF].

Once acknowledged, this feature provides important insights with regard to Wolbachia compatibility types. First, they are numerous: the number of mutually incompatible strains in C. pipiens exceeds by far those reported in D. simulans. For instance, by reciprocally crossing 19 different C. pipiens lines, 15 distinct compatibility profiles have been revealed [START_REF] Duron | Tracking factors modulating cytoplasmic incompatibilities in the mosquito Culex pipiens[END_REF][START_REF] Nor | On the genetic architecture of cytoplasmic incompatibility: Inference from phenotypic data[END_REF]. Second, they can diverge very rapidly, as indicated by the study of intra-population CI variations [START_REF] Duron | Rapid evolution of Wolbachia incompatibility types[END_REF][42][43][44]. Third, they can be asymmetrical: strain A may rescue the CI induced by strain B without the reciprocal being true. Fourth, they can be non-transitive: strains A and B may be mutually compatible and yet display distinct compatibility patterns with strain C [START_REF] Atyame | Multiple Wolbachia determinants control the evolution of cytoplasmic incompatibilities in Culex pipiens mosquito populations[END_REF]45]. On the basis of these properties, it was possible to predict not only that the CI genes should display an important diversity in C. pipiens, but also that a single Wolbachia genome should carry more than one mod/resc pair [START_REF] Atyame | Multiple Wolbachia determinants control the evolution of cytoplasmic incompatibilities in Culex pipiens mosquito populations[END_REF][START_REF] Nor | On the genetic architecture of cytoplasmic incompatibility: Inference from phenotypic data[END_REF].

The recently uncovered diversity of the cif genes in wPip matched these predictions rather neatly [43]. Both cin and cid genes are found in wPip genomes but all strains carry these two enzymatic types of CI effectors, ruling out the hypothesis that incompatibilities could result from carrying either one or the other. Most critically, while the cin genes are monomorphic in all wPip genomes studied so far, the cid genes are amplified and diversified, with up to six different copies of a given cid gene within a single genome [43]. Furthermore, cid diversity strongly correlates with CI patterns [43,44]. Notably, mutual incompatibilities have nothing to do with differences in the Cid deubiquitinase functional domain, which is monomorphic, but rather with polymorphism in CidA/CidB interaction sites [43,44], Hochstrasser 2021.

It may be obvious already to the reader that these latter results better fit a TA than an HM framework. The large diversity of compatibility types constitutes a serious difficulty for the HM models, because it unreasonably multiplies the number of distinct host targets, both on the male and female sides, through which Wolbachia strains would produce the same cellular phenotype. The goalkeeper model [46] may in theory offer a partial solution to this problem in suggesting that within an HM framework, the high number of targets required could be attenuated by quantitative differences in cif genes expressions. Yet, such quantitative variations would fail to explain the non-transitive CI relationships observed in Culex pipiens, unless one assumes host genetic variation is also involved, an hypothesis that, once again, was repeatedly ruled out in Culex pipiens, e.g. lastly in [START_REF] Sicard | Cytoplasmic incompatibility variations in relation with Wolbachia cid genes divergence in Culex pipiens[END_REF]. Furthermore, if changing the gene expression is indeed presumably "easier than evolving a new functional protein" [46], the existence of multiple polymorphic copies of cid genes is demonstrated in Culex [43]. Much more parsimoniously, within the TA framework, mutual incompatibilities can simply result from differences in CidA/CidB interaction sites. Intermediate affinities among Cid proteins and quantitative variations in their production could also readily explain the intermediate levels of rescue recently described [START_REF] Sicard | Cytoplasmic incompatibility variations in relation with Wolbachia cid genes divergence in Culex pipiens[END_REF].

Conclusion

It is arguably premature to conclude that the HM versus TA debate is over. A better understanding of what the Cif proteins actually do will be required to reach that point.

However, unless one assumes that the Culex case represents an oddity relying on peculiar mechanisms, the multiplicity of Wolbachia compatibility types weighs heavy in the balance and has perhaps received too limited attention in previous discussions of CI mechanistic models [17,19,[START_REF] Shropshire | Symbiont-mediated cytoplasmic incompatibility: What have we learned in 50 years?[END_REF]. Among other approaches, further studies focusing on the binding affinities between different versions of the wPip Cid proteins along with a deeper assessment of cid genes diversity in this system will certainly contribute to deciphering the mechanisms and evolution of CI.

Box 1: A quick guide to the CI models jungle

Mod/resc, two-by-one, lock/key, toxin/antidote, host-modification, mistiming, goalkeeper, titration-restitution… Many terms are used in the current literature to describe molecular models of CI, from which a certain confusion may arise to the newcomer. To complement the main text, we provide here brief definitions of these various terms.

• Mod/resc stands for "modification / rescue". The mod/resc model is currently not considered as a specific molecular model, but rather a summary formulation of the general idea that CI has both a paternal and maternal side: it is induced by some "modification" occurring in infected males, and it is "rescued" by something occurring in the egg. All the following models can fit in this general framework.

• Similarly, the two-by-one model should not be considered as a specific molecular model, but rather a synthetic reminder that, at least in the Drosophila transgenic system, expressing cifA alone is sufficient to rescue the eggs, while both cifA and cifB are required in males to recapitulate the induction of CI. Debates are taking place with regard to the generality of this pattern, or its explanation, which may or not involve cifA being considered as a dual factor, involved not only in the rescue of CI, but also its induction.

• The toxin-antidote (TA) and host-modification (HM) terms are currently used to describe and distinguish two main families of explanations.

• In TA models, a toxin is deposited in sperm and transmitted to the egg upon fertilization, leading to embryonic defaults unless an antidote is there, transmitted from the maternal Wolbachia, that binds the toxin, which somehow prevents its effetcs. Compatibility patterns among strains then result from the (in)ability of male-borne toxins and femaleborne antidotes to bind. In that sense, the TA model, in the current usage of this term, is strictly equivalent to a lock-key model, although in principle the toxin-antidote wording could be seen as more general.

• In HM models, CI is due to some modifications of host factors, ultimately affecting the paternal DNA and leading to embryonic defaults unless they are either compensated or reversed in the egg by the maternal Wolbachia. HM models thus do not posit any direct interaction between the cif proteins. (In)compatibiliy between strains may then arise through qualitative and/or quantitative differences in the targeted host factors. Several more specific models fit within the HM family. The mistiming model proposes that CI stems from asynchrony in the timing of events taking place in the paternal and maternal pronuclei after fertilization. The titration-restitution model posits that some host factor is removed (titrated) in sperm and conversely made more abundant (restituted) in the egg. Finally, the goalkeeper model is slightly more abstract. Its name metaphorically refers to the idea that, just as a goalkeeper moves both horizontally and vertically, CI induction may imply a qualitative and a quantitative dimension: not only the hostmodification but also its degree would matter.

Figure 1 :

 1 Figure 1: CI in its simplest form: compatibility relationships between infected and

Figure 2 :

 2 Figure 2: Schematic views of the two main families of CI mechanistic models. In both

Figure 3 :

 3 Figure 3: A simple case of mutual incompatibilities explained under the TA and HM

  . The existence of incompatibilities within functional families is further suggested by pull down assays indicating that CidB from wPip only interacts with its cognate CidA, and not with either CinA from wPip or CidA from the wMel operon [15]. Recent transgenesis experiments in Drosophila provide mixed results regarding the degree of specificity associated with different Cif factors from distinct strains carrying clade I and/or clade II operons, presumably differing in their enzymatic activities [17]. The CI induced by two distinct clade I operons was interchangeably rescued by their CifA factors, suggesting compatibility would be preserved within the deubiquitinase functional family, at least between the two Wolbachia strains tested. Yet, surprisingly, a clade I CifA protein was found able to rescue the CI induced by a clade II operon, although the reciprocal was not true.

  induces cytoplasmic incompatibility. Nat Microbiol 2017, 2:1-7. 16. LePage DP, Metcalf JA, Bordenstein SR, On J, Perlmutter JI, Shropshire JD, Layton EM, Funkhouser-Jones LJ, Beckmann JF, Bordenstein SR: Prophage WO genes recapitulate and enhance Wolbachia-induced cytoplasmic incompatibility. Nature 2017, 543:243-247. 17. Shropshire JD, Rosenberg R, Bordenstein SR: The impacts of cytoplasmic incompatibility factor (cifA and cifB) genetic variation on phenotypes. Genetics 2021, 217. 18. Shropshire JD, Bordenstein SR: Two-By-One model of cytoplasmic incompatibility: Synthetic recapitulation by transgenic expression of cifA and cifB in Drosophila. PLoS Genet 2019, 15:1-20. ** This transgenesis study in Drosophila melanogaster suggests that wMel cifB alone does not suffice to recapitulate cytoplasmic incompatibility, which requires expression of both cifA and cifB. This result forms the basis of "Two-by-One" pattern in wMelinfected D. melanogaster. 19. Beckmann JF, Bonneau M, Chen H, Hochstrasser M, Poinsot D, Merçot H, Weill M, Sicard M, Charlat S: The Toxin-Antidote Model of Cytoplasmic Incompatibility: Genetics and Evolutionary Implications. Trends Genet 2019, 35:175-185. 20. Adams KL, Abernathy DG, Willett BC, Selland EK: Wolbachia cifB induces cytoplasmic incompatibility in the malaria mosquito. bioRxiv 2021, ** One of the very first studies on an artificial Wolbachia infection in anopheline species. It investigates whether cifA and cifB from the wPip strain are able to induce CI in Anopheles gambiae, and shows that cifB expression alone in males is sufficient to induce CI. Furthermore, high cifA expression in males appears to alleviate the CI effects. These results are supporting a typical toxin-antidote model and do not fit the two-by-one pattern. 21. Lindsey ARI, Rice DW, Bordenstein SR, Brooks AW, Bordenstein SR, Newton ILG: Evolutionary Genetics of Cytoplasmic Incompatibility Genes cifA and cifB in Prophage WO of Wolbachia. Genome Biol Evol 2018, 10:434-451. 22. Martinez J, Klasson L, Welch JJ, Jiggins FM: Life and Death of Selfish Genes: Comparative Genomics Reveals the Dynamic Evolution of Cytoplasmic Incompatibility. Mol Biol Evol 2021, 38:2-15. * This paper, following up the work from Lindsey et al. (GBE, 2018) is a wide study of cif genes homologs as found in many species. It provides insights into their evolution and shows notably (i) that protein domains required for CI are widely conserved, (ii) that the modification function is usually lost before the rescue function, and (iii) that cif genes often recombine, although recombination is mostly restricted to close homologs. 23. Bing XL, Zhao DS, Sun JT, Zhang KJ, Hong XY, Sloan D: Genomic Analysis of Wolbachia from Laodelphax striatellus (Delphacidae, Hemiptera) Reveals Insights into Its "jekyll and Hyde" Mode of Infection Pattern. Genome Biol Evol 2020, 12:3818-3831. 24. Chen H, Ronau JA, Beckmann JF, Hochstrasser M: A Wolbachia nuclease and its binding partner provide a distinct mechanism for cytoplasmic incompatibility. Proc Natl Acad Sci U S A 2019, 116:22314-22321. ** This study provides the first evidence that the cin operon from wPip can induce CI. After validating its nuclease activity, it shows that transgenesis of the wPip cin operon in D. melanogaster induces male sterility and produces typical CI defects. The cinB gene narrow contact zone in mosquito field populations. Mol Ecol 2015, 24:508-521. 43. Bonneau M, Atyame CM, Beji M, Justy F, Cohen-Gonsaud M, Sicard M, Weill M: Culex pipiens crossing type diversity is governed by an amplified and polymorphic operon of Wolbachia. Nat Commun 2018, 9:319. ** This study demonstrates the amplification and diversification of the cid operons in wPip (with up to 6 different copies of cid genes found within a single Wolbachia genome), and shows some correlation between cid gentic variations and CI phenotypes, thus supporting the role of these genes in the unparalleled diversity of CI patterns in C. pipiens. 44. Bonneau M, Caputo B, Ligier A, Caparros R, Unal S, Perriat-Sanguinet M, Arnoldi D, Sicard M, Weill M: Variation in Wolbachia cidB gene, but not cidA, is associated with cytoplasmic incompatibility mod phenotype diversity in Culex pipiens. Mol Ecol 2019, 28:4725-4736. 45. Atyame CM, Labbé P, Dumas E, Milesi P, Charlat S, Fort P, Weill M: Wolbachia divergence and the evolution of cytoplasmic incompatibility in Culex pipiens. PLoS One 2014, 9:e87336. 46. Bossan B, Koehncke A, Hammerstein P: A new model and method for understanding Wolbachia-induced cytoplasmic incompatibility. PLoS One 2011, 6.

Acknowledgments

This work was funded by the French ANR (project "CIAWOL" ANR-16-CE02-0006-01).


Franck Plouraboue

Pascal Kurowski

Jean-Marc Boffa

Jean-Pierre Hulin
email: hulin@ariane.fast.u-psud.fr

Stephane Roux

Experimental study of the transport properties of rough self-affine fractures

Keywords: Fracture, Self-affine, Roughness, Flow, Permeability, Conductivity )

An experimental study of the transport properties of fluid-saturated joints composed of two complementary rough fracture surfaces, translated with respect to each other and brought in contact, is reported. Quantitative roughness measurements on different fractured granite samples show that the surfaces have a self-affine geometry from which the dependence of the mean aperture on the relative displacement of fracture surfaces kept in contact can be predicted. Variations of the hydraulic and electrical conductances of the joint are measured as functions of its mean aperture. A simple parallel plane model accounts for the global trend of the measurements, but significant deviations are observed when a relative lateral displacement of the surfaces is introduced. A theoretical analysis of their origin shows that they are due both to the randomness of the aperture field and to a nonzero local slope of the surface near the injection hole; the corresponding conductivity fluctuation amplitudes have power law and linear variations with the lateral displacement, and are enhanced by the radial injection geometry.

Introduction

Interest and geometrical characteristics of natural fractures

Fractures are very common features of earth's upper crust and range from large-scale faults to microcracks. Fluid flow through fractures is a very important mechanism in porous media with small intrinsic permeability; it greatly influences transport processes Ž both in natural crystalline and tight sedimentary rocks granite, basalt, low-permeability . Ž . sandstone, etc. and in artificial structures concrete, cement . In addition, fluid transport in fractured rocks is relevant to a number of geophysical applications, such as hazardous Ž waste storage, hydrogeology, and geothermal and petroleum engineering [START_REF] Evans | Technical report[END_REF][START_REF] Petitjean | Modelling of fracture propagation during overbalanced perforating[END_REF] Understanding transport properties of natural frac-™ Ž . Ž . tures, particularly their relation to variations of the local aperture a x cf. Fig. 1 , is, thus, crucial for numerous applications.

The objective of the present work is to experimentally analyze the dependence of these transport properties on fracture roughness. An important step in modelling the topography of fracture surfaces has been the observation that these surfaces are Ž . well-described as self-affine structures [START_REF] Bouchaud | Models of fractal crack[END_REF] in various natural and Ž . Ž . man-made materials, such as granite [START_REF] Brown | Correction to AA note on the description of surface roughness using fractal dimensionB[END_REF], marble Poon et al., 1992 , Ž . Ž . sandstone [START_REF] Plouraboue | Experimental study of the influence of ´Ž . crack velocity on fracture surface roughness of rocks[END_REF]concrete Schmittbuhl et al., 1995 , as well as in Ž .

Ž natural faults [START_REF] Power | Roughness of natural fault surface[END_REF] . Such surfaces satisfy a scale invariance Feder, 1988;

. [START_REF] Falconer | Fractal Geometry: Mathematical Foundations and Applications[END_REF]

Ž

. Ž . roughness exponent or Hurst exponent [START_REF] Feder | Fractals[END_REF] . Experimental z values are close Ž to 0.8 " 0.05 for many materials and fracture mechanisms lower values of the exponent have, however, been obtained at small scales or for materials fractured at a finite fraction of the acoustic Rayleigh wave speed, or yet, for some loosely consolidated . granular materials . These geometrical properties imply the existence of long-range Ž . infinite correlations in the topography of the fracture surface; these result in scale effects that must be considered when extrapolating laboratory-scale measurements to field applications. The question, then, is: how can scaling laws relating length scales parallel and perpendicular to the fracture plane be translated into scaling laws for Ž transport properties, such as the hydraulic and electrical conductances of fractures [START_REF] Roux | Some physical properties of self-affine rough[END_REF] For such applications, it is essential to take into account the correlation Ž . between the two faces of the crack [START_REF] Brown | Correlation between the surfaces of natural joint rocks[END_REF] , and simple phenomenological models Ž have been developed to introduce such geometrical correlations [START_REF] Brown | Simple mathematical model of a rough fracture[END_REF] . [START_REF] Mourzenko | Geometry of simulated fractures[END_REF] The present paper considers a crack bounded by two complementary self-affine Ž . surfaces the two halves of the same fracture block translated with respect to each other Ž . Ž . distance between each block is thus h u s h y h , which is also the mean aperture. The fracture height z x 1 2

Ž . is represented as well as the local aperture a x .

Ž .

Ž . and brought into contact Fig. 1d , or were left separated Fig. 1c . A model describing the local thickness distribution in such a configuration had been previously suggested Ž . [START_REF] Plouraboue | Aperture of rough cracks[END_REF] Hereafter, the in-plane component of the relative displacement ís referred to as the AshiftB u. When both surfaces are in contact, the normal component Ž . h is a function of u, i.e. h ' h u for a given direction of the displacement. The local Ž . Ž . Ž . aperture a x, y is then related to h u and to the height profile z x, y of the fracture surface by a x , y s z x q u, y y z x , y q h u . 2

Ž . Ž . Ž . Ž . Ž .
Ž . This assumes that the two surfaces are perfectly mated for a zero lateral shift u s 0 , Ž . which amounts to neglect the effect of surface wear. The results from Eq. 2 show that Ž . Ž . ² Ž .: h u coincides with the spatial average of the aperture: h u s a x, y . In a previous Ž . Ž . work [START_REF] Plouraboue | Aperture of rough cracks[END_REF]

Ž .

1995 that h u also depends on the length L of the crack and varies as

X z ( h u s u K 1 q K log Lru 3 Ž . Ž . Ž .
where urL is the relative amplitude of the shift of the joint, z is the roughness exponent of the crack surface, and K and K X are two constants. This result has been verified both experimentally on a fractured granite block and through numerical simulations on Ž . computer-generated self-affine surfaces [START_REF] Plouraboue | Aperture of rough cracks[END_REF] This scaling law is ímportant for extrapolating experimental laboratory results to the field scale.

Models of flow in fractures

Many authors have studied the influence of variations of the local fracture aperture on Ž . fluid transport properties experimentally [START_REF] Witherspoon | Validity of cubic law for fluid flow in a deformable rock fracture[END_REF][START_REF] Brown | Transport of fluid and electric current through a single fracture[END_REF] , Ž . Ž theoretically [START_REF] Ge | A governing equation for fluid flow in rough fractures[END_REF]numerically Tsang et al., 1988;[START_REF] Rasmuson | Chemical transport in a fissure rock: verification of a numerical model[END_REF][START_REF] Brown | Fluid flow through rock joints: the effect of surface roughness[END_REF][START_REF] Gutfraind | Study of hydrodynamic permeability of fractures using gas automata[END_REF]Mourzenko et al., 1995;[START_REF] Zimmerman | Lubrification theory analysis of the permeability of Ž . rough-walled fractures[END_REF] . A lowest-order approach to model flow properties is to represent the fracture as two parallel surfaces separated by a distance h. Then, the fluid flow rate across a section of the fracture is proportional to the applied hydraulic pressure gradient and to the cube of the separation distance h. Such a simple law is, however, of limited Ž . use [START_REF] Evans | Technical report[END_REF] , particularly when h is of the order of magnitude of the surface roughness.

™ ™

V s y = P . 4

Ž .

m Ž The Reynolds approximation is valid when the aperture gradient and thus, the local . < < slope of the surface is much smaller than unity, i.e. = a < 1. For self-affine cracks, the averaged slope estimated at a scale l varies as l zy1 , l y0 .2 , and therefore increases at smaller length scales. However, experimental profile measurements on all tested samples have shown that local slopes are small above a scale l of 50 mm. Hence, the Reynolds approximation should be valid for fractures with mean apertures above a few tens of microns. Experimental data will be analyzed in Section 5 in the context of the Reynolds approximation.

v

In the inertial regime, even with small local slopes, the Reynolds approximation fails Ž because nonlinear terms of Navier-Stokes equation are not negligible anymore [START_REF] Gutfraind | Study of hydrodynamic permeability of fractures using gas automata[END_REF][START_REF] Brown | Simple mathematical model of a rough fracture[END_REF] Such effects are not considered in the present work since all measurements are carried out in the range of linear variation of the total flux with the applied pressure.

Although a large number of theoretical and numerical studies have been devoted to the transport properties of rough cracks, few have experimentally addressed the effect of Ž roughness in natural fractured rocks [START_REF] Witherspoon | Validity of cubic law for fluid flow in a deformable rock fracture[END_REF][START_REF] Brown | Transport of fluid and electric current through a single fracture[END_REF] . [START_REF] Brown | Transport of fluid and electric current through a single fracture[END_REF] . The present work first confirms and experimentally extends previous geometrical results obtained on self-affine rough fractures to the different samples for which transport properties will be analyzed. After describing the experimental set-up and the characterization of fracture surface roughness, measurements of the influence of the sample size on the geometrical aperture of different granite joints are discussed. More Ž . Ž . precisely, fluctuations standard deviation s of the mean aperture h u for two blocks h ™ in contact are studied for a varying orientation of a relative displacement u of constant amplitude in the mean plane of the fracture. This experimental study is complemented by numerical simulations of geometrical properties of the joint. Then, the transport properties of rough joints saturated with a fluid are investigated experimentally. Flow measurements are realised in the domain of linear variation of the pressure gradient with the flux. The variation of both electrical and hydraulic conductances is analyzed as a Ž . function of the normal, h u , and lateral, u, relative displacements between joint faces. Large fluctuations of the electrical conductance, only partly correlated to aperture variations, are observed and studied theoretically.

Experimental procedure

Rock samples

Experiments were performed on fractured granite and basalt blocks of initial size 25 = 25 = 15 cm 3 . Two parallel notches were first carved in the middle of two opposite Ž . faces, and a mode I crack [START_REF] Scholz | The Mechanics of Earthquakes and Faulting[END_REF] was propagated from one notch to the other. Two types of rocks have been used to test the sensitivity of the results to the nature of the material. The first two blocks are made of grey granite with coarse millimetric-size Ž . heterogeneities sample nos. 1 and 3, originating from Brittany, France . The basalt Ž . block is black and finely grained sample no. 2, originating from Africa . Both materials are very strongly consolidated so that fractures propagate through the heterogeneities. In both cases, the fracture roughness displays a very broad range of characteristic length Ž . scales from 100 mm to a few cm , as expected for a self-affine surface.

Profilometry measurements

Fracture surfaces are characterized using a mechanical profilometer of 2 mm resolution perpendicular to the mean fracture plane and 25 mm parallel to it. Height profiles are measured on the fracture surface along parallel straight segments, and the repeatabil- Ž . Dashed lines correspond to the slope of the estimated rough exponent z of Table 1. a Superposition of two Ž . power spectra obtained from measurements carried out on two perpendicular directions on sample no. 2. b Superposition of two power spectra obtained from measurements carried out on samples nos. 1 and 2. ity on the height values is of order 10 mm. The typical length of the segments is 100 mm and the distance between neighboring profiles is 400 mm. One of the halves of the fractured block is put horizontally on the sensing assembly with the fracture surface oriented upward. The local height is measured by lowering a fine-tipped sensor onto the surface and recording the position at which contact is detected; the cycle is then repeated for all points over the profile with a horizontal spacing of 50 mm. The roughness Ž . exponent z can be determined, for instance, by computing the Fourier spectrum P k of Ž .

Ž . the height z x of the surface along a profile. For self-affine profiles, P k can be Ž shown to obey Brown, 1988;[START_REF] Poon | Surface measurement and fractal characterisation of naturally fractured rocks[END_REF][START_REF] Plouraboue | Experimental study of the influence of ´Ž . crack velocity on fracture surface roughness of rocks[END_REF]Schmittbuhl et . al., 1995

< < 2 y1y2 z P k s z k A k . 5
Ž . Ž . Ž . Fig.

2a displays two average Fourier power spectra obtained from block no. 2 using profiles measured in two perpendicular directions. Each spectrum corresponds to the average of ten 2048-point parallel profiles. The near coincidence of the spectra suggests an isotropy of the self-affine properties of the surface. Fig. 2b displays superimposed Ž . average spectra corresponding to a granite and a basalt sample block nos. 1 and 2 .

Ž . Roughness exponents z are estimated from the slopes of Fig. 2a and b using Eq. 5 . Both exponents are very close to the 0.8 value reported by other authors on granite and Ž . other types of materials. The standard deviation s L of the height of the sample z surface of a profile segment of length L is finally expected to scale according to

s L s ML z 6 Ž . Ž .
z

The parameter M characterizes the amplitude of the roughness and is lower for the finer grained basalt sample than for the coarser granite ones. Values obtained by Ž . computing s L over segments of length L s 10 cm are listed in Table 1.

Transport properties measurement set-up

Transport properties are analyzed on a second set-up allowing both hydraulic Ž . permeability and electrical conductance measurements Fig. 3 . The lower block can be Ž . moved vertically with an amplitude of 12 mm and in both horizontal directions with an amplitude of "25 mm by micrometric translation devices. The upper block is allowed to move only upwards from a rest position adjusted initially with both halves of the fractured block in contact. Sideway motions are blocked by ballpoint devices pressed against the edges. Downward motion from the initial position is prevented by adjusting the propping system with the two blocks in contact. Displacements of both blocks in all directions are permanently monitored with a resolution of "5 mm by mechanical sensors. After each set of measurements, one verifies that the blocks come back in contact over the whole fracture surface when the lower block is brought back to its initial position. This is checked by raising the lower block after it has touched the upper one and by verifying that upward displacements measured by three sensors at different points on the upper block are identical. In situations where the blocks have a relative horizontal displacement, they are considered to be in contact when any of the sensors located on the upper block detects an upward motion larger than 10 mm.

The sample is then plunged inside a parallelepipedic basin filled with water so that the liquid surface is always above the uppermost point of the fracture. In this way, it always remains saturated with water once trapped air has been sucked away. For a given relative position of the blocks, the hydraulic conductance is measured by injecting fluid at a constant flow rate through a vertical 5-mm borehole drilled at the center of the upperblock. Injected fluid flows out on the sides. The injected flow field is only approximated to be radial because of the square outer shape of the blocks. This does not have much influence on the transport properties since most of the pressure drop is localized near the injection. The effective hydraulic conductance of the fracture joint is determined by measuring the pressure difference D P between the injection point and the rim of the fracture at five different flow rates spanning one order of magnitude. Finally, a linear regression is performed on each data set. A 40%:60% glycerol-water mixture is used as the flowing fluid in order to increase the measured pressure and Ž . reduce the Reynolds number resulting in more linear pressure vs. flow rate relations . D P is measured through a sensing tube opening close to the fracture surface in the injection borehole. This allows the elimination of the effect of the pressure drops in the injection circuits. The other end of the tube is connected to one of the ports of a Validyne DP15 differential pressure sensor. The other port is connected to another measurement point located at the rim of the fracture. The range of the sensor is "10 4 Pa and its precision is 0.1%. Flow rates range between 4 = 10 y8 and 2.5 = 10 y6 m 3 rs, and can be adjusted by steps of 2 = 10 y9 m 3 rs. The electrical admittance is determined with an HP 4192A computer-controlled bridge after saturating the fracture with a 1 grl NaNO -water solution. The measurement frequency is 20 kHz, which minimizes the 3 reactive component of the impedance and allows a precision of 0.01% to be obtained. Measurements are performed between a cylindrical movable electrode inserted into the injection borehole and a set of brass plates located in the fluid bath outside the blocks. The center electrode is lowered down to the lower block surface in order to extend vertically across the full fracture aperture.

Geometrical measurements

The geometrical characterization of the samples by profilometry has been comple- 2 for sample nos. 2 and 3. The relative deviation s u rh u clearly increases h < < with the lateral displacement u . This variation is a relevant additional information on the structure of the fracture joint. A numerical investigation of this problem is reported in Section 5.1.

Ž .

2 Ž . 2 z Finally, the scaling law of Eq. 3 is checked by plotting the ratio h u ru as a Ž . function of log urL . The values of z used are those obtained from mechanical Ž . Ž . profilometry measurements Table 1 . In Fig. 6a-b, h u is the average of the four Ž ) . curves displayed in Fig. 4a-b. Data points symbol corresponding to an average over all u values in Fig. 5 have also been plotted in Fig. 6b, and coincide well with the to the smallest scales. This vertical shift represents a residual aperture, and will also appear in transport measurements discussed in Section 4. Small secondary fractures branching away from the surface at a small angle have been observed and may account, in part, for the residual aperture. Another possible origin is the removal of some asperities during fracture propagation. It is also worth noting that the plots of Fig. 6a-b 2 Ž . 2 z are very sensitive to the choice of z , and that the variation of h u ru is no longer X Ž . linear if z varies by 15% or so. Values of the coefficients K and K in Eq. 3 deduced 3. K is a useful characterization of the roughness amplitude. Experiments on two very different types of samples reported in this section therefore confirm previous experimental results and theoretical predictions [START_REF] Plouraboue | Aperture of rough cracks[END_REF] as Eq. 3 relating the horizontal and normal displace-ḿents of two complementary self-affine surfaces in contact. Significant relative varia-™ tions of the aperture as a function of the orientation of the displacement vector u have, however, been observed, particularly at large amplitudes. In the next part, the depen-™ dence of transport parameters on the amplitude and orientation of u is investigated on the same set of samples.

Transport measurements

Electrical conductance measurements

The measurement technique is calibrated by moving vertically apart the two halves of the fractured block from the position of optimal contact; the blocks are kept parallel with a zero relative horizontal shift so that their local spacing is constant over the fracture surface. Then, the effective electrical admittance should be similar to that of two parallel planes with a distance equal to the relative normal displacement h. One may only expect a second order effect due to the roughness that increases the effective distance between Ž ) . the center electrode and the outside ones. The dotted lines in Fig. 7a-b symbols represent the variation of the electrical admittance with the vertical relative displacement h; this variation is nearly linear and its slope can be estimated by modelling the fracture as a set of two parallel disks of radius R and separation h with an inside electrode of Ž . radius r . The corresponding electrical admittance A h is which, for simplicity, will be referred to in the following sections as the electrical conductance. Actually, the complex admittance also includes reactive components of Ž . various origins double layer capacitance on electrodes, for example . The measurement frequency is chosen so that the capacitive component only represents a small fraction Ž . Ž . 0.1% of the overall admittance and can be neglected. In Eq. 7 , the shape of the outer boundary of the domain is approximated by a circle instead of the actual square boundary. However, the conductance only weakly depends on the actual shape of the Ž . outer boundary due to the logarithmic dependence on R in Eq. 7 . R is, in the following equations, taken to be equal to half the diagonal of the boundary of the granite samples.

Using s values measured on the bulk solution, one obtains theoretical values for 0 Ž . the slope 2ps rlog Rrr of the dotted curves of Fig. 7b, about 5-7% higher than the 0 0 experimental ones. This small difference may be due to the roughness of the surface, increasing the effective distance between the inner and outer electrode. Extrapolating the experimental linear variation to zero admittance gives a residual aperture estimate of 120 Ž . mm, which confirms geometrical measurements Fig. 7b . displacement are compared in Fig. 7a-b to the previous results obtained for u s 0. In all cases, and particularly for h -4 mm, the admittance A roughly increases linearly with Ž . h u ; this confirms that the distance h is the adequate scaling variable to describe the Ž Ž .. electrical conductance, as suggested by the parallel plane model Eq. 7 . Strong Ž . relative deviations from the linear law up to 20-30% are, however, observed, Ž . Ž . particularly for large values of h u and in the x q direction . These variations were investigated in more detail by measuring the electrical admittance A in the same Ž < < geometrical configuration as in Fig. 5 constant amplitude u and variable orientation of ™ . the displacement u . Fig. 8a displays the variations of A with the angle u defined as in Fig. 5. The relative amplitude of the variations of A are of orders 14%, 13.5% and 25% for u s 1, 3 and 10 mm, respectively. However, a significant fraction of these variations Ž . Ž may be due to the fact that h u varies significantly with the angle u as shown in Fig.

.

.

In order to separate out this effect, the variations of A corresponding to the data < < points of Figs. 5 and8a, and to the three displacement amplitudes u s 1, 3 and 10 mm, Ž . are plotted in Fig. 8b as functions of h u . As expected, the mean trend of the variation is clearly linear, but significant variations of up to "15% are still observed, and their Ž . amplitude increases markedly with h u .

Ž . In all cases, the normal displacement h u is, therefore, the adequate scaling variable controlling the overall variation of the admittance A. The amplitude of the deviations of Ž . A from this global linear dependence on h u is, however, quite important and tends to increase with h. As will be seen in the next part, similar effects can be observed experimentally for the hydraulic conductance.

Hydraulic conductance measurements

As for the admittance, hydraulic conductance measurements are first calibrated in a Ž . < < constant aperture geometry Fig. 1b , with a zero lateral relative shift u of the blocks. Practically, the variation of the pressure difference D P is measured as a function of the flow rate Q. For a parallel plane geometry and a radial injection, the ratio QrD P should verify

Q p h 3 1 s s K h 8 Ž . Ž . hy D P 6h log Rrr h Ž . 0 Ž
. in which h is the dynamic viscosity of the fluid kept constant for these experiments , r 0 the radius of the central hole for injection, and K the hydraulic conductance. The hy dotted line in Fig. 9 represents, for sample no. 2, the experimental variations of Ž . 1r3 QrD P with the relative displacement h of the blocks normal to the fracture plane Ž < < . the lateral displacement u is zero . A nearly linear variation is observed as expected Ž . from Eq. 8 . By applying this equation and including the measured fluid viscosity, one Ž . 1r3 obtains for each value of QrD P an effective hydraulic aperture h equal to within hy 8% to the mechanical displacement. As in the case of the admittance, h reaches for hy h s 0 a finite limit h ; it is the difference h y h that is actually compared to h. Hydraulic conductance variations are measured by moving the lower block horizon-Ž . tally over a distance u and putting it in contact again with the upper one Fig. 1c andd . As shown in Fig. 9, the overall variation of the hydraulic conductance with the vertical ™ Ž . displacement h u is roughly cubic for all four orientations of the displacement u investigated. Observed relative deviations from this mean law are of order 50%. These variations are not random fluctuations around the mean cubic law, but are systematic ™ deviations which depends on the orientation of the horizontal shift u. These deviations are expected to be enhanced by the radial geometry; the hydraulic conductance K is hy indeed strongly influenced by the structure of the flow field in the immediate vicinity of the injection hole, where pressure gradients are largest. Let us note that, within Reynolds approximation, K is much more sensitive to aperture variations than the electrical hy ™ y1 y3 Ž ² Ž . : . conductance, due to its cubic dependence on the local aperture K A a x .

hy

In this section, it has been shown that the parallel plane model provides a reasonable first-order description of the transport properties of natural fracture joints. Electrical and hydraulic conductance values in the radial geometry display globally linear and cubic variations with the mean aperture. On the other hand, lateral displacements of the fracture walls -at constant normal spacing h -induce significant variations in the electrical and hydraulic conductance. These variations strongly depend on the orientation ™ of the displacement u. Their origin, and particularly their dependence on the radial geometry, will now be discussed theoretically.

Theoretical analysis of experimental data

In this section, spatial mean aperture fluctuations of a self-affine rough joint are first investigated with the help of numerical simulations. Then, the influence of the orienta-™ tion of the displacement u on deviations from the parallel plane model found for the electrical and hydraulic conductances is discussed.

Numerical analysis of aperture fluctuations

Ž

. Numerical simulations were previously used [START_REF] Plouraboue | Aperture of rough cracks[END_REF] to analyze áperture variations with the relative displacement of rigidly translated self-affine surfaces kept in contact. Random two-dimensional self-affine surfaces are first generated Ž . numerically; they are defined by their height z x, y , and have the same statistical properties as the experimental fracture surfaces. A dichotomy algorithm introduced by Ž . Voss 1985 is used to generate self-affine surfaces with a roughness exponent close to Ž . the experimentally observed value z s 0.8 , and with amplitudes similar to those of Table 1. Using normalization, such that the small-scale roughness amplitude remains Ž . constant, the standard deviation s L of the height of the surface increases with the z Ž . sample size L following Eq. 6 . The simulations reported here are purely geometrical Ž . no fluid flow is computed ; they analyze statistical fluctuations of the spatial mean Ž . aperture h u from one realisation to another and, more specifically, deviations from the Ž . Ž . predictions of Eq. 3 . The generated self-affine crack surfaces z x, y are first ™ translated by a vector u, as displayed in Fig. 1. This allows the generation of the Ž . aperture field a x, y , and to determine numerically the normal relative displacement Ž .

² Ž .: qualitatively in experiments. The ratio between s u and h u , which represents the h relative fluctuation of the aperture, has also been plotted and is roughly constant with u; this implies that both the aperture fluctuations and its mean value have the same Ž ² Ž .:. dependence on the displacement i.e. s A h u . This result is consistent with h Ž . experimental geometrical measurements Table 2 , but the small number of available data points does not allow quantitative tests of the relation. Large fluctuations are also observed on the values of the electrical admittance in the same configuration of the < < Ž . surfaces for a given lateral displacement amplitude u Fig. 8a . As discussed above, a large fraction of these fluctuations is directly associated to the parameter s through the h linear relation between the normal displacement h and the electrical admittance A in the parallel plane model. The corresponding component s of the fluctuations of the A h ² Ž .: admittance will, thus, verify s A s or s A h u . However, there still remains a

A h h A h
significant scatter of the admittance for a fixed aperture once the admittanceraperture correlation has been taken into account, as can be seen in Fig. 8b. The corresponding fluctuation amplitude s also increases with u. In the following, the origin of such Af observed deviations from the parallel plane model is analyzed theoretically.

DeÕiations from parallel plane model

These deviations are studied here within the Reynolds approximation and in a radial flow geometry. Three possible origins of the fluctuations s will now be examined, Af namely, the influences of the finite sample size, of the local variability of the fracture aperture, and of a small angle of the fracture surface with respect to the horizontal plane. The first effect -purely geometrical -corresponds to the decrease of the length ™ of the interval between the fracture walls from R to R y u for a lateral translation u.

Ž . This effect would lead, from Eq. 7 , to a systematic correction independent of the ™ orientation of u for the value of the electrical admittance. This correction can be Ž . estimated through a first order Taylor expansion of Eq. 7 , with respect to the parameter < < Ž Ž .. urL, which leads to a term of order u r Rlog Rrr . Numerical estimates of this 0 correction using experimental parameters lead to a maximum deviation of 2.5% from linearity for the maximum relative displacement urL s 0.1. This effect is one order of magnitude smaller than that observed in experiments, and can be neglected here for the electrical as well as for the hydraulic conductances.

v

The local aperture variability may also modify the transport parameters. Previous Ž . analytical work Dewitt, 1995 has analyzed the influence of permeability heterogeneities on the macroscopic permeability for parallel flows; the approach of this author Ž . is based on a perturbation expansion [START_REF] Matheron | Elements Pour Une Theorie Des Milieux Poreux[END_REF] with respect to the relative Ž . amplitude of permeability fluctuations assuming a weak disorder . The heterogeneous Darcy equation is expanded into a homogeneous one, assuming a small ratio between ² : ² : the standard deviation s of the permeability and its mean value k , i.e. s r k -

1 k k
Žthe characteristic small parameter chosen is often the amplitude of log-permeability . fluctuations; [START_REF] Matheron | Elements Pour Une Theorie Des Milieux Poreux[END_REF] This same approach can be applied to radial flows Ž . [START_REF] Schwydler | Sur la precision du debit des puits[END_REF] expansions are also easily adapted to relative fluctuations Ž . of the aperture. The aperture field a x, y can be written as a sum of the spatial mean u Ž .

Ž

. aperture h u and of local fluctuations related to the height z x, y of the individual surfaces ² : a x , y s a x , y q d z s h u q z x q u, y y z x , y . 9

Ž . Ž . Ž . Ž . Ž . Ž . u u
The aperture field is assumed to be random but statistically homogeneous. The expected value of the aperture is independent of the point considered, and its Ž . Ž . standard mean deviation s a may be related to the translation u by s a s ( Mr K 1 q K log Lru . This correction is much smaller than experimental devia-Ž .

™ < < tions. Moreover, it depends only on the modulus u of the translation vector and not on its orientation, in contrast with the experimental results. where S is the joint surface area, and f is the average tilt angle of the surface with Ž . respect to the horizontal displacement plane. Eq. 10 predicts that the measured admittance will vary linearly with the tilt angle f; quantitatively, for f s 48, the estimated relative variation is of order 20%. Although actual tilt angles are somewhat smaller, the effect is strongly amplified by the radial flow geometry. In this case, the Ž . current increases close to the injection point as the inverse ; 1rr of the distance r.

Ž . The effective slope f in Eq. 10 then corresponds to a weighted average of local slopes with a particularly large influence of zones close to the injection point. Since the local slope in this region rarely coincides with the global plane of the fracture, f may take values differing significantly from zero. In a more precise analysis, one must introduce the radius r of the inner electrode as a new characteristic length scale. Two different 0 types of behaviours are then expected, depending on the relative values of u and r . Ž .

ii u 4 r In this case, the edge effect near the electrode is negligible and the main 0 contribution is due to a dipolar field

™ ™™ ™ K y K u 2 u r r Ž . el ref ™ ™ A y z r d r . 13 Ž . Ž . H 2 4 ž / K r r ref ™
5 5 This integral is computed over the domain rr -R. In order to estimate this 0 effect, a two-dimensional height map of nearly 1 cm 2 close to the injection point has been analyzed using a mechanical profilometer to estimate the order of magnitude of Ž . Ž Ž . expected corrections in the first regime i the second regime ii would require a much . Ž . larger map and will not be examined here . Fig. 11a shows the measured height z r in 0 a r s 3 mm radius circle. Corrections due to the radial geometry were estimated both 0 from a direct calculation on the two-dimensional map, and by computation of the tilt Ž . parameter f from Eq. 12 . These two methods give the same order of magnitude for the corrections. provides an equally good coincidence and supports the observation that D KrK is Ž . independent of u in regime i .

To conclude, the deviation of the electrical conductance from the parallel plane behaviour observed in a radial geometry can be largely accounted for by a finite tilt angle of the fracture surface near the injection point. This result can be easily generalized to the hydraulic conductance, which leads to a relative correction term about Ž three times larger than for the electrical conductance the factor 3 results from the cubic . law dependence of the hydraulic conductance on the aperture . Radial flow geometries are frequently encountered in configurations of practical interest, such as flow in fractures intersecting a borehole in oil reservoir or geothermal fields. The above analysis underlines the well-known fact that the structure of the rock near the borehole influences very much the flow properties in such geometries.

Conclusion

In the present paper, previous theoretical suggestions of scaling laws relating the electrical and hydraulic conductances of rough granite and basalt fracture joints, the relative normal and lateral displacements of the fracture surfaces, the total crack size, and the roughness exponent have been evaluated experimentally.

For purely geometrical aspects, finite-size effects on the relation between the lateral u Ž . and normal h u displacements of complementary fracture surfaces have been confirmed with these self-affine surfaces kept in contact.

Concerning transport properties, it has first been observed that a simple parallel plane model gives a good first order approximation of the variation of electrical and hydraulic conductivities with the relative normal displacement h of crack surfaces. A more detailed analysis has confirmed that, as expected, the electrical conductance varies linearly with the normal spacing h for complementary surfaces with no relative lateral Ž . displacement this corresponds to a constant local distance between surfaces . On the contrary, when the two surfaces are displaced laterally, large relative deviations of few tens of percents from the linear dependence of the conductance on h are measured. These deviations increase with the lateral displacement and depend strongly on its direction. Such features are observed both for the hydraulic and the electrical conductances.

These observations have been analyzed in the frame of the Reynolds approximation using a first order perturbation theory. Our analysis shows that these fluctuations are associated with the radial flow geometry used in our experiments. The transport properties are strongly influenced by the geometry of the flow paths near the injection, and weakly perturbed by the flow field heterogeneity. In particular, we have shown that a small tilt angle between the global surface plane and the local orientation of the surface in the vicinity of the injection point may strongly influence the transport properties. Significant anisotropies of the conductance variations with respect to the direction of the relative displacement are also expected. Such a hypothesis can be verified by comparing the results obtained in the present work to experiments realised in a parallel flow geometry. In such configuration, deviation from the parallel plane approximation should be mainly due to the heterogeneity of the flow field, which will be well described by a first order weak disorder perturbation approximation.

H

  in which z is perpendicular to the mean fracture plane, and x and y are parallel to it. The dilation ratios l and l are related through l s l z , where z is called the

Fig. 1 .

 1 Fig. 1. Schematic representation of the experimental rough crack. Rough faces are complementary and in Ž . Ž . perfect contact at the begining of the experiment a . Then, the lower block is moved vertically b by a Ž . Ž . distance h , laterally by a distance u c , and back vertically to contact by a distance h d . The vertical 1 2

v

  In the laminar flow regime, variations of the local spacing between fracture surfaces Ž . as well as the tortuosity of the flow paths may influence the permeability[START_REF] Ge | A governing equation for fluid flow in rough fractures[END_REF] . Ž The simplest way to take these effects into account is to use the Reynolds or . Ž . lubrication approximation Pinkus and Sternlicht, 1961 . In this case, one considers Ž . Ž only the in-plane x, y component of the velocity field the vertical component is . ignored , whose distribution is approximated by a local parabolic Poiseuille profile. The Ž . 2 square a x, y of the local aperture represents a local Darcy permeability, and the macroscopic permeability depends on the distribution of the local apertures. The local ™ ™ velocity V should then be related to the viscosity m and to the pressure gradient = P by

Fig. 2 .

 2 Fig. 2. Mean Fourier power spectra of height profiles measured on the surface of fractured granite blocks.Ž . Dashed lines correspond to the slope of the estimated rough exponent z of Table1. a Superposition of two Ž . power spectra obtained from measurements carried out on two perpendicular directions on sample no. 2. b Superposition of two power spectra obtained from measurements carried out on samples nos. 1 and 2.

  z

Fig. 3 .

 3 Fig. 3. Experimental set-up for measuring transport properties on an artificial fracture made of two fractured Ž granite or basalt blocks translated and brought back to contact along the fracture surface horizontal and top . view .

  mented by measuring the vertical spacing h u s u between the two blocks kept in ™ Ž . contact as a function of the horizontal displacement u Plouraboue et al., 1995 . Éxperimental results obtained with sample nos. 2 and 3 for four different orientations of ™ u at 908 intervals are shown in Fig. 4a-b. The four curves coincide quite well up to relative displacements urL ) 0.05 for which "30% deviations may be observed. Let us ™ Ž . now investigate in more details these variations of h u as a function of the orientation of the displacement. ™ Ž . Fig. 5 displays for sample no. 2 the variation of h u as a function of the angle u of ™ the vector u with respect to a reference axis. The data correspond to three displacement < < amplitudes u equal to 1, 3 and 10 mm, and to 32 different u values. The corresponding Ž . Ž . ² Ž .: values of the standard deviation s u of h u from its mean h u are shown in u h Ž . Ž . Table

  . Mean aperture h u in mm vs. horizontal relative displacement u r L Ls 25 cm for the two parts of fractured granite and basalt blocks displaced relative to each other, with respect to their initial adjustment Ž position while remaining in contact at least on one point. The four symbols open and shaded triangles and . circles correspond to orientations of the displacement in four different directions noted x q, x y, y q and Ž . Ž . y y in Fig.3. a Sample no. 3. b Sample no. 2. overall variation. This suggests that averages performed over four different orientations u are sufficiently representative. As expected, h u ru decreases Ž . linearly with log urL , except at very small displacements. One possible origin of the discrepancy is the fact that, due to damage associated with the fracturing process at small scales, the two halves of a fractured block do not coincide perfectly when put back together. Introducing a small constant shift of roughly 100 mm on the aperture values Ž . Plouraboue et al., 1995 improves the agreement with the theoretical predictions down´ Spatial mean aperture h u,u y h u , as a function of the orientation us ue u with u g 0,2p u of the horizontal displacement for constant displacement moduli equal to us1, 3 and 10 mm. Respective ² Ž .: values of aperture h u are: 0.794, 1.78, 4.47 mm.

Fig. 6 .

 6 Variations of the normalized squared mean vertical displacement h u r u as a function of ln u r L Ž . for the granite and basalt samples used in the present study. e corresponds to average of h u over four different displacement orientations at right angle from each other, ) corresponds to average over 32 different ™ Ž . Ž . orientations of u for the three data series of Fig. 4a. a Granite sample no. 3. b Basalt sample no. 2. from the data of Fig. 6a-b are listed in Table

  the electrical conductivity of the fluid and K is a geometrical factor 0 el

Fig. 7 .

 7 Fig. 7. Variations of the electrical admittance with the relative vertical displacement of two fractured blocks in ) Ž a radial measurement geometry. Dotted line with symbols corresponds to a pure vertical displacement no . horizontal shift . Long-dashed line corresponds to theoretical expectation with a parallel plane model. The four ™ Ž . Ž . symbol sets correspond to different directions of displacement u. a Granite sample no. 3. b Basalt sample no. 2.

  Variations of QrD P vs. the mean aperture h u for the basalt sample no. 2.

  which the two translated surfaces are in contact. The procedure is repeated for a ™ large number of realisations of the surface and of orientations of the displacement u. ² Ž .: Ž . The mean value h u and the standard deviation s u of the normal displacement h h are finally computed. Fig. 10 displays variations of these two parameters with urL obtained from a large number of simulations in which the surface size increases from 128 = 128 to 1024 = 1024 points. These quantities increase with the relative displacement urL, as observed

Fig. 10 .

 10 Fig. 10. Numerical simulations of the spatial mean aperture fluctuations. We plot in the log-log coordinates ² : Ž . the variations of the mean aperture h averaged over 200 realisations solid lines , the mean square deviation Ž . ² : Ž . of the aperture s dashed lines , and the relative fluctuations of the aperture s r h dotted lines vs. the h h relative displacement of the joint.

  v

  s d z s Mu . Using Eq. 3 , one sees that the expansion parameter s a rh u u Ž . varies only slowly with u logarithmic dependence . In the Reynolds approximation, by developing the effective conductance K with respect to the displacement u, K s h u 1 y c urL . The function c urL characterizes the deviael Ž . tion from the parallel plane model and varies only slowly with u: c urL s X Ž .

v

  A third possible origin of the deviation from the parallel plane model is a small angle of the plane of the joint with respect to the reference horizontal plane. Let us first estimate a simple correction in the case of a parallel flow. The electrical conductance is classically bounded by the harmonic and arithmetic spatial averages of the local Ž . conductance, which is proportional to the local aperture a x, y . In the limit of small u Ž . fluctuations of a x, y and small relative displacements urL < 1, both upper and u lower boundaries of the conductance have the same expression to first order uf 1rS a x , y d xd yh u 1 q 10

Ž

  One recovers the previous effect due to the local tilt slope of the surface 0 near the inner electrode, and the relative change of the conductance K compared to the el case of parallel plates that K obeys ref . a Height z r,u of the surface along the boundary of the injection hole obtained from a mechanical ™ Ž . Ž . profilometry. b Correction, K y K rs in mm , vs. the angle u between u and reference axis, us u : e us 200 mm, us 400 mm, v us1000 mm. with f s z r ,u cos u du rr 12 polar coordinates r,u are centered at the injection point.

  Fig. 11b displays the correction deduced from a complete calculation of the radial ™ 5 5 geometry correction on conductances for different displacement amplitudes u s u and ™ different angle u between u and the reference axis Ox. This correction has to be Ž . Ž .Ž . compared to 2p h u rlog Rrr , 1.43h u , which can be estimated from Fig.8a, and 0 indicates that the mean plane tilt near the injection hole is responsible for an 18% deviation in conductance. This order of magnitude is in agreement with experimental Ž . measurements. A linear rescaling of the corrective term K y K ru with u gives a el ref satisfactory superposition supporting the linear dependence of the correction K y K in the first regime i of Eq. 11 . However, a rescaling K y K ru el ref

  statistical analysis of h u for a self-affine surface was ´Ž suggested. It was demonstrated theoretically and experimentally Plouraboue et al., .

Table 1

 1

	Sample number	Roughness exponent z	Ž . s cm z	Roughness magnitude M
	1	0.83"0.05	1.33	0.19 cm	0.17
	2	0.79"0.06	0.94	0.15 cm	0.21
	3	0.80"0.06	1.39	0.22 cm	0.2

Acknowledgements

We wish to gratefully acknowledge Y. Berthaud and J. Schmittbuhl for their help with the profilometry measurements, O. Pinet, J. Lacoste and S. Valette for their participation in the transport properties experiments, and O. Brouard and D. Vallet for the preparation of the experimental set-up. We also thank E.G. Flekkøy for the stimulating discussions and S. Tyler for his remarks and a careful reading of the Ž manuscript. This work has been supported by the ECODEV AARC Geothermie des


Improved Fréchet Distance for Time Series

Ahlame Chouakria Douzal, P.