Fast geomagnetic field intensity variations between 1400 and 400 BCE: New archaeointensity data from Germany

Physics of the Earth and Planetary Interiors 270 (2017) 143–156

Contents lists available at ScienceDirect

Physics of the Earth and Planetary Interiors journal homepage: www.elsevier.com/locate/pepi

Fast geomagnetic field intensity variations between 1400 and 400 BCE: New archaeointensity data from Germany Gwenaël Hervé a,⇑, Jörg Fabbinder a,b, Stuart A. Gilder a, Carola Metzner-Nebelsick c, Yves Gallet d, Agnès Genevey e, Elisabeth Schnepp f, Leonhard Geisweid c, Anja Pütz g, Simone Reub c, Fabian Wittenborn h, Antonia Flontas c, Rainer Linke i, Gerd Riedel j, Florian Walter c, Imke Westhausen c a

Department für Geo- und Umweltwissenschaften, Ludwig-Maximilians Universität, Theresienstrasse 41, 80333 München, Germany Bayerisches Landesamt für Denkmalpflege, Hofgraben 4, 80539 München, Germany Institut für Vor- und Frühgeschichtliche Archäologie und Provinzialrömische Archäologie, Ludwig-Maximilians Universität, Geschwister-Scholl-Platz 1, 80539 München, Germany d Institut de Physique du Globe de Paris, UMR 7154 CNRS, Paris, France e Laboratoire d’archéologie moléculaire et structurale (LAMS), UMR 8220 CNRS, Paris, France f Paleomagnetic Laboratory Gams, Montanuniversität, Leoben, Austria g Aschheim Museum, Germany h Institut für Ur- und Frühgeschichte, Heidelberg Universität, Germany i Museum Königsbrunn, Germany j Stadtmuseum Ingolstadt, Germany b c

a r t i c l e

i n f o

Article history: Received 10 May 2017 Received in revised form 5 July 2017 Accepted 18 July 2017 Available online 19 July 2017 Keywords: Archaeointensity Geomagnetic secular variation Central western Europe Bronze age Iron age Ceramics

a b s t r a c t Thirty-five mean archaeointensity data were obtained on ceramic sherds dated between 1400 and 400 BCE from sites located near Munich, Germany. The 453 sherds were collected from 52 graves, pits and wells dated by archaeological correlation, radiocarbon and/or dendrochronology. Rock magnetic analyses indicate that the remanent magnetization was mainly carried by magnetite. Data from Thellier-Thellier experiments were corrected for anisotropy and cooling rate effects. Triaxe and multispecimen (MSPDSC) protocols were also measured on a subset of specimens. Around 60% of the samples provide reliable results when using stringent criteria selection. The 35 average archaeointensity values based on 154 pots are consistent with previous data and triple the Western Europe database between 1400 and 400 BCE. A secular variation curve for central-western Europe, built using a Bayesian approach, shows a double oscillation in geomagnetic field strength with intensity maxima of 70 mT around 1000-900 BCE and another up to 90 mT around 600-500 BCE. The maximum rate of variation was 0.25 mT/yr circa 700 BCE. The secular variation trend in Western Europe is similar to that observed in the Middle East and the Caucasus except that we find no evidence for hyper-rapid field variations (i.e. geomagnetic spikes). Virtual Axial Dipole Moments from Western Europe, the Middle East and central Asia differ by more than 2
1022 A
m2 prior to 600 BCE, which signifies a departure from an axial dipole field especially between 1000 and 600 BCE. Our observations suggest that the regional Levantine Iron Age anomaly has been accompanied by an increase of the axial dipole moment together with a tilt of the dipole. Ó 2017 Elsevier B.V. All rights reserved.

1. Introduction Knowing how geomagnetic field intensity has changed through time represents an important challenge in palaeomagnetic research. Such information is crucial, for example, to understand how the geodynamo functions (e.g. Aubert et al., 2013). Defining the geomagnetic dipole moment over centennial timescales is also ⇑ Corresponding author at: CEREGE, UMR 7330 CNRS, AMU, IRD, Europôle de l’Arbois, BP80, 13545 Aix-en-Provence, cedex 04, France. E-mail address: [email protected] (G. Hervé). http://dx.doi.org/10.1016/j.pepi.2017.07.002 0031-9201/Ó 2017 Elsevier B.V. All rights reserved.

essential to reconstruct the past solar activity through the production rate of cosmogenic nuclides (e.g. Usoskin et al., 2016). Finally, once a time-dependent palaeointensity curve is established for a region, one can use it to date archaeological sites. An absolute estimate of field strength relies on the study of the thermoremanent magnetization of ferrimagnets hosted within a material that cooled through their Curie temperatures. The growth of the worldwide reference database in the past ten years (e.g. Brown et al., 2015; Donadini et al., 2006; Genevey et al., 2008) has facilitated the construction of increasingly sophisticated global geomagnetic models over the last three to fourteen thousand years

144

G. Hervé et al. / Physics of the Earth and Planetary Interiors 270 (2017) 143–156

(e.g., Constable et al., 2016; Licht et al., 2013; Nilsson et al., 2014; Pavón-Carrasco et al., 2014a for the most recent examples). Despite these progresses, the models are still limited by an inhomogeneously-distributed data in space and time and by an unequal data quality, especially for the BCE period (e.g., PavónCarrasco et al., 2014b). The period between 1000 and 500 BCE contains the highest virtual axial dipole moments in the Holocene, with a peak around 40% higher than at present (Genevey et al., 2008; Usoskin et al., 2016). This is a key period where experimental constraints are needed to improve our understanding of the geodynamo. Moreover, in the Levantine area, Shaar et al. (2011, 2016) found that the fastest secular variation rates occurred at that time. They found two geomagnetic spikes at the beginning of the tenth and eighth century BCE. However, these sharp increases in field intensity at the subcentennial scale are debated because they would require variation rates of several mT/yr, inconsistent with our present understanding of core flow patterns at the surface of the Earth’s outer core (Livermore et al., 2014). The spikes could be associated with a regional geomagnetic anomaly centred on the Middle East between 1050 and 700 BCE (Shaar et al., 2016, 2017). These authors called it the ‘‘Levantine Iron Age Anomaly” (LIAA), which is characterized by a high average geomagnetic field with VADM values close to 14
1022 A
m2 and by a 20° deviation from the Geocentric Axial Dipole (GAD) direction. The LIAA extended to Caucasus and Central Asia but maybe not to Europe (Shaar et al., 2017). In Western Europe, the direction of the geomagnetic field during the LIAA period deviated 12–14° from the GAD direction (Hervé et al., 2013a), but intensity secular variation remains poorly constrained between the 15th and the 4th century BCE, as only twenty high-quality archaeointensity values are available (e.g., Hervé et al., 2013b; Kapper et al., 2015). These data hint toward fast secular variation with possibly two successive intensity maxima, yet no evidence for the occurrence of geomagnetic spikes (Hervé et al., 2016). Better temporal resolution and better knowledge of the amplitude and the rate of secular variation are clearly needed. For these reasons, we collected samples from Late Bronze Age and Early Iron Age archaeological sites in Bavaria. After presenting the archaeological context and sampling, we describe the laboratory procedures using three complementary protocols. The results shed new light on secular variation in central-western Europe from 1400 BCE to the beginning of the CE period. We discuss the speed and amplitude of secular variation from Europe to Central Asia and their geomagnetic implications.

2. Archaeological context The archaeological sites were located close to Munich in Bavaria (Germany) (Fig. 1a), with each set of ceramics coming from a distinct feature (grave, pit or well). The 52 sampled features were dated between 1400 and 400 BCE within the end of the Middle Bronze Age, the Late Bronze Age (also called the Urnfield period in South Germany), the Early Iron Age and the beginning of the Late Iron Age (Fig. 1b). Dating methods include radiocarbon, dendrochronology and, in most cases, archaeological constraints, stratigraphy and typochronology. Typochronology of centralwestern Europe is divided in stages defined by a characteristic assemblage of ceramics and metallic objects (e.g. Müller-Karpe, 1959; Sommer, 2006). The three main phases are called Bronze (Bz), Hallstatt (Ha) and La Tène (Lt), themselves being divided in four subphases (A to D) (Fig. 1b). Each sampled feature was dated by comparing its artefacts against the archaeological chronology. Because the style of the artefacts did not change rapidly, the resolution of the reference chronology is generally not better than 150– 200 years.

All the graves that we sampled were dated using typochronology. They were excavated in the cemeteries of Ingolstadt – Mailing Schindergrubäcker (Malcher and Weinig, 1995), München – Obermenzing (Winghart, 1984), Aschheim – Ostumgehung (Pütz, 2008) and, in the Lech valley, of Kleinaitingen – Gewerbegebiet Nord, Kleinaitingen – Kiesgrube Weigl, Königsbrunn – Firma Ampack, Oberottmarshausen – Kiesgrube Lauter and Oberottmarshausen – Leberbichl (Büttner et al., 2006). In the Grünwald – Gymnasium settlement (Metzner-Nebelsick et al., 2016), we sampled an Early Iron Age grave and two Late Bronze Age pits, called grgy661 (Fig. 1c) and grgy665. The stratigraphy in the grgy661 pit clearly differentiated ritual deposits in the bottom layers and filling of the pit in the upper layers. The ceramics we sampled were distinguished in time according to this stratigraphy. Two AMS radiocarbon dates on animal bones associated with the sampled ceramics dated pit grgy661 (Poz-80712, 2955 ± 35 uncal. years BP and Poz-80713, 2950 ± 35 uncal. years BP). The synthesis of these dates using Bayesian Chronomodel software (Lanos and Philippe, 2017) constrained the interval [1270; 1045] BCE at 95% of confidence. We assumed the same date interval for pit grgy665 regarding the stratigraphical context. Besides the graves, we sampled two wells with wooden walls in Aschheim (Zach et al., 2010). Dendrochronological dates on the wood are 762 ± 10 BCE for well 692 of Aschheim – Ostumgehung (asou692) (Fig. 1d) and 700 ± 10 BCE for well 390 of Aschheim – XXXLutz (asxx390) (Herzig 2008a,b). However, we do not think that the ceramics filling the wells are dated so precisely. The date of the ceramic assemblage must take in consideration that older and younger ceramics may have been dumped into the well during and after its use. According to the position of the sampled ceramics in the filling layers and to the assumed lifetime inferred from the dendrochronological study, we respectively assigned [780; 720] BCE and [760; 600] BCE for the ceramic sets of asou692 and asxx390 (Table 1). The youngest of the studied sites, Haffstrabe in München – Trudering, was a settlement dated at the end of the Early Iron Age and at the beginning of the Late Iron Age (Bagley et al., 2010). There, we sampled ceramics from two wells and three pits.

3. Sampling and magnetic characterization In total, 453 ceramic sherds from 271 different pots were collected. In Oberottmarshausen – Kiesgrube Lauter, we also sampled ten fragments of burnt daub (obla246e to obla246n). The pots were generally common wares that, for the period, were not produced farther than a few tens of kilometres from the sampled sites. Sampling preferentially focused on red-coloured sherds because they were less sensitive to mineralogical alteration during the archaeointensity protocol than the grey-black sherds (e.g. Osete et al., 2016). Thermomagnetic curves measured with a Petersen Instruments variable field translation balance (VFTB) using a 200 mT field on chips from 109 pots confirmed this observation (Fig. 2a). Almost all heating-cooling cycles of red sherds were reversible up to 600 °C. Some black sherds also showed a high degree of reversibility, but most did not, even at lower temperatures (400 °C), so they were rejected from further experimentation. Most sherds from our collection were red only within a few millimeters from the surface, whereas their cores were greyish (Fig. 1e). This two-colour appearance is related to variations in oxygen conditions in the kiln during baking. In order to increase the success of the archaeointensity experiments, 100 mg specimens were cut from the red-coloured part of the sherds. They were fixed in 8 mm diameter quartz holders filled with quartz wool. In all 109 samples, the thermomagnetic curves highlighted a dominant ferromagnetic phase with Curie temperatures between

G. Hervé et al. / Physics of the Earth and Planetary Interiors 270 (2017) 143–156

145

Fig. 1. (a) Location map of the studied sites around Munich. (b) Archaeological chronology and temporal distribution of the sampled features. (c) Photograph of pit 661 from Grünwald – Gymnasium. The base layers, associated with ritual deposits, include numerous ceramic ware pots together with animal bones. (d) Photograph of well 692 from Aschheim – Ostumgehung (asou692). (e) Photograph illustrating two representative sherds and extracted specimens that were then packed into 8 mm diameter quartz holders.

540 and 580 °C, which indicated the presence of pure magnetite to magnetite <10% cation substitution (Fig. 2a). In 70% of the samples, backfield curves of an isothermal remanent magnetization (IRM) revealed a minor contribution of a higher coercivity phase, as the curves were not saturated by 300 mT (Fig. 2b). Hysteresis curves were performed with a Lakeshore-PMC vibrating sample magnetometer on the same 109 chips, which were used for the thermomagnetic curves. The samples lay in the pseudo-single domain (PSD) region of the Day plot along the single domain (SD) to multidomain (MD) mixing curve (Dunlop, 2002) (Fig. 2c). The presence of a high coercivity phase shifted them to the right of the mixing curve. No clear differences in rock magnetic properties were visible among the sites. For six samples, relatively small magnetic grain sizes (SD-PSD) were confirmed by first-order reversal curves (FORCs) (Fig. 2d). The 10 to +20 mT spread on the vertical axis indicated moderate interactions, which make our collected pots well suited for archaeointensity experiments.

4. Archaeointensity protocols Prior to the archaeointensity experiments, one specimen per pot was thermally demagnetized in a Schoenstedt furnace in order to identify the unblocking temperature range of the thermoremanent magnetization. After removing a minor magnetization component below 100–300 °C, 70% of the pots defined a linear magnetization component that trended toward the origin on orthogonal diagrams. This component was interpreted as a thermoremanent magnetization (TRM) acquired during the baking of the pots (blue in Fig. 3a and b). For 20% of the pots, the low temperature magnetization component dominated the unblocking spectrum, while the TRM was unblocked only above 450–500 °C (Fig. 3c). In such cases, the secondary component was interpreted as a partial TRM acquired during the use of the pot in a domestic

or funeral (pyre) context. When the unblocking spectra of the two components of magnetization could be clearly separated, we also tested the archaeointensity protocol on them, in order to have more successful archaeointensity results. The 10% remaining pots had more complex behaviour with two or three unclearly separated magnetization components; they were rejected from further analyses. Archaeointensities were determined using the classical Thellier-Thellier method (Thellier and Thellier, 1959) with partial thermoremanent magnetization (pTRM) checks and TRM anisotropy and cooling rate corrections. All experiments were conducted in the palaeomagnetic laboratory of LMU Munich in Niederlippach. Remanent magnetizations were measured with a three-axis, 2G Enterprises, cryogenic magnetometer. Heatings were performed using 9–13 temperature steps up to 575 °C in a Magnetic Measurements Thermal Demagnetizer (MMTD) furnace. Four different laboratory field values were used: 50, 60, 70 and 80 mT. We studied 548 specimens in total. At each temperature step, specimens were heated and cooled with the laboratory field applied first along their Z-axes and then in the opposite sense. Low-field susceptibility was measured after each temperature step in order to monitor possible changes in magnetic mineralogy. The anisotropy of TRM (ATRM) was determined at the step around 500 °C using 6 positions followed by a stability check (Chauvin et al., 2000). The anisotropy correction was not applied if the moments of the pTRM between the first and seventh (repeated) steps exceeded 10%. The cooling of the minispecimens in the MMTD furnace lasted about 15 minutes, which was much faster than the archaeological cooling. To take into account the cooling rate effect on TRM intensity, the procedure of Gómez-Paccard et al. (2006) was applied with four heating steps at one of the last temperature steps when 5–10% of the NRM remained (usually between 480 and 550 °C). The cooling decreased exponentially over 8 h. The ATRM and the cooling rate effects were corrected at the specimen level.

146

G. Hervé et al. / Physics of the Earth and Planetary Interiors 270 (2017) 143–156

Table 1 Archaeointensity data from this study. Abbreviations are: Age (yrs BCE), dating intervals in years BCE; AI, archaeointensity protocol: Thellier-Thellier (Th-Th), Triaxe and Multispecimen protocol – domain state corrected (MSP-DSC); Number of accepted pots (N) and specimens (n) with the number measured in brackets; FATRM+CR, average archaeointensity with standard deviation (SD) corrected for TRM anisotropy and cooling rate; F46°N,7°E, average archaeointensity relocated to 46°N, 7°E; VADM, Virtual Axial Dipole Moment. Feature

Age (yrs BCE)

Dating method

AI method

Npots/ nspecimens

FATRM+CR (mT) SDF F46°N,7°E VADM (mT) (mT) (1022 A
m2)

Aschheim, Ostumgehung (48.1°N; 11.7°E) Grave 47 (asou47) Grave 48 (asou48) Grave 51 (asou51) Well 692 (asou692)

[800; [800; [800; [780;

620] 620] 620] 720]

Arch. (Ha C) Arch. (Ha C) Arch. (Ha C) Dendrochronology

Th-Th Th-Th Th-Th Th-Th

2/4 3/8 3/5 4/6

72.0 63.0 63.5 59.2

7.2 17.4 13.3 2.6

Aschheim, XXXLutz (48.1°N; 11.7°E) Well 390 (asxx390)

[760; 600]

Dendrochronology

Th-Th

6/12 (8/19)

59.4

17.4 58.2

9.4

Th-Th + Triaxe 16/19 (20/34) 54.8 Th-Th + Triaxe 12/17 (16/25) 47.9 Th-Th + Triaxe 15/27 (17/39) 54.1

4.1 3.5 2.8

53.7 46.9 53.0

8.7 7.6 8.6

Th-Th Th-Th Th-Th Th-Th

1/4 2/4 2/5 3/4

93.9 68.9 48.1 53.7

6.2 12.2 1.7 6.5

91.3 67.0 46.8 52.2

14.8 10.8 7.6 8.5

Lechtal – Kleinaitingen, Gewerbegebiet Nord (48.2°N; 10.9°E) Grave 215 (klgn215) [1360; 1200] Arch. (Bz D)

Th-Th

2/6 (3/8)

54.2

2.9

53.0

8.6

Lechtal – Kleinaitingen, Kiesgrube Weigl (48.2°N; 10.9°E) Grave 24 (klwe24) [1360; 1200] Grave 68 (klwe68) [1100; 1020] Grave 71 (klwe71) [1360; 1200] Grave 162 (klwe162) [1200; 1100] Grave 186 (klwe186) [1360; 1200]

Th-Th Th-Th Th-Th + Triaxe Th-Th Th-Th

3/8 1/3 5/8 3/8 3/3

(3/11) (3/9) (5/15) (3/12) (5/12)

51.9 64.3 55.6 59.0 48.5

4.3 1.3 2.7 6.1 1.8

50.8 62.9 54.4 57.7 47.4

8.2 10.2 8.8 9.3 7.7

Th-Th Th-Th

2/5 (2/11) 2/4 (2/6)

51.6 61.5

0.9 3.7

50.5 60.2

8.2 9.7

Th-Th Th-Th Th-Th + Triaxe Th-Th + Triaxe + MSP-DSC

3/7 (3/9) 6/9 (7/17) 4/10 (9/17) 16/42 (16/57)

53.7 52.9 52.9 62.1

5.0 2.1 3.6 9.6

52.5 51.7 51.7 60.7

8.5 8.4 8.4 9.8

Lechtal – Oberottmarshausen, Leberbichl (48.2°N; 10.9°E) Grave 13b (oble13b) [1400; 1300] Arch. (Transition Bz C2 – Bz D) Th-Th Grave 16 (oble16) [1400; 1300] Arch. (Transition Bz C2 – Bz D) Th-Th Grave 17 (oble17) [1360; 1280] Arch. (first half Bz D) Th-Th

2/6 (2/8) 2/6 (2/6) 4/10 (5/16)

54.1 53.4 51.9

0.2 4.1 3.4

52.9 52.2 50.8

8.6 8.5 8.2

München – Obermenzing, Ausgrabungen 1984 (48.1°N; 11.4°E) Grave 22 (ober22) [1050; 800] Arch. (Ha B) Grave 26 (ober26) [950; 800] Arch. (Ha B2-3) Grave 50 (ober50) [950; 800] Arch. (Ha B2-3)

Th-Th Th-Th Th-Th

2/6 (2/6) 3/4 (7/11) 3/7 (6/13)

66.0 57.1 49.0

6.3 3.8 3.8

64.6 55.9 48.0

10.5 9.1 7.8

München – Trudering, Haffstrabe (48.1°N; 11.6°E) Well 74 (trha74) [460; Well 108 (trha108) [460; Pit 541 (trha541) [420; Pit 702A (trha702) [460; Pit 922 (trha922) [540;

Th-Th Th-Th Th-Th Th-Th Th-Th

3/8 (5/11) 5/11 (5/13) 2/7 (3/12) 6/16 (6/19) 3/8 (6/14)

72.4 79.2 60.3 74.9 86.9

7.9 6.5 0.8 4.9 2.5

70.9 77.6 59.0 73.3 85.1

11.5 12.6 9.6 11.9 13.8

Grünwald, Gymnasium (48.0°N; 11.5°E) Pit 661, Deposition layers (grgy661-3/4/6) [1270; 1045] Calibrated Pit 661, Filling layers (grgy661-1/2) [1270; 1045] Calibrated Pit 665 (grgy665) [1270; 1045] Calibrated Ingolstadt – Mailing, Schindergrubäcker (48.8°N; 11.5°E) Grave 4 (inma4) [660; 540] Grave 57 (inma57) [800; 620] Grave 60A (inma60) [800; 620] Grave 67 (inma67) [800; 620]

Arch. Arch. Arch. Arch.

Arch. Arch. Arch. Arch. Arch.

(Ha (Ha (Ha (Ha

14

C (95%) C (95%) 14 C (95%) 14

C2-D1) C) C) C)

(Bz D) (Ha A2) (Bz D) (Ha A1) (Bz D)

Lechtal – Königsbrunn, Firma Ampack (48.2°N; 10.9°E) Grave 4 (koam 4) [1360; 1200] Arch. (Bz D) Grave 8 (koam 8) [1220; 1120] Arch. (end Bz D – Ha A1) Lechtal – Oberottmarshausen, Kiesgrube Lauter (48.2°N; 10.9°E) Grave 8 (obla8) [1360; 1200] Arch. (Bz D) Grave 129 (obla129) [1360; 1280] Arch. (first half Bz D) Grave 130 (obla130) [1360; 1200] Arch. (Bz D) Grave 246 (obla246) [1200; 900] Arch. (Ha A) + 14C

360] 360] 320] 360] 440]

Arch. Arch. Arch. Arch. Arch.

(Lt A – begin Lt B1) (Lt A – begin Lt B1) (end Lt A – Lt B1) (Lt A – begin Lt B1) (Late Ha – Lt A1)

The Thellier-Thellier experiments were complemented using the Triaxe protocol of Le Goff and Gallet (2004) at the Institut de Physique du Globe de Paris. This protocol consists of continuous high-temperature magnetization measurements through a series of five heating and cooling steps with or without applied laboratory field (Le Goff and Gallet, 2004, see also Genevey et al., 2009, 2016). The magnetization was measured every 5 °C from 150 °C up to 465–550 °C depending on the unblocking temperature of the specimen. Experiments were conducted on 60, 1 cm-wide chips from 41 pots. Only pots having a large primary TRM and strictly reversible thermomagnetic curves were selected. The intensity of the laboratory field was within 5 mT of the value determined with the Thellier-Thellier protocol. The direction of the field

(4/9) (3/9) (3/8) (7/20)

(1/6) (6/12) (5/13) (4/13)

70.5 61.7 62.2 58.0

11.4 10.0 10.1 9.4

was adjusted in such a manner that the laboratory TRM was parallel to the NRM direction, which thus avoided TRM anisotropy effects. It was experimentally shown that the Triaxe protocol was largely immune to cooling rate effects (e.g. Le Goff and Gallet, 2004; Genevey et al., 2009). On four features (burnt daub fragments of obla246, pottery sherds of grgy661-3/4/6, grgy661-1/2 and grgy665), we performed the domain state corrected multi-specimen protocol (MSP-DSC) at the Gams paleomagnetic laboratory (Montanuniversität, Leoben, Austria). This method compares the NRM intensity with the one of a partial-TRM (pTRM) acquired in a laboratory field Hlab along the NRM direction (Dekkers and Böhnel, 2006). Different specimens from the same site were heated at the same temperature

G. Hervé et al. / Physics of the Earth and Planetary Interiors 270 (2017) 143–156

147

Fig. 2. Rock magnetic results from the Bavarian ceramics. (a) Representative thermomagnetic curves of red and black-colored ceramics. Heating curves are shown in red, cooling curves in blue. (b) Backfield curves of isothermal remanent magnetization. (c) Day plot of the remanence ratio (remanent saturation magnetization [Mrs]/saturation magnetization [Ms]) versus the coercivity ratio (coercivity of remanence [Hcr]/bulk coercive force [Hc]) with SD–MD and SP–SD mixing curves for magnetite (Dunlop, 2002). (d) FORC diagrams of two selected sherds. The FORCs were calculated from 320 loops and processed using VARIFORC software (Egli, 2013). (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

148

G. Hervé et al. / Physics of the Earth and Planetary Interiors 270 (2017) 143–156

Fig. 3. Accepted (a) and rejected (b) archaeointensity results estimated on the primary TRM. (c) Examples of accepted archaeointensity estimations on the secondary partial TRM component. The NRM-TRM diagrams were corrected by subtracting the remaining NRM at the junction temperature from the NRM at each previous heating step Solid circles on NRM-TRM diagrams indicate the temperature steps used in the intensity determination. The orthogonal diagrams compare the variation of the direction during the Thellier-Thellier protocol in black and during the thermal demagnetization in blue. The two sister-specimens do not have the same orientation. Open (solid) circles denote the projection on the vertical (horizontal) plane. Temperatures are in degrees Celsius. The arrow marks the direction of the applied laboratory field. (d) Archaeointensity comparisons of specimens from five pots with and without TRM anisotropy correction. (e) Histogram of the cooling rate correction factor of the 306 accepted specimens. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

G. Hervé et al. / Physics of the Earth and Planetary Interiors 270 (2017) 143–156

using different Hlab intensity values. The advantages of the MSP are to reduce the risk of alteration by heating at lower temperatures, to take account of the grain size effects and to avoid the need for anisotropy correction. Moreover, Tema et al. (2016) argued that the cooling rate correction was not required for this protocol. Fabian and Leonhardt (2010) expanded the original Dekkers and Böhnel (2006) protocol by correcting for the effects of NRM fraction and domain state. Four heating-cooling steps were performed at each Hlab: (1) +Hlab applied parallel to the NRM during heating and cooling; (2) Hlab applied antiparallel to the NRM during heating and cooling; (3) +Hlab applied parallel to the NRM only during cooling; (4) same as (1) to test for alteration. The magnetizations after each step were labelled m1, m2, m3 and m4, respectively. The first two steps provided the fraction of demagnetized NRM. The third step corrected for domain state and also served as an alteration check (Schnepp et al., 2016). We selected 25 fragments with no or a slight secondary component. From these, 66 specimens were embedded in plaster to make 2.2 cm3 cubic specimens. Heating experiments used a MMTD furnace at 400 °C or 450 °C with ten fields between 25 mT and 87 mT. Magnetizations were measured with a three-axis, 2G Enterprises, cryogenic magnetometer. Before each measurement, nine specimens were partly demagnetized with a 10 mT peak alternating field to remove the weak low-temperature magnetization component. The adjustment of the position of the specimens was achieved by home-made sample holders. The pTRM acquired in Hlab has to be parallel to the archaeological TRM. As the specimens were anisotropic, Hlab was not applied parallel to the NRM but in a direction determined using the TRM anisotropy results of the Thellier-Thellier protocol.

5. Archaeointensity results 5.1. Thellier-Thellier protocol Fig. 3a shows two examples of accepted specimens with linear NRM-TRM diagrams, positive pTRM-checks and no deviation of the NRM direction towards the laboratory field direction. Intensities were computed with temperature steps that displayed no evidence for alteration using the following criteria: number of temperature steps higher than 4, NRM fraction (f) higher than 0.3, quality factor (q) higher than 5, maximum angular deviation (MAD) lower than 5°, deviation angle (DANG) lower than 5°, ratio of the standard error of the slope to the absolute value of the slope (b) lower than 0.06 and difference ratio sum (DRATS) lower than 10. Fig. 3b shows rejected specimens based on deviation of the magnetization direction towards the laboratory field direction and on concave-up behaviour on the NRM-TRM diagram. For specimens possessing secondary TRMs, the NRM-TRM diagrams were corrected by subtracting the remaining NRM at the junction temperature between the two components from the NRM at each previous heating step (Hervé et al., 2013b). Archaeointensity determinations on these specimens required additional selection criteria (Fig. 3c). First, the samples must show no evidence for mineralogical changes, i.e. negative pTRM-checks, until at least the temperature between the two components determined in the thermal demagnetization. Second, the intensity value must be computed on the same temperature interval as the magnetization component isolated from thermal demagnetization. Specimens from 24 pots fulfilled these criteria. At the feature scale, pots with one or two magnetization components provided consistent archaeointensity values (Table 1S, Supplementary Material). 306 specimens from 154 pots were accepted for an acceptance rate of 56%. Table 1S (Supplementary Material) lists the data from all accepted specimens. The use of minispecimens raises the ques-

149

tion of accuracy, especially to know if the TRM anisotropy was adequately corrected. They exhibited high corrected degrees of TRM anisotropy (P’) generally between 1.1 and 1.8. Some pots had P’ up to 2.40. The burnt daub fragments of obla246 had P’ less than 1.3. During the manufacturing process, the clay paste was less deformed for this material than for ceramics, which resulted in the lower anisotropy. The consistency of the intensity results at the pot level with and without the ATRM correction allowed us to check the reliability of the ATRM tensor determinations (Fig. 3d). Whereas the raw values can vary by a factor of two, the ATRM correction generally reduced the dispersion by 10–95% (Table 1S, columns Fraw and FATRM). As the TRM difference between the fast and slow cooling was always higher than the alteration check, a cooling-rate correction was applied to all specimens (Gomez-Paccard et al., 2006). The correction factors varied between 5% and 21.5%, which are higher in comparison to literature values (e.g. Genevey et al., 2008) (Fig. 3e). This difference can be explained by a faster cooling for the 100 mg minispecimens than for the 10 cm3 specimens that we commonly used. For 75% of the specimens, the cooling rate correction further reduced the within-pot dispersion up to 90% (Table 1S). 5.2. Triaxe protocol Triaxe archaeointensity measurements serve as a check against those determined from Thellier-Thellier. As developed by Le Goff and Gallet (2004), the Triaxe protocol has an advantage over the classical Thellier-Thellier technique in that the former accounts for the anisotropy and cooling rate effects on the TRM acquisition (Genevey et al., 2009, 2016). The basic principles of the Triaxe protocol mimic those of Thellier-Thellier by comparing the NRM lost with a TRM acquired in a known laboratory field. Intensity experiments were conducted between T1, fixed to 150 °C, and T2 chosen close to 500 °C. At each running temperature Ti between T1 and T2, the ratio of the TRM and NRM fractions demagnetized between T1 and Ti, multiplied by the intensity of the laboratory field, gives R’(Ti) (Le Goff and Gallet, 2004). In case of a secondary magnetization component, the temperature T1 can be adjusted to T’1, so that the intensity values at the specimen level are only derived from the primary univectorial NRM. The R’(Ti) data are averaged over the T’1-T2 temperature range. The same set of quality criteria as Genevey et al. (2009, 2016) was used to select the Triaxe data. In particular, the R’(Ti) data must be constant between T’1 and T2, to within 10%. 22 specimens from 11 fragments respected this criterion. Fig. 4a shows examples of accepted data at the specimen level from one pot (grgy6652d) and one daub fragment (obla246e). The archaeointensity values were very consistent at the pot level, except for pot trha702b1, which also showed a higher dispersion in the Thellier-Thellier results (Table 1S, Supplementary Material). Without the cooling rate correction, archaeointensity values determined using the Thellier-Thellier protocol were overestimated compared to the Triaxe determinations (Fig. 4b). After correction, both methods gave consistent results within 5% uncertainty. This provided confidence that the Thellier-Thellier archaeointensities were adequately corrected for the cooling rate effect. 5.3. MSP-DSC protocol The MSP-DSC experiments were interpreted using MSPTool software (https://github.com/leonro/MSPTool). The first step was to evaluate the results using simple selection criteria. The fraction of NRM used for the archaeointensity calculation must range between 25 and 75%. Almost all samples had NRM fractions

150

G. Hervé et al. / Physics of the Earth and Planetary Interiors 270 (2017) 143–156

a = 0.2, 0.5 or 0.8 were similar (red, black and green lines,

Fig. 4. (a) Triaxe archaeointensity results of a burnt daub fragment (obla246e) and a pot fragment (grgy6652d). Each curve represents the intensity data obtained for one specimen over the temperature interval used for the archaeointensity calculation. (b) Comparison of Thellier-Thellier and Triaxe archaeointensity values with and without the cooling rate (CR) correction.

between 30 and 60% (Fig. 1S b and f, Supplementary Material). Specimens were also selected according to the NRM and TRM differences. The NRM difference was the angle between the full NRM and the remaining NRM directions. The TRM difference was the angle between the acquired pTRM and the remaining NRM directions. Both the TRM and NRM differences should be lower than 10°, otherwise the specimens were rejected (N = 12 in our study) (Fig. 1S c-d and g-h, Supplementary Material). Deviations of more than 10° can be explained by un- or mis-corrected anisotropy or by insufficient removal of the secondary magnetization component. The third commonly-used selection criteria was the relative alteration error ealt calculated from the difference between m1 and m4 (e.g. Tema et al., 2016). Here, ealt was lower than 5%, except for one specimen. On that one, we used two more criteria based on the third heating-cooling step (Schnepp and Brüggler, 2016; Schnepp et al., 2016). The intensity of a pTRM should be lower when Hlab was applied only during cooling (third step) than when it was applied both during heating and cooling (first and fourth steps). We discarded five other specimens having m3 higher than m1 and/or m4. In the end, 17, 11, 6 and 10 specimens were accepted for obla246, grgy665, grgy661-1/2 and grgy661-3/4/6. The number of field steps with accepted specimens was respectively 10, 8, 5 and 6. We considered that the number of accepted specimens should be higher than 10, which led us to discard grgy661-1/2 and grgy661-3/4/6. In the case of obla246 and grgy665, the fraction and domain-state corrections reduced the scatter (Fig. 5a–c and Fig. 1S a and e, Supplementary Material). The line fit to the data was fixed at 1 for zero-field. Following Fabian and Leonhardt (2010), the domain-state proxy (a-parameter) was set to 0.5. The a value did not play a significant role as the linear fits using

Fig. 5a-c). The reliability of the linear regression was tested with a jackknife technique, by randomly removing 20% of the collection (Schnepp and Brüggler, 2016; Schnepp et al., 2016). The intensity result is the median value of the jackknife distribution and the mean of the lower and upper quantiles as 2r uncertainties. In the case of the obla246 set of burnt daub fragments, the unimodal distribution of archaeointensity values highlights the reliability of the MSP-DSC result (Fig. 5b). The archaeointensity (65.0 ± 1.7 mT, standard error at 2r) was within 2r uncertainty limits of the mean Thellier-Thellier and Triaxe determination corrected for cooling rate effect (69.2 ± 5.2 mT, standard error at 2r). This supports the observation of Tema et al. (2016) that the cooling rate correction was not necessary for the MSP-DSC protocol; this point however requires further investigation. Site grgy665 had a bimodal distribution indicating a high scatter between specimens (Fig. 5d). The comparison of obla246 and grgy665 showed that the jackknife resampling distribution was a better indicator of the scatter than the correlation coefficient (R2). This is related to the uncertainty generated by the domainstate correction on each specimen. The archaeointensity values determined using Thellier-Thellier and Triaxe methods did not present a non-ideal bimodal distribution, which led us to reject the MSP-DSC results from grgy665. Besides the MSP-DSC intensity estimation would be on average 20% lower than those of Thellier-Thellier and Triaxe protocols. We think that undetected alteration, an under-corrected TRM anisotropy and/or slightly different times of TRM acquisition between pots, played a role in the scatter and the underestimation. For these reasons, a set of anisotropic ceramics from a pit appears to be an inappropriate archaeological material for the MSP protocol; when used, it is necessary to insure the laboratory field is well aligned parallel to the anisotropy corrected NRM direction of each specimen. The MSP method appears better suited for burnt features being less anisotropic material with a single TRM acquisition (Schnepp and Brüggler, 2016; Schnepp et al., 2016). The MSP protocol gave reliable results for homogeneous sets of baked clays fragments, for example brick or burnt daub fragments originating from a single building (e.g., obla246 feature).

6. Western Europe intensity secular variation from 1500 BCE to 200 CE 6.1. New archaeointensity data The individual determinations of the Thellier-Thellier and Triaxe protocols were averaged at the pot level. All burnt daub fragments from pit obla246 were treated as a single set. The final archaeointensity value of this pit was the average of the ThellierThellier/Triaxe and MSP-DSC averages weighted by the inverse of the standard error (Coe et al., 1978). The standard deviation at the pot level was generally lower than 5% or 5 mT. We rejected two pots with standard deviations >10%. Fig. 6a plots all the other archaeointensities results, which range between 40 mT and 95 mT. Archaeointensities were consistent within 10% between pots of the same feature at most periods (Table 1S, Supplementary material). The consistency was better for graves than pits and wells because mixing of domestic-used pots of slightly different ages was more likely to occur in domestic features. The dispersion between pots was also higher at certain periods, particularly in the Ha C period from 800/780 to 620 BCE and to lesser extent around 400 BCE. For example, the six accepted pots from well 390 of Aschheim – XXXLutz (asxx390) had archaeointensities from 40 to 87 mT. The fact that a high degree of scatter was observed in

G. Hervé et al. / Physics of the Earth and Planetary Interiors 270 (2017) 143–156

151

Fig. 5. MSP-DSC results for (a-b) obla246 and (c-d) grgy665 sites from Oberottmarshausen – Kiesgrube Lauter and Grünwald – Gymnasium. (a,c) Intensity difference ratios plotted versus applied field. Each colour refers to a different burnt daub fragment or sherd. The domain-state correction (DSC) data are shown with one sigma uncertainties. The red, black and green lines are the linear fit lines with the a-parameter set to 0.2, 0.5 and 0.8. Archaeointensity values are the median of the jackknife distribution and the mean of the lower and upper quantiles as 2r uncertainties. (b,d) Archaeointensity distributions obtained from the jackknife resampling technique. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

contemporaneous graves (e.g. asou48, asou51 and inma57) suggests that secular variation was faster during this period. Faster secular variation was also indicated by the archaeointensity results from kilns and fireplaces, where the samples acquired their thermoremanent magnetization at the same time. The Late Bronze Age kiln 900 BCE of Bevaix (Switzerland, Kovacheva et al., 2009 database) and fireplaces of Corent (France, Hervé et al., 2016) have archaeointensity values close to 50–60 mT, whereas those of the Early Iron Age 700-600 BCE are around 80 mT (France, Hervé et al., 2011, 2013b). In the present case, as we do not know precisely the age range of the individual pots within the age interval of the feature, we averaged the archaeointensity values at the feature level, therefore yielding larger experimental uncertainties up to 17 mT. The 35 new average data are listed in the Table 1 and plotted on Fig. 6b, relocated at 46°N-7°E using the Virtual Axial Dipole Moment (VADM) correction. This reference location was chosen as it was central to the central-western Europe database between 1400 BCE and 200 CE (Fig. 2S). Our study triples the number of reference data in Western Europe for the time interval from 1400 to 400 BCE. They are consistent with published results obtained within 600 km from the reference location (Aidona et al., 2006; Chauvin et al., 2000; Gallet et al., 2009a; Hervé et al., 2011, 2013b, 2016; Hill et al., 2007; Kapper et al., 2015; Kovacheva et al., 2004). According to Hervé and Lanos (2017), the Bevaix kiln should be dated at [970; 800] BCE using dendrochronology and archaeological artefacts. It is worth pointing out that the pottery set CES-BUR of Kapper et al. (2015), dated in the 8th and 7th centuries BCE, has a standard deviation of 15.5 mT similar to our contemporaneous data from Aschheim and Ingolstadt.

6.2. Bayesian secular variation curve The mean secular variation curve with its 95% confidence envelope was computed using a Bayesian framework implemented in Rencurve software (Lanos, 2004; Schnepp et al., 2015; Hervé and Lanos, 2017). The secular variation curve is given as a smooth continuous curve obtained by averaging cubic splines. The inversion process minimizes the misfit of each data with respect to the curve by exploring the multi-dimensional space of probability densities using Monte Carlo Markov Chains (MCMC). It takes into account all prior information on the selected reference data, including the uncertainties on the average intensity and its age. The prior probability for the uncertainty on age is uniform for data dated by archaeological techniques and dendrochronology, or irregular, resulting from the calibration process, for radiocarbon ages. The modelling takes into account stratigraphical relationships between features: oble13b is older than oble16 and oble17, grgy661-3/4/6 is older than grgy661-1/2). The curve in Fig. 6b was computed from 100,000 MCMC iterations (see also, Table 2S, Supplementary Material). The curve shows low variability in intensity around 50–55 mT from 1400 to 1200 BCE. The intensity then increased to 65 mT at 1000-900 BCE, before decreasing to 50 mT, close to the present day value, around 800 BCE. There was a subsequent increase up to 90 mT at 600 BCE. The geomagnetic field strength remained high during the 6th century BCE, before dropping between 450 and 350 BCE. Intensity secular variation is poorly constrained from the 4th to 2nd century BCE. It seems possible that a minimum occurred around 200 BCE but new data are required to confirm its existence. In Spain, Osete et al. (2016) noted lower intensity values in the 2nd

152

G. Hervé et al. / Physics of the Earth and Planetary Interiors 270 (2017) 143–156

Fig. 6. (a) Archaeointensity results at the pot level from this study. The 43 accepted pots from Grünwald – Gymnasium are not plotted by a coloured symbol like those of the other features. For clarity, the light grey box represents their range in intensity. All values were corrected for TRM anisotropy and cooling rate and relocated to Munich using the VADM assumption. (b) Secular variation of geomagnetic field strength in Western Europe from 1400 BCE to 200 CE. All data are relocated to 46°N, 7°E using the VADM assumption. The new data are plotted in colour and published data from within 600 km of the reference location in black. The new Bayesian secular variation is plotted in grey with its 95% confidence envelope.

than in the 3rd century BCE. Finally, the geomagnetic field intensity was fairly stable between the 1st century BCE and the 4th century CE (i.e. during Roman times). 7. Discussion 7.1. Secular variation rate The double oscillation of the intensity variations observed between 1100 and 300 BCE provides further evidence that the field can change by 50% over one or two centuries (e.g. Ben-Yosef et al., 2017; Shaar et al., 2016). The faster rate is estimated at 0.25 mT/yr at circa 700 BCE when intensity sharply increased. We note that this high rate of change is associated with higher standard deviations on the mean archaeointensity values around that period. This rate is three times higher than the average variation rate of 0.07– 0.1 mT/yr, in Western Europe over the last two millennia but similar to the rate at 1000 CE (Genevey et al., 2016; Gómez-Paccard et al., 2016). 0.25 mT/yr is less than half of the maximum possible

value (0.62 mT/yr) predicted by a magnetohydrodynamic toroidal flow regime at the core-mantle boundary (Livermore et al., 2014). In the Middle East, at the beginning of the 10th and 8th century BCE, Shaar et al. (2011, 2016) identified two geomagnetic spikes with VADMs higher than 16
1022 A
m2. There is no evidence of such sharp intensity increases at the sub-centennial scale in Western Europe considering the intensity data both at the feature (Fig. 6b) and pot levels (Fig. 6a). 7.2. Link with directional variations The directional secular variation was also fast in Western Europe between 1200 and 500 BCE, with an abrupt change in direction (cusp) around 800 BCE when declination values were above 30° (Gallet et al., 2002; Hervé et al., 2013a). When comparing intensity and directional secular variation curves from Western Europe over the two CE millennia, one of the most intriguing observations is that times of directional cusps coincide with times of field intensity maxima. These features were referred to as archaeomagnetic jerks

G. Hervé et al. / Physics of the Earth and Planetary Interiors 270 (2017) 143–156

(Gallet et al., 2003). For earlier periods, an archaeomagnetic jerk was proposed at circa 800 BCE (Gallet et al., 2003; Hervé et al., 2013b), although the database for Western Europe was rather poor and only an intensity decrease after 800 BCE was supported by the data. However, our new intensity data showing an intensity minimum challenge the existence of an archaeomagnetic jerk at 800 BCE (Fig. 6b). Even accounting for the uncertainties on the secular variation curves, it seems more probable that an age of 800 BCE for the directional cusp (Fig. 10a in Hervé et al., 2013a) does not coincide with an intensity maximum. This discrepancy may arise from an insufficient number of directional reference data, which would result in an over-smoothing of the corresponding secular variation curve. Furthermore, the fact that we observe two intensity maxima at 950 and 600 BCE could indicate the occurrence of two directional cusps—something the present database does not allow to test. The acquisition of new directional and/or full vector data from the first half of the first millennium BCE would be particularly important to decipher a possible relationship between directional and intensity variations. 7.3. Comparison with geomagnetic models Fig. 7a compares the Western Europe intensity curve with global geomagnetic models CALS10k.2 (Constable et al., 2016), pfm9k.1b (Nilsson et al., 2014) and SHA.DIF.14k (Pavón-Carrasco et al., 2014a) predicted for 46°N, 7°E. All models were developed by data inversion using spherical harmonic analysis in space and cubic B-splines in time. The input data of the pfm9k.1b and CALS10k.2 models come from archaeomagnetic, volcanic and sedimentary sources. Model pfm9k.1b identifies a single maximum during the sixth and fifth centuries BCE, while CALS10k.2 indicates a weak second maximum in the 10th century BCE. Both models smooth secular variation and do not reproduce the intensity minimum at 800 BCE. Inhomogeneous quality in archaeointensity and inhomogeneous geographical site distribution probably explains the inconsistencies, as well as the smoothing induced from the relative sedimentary records. SHA.DIF.14k built using only archaeomagnetic and volcanic data is essentially valid in the northern hemisphere as the data coverage is sparse in the southern hemisphere. This model predicts the two maxima in geomagnetic field strength with slight shifts in time and amplitude. It is worth pointing out that only 20% of the data used in our Western Europe curve were included in SHA.DIF.14k.

153

7.4. A strong departure from an axial dipole field in 1000-700 BCE SHA.DIF.14k predicts an increase in the axial dipole moment from 9
1022 A
m2 in 1400 BCE to 11
1022 A
m2 at circa 500 BCE before slightly decreasing to 10
1022 A
m2 at the transition from BCE to CE (Fig. 7b). This prediction is consistent in time and in amplitude with the GMAG.9k global VADM reconstruction from temporally and spatially weighted archaeomagnetic and volcanic data (Usoskin et al., 2016). The maximum intensity in Western Europe is approximately synchronous with the maximum global axial dipole moment. However, the Western Europe curve clearly departs from predictions provided by the SHA.DIF.14k and GMAG.9k models between 1000 and 400 BCE. This suggests that the VADM minimum 8
1022 A
m2 at circa 800 BCE and the VADM maximum 14
1022 A
m2 at 600-500 BCE are related to stronger non-axial dipole fields. According to model SHA.DIF.14k, the dipole axis was tilted 10° toward 20-40°E longitude between 1000 and 600 BCE, which suggests a stronger equatorial dipole during this period (Fig. 8a). Virtual Geomagnetic Poles (VGP) at 800 BCE computed from Western Europe archaeomagnetic directions deviated 20° from the axis toward 80°E (Hervé et al., 2013b) (Fig. 8a). Although this discrepancy may partly arise from the limited resolution of the geomagnetic field model, it seems more probable that, besides the dipole tilt, higher spherical harmonic degrees are involved in the rapid secular variation. The non-dipole field was investigated by comparing the Western Europe VADMs with other regional records of the geomagnetic field strength (Fig. 8). We compiled all intensity data with alteration checks from five regions close to Europe (data in Table 3S): Canary Islands and West Africa (Fig. 8b), Middle East and Caucasus (Fig. 8c) and Central Asia (Fig. 8d). The Middle East data set gathers all data from Israel, Jordan, Syria, Turkey and Cyprus. VADMs from all these areas are consistent within uncertainty from 700-600 BCE to at least 100 CE, whereas significant differences are observed for earlier centuries. As in Western Europe, the intensity profile presents two maxima in the Canary Islands, the Middle East, the Caucasus and, to a lesser extent, Central Asia. VADMs from the Canary Islands (Kissel et al., 2015) and West Africa (Mitra et al., 2013) are consistent both in time and in amplitude with those of Western Europe (Fig. 8b). The scatter of the Canary Islands data between 800 and 400 BCE is due to the Hallstattian radiocarbon plateau effect but

Fig. 7. (a) Comparison of the Western Europe secular variation curve with the prediction at 46°N, 7°E of the geomagnetic global models CALS10k.2 (green), pfm9k.1b (orange) and SHA.DIF.14k (blue). (b) Comparison of the Virtual Axial Dipole Moment (VADM) variation curve inferred from the Western Europe curve with the averaged VADM global curve (pink, Usoskin et al., 2016) and the axial dipole moment (ADM) predicted by SHA.DIF.14k (blue, Pavón-Carrasco et al., 2014a). (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

154

G. Hervé et al. / Physics of the Earth and Planetary Interiors 270 (2017) 143–156

Fig. 8. (a) Mean VGP curve for Western Europe derived from the directional data (red, Hervé et al., 2013b) and the curve predicted by SHA.DIF.14 k (blue) of the position of the north geomagnetic pole between 1400 BCE and 100 CE. Also plotted in (a) are the sampling locations of the VADM data from (b) the Canary Islands and West Africa, (c) the Middle East and Georgia and (d) Central Asia, which are compared to those from Western Europe. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

the data suggest a minimum at the beginning of this period followed by a fast increase up to the 5th century BCE. In the Middle East and the Caucasus, the intensity maximum at 700 BCE has the same amplitude 14–16
1022 A
m2 but is apparently shifted one century earlier than in Western Europe. The maximum in the 10th century BCE is contemporaneous between the two areas taking into account the age uncertainties, but the amplitude is 14–16
1022 A
m2 in Middle East yet much lower 10–12
1022 A
m2 in Western Europe and the Canary Islands. The intensity minimum circa 800 BCE is observed in the Georgian data with VADM values close to those of Western Europe. The Middle East data do not show evidence of the intensity minimum. Variations from 1400 to 1000 BCE are similar in the three regional data sets, but with higher intensities in the Middle East and lower intensities in the Caucasus in comparison to Western Europe (Fig. 8c). The data from each of these areas were published by different researchers, acquired from different archaeomagnetic materials and sometimes used different laboratory protocols. It can therefore be reasonably concluded that the intensity differences, especially between Western Europe and the Middle East, are not dominated by experimental error and have rather a geomagnetic origin. Contrary to other data sets, VADMs from Central Asia (Uzbekistan and Turkmenistan) were not corrected for the cooling rate effect and may be overestimated by 5–10% (Genevey et al.,

2008). They are consistent with the Western Europe data between 1400 and 1000 BCE and after the 7th century BCE (Fig. 8d). Differences are observed between 1000 and 700 BCE with no clear intensity minimum in Central Asia and higher VADMs (14
1022 A
m2) than in Western Europe. The 20–30% difference between the VADMs seems too high to be solely explained by the absence of the cooling rate correction. To explain the high 14
1022 A
m2 intensities and the fast variation between 1050 and 700 BCE in the Middle East, Caucasus and Central Asia, Shaar et al. (2016, 2017) proposed the existence of a regional geomagnetic anomaly, known as the ‘‘Levantine Iron Age Anomaly” (LIAA). Our new Western Europe data compilation fits well with this hypothesis, as the greatest differences in VADM between Western Europe, the Middle East and Central Asia lie within the LIAA interval between 1000 and 700 BCE. Between 1400 and 1000 BCE, the Western Europe VADMs are consistent with those of Central Asia but deviate from those of the Middle East and the Caucasus. An explanation for this could be that the LIAA started earlier, at least from 1400 BCE, and then grew toward Central Asia circa 1000 BCE. However, the lack of intensity data farther north of the Caucasus and of Central Asia prohibits investigation whether the proposed LIAA could be related to the growth, the migration of a flux lobe from south to north, or to the coalescence of several lobes. Nevertheless, the strengthening of a major flux

G. Hervé et al. / Physics of the Earth and Planetary Interiors 270 (2017) 143–156

lobe at the beginning of the first millennium BCE may be consistent with the hypothesis of an eccentric dipole during this period (Gallet et al., 2009b). It is worth noting that the tilt of the dipole predicted by global models between 1000 and 800 BCE is toward the Middle East (Fig. 8a). In sum, a non-dipole field component that superimposes an increase of the axial dipole and a tilt of the dipole (or an eccentric dipole) could explain the geomagnetic intensity variations that occurred in Western Eurasia between 1000 and 500 BCE. At this stage, however, deciphering the respective contribution of these processes remains premature.

8. Conclusion Thirty-five new mean archaeointensity data were obtained from archaeological sites in Bavaria, Germany. The new data were dated between 1400 and 400 BCE and triple the database for Western Europe during this period. The reliability of the ThellierThellier archaeointensity values, especially regarding the anisotropy and cooling rate corrections, has been strengthened by Triaxe measurements. The new secular variation curve of Western Europe from 1400 BCE to the beginning of the CE period shows a double oscillation of the geomagnetic field strength with two maxima, one during the 10th century BCE reaching 70 mT and another at circa 600 BCE of 90 mT. The two are separated by a minimum of 50 mT near 800 BCE. This evolution resulted in fast secular variation rates up to 0.25 mT/yr circa 700 BCE. Our results confirm that the geomagnetic field strength in Europe had varied over a large range of amplitudes and rates over the past few millennia. One can question whether the configuration of the period between 1000 and 500 BCE is a typical feature of the geomagnetic field or whether it is somehow exceptional. Interestingly, secular variation in Western Europe during the first half of the first millennium BCE presents similarities with that having prevailed between 600 and 1000 CE, including two successive intensity maxima 80– 90 mT, rapid variation rates of the same order of magnitude (Genevey et al., 2016; Gómez-Paccard et al., 2016), and directions that deviated 11° away from the GAD direction (Bucur, 1994). In this respect, the geomagnetic field configuration in Western Europe between 1000 and 500 BCE would not be exceptional. Besides improving our understanding of the geomagnetic field, the better definition of secular variation will improve the archaeomagnetic dating technique during the Late Bronze and Iron Ages in Western Europe. This is especially true for the Hallstatt C period between 800 and 600 BCE, where intensity alone can now provide precise dating constraints thanks to high rates of secular variation.

Acknowledgements Funding was provided to GH by Deutsche Forschungsgemeinschaft project HE7343/1-1. This is IPGP contribution no. 3863. Jim Pincini and Sandra Ostner are kindly acknowledged for their contribution to the laboratory experiments. We are grateful to Philippe Lanos for his help in the computation of the Bayesian curve, to Maxime Le Goff for his laboratory support at IPG Paris and to Roman Leonhard for providing the MSPTool python software. The article also benefits from fruitful discussions with Florian Lhuillier and from the comments of two anonymous reviewers.

Appendix A. Supplementary data Supplementary data associated with this article can be found, in the online version, at http://dx.doi.org/10.1016/j.pepi.2017.07.002.

155

References Aidona, E., Scholger, R., Mauritsch, H.J., Schnepp, E., Klemm, S., 2006. Spatial distribution of archaeomagnetic vectors within archaeological samples from Eisenerz (Austria). Geophys. J. Int. 166 (1), 46–58. Aubert, J., Finlay, C.C., Fournier, A., 2013. Bottom-up control of geomagnetic secular variation by the Earth’s inner core. Nature 502 (7470), 219–223. Bagley, J.M., Eggl, C., Neumann, D., Schumann, R., 2010. Die späthallstatt-/ frühlatènezeitliche Siedlung an der Haffstraße in München-Trudering. Bericht Der Bayerischen Bodendenkmalpflege 51, 67–125. Ben-Yosef, E., Millman, M., Shaar, R., Tauxe, L., Lipschits, O., 2017. Six centuries of geomagnetic intensity variations recorded by royal Judean stamped jar handles. PNAS 114 (9), 2160–2165. Brown, M.C., Donadini, F., Korte, M., Nilsson, A., Korhonen, K., Lodge, A., Lengyel, S. N., Constable, C.G., 2015. GEOMAGIA50.v3: 1. general structure and modifications to the archeological and volcanic database. Earth Planets Space 67, 83. http://dx.doi.org/10.1186/s40623-015-0232-0. Bucur, I., 1994. The direction of the terrestrial magnetic field in France, during the last 21 centuries. Recent progress. Phys. Earth Planet. Inter. 87 (1–2), 95–109. Büttner, A., Linke, R., Wirth, S., 2006. Reiche Gräber der Bronze- und Urnenfelderzeit und ein Depotfund mit Wagenbronzen von Oberottmarshausen. Das archäologische Jahr in Bayern, 53–56. Chauvin, A., Garcia, Y., Lanos, P., Laubenheimer, F., 2000. Paleointensity of the geomagnetic field recovered on archaeomagnetic sites from France. Phys. Earth Planet. Inter. 120, 111–136. Coe, R., Grommé, S., Mankinen, E., 1978. Geomagnetic paleointensities from radiocarbon-dated lava flows on Hawaii and the question of the Pacific nondipole low. J. Geophys. Res. B 4, 1740–1756. Constable, C., Korte, M., Panovska, S., 2016. Persistent high paleosecular variation activity in southern hemisphere for at least 10 000 years. Earth Planet. Sci. Lett. 453, 78–86. Dekkers, M.J., Böhnel, H.N., 2006. Reliable absolute palaeointensities independent of magnetic domain state. Earth Planet. Sci. Lett. 248, 508–517. Donadini, F., Korhonen, K., Riisager, P., Pesonen, L., 2006. Database for Holocene geomagnetic intensity information. Eos Trans. AGU 87 (14), 137. Dunlop, D., 2002. Theory and application of the Day plot (Mrs/Ms versus Hcr/Hc) 1. Theoretical curves and tests using titanomagnetite data. J. Geophys. Res 107 (B3). http://dx.doi.org/10.1029/2001JB000486. Egli, R., 2013. VARIFORC: an optimized protocol for calculating non-regular firstorder reversal curve (FORC) diagrams. Global Planet. Change 110, 302–320. Fabian, K., Leonhardt, R., 2010. Multiple-specimen absolute paleointensity determination : an optimal protocol including pTRM normalization, domainstate correction, and alteration test. Earth Planet. Sci. Lett. 297, 84–94. Gallet, Y., Genevey, A., Le Goff, M., 2002. Three millennia of directional variation of the Earth’s magnetic field in Western Europe as revealed by archaeological artefacts. Phys. Earth planet. Int. 131, 81–89. Gallet, Y., Genevey, A., Courtillot, V., 2003. On the possible occurrence of ‘‘archaeomagnetic jerks” in the geomagnetic field over the past three millennia. Earth Planet. Sci. Lett. 214 (1–2), 237–242. Gallet, Y., Genevey, A., Le Goff, M., Warmé, N., Gran-Aymerich, J., Lefèvre, A., 2009a. On the use of archeology in geomagnetism, and vice-versa: recent developments in archeomagnetism. C.R. Phys. 10 (7), 630–648. Gallet, Y., Hulot, G., Chulliat, A., Genevey, A., 2009b. Geomagnetic field hemispheric asymmetry and archeomagnetic jerks. Earth Planet. Sci. Lett. 284 (1–2), 179– 186. Genevey, A., Gallet, Y., Constable, C.G., Korte, M., Hulot, G., 2008. ArcheoInt: an upgraded compilation of geomagnetic field intensity data for the past ten millennia and its application to the recovery of the past dipole moment. Geochem. Geophys. Geosyst. 9 (4). http://dx.doi.org/10.1029/2007GC001881. Genevey, A., Gallet, Y., Rosen, J., Le Goff, M., 2009. Evidence for rapid geomagnetic field intensity variations in Western Europe over the past 800 years from new French archeointensity data. Earth Planet. Sci. Lett. 284 (1–2), 132–143. Genevey, A., Gallet, Y., Jesset, S., Thébault, E., Bouillon, J., Lefèvre, A., Le Goff, M., 2016. New archeointensity data from French Early Medieval pottery production (6th–10th century AD). Tracing 1500 years of geomagnetic field intensity variations in Western Europe. Phys. Earth Planet. Inter. 257, 205–219. Gómez-Paccard, M., Chauvin, A., Lanos, P., Thiriot, J., Jiménez-Castillo, P., 2006. Archeomagnetic study of seven contemporaneous kilns from Murcia (Spain). Phys. Earth Planet. Inter. 157 (1–2), 16–32. Gómez-Paccard, M., Osete, M.L., Chauvin, A., Pavón-Carrasco, F.J., Pérez-Asensio, M., Jiménez, P., Lanos, P., 2016. New constraints on the most significant paleointensity change in Western Europe over the last two millennia. A nondipolar origin? Earth Planet. Sci. Lett. 454, 55–64. Hervé, G., Schnepp, E., Chauvin, A., Lanos, P., Nowaczyk, N., 2011. Archaeomagnetic results on three Early Iron Age salt-kilns from Moyenvic (France). Geophys. J. Int. 185 (1), 144–156. Hervé, G., Chauvin, A., Lanos, P., 2013a. Geomagnetic field variations in Western Europe from 1500BC to 200AD. Part I: directional secular variation curve. Phys. Earth Planet. Inter. 218, 1–13. Hervé, G., Chauvin, A., Lanos, P., 2013b. Geomagnetic field variations in Western Europe from 1500 BC to 200 AD. Part II: new intensity secular variation curve. Phys. Earth Planet. Inter. 218, 51–65. Hervé, G., Chauvin, A., Milcent, P.-Y., Tramon, A., 2016. Archaeointensity study of five Late Bronze Age fireplaces from Corent. J. Arch. Sci.: Rep. 7, 414–419.

156

G. Hervé et al. / Physics of the Earth and Planetary Interiors 270 (2017) 143–156

Hervé, G., Lanos, P., 2017. Improvements of archaeomagnetic dating in Western Europe from the Late Bronze to the Late Iron Ages: an alternative to the problem of the Hallstattian radiocarbon plateau, accepted. Archaeometry. http://dx.doi. org/10.1111/arcm.12344. Herzig, F., 2008a. Dendroarchäologische Untersuchungen Aschheim-Ostumgehung 2008. Bayerisches Landesmat für Denkmalpflege, Referat BV – Restaurierung Archäologie und Dendrolabor Thierhaupten, pp. 3. Herzig, F., 2008b. Dendroarchäologische Untersuchungen Aschheim-XXX-Lutz 2008. Bayerisches Landesmat für Denkmalpflege, Referat BV – Restaurierung Archäologie und Dendrolabor Thierhaupten, pp. 9. Hill, M.J., Lanos, P., Chauvin, A., Vitali, D., Laubenheimer, F., 2007. An archaeomagnetic investigation of a Roman amphorae workshop in Albinia (Italy). Geophys. J. Int. 169 (2), 471–782. Kapper, K.L., Donadini, F., Hirt, A.M., 2015. Holocene archeointensities from mid European ceramics, slags, burned sediments and cherts. Phys. Earth Planet. Inter. 241, 21–36. Kissel, C., Laj, C., Rodriguez-Gonzalez, A., Perez-torrado, F., Carracedo, J.C., Wandres, C., 2015. Holocene geomagnetic field intensity variations: contribution from the low latitude Canary Islands site. Earth Planet. Sci. Lett. 430, 178–190. Kovacheva, M., Boyadziev, Y., Kostadinova-Avramova, M., Jordanova, N., 2009. Updated archeomagnetic data set of the past eight millenia from the Sofia laboratory, Bulgaria. Geochem. Geophys. Geosyst. 10 (5). http://dx.doi.org/ 10.1029/2008GC002347. Kovacheva, M., Hedley, I., Jordanova, N., Kostadinova, M., Gigov, V., 2004. Archaeomagnetic dating of archaeological sites from Switzerland and Bulgaria. J. Archaeol. Sci. 31 (10), 1463–1479. Lanos, P., 2004. Bayesian inference of calibration curves: application to archaeomagnetism. In: Buck, C.E., Millard, A.R. (Eds.), Tools for Constructing Chronologies: Crossing Disciplinary Boundaries. Springer, London, pp. 43–82. Lecture No. Lanos, P., Philippe, A., 2017. Hierarchical Bayesian modeling for combining dates in archaeological context. Journal de La Société Française de Statistique. Le Goff, M., Gallet, Y., 2004. A new three-axis vibrating sample magnetometer for continuous high-temperature magnetization measurements : applications to paleo- and archeo-intensity determinations. Earth Planet. Sci. Lett. 229, 31–43. Licht, A., Hulot, G., Gallet, Y., Thébault, E., 2013. Ensembles of low degree archeomagnetic field models for the past three millennia. Phys. Earth Planet. Inter. 224, 38–67. Livermore, P.W., Fournier, A., Gallet, Y., 2014. Core-flow constraints on extreme archeomagnetic intensity changes. Earth Planet. Sci. Lett. 387, 145–156. Malcher, G. J., Weinig, J. G., 1995. Der urnenfelder- und hallstattzeitliche Friedhof von Mailing, ‘‘Schindergrubäcker”. Das Archäologische Jahr in Bayern, 63. Metzner-Nebelsick, C. (dir.), 2016. Einblicke in die Vorgeschichte Grünwalds, Ergebnisse eines archäologischen Pilotprojekts, Institut für Vor- und Frügeschichtliche Archäologie und Provinzialrömische Archäologie, München, pp. 97. Mitra, R., Tauxe, L., Keech McIntosh, S., 2013. Two thousand years of archeointensity from West Africa. Earth Planet. Sci. Lett. 364, 123–133. Müller-Karpe, H., 1959. Beiträge zur Chronologie der Urnenfelderzeit nördlich und südlich der Alpen, Römich-Germanischen Forschungen, Berlin, W. De Gruyter, pp. 334.

Nilsson, A., Holme, R., Korte, M., Suttie, N., Hill, M., 2014. Reconstructing holocene geomagnetic field variation: new methods, models and implications. Geophys. J. Int. 198 (1), 229–248. http://dx.doi.org/10.1093/gji/ggu120. Osete, M.-L., Chauvin, A., Catanzariti, G., Jimeno, A., Campuzano, S.A., BenitoBatanero, J.P., Tabernero-Galan, C., Roperch, P., 2016. New archaeomagnetic data recovered from the study of celtiberic remains from central Spain (Numantia and Ciadueña, 3rd-1st centuries BC). Implications on the fidelity of the Iberian paleointensity database. Phys. Earth Planet. Inter. 260, 74–86. Pavón-Carrasco, F.J., Osete, M.L., Torta, J.M., De Santis, A., 2014a. A geomagnetic field model for the Holocene based on archaeomagnetic and lava flow data. Earth Planet. Sci. Lett. 388, 98–109. Pavón-Carrasco, F.J., Gómez-Paccard, M., Hervé, G., Osete, M.L., Chauvin, A., 2014b. Intensity of the geomagnetic field in Europe for the last 3 ka: influence of data quality on geomagnetic field modeling. Geochem. Geophys. Geosyst. 15. doi. org/10.1002/2014GC005311. Pütz, A., 2008. Aschheim-Ostumgehung, Grabungsbericht Firma ARDI. Schnepp, E., Brüggler, M., 2016. Archaeomagnetic investigation of a Roman glass workshop in Goch-Asperden, Germany. J. Arch. Sci.: Rep. 10, 322–330. Schnepp, E., Obenaus, M., Lanos, P., 2015. Posterior archaeomagnetic dating: an example from the Early Medieval site Thunau am Kamp, Austria. J. Arch. Sci.: Rep. 2, 688–698. Schnepp, E., Leonhardt, R., Korte, M., Klett-Drechsel, J., 2016. Validity of archaeomagnetic field recording: an experimental pottery kiln at Coppengrave, Germany. Geophys. J. Int. 205, 622–635. Shaar, R., Ben-Yosef, E., Ron, H., Tauxe, L., Agnon, A., Kessel, R., 2011. Geomagnetic field intensity: how high can it get? How fast can it change? Constraints from Iron Age copper slag. Earth Planet. Sci. Lett. 301 (1–2), 297–306. Shaar, R., Tauxe, L., Ron, H., Ebert, Y., Zuckerman, S., Finkelstein, I., Agnon, A., 2016. Large geomagnetic field anomalies revealed in Bronze to Iron Age archeomagnetic data from Tel Megiddo and Tel Hazor, Israel. Earth Planet. Sci. Lett. 442, 173–185. Shaar, R., Tauxe, L., Goguitchaichvili, A., Devidze, M., Licheli, V., 2017. Further evidence of the Levantine Iron Age geomagnetic anomaly from Georgian pottery. Geophys. Res. Lett. 5. http://dx.doi.org/10.1002/2016GL071494. Sommer, S. (Ed.), 2006. Archäologie in Bayern: Fenster zur Vergangenheit, Gesellschaft für Archäologie in Bayern, Regensburg: Pustet Verlag, pp. 336. Tema, E., Ferrara, E., Camps, P., Conati Barbaro, C., Spatafora, S., Carvallo, C., Poidras, T., 2016. The Earth’s magnetic field in Italy during the Neolithic period: new data from the Early Neolithic site of Portonovo (Marche, Italy). Earth Planet Sci. Lett. 448, 49–61. Thellier, E., Thellier, O., 1959. Sur l’intensité du champ magnétique terrestre dans le passé historique et géologique. Annales de Géophysique 15, 285–376. Usoskin, I.G., Gallet, Y., Lopes, F., Kovaltsov, G.A., Hulot, G., 2016. Solar activity during the Holocene: the Hallstatt cycle and its consequence for grand minima and maxima. Astron. Astrophys. 587 (A150). http://dx.doi.org/10.1051/00046361/201527295. Winghart, S., 1984. Neue Ausgrabungen im urnenfelderzeitlichen Gräberfeld von München – Obermenzing, Landeshauptstadt München, Oberbayern. Das Archäologische Jahr in Bayern, 60–61. Zach, B., Knipping, M., Pütz, A., 2010. The well Ostumgehung 692, Hallstattzeitliche Pflanzenreste aus einem Brunnen bei Aschheim. Bericht der bayerischen Bodendenkmalpflege 51, 59–66.