Individual skyrmion manipulation by local magnetic field gradients

Magnetic skyrmions are topologically protected spin textures, stabilised in systems with strong Dzyaloshinskii-Moriya interaction (DMI). Several studies have shown that electrical currents can move skyrmions efficiently through spin-orbit torques. While promising for technological applications, current-driven skyrmion motion is intrinsically collective and accompanied by undesired heating effects. Here we demonstrate a new approach to control individual skyrmion positions precisely, which relies on the magnetic interaction between sample and a magnetic force microscopy (MFM) probe. We investigate perpendicularly magnetised X/CoFeB/MgO multilayers, where for X = W or Pt the DMI is sufficiently strong to allow for skyrmion nucleation in an applied field. We show that these skyrmions can be manipulated individually through the local field gradient generated by the scanning MFM probe with an unprecedented level of accuracy. Furthermore, we show that the probe stray field can assist skyrmion nucleation. Our proof-of-concepts results pave the way towards achieving current-free skyrmion control. Skyrmions are topologically nontrivial spin textures which could be used as energy efficient carriers of information in future memory devices, but first reliable and efficient control of their movement is required. Here, the authors demonstrate a method to control the movement of individual skyrmions by making use of the magnetic interaction between a sample and a magnetic force microscopy probe.

M ost ferromagnetic materials exhibit collinear magnetic ordering due to the Heisenberg exchange interaction between neighbouring spins. However, in magnetic systems with large spin-orbit coupling and lack of structural inversion symmetry, the Dzyaloshinskii-Moriya interaction (DMI) 1,2 can be strong enough to stabilise non-collinear spin textures, most notably magnetic skyrmions. Since their first discovery in chiral magnets 3,4 , magnetic skyrmions have attracted significant attention within the spintronics community due to their fascinating properties and potential applications 5,6 . Indeed, their nanoscale dimensions combined with topological stabilisation 7 and low threshold current densities required for motion [8][9][10][11][12][13] , make magnetic skyrmions promising candidates for storage and logic devices. In this technological context, particularly interesting are skyrmions found in systems composed of a heavy metal layer and an ultra-thin magnetic film with perpendicular magnetic anisotropy (PMA) [13][14][15][16][17][18][19] . Here the DMI arises due to the broken inversion symmetry at the interface and the large spin-orbit coupling in the heavy metal atoms, which mediate the interaction between spins in the ferromagnet. The interfacial nature of the DMI in these material stacks promotes the stabilisation of skyrmions of the Néel type 15 . Typical systems investigated utilise Co or CoFeB as magnetic films and Pt as heavy metal layer, while either an oxide or another heavy metal are used as top layers. In these material stacks, current pulses are usually employed to drive skyrmion motion, exploiting the spin Hall effect in the heavy metal layer and consequent spin-orbit torque [16][17][18][19] . However, despite the low values of depinning currents required, fast skyrmion motion is only attainable using large current densities, which are accompanied by undesired Joule heating effects. Several alternative methods have been proposed to drive skyrmions, such as magnetic field gradients [20][21][22] , temperature gradients 20,23 and magnons 24,25 . Among these, using magnetic field gradients has been very recently shown to be an effective, current-free method to move skyrmion lattices in a chiral magnet 22 . While a few theoretical studies have investigated this approach also for isolated skyrmions in thin films 21,26,27 , to the best of our knowledge no experimental work has yet been performed in these systems.
Here we show that isolated skyrmions can be manipulated individually and precisely using the local field gradient generated by a magnetic force microscopy (MFM) probe in ultra-thin films with PMA. The samples investigated consist of [X/CoFeB/ MgO] 15 multilayers, with X = Ta, W or Pt. CoFeB/MgO was chosen due to its excellent homogeneity in terms of both thickness and PMA strength, which results in a low pinning energy landscape 16,28 , making it a good candidate for reliable skyrmion motion. CoFeB is, however, a very soft magnetic material, meaning that imaging via MFM can be challenging due to potential perturbations caused by the probe stray field on the sample magnetic state. Indeed, magnetic stripe domains and their evolution into skyrmions, from which the DMI strength can also be extracted, have previously only been imaged at synchrotron facilities for these ultra-thin CoFeB layers, using non-interfering X-ray based methods 16,18,19,[29][30][31] . On the other hand, MFM is a more established tool for skyrmion investigations in harder materials 17,[32][33][34][35][36][37][38] . In this work we show that a technique as straightforward and readily available as MFM can be used on CoFeB-based stacks at room temperature and in ambient conditions to: (i) determine the DMI strength through measurements of magnetic stripe-domain widths, (ii) image magnetic skyrmions and manipulate their motion with a high level of accuracy, and (iii) assist skyrmion nucleation in field. The labyrinth-like demagnetised state and its evolution with perpendicular field is imaged for all sample stacks, and the DMI constant D is extracted from the average equilibrium domain width of the demagnetised state. Individual skyrmions are nucleated by the controllable shrinking of stripe domains in field. Since the probe-sample interaction is detrimental for imaging, we show that skyrmions can only be imaged in single pass MFM mode, with topography acquired prior to skyrmion nucleation. On the other hand, by exploiting the probe-sample interaction skyrmions can be moved in standard two pass MFM mode, through the local field gradient generated by the probe. This technique allows for the reversible and controllable manipulation of individual skyrmions, contrary to more conventional current-driven motion where skyrmions move in a collective manner. Finally, we show that the probe-sample interaction can also be used to create individual skyrmions at magnetic fields lower than their intrinsic nucleation field, by either collapsing stripe domains into skyrmions or by directly cutting through the stripe domains to generate skyrmions.

Results
Probe-sample interaction. Since MFM was used to investigate the magnetic configuration for all sample stacks, as well as to image, nucleate, and manipulate individual skyrmions, it is important to assess the effects of probe-sample interaction during MFM scanning. The principle behind standard MFM imaging relies on detecting the probe phase shift resulting from its interaction with the magnetic stray field generated by the sample. These phase images show field gradient, as the probe can be considered an harmonic oscillator and thus the change in phase is directly proportional to the first derivative of the force experienced by the probe during scanning. It follows that image perturbations caused by the probe-sample interaction are a common artefact in MFM imaging, particularly for relatively soft magnetic materials like CoFeB.
Probe-induced modifications are represented clearly in Fig. 1, which shows the demagnetised state for Ta/CoFeB/MgO imaged at different tip lift heights Δz and the corresponding image histograms, respectively. While the contrast of domains magnetised up and down increases significantly upon reducing Δz, the influence of the tip, which has a stray field of $21 mT at its apex (see Supplementary Note 1 and Supplementary Fig. 1), becomes more evident since the domains magnetised in the same direction as the local tip stray field expand at the expense of the others. This is reflected by an imbalance in the histogram peaks corresponding to the two domain types, which is visible already below Δz as high as 245 nm. Significant domain distortions caused by the probe were also found in W/CoFeB/MgO, albeit of a reduced amount, while a much smaller effect was visible for the harder Pt/CoFeB/ MgO sample even down to Δz = 70 nm (see Supplementary  Figs. 2 and 3). This is consistent with magnetometry measurements which indicate that the perpendicular saturation field is lowest for Ta/CoFeB/MgO ($27 mT), intermediate for W/ CoFeB/MgO ($36 mT) and highest for Pt/CoFeB/MgO ($94 mT) (see Supplementary Fig. 4).
Distortions of the magnetic configuration caused by the probe-sample interaction are evidently detrimental for imaging, and precautions have to be taken in order to image the samples in the demagnetised state as well as to image individual skyrmions. Later we show that the same probe-sample interaction can actually be employed advantageously to control skyrmion motion.
DMI measurements. The DMI constant D is known to reduce the domain wall energy density of Néel walls σ DW , following the relationship 39,40 : where A is the exchange stiffness and K eff is the effective ARTICLE COMMUNICATIONS PHYSICS | https://doi.org/10.1038/s42005-019-0242-5 perpendicular anisotropy constant. Accordingly, domain formation becomes more favourable in the presence of a finite DMI and measurements of the domain width can then be used to directly quantify D 16,41 . In particular, following the approach implemented by Woo et al. 16 , we extracted D from the equilibrium domain width of the demagnetised state, where up and down domains have the same average size.
Since the MFM probe alters significantly the relative domain population when scanning at Δz low enough to grant a good magnetic contrast (see Fig. 1), the demagnetised state of the samples was imaged while simultaneously applying an out-ofplane external bias field to counteract the perturbing effect of the tip. Only in this way it became possible to achieve an equal balance between up and down magnetised domains and image the demagnetised configuration as close as possible to an unperturbed state. This is schematically illustrated in Fig. 2a, b for the Ta/CoFeB/MgO multilayer: by scanning at Δz = 145 nm, an equal domain distribution is restored upon application of a perpendicular field H z $ 7 mT. In general, it was found that varying strengths of this bias field H z $ 7 ± 2 mT were needed to offset the local probe stray field at Δz = 145 nm, which was attributed to variations among tips or their ageing stage. This value is in good agreement with the stray field value calculated through the tip transfer method for the field seen by the sample at the bottom of the probe oscillations when the probe is at 145 nm, as explained in Supplementary Note 1. One of the MFM images for the demagnetised state of the Ta/CoFeB/MgO multilayer, that was analysed to quantify D, is shown in Fig. 2c. Similar images for the W-based and Pt-based stacks are displayed in Supplementary Fig. 5. The labyrinth domain structure observed in all samples is typical of multilayers with PMA and is also consistent with the slanted perpendicular hysteresis loops (see Supplementary Fig. 4).
Two different methods are used to extract the equilibrium domain width from the demagnetised images. The first approach relies on performing a 2-D fast Fourier transform (FFT) which, averaged across its radial profile, yields directly an amplitude peak at the reciprocal wavelength corresponding to the average domain periodicity 42 . The second approach uses an image segmentation process based on the random walker algorithm 43 to accurately separate pixels belonging to domains magnetised up or down. From this, the average domain width for each domain type is calculated as the ratio between the total area covered by one domain type and half of the total length of all domain edges detected (see Supplementary Note 2 and Supplementary Fig. 6). The average equilibrium domain widths extracted across several images via these two different methods are in excellent agreement. Using the approach implemented by Woo et al. 16 , the domain wall energy density σ DW can then be extracted from the average equilibrium domain width according to the analytical domain-spacing model applicable to multilayers 44 . Finally, once σ DW is known, D is easily calculated according to Eq. (1). For the exchange stiffness A a value of 10 pJ m −1 was assumed. The DMI magnitude extracted for the different samples is summarised in Table 1, together with the average equilibrium domain width, and the relevant magnetic parameters, which were measured through vibrating sample magnetometry (VSM).
The magnitude of D s extracted for Ta/CoFeB/MgO and W/ CoFeB/MgO multilayers is mostly in line with values reported for similar systems, in experiments that quantify the DMI through current-driven or field-driven domain wall motion, as well as from measurements of equilibrium domain patterns 30,[45][46][47] . On the other hand, lower estimations of D s in these systems have been extracted through Brillouin light scattering (BLS) measurements 48,51 . Regarding the Pt/CoFeB/MgO multilayer our measured D s is lower than what observed in similar multilayers through measurements of equilibrium domain patterns 18,29,31 , but in good agreement with BLS experiments performed on individual layers 48,52,53 . In any case, only for W/ CoFeB/MgO and Pt/CoFeB/MgO multilayers the DMI magnitude is significantly larger than the critical amount required to stabilise Néel domain walls 40 , with jD c j ¼ 2lnð2Þμ 0 M 2 s t=π 2 being 0.18, 0.15 and 0.33 mJ m À2 for Ta, W and Pt-based stacks, respectively. Furthermore, the DMI magnitude is always lower than the threshold value for spontaneous skyrmion generation 54 , indicating that skymrion lattices or individual skyrmions cannot be generated in these systems without some form of external excitation. Indeed, as will be shown in the next section, skyrmions could only be stabilised in W/CoFeB/ MgO and Pt/CoFeB/MgO multilayers by applying large perpendicular magnetic fields.
Skyrmion imaging and manipulation. For all investigated samples, application of a perpendicular magnetic field results in the contraction of stripe domains that are magnetised against the field direction. In systems with a large enough DMI these shrinking stripe domains can collapse into stable individual skyrmions upon increasing the perpendicular field towards saturation 14,16,19,33 . This was observed in W/CoFeB/MgO and Pt/ CoFeB/MgO multilayers, while the small DMI in the Ta-based stack prevented skyrmions from being stabilised in field. A sequence of images illustrating the process of skyrmion In order to image the nucleation of skyrmions, standard two pass MFM imaging was initially attempted. However, the magnetic interaction between probe and sample during first pass topography line scans, when the probe is closer to the sample, was found to be too strong to allow for skyrmion imaging using this MFM mode, even when using commercial low moment field probes. As a matter of fact, individual skyrmions could only be imaged in single pass MFM mode, having acquired the topography of the whole region of interest separately (for an evaluation of the stray field at the probe apex and its decay with Δz see Supplementary Note 1 and Supplementary Fig. 1, while for the difference between two pass and single pass mode see Methods and Supplementary Fig. 8). In particular, the topography was recorded at a field below skyrmion nucleation, where stripe domains are still present and able to better withstand strong perturbations caused by probe-sample interactions. The probe is magnetised in the same direction as the external field, as shown by the sketch in Fig. 3a. Upon increasing the field up to skyrmion nucleation, it becomes then possible to image unperturbed individual skyrmions by using the MFM in single pass mode, to record the phase signal only, at a set lift height. This is shown for Pt/CoFeB/MgO in Fig. 3b, where two skyrmions and a long stripe domain are imaged in single pass MFM at 71 mT, after having recorded the topography of the corresponding area at 50 mT.
Skyrmion nucleation occurred at 32 and 57 mT for W/CoFeB/ MgO and Pt/CoFeB/MgO multilayers, respectively. Measurements of an upper bound for skyrmions diameter in both samples, at their respective nucleation fields, alongside with results from micromagnetic simulations, are discussed in Supplementary Note 3 and Supplementary Figs. 9 and 11. The skyrmion nucleation field was found to be highly reproducible, also across different sample regions. Furthermore, the position of the individual skyrmions was observed to vary from one nucleation cycle to another, indicating that pinning does not play a significant role in the skyrmion nucleation process. This is likely to be enabled by the energy landscape typical of CoFeB which features lower pinning than Co, as shown for both domain wall motion 28 and skyrmion motion 16 .
Although detrimental for imaging skyrmions, the probe-sample interaction can actually be exploited to move them. This is illustrated for Pt/CoFeB/MgO in Fig. 3c, which was acquired in standard two pass MFM, with the probe magnetised along the external field direction. Here, skyrmions are moved downwards during topography line scans (in first pass), when the tip is in close proximity to the sample surface and the field gradient experienced by the skyrmions is larger. It is this field gradient that causes skyrmions to move further down at each subsequent topography line scan. In the meantime, magnetic contrast is recorded during each line scan in second pass, but what can be visualised is only the beginning portion of the skyrmions, before skyrmions are pushed further down by the following topography line scan. So, what eventually appears from the second pass image, can effectively been interpreted as a skyrmoin trajectory or 'footprint'. A comparison between Fig. 3b and Fig. 3d, which image the skyrmions before and after the two pass MFM scan, respectively, indicates that skyrmions have moved downwards by a few micrometres. Several repetitions of the single pass scan in Fig. 3d were acquired to ensure the stability of the new skyrmion position. It is interesting to observe that only skyrmions are moved by the MFM probe, while the position of the stripe domain visible at the top of the images in Fig. 3b-h is   Fig. 3b-h. Furthermore, it is important to note that the spacing δ between consecutive line scans along the slow scan axis should be kept small in order to attain reliable skyrmion motion. In our case, skyrmions could be moved effectively for δ $ 10-20 nm, while increasing the spacing to δ $ 40 nm or above resulted in stationary skyrmions or skyrmions that could only be moved a little before decoupling from the scanning tip. Arguably, a reliable motion cannot be expected if δ becomes comparable to the size of the skyrmion core. By using δ $ 10-20 nm in two pass MFM scans, individual skyrmion motion can be controlled with a high level of accuracy and the MFM tip can effectively become a tool to manipulate skyrmions in any desired manner. This is demonstrated in Fig. 4, where arbitrary geometrical patterns have been drawn with skyrmions, by moving them from random positions through a series of two pass MFM scans with different slow scan axis directions. It is worth mentioning that some recent works have reported on skyrmion creation by MFM 35,36,55 , but have not shown the possibility of achieving skyrmion manipulation as presented here. Skyrmion manipulation was instead demonstrated in PdFe atomic bilayers on Ir(111) using a scanning tunnelling microscopy (STM) tip 56 . In this work, lateral skyrmion motion was achieved by moving defects, to which the skyrmions are pinned, with the STM tip. This is different from our manipulation approach, where skyrmions interact directly with the MFM tip. Finally, we investigate the dynamics of two adjacent skyrmions when being pushed against each other by the MFM probe. As illustrated in Fig. 5a, two individual skyrmions are initially positioned by the MFM tip next to each other, with skyrmion S 1 being directly above skyrmion S 2 . When performing a two pass MFM scan from top to bottom, the local field gradient generated by the probe drives S 1 towards S 2 . The two skyrmions footprints are visible in Fig. 5b and indicate that S 1 , when being pushed downwards, avoids interacting with S 2 by suddenly jumping towards the right. Analogously, S 2 moves slightly towards the left. The final skyrmions positions are illustrated in Fig. 5c. The interaction between the two skyrmions seems therefore to be repulsive in nature, similarly to those reported elsewhere 57 .
Skyrmion nucleation at lower magnetic fields. We already mentioned that the skyrmion nucleation field for Pt/CoFeB/MgO was reproducibly found to be 57 mT. At this field, a few individual skyrmions are stabilised, while many long stripe domains are still present in a typical 10 10 μm 2 image (see Fig. 6a). In order to convert all stripe domains into skyrmions, fields higher than 75 mT need to be applied. Here we show that the probe stray field can be used to promote skyrmion nucleation by assisting domain contraction without increasing the external field. This is illustrated through the sequence of MFM single pass images shown in Fig. 6. Between each image a topography scan is performed (not shown) with the direction of the slow scan axis as indicated. Many elongated domains are contracted into skyrmions after only one topography scan (see Fig. 6b), but repeating the topography scan in the same direction does not change significantly further the magnetic configuration, as several long domains remain unvaried (see Fig. 6c, d). However, it is found here that even these more persistent domains can collapse into skyrmions by changing the slow scan axis direction in such a way that the scanning tip moves against one of the domain extremities. For instance, the topography between Fig. 6e and f is scanned from left to right and leads to the contraction of stripe domain D 2 into skyrmion S 2 , while leaving stripe domain D 1 almost unaffected. To convert D 1 into a skyrmion, the topography has to be scanned firstly from right to left, with the tip moving against the extremities of D 1 , leading to a shorter vertical stripe domain (see Fig. 6g), and then from bottom to top, finally resulting in skyrmion S 1 (see Fig. 6h). Thus, this method allows to convert all stripe domains into skyrmions without increasing the external field, and more efficiently than what would otherwise be possible by keeping the same direction of the slow scan axis. It is c, e, g Skyrmions are moved by the probe stray field during the topography line scans of a standard two pass MFM scan (also with Δz = 145 nm), which reveal the trajectory followed by skyrmions in the form of actual skyrmion footprints. These scans were stopped before reaching the end only to prevent skyrmions from moving outside of the selected region. The directions of fast (with trace and retrace line scans) and slow MFM scan axes are marked in the individual images. The spacing between consecutive line scans along the slow scan axis is δ $ 20 nm. As indicated by the sketch in (a), the tip is magnetised along the external field direction for all MFM images, thus opposing the skyrmions core magnetisation important to note that skyrmions in Fig. 6 do not move during topography scans in the same reliable and controllable manner that was shown in Fig. 3. This is due to the fact that a larger line spacing δ $ 40 nm was used in Fig. 6, resulting in small or negligible skyrmion motion, as already explained.
Finally, we discuss the possibility of nucleating skyrmions by 'cutting' the stripe domains with the tip stray field. To this end, it is important to note that the switching field of the tip magnetisation was measured to be $40-50 mT (see Supplementary Fig. 10), which means that when applying external fields larger than this value, the tip is always magnetised along the field direction. This is the case for all the MFM images presented for Pt/CoFeB/MgO thus far (i.e. Figs. 3-6), whereby external field and tip stray field are aligned in the same direction, opposing the magnetisation of the skyrmions core. Skyrmion nucleation through domain slicing was first attempted in these conditions for Pt/CoFeB/MgO, by scanning the tip several times in topography mode across a few stripe domains, under an applied perpendicular field H z = 57 mT. Despite external field and tip stray field being aligned with each other, no visible modification to the sample magnetic configuration was obtained. However, for the softer W/CoFeB/MgO only one topography scan in a modest external field ($10 mT) was sufficient to significantly alter the labyrinth-like magnetic state. This is shown in Fig. 7 for two possible magnetic configurations of the tip. With the tip magnetisation aligned against the external field, the topography scan mainly leads to a rearrangement of the labyrinth domains (see Fig. 7a), with no significant skyrmion nucleation. On the other hand, when the tip magnetisation is in the same direction as the field, so that the tip stray field adds to the external field, some  The dashed arrows are guides for the eye indicating the trajectories followed by the skyrmions during the probe-driven interaction. b Skyrmion S 1 is pushed against skyrmion S 2 by the probe stray field during the topography line scans of a standard two pass MFM scan (also with Δz = 145 nm). Skyrmion S 1 avoids the collision with skyrmion S 2 by jumping abruptly towards the right. The directions of fast and slow scan axes are marked by the arrows. The tip is magnetised along the external field direction, thus opposing the skyrmions core magnetisation of the dark stripe domains are effectively sliced through by the scanning tip and individual skyrmions are created (see Fig. 7b), which is also consistent with previous work 35 .

Discussion
While MFM is well established as a tool to measure DMI and investigate nanoscale magnetic skyrmions in Co-based multilayer stacks 17,[32][33][34][35]37 , similar studies for CoFeB-based systems have involved using less widely available techniques, such as scanning transmission X-ray microscopy (STXM) 16,18,30,31 , magnetic transmission soft X-ray microscopy (MTXM) 16,19 or X-ray holography 29 . MFM imaging of soft magnetic materials like CoFeB can actually be very challenging due to image deformations caused by the probe-sample magnetic interaction. However, in this work we performed extensive MFM imaging of CoFeB/ MgO multilayers, and we demonstrated that MFM can be used as a versatile technique to ascertain DMI and study magnetic skyrmions also in these systems.
We determined the DMI from measurements of equilibrium domain widths in the demagnetised state. To minimise the effect of tip-induced perturbations, commercial low moment tips were used, in combination with large lift heights. A small bias field was applied during imaging of the demagnetised state to compensate for the imbalance between up and down magnetised domains due to the probe stray field. Image analysis was then performed to extract the average equilibrium domain width with two different approaches, which yielded results in excellent agreement, and the DMI was then quantified through analytical domain-spacing model for multilayers.
Furthermore, we explored how the probe-sample interaction can be beneficially exploited to attain individual skyrmion manipulation, as well as to assist skyrmion nucleation in field. In both cases, the local field gradient generated by the scanning MFM tip represents a tool not only to image skyrmions, but also to actively interact with them and manipulate their position controllably. Interestingly, while skyrmion imaging is only possible in single pass MFM scans, with the topography of the corresponding area acquired separately, standard two pass MFM provides a means to drive skyrmions (during topography scan, in first pass) and image their trajectory (during lift scan, in second pass). Skyrmion motion by this technique is reversible and highly controllable, also enabled by the low pinning landscape of CoFeB. Since the probe stray field opposes the direction of the skyrmions  core magnetisation, we found that skyrmions move away from the scanning probe, in a direction approximately perpendicular to the probe fast scan axis (i.e. the trace/retrace axis). Furthermore, we show that the tip stray field in topography scans can also be used to promote the creation of skyrmions at fields lower than their intrinsic nucleation field by either collapsing labyrinth-like domains or directly cutting through them. Overall, our proof-ofconcept results open alternative ways for the creation and manipulation of individual skyrmions in multilayers with perpendicular magnetic anisotropy.  15 /Ta (5). They were grown at room temperature by d.c. magnetron sputtering in a sputter system with a base pressure of 5 × 10 À8 mbar. After growth, the stacks were annealed in vacuum for 2 h at a temperature of 300°C, for sample (i), and 400°C for samples (ii) and (iii). The 5 nm Ta capping layer was used to prevent sample oxidation. For sample (iii), the 5.7-nm thick Ta seed layer was used to improve adhesion between film and substrate. The hybridization of Fe and Co is different with Pt, W and Ta and the DMI and perpendicular magnetic anisotropy depend on the CoFeB composition and thickness 30 . We have developed different stacks with the correct properties to host skyrmions 58,59 .

Methods
Magnetometry. The magnetic properties of the multilayers were studied using vibrating sample magnetometry (VSM) at room temperature. Both in-plane and out-of-plane magnetic hysteresis loops were measured (see Supplementary Fig. 4). All investigated samples exhibit PMA, with a saturation field that increases in going from sample (i) to sample (iii). The saturation magnetisation M s and the effective perpendicular anisotropy constant K eff , which was extracted from the in-plane saturation field, are listed in Table 1.
Magnetic force microscopy imaging. MFM imaging of the multilayers was performed at room temperature with a NT-MDT Ntegra Aura scanning probe microscope (SPM). The system is fitted with an electromagnet which allows to apply out-of-plane magnetic fields up to $115 mT during scanning. Low moment tips (NT-MDT MFM-LM) were chosen in order to minimise probe-sample interaction. MFM imaging of the demagnetised state of the samples was carried out in standard two pass tapping mode, whereby for each scanned line topography is acquired during the first pass, while the MFM phase signal of that line is acquired during the second pass at a set lift height Δz. MFM imaging of individual skyrmions was instead performed in single pass tapping mode, which allows to acquire topography for the whole scanned area first, followed by the MFM phase signal for the same area at a set lift height Δz. By working in single pass mode, it is then possible to measure topography and phase separately, which was essential to image skyrmions. Furthermore, single pass mode allows to acquire the phase signal at different magnetic fields without the need to lower the probe close to the sample to retake topography. For the difference between two pass and single pass mode see also Supplementary Fig. 8. Finally, individual skyrmion manipulation was attained using standard two pass tapping mode: skyrmions move due to the field gradient generated by the probe during first pass, when topography is recorded, while their trajectory is imaged during second pass.
Micromagnetic simulations. We simulated the skyrmion dynamics with the open source micromagnetic code mumax3 60 . This software solves the Landau-Lifshitz-Gilbert equation using a finite difference approach and includes the DMI interaction. We performed simulations for a multilayer of Pt/CoFeB/MgO with the parameters shown in Table 1. To simplify the magnetic stack, we used the effective medium approach 16,41 to simulate the magnetic and nonmagnetic layers as a single effective magnetic layer. The simulated sample had dimensions 1024 × 1024 × 84 nm 3 and a discretization cell of 2 × 2 × 5.6 nm 3 . Using this approach, we calculated the skyrmion diameter as a function of applied field (see Supplementary Note 3 and Supplementary Fig. 11) and the several realizations of the skyrmion motion under the field created by the tip dipole (see Supplementary Note 4 and Supplementary  Fig. 12).

Data availability
The datasets generated during and/or analysed during the current study are available from the corresponding author upon reasonable request.