La région d’étude se caractérise par les dépôts d'une épaisse série sédimentaire couvrant l'intervalle stratigraphique Mésozoïque-Cénozoïque. La nomenclature des différentes formations étudiées ici est adaptée de celle de Rerbal (2008). Ainsi, on distingue dans l'ordre ascendant la succession lithostratigraphique suivante:

La Formation d'El Rhelida: incomplète dans la région d’étude (18 m d’épaisseurs), cette formation est attribuée à l'Albien supérieur (Vraconien). Elle est constituée essentiellement d'une série argilo-gréseuse admettant quelques niveaux carbonatés (calcaires et dolomies). Stratigraphiquement, la série repose sur les dépôts deltaïques de la formation de Tiout d’âge Néocomien-Albien (Rerbal, 2008).

Enfin les formations rapportées au Paléogène et au Néogène sont représentés par des dépôts détritiques essentiellement continentaux ne dépassant pas 1000 m. (Guillemot & Estorges, 1981; Rerbal, 2008).

Les observations et les mesures des fractures et des bandes de déformation ont été effectuées à travers 4 stations de mesures, 2 stations sur le flanc court ou avant (forelimb) de la structure (1 et 2) et 2 stations sur le flanc long ou arrière (backlimb) (3 et 4) (Fig. 1B; Tableau 1). La méthodologie d'approche consiste à:

A analyser la fractalité du réseau de failles à partir de la carte géologique (Fig. 1B), ainsi que la géométrie fractale en 2 dimensions des mesures de réseaux de fracturation et des bandes de déformation relevées sur le terrain, afin de calculer la dimension fractale à l'aide du logiciel Fractalyse 2.3.2 développé par Frankhauser et al. (2005) selon la méthode “The Two Point Correlation Function Method” de Hentschel & Procaccia (1983) et Grassberger & Procaccia (1983). (Equation 1),

θ étant l'angle en degré entre le profil d’échantillonnage et les fractures d'une même famille subparallèles. Le traitement des données de fracturation à l'affleurement ont été analysés à l'aide du logiciel Matter Cliff 1.0 développé par Deleze et al. (2005) pour l’étude statistique et la distribution des familles de discontinuités, et l’établissement des types de lois des espacements.

A calculer les paramètres d'intensité et de densité (Pxy). A savoir, P0, P10, P20, P21 et P32 des fractures sur un profil d’échantillonnage selon le modèle de Dershowitz & Herda (1992).

L'analyse fractale effectuée par la méthode “The Two Point Correlation Function Method” pour la carte de la région d'El Kohol (Fig. 3A) montre un caractère fractal de dimension 2 pour la répartition des failles. En effet, la courbe estimée du graphe bi-logarithmique de la fonction de corrélation [C(r)] en fonction de la distance entre fractures (r) peut être ajustée par une loi de puissance (Fig. 3B). La dimension fractale (Dc) obtenue est de 1,42; donc comprise entre 1 et 2, avec un coefficient de corrélation (R²) de 0,99. Cette valeur Dc calculée est cohérente avec les résultats publiés dans la littérature sur les études en 2D des réseaux de fracturation (Bonnet et al., 2001).

L'interprétation de la rosace des orientations des failles dans la région d'El Kohol (Fig. 3C) a permis de distinguer une famille N-S (dextre), NO-SE (dextre) et NE-SO (senestre) avec des proportions respectivement de l'ordre de 17,54%, 14,02% et 10,51%.

A partir du fichier des données de la carte numérisée, la longueur des failles à été calculée puis analysée. Les longueurs des fractures sont comprises dans l'intervalle [220,34 m: 8161,50 m] et s’étendent sur plusieurs ordres de grandeurs. Le diagramme log-log des données de fréquences cumulées des longueurs (classes de longueur = 50 m) est représenté sur la figure 3D. Seule une partie de la courbe entre 695 m et 3145 m à un comportement linéaire et peut être ajustée par une loi de puissance présentant une droite de pente (a) égal à 0,59 avec un coefficient de détermination de 0,97. La loi de puissance (Fig. 3D) déduite à partir de la distribution des fréquences, ne rend pas compte des valeurs inférieures à 695 m et supérieures à 3145 m. Ces biais sont causés par les effets de taille des failles. Les grandes failles présentent une probabilité plus élevée d’être échantillonnées que les plus petites et les failles plus courtes qu'une longueur donnée ne sont pas échantillonnées. D'après Odling et al. (1999), la distribution des longueurs de failles selon une loi de puissance présentent des seuils inférieurs (troncature) et des seuils supérieurs (censure), au deçà et au-delà desquels la loi de puissance n'est pas valide. La distribution en loi de puissance se situe entre 695 m et 3145 m et couvrent seulement 2 ordres d’échelle après le retrait des artefacts de censure et de troncature.

Dans notre approche, nous avons considéré sur un profil d’échantillonnage l'espacement entre fractures appartenant à une même famille par stations de mesures (Fig. 1B). Selon Bodin & Razack (1999), les travaux sur les lois de distribution des espacements entre les fractures d'un réseau remontent à plusieurs dizaines d'années (Priest & Hudson, 1976), et demeurent encore sujet de controverse (Dershowitz & Einstein, 1988). Dans cette optique, des différentes lois de distribution ont été ajustées aux données d'espacement: loi log-normale, négative exponentielle, loi normale et loi de puissance (Bonnet et al., 2001; Peacock & Mann, 2005; Zazoun, 2008; Kruhl, 2013).

Nous avons calculées la distribution des paramètres P0, P10, P20, P21 et P32 pour l'ensemble des familles de fractures dans le plan XY. Les résultats sont présentés dans le Tableau 2.

P0 étant le nombre total de fractures pour une même famille,

avec P0=N   [3]

P10 étant la densité linéique ou le nombre de fractures observées pour une même famille interceptées sur un profil d’échantillonnage de longueur L. L'inverse de P10 représente l'espacement moyen des fractures d'une même famille (Ortega et al., 2006).

avec P10=N/L   [4]

P20 est le nombre de fractures d'une même famille par unité de surface (m²),

avec P20=N/m²   [5]

P21 représente la longueur totale des fractures d'une même famille (mètre) dans une zone spécifiée (m²),

avec P21=L totale par m²   [6]

P32 correspond à la densité tridimensionnelle de fracturation et elle est égale à la somme des surfaces des fractures par unité de volume.

avec P32=S(m²)/m³   [7]

Le flanc court ou avant de la structure montre un P0 égal à 54 (familles de fractures F1, F2 et F3), alors que le flanc long ou arrière, le P0 est égal à 49 (familles de fractures F1, F2, F4 et F5). Au niveau du flanc court, nous avons remarqué que la valeur importante de P10 pour la famille F3 (P10=6,09), est liée probablement au faible espacement moyen des fractures (5,06 cm) à la proximité du plan axial du pli. Pour le flanc long ou arrière la forte valeur de P10 est observée pour la famille F4 (P10=3,68) et F5 (P10=8,76). Cette forte valeur constitue une réponse à l'influence de la longueur de la fracture sur la distribution de P10 (Bisdom, 2011). La valeur du P20 est généralement inférieure au P10. Ceci est le cas pour toutes les familles observées aux niveaux des stations et cela signifie que les fractures courtes sont prises en considération et elles sont mesurées. Toutefois, la famille F2 sur les deux flancs présente une valeur de P20 identique et égale à 0,04. La famille de fracture F1 montre des valeurs de P20 différentes pour les stations 2 et 3. La forte valeur (0,27) de cette dernière est observée sur le flanc court et elle est liée également au faible espacement moyen des fractures (2,54 cm) et de la proximité du plan axial du pli. Quant aux familles F4 et F5, elles présentent la même valeur de P20, égale à 0,07. Sans exceptions, la distribution de P21 est très similaire à la distribution du P10. Sa distribution a l'avantage sur le P20 qu'elle ne dépend pas uniquement du nombre de fractures, mais aussi sur la longueur des fractures (Bisdom, 2011). L'utilisation du paramètre P21 en 2D calculé comme une fonction de la distribution des longueurs de fractures peut être utilisé pour prédire la distribution de l'espacement d'une famille de fractures (Fig. 5A).

Le graphe bi-logarithmique du paramètre P32 en fonction du paramètre P10 d'une part et le graphe bi-logarithmique du paramètre P32 en fonction du paramètre P20 d'autre part montrent une relation proportionnelle (Fig. 5B et 5C).

Le flanc long ou arrière de la structure présente un réseau ou les fractures sont plus longues et moins denses. A l'inverse, le flanc court ou avant présente les fractures les moins longues et les plus denses. Il est à noter que la connectivité d'un réseau de fracturation peut être estimée de manière qualitative en regardant la densité de fractures du P21. Lorsque le P21 est élevé, cela signifie que le potentiel d’écoulement est plus favorable du fait des grandes longueurs de fractures et lorsque les paramètres de densité P10 et P20 sont presque identiques, il est attendu à ce que le réseau de fractures présente des caractéristiques d’écoulement de fluide (connectivité et capacité de stockage) similaires (Bekkers, 2013). Dans cette analyse, nous distinguons que la famille de fracture F1 présente le plus grand potentiel d’écoulement (1,44≤P21≤1,96) sur les deux flancs de la structure et ceci malgré que le réseau présente des caractéristiques d’écoulement de fluide non similaires (P10 différent de P20) ( Tableau 2 et Fig. 5D).

L'unité mécanique (Gross, 1993; Laubach et al., 2009) correspond à une ou à un groupe de couches qui affichent des motifs de ruptures similaires, selon des attributs géométriques et des propriétés physiques de la roche hôte (Lavenu et al., 2013). Nous avons adopté la méthode décrite par Underwood et al. (2003), qui consiste à digitaliser les fractures et les couches à partir d'une photographie (Bertotti et al., 2007). L'analyse mécanostratigraphique (Fig. 6A) montre une subdivision des trois formations en 12 unités mécaniques. De cette analyse, on distingue de bas en haut une unité mécanique pour la formation d'El Rhelida, sept unités mécaniques pour la formation de M'daouer et 4 unités mécaniques pour la formation de Rhoundjaïa. Le coefficient de variation calculé pour chaque unité montre des valeurs inférieures à 1, ce qui permet de souligner le contrôle de la fracturation par l’épaisseur des bancs ( stratabound). Il est communément admis que l'intensité de la fracturation dépend de l’épaisseur des bancs ou des unités mécaniques (Wu & Pollard, 1995) et les bancs les plus minces ont tendance à être plus fracturés que les bancs épais (Hobbs, 1967). En effet, la figure 6B illustre bien cette relation où l'intensité de la fracturation est inversement proportionnelle à l’épaisseur de l'unité mécanique (pente négative avec R²=0,78). Par ailleurs, un raisonnement mettant en relation l’épaisseur des formations et l'intensité de la fracturation ne permet pas d'arriver à la même conclusion (Fig. 6C). Les valeurs du coefficient d'aplatissement (kurtosis) calculées à partir des espacements des fractures sur les unités mécaniques montrent des valeurs comprises entre −1,39 (UmM-5) et 1,03 (UmR-2b). D'après Ruf et al. (1998), les fractures perpendiculaires à la stratification présentent généralement de faibles valeurs du coefficient d'aplatissement, et leur nombre est fonction de l’épaisseur du banc (stratabound). Notre analyse corrobore ces observations et l'essentiel des résultats est résumé dans le Tableau 3. La relation entre la longueur et l'espacement des fractures peut être quantifiée sous l'appellation de l'index de l'espacement des fractures (FSI pour Fracture Spacing Index) (Narr & Suppe, 1991; Ruf et al., 1998), qui est égal à la pente de la droite de régression du diagramme mettant en relation l’épaisseur des unités mécaniques en fonction de l'espacement médian des fractures. La courbe tracée ( Fig. 6D), peut être ajustée par une loi de puissance dont la pente est égale à 1,14. Les valeurs typiques pour le FSI sont comprises généralement entre 0,8 et 1,5 (Narr & Suppe, 1991; Groos, 1993; Engelder et al., 1997), ce qui démontre que la valeur calculé dans notre cas reste dans la fourchette au regard des limites mentionnées dans la littérature. Toutefois, la pente négative soutient que l'espacement le plus faible entre les fractures est observé dans les unités mécaniques les plus épaisses, ce qui est contraire aux résultats publiés dans la littérature (Ladeira & Priece,1991; Engelder et al., 1997). Selon Eyal et al., ( 2001), une pente négative indique que l’épaisseur de l'unité mécanique n'est pas seule responsable dans le contrôle des espacements des fractures et une possible explication est à rechercher probablement dans le contraste lithologique. Nous avons calculé pour chaque unité mécanique le coefficient FSR (Fracture Spacing Ratio) défini par Gross (1993) qui est égal à l’épaisseur de l'unité mécanique divisé par l'espacement médian des fractures (Tableau 3). Le FSR permet de montrer qu'au sein d'une même formation, l'intensité de la fracturation peut varier d'une unité mécanique à une autre (Eyal et al., 2001). Le FSR calculé dans notre cas ( Tableau 3 et Fig. 6A) montre cette variation avec un gradient moyen décroissant de la base (Formation de Rhelida) vers le sommet (Formation de Rhoundjaïa). Dans le détail, le FSR peut être subdivisé en 3 niveaux (FSR1, FSR2 et FSR3), toutes trois décroissants de la base vers le sommet. Les fortes valeurs matérialisant les bases des niveaux sont liées aux unités mécaniques UmR-2b, UmM-4 et UmE-1 qui correspondent respectivement aux niveaux francs de calcaires bioturbés, calcaires et grès. Cette décroissance pourrait probablement être mise en relation avec le gradient décroissant de la compaction de la base vers le sommet des formations, donc une diminution probable de la porosité du sommet vers la base. En effet, Engelder et al. (1997), dans leur étude de l'anticlinal de Elk Basin (Montana-Wyoming) relient les valeurs de grande porosité avec les valeurs faibles de FSR. La décroissance du FSR est à mettre aussi avec la diminution de la distance des unités mécaniques par rapport au décro-chevauchement principale d'El Kohol (ADC) ( Fig. 1C), sans toutefois écarter l'hypothèse de l'existence de plusieurs événements de déformation ayant aboutit au réseau de fracturation observé actuellement. En plus des familles de fractures F1 et F2, décrites précédemment, l'affleurement montre l'existence de deux familles de fractures F4 et F5, de direction et de pendage respectivement N010°-76°O et N130°-55°NE (Fig. 6A). Ces dernières, à l'inverse des autres familles de fractures, présentent une grande extension verticale indépendamment des unités mécaniques qu'elles traversent. L'analyse fractale effectuée sur les fractures digitalisées en figure 6A montre un caractère fractal de dimension 2 pour la répartition des fractures. La courbe estimée du graphe bi-logarithmique de la fonction de corrélation [ C(r)] en fonction de la distance entre fractures (r) peut être ajustée par une loi de puissance. La dimension fractale (Dc) obtenue est de 1,34; donc comprise entre 1 et 2, avec un coefficient de corrélation (R²) de 0,99 (Fig. 6E).

Les roches poreuses tels le grès, le calcaire, le tuff pyroclastique et la brèche, organisés en réseaux cataclasique présentent parfois des bandes de déformation (BD) en réponse à une déformation par cisaillement, conduisant finalement à des zones de failles géométriquement complexes (Aydin, 1978; Aydin & Johnson, 1978; Shipton & Cowie, 2001; Fossen et al., 2007). Les (BD) sont des structures d’épaisseur millimétrique, caractérisés par une réorganisation de grain (glissement de grain par rotation et/ou à de la fracturation associée à de la dilatation, du cisaillement et/ou des mécanismes de compactage). Parmi les différents types de bandes classées en fonction des mécanismes de génèse (cataclastiques, désagrégation, phyllosilicates, solution et cimentation) (Fossen et al., 2007; Schueller et al., 2013), nous nous intéresserons uniquement à un seul type de bandes de déformation présente dans la région. Il s'agit de bandes de déformation compactantes (BDc) (Fig. 7A et 7B) (Mollema & Antonellini, 1996) et cisaillantes (Fig. 7A) (Du Bernard et al., 2002). Les figures 7A et 7C, illustrent les bandes de déformation observées dans les grès de la formation d'El Rhelida et leurs orientations (Fig. 1B). Ces dernières montrent une augmentation de leurs nombres, ainsi qu'une réduction de leurs espacements à l'approche des accidents. Ainsi, la famille BD1 se concentre à l'approche de l'accident décro-chevauchant d'El Kohol (ADC), et les familles BD2a et BD2b sont à mettre en relation avec les décrochements NO-SE dextres (FEK) et N-S senestres ( Fig. 1B). L'analyse fractale effectuée par la méthode “The Two Point Correlation Function Method” montre un caractère fractal de dimension 2 pour la répartition des bandes de déformation. La courbe estimée du graphe bi-logarithmique de la fonction de corrélation [C(r)] en fonction de la distance entre fractures (r) peut être ajustée par une loi de puissance. La dimension fractale (Dc) obtenue est de 1,57 (R²=0,99) pour les bandes de déformation (BD1 et BD2a) et de 1,67 (R²=0,99) pour les bandes de déformation (BD1+BD2a+BD2b) (Fig. 7E et 7F). Le graphe mettant en relation les longueurs des bandes de déformation en fonction de l'orientation de ces dernières (Fig. 7 G), montre que les bandes de déformation les plus longues dont les orientations sont comprises entre 130° et 150°, sont liées probablement aux décrochements NO-SE dextres (FEK). La distribution des fréquences cumulées des longueurs de bande déformation montre une loi de puissance entre 272 mm et 542 mm (Fig. 7H) caractérisée par une droite de pente −0,05 (R²=0,97). Notons que la famille F1 est parallèle à la bande de déformation BD1 et au déco-chevauchement principal d'El Kohol et les bandes de déformation BD2a et BD2b sont respectivement plus ou moins parallèles aux familles F5 et F2+F4.

A partir des données de terrain (Tableau 2, 3 et 4), nous avons tenté une simulation stochastique du réseau de fracturation par l'approche de la modélisation des réseaux de fractures discrètes (DFN) pour l'unité mécanique UmR-1 (Fig. 6A). Cette méthode permet de définir des modèles de distribution des paramètres relatifs à la densité, à l'orientation et à la taille des fractures à partir des données de terrain (Henrion, 2011). Le nombre et la répartition spatiale des centres des fractures sont générés stochastiquement par le processus de Poisson en fonction d'un modèle de densité de fractures tandis que les caractéristiques géométriques (dimensions, orientation) de chaque fracture sont générées aléatoirement par le tirage de Monte-Carlo dans leurs lois de distribution respectives (Stoyan et al., 1995). Les résultats obtenus pour le flanc court et le flanc long de la structure sont observables en figure 8. Le réseau 2D (Fig. 8B et 8D) fourni résulte de l'intersection du réseau 3D (Fig. 8A et 8C) avec un plan de coupe horizontal représenté par le toit de l'unité mécanique UmR-1.