Shear Elasticity of Magnetic Gels with Internal Structures

We present the results of the theoretical modeling of the elastic shear properties of a magnetic gel, consisting of soft matrix and embedded, fine magnetizable particles, which are united in linear chain-like structures. We suppose that the composite is placed in a magnetic field, perpendicular to the direction of the sample shear. Our results show that the field can significantly enhance the mechanical rigidity of the soft composite. Theoretical results are in quantitative agreement with the experiments.


Introduction
Magnetic gels are composites of nano-and micron-sized magnetic particles in soft polymer matrixes. The combination of a rich set of physical properties of the polymer and magnetic materials is very valuable for many progressive industrial, bioengineering, and biomedical applications [1][2][3][4][5][6][7]. In part, for address drug delivery, for industrial and biological sensors [8][9][10][11][12][13][14], for the construction of soft actuators and artificial muscles [2,15], and for regenerative medicine and tissue engineering . An overview of the works on magnetic gel synthesis and their biomedical applications can be found in [23].
One of the remarkable properties of magnetic gels is their ability to change, under the action of an external magnetic field, their microstructure, magnetic, mechanical, and other microscopic properties, size, and shape. This provides the opportunity to control, with the help of the field, mechanic behavior, transport, and electrical processes in these systems, and this possibility presents a significant advantage for biosensoric, tissue engineering, and other biological applications [8,9,20,23,39,41].
During the magnetic gels' synthesis, the particles are usually embedded in the liquid polymer, and their spatial distribution is fixed after the composite gelation. If the host polymer is cured without an external magnetic field, the particles, as a rule, are distributed more or less homogeneously and isotropically. If the composite is polymerized under the field (field of polymerization), the particles form various anisotropic structures, elongated in the field direction. The appearance of these internal structures significantly changes the sensitivity of the gels to mechanic, electrical, magnetic, and other external impacts, and changes the kinetics of the internal transport phenomena and chemical reactions, the rate of cell proliferation, and other phenomena in these systems. This opens perspectives of the tunable synthesis of the magnetically controlled sensors, and scaffolds for the growth of cell tissues with willing structure and properties, artificial muscles, and other materials for biological and industrial applications.
The simplest kind of the internal structures that are formed by the magnetic particles in liquid media are liner chains, where particles are bounded as "head to tails" by the forces of magnetic attraction. These structures appear in the systems with a low and moderate volume concentrations of magnetic particles (usually in the range of 10-15%). Some photos of these chains can be found, for example in [42,43]. In the case of the higher concentrations, the particles can form topologically more complicated branched, net-like, bulk and other structures (see, for example, [44][45][46][47][48][49][50]).
The effect of the magnetic field on the elastic shear modulus of the gels with a homogeneous and isotropic (gas-like) spatial distribution of magnetic particles has been theoretically studied in refs. [51,52]. It was shown that a magnetic field that is applied perpendicular to the macroscopic shear of the composite, enhances the composite elastic modulus. At the same time, the experiments [53] show that the rigidity of the composites with the chains aligned perpendicularly to the shear, is significantly more than that of the systems with a chaotic distribution of the particles. Since the micromechanic (on the level of the particles and their aggregates) situation in the composites with the heterogeneous aggregates is significantly different from that in the homogeneous systems, the microscopic analysis of macroscopic properties of the composites with internal heterogeneous structures requires the development of a special theoretical approach.
The aim of this work is to theoretically study the effect of an external magnetic field on the shear elastic modulus of magnetic gels with internal chain-like structures, that are formed by magnetizable non-Brownian spherical particles. Physically, this means that the size of the particles are supposed to range from several tens of nanometers to microns. The particles of this size are very often used for the preparation of magnetic gels for bio-medical applications.
We take into account that the chains appear at the stage of the matrix polymerization under the action of an external magnetic field. That is why all chains are parallel to this field. We suppose that the actual magnetic field H has the same direction as the field of polimerization, i.e., that the field H is parallel to the chains. We consider the situation when the macroscopic shear of the sample is perpendicular, whereas the gradient of the shear is parallel to the chains. The length of the chains is supposed to be much less than the size of the sample.

Physical and Mathematical Model
For maximal simplification of the mathematical part of the problem, we will suppose that the particles are identical. Like in [54,55], we will neglect fluctuations of the chain's shape and will consider them as ideally straight aggregates, aligned along the applied magnetic field H. This model of the chain is illustrated in Figure 1.
The typical size of the cell of the polymer matrix in ferrogels is several nanometers, and the size of the particles vary from several tens of nanometers to microns. Thus, the particles are much larger than the gel cell. That is why we will consider the host polymer as a continuous medium with respect to the particles.
We will restrict ourselves by the analysis of small deformations of the composites and will suppose the linear relations between the mechanic stress and deformations in the matrix.
We will also neglect any interactions between the chains. This approximation is based on the results of ref. [56], which shows that the effects inside the chains play a dominant role in the formation of macroscopic properties of the composites, compared with the effects of the interchain interaction.
For mathematical definiteness, we will suppose that the chain consists of an odd number of particles. This assumption is not of principle for the physical analysis. We will restrict ourselves by the analysis of small deformations of the composites and will suppose the linear relations between the mechanic stress and deformations in the matrix.
We will also neglect any interactions between the chains. This approximation is based on the results of ref. [56], which shows that the effects inside the chains play a dominant role in the formation of macroscopic properties of the composites, compared with the effects of the interchain interaction.
For mathematical definiteness, we will suppose that the chain consists of an odd number of particles. This assumption is not of principle for the physical analysis.
Let us denote the mean vector of a material point displacement in the composite as u. In the coordinate system, shown in Figure 2, the vector u is the mean displacement in the composite (i.e., the displacement at the infinitive distance from the chain). The vector u can be presented as: = where γ is the mean shear of the system. We will consider small shear deformations of the composite, which means that the strong inequality ≪ 1 is held. Let the total number of the particles N in a chain be = 2 + 1, where n is an integer. In the framework of the used approximations, the equations of the stationary displacement of the particles in the chain can be presented in the following form [54]: Let us denote the mean vector of a material point displacement in the composite as u. In the coordinate system, shown in Figure 2, the vector u is the mean displacement in the composite (i.e., the displacement at the infinitive distance from the chain). The vector u can be presented as: u x = γx where γ is the mean shear of the system. We will consider small shear deformations of the composite, which means that the strong inequality γ 1 is held. We will restrict ourselves by the analysis of small deformations of the composites and will suppose the linear relations between the mechanic stress and deformations in the matrix.
We will also neglect any interactions between the chains. This approximation is based on the results of ref. [56], which shows that the effects inside the chains play a dominant role in the formation of macroscopic properties of the composites, compared with the effects of the interchain interaction.
For mathematical definiteness, we will suppose that the chain consists of an odd number of particles. This assumption is not of principle for the physical analysis.
Let us denote the mean vector of a material point displacement in the composite as u. In the coordinate system, shown in Figure 2, the vector u is the mean displacement in the composite (i.e., the displacement at the infinitive distance from the chain). The vector u can be presented as: = where γ is the mean shear of the system. We will consider small shear deformations of the composite, which means that the strong inequality ≪ 1 is held. Let the total number of the particles N in a chain be = 2 + 1, where n is an integer. In the framework of the used approximations, the equations of the stationary displacement of the particles in the chain can be presented in the following form [54]: Let the total number of the particles N in a chain be N = 2n + 1, where n is an integer. In the framework of the used approximations, the equations of the stationary displacement of the particles in the chain can be presented in the following form [54]: (1) Here, G 0 is the shear modulus of the matrix, d is the diameter of the particles, f (m) i is the force of the magnetic interaction of the i-th particle with the neighboring particles of the chain, and for the central particle i = 0. The equality u i = −u −i for 0 < i ≤ n follows from the symmetry of the problem.
For convenience, we introduce the dimensionless magnetic force f (m) i and dimensionless displacement u i : By using these notations, one can rewrite Equation (1) as: Under the assumption that γ 1, the inequality u i 1 is held.
The magnetic force f , we need to determine the magnetic moments of the particles in the chain. Strictly speaking, the value of a particle moment depends on the number i of the particle position in the chain. However, analysis [55] shows that approximation, where magnetic moments of all of the particles in the chain are supposed identical, leads to not significant deviations from the strict approach. That is why we will use the simplest approximation of identity of the particle moments in the chain.
When the composite experiences the macroscopic shear deformation, the axis of the chain deviates from the z-axis, as is illustrated in Figure 2. Because of the mutual magnetization of particles, the vector of the particle magnetic moment will also be deviated from this axis. Therefore, both components M x and M z of the vector M of the particle magnetization will take place in the deformed composite.
We estimate the magnetic force f (m) i by using the nearest-neighbor approximation, taking into account the magnetic dipole-dipole interaction only between the neighbor particles in the chain. By using the well-known relation for the force of the dipole-dipole interaction (see, for example, [57]), after simple but cumbersome transformations, in the linear approximation with respect to the displacements u i , one can get: , Here, f i,i±1 is the force of the magnetic interaction between the i-th particle in the chain with the neighbor particle, M s is the saturated magnetization of the particle's material, β is the parameter which defines the ratio of the energy of magnetic interaction between two magnetically saturated particles to the energy of elastic deformation of the matrix, and µ 0 is the vacuum magnetic permeability.
The main problem now is to estimate the dimensionless components M x and M z . The strict solution of a problem of determination of magnetic moments of two closely situated magnetizable particles has not been obtained in literature because of not overcoming mathematical complexity. Here, we use the approach [58], where each particle is considered to be situated in a uniform magnetic field H (e) , consisting of the external, with respect to the particle's field H (i.e., the mean field in the sample), and the field created by the other particle in the center of the first one. It is supposed that the magnetization of the particle obeys to the nonlinear Frolich-Kennelly relation [59]: Here, χ p is the initial magnetic susceptibility of the particle material, and H (i) is the magnetic field inside the particle. For the spherical particle the last field, this can be found from the general relation [60]: Taking into account the magnetic interaction between the neighbor particles in the chain, in the dipole-dipole approximations we get: n = 2(a n−1 u n−1 + a n u n ), Combining Equations (5)- (7) and Notations (4), after some transformations we come to the following system of nonlinear algebraic equations with respect to the dimensionless components M x and M z of the particle magnetization: This system can be solved analytically only under condition M x M z . The last is true when the sample deformation is small (i.e., when γ 1). By using the linear approximation with respect to M x , one gets from (8): M x = M z (a n−1 u n−1 +a n u n )

2NE
, Parameter a n can be found from the definition of a i in Equation (7). Let us remind that the result (9) is obtained in the approximation [55] of the identity of magnetic moments and therefore, of is the identity of magnetization M of all of the particles in the chain.
Combining the relations (3), (4), and (9), we come to the system of the linear algebraic equations with respect to the dimensionless displacement u i : This system can be solved either analytically, or numerically. The total shear stress σ in the composite that is placed in magnetic field can be presented as: Here, σ (0) and σ (m) are the nonmagnetic and magnetic parts of the stress that are produced by the aggregates. The nonmagnetic part of the stress appears because of the local inhomogeneous deformations of the elastic matrix, caused by the presence of the chains.
Nowadays, there is no a strict theoretical description of the elastic interaction of the chain with the environment. That is why we estimate σ (0) by using the approach [61], modeling the N-particles chain as a prolate ellipsoid of revolution with the minor and major axes equal to the particle diameter d and Nd, respectively. It is of fundamental importance that the volume of this ellipsoid is equal to the total volume of all of the particles in the chain. Therefore, the volume concentration of these ellipsoids is equal to the volume concentration of the particles in the ferrogels.
By using the results of the mechanics of suspensions of ellipsoidal particles [62], in the linear approximation with respect to the shear, we get: Here, ϕ is the volume concentration of the particles; α N , β N , λ N , and ζ N are some coefficients, which depend on the length of the chain. The explicit forms of these coefficients are given in the Appendix A; F N is a function of distribution over the number N of the particles in the chains, normalized so that ∑ N F N = 1. This function depends on many factors and features of the composite synthesis (size and concentration of the particles, viscosity and kinetics of the host polymer curing, the strength of the field of polymerization, etc.). The determination of the function F N presents a separate problem. Theoretical study of evolution over time of this function in the magnetic suspensions with a permanent viscosity of the currier liquid has been done in [63]. This model is based on the analysis of a system of the Smoluchowski equations, which describes the kinetic of the aggregation of the chains with a various number of particles.
In the case of cured magnetic gel, this evolution must be studied by taking into account change and time of the rheological properties of the host polymer. Depending on the molecular structure of the polymer, the concentration, and the chemical properties of the curing agent, these properties of evolution can obey to different laws, which can hardly be presented in a general form. That is why here, we suppose that the function F N is known from either independent experiments or theoretical analysis.
By using the results [64,65] for the macroscopic stress in a system of chain-like polymer macromolecules, we get the following estimate for the magnetic stress σ (m) : Here, v p is the volume of the particle. By definition, the shear modulus of the composite G: Let us introduce the dimensionless stresses: z N u n − n(a n−1 u n−1 +a n u n ) Sensors 2018, 18, 2054 7 of 13 and the dimensionless shear modulus of the composite G: The shear modulus G (0) of the composite without a magnetic field can be found from Equation (12): The magnetically induced part ∆ G of the shear modulus can be calculated as (14): For the chain with each number N of the particles, the dimensionless displacements x i are determined from the system of Equation (10), which can be presented in the form: Substituting the solution of Equations (19) into (18), we determine the dimensionless modulus ∆ G.

Results
In this part, we compare the results of our calculations with the experiments of [44]. It should be noted that the distribution function F N over number N of particles in the chains has not been determined in [44]. That is why we used the simplest approximation that all chains consist of an identical number N of particles. This number has been determined from the condition of the best agreement between the calculated and the measured [44] values of the "no field" modulus G (0) . Some of the results of this comparison and the estimated magnitudes of N are given in Table 1. With respect to the real systems, the estimated N can be considered as a characteristic number of the particles per chain.
A comparison of our calculations of the magnetically induced part ∆G of the shear modulus G with the experiments of [44] are shown in Figure 3. With respect to the real systems, the estimated N can be considered as a characteristic number of the particles per chain.
A comparison of our calculations of the magnetically induced part ∆ of the shear modulus G with the experiments of [44] are shown in Figure 3. For the systems with the relatively low volume concentration of the particles (Figure 3a), our results are in good agreement with experiments. Note that any unjustly fit parameters have not been used in our calculation of the modulus ∆ . For the higher concentrations (Figure 3b), the agreement is worse and rather, is only in the frame of the order of magnitude. The physical reason of the worsening of the agreement between the theory and the experiment lies in the fact that besides the linear chains, more topologically complicated branched, net-like, and bulk structures appear in magnetic suspensions with a high concentration of particles [48][49][50][51][52][53][54]. The spatial disposition of the particles is fixed with the host polymer gelation and determines the experimental results for the cured composite. The analysis of morphology of these structures and their effect on the macroscopic properties of the magnetic gels requires a special study. Note that, as a rule, the volume concentration of the particles in ferrogels prepared for biological applications is in the frame of several per cent, or even less than one per cent [9,39,41,66,67]. Figure 3a demonstrates that the present model leads to appropriate results for the low concentrated systems with the internal chains.
Our analysis shows that the elastic modulus of the composite significantly depends on the characteristic number N of the particles in the chains. The calculated dependencies of G on N are shown in Figure 4 for the gels, with two different magnitudes of the elastic modulus of the polymer matrix.
These results demonstrate that by varying the number N of particles in the chains, one can vary in a wide range of magnitudes, the mechanic modulus of the ferrogel. The relative increase of the modulus under the field action is high in soft gels and is less pronounced in the rigid ones. Biological ferrogels, which are used in various applications, are usually soft, with the modulus less than 10 kPa. For the systems with the relatively low volume concentration of the particles (Figure 3a), our results are in good agreement with experiments. Note that any unjustly fit parameters have not been used in our calculation of the modulus ∆G. For the higher concentrations (Figure 3b), the agreement is worse and rather, is only in the frame of the order of magnitude. The physical reason of the worsening of the agreement between the theory and the experiment lies in the fact that besides the linear chains, more topologically complicated branched, net-like, and bulk structures appear in magnetic suspensions with a high concentration of particles [48][49][50][51][52][53][54]. The spatial disposition of the particles is fixed with the host polymer gelation and determines the experimental results for the cured composite. The analysis of morphology of these structures and their effect on the macroscopic properties of the magnetic gels requires a special study. Note that, as a rule, the volume concentration of the particles in ferrogels prepared for biological applications is in the frame of several per cent, or even less than one per cent [9,39,41,66,67]. Figure 3a demonstrates that the present model leads to appropriate results for the low concentrated systems with the internal chains.
Our analysis shows that the elastic modulus of the composite significantly depends on the characteristic number N of the particles in the chains. The calculated dependencies of G on N are shown in Figure 4 for the gels, with two different magnitudes of the elastic modulus G 0 of the polymer matrix.
These results demonstrate that by varying the number N of particles in the chains, one can vary in a wide range of magnitudes, the mechanic modulus of the ferrogel. The relative increase of the modulus under the field action is high in soft gels and is less pronounced in the rigid ones. Biological ferrogels, which are used in various applications, are usually soft, with the modulus less than 10 kPa. Therefore, their mechanical properties and behavior can be effectively controlled with the help of an applied magnetic field.
The characteristic length of the chains is determined by the condition of the gel polymerization, which is in part determined by the ratio between the kinetics of the particle's aggregation and the rate of the host polymer's curing. Thus, by changing the condition of the system gelation, one can tune in a wide range of the magnitudes of the macroscopic properties of the composite material.
Therefore, their mechanical properties and behavior can be effectively controlled with the help of an applied magnetic field.
The characteristic length of the chains is determined by the condition of the gel polymerization, which is in part determined by the ratio between the kinetics of the particle's aggregation and the rate of the host polymer's curing. Thus, by changing the condition of the system gelation, one can tune in a wide range of the magnitudes of the macroscopic properties of the composite material.

Discussion
We present the results of the theoretical modeling of magnetorheological effects in magnetic gels with chain-like aggregates. Unlike the previous theoretical models suggested in [44], this model does not contain any unjustly fit parameters. In the frames of applicability of the hypothesis that only linear chains appear at the stage of the composite synthesis (i.e., that only materials with low or moderate concentrations of particles are considered), the model is in good agreement with the experiments of [44] (see Figure 3a). This agreement indicates that the proposed model leads to adequate results for the composite with the concentration of particles, at least in the frames of ten per cent. For the higher concentrations, the appearance of topologically complicated structures is quite probable, and that is why these concentrations are out of the scope of this model.
Note that the volume concentration of the particles in the magnetic gels that are synthesized for biomedical applications, as a rule, is low and in the frames of several per cent (see, for example [9,39,41,66,67] and the references therein). The systems with higher concentrations are usually synthesized for various mechanical systems (dampers, actuators, etc.).
The obtained results can be considered as a theoretical background for the development of technologies of magnetically controllable biosensors, scaffolds with tunable properties, and for engineering and the regeneration of biological tissues.
In real magnetic gels, the distribution over chain size can be quite broad. At the same time, the assumption of the identity of the chains in the composite, which was used in part 2 of this work, is not necessary for the present model. Indeed, the relation 17 and relation 18 can be used to estimate the elastic modulus if the distribution function FN is known.
Unfortunately, the law of the size of distribution has not been studied in [44]. That is why, and only because of that, that we have used the model of identical chains to compare our results with the experiments [44]. The determined number N can be considered as an estimate of the characteristic size of the chains in the real composites.

Discussion
We present the results of the theoretical modeling of magnetorheological effects in magnetic gels with chain-like aggregates. Unlike the previous theoretical models suggested in [44], this model does not contain any unjustly fit parameters. In the frames of applicability of the hypothesis that only linear chains appear at the stage of the composite synthesis (i.e., that only materials with low or moderate concentrations of particles are considered), the model is in good agreement with the experiments of [44] (see Figure 3a). This agreement indicates that the proposed model leads to adequate results for the composite with the concentration of particles, at least in the frames of ten per cent. For the higher concentrations, the appearance of topologically complicated structures is quite probable, and that is why these concentrations are out of the scope of this model.
Note that the volume concentration of the particles in the magnetic gels that are synthesized for biomedical applications, as a rule, is low and in the frames of several per cent (see, for example [9,39,41,66,67] and the references therein). The systems with higher concentrations are usually synthesized for various mechanical systems (dampers, actuators, etc.).
The obtained results can be considered as a theoretical background for the development of technologies of magnetically controllable biosensors, scaffolds with tunable properties, and for engineering and the regeneration of biological tissues.
In real magnetic gels, the distribution over chain size can be quite broad. At the same time, the assumption of the identity of the chains in the composite, which was used in part 2 of this work, is not necessary for the present model. Indeed, the relation 17 and relation 18 can be used to estimate the elastic modulus if the distribution function F N is known.
Unfortunately, the law of the size of distribution has not been studied in [44]. That is why, and only because of that, that we have used the model of identical chains to compare our results with the experiments [44]. The determined number N can be considered as an estimate of the characteristic size of the chains in the real composites.
In principle, our approach allows the studying of the large shear deformations, including the rupture of the chains, and Equations (1)-(3) can be solved numerically in nonlinear approximation with respect to the particles' displacement u i . This can be a natural continuation of the present work.

Conflicts of Interest:
The authors declare no conflict of interest.