Rapid Intensity Decrease During the Second Half of the First Millennium BCE in Central Asia and Global Implications

Recent paleomagnetic studies have shown that important short‐lived intensity fluctuations occurred during the first millennium BCE. However, the knowledge of the spatial and temporal extension of these features is still limited by the scarce availability of robust data. In this study we focus on the study of the intensity decrease that took place in Central Asia during the second half of the 1st millennium BCE after the high intensities that characterized the Levantine Iron Age Anomaly. Since previous archeointensities available for this period and region were obtained without accomplishing modern standards of quality, we present here new archeointensities that are derived from classical Thellier and Thellier experiments, including partial thermoremanent magnetization (pTRM) checks, thermoremanent magnetization (TRM) anisotropy and cooling rate corrections at the specimen level. The new 51 archeointensities, together with previous archeointensities, have been used to present a new local paleosecular variation curve for Central Asia. The results confirm the existence of an important geomagnetic field intensity decrease in South Uzbekistan from the 4th century BCE to the end of the 1st century BCE associated with rates of changes up to −15 μT/century. A critical analysis of the archeointensity global database indicates that this feature was present at continental scale, from Western Europe to Central Asia. However, this trend is not identified in other regions such as Japan or Mexico. Finally, the comparison with the dipole moment derived from recent global geomagnetic field reconstructions suggests a strong influence of non‐dipolar sources upon this continental intensity feature.

Since the middle of the 19th century, the spatial and temporal evolution of the main geomagnetic field is well known through the instrumental data. However, the only way to decipher past geomagnetic field changes before the era of direct measurements comes from the paleomagnetic study of the remanent magnetization acquired by the iron oxides present in certain rocks and baked archeological materials and retained over historical or geological time scales. This magnetization can be analyzed in the paleomagnetic laboratories to infer the direction and strength of the ancient geomagnetic field. Since past changes of the geomagnetic field observed at the Earth's surface are not equal throughout the world due to the different contributions of dipolar and non-dipolar sources at each location, in order to have a global view of the past geomagnetic field it is necessary to investigate its behavior at different regions.
Strong variations in geomagnetic field intensity occurred during the first millennium BCE. Extremely high paleointensities were first reported by Ben-Yosef et al. (2009) and Shaar et al. (2011) from the archeointensity study of Jordanian and Israelian copper slags dated around 1000 BCE. Later studies confirmed these anomalous high strength field values in the Levant and surrounding regions (e.g., Ben-Yosef et al., 2017;Ertepinar et al., 2012Ertepinar et al., , 2020Shaar et al., 2016Shaar et al., , 2017. Shaar et al. (2016) named this high intensity feature observed in the Levant from about 1050 BCE to 700 BCE as the Levantine Iron Age Anomaly (LIAA). Rapid changes in geomagnetic field strength during the first millennium BCE have also been reported in other parts of the world: Canary Islands (de Groot et al., 2015;Kissel et al., 2015), the Azores Archipelago (Di Chiara et al., 2014), Western and Central Europe (Hervé et al., 2013(Hervé et al., , 2017Molina-Cardín et al., 2018;Osete et al., 2020;Schnepp et al., 2020), the Mediterranean area (Béguin et al., 2019;Rivero-Montero, Gómez-Paccard, Kondopoulou, et al., 2021), eastern and central China (Cai et al., 2017(Cai et al., , 2020, Korea (Hong et al., 2013) and Texas, USA (Bourne et al., 2016). The youngest intensity maxima in the European region for the first millennium BCE occurred around 500 BCE (Rivero-Montero, Gómez-Paccard, Kondopoulou, et al., 2021). After 500 BCE, the strong field vanished and a subsequent time interval characterized by a decreasing trend of geomagnetic field strength was observed over the Mediterranean area (Rivero-Montero, Gómez-Paccard, Kondopoulou, et al., 2021). To further investigate if this behavior is also observed eastwards, we present here the archeointensity and rock-magnetic study of 70 pottery fragments from Southern Uzbekistan (Central Asia). The new data, together with a critical revision of previous results, provide an improved vision of the geomagnetic field intensity changes in the region between 600 BCE and 600 CE. The present study is the first step of a long-term project which aims to decipher geomagnetic field strength changes in Central Asia over the past few millennia and that will complement the recently published archeointensity results from Bukhara (Uzbekistan) dated between the end of the 16th century and the beginning of the 19th century (Troyano et al., 2021).

Archeological Contexts
The ceramic assemblages studied were collected in three archeological sites located in Bactria: Kurganzol (KGZ), Kampyr Tepe (KPT) and Termez (TRZ). This historical region is comprised by the current territories of southern Uzbekistan, Tajikistan, and northern Afghanistan (see Figure 1 and Supporting Information S1 for further details).

Archeological Sites and Dating Constraints
Archaeological data suggest that the three investigated sites were founded during or shortly after Alexander the Great's campaigns in ancient Bactria (late 4th century BCE). The relative chronologies proposed for each site are based on numismatics and the techno-typological analysis of the ceramic wares, tracing parallels with other well-dated coeval sites in Bactria. All 14 C and techno-typological dates available, with calibrated calendar years and probability percentages for 1σ and 2σ, are listed in Table S1 in Supporting Information S1. In the present study, the dating assigned to each pottery sample was based on the typological chronology and the 14 C data when available, and was specified by considering the presence of coins, the chronologies determined to each unit of the stratigraphic sequence and the chronologies proposed to the same ceramic prototypes and types in other regional archaeological sites.
The first site, Kurganzol, was excavated by the Institute of Fine Arts of the Academy of Sciences of Uzbekistan in 2008. It lies at 924 m a.s.l. (above sea level), in the southern foothills of the Hissar Range (current Uzbekistan), and on the southern bank of the Uchkul river, a tributary of the Surkhan Darya. Kurganzol consists of a circular fortress, about 35 m in diameter, surrounded by a thick wall made with mud bricks, with several defensive towers attached to its north-eastern side. The archaeological context and absolute dating (dendrochronology) suggest a date for the occupation sequence which should be established in the late 4th century BCE-beginning of the 3rd century BCE (Sverchkov, 2013; see Supporting Information S1 for details).
The second site, Kampyr Tepe, comprises a citadel and a dwelling area in the fluvial terrace, which included a river moat. The materials investigated here were recovered from different archaeological contexts: the citadel, the low city, and a pottery workshop (see Martínez Ferreras et al., 2016 and Supporting Information S1 for details). The pottery workshop includes two circular kilns that were excavated in the eastern sector outside the walls which have been dated by typological dating between late 4th BCE and 1st century BCE. Several small firing structures were also found in the eastern side of the citadel (Bolelov, 2001(Bolelov, , 2011. 14 C analysis carried out at the Radiocarbon Dating Laboratory of the University of Barcelona on charcoal from one of these firing structures from citadel dates its use between the late 4th and mid-2nd century BCE. Also, 14 C absolute dating has been performed on materials from the low city which dated from 2nd to first half of the 1st century BCE. Termez was settled close to the mouth of the Surkhan Darya, a tributary of the Amu Darya. It has systematically been excavated since the 1990s; one of the most active archaeological research teams from the last decades is the Uzbek-Spanish team IPAEB, led by the University of Barcelona and the Institute of Fine Arts of Uzbekistan. The city developed from the original fortress, the citadel, was built probably by the successors of Alexander in Seleucid times (first half of the 3rd century BCE) (Martínez Ferreras et al., 2019). It was extended during the Greco-Bactrian (mid-3rd to mid/late-2nd century BCE), Yuezhi (mid/late 2nd century and 1st century BCE) and Kushan periods (1st to mid-3rd centuries CE), becoming an important Buddhist center. The materials for this study were recovered from two different areas. The sherds dated in the Yuezhi/early Kushan nomadic period come from the sector AC located in the alluvial plain. The fortress of Tchingiz Tepe provides materials dated in the Kushan and Kushano-Sasanian sequence. Radiocarbon analysis has been performed at the Radiocarbon Dating Laboratory of the University of Barcelona and the National Center of Accelerators in Seville on vegetable charcoal and bone found in several strata of the two settlements. The results suggest an occupation between the late-3rd century BCE and the 1st century CE for the alluvial plain. The materials from Tchingiz Tepe here studied are dated between the 2nd and the 4th centuries CE (Martínez Ferreras et al., 2014 and Supporting Information S1 for details).

Materials: The Tableware of Ancient Bactria
The chemical, mineralogical, and petrographic composition of the studied collection was characterized in previous studies (Gurt Esparraguera et al., 2015;Martínez Ferreras et al., 2016 by using an ensemble of archeometric techniques, including Wavelength Dispersive X-Ray Fluorescence, X-Ray Diffraction, Scanning Electron Microscope with Energy Dispersive X-Ray Spectroscopy, and thin section optical microscopy. These studies prove that the chemical variability is very low between the ceramic vessels found in the same site, suggesting the use of very similar raw materials and therefore, a common local/regional origin. The variability detected is frequently explained by the use of different technological processes in the pottery manufacture (mainly in the processing of the clayey paste) and the development of secondary minerals due to weathering processes.
The materials used in this study (see Figure 1) are fine-to-medium-grained calcareous pottery fragments, fairly homogeneous in their composition but with distinctive inclusions depending on their provenance and their intended use. The ceramic wares from the alluvial plain (Kampyr Tepe and Termez) show non-plastic inclusions consisting of fragments of igneous and metamorphic rocks, together with crystals derived from granitoids (mainly quartz, feldspar and phyllosilicates); fragments and crystals of volcanic effusive rocks are scarce. By contrast, pottery from Kurganzol exhibits finer non-plastic inclusions consisting of the same lithology, although crystals and fragments of volcanic rocks are rare or absent, being a distinctive finding in the vessels from the Amu Darya basin. Since Kampyr Tepe comprised at least two pottery workshops dated to the Hellenistic and probably Greco-Bactrian periods, the provenance of most of the tableware analyzed from this site has been determined as local . In contrast, no workshop has been excavated from either in Kurganzol nor in Termez. It is highly unlikely that pottery was produced in the former military fortress, and therefore the vessels from Kurganzol have been considered as regional products. Moreover, the compositional differences detected among these samples and those found at Kampyr Tepe and Termez reveal that they were not produced in the Amu Darya floodplain. At Termez, all the samples analyzed were recovered from the alluvial plain (sectors AC1 and AC2) and have been considered as local/regional products (Martínez Ferreras et al., 2019). In contrast, for the tableware from the Tchingiz Tepe, dated to the Kushan and Kushano-Sasanian periods, a local provenance has been established (Gurt Esparraguera et al., 2015). It should be noted that Termez was an important ceramic production center during the Kushan and probably the Kushano-Sasanian period, as attested by several pottery kilns associated with the Buddhist communities (Martínez Ferreras et al., 2014;Tsantini et al., 2016).

Methods
Each pottery sherd (sample) selected for archeomagnetic analysis was sub-divided into different specimens-per-sample, yielding a total of 175 specimens that underwent experimental procedures. Specimens were named and marked, weighed, and encased in 2.3 cm each side quartz cubes, using a mix of glass wool and sodium silicate to prevent their moving during the heating and measuring stages. The rock-magnetism and paleointensity experiments were carried out at the paleomagnetic laboratories of the Complutense University of Madrid (Spain) and Géosciences-Rennes (Univ. Rennes, CNRS, UMR 6118, France).

Rock-Magnetic Experiments
Susceptibility versus temperature (χ-T) curves were carried out for a total of 29 samples using a Kappabridge KLY-4 susceptibility meter coupled to a CS4 furnace (AGICO). This device increases the temperature progressively and measures the bulk susceptibility of the powdered sample in each step. In this study, the temperature was increased from room temperature to 700°C, and then cooled back to room temperature in an Ar atmosphere. In addition, 28 saturation magnetization versus temperature (Ms-T) curves were done in a Variable Field Translation Balance (VFTB) (Mag-Instruments, Petersen & Petersen, 2008). Magnetization was measured in an applied field of 1 T warming up to around 650°C and further cooling back to 40°C in air. These experiments can be used to estimate the Curie/Néel temperature of the magnetic minerals contained in our samples, and to check for thermally induced magnetic alteration by comparing the heating and cooling branches of the plots. From the Ms-T curves, we calculated the corresponding Curie/Néel temperatures following the second derivative method proposed by Tauxe (1998) to find the inflection points and convexities of the curves. This was done with the RockMagAnalyzer 1.1 software (Leonhardt, 2006).
Hysteresis loops, isothermal remanent magnetization (IRM) acquisition curves and further direct current IRM demagnetization or back field saturation IRM (SIRM) curves were performed in 53 samples using a J-Meter Coercivity-Spectrometer, originally designed at Kazan federal University's paleomagnetism laboratory and improved by Jasonov et al. (1998). This instrument consists of a non-magnetic wheel, 50 cm in diameter, which rotates 18 times per second through an electromagnet. The in-field magnetic remanence of the sample, which is secured in a holder on the outer rim of the wheel, is measured while the electromagnetic field increases up to 500 mT and then decreases to the opposite polarity (−500 mT). The hysteresis loops, after calculation and subtraction of the paramagnetic contribution, were used to calculate the saturation magnetization (M s ), the saturation remanent magnetization (M rs ), the coercive field (H c ). The remanent coercive field (H cr ) associated with the samples is derived from the back-field IRM curves. The hysteresis parameters have been classically plotted in the so-called Day diagram, that was specifically designed to estimate the domain state of magnetite/titanomagnetite bearing rocks (Day et al., 1977) and where the magnetization parameters ratio (M rs /M s ) is plotted against the coercivity parameters ratio (H cr / H c ). We also plotted in this diagram the SD-MD mixing curves for magnetite (Dunlop, 2002). However, as noted by Roberts et al. (2018), the use of the Day plot to infer the domain state of the samples should be discarded. We have therefore only used the Day plot to identify anomalous samples in our collection. The shape of the hysteresis cycles can also be used to identify samples containing two or more magnetic populations with different coercivities, usually associated with wasp-waisted hysteresis loops (Tauxe, 1998). For this purpose, coercivity spectra of the IRM acquisition curves were analyzed using the model proposed by Kruiver et al. (2001) to determine the number and properties of the magnetic components. We also calculate the saturation IRM and the median destructive field at which half of the SIRM is reached. To do this, a curve fitting of the IRM versus the logarithm of the applied field is done using three different approaches: linear scale (LAP), gradient (GAP), and probability scale (SAP). Here we obtained, for the majority of the samples, the best fit of the IRM curves considering three magnetic components but for some specific cases it was necessary to also consider a 4th and 5th components, but being these components related with less than 10% of the total content. The component related to the lowest coercivity has been interpreted as a thermal activation of the magnetic particles and therefore neglected for further interpretation (Heslop et al., 2004).
Finally, in order to determine the existence of a single or multiple natural remanent magnetization (NRM) components we performed alternating field (AF) demagnetization experiments of the initial NRM over 30 pottery fragments, using a 2G cryogenic Superconducting Quantum Interference Device (SQUID) magnetometer. These experiments have been performed up to maximum fields of 120 mT in two different cycles, the first one from 0 to 40 mT with progressive steps of 5 mT; and the second one with steps of 20 mT, from 40 to 120 mT. 6 of 18

Paleointensity Experiments
To estimate the past geomagnetic field strength, we applied the classical Thellier-Thellier method (Thellier & Thellier, 1959), including partial thermoremanent magnetization (pTRM) checks, and the thermoremanent magnetization (TRM) anisotropy and cooling rate corrections upon intensity estimates. The magnetization was measured with a 2G SQUID magnetometer. Sample heating was performed in a Thermal Demagnetizer MMTD-24 (Magnetic Measurements) and the locally-manufactured "Ramsés" furnace at Géosciences-Rennes. The temperature was progressively raised from 100°C to 530-560°C in 10-11 steps, each step consisting of two heating and cooling cycles with an applied field of 50 µT along the Z axis of  (Day et al., 1977) with SD-MD mixing curves for magnetite (Dunlop, 2002) showing the results obtained for all the studied samples. (e) Isothermal remanent magnetization (IRM) acquisition curve and separation into three components to obtain the best fit (Kruiver et al., 2001). (f) Typical example of an alternating field (AF) demagnetization experiment from 5 to 120 mT. the samples: first in +Z, then in the opposite way. In order to check the thermal stability of the samples we measured the susceptibility after each temperature step with a MS3 susceptibility meter (Bartington Instruments), and performed pTRM checks every two steps of temperature. The effect of the TRM anisotropy on paleointensity determinations was calculated by measuring the acquisition of TRM in six different directions (±X, ±Y, ±Z) at a heating step for which around the 70% of the NRM is lost. In addition, and in order to check the TRM anisotropy protocol applied we decided to place sister specimens from one fragment (TRZ302) in three different faces of the cubes (denoted as specimens B, C, and L in Table S2 in Supporting Information S1; see also Osete et al., 2016 for a similar experiment). Finally, the cooling rate correction factors were calculated following the method described in Gómez-Paccard et al. (2006 after the TRM anisotropy determination or at the end of paleointensity experiments (when around 85%-95% of the NRM was lost). For this, we heated the samples at the same temperature in four different steps: the first two steps were made in the usual way, that is to say, during the first step the sample was heated during 45 min with the magnetic field applied in +Z direction and after cooled in about 45 min with the same magnetic field; in the second step the duration of the heating and cooling is the same as for the first step, but the magnetic field was applied in −Z direction. In order to determine the variation of the TRM acquisition, on the third step the sample was heated during 45 min with the magnetic field applied in +Z direction but the time to cool off was 24 hr; at the end of this additional cycle the samples are heated in −Z direction and the cooling is done in about 45 min. Comparing the third and first step we can calculate the effect of the cooling rate upon TRM acquisition capacity. Thanks to the second and fourth steps we can control potential changes in the TRM acquisition capacity related to thermally induced mineralogical changes. If the mineralogical evolution during the cooling rate is higher than the effect the cooling rate, this correction is not considered.

Rock-Magnetic Characterization
In order to identify the principal magnetic minerals responsible for the magnetization different rock-magnetism experiments have been performed (see Section 3.1). The majority of the thermomagnetic curves performed are almost reversible, which points to no significant magnetic alteration during the thermal treatment of the samples (Figures 2a and 2b). This indicates that the studied samples might be promising for successful paleointensity determination. The Curie temperatures (Figure 2a) obtained range between 269 and 640°C corresponding to different magnetic carriers. According to the results we divided our collection in two groups. The first group (including 15 samples, 55.6%) presents two different Curie temperatures: one at low temperatures, ranging from 269 to 380°C, which could be related to maghemite shells around magnetite cores (O'Reilly, 1984); another one is observed at high temperatures between 553 and 635°C which can be related to either magnetite or maghemite. At this stage, we cannot distinguish between these minerals, with quite similar values for the Curie/Néel temperatures 580°C (magnetite) and 610°C (maghemite) due to the large error in Curie temperature determination and linked to the smoothed Ms-T curves obtained. The second group (11 samples, 40.7%) corresponds to a single Curie temperature, ranging between 545 and 604°C, pointing to magnetite or maghemite as the principal magnetic carriers. Also, one exception has been identified, the sample TRZ595 shows two phases corresponding to titanomagnetites with different Ti content (T c = 481°C) and magnetite (T c = 561°C). In Figure 2b, one representative χ-T thermomagnetic curve per studied city is represented. In these plots we can visually ascertain a peak between 300 and 400°C, which could be associated with maghemite shells around magnetite cores. This observation seems to be consistent with the phenomenon and shapes observed by O'Reilly (1984) and mentioned earlier.
Hysteresis loops were used to calculate the hysteresis parameters (Figure 2c). M s varies between 0.024−0.520 Am 2 /kg and M r between 0.190-0.003 Am 2 /kg. H c ranges between 13.5-2.2 mT and saturation is reached at 150 mT, suggesting the presence of low-coercivity minerals. H cr varies between 47.0-14.0 mT. The hysteresis ratios obtained are presented in a Day plot (Day et al., 1977), modified as proposed by Dunlop (2002) (Nagata et al., 1963) (x-axis normalized by total partial thermoremanent magnetization) together with Zijderveld (Zijderveld, 1967) plots for representative specimens for (a) Kurganzol, (b) Kampyr Tepe's pottery workshop, (c) low city, (d) citadel, and (e-l) Termez. The initial natural remanent magnetization (NRM), the archeointensity obtained before (H) and after (H a ) the thermoremanent magnetization (TRM) anisotropy correction and the quality factor q defined by Coe et al. (1978) are indicated.
using J-coercivity measurements, which can be used to estimate SP existence with the decay of the SIRM when reversing the field (Jasonov et al., 1998). Therefore, it is very likely that the displacement towards higher values of the H cr /H c ratio could be related to some maghemitization of the samples. This effect is also observed for the two outliers from Kampyr Tepe which fall close to the MD (multi-domain) range. These two samples (TRZ554 and TRZ230) present high coercivity values which might indicate the presence of maghemite or hematite.
IRM acquisition up to 500 mT (Figures 2e and S1 from Supporting Information S1) were measured and the coercivity component analysis proposed by Kruiver et al. (2001) was performed. These experiments confirm the presence of a main low-coercivity fraction (1st component) with a median destructive field of about 30 mT, which contributes to 70%-80% of the total magnetization. The second component contributes to 10%-20% and corresponds to a median destructive field of about 70 mT. Finally, the third component, corresponding to around 5%-10% shows a very low median destructive field (of about 5 mT) that can be associated with thermal activation of the magnetic particles (see Heslop et al., 2004). These results suggest that the magnetite or titanomagnetite particles contained in our samples are close to the single domain state, ideal conditions to carry out successful paleointensity experiments. Finally, AF demagnetization of the NRM had been performed using sister specimens. Almost all the samples have a single component of magnetization which is completely demagnetized at 120 mT (see a typical example in Figure 2f). These results confirm the rock magnetic results previously described. In summary, rock magnetic analysis indicates that the principal magnetic carriers in our ceramic assemblage probably are low Ti-content magnetite, magnetite and/or maghemite. . Quality parameters at the specimen level obtained after the thermoremanent magnetization anisotropy correction. N, number of steps for the slope calculation; f, fraction of natural remanent magnetization used for the slope calculation; g, gap factor; q, quality factor; MAD, maximum angle of deviation (in °); DANG, deviation angle (in °); β, ratio of the standard error of the slope to the absolute value of the slope; k, curvature parameter (see main text for details).

New Paleointensity Data for Uzbekistan
A total of 175 specimens from 70 different fragments were prepared and measured. In Figure 3, representative Arai plots are presented. For paleointensity determination, in some cases we discarded the first steps of temperature (100−200°C) which correspond to a viscous component. Samples were completely demagnetized at temperatures between 480 and 560°C. In order to ensure the reliability of paleointensity determinations, we used common quality criteria similar than those used in our previous publications (e.g., Gómez-Paccard et al., 2016). The paleointensities have been only calculated for specimens that corresponded to a single component of magnetization with a straight line pointing to the origin in the Zijderveld plots and linear Arai diagrams (Figures 3a-3j). The maximum difference accepted between pTRM checks ant the original TRMs obtained at the same temperature step was 10%, this criterion was also applied for pTRM anisotropy checks. The fraction of the TRM used to determinate the paleointensity data, f parameter (Coe et al., 1978) should be at least 50% and 5 or more temperature steps must be used to determinate the slope. Following Shaar et al. (2016) we have considered a maximum value of 5° and 10° for the angular deviation MAD (Kirschvink, 1980) and the deviation angle DANG (Tauxe & Staudigel, 2004), respectively. It is worth noting that we have accepted a single specimen (KPT2229) with a MAD value slightly higher than 5° because all the other parameters are correct, and the intensity value obtained is very similar to the values obtained for the other two sister specimens studied (see Table S2 in Supporting Information S1). In addition, this high MAD value seems to be related to an anomalous high noise level in one single step rather than to a bad magnetic behavior. Around 74% (129 from 175) of the studied specimens have been retained for paleointensity determination.
In Figure 4, the different parameters used to evaluate the quality of the paleointensity estimates are showed. The fraction of the NRM used to calculate the paleointensity data ranges between 0.54 and 0.95. The quality factor q (Coe et al., 1978) varies between 10.2 and 115.9. The curvature parameter k (Paterson, 2011) and the ratio of the standard error of the slope to the absolute value f the slope (β ratio) range between 0.01 and 0.836 and between 0.0063 and 0.0466, respectively, although these parameters were not used as rejection criteria. As indicated before the highest MAD value obtained is 5.3° for one single specimen, being all the other values lower than 5° with a minimum of 0.8°. DANG values range between 0.2° and 3.3°. Thanks to the applied criteria we consider that the paleointensity estimates proposed can be considered as highly reliable.
The influence of TRM anisotropy and cooling rate on paleointensity determinations has been largely recognized. Here we analyze these effects at the specimen level. Important differences (up to 30%) are found between the uncorrected and TRM anisotropy corrected intensities (see Supporting Information S1 online). These results agree with previous archeomagnetic studies where very important differences have been obtained between the uncorrected and corrected intensities when pottery fragments are studied (e.g., Gómez-Paccard et al., 2019;Osete et al., 2016). In order to check the TRM anisotropy protocol applied we decided to place sister specimens from the Termez collection in different faces of the cubes. In general, specimens named with the letters A or B were placed in the central part of the base of the quartz cubes. Specimens denoted as C were placed perpendicular to the base of the cube and, finally, in some cases a third sister specimen, denoted as L, was placed in one of the lateral faces of the quartz cubes. It can be seen in Table S2 in Supporting Information S1 how very different intensity values are obtained before and after the TRM correction, especially for the Termez collection. The standard deviations (s.d) of the means are substantially reduced after the TRM anisotropy correction (see Table S2 in Supporting Information S1). For example, the s.d. for the fragment TRZ364 was 10.5 µT before the correction and 0.7 µT after. Concerning the cooling rate experiments, the corrections factors obtained (see Supporting Information S1) range between 0.7% and 11.3%, with a mean value of about 5%. This is in agreement with recent results obtained for archeomagnetic materials (Hervé, Chauvin, et al., 2019).
From the paleointensity determinations obtained at the specimen level (Table S2 in Supporting Information S1) we have derived 51 mean paleointensities: 11 calculated from 4 specimens, 13 from 3 specimens, 19 from 2 and 9 from just one single specimen (see Table 1). The more specimens we use, the greater statistical confidence on the mean value we have. Therefore, and following our recent publications (Osete et al., 2020;Rivero-Montero, Gómez-Paccard, Kondopoulou, et al., 2021), we classified our new data in two categories. The first one (called high-quality paleointensities) corresponds to data derived from three specimens and the lower quality group includes the other data obtained from 1 or 2 specimens. It is worth noting that since the main objective of this paper is to detect the rapid variations of the geomagnetic field, we decided to not calculate a mean value for the different age groups given in Table 1 since some of them are associated with large age uncertainties.
To have a global view of the geomagnetic field strength variations in our targeted area, the new mean intensities were relocated to Termez coordinates (37.3°N, 67.2°E) using the virtual axial geomagnetic pole ( Figure 5). The general trend described by our new collection suggests a V-shaped behavior of geomagnetic field strength in South Uzbekistan between 300 BCE and 200 CE. The new high-quality paleointensity results obtained in this study (calculated from 3 or more specimens, Table 1) suggest that the intensity of the geomagnetic field was much higher around 300 BCE than around 200 BCE where the lowest high-quality values are observed (see big dots in Figure 5) indicating fast intensity changes at multidecadal scales. Moreover, the new high-quality determinations derived from the well-dated site of Kurganzol (310-290 BCE, orange dots in Figure 5) suggest a large geomagnetic field strength variation (between 66.6 and 51.5 µT, in situ coordinates) at decadal scales. Therefore, the new data evidence a rapid intensity decrease in South Uzbekistan around 300 BCE. In Kampyr Tepe (blue dots in Figure 5), a total of 22 new intensities have been obtained dated between the end of 310 BCE and 100 BCE, with in general lower values than those obtained for Kurganzol, specially for the group of fragments from the pottery workshop (2d in Table 1). The high-quality paleointensities obtained for Kampyr Tepe range between 62.1 and 38.4 µT (in situ coordinates). However, the large age uncertainties associated with this collection makes it difficult to infer a detailed description about multidecadal intensity changes during this period. In Termez (purple dots in Figure 5), samples from 6 different age groups have given successful results, but only one of them corresponds to a high-quality intensity (fragment TRZ302 derived from 3 specimens). The 2 results obtained for the oldest group (3a in Table 1) indicate a variation between 47.6 and 41.9 µT during the 3rd century BCE. The results obtained for the stratigraphic unit AC1-SU 12 (group 3b) indicate a range of the intensity between 45.6 µT and 38.7 between 100 BCE and 100 CE. Finally, the mean value obtained for Tchingiz Tepe samples (groups 3c-f, 100-350 CE) of about 50 ± 4 µT further support a higher geomagnetic field intensity around this time, in comparison with our BCE collection. However, it is necessary to investigate more well-dated fragments in order to better describe geomagnetic field intensity features around the targeted periods.

A New Paleointensity Secular Variation Curve for Uzbekistan Between 600 BCE and 600 CE
The new mean intensities are compared in Figure 5 with previous paleointensity results from archeological sites located within a radius of 1,000 km around Termez and corresponding to the period 600 BCE-600 CE. As with the new data, the previous dataset was also relocated to the Termez coordinates. The number of the previous studies is significant and 72 paleointensities derived from the study of bricks, potteries or baked materials (Burlatskaya et al., 1995;Nachasova & Burakov, 1994, 1996, 1997. These older data were obtained with the Thellier and Thellier intensity method (1959) but, in case of Burlatskaya et al. (1995), without considering the TRM anisotropy effect upon paleointensity estimations. The importance of this effect on potteries has been recognized since several decades (e.g., Rogers et al., 1979) and it is now accepted that this correction needs to be performed at the specimen level. If not considered the real past intensity field values can be erroneously estimated with potential errors up to 90% (e.g., Genevey et al., 2008;Osete et al., 2020). In this study, for example, differences up to 30% between the uncorrected and corrected values are obtained (see Supporting Information S1). Moreover, the cooling rate effect upon intensity estimates was also not considered in these studies. In a recent work, Hervé, Chauvin, et al. (2019) have studied the influence of cooling rate in a set of 35 baked bricks and concluded that, on average, the non-corrected results overestimated the expected value by 5%-6%. These values are of the same order of magnitude as the mean correction proposed by Genevey et al. (2008) from an analysis of the global database. In our studied collection cooling rate factors up to 12% have been obtained with a mean value for the overall collection of 5% in agreement with the results obtained by Genevey et al. (2008) and Hervé, Chauvin, et al. (2019). Therefore, another potential source of error can be identified in the ancient data available for the targeted region. In order to consider this effect and in agreement with the studies cited before, a 5% reduction in the paleointensity estimates have been applied for the previous intensities available.
As it can be seen in Figure 5 new and previous data are quite consistent, although some discrepancies can also be observed. As discussed before, the reliability of the data cannot be considered as equal, and we decided to classify the data in three different categories using an appropriate weighting scheme designed for this purpose. Data obtained with 3 or more specimens and including the TRM and cooling rate corrections and pTRM checks were considered as the most reliable group and a weight of 10 was used (Campuzano et al., 2019). Data with TRM anisotropy and cooling rate corrections and pTRM checks but obtained with less than 3 specimens were assigned to an intermediate weight of 5 and, finally, the other data were considered as the less reliable group with an associated weight of 1. We calculated a new paleointensity paleosecular variation (PSV) curve for Central Asia using the compiled dataset together with the new results (shown in Figure 5) and using the mentioned weights. To fit the paleointensity data to a time-continuous curve we applied the same approach than Rivero-Montero, Gómez-Paccard, Kondopoulou, et al. (2021). In summary, this method provides an ensemble of 5,000 PSV curves by bootstrapping the paleointensity data in both age and intensity domains (including the error estimations). In addition, to obtain a more robust PSV curve, we eliminate possible outliers by detecting the data which residuals (the difference between data and curve) are three times higher than the root mean square error of the curve fitting. Only one data point has been identified as an outlier, data around 485 BCE (Burlatskaya et al., 1995) (see Table S3 from Supporting Information S1 for more Figure 5. New archeointensity results obtained for South Uzbekistan (Kurganzol, orange dots; Kampyr Tepe, blue dots; Termez, purple dots) and the new paleosecular variation (PSV) curves proposed here for central Asia between 600 BCE and 600 CE (orange line), plotted with its 2-sigma error band (shaded area). Also, the prediction of different global geomagnetic field models has been represented: SHA.DIF.14k (Pavón-Carrasco et al., 2014), SHAWQ Iron Age (Osete et al., 2020), SHAWQ.2k (Campuzano et al., 2019), ARCH10k (Constable et al., 2016), CALS10k.2 (Constable et al., 2016). The data, curves, and models are all relocated to Termez coordinates (37.3°N, 67.2°E). details). The mean and the standard deviation (1σ or 2σ) of the ensemble of 5,000 PSV curves provide the final PSV curve for South Uzbekistan. The new local PSV curve obtained ( Figure 5) shows an intensity maximum around 400 BCE followed by a rapid decrease of about −15 µT/century until 100 BCE. However, we must notice that this maximum is constrained by some archeointensity values higher than 80 µT that are included in the less reliable group, suggesting that this intensity bump might not reflect a real feature of the past geomagnetic field. New quality data covering this time interval are needed to confirm this hypothesis. The curve also shows an intensity minimum value around 100 BCE followed by a slowly growing trend until 500 CE. Finally, the compiled dataset has been compared with global paleomagnetic reconstructions: SHA.DIF.14k (Pavón-Carrasco et al., 2014), ARCH10k.1 and CALS10k.2 (Constable et al., 2016), and the SHAWQ family models (Campuzano et al., 2019;Osete et al., 2020). Results show notable discrepancies between data and model estimations. In detail, the models do not reproduce the minima and maxima values reported by the data and the intensity decrease observed in the 1st millennium BCE is much less pronounced in all these paleoreconstructions. These discrepancies (quantified for some time periods up to 15 µT) point out the need to obtain more quality data in the region to develop new and improved paleomagnetic reconstructions, and also emphasizing the importance of well-established local paleosecular variation curves.

Comparison With Neighboring Regions and Global Implications
To analyze if intensity changes in Uzbekistan are reflecting global or regional patterns of the past geomagnetic field, we have compared the virtual axial dipole moment (VADM) directly calculated from the data and the dipole moment (DM) and Axial Dipole Moment (ADM) derived from the global reconstructions previously mentioned. As shown in Figure 6, the VADMs exhibit much stronger fluctuations than the global ADM-DM, showing higher values between 600 and 300 BCE and lower estimations from 200 BCE to 300 CE. This suggests an important regional geomagnetic field behavior due to the non-dipolar contributions.
In order to investigate if similar geomagnetic field intensity trends are observed in other regions, we used the complete global database from 600 BCE-600 CE (GEOMAGIA50V3.4 database, Brown et al., 2015updated with recent publications: García Redondo et al., 2020Rivero-Montero, Gómez-Paccard, Kondopoulou, et al., 2021) to compile the VADM values available for different regions (1,000 km of radius, see Figure 7). The results are plotted from west to east in Figure 7. In Mexico, a VADM maximum of 16·10 22 Am 2 can be seen around 200 BCE but this feature is only based on two high-quality values . Nevertheless, this maximum is not observed in Peru, where high-quality data seem to present a minimum of VADM around this age. Moving to Europe, a maximum of 17·10 22 Am 2 is observed around 500-450 BCE followed by a decreasing trend until 100 BCE in Western Europe (Osete et al., 2020), Germany (Hervé et al., 2013(Hervé et al., , 2017Schnepp et al., 2020), Italy and Greece (Rivero-Montero, Gómez-Paccard, Kondopoulou, et al., 2021;Tema & Lanos, 2021), the Caucasian region and the Levant . Moving eastwards, in Central Asia, we also observe a VADM maximum around 400 BCE in Uzbekistan, around one century later than in Western Europe, but this maximum is based on low quality data that do not accomplish modern standard criteria of quality. In India, and due to the low number of data, a maximum cannot be clearly identified, but a decreasing tendency is also observed during the second half of the 1st millennium BCE. A similar V-shape trend is observed from Western Europe to India, with the minima values between 100 BCE and 200 CE. In Eastern Asia and Japan (Cai et al., 2017(Cai et al., , 2020Hong et al., 2013), there is no evidence of a maximum around 500-400 BCE if the high-quality dataset is considered, but some high VADM values of low quality are present in China. It is, however, hard to define a clear tendency in these areas due to the important dispersion of the corresponding regional databases. We also calculated the PSV curves obtained from the different regional data sets (see Supporting Information S1) following the methodology described before, using the mentioned categories of high-and low-quality data and different weights for each category, giving weights from 1 to 10 depending on the quality of the data. From this, the intensity secular variation rates (μT/century) have been calculated in 100-year time windows. The results, presented in Figure 8, showed a very similar trend all over Europe. The maximum variation (9.6 µT/century) is observed in Western Europe (Iberia and France-Germany) between 400 BCE and 300 BCE. This maximum variation is lower eastwards where values of 3.3 µT/century are obtained for the Caucasus between 200 BCE and 100 BCE. In Central Asia, the maximum values (14.7 µT/century) are achieved in Uzbekistan between 400 and 300 BCE, at the same time as in Western Europe. However, in the Levant and India the maximum variation is around three times lower than in Uzbekistan and appears between 300-100 BCE. As is mentioned previously, these rapid changes are not observed in Eastern Asia and the Americas.
The maximum rate of change obtained for Uzbekistan is similar to the one associated with the sharp intensity decrease recently identified between 1560-1725 CE in the same region (Troyano et al., 2021). It is, however, significantly lower than the rates of change of about 75-150 µT/century proposed for the LIAA (Ben-Yosef et al., 2017;Shaar et al., 2016). It is worth noting that the rates of changes proposed here are derived from PSV smoothed curves and therefore the maximum rates of change are most probably underestimated. This highlights the need to obtain a more detailed description of the past changes of geomagnetic field strength at the decadal scale. Finally, our results suggest a similar V-shape behavior for the intensity element (i.e., a maximum around 500-400 BCE followed by a rapid decrease) in a large region of more than 7,000 km from Western Europe (5°W longitude) to central Asia (80°E longitude). This wavelength on the surface of the Earth is equivalent to a spherical harmonic of degree n = 5-6, which is compatible with a non-dipolar source of the Earth's magnetic field generated in the outer core.

Conclusions
Fifty-one new paleointensity data, at the sample level, have been obtained for South Uzbekistan with ages ranging between 310 BCE and 300 CE, a period for which no high-quality paleointensity data were previously available. Our new results provide evidence of a rapid drop in intensity during the second half of the 1st millennium BCE showing a minimum of about 30 µT around 200 BCE-1 CE. Together with previous archeointensities our results allow us to present a new local PSV curve for Central Asia. The results suggest a V-shaped intensity behavior in Central Asia during the 600 BCE up to 600 CE period. High-intensity field values around 400 BCE were achieved in Uzbekistan followed by a sharp decrease until 100 BCE-1 CE associated with rates of changes of about −15 µT/century. A progressive recovery of geomagnetic field strength is then observed during the first centuries of the CE. This behavior is also observed at the continental scale, from Western and Central Europe to the Levant and Central Asia but cannot be detected eastwards in areas such as China and Japan or the Americas, although the low number of high-quality intensity data from this last continent hamper some more meaningful analysis. The important differences observed between the VADMs values available for Uzbekistan and the DM (ADM) estimations derived from different global geomagnetic field paleoreconstructions suggest the strong influence of non-dipolar sources upon this continental feature.