Linear Dipole Model

The linear dipole model is mechanically the same as the Linear Model, comprised of both linear and dashpot components. In addition to unbonded mechanical interactions, forces and moments based on discrete magnetic dipoles may be included to explicitly model permanent and induced magnetization.


The linear dipole model provides mechanical linear and dashpot components that act in parallel with one another. The linear component provides linear elastic (no-tension) frictional behavior, while the dashpot component provides viscous behavior. Both components act over a vanishingly small area, and thus transmit only a force. The mechanical model is thoroughly elucidated in the Linear Model document, and the user is referred there for further details of the mechanical characteristics of the model.

In order to model magnetic interactions, each piece is abstracted as a magnetic dipole based on [Thomaszewski2008]. Forces and moments due to the dipole interactions are computed and also applied to the pieces. The moments applied to each piece due to the dipole interactions are, in general, not the same. One must also use the proximity concept (see the contact cmat proximity command) so that contacts may be deemed active despite the lack of physical overlap. Inheritance of the dipole_mo surface property is also required along with orientation tracking (see the model orientation-tracking command).

Behavior Summary

The linear dipole model provides the behavior of an infinitesimal interface that does not resist relative rotation. The mechanical contact force is resolved into linear and dashpot components. The linear component provides linear elastic (no-tension), frictional behavior. The dashpot component provides viscous behavior. Review the Linear Model document for further details of the mechanical response. Additional forces and moments due to interacting magnetic dipoles are also added. The total contact force is \(\mathbf{F_{c}} = \mathbf{F^{l}} + \mathbf{F^{d}} + \mathbf{F^{di}}\), where \(\mathbf{F^{l}}\) is the linear force, \(\mathbf{F^{d}}\) is the dashpot force, and \(\mathbf{F^{di}}\) is the magnetic dipole force. The moments applied to each piece from the magnetic dipoles are generally different from each other and different from zero (i.e., \(\mathbf{M_{1}} \neq \mathbf{M_{2}} \neq 0\), where \(\mathbf{M_{1}}\) is the moment applied to piece 1 and \(\mathbf{M_{2}}\) is the moment applied to piece 2). These forces and moments are based on magnetostatic computations assuming dipoles as outlined below.

Activity-Deletion Criteria

A contact with the linear dipole model is active if the surface gap is less than or equal to the dipole distance \(d_d\) (0 by default). The force-displacement law is skipped for inactive contacts. Mechanical forces may not arise even if dipole forces/moments exist.

Force-Displacement Law

In addition to the force displacements computed in the linear model, this model computes forces and moments due to magnetic dipoles interacting at a distance. These forces are added to the linear and dashpot forces.

Given the magnetic field \(\mathbf{B}\), the magnetizing field \(\mathbf{H}\), the current density \(\mathbf{J}\), and the magnetization \(\mathbf{M}\), Maxwell’s equations can be written as

(1)\[\nabla \cdot \mathbf{B} = 0\]
(2)\[\nabla \times \mathbf{H} = \mathbf{J}\]
(3)\[\mathbf{H} = {1 \over{\mu_0}}\mathbf{B} - \mathbf{M}\]

where \(\mu_0\) is the permeability of free space and \(\nabla\) is the gradient.

Suppose there is a magnetic moment \(\mathbf{m}\) positioned at the center of a rigid piece given a specified current density. Assuming this can be represented as a vanishingly small magnetic dipole, the force exerted on the rigid body in the presence of a magnetic field is given by

(4)\[\mathbf{F^{di}} = \nabla (\mathbf{m} \cdot \mathbf{B})\]

while the moment is

(5)\[\mathbf{M^{di}} = \mathbf{m} \times \mathbf{B}\]

The magnetic field as a function of distance from a vanishingly small magnetic dipole is

(6)\[\mathbf{B}(r) = {\mu_0 \over {4 \pi}} \left[{{{3\mathbf{n}(\mathbf{n} \cdot \mathbf{m}) - \mathbf{m}}} \over {{(|\mathbf{r} - \mathbf{r_0}|)}^3}}\right]\]

where \(\mathbf{n} = (\mathbf{r} - \mathbf{r_0}) / |\mathbf{r} - \mathbf{r_0}|\) and \(\mathbf{r_0}\) is the center of the rigid piece.

As noted in [Thomaszewski2008], for an arbitrarily shaped object with uniform magnetization, one can decompose it into \(N\) magnetic dipoles so that the total magnetic field is

(7)\[\mathbf{B}(r) = {\mu_0 \over {4 \pi}} \sum_{i=1}^N \left[{{{3\mathbf{n_i}(\mathbf{n_i} \cdot \mathbf{m_i}) - \mathbf{m_i}}} \over {{(|\mathbf{r} - \mathbf{r_i}|)}^3}}\right]\]

For an infinite number of dipoles, the solution converges to the total magnetic field if their spatial distribution discretizes the piece shape. As the shapes of the pieces/bodies are not incorporated explicitly in the model, one can discretize a body into multiple clump pebbles to simulate non-spherical magnet geometries using this contact model.

The forces/moments represent the interactions of two magnetic dipoles with magnetic moments \(\mathbf{m_1}\) and \(\mathbf{m_2}\) at positions \(\mathbf{r_1}\) and \(\mathbf{r_2}\); these are the piece positions. Given the definition of the magnetic field due to a dipole, the resulting force is

(8)\[\begin{split}\begin{array}{l} \mathbf{F^{di}} = {3 \mu_0 \over {4 \pi {|\mathbf{r_1} - \mathbf{r_2}|}^4}} \left[ {\mathbf{n_{12}}(\mathbf{m_1} \cdot \mathbf{m_2}) + \mathbf{m_1}(\mathbf{m_2} \cdot \mathbf{n_{12}}) + \\ \mathbf{m_2}(\mathbf{m_1} \cdot \mathbf{n_{12}}) - 5\mathbf{n_{12}}(\mathbf{m_1} \cdot \mathbf{n_{12}})(\mathbf{m_2} \cdot \mathbf{n_{12}})}\right] \end{array}\end{split}\]

where \(\mathbf{n_{12}} = (\mathbf{r_1} - \mathbf{r_2}) / |\mathbf{r_1} - \mathbf{r_2}|\). This is the force applied to piece 2. The force applied to piece 1 is equal and opposite, replacing \(\mathbf{n_{12}}\) with \(\mathbf{n_{21}}\). As the distance between pieces shrinks to a vanishingly small value the force becomes infinite. The \(d_{cap}\) is used to place a cap on the minimum distance used for force and moment computations.

The moments are

(9)\[\mathbf{M_{2}} = {\mu_0 \over {4 \pi {|\mathbf{r_1} - \mathbf{r_2}|}^3}} \left[ {3 (\mathbf{m_2} \times \mathbf{n_{12}})(\mathbf{m_1} \cdot \mathbf{n_{12}}) - \mathbf{m_2} \times \mathbf{m_1}}\right]\]


(10)\[\mathbf{M_{1}} = {\mu_0 \over {4 \pi {|\mathbf{r_1} - \mathbf{r_2}|}^3}} \left[ {3 (\mathbf{m_1} \times \mathbf{n_{21}})(\mathbf{m_2} \cdot \mathbf{n_{21}}) - \mathbf{m_1} \times \mathbf{m_2}}\right]\]

These moments are not equal and opposite. Also notice that the moments decay less rapidly than the forces, meaning that moments dominate at larger distances.

The above forces and moments presume the magnetic moments are in the local coordinate system. They should be specified as a piece surface property and inherited. Since these are vectors, though, they must be rotated to account for the rigid body rotation of the pieces. Thus, at the point in the simulation when they are specified as piece surface properties, it is presumed no rotations have occurred and orientation tracking is active. Subsequently they are rotated as the pieces spin.

Energy Partitions

The energies computed are the same as for the linear energies.


The properties defined by the linear dipole contact model are listed in the table below as a concise reference. See the Contact Properties section for a description of the information in the table columns. The mapping from the surface inheritable properties to the contact model properties is also discussed below.

Table 1: Linear Dipole Model Properties
Keyword Symbol Description Type Range Default Modifiable Inheritable
linear Model name
Linear Group:
kn \(k_n\) Normal stiffness [force/length] FLT \([0.0,+\infty)\) 0.0 YES YES
ks \(k_s\) Shear stiffness [force/length] FLT \([0.0,+\infty)\) 0.0 YES YES
fric \(\mu\) Friction coefficient [-] FLT \([0.0,+\infty)\) 0.0 YES YES
rgap \(g_r\) Reference gap [length] FLT \(\mathbb{R}\) 0.0 YES NO
lin_mode \(M_l\) Normal-force update mode [-] INT \(\{0,1\}\) 0 YES NO
    \(\;\;\;\;\;\;\begin{cases} \mbox{0: update is absolute} \\ \mbox{1: update is incremental} \end{cases}\)          
emod \(E^*\) Effective modulus [force/area] FLT \([0.0,+\infty)\) 0.0 NO N/A
kratio \(\kappa^*\) Normal-to-shear stiffness ratio [-] FLT \([0.0,+\infty)^*\) 0.0\(^*\) NO N/A
    \(\kappa^* \equiv \frac{k_n}{k_s}\)          
lin_slip \(s\) Slip state [-] BOOL {false,true} false NO N/A
    \(\;\;\;\;\;\;\begin{cases} \mbox{true: slipping} \\ \mbox{false: not slipping} \end{cases}\)          
lin_force \(\mathbf{F^l}\) Linear force (contact plane coord. system) VEC \(\mathbb{R}^3\) \(\mathbf{0}\) YES NO
    \(\left( -F_n^l,F_{ss}^l,F_{st}^l \right) \quad \left(\mbox{2D model: } F_{ss}^l \equiv 0 \right)\)          
user_area \(A\) Constant area [length*length] FLT \((0.0,+\infty)\) 0.0 YES NO
Dashpot Group:
dp_nratio \(\beta_n\) Normal critical damping ratio [-] FLT \([0.0,1.0]\) 0.0 YES NO
dp_sratio \(\beta_s\) Shear critical damping ratio [-] FLT \([0.0,1.0]\) 0.0 YES NO
dp_mode \(M_d\) Dashpot mode [-] INT {0,1,2,3} 0 YES NO
    \(\;\;\;\;\;\;\begin{cases} \mbox{0: full normal & full shear} \\ \mbox{1: no-tension normal & full shear} \\ \mbox{2: full normal & slip-cut shear} \\ \mbox{3: no-tension normal & slip-cut shear} \end{cases}\)          
dp_force \(\mathbf{F^d}\) Dashpot force (contact plane coord. system) VEC \(\mathbb{R}^3\) \(\mathbf{0}\) NO NO
    \(\left( -F_n^d,F_{ss}^d,F_{st}^d \right) \quad \left(\mbox{2D model: } F_{ss}^d \equiv 0 \right)\)          
Dipole Group:
dipole_cap \(d_{cap}\) Capping distance for force/moment [length] FLT \((0.0,+\infty)\) 0.0 YES NO
dipole_d \(d_d\) Activity distance [length] FLT \((0.0,+\infty)\) 0.0 YES NO
dipole_f \(\mathbf{F^{di}}\) Dipole force (contact plane coord. system) VEC \(\mathbb{R}^3\) \(\mathbf{0}\) YES NO
dipole_m1 \(\mathbf{M_1}\) Moment applied to piece 1 (contact plane coord. system) VEC \(\mathbb{R}^3\) \(\mathbf{0}\) YES NO
dipole_m2 \(\mathbf{M_2}\) Moment applied to piece 2 (contact plane coord. system) VEC \(\mathbb{R}^3\) \(\mathbf{0}\) YES NO
dipole_mo1 \(\mathbf{m_1}\) Dipole moment of piece 1 VEC \(\mathbb{R}^3\) \(\mathbf{0}\) YES YES
dipole_mo2 \(\mathbf{m_2}\) Dipole moment of piece 2 VEC \(\mathbb{R}^3\) \(\mathbf{0}\) YES YES
\(^*\) By convention, \(\kappa^* = 0.0\) if either the normal or the shear stiffness is 0.


Out-of-balance forces acting on bodies are accumulated during force-displacement calculations. Modifying the forces stored in the contact models will not alter them instantaneously. Therefore, any change to \(\mathbf{F^l}\) may only be effective during the next force-displacement calculation.

Surface Property Inheritance

The linear stiffnesses, \(k_n\) and \(k_s\), and the friction coefficient, \(\mu\), may be inherited from the contacting pieces, along with the original dipole moments \(\mathbf{m_1}\) and \(\mathbf{m_2}\). Remember that for a property to be inherited, the inheritance flag for this property must be set to true (default), and both contacting pieces must hold a property with the exact same name.

The linear stiffnesses are inherited, assuming that both pieces’ stiffnesses act in series (Figure 5):

(11)\[\begin{split}\frac{1}{k_n} = \frac{1}{k^{(1)}_n} + \frac{1}{k^{(2)}_n} \\ \frac{1}{k_s} = \frac{1}{k^{(1)}_s} + \frac{1}{k^{(2)}_s}\end{split}\]

where \((1)\) and \((2)\) denote the properties of piece 1 and 2, respectively.

The friction coefficient is inherited using the minimum of the values set for the pieces:

(12)\[\mu = \min(\mu^{(1)},\mu^{(2)}).\]

The dipole moments are taken from each piece and stored.


Table 2: Linear Dipole Model Methods
Method Arguments Symbol Type Range Default Description
area Set user_area to the area
deformability Set deformability
  emod \(E^*\) FLT \([0.0,+\infty)\) N/A Effective modulus
  kratio \(\kappa^*\) FLT \([0.0,+\infty)^*\) N/A Normal-to-shear stiffness ratio
\(^*\) By convention, setting \(\kappa^* = 0.0\) sets \(k_s = 0.0\) but does not alter \(k_n\).

See here for a description of the methods.

Callback Events

Table 3: Linear Dipole Model Callback Events
Event Array Slot Value Type Range Description
contact_activated Contact has become active
  1 C_PNT N/A Contact pointer
slip_change Slip state has changed
  1 C_PNT N/A Contact pointer
  2 INT {0,1} Slip change mode
        \(\;\;\;\;\;\;\begin{cases} \mbox{0: slip has initiated} \\ \mbox{1: slip has ended} \end{cases}\)

Usage and Verification Examples

The verification “Simulating Square Permanent Magnets” compares the analytic response of the attractive forces/moments of a square magnet with increasing levels of discretization of a cube with dipoles.

Model Summary

An alphabetical list of the linear dipole model methods is given here. An alphabetical list of the linear dipole model properties is given here.


[Thomaszewski2008](1, 2) Thomaszewski, B., A. Gunmann, S. Pabst and W. Straßer. “Magnets in Motion,” in ACM Transactions on Graphics, 27, pp. 1-9, 2008.