Simultaneous magnetic field and field gradient mapping of hexagonal MnNiGa by quantitative magnetic force microscopy

Magnetic force microscopy (MFM) is a scanning microscopy technique that is commonly employed to probe the sample’s magnetostatic stray fields via their interaction with a magnetic probe tip. In this work, a quantitative, single-pass MFM technique is presented that maps one magnetic stray-field component and its spatial derivative at the same time. This technique uses a special cantilever design and a special high-aspect-ratio magnetic interaction tip that approximates a monopole-like moment. Experimental details, such as the control scheme, the sensor design, which enables simultaneous force and force gradient measurements, as well as the potential and limits of the monopole description of the tip moment are thoroughly discussed. To demonstrate the merit of this technique for studying complex magnetic samples it is applied to the examination of polycrystalline MnNiGa bulk samples. In these experiments, the focus lies on mapping and analyzing the stray-field distribution of individual bubble-like magnetization patterns in a centrosymmetric [001] MnNiGa phase. The experimental data is compared to calculated and simulated stray-field distributions of 3D magnetization textures, and, furthermore, bubble dimensions including diameters are evaluated. The results indicate that the magnetic bubbles have a significant spatial extent in depth and a buried bubble top base. A quantitative magnetic force microscopy technique is presented that maps one magnetic stray-field component and its spatial derivative at the same time. Furthermore, this technique is applied to investigate individual circular magnetic nano-domains in MnNiGa bulk samples providing bubble diameters and the spatial extent in depth.

C urrent developments in the area of nanomagnetism, including studies on very small non-collinear magnetic textures, call for more advanced, quantitative magnetic microscopy solutions. Nanoscopic magnetization patterns can be accessed indirectly via their emanating stray-field distributions by a variety of microscopy techniques such as scanning probe microscopy based on single nitrogen-vacancy defects in diamond 1 , superconducting quantum interference devices 2 and force sensors equipped with ferromagnetic tips. The latter, magnetic force microscopy (MFM), provides spatially resolved maps of magnetostatic probe-sample interactions at the nanometer scale [3][4][5][6] . Basic MFM is a relatively straightforward method being widely available in labs worldwide. Therefore, it is extensively used for nanomagnetic characterization 7 and, in particular, for skyrmion research 8,9 . Here, MFM is not only beneficial for providing stray-field maps; the magnetic tip may also act as a tool for localized skyrmion creation and manipulation 10 . Sophisticated quantitative MFM (qMFM) data evaluation methods even enable the analysis of the skyrmion type (Néel vs. Bloch), helicity and size of skyrmions, and the force required to move individual skyrmions [11][12][13] .
The majority of commercially available MFM probes is based on Si cantilevers with magnetic thin-film coated pyramidal or conical tips. It turns out that their magnetostatic behavior differs from that of a magnetic point dipole or monopole 27 .
A number of reports describe the application of magnetic nanowires as major components of MFM probes. Such nanowire probes offer a variety of promising features including high spatial topographic and magnetic resolution 36 , predictable single domain magnetization stabilized by shape anisotropy [37][38][39] , superior mechanical stability in case of magnetic nanowires encapsulated in multi-wall carbon nanotubes 40,41 , and in addition, monopoletype stray-field characteristics for high aspect ratio single-domain magnetic nanowires [31][32][33]42,43 . Furthermore, in addition to their magnetic functionalities, nanowires may directly constitute mechanical force transducers for MFM operation [43][44][45] . Mattiat et al. 43 used a nanowire-based MFM approach to quantitatively measure the amplitude of an oscillating magnetic field with nT/ √Hz sensitivity at low temperatures.
Here, we report on a single pass multi-frequency MFM approach that provides immediately accessible quantitative maps of both the magnetic field component B z and its spatial derivative B 0 z ¼ ∂B z =∂z. B 0 z maps show much sharper magnetic features compared to B z maps. This is not a consequence of spatial resolution or sensitivity, but it is rather a demonstration that field derivatives simply amplify higher spatial frequencies compared to their respective field distributions 22 . Our simultaneously measured B z and B 0 z maps immediately illustrate these characteristic differences. Furthermore, recording both quantities may provide information on the magnetization distribution in the depth of the sample in specific cases as shown herein.
Our magnetostatically active probe element is an iron-filled carbon nanotube (FeCNT), containing a single domain iron nanowire (FeNW) with an aspect ratio as high as 250. We simultaneously measure the deflection of a tailored low-stiffness cantilever and the frequency shift of one of its flexural vibration modes in fully noncontact operation. Snap-in events are avoided by electrostatic distance control employing a higher order flexural vibration mode 46 . The monopole-type characteristics of the FeNW tip allows for a straightforward translation of the cantilever deflection and the frequency shift into B z and B 0 z , respectively. The validity of the pointmonopole description with respect to particular spatial wavelength intervals of magnetic stray-field patterns will be thoroughly discussed herein.
Skyrmionic spin textures constitute a promising class of nanomagnetic configurations with particle-like properties. Magnetic skyrmions are associated with topological protection and attract a lot of interest both in fundamental magnetism research and spintronics such as advanced memory technologies and neuromorphic computing 47 . There are various types of magnetic skyrmions with different topological charge Q t including Bloch skyrmions (Q t = 1), Néel skyrmions (Q t = 1), and biskyrmions (Q t = 2). Biskyrmions, i.e. bound pairs of Bloch skyrmions, have been observed by Lorentz transmission electron microscopy in a centrosymmetric bilayered manganite 48 and in hexagonal MnNiGa [49][50][51][52] . Skyrmions in MnNiGa were observed at various temperatures including room temperature and were reported to survive even at zero external field 50 . In other studies, however, it was shown that Lorentz transmission electron microscopy data of magnetic hard bubbles, i.e., type-II bubbles with Q t = 0, might be misinterpreted under certain conditions as Q t = 2 biskyrmions 53,54 . Also, a recent soft X-ray holography investigation found single-core type-II bubbles in MnNiGa 55 . Nevertheless, theoretical studies present the prospect of a biskyrmion formation in centrosymmetric materials with uniaxial magnetic anisotropy [56][57][58] . The above mentioned experimental MnNiGa studies were limited to thin lamella samples. Here, we investigate magnetic nano-domains in bulk MnNiGa.
After detailing the fundamentals of the approach for simultaneous B z and B 0 z measurements in the "Concept of combined field and field gradient microscopy" section of the Results and Discussion part, we report on the experimental realization of a corresponding MFM probe design in "Sensor design and fabrication" in the Results and Discussion. Finally, we apply our technique to the study of magnetic bubble-like domains in bulk [001] MnNiGa ("Quantitative B z and B 0 z measurements on bulk [001] MnNiGa" in the Results and Discussion). We focus on the specific question of how the investigated bubbles continue along the third dimension, as discussed for other materials 59,60 , and determine their diameters ranging from 120 nm to 200 nm. Our results indicate that the magnetic bubbles have a significant spatial extent in depth and a buried bubble top base.
The presented method will boost the capabilities of commonly used magnetic force microscopy equipment towards quantitative nano-magnetic microscopy.

Results and discussion
Concept of combined field and field gradient microscopy. Our MFM measurement concept allows for combined quantitative magnetic field B z and field derivative B 0 z scanning with simultaneous sample topography tracking in a single-pass, true non-contact operation mode (see Supplementary Note 1 and  Supplementary Table S1). The presented microscopy modes can be performed in standard AFM/MFM setups but require a specific probe design as introduced in "Sensor design and fabrication" in the Results and Discussion. Such probes are not yet commercially available.
As a key feature, once the sensor is calibrated, our MFM approach provides direct B z and B 0 z data at any location in real space without any postprocessing and without the need for carrying out additional data evaluation in Fourier space.
The concept map ( Fig. 1) illustrates an example on how qMFM for both magnetic fields and magnetic field gradients can be achieved. A magnetic point charge exposed to a magnetic field or field gradient is subject to a force or a force gradient, respectively. Forces and force gradients can easily be measured via deflections and frequency shifts of micro-cantilevers. In terms of magnetostatics, the end of an ideal infinitely long and infinitesimally thin ferromagnetic wire, i.e. a Dirac string, behaves like a point monopole (probe nanowire shown in the upper part of Fig. 2a) and generates a radially symmetric magnetostatic stray-field. Thus, in our approach, high-aspect ratio iron nanowires are used to approximate such magnetic point monopole probes. A quantitative evaluation of the validity and the limits of the point monopole model when describing such magnetic nanowire probes is discussed in what follows (see Fig. 2b-f).
To enable sensitive static force measurements, low spring constant cantilevers are necessary. Generally, they work well for contact mode measurements or in pendulum geometry. However, they are not easily usable for dynamic mode measurements in conventional cantilever geometry since attractive contributions to the overall tip-sample interaction may lead to snap-in events that can damage or change the tip's characteristics. Thus, intermittent contact or tapping mode is very difficult to implement especially when using soft cantilevers. To enable the use of low static spring constant k stat cantilevers in dynamic operation modes, an independent tip-sample distance control is required that prevents any direct tip-sample contact. Therefore, we employ an electrostatic distance control 17 that ensures a safe tip-sample distance for true non-contact operation.
Once an appropriate tip-sample distance z ts is set and maintained and, additionally, the tip monopole moment q t is known, measurements of the cantilever deflection Δz and of the i-th order flexural mode resonance frequency shift Δf i yield B zẑ Ã ð Þ and B 0 zẑ Ã ð Þ, respectively (see Fig. 1). Here,ẑ Ã refers to the z-position of the probe's magnetic pole q t that is closest to the sample (Fig. 2a).
Prospects of simultaneous quantitative field and gradient measurements. Simultaneous quantitative magnetic field B z x; y;ẑ Ã À Á and field derivative B 0 z x; y;ẑ Ã À Á mapping allow for direct linear Taylor approximations of B z at z-locations in the vicinity of the Fig. 1 Concept map. Concept map of combined field and field gradient microscopy as proposed and employed in this publication. Here, k stat , 4z, k i , f i , 4f i , and q t are the cantilever's static spring constant, the cantilever deflection, the spring constant of the i-th flexural vibration mode, the resonance frequency of the i-th flexural mode, the corresponding resonance frequency shift, and the monopole-like magnetic tip moment. contained in a protecting carbon nanotube (in grey); the lower part shows a magnetic sample (sample surface at z ¼ 0) containing a bubble-like magnetic single-domain configuration (in blue), which is embedded in an oppositely magnetized environment (red). q s indicates the upper magnetic pole associated with the magnetic charge of the bubble base. b-f indicate several types of nanowire (NW) geometries and magnetization configurations. b The end of an ideal infinitely long and infinitely thin nanowire probe generates a stray-field resembling the field of a point monopole q t . Deviations from point monopole behaviour can be caused by: (c) finite NW length, (d) non-zero NW diameter, (e) non-homogeneous internal magnetization, and (f) a NW magnetization direction that deviates from the long (easy) axis. measurement location ðx; y;ẑ Ã Þ, for which: We illustrate the potential of simultaneous B z and B 0 z measurements for the case of three simple and idealized strayfield distributions (Fig. 3). In case of buried point-like stray-field sources, both their depth z* below the sample surface and their magnetic monopole moment (i.e., magnetic charge) q s or dipole moment m s can be directly inferred from simultaneously available B z and B 0 z data ( Fig. 3a, b, and Supplementary Note 2). Furthermore, in case of a sample that generates a B z distribution that is sinusoidal along in-plane directions with wavenumber , the ratio of B 0 z to B z immediately provides k xy without the necessity of performing a full 2D surface MFM scan ( Fig. 3c and Supplementary Note 2). k xy ¼ k x ; k y is the corresponding in-plane wavevector. Figure 3 sketches the potential of a simultaneous B z and B 0 z measurement. Inferred quantities like depth or width of the magnetization textures can be directly accessed, if the measured magnetization textures are well approximated by idealized monopoles or dipoles, or if being periodic in a lateral direction (i.e. parallel to the surface). B z and B 0 z complement each other. On the one hand, B z might be considered as more fundamental data. On the other hand, as mentioned in the Introduction, B 0 z maps generally show sharper features since higher spatial frequencies are amplified by a factor of k xy as compared to B z mapping 22 .
A perfect magnetic point monopole probe enables us to measure B z and B 0 z . A single domain magnetic nanowire exhibits several aspects that lead to deviations from the magnetic behaviour of a magnetic point monopole: (i) a finite nanowire (NW) length, (ii) a non-zero NW diameter, and deviations from a homogenous NW magnetization parallel to the NW axis due to (iii) demagnetization effects, and (iv) external magnetic fields ( Fig. 2b-f).
Validity of the point monopole approximation: finite length and non-zero diameter effects of nanowire probes. In the experimental part: "Sensor design and fabrication", our MFM sensor is presented which features a d Fe ¼ 20 nm and l Fe ¼ 5 μm iron NW as magnetic probe. Compared to the hypothetical singlepoint monopole probe, this real NW comprises two magnetic poles with the magnetic charge of each pole distributed over the NW's cross section. We quantitatively evaluate the effects of non-zero diameter and finite length that constitute deviations from a point monopole probe behaviour. According to Hug et al. 22 , the relation between the force at a magnetic NW probe and the field can be described by a force transfer function (FTF) that depends on the wavenumbers k x , k y of the field distribution.
For a NW probe with square-shaped cross-section of side length d, the FTF reads 22,38 : in which M S and l are the NW probe saturation magnetization and the NW length, respectively. Note that for l → ∞ and k x ; k y ( 2 d the FTF NW approaches the NW's magnetic charge q t ¼ d 2 M S . Under these conditions, the probe shows the same characteristics as a point monopole.
Analogously, the FTF of a NW probe with circular NW crosssection can be derived by taking advantage of an expression for the Fourier transformation of a unity circle, described by spherical Bessel functions of the first kind. After scale conversion and normalization with respect to the magnetic NW diameter d, we obtain: . Furthermore, for simplicity, here we use the same symbol d as above for the side length of a square-shaped cross-section. Again, for l → ∞ and k x ; k y ( 2 d the FTF NW approaches the magnetic charge q t ¼ π 4 d 2 M S . In Fig. 4, calculated FTF(λ) normalized by the FTF of a point monopole probe, i.e., by q t , are shown. Instead of k xy we now use the wavelength λ ¼ 2π=k xy to describe the spatial variation of the sample's stray-field. The FTFðλÞ=q t ratio for a magnetic probe with NW dimensions similar to that of our fabricated MFM probe introduced below is shown in Fig. 4a. For a broad band of λ values the FTF(λ) almost perfectly matches q t . It implies that in this range the simple point monopole description with magnetic charge q t is fully sufficient for high precision B z and B 0 z mapping. The FTF(λ) reduction for very small and very large λ can be attributed to the finite dimensions of d and l, respectively. Figure 4b-e reveal corresponding calculations for a variety of dl pairs.
In Fig. 4c, e, λ-l and λ-d phase diagrams depict calculated areas with FTFðλÞ=q t ratios exceeding 99%, 95%, and 90%. These calculated areas can be approximated by triangles even though there are no sharp vertices on the left-hand side of Fig. 4c and the right-hand side of Fig. 4e. The upper area boundaries are mainly predefined by l. For large l (in Fig. 4c) and small d (in Fig. 4e), these upper boundaries λ max are proportional to the NW length l. They can be approximated by equating the ð1 À e Àk xy l Þ term of equation 3 with 99%, 95% or 90% and substituting k xy by λ ¼ 2π=k xy . For the 99% case λ max % 1:36l and for the 90% case λ max % 2:73l. As an example, an l=1 μm NW MFM probe can be used to map B z and B 0 z of a periodic λ ≤ 1.36 μm strayfield distribution with ≥99% accuracy using the simple point monopole model.
The lower boundaries λ min of the triangles in Fig. 4c, e are mainly defined by the NW diameter d. For large l (in Fig. 4c) and small d (in Fig. 4e) we get λ min % 11:0d (or λ min % 3:44d) corresponding to 99% (or 90%). As an example, a NW MFM probe with d = 10 nm can be used to map B z and B 0 z of a periodic λ ≥ 110 nm stray-field distribution with ≥ 99% precision using the simple point monopole approach.
Validity of the point monopole approximation: nanowire probes with non-homogeneous and non-aligned magnetization. Both non-homogeneous magnetization of the NW and a NW magnetization direction deviating from the long NW axis lead to deviations of a real NW's behaviour from that of a magnetic point monopole. Apart from non-trivial internal magnetization effects deviations can additionally be promoted by local sample stray-fields or an external applied field.
Therefore, the internal magnetization structure of a simple l = 5 μm Fe cylinder with its axis oriented along the z-direction was modelled using the MuMax3 micromagnetic simulation software (version 3.9.1c) 61,62 . For an individual cell, the spatial extent in z-direction amounts to 2.4 nm, not exceeding the exchange length of iron of approximately 3.5 nm. For the exchange constant, we assume 22 pJ/m being within the reported range of values for pure iron [63][64][65][66] .
For the saturation magnetization, M S ðFeÞ ¼ 1710 kA=m was used 67 . Here, we neglected the effect of magneto-crystalline anisotropy as the overall magnetic anisotropy is mainly governed by the high-aspect ratio NW's shape anisotropy 37 .
In the simulation results of Fig. 5, M z /M S along the long NW axis (z-axis) is shown. For a d = 21 nm NW, M z =M s > 0:97 at the free NW end even in case of a 300 mT (−200 mT) applied field nearly parallel (antiparallel) to the z-axis.
Here, we assumed that the applied fields do not exceed the switching field of the NW. For similar FeNWs and external fields parallel to the NW axis, switching fields between 100 mT and 400 mT are reported [68][69][70] . Small deviations from M z =M S ¼ 1 are present for the last 10 nm close to the NW tip apex. The application of a perpendicular 100 mT field leads to a small overall reduction of M z , i.e., M z =M S % 0:99 along the entire NW and a further decrease towards the NW end. For larger NW diameters, vortex states form close to the NW free end which are accompanied by a local reduction of M z =M S . Altogether, small NW diameters and external fields much smaller than both NW anisotropy and switching fields support the application of a point monopole model for the magnetostatic behaviour of magnetic NWs.
Sensor design and fabrication. We developed a cantilever design based on the following demands: (i) low static spring constant k stat ; (ii) higher order flexural modes with high resonant frequencies, i.e., the frequency ratios f i /f 1 should be much larger as compared to the case of a simple Euler-Bernoulli cantilever; (iii) efficient fabrication routine by focused ion beam (FIB) preparation and micromanipulation.
Note that requirement (ii) is related to the i-th mode's settling time τ i , which represents the bandwidth-limiting factor for a distance control relying on the oscillation amplitude A i as a feedback loop input parameter. Further details are introduced in section: "Magnetic force microscopy setup and modes of operation". The sensor preparation according to our MFM sensor concept is based on a commercial picket-shaped cantilever (Nanoworld TL1) consisting of highly doped (001) single-crystalline Si with nominal k stat ¼ 30 mN=m. The dimensions of the cantilever geometry as provided by the manufacturer read: length l ¼ 500 μm, thickness t ¼ 1 μm, width w = 100 μm. To support low spring constants, we selected a cantilever with a particularly low thickness. The cantilever dimensions, as measured by scanning electron microscopy (SEM) read: total cantilever length l = (503.8 ± 2.2) μm, the length of the triangular part l t ¼ 93:0 ± 1:7 ð Þμm, the total width w = (93.6 ± 1.4) μm, and the thickness t = (0.795 ± 0.015) μm.
To fulfil the requirements mentioned above, we applied FIBbased line cutting along a pre-defined Π-shape using a combined SEM-FIB (ZEISS 1540 XB). The result is a double-cantilever beam structure and an inner paddle as shown in the SEM image of Fig. 6a. The desired low k stat is mainly determined by the combined Fig. 4 Validation of the point monopole approximation. Normalized force transfer function FTF=q t for nanowire (NW) probes with circular cross-sections. a NW dimensions are: diameter d ¼ 20 nm and length ¼ 5 μm. The red line refers to a calculation based on a circular NW cross-section. For comparison, results for a square-shaped NW are also shown (dotted line; note that a side of the square is aligned parallel to k xy ). Purple and reddish-colored areas indicate the λ ranges in which FTFðλÞ=q t is >99% or >90%, respectively. k xy and λ ¼ 2π=k xy are the wavevector and the wavelength of the spatial magnetic field distribution, respectively. b, c The impact of varying l is shown while d ¼ 20 nm is kept constant. d, e The impact of varying d is shown while l ¼ 5 μm is kept constant. c, e λ À l and λ À d phase diagrams showing calculated areas with FTF=q t ratios exceeding 99%, 95%, and 90%.
width 2 · b of the two remaining high aspect ratio beams. Those beams with b = (10 ± 1.0) μm connect the paddle structure to the carrier chip as shown in Fig. 6a. The width b is much smaller than the pristine cantilever width w. In fact, the rectangular inner paddle has a huge impact on the cantilever dynamics but negligible influence on the static spring constant.
For optimizing the sensor design and to characterize its mechanical properties, we employed simulations based on finite  In the inset, the iron nanowire (FeNW) contained in the carbon nanotube is visible. The diameter d Fe of the FeNW and the non-magnetic cap thickness t cap are indicated. c SEM side view of the sensor's free end disclosing the full NW length. The 275 nm gap splits the FeNW into two parts with corresponding lengths l Fe;1 and l Fe;2 . d COMSOL beam shape simulations of static bending caused by an applied tip load of F z ¼ À100 pN and for the 2nd and 5th flexural oscillation modes. The colors describe the absolute value of the deflection. e Magnetic induction (B-field) mapping of a typical iron-filled carbon nanotube (FeCNT) obtained from a phase image reconstructed by off-axis electron holography in the transmission electron microscope. The red iso-lines in steps of π⁄20 display the projected stray-field lines around the FeCNT, which is visualized as greyscale phase image. element analysis (Fig. 6d) with the structural mechanics module of COMSOL Multiphysics® v. 5.4 (COMSOL AB, Stockholm, Sweden). For the simulations, we neglected non-linearities in our study since the cantilever deflection is small compared to the total sensor length. The material properties of the simulated cantilever read: Silicon mass density ρ = 2.329 g/cm³, elastic modulus E = 170 GPa and Poisson ratio ν = 0.28. Similar values for E and ν are reported in literature 71 .
The selected FeCNT with d Fe ¼ 20 ± 2 ð Þnm was attached to the cantilever end using a Kleindiek nanomanipulator by electronbeam assisted deposition of carbon (Fig. 6b, c). To experimentally confirm the monopole-type magnetic stray-field near the apex of FeCNTs, we performed transmission electron microscopy (TEM) based off-axis electron holography 77 to an about 8 μm long FeCNT stemming from the same sample as the FeCNT that was used for assembling our sensor tip (see Supplementary Note 6, Supplementary Fig. S4). Figure 6e depicts the reconstructed phase image of an FeCNT tip in which the magnetic phase shift in vacuum is emphasized by red isolines. The latter are proportional to the field lines of the projected magnetic induction revealing indeed a monopole-like behaviour except at the lower side. These perturbations can be explained by stray fields generated by contamination parts attached to the FeCNT in the vicinity of the field of view (see Supplementary Note 6, Supplementary Fig. S5). We also supported these TEM observations by magnetostatic and image simulations in combination with tomographic methods (see Supplementary Note 6, Supplementary Fig. S5).
Using a piezo-driven oscillation stage within the SEM-FIB, we comprehensively investigated the fundamental and higher order flexural oscillation modes. Imaging the oscillation envelopes allows for a visual analysis of flexural mode shapes. Moreover, amplitude-frequency curves of the sensor were measured and evaluated using the model of a damped harmonic oscillator 78 providing mode-dependent resonant frequencies and Q-factors. These and all further MFM calibration and application measurements described in the next section were performed in our highvacuum MFM (NanoScan AG hrMFM) supported by a lock-in amplifier (Zurich Instruments HF2LI). Our experimental setup is sketched in Fig. 7. All modes of operation presented in this work can also be performed in ambient conditions. Magnetic force microscopy setup and modes of operation. The position-sensitive detector (PSD) signal, which contains information about the measurement quantities f i , A i , and 4z, is fed into and demodulated by both the MFM main controller and the additional lock-in amplifier. The former provides the static cantilever deflection readout and controls the second flexural mode with the first phase-locked loop (PLL 1) tracking changes of the resonance frequency f 2 . Moreover, an integrated proportionalintegral-(derivative) controller (PID 1) ensures a constant second mode vibration amplitude A 2 . The second flexural mode is selected to measure force gradients via frequency shifts as it is a stable low order vibration mode associated with low k 2 .
Obviously, the fundamental mode comes with an even lower spring constant, but its low resonant frequency was not compatible with our MFM setup. The second mode vibration is excited by piezo-driven mechanical actuation of the cantilever base.
An additional phase-locked loop (PLL 2) provided by the external lock-in amplifier tracks the resonance frequency f 5 of the 5th flexural mode and drives this mode by applying a tip-sample AC bias voltage at f AC ¼ f 5 =2, i.e., at half of the mode's resonance frequency. A 5 is proportional to V 2 AC , Q 5 , and to the tip-sample capacitance gradient which is a measure of the tip-sample distance. Keeping A 5 constant by a properly set PI(D)-regulated tip-sample distance control enables tracking of the sample surface topography. Furthermore, a well-defined tip-sample distance prevents snap-in events that may arise from a large sample surface roughness and the usage of low-k probes. Selecting the 5th mode provides a sufficiently small settling time τ 5 (see Supplementary Note 3) which guarantees an acceptable bandwidth for the distance control feedback loop.
Spring constant calibration. The static spring constant k stat of our sensor shown in Fig. 6 was determined by means of the above mentioned COMSOL simulations by virtually applying a small force (e.g., 100 pN) to the sensor tip, with the force vector pointing perpendicular to the cantilever basal plane, resulting in k stat ¼ 3:3 mN=m. We compared this simulation result with an analytical approximation describing an arrow-shape cantilever with a shaft of width 2 b and a triangle of base width w and height (i.e. length of triangular free end) l t 79 : resulting in k stat ¼ 3:4 ± 0:6 ð ÞmN=m, which is consistent with the simulation. The dynamic spring constants k i were determined by COMSOL simulations by virtually attaching a small spring (4k ( k i , e.g., 4k ¼ 1 μN=m) to the sensor tip resulting in a small eigenfrequency alteration 4f i of the respective mode. Based on the harmonic oscillator model, the dynamic spring constant k i Fig. 7 Quantitative magnetic force microscopy setup. Schematic of the experimental setup and the multi-mode approach that allows for simultaneous magnetic field and field gradient microscopy supported by non-contact topography tracking. The magnetic force microscopy (MFM) setup (based on a Nanoscan AG hrMFM high-vacuum MFM) comprises a laser-beam deflection detection based position-sensitive detector (PSD). The MFM main controller features one integrated phase-locked loop (PLL 1) used for tracking the resonance frequency f 2 of the mechanically driven second mode (dither). An additional PLL 2 realized by an external lock-in amplifier (Zurich Instruments HF2LI) tracks the AC-driven high frequency mode resonance f 5 . A 5 is kept constant by a proportional-integral-(derivative) distance control feedback-loop (PID 2) adjusting and maintaining a constant tip-sample distance z ts during scanning. Cantilever mode shapes including the static deflection Δz and two dynamic flexural modes (f 2 and f 5 ) with their respective amplitudes (A 2 and A 5 ) are indicated separately (dotted lines) and superimposed (solid black line).
can be calculated according to: As far as the important second flexural mode is concerned, we obtain k 2 ¼ 0:47 N/m. Further k i results are presented in Supplementary Note 3 (Supplementary Table S2).
To experimentally confirm the simulation results, the second order dynamic spring constant k 2 was calibrated by comparing F z ðz ts Þ and 4f 2 ðz ts Þ sweeps (Fig. 8). For these sweeps a sample location was selected that provides attractive magnetostatic interactions and thus monotonic tip-sample forces. The F z ðz ts Þ sweep was taken without any excitation of dynamic modes.
Considering a non-linear F z ðz ts Þ, the following relation based on an analysis of F. J. Giessibl 80 is used: Solving this equation by numerical integration yields k 2 ¼ 0:45 ± 0:07 ð ÞN=m. Note, that for the calibration measurement shown in Fig. 8, all z ts data are based on the height control piezo data corrected by the static cantilever deflection. Hence, z ts reflects the real distance between the tip apex and the sample surface.
Calibration of the magnetic tip moment. The knowledge of k stat , k i , and f i allows us to quantitatively measure magnetostatic forces and force gradients, i.e., F z and F 0 z , respectively. These quantities can be translated into B z and B 0 z using the magnetic charge q t according to the point monopole description of our NW tip. The calibration procedure is based on MFM scans of a known magnetic reference sample 33,81 (see "Magnetic reference sample" in the "Methods" section) and a quantitative analysis of the MFM data. Furthermore, the results are compared with a simple q t determination using the NW's geometry and saturation magnetization M s ðFeÞ.
By calculating an effective magnetic surface charge map of the reference sample based on its materials parameters (given in the "Magnetic reference sample" section) and the known sample thickness of 130 nm, the MFM tip's transfer function is derived and fitted to a point monopole model. The details of this calibration procedure are described elsewhere 22,33,81 . The calibration results in the effective monopole moment q t ¼ 0:50 ± 0:10 ð ÞÁ10 À9 Am and the position of the effective point monopole with respect to the physical tip apex. This distance is mainly defined by the FeCNT carbon cap layer thickness t cap , even if, technically, both are not the same. For clarity, we use the term t cap for both quantities and obtain t cap ¼ 33 nm in agreement with SEM analysis. This value is also used for the evaluation of the MFM scan images shown in section: "Quantitative B z and B 0 z measurements on bulk [001] MnNiGa". Note that due to electron-beam induced local contamination during repeated SEM inspection of our MFM probe, the thickness of the non-magnetic cap layer increased gradually by a factor of two and finally reached t cap % 70 nm (Fig. 6).
For comparison, a q t determination via direct measurement of the FeNW diameter d Fe ¼ 20 ± 2 ð Þnm and taking account of the saturation magnetization M s ðFeÞ ¼ 1710 kA=m would result in a tip moment q t ¼ π 4 d 2 Fe M s ¼ 0:54 ± 0:11 ð ÞÁ10 À9 Am.
Detection sensitivities. In all types of cantilever-based scanning force microscopy measurements, including those presented herein, thermal noise and other noise sources limit the detection sensitivity. The minimal detectable force limited by thermal noise F min = ffiffiffiffiffiffiffi ffi BW p is given by 82 : in which BW denotes the detection bandwidth. Considering a force that is periodic with the fundamental mode frequency f 1 , the minimal detectable F min = ffiffiffiffiffiffiffi ffi BW p of our sensor as introduced above is 0:93 fN= ffiffiffiffiffiffi Hz p at T ¼ 300 K and 2:72 fN= ffiffiffiffiffiffi Hz p for f 2 . Given the monopole-like moment q t ¼ 5 Á 10 À10 Am, the minimal detectable magnetic fields at T ¼ 300 K read B z ¼ 1:9 μT= ffiffiffiffiffiffi Hz p and B z ¼ 5:44 μT= ffiffiffiffiffiffi Hz p for f 1 and f 2 , respectively. The minimal detectable force gradient limited by thermal noise F 0 min = ffiffiffiffiffiffiffi ffi BW p is given by 83 : In the present case, F 0 min = ffiffiffiffiffiffiffi ffi BW p at f 2 ¼ 34130 Hz is 270 aN nm = ffiffiffiffiffiffi Hz p at T ¼ 300 K when using the second flexural mode. Considering q t ¼ 5 Á 10 À10 Am, the minimal detectable field derivative B 0 z at T ¼ 300 K is 545 nT nm = ffiffiffiffiffiffi Hz p . For these calculations, spring constants k i , quality factors Q i , and frequencies f i for the fundamental and second flexural modes are applied according to Supplementary Note 3 (see Supplementary Table S2). The vibration amplitude is A 2 ¼ 10 nm.
Please note that the mean square deflection induced by thermal noise of an Euler-Bernoulli cantilever appears to be enhanced by a factor of ffiffiffiffiffiffiffi ffi 4=3 p when the inclination at the cantilever free end is measured instead of a direct deflection measurement by, e.g., interferometry 84 . As an example, the optical beam deflection method measures the cantilever inclination rather than the deflection. If the cantilever geometry deviates from that of an Euler-Bernoulli beam as in our case, this factor needs to be calculated individually.
Please note that the considerations above address thermal noise limitations only. The investigation of other contributions to the total noise was beyond the scope of this study.
MnNiGa sample fabrication and surface characterization. A polycrystalline MnNiGa sample was fabricated using the arc melting technique (for details see "MnNiGa sample" in the Methods section). To identify a suitable homogeneous sample location with the desired [001] crystal orientation (Fig. 9a, b) we performed electron backscatter diffraction (EBSD) at 15 kV Fig. 8 Tip-sample distance dependent perpendicular force and frequency shift. Perpendicular force F z and Δf 2 (frequency shift of the second flexural mode) graphs measured by sequential distance sweeps with the arrows pointing to the respective axes (F z measurement without any vibration excitation; second mode amplitude A 2 ¼ 4:8 nm during Δf 2 data acquisition with free vibration resonance frequency f 2 ¼ 33740 Hz). The inset covers the full z ts sweep range.
acceleration voltage (Fig. 9c). Furthermore, the chemical composition and the crystal homogeneity were investigated by surface-sensitive energy dispersive X-ray spectroscopy (EDX) L-edge mapping at 5 kV and by employing backscatter electron (BSE) micrographs, respectively. For MFM investigations, we selected a region with homogeneous (red-colored) EBSD [001] surface pattern, located next to a region with [110] crystal orientation, see Fig. 9c. This ensures that the easy magnetization direction in the studied [001] crystallite is perpendicular to the sample surface.
Quantitative B z and B 0 z measurements on bulk [001] MnNiGa. The MnNiGa sample was introduced into the MFM vacuum chamber 48 h prior to the measurement to ensure stable high vacuum conditions and to approach thermal equilibrium at room temperature. After calibration measurements and compensation of the contact-potential difference V cpd by application of an appropriate constant tip-sample DC bias, the main measurement was performed as a single-pass ð10 μmÞ 2 /ð512Þ 2 pixel scan at a controlled scanning distance z ts ¼ 70 nm and a tip velocity of 1:67 μm=s, resulting in a full frame completion time of 103 min. To guarantee a well-defined magnetic state, the MnNiGa sample and the magnetic probe tip were initially magnetized by a perpendicular À400 mT field. Note that the applied field is larger than the switching field of the magnetic nanowire probe (275 mT ± 10 mT). Figure 10a-d show simultaneously recorded maps of the [001] oriented uniaxial MnNiGa crystallite in an out-of-plane field of À270 mT pointing along the magnetic easy axis. Included are the raw z-piezo position data (Fig. 10a), the derived topography map ( Fig. 10b; representing the sum of the z-piezo output and the static cantilever deflection 4z data), the magnetic B z map (Fig. 10c), and the B 0 z map (Fig. 10d). Note, B z / 4z and B 0 z / 4f 2 . Several bubble-like features can be observed, accompanied by other types of magnetization patterns including magnetic stripe domains and field background generated by adjacent crystallites. The spatial resolution of the probe-sample interaction is identical for the B z and the B 0 z channel but in B 0 z data higher spatial frequencies are amplified compared to corresponding B z data 22 . As illustrated in Fig. 4a, our MFM approach is not sensitive to homogeneous field contributions, since FTF λ ð Þ approaches zero for large λ. Nevertheless, the finite size of the studied crystallite causes a small but measurable field variation associated with a large λ. Such contributions, however, only constitute a field background that is nearly constant within the length scale of the bubble diameters. Therefore, it does not affect our detailed bubble analysis discussed in the next section.
Evaluation of individual magnetic nano-domains. Finally, we analyze the observed circular nano-domains. Typical magnetic bubble domains form in a material with perpendicular anisotropy as minority domains in an applied external magnetic field. To determine the diameters d qs and the vertical extent of the considered magnetic bubbles, we evaluate angularly averaged B z and B 0 z profiles. An investigation of further properties such as depthdependent helicity 85 , bubble tube bending and merging 59 is not reported in this paper.
We limit our detailed evaluation to radially symmetric bubbles. In our analysis, we compare the measured and angularly averaged stray-field distributions of ten individual MnNiGa bubbles with stray-field characteristics of four idealized magnetization distributions: (a) a magnetic point dipole pointing in the direction perpendicular to the sample surface and exhibiting two degrees of freedom: the magnetic moment m s and its z Ã position, i.e., its depth below the surface; (b) a high aspect ratio magnetic type-I bubble consisting of a homogenously magnetized core that is separated from the magnetic environment by a Bloch wall, again with two degrees of freedom: diameter d qs and z Ã position of the bubble's top base; (c) a double monopole arrangement (Dirac string) with three degrees of freedom: the top monopole with moment q s and its z Ã position and the second monopole Àq s buried at a depth of z Ã þ h; (d) a single magnetic point monopole with two degrees of freedom: the magnetic monopole moment q s and its z Ã position.
In all cases, the considered magnetic structures are embedded in a homogenously magnetized sample with perpendicular magnetization. The magnetization direction of the considered magnetic structures is opposite to that of the homogeneously magnetized environment.
To assess the applicability of the respective magnetostatic models to our measured MFM data, we make use of the coefficient of determination R 2 (see Supplementary Note 5 and Supplementary  Fig. S1 for illustration). For the simulations and fitting functions, a MnNiGa saturation magnetization of M S ðMnNiGaÞ ¼ 470 kA=m was used which was determined by SQUID magnetometry measurements at room temperature. For the Bloch domain wall width, a value of 47 nm, as reported in literature 54 , was assumed.
Magnified images of an individual bubble (in Supplementary  Fig. S2 named bubble #3) are shown in the insets of Fig. 10c, d. In Table 1, the angularly averaged B z profile of bubble #3 (Fig. 10e) is compared to the stray-field distributions of the above introduced models. The fitting functions are given in Supplementary Note 4 (Supplementary Equations S13, S14, S16, and S17). The data evaluation indicates that the dipole-like description is inferior to the bubble domain and monopole models. The significance of the resulting R 2 differences is proven in Supplementary Note 5 by comparing to FEMMsimulations of related nanomagnetic configurations 86 . As far as the mutual, vertical pole distance h is concerned, the double pole fitting results in h ≥ 10 μm, which essentially corresponds to the monopole model with identical fitting results. For the bubble and monopole-type models, q s , d qs ;and z Ã of bubble #3 are given in Table 1.
Although both the magnetic bubble and the monopole models describe the measured magnetic field distribution well, there is no microscopic magnetization configuration with true point-like character. In all cases a non-vanishing Heisenberg exchange interaction forces a magnetization structure to be spatially extended. This implies that the measured B z profiles of the localized domains in [001] MnNiGa are well described by circular bubble domains.
An interesting aspect of the analysis is the determination of a non-zero depth z Ã of the bubble's top surface. This can be explained by a possible formation of flux closure domains near the sample surface. An additional source of indeterminate error  The coefficient of determination R 2 is estimated by comparing angularly averaged Bz data of bubble #3 to stray-field distributions of a simple point dipole model, a simulated type-I bubble model, a double point-monopole model, a single point monopole model, and an approach based on a pair of Bz and B 0 z values to evaluate properties of buried magnetic monopoles (as sketched in Fig. 3a). For the bubble and monopole models, the magnetic charge qs, the diameter d qs and the z-position of the bubble's top base, i.e., the depth of the (top) monopole, z Ã are given. Note that q s ¼ π 4 d 2 qs M S ðMnNiGaÞ.
in the specification of z Ã might be the possible existence of a magnetic dead layer on the sample surface. The magnetic point pole evaluation of the remaining bubbles based on angularly averaged B z and B 0 z MFM data is included in Supplementary Note 5 ( Supplementary Fig. S2 and Supplementary  Table S3). In any case, the fitting procedures using B 0 z data confirm magnetic bubbles with monopole-like stray-field characteristics and a significant vertical extent, presumably proliferating through the depth of the entire crystallite in agreement with a small-angle neutron scattering study 52 . Furthermore, all bubbles are buried as described by a non-zero z Ã . The simple monopole fits result in bubble diameters d qs and depths z Ã covering ranges from 120 nm to 200 nm and 70 nm to 110 nm, respectively. Thus, our study indicates that the bubble diameters of bulk MnNiGa bubbles are somewhat larger than the diameters of the assumed biskyrmions in thin MnNiGa lamellae prepared for transmission electron microscopy investigations 49,50 .

Conclusions
In this work we demonstrated a unique magnetic force microscopy approach that employs the combination of specially structured low-stiffness cantilevers, high-aspect-ratio magnetic nanowire tips with excitation and control of multiple vibration modes. This technique enables simultaneous single-pass distancecontrolled non-contact quantitative measurements of both a magnetic sample's stray-field and its stray-field gradient. As proof of principle we used this technique to investigate magnetic bubbles of MnNiGa bulk samples. We presented data that suggest that the magnetic bubbles have significant vertical extent and a buried bubble top base.
The presented method is well suited to study complex magnetic structures and can be easily implemented to enhance the capabilities of commonly used magnetic force microscopy equipment. This further extends the set of available microscopy methods for studying magnetic systems and may be applied to contribute to the characterization and understanding of a large variety of nano-and micron-sized magnetic stray-field patterns generated by ferromagnetic, superparamagnetic, and superconducting samples as well as electric current distributions.

Methods
Magnetic reference sample. A Pt(5 nm)[Pt(0.9 nm)/Co(0.4 nm)] 100 Pt(2 nm) multilayer with perpendicular (easy-axis) anisotropy is used as magnetic reference sample for the tip calibration procedure. Vibrating sample magnetometry reveals an effective multilayer magnetization M s ðCoPtÞ ¼ 510 kA=m and a perpendicular anisotropy constant K u ¼ 0:52 MJ=m 3 . The large ratio K u =K d ¼ 3:2 with the shape anisotropy constant K d ¼ 0:5 μ 0 M 2 s promotes the required perpendicular magnetization orientation within individual domains 33,81 .
MnNiGa sample. A polycrystalline MnNiGa sample was fabricated using arc melting under Ar atmosphere using an ingot with 99.99% nominal composition of its constituent elements Ni, Mn, and Ga. During the preparation process, the ingot was molten and recrystallized several times to ensure good homogeneity. The synthesis was accompanied by a final weight loss of less than 1%. Afterwards, the as-cast button shaped ingot was annealed at a temperature of 1073 K for 6 days in a sealed quartz tube under vacuum to further improve its homogeneity. Finally, the ingot was quenched in ice water to sustain the hexagonal phase. For the finished MnNiGa alloy sample M s ðMnNiGaÞ ¼ 470 kA=m was determined by SQUID magnetometry.

Data availability
The data analyzed and presented in this study are available from the corresponding authors upon reasonable request.

Code availability
All code developed in this study is available from the corresponding authors upon reasonable request.