In order to design magnetoresistive sensor devices, one needs detailed insight into microstructure and mechanisms that are involved. For this purpose, a strong entanglement of experimental investigation and simulation is required. Here, we present numerical achievements alongside experimental methods with the objective of improving the comprehension of fundamental mechanisms leading to beneficial magnetoresistive properties.
2.1. Towards a Rational Design Tool for Magnetoresistive Sensor Devices by Means of Hybrid Molecular Dynamics and Spin Dynamics Simulations
The spatial arrangement of the magnetic nanoparticles in the gel matrix is mainly responsible for the magnetoresistive properties of the sensor device. Hence, detailed information about the structuring process driven by external magnetic fields is necessary. Classical spin dynamics (SD) is a numerical method which allows to compute the static and dynamic magnetic properties of a spin system as a function of temperature. It can be applied to microscopic and mesoscopic ensembles of classical magnetic moments. In the following we assume single domain nanoparticles,
i.e., an effective magnetic moment is assigned to each particle, related to its diameter and saturation magnetization [
6]. In addition to that, we assume that exchange interactions between particles do not play a role and therefore only magnetic dipole dipole interaction is taken into account. We used a SD algorithm where the heat bath contact is modeled using a Langevin approach. Details concerning the algorithm can be found in [
7]. Such SD simulations assume fixed positions of the magnetic moments and determine the orientations of the magnetic moments by solving the spin equations of motion. However, in the liquid gel, the magnetic forces due to inter-particle interactions and external fields lead to particle movement that has to be addressed by a second simulation method, namely molecular dynamics (MD). MD is a numerical method that computes trajectories of classical point particles,
i.e., dimensionless particles which are set in motion by forces that are derived from classical potentials. The simulation of nanoparticles that interact by magnetic dipole dipole interaction and which are free to move due to magnetic forces is only possible by following a hybrid molecular and spin dynamics approach. For the MD part of our simulation, the highly-specialized software package HOOMD-blue [
8] has been chosen because of its remarkable range of available functions. In order to include the magnetic dipole dipole interaction in the MD code, a potential energy contribution according to Equation (1) is evaluated for every particle pair in the simulated system for every step in time:
Within this equation,
denotes the distance between the magnetic moments
and
of the particles 1 and 2, while
represents the vector that joins the centers of particles 1 and 2. The usage of point particles in HOOMD-blue in combination with magnetic dipole dipole interaction poses a general problem for our hybrid approach that has to be addressed first. The magnetic dipolar forces exhibit a singularity that occurs when the particle distance approaches zero. In order to avoid such singularities we have implemented the so-called Weeks-Chandler-Anderson potential [
9] that provides a short-range repulsion so that the simulated particles behave like hard spheres that interact by their magnetic dipolar fields:
Therein, ε represents the depth of the attractive potential well, given in energy units and σ, given in units of distance, denotes the distance at which the potential becomes zero. The potential is cut off at a distance of . By choosing , the simulated particles behave like hard spheres with given diameters. Similarly to the magnetic dipole dipole interaction, the WCA potential is evaluated for every particle pair in every step in time.
In general, a hybrid MD and SD approach requires to solve the equations of motion simultaneously. For our hybrid approach we assumed that the translational degrees of freedom for each nanoparticle are much slower than its magnetic degrees of freedom,
i.e., the change of the orientation of the magnetic moment. This allows us to decouple the calculation of the degrees of freedom reminiscent of the Born-Oppenheimer approximation [
10]. Starting with an initial configuration consisting of the positions and magnetic moment orientations for each nanoparticle, we first use SD in order to compute the magnetic ground state for that configuration. After that, the new magnetic moment orientations are used in the MD part to compute the corresponding forces and integrate the equations of motion for a predefined number of time steps, resulting in new positions which are again passed over to the SD calculation. These steps are repeated alternately until the maximum number of time steps is reached. We applied our method to a sample of 192 cobalt nanoparticles with diameters of 20 nm which is chosen to closely mimic the experimental situation which is depictured in
Figure 1. In our model, the particles are placed in an area of 350 nm × 350 nm. Therefore, 2 × 10
7 hybrid MD/SD time steps have been performed with one full SD run every 1 × 10
5 time steps. In zero external magnetic fields, the formation of chain fragments and small islands containing magnetic vortices can be observed. During the structuring process, the magnetic dipolar energy and the total magnetization are gradually reduced until a minimum has been reached.
Figure 1.
GMR effect of cobalt particles embedded in agarose. (a) Particle size distribution of Co core/conductive C-shell nanoparticles used for magnetoresistive measurements in agarose matrix utilizing an AC mode; (b) Experimentally achieved GMR effect as a function of particle concentration c which is defined by the ratio of particle mass over gel mass. The GMR effect amplitude is given in the ratio of ΔR to R with R, the electrical resistance, which relates to an effect amplitude of about 10% to 20%. The reliability of an effect amplitude of above 20% is given in (c) where the GMR effect of a sample has been repeated for 20 times.
Figure 1.
GMR effect of cobalt particles embedded in agarose. (a) Particle size distribution of Co core/conductive C-shell nanoparticles used for magnetoresistive measurements in agarose matrix utilizing an AC mode; (b) Experimentally achieved GMR effect as a function of particle concentration c which is defined by the ratio of particle mass over gel mass. The GMR effect amplitude is given in the ratio of ΔR to R with R, the electrical resistance, which relates to an effect amplitude of about 10% to 20%. The reliability of an effect amplitude of above 20% is given in (c) where the GMR effect of a sample has been repeated for 20 times.
Figure 2.
Result of hybrid SD-MD simulation in zero external fields for a test structure containing 192 Co particles. The formation of chains and islands containing magnetic vortices can be found which is in agreement with the experimental observations. The total magnetic energy decreases from Emag = −1.6512E−13 J for the initial state (a) to Emag = −3.34852E−15 J for the final state (b).
Figure 2.
Result of hybrid SD-MD simulation in zero external fields for a test structure containing 192 Co particles. The formation of chains and islands containing magnetic vortices can be found which is in agreement with the experimental observations. The total magnetic energy decreases from Emag = −1.6512E−13 J for the initial state (a) to Emag = −3.34852E−15 J for the final state (b).
Figure 3.
Calculated GMR curve for our SD/MD-optimized sample of 192 Co particles. The magnetization curve (black dots) has been obtained by means of spin dynamics simulations after performing a complete hybrid SD-MD simulation. From the magnetization data, the qualitative GMR curve (blue squares) can be extracted according to Equation (3). For the sake of simplicity, the GMR effect amplitude is assumed to be 20%, according to the experimental results shown in
Figure 1.
Figure 3.
Calculated GMR curve for our SD/MD-optimized sample of 192 Co particles. The magnetization curve (black dots) has been obtained by means of spin dynamics simulations after performing a complete hybrid SD-MD simulation. From the magnetization data, the qualitative GMR curve (blue squares) can be extracted according to Equation (3). For the sake of simplicity, the GMR effect amplitude is assumed to be 20%, according to the experimental results shown in
Figure 1.
This situation corresponds to a “freezing” with respect to the positions and the orientations of the magnetic moments as well. The initial configuration and the final configuration of the particles are shown in
Figure 2. In analogy to experiments, the next step is to dry out the gel matrix in order to fix the positions of the nanoparticles. In the simulation, this can be done by simply switching off the MD part and just perform SD simulations. In order to compare to experimental results [
5] we have to determine the magnetoresistive properties. Experimentally, the electrical resistance of such structures can be measured with a four-point-probe setup [
5] and the GMR is determined by measuring the resistance for different external magnetic fields. In the framework of SD simulations the GMR characteristics can be derived from the total magnetization according to [
11,
12]:
with
M, the magnetization and
MS, the saturation magnetization. The GMR effect amplitude
can be calculated by means of quantum mechanical techniques [
13] or it can be determined experimentally. We matched
to our experimental data.
Figure 3 shows a calculated magnetization curve for the structure under consideration obtained by means of classical spin dynamics simulations and the resulting GMR curve according to Equation (3).
2.2. Efficient Calculation of Low Energy Configurations of Frustrated Dipolar Systems
After gelatinization of the gel matrix, the particle positions are fixed in space. The magnetic state of a rigid configuration can be investigated by means of SD. Simulations at T = 0 K reveal a sophisticated behavior due to the dipolar coupling between the particles and their spatial disorder. Instead of one unique ground state, a multitude of low energy configurations has been found. Here, the inherent magnetic frustration of the particle arrangements due to the dipole dipole interaction poses computational challenges reminiscent to those of calculating the properties of large-scale spin glass systems [
14]. A standard relaxation algorithm fails to predict a global energy minimum because the system gets trapped in local energy minima that are separated by energy barriers. In order to solve that problem we have implemented a demagnetization routine [
6] that has been successfully applied to the numerical investigation of artificial spin ice [
15] which belongs to a class of spin systems that is very similar to our system. In our demagnetization routine the sample is exposed to a sinusoidally varying, rotating, and damped magnetic field. Using this technique one can overcome the energy barriers which occur due to strong inter-particle interactions with the result that lowest energy configurations can be found much more efficiently as compared to the standard relaxation algorithm. Results of simulations obtained by a standard relaxation and simulations using our demagnetization routine are shown in
Figure 4. For the demagnetization routine, 30 sinusoidal oscillations, 7 azimuthal, and 3 polar turns have been applied during a simulation run of 1 × 10
−7 s while damping the magnetic field from 0.05 T to 0 T.
Figure 4.
Evaluation of a total of 1000 simulations following a simple relaxation procedure and 1000 simulations using our demagnetization protocol resulting in lower total magnetic energies and a narrower distribution (a); Application of a demagnetization protocol leads to energetically favorable configurations with larger vortices (c) compared to the simple relaxation procedure (b); The colors of the arrows in (b,c) indicate the x-component of the magnetic moment of each particle.
Figure 4.
Evaluation of a total of 1000 simulations following a simple relaxation procedure and 1000 simulations using our demagnetization protocol resulting in lower total magnetic energies and a narrower distribution (a); Application of a demagnetization protocol leads to energetically favorable configurations with larger vortices (c) compared to the simple relaxation procedure (b); The colors of the arrows in (b,c) indicate the x-component of the magnetic moment of each particle.
2.3. Combining Experiment and Simulation: Spin Dynamics Simulations of 3D Reconstructed Sensor Systems
The GMR effect amplitude depends on several parameters including the volume fraction of particles along the current path and the overall particle concentration. The microstructural arrangement of magnetic nanoparticles can be investigated using a dual system composed of a focused ion beam (FIB) and a scanning electron microscope (SEM). For a detailed study, the dual method was applied along with SD simulations. With the dual system, SEM images are taken while cutting away thin slices of the sample successively. Afterwards a three-dimensional reconstruction of the SEM images provides detailed topological information of the magnetic nanoparticle distribution, whereas the magnetoresistive properties can be determined by a four-point probe method with the contacting needles configured in-line as proposed in [
5]. To gain insight into the mechanisms that control the magnetoresistive properties, detailed information about the microscopic magnetic structure is necessary, which can be obtained by performing SD simulations. This method has been applied to two samples consisting of cobalt nanoparticles in an agarose matrix. In order to perform the dual FIB-SEM experiments, the sample was freeze dried whereby the agarose matrix is losing its conductivity. Thus, the sample was made conductive by placing a thin gold layer on top. By extracting the particles’ positions and sizes from 3D reconstructions and using the information as input for SD simulations, it can be shown that potential sensor configuration are composed of particle chains, vortices and ferromagnetically ordered areas, as shown in
Figure 5. Because the structure exhibits strong geometric frustration, the numerical technique described in
Section 2.2 has been applied. By comparing the magnetoresistive and magneto-structural information, it can be shown that the strong ferromagnetic order of closely packed particles decreases the GMR effect amplitude, whereas antiferromagnetically ordered chains increase the GMR effect (details on magnetoresistive characteristics will be published elsewhere). Moreover, information about the GMR characteristics can be extracted from simulation data, as described in
Section 2.1. As a result, the combination of the experimental dual FIB-SEM technique and numerical spin dynamics simulations provides a powerful tool for the investigation of magnetoresistive effects in detail and it is crucial for the rational development of future sensor devices. By combining SD and MD simulations to one hybrid method and making use of improved SD for frustrated systems, together with experiments done by the dual FIB-SEM method, a very elaborate tool chain has been created. With this, one is able to gain detailed insight in the formation of the sensor’s microstructure and can investigate the role of different parameters on its magnetoresistive properties.
Figure 5.
Results of SD simulations of two experimentally investigated structures, one prepared in a homogeneous external magnetic field (a) and one prepared in a rotational field (b). The particle positions and sizes were obtained by 3D reconstruction of images taken by a dual FIB-SEM system. The length and thickness of the arrows indicates the size of the magnetic moment, whereas the color indicates the x-component of the magnetic moment. The simulations reveal vortex-like structures (blue), ferromagnetically ordered (green), and chain-like (red) areas.
Figure 5.
Results of SD simulations of two experimentally investigated structures, one prepared in a homogeneous external magnetic field (a) and one prepared in a rotational field (b). The particle positions and sizes were obtained by 3D reconstruction of images taken by a dual FIB-SEM system. The length and thickness of the arrows indicates the size of the magnetic moment, whereas the color indicates the x-component of the magnetic moment. The simulations reveal vortex-like structures (blue), ferromagnetically ordered (green), and chain-like (red) areas.