---
title: Discrete Configurational Forces
url: https://www.emergentmind.com/topics/discrete-configurational-forces
type: topic
---

# Discrete Configurational Forces

Searching arXiv for the cited papers to ground the article in current research.
Discrete configurational forces are energy-based material-force constructs defined on discrete lattices, pixel grids, finite-element meshes, defect lines, or discretized configuration spaces, and they measure how a system’s total energy changes under changes of material configuration rather than ordinary spatial displacement. In recent work, the concept appears in several technically distinct but related forms: equivalent nodal forces assembled from the Eshelby stress in nonlinear fracture and topology optimization, J-type equivalent domain integrals evaluated from measured elastic fields near blocked slip bands, generalized inner-variation forces and stresses in Kohn–Sham density functional theory, Peierls–Nabarro-modulated driving forces for domain walls in multistable media, Noether-current force densities in mesoscopic Cosserat continua with distributed defects, and a discrete counting of admissible force configurations in jammed granular matter [2507.12247] [2507.15570] [2603.24129] [1712.05535] [2512.20453] [2604.12984] [1912.03504].

## 1. Variational basis and material-force definitions

In configurational mechanics, the central object is an energy–momentum tensor of Eshelby type. For finite-strain hyperelasticity, the mixed-mode fracture study defines
\[
\Sigma(F)=\Psi\,I-F^{T}\frac{\partial\Psi}{\partial F},
\]
with stored energy density \(\Psi(F)\) and first Piola–Kirchhoff stress \(P=\partial\Psi/\partial F\). In the topology-optimization formulation, the analogous quantity is written
\[
\Sigma=U_0\,I-F^T P,
\qquad
U_0=W_0(F;X)+V_0(\varphi;X),
\]
and the local configurational traction is \(T_0=\Sigma\cdot N\) on \(\partial B_0\). In small-strain blocked-slip analysis, the corresponding material tensor is
\[
P_{kj}=W\,\delta_{kj}-\sigma_{ij}\,u_{i,k},
\qquad
W=\tfrac12\,\sigma_{mn}\,\varepsilon_{mn},
\]
which yields configurational forces by contour or domain integration [2507.12247] [2507.15570] [2603.24129].

The defining principle is that configurational forces act in material space. In fracture mechanics, the energy release rate \(G\) is obtained from the flux of \(\Sigma\) through a surface around the crack tip, with \(G=-\,d\Pi/da\). In blocked-slip problems, Rice’s \(J\)-integral appears as the \(k=1\) component of the configurational force, with the convention \(J=-F_1\). In Kohn–Sham density functional theory, the same variational idea is implemented through inner variations \(x\mapsto x+\varepsilon\Theta(x)\), so that the derivative of the constrained ground-state energy is written in terms of a generalized Eshelby tensor, nuclear self terms, and, in pseudopotential calculations, additional non-local contributions [1712.05535].

A recurring distinction in these formulations is between spatial equilibrium and configurational balance. The blocked-slip formulation states the usual equilibrium equations \(\sigma_{ij,j}=0\), while the topology-optimization formulation derives the conventional force balance \(\mathrm{Div}\,P+b_0=0\). Configurational forces are additional energetic descriptors: they govern where cracks advance, where defects are driven, where a mesh should refine, or which extension directions are energetically favored [2603.24129] [2507.15570].

## 2. Discrete realizations

Recent work uses several discrete realizations of configurational-force concepts.

| Setting | Discrete quantity | Construction |
|---|---|---|
| Nonlinear FE fracture | \(\mathbf{F}^A_{\rm CNF}\) | \(\sum_e \int_{\mathcal B_e}\Sigma:\nabla_0N^A\,dV\) |
| SIMP topology optimization | \(F^I_{\rm CNF}\) | \(\sum_e \int_{\Omega_e}\Sigma_\rho(\hat\rho_e):\nabla_XN^i\,dV\) |
| HR-EBSD blocked slip | \(F_k\), \(J=-F_1\) | J-type equivalent domain integral with \(q(x)\) |
| Axisymmetric domain wall | \(F_{\rm conf}(R)\) | \(-\,dU_{\rm tot}/dR\) with periodic modulation |
| Cosserat line defect | \(f_A\) | \(f_A=b^i\Sigma_{iA}+K^{ij}M_{jiA}\) |
| Granular force ensemble | \(\Omega(P)/h_f^{N_f}\) | discretized force-space multiplicity |

In finite-element fracture, the discrete configurational-force vector at node \(A\) is assembled elementwise as
\[
\mathbf{F}^{A}_\text{CNF}
=\sum_{e=1}^{n_\text{el}}
\int_{\mathcal{B}_e}
\Sigma : \nabla_0 N^A
\,dV.
\]
The algorithm is explicitly post-processing based: after the nonlinear FE solve, one computes \(F=\nabla_0u+I\), \(\Psi(F)\), \(P\), and \(\Sigma\) at Gauss points; integrates element contributions \(f_e^A\); assembles the global nodal vector; identifies nodes on the crack-surface edge plus a narrow “spurious” neighbor band; transforms forces to local cylindrical or Frenet frames; and sums the forces over crack-front segments per unit arc length \(s\) [2507.12247].

In topology optimization, the discrete force is modified by the design field. Each element carries a pseudo-density \(\rho_e\), from which one defines
\[
\hat\rho_e=\frac{\tilde\rho_e}{\varepsilon(1-\tilde\rho_e)+\tilde\rho_e},
\qquad
\Sigma_\rho(\hat\rho_e)=f_\varepsilon(\hat\rho_e)\Sigma,
\]
with
\[
f_\varepsilon(\hat\rho)=\frac{\hat\rho}{\varepsilon(1-\hat\rho)+\hat\rho}.
\]
The nodal configurational force then becomes
\[
F_{\rm CNF}^I
=
\sum_{e=1}^{n_{\rm el}}
\int_{\Omega_e}
\Sigma_\rho(\hat\rho_e):\nabla_XN^i(X)\,dV.
\]
The scalar magnitude \(\|F_{\rm CNF}^I\|\) is used directly as an adaptive-mesh indicator [2507.15570].

In HR-EBSD slip-band analysis, a contour integral is converted to the J-type equivalent domain integral
\[
J
=
\int_{\Omega}
\Bigl(\sigma_{ij}\,\frac{\partial u_i}{\partial x_1}
-W\,\delta_{j1}\Bigr)\,
\frac{\partial q}{\partial x_j}\,dA,
\]
where \(q(x)\) is a plateau function equal to \(1\) inside an inner domain, \(0\) outside an outer domain, and varying smoothly between them. On the pixel grid, the integral is evaluated by a Riemann sum using HR-EBSD displacement-gradient data, anisotropic elasticity for \(\alpha\)-Ti, and finite differences of \(q\) [2603.24129].

These constructions share a discrete evaluation strategy, but the source of discreteness differs by field. In FE applications, discreteness enters through the mesh and shape functions. In blocked-slip analysis, it enters through the measured pixel grid. In multistable lattices, it enters through the lattice spacing \(d\) and the associated Peierls–Nabarro modulation. In granular materials, it enters through the subdivision of continuous force space into \(h_f\)-cells [2512.20453] [1912.03504].

## 3. Crack-front driving forces in soft fracture

For mixed-mode I + III fracture in soft, highly deformable solids, the discrete configurational-force method is implemented as a post-processing algorithm on finite-strain, neo-Hookean FE simulations with no body forces. The material model uses
\[
\Psi(F)=\Psi_{\rm iso}(\overline F)+\Psi_{\rm vol}(\det F),
\]
with \(\Psi_{\rm iso}=(G/2)[\overline F^T:\overline F-3]\) and \(\Psi_{\rm vol}=(\kappa/2)[\det F-1]^2\), together with \(\nu=0.451\), \(G=9\,\mathrm{kPa}\), and \(\kappa=2G(1+\nu)/[3(1-2\nu)]\). Loading is imposed by Dirichlet conditions: first an increasing vertical separation \(\Delta L\) for Mode I, then an increasing bottom-to-top twist angle \(\alpha\) for Mode I + III [2507.12247].

The front is partitioned into planar-segment nodes and facet-edge nodes. For planar segments, the force is expressed in a cylindrical basis \(\{\hat r,\hat\theta,\hat z\}\); for facets, in the local Frenet triad \(\{\hat T,\hat N,\hat B\}\). The segmentwise measure is
\[
\mathbf{F}_{\rm CNF}/s
=
\frac{1}{s}
\sum_{A\in \mathcal N_{\rm segment}}
\mathbf{F}^A_{\rm CNF},
\]
and its magnitude and orientation are then analyzed [2507.12247].

The resulting interpretation is explicitly predictive. The scalar magnitude \(|F_{\rm CNF}/s|\) is proportional to the local energy-release rate \(G\), and crack advance is expected when \(|F_{\rm CNF}/s|\) exceeds a material-specific threshold \(F_c\), equivalently \(G\ge G_c\). The orientation of \(\mathbf{F}_{\rm CNF}/s\) indicates the preferred growth direction. On tilted facets, the dominant component lies along the binormal \(\hat B\), showing that cracks advance while maintaining their tilt angle \(\phi\). A nonzero component along the facet normal \(\hat N\) indicates a slight tendency for facet rotation, and under mixed-mode I + III the angle \(\beta\) between \(\mathbf{F}_{\rm CNF}/s\) and the facet-plane normal can change sign for small \(\phi\), revealing shear-induced reorientation away from the original tilt [2507.12247].

The same discrete force field clarifies elastic interactions along complex crack fronts. On planar segments, \(|F_{\rm CNF}/s|\) often exceeds the values on facets, especially under mixed-mode loading, which predicts the later growth of “type B” bridging cracks. Varying the facet spacing \(\Lambda\) shows that close spacing amplifies \(|F_{\rm CNF}/s|\) on planar regions through constructive elastic interaction, whereas very small \(\Lambda\) leads to shielding. In facet-coalescence studies, mesh relaxation by DynaMesh-R is used to join facets smoothly, after which configurational-force redistribution shows a drop near the new junction and amplification farther away, producing an asymmetric driving-force distribution that controls non-uniform growth and the final echelon morphology [2507.12247].

## 4. Interface motion, pinning, and blocked-slip energetics

In elastically coupled multistable metamaterials, an axisymmetric domain wall of normalized radius \(R\) is described by a reduced-order model in which the total potential energy is the sum of plate bending and foundation contributions. In normalized form,
\[
\bar U_{\rm plate}(R)
=
\frac{1}{15\,W^2}
\Biggl(
256\,\frac{R}{W}
+
20\,\frac{W}{R}
\Biggr),
\]
\[
\bar U_{\rm fond}(R)
=
\frac{W^2}{384}
\Bigl[
-\Delta
+
4\,\frac{R}{W}
-
32\,\Delta
\Bigl(\frac{R}{W}\Bigr)^2
\Bigr],
\]
and
\[
\bar U_{\rm cont}(R)=\bar U_{\rm plate}(R)+\bar U_{\rm fond}(R).
\]
Here \(\Delta=1-2\alpha\) is the asymmetry parameter and \(\varepsilon=\gamma=d/W\) is the discreteness parameter, with \(W\approx 8L\) and \(\Delta>0\) favoring expansion [2512.20453].

The configurational driving force is
\[
F_{\rm conf}(R)=-\,\frac{dU_{\rm tot}}{dR}.
\]
In the continuum limit \(\varepsilon\to 0\), the force follows from differentiating \(\bar U_{\rm cont}\). For finite discreteness, the model superposes a periodic Peierls–Nabarro correction of period \(\lambda\approx \sqrt3\,d\) and amplitude \(A(\Delta,\varepsilon)\),
\[
\bar U(R)\approx \bar U_{\rm cont}(R)+A(\Delta,\varepsilon)\cos\!\Bigl(\tfrac{2\pi}{\lambda}R\Bigr),
\]
which yields a modulated configurational force with local minima whenever \(\bar F_{\rm conf}(R^*)=0\) and \(d\bar F/dR(R^*)>0\). These minima create pinning wells. The paper identifies three regimes—expansion, shrinking, and metastable pinning—and states that local minima exist only for \(\varepsilon>\varepsilon_{\rm cr}(\Delta)\), within a finite band \(R_{\min}\le R^*\le R_{\max}\) [2512.20453].

The continuum nucleation radius \(R_{\rm cr}(\Delta)\) is determined by the global maximum of \(\bar U_{\rm cont}(R)\). If the initial wall radius satisfies \(R_0<R_{\rm cr}\), then \(F_{\rm conf}<0\) and the wall shrinks; if \(R_0>R_{\rm cr}\), then \(F_{\rm conf}>0\) and the wall expands. The same axisymmetric reduced-order model extends to non-axisymmetric polygonal walls. For convex \(120^\circ\) polygons, a regular hexagon of side length \(\ell=\sqrt3\,d\,n\) is metastable whenever \(n\) lies in the discrete stability set \(\{n_{\min}(\Delta,\varepsilon),\dots,n_{\max}(\Delta,\varepsilon)\}\), and an irregular convex polygon is fully pinned if every edge length \(\ell_i=\sqrt3\,d\,n_i\) satisfies \(n_{\min}\le n_i\le n_{\max}\). For concave \(240^\circ\) corners, pinning requires that the complementary regular hexagons be stable under the complementary asymmetry \(\Delta'=-\Delta\) [2512.20453].

A different interface-localization problem appears in blocked slip bands at grain boundaries in \(\alpha\)-Ti. There, the configurational force vector
\[
F_k
=
\int_{\Omega}
\bigl[
\sigma_{ij}\,\partial_k u_i
-
W\,\delta_{jk}
\bigr]\,
\partial_j q\,dA
\qquad (k=1,2)
\]
is evaluated from HR-EBSD measurements. Directionality is resolved by defining a Virtual Extension Direction and rotating the displacement-gradient field into the frame of each candidate slip-system trace in the neighboring grain. The polar function \(J(\theta)=-F_1(\theta)\) exhibits peaks at the most energetically favorable extension directions. In the reported example, the maximum \(J\) occurs on a pyramidal \((c+a)\) trace that is geometrically admissible but not the one with the highest Schmid factor; the highest-Schmid prismatic variant has only moderate \(J\), and a basal variant with \(SF\approx 0\) still produces non-zero \(J\). The paper therefore states that Schmid factor, \(m'\), and residual Burgers vector do not uniquely predict the local energetic driving force \(J\) [2603.24129].

## 5. Adaptive computation and electronic-structure forces

In SIMP-based topology optimization, discrete configurational forces are used as a mesh-adaptivity criterion because, with the relaxed Eshelby-stress interpolation, they localize both in grey transition regions and in highly stressed regions. The localization mechanism is explicit in the paper: fully void elements have \(\hat\rho_e\approx 0\) and therefore \(\Sigma_\rho\approx 0\); fully solid elements have \(\hat\rho_e\approx 1\), but \(\Sigma\) is often only moderate if energy density and stress vary smoothly; by contrast, grey transition elements generate large values because the density gradients align with \(\nabla_XN^i\), while high-stress solid regions also generate large values through the local Eshelby stress [2507.15570].

The marking rule uses the normalized maximum
\[
F^{\rm max}_{\rm CNF}=\max_I \|F^I_{\rm CNF}\|,
\]
and thresholds
\[
\|F^I_{\rm CNF}\|\ge c_r\,F^{\rm max}_{\rm CNF}
\quad\text{for refinement},\qquad
\|F^I_{\rm CNF}\|\le c_c\,F^{\rm max}_{\rm CNF}
\quad\text{for coarsening},
\]
with the common choice
\[
c_r=0.25,\qquad c_c=0.01.
\]
Refinement is aggressive, whereas coarsening is conservative: an element is coarsened only if all its vertices satisfy the coarsening threshold. The algorithm stores the refinement tree, allowing multilevel coarsening up to a user-defined maximum refinement level \(L_{\max}\) [2507.15570].

The numerical results reported in the paper are problem-specific. For a cantilever compliance-minimization example with a \(40\times 20\) base mesh and max-level \(4\), configurational-force-based adaptivity achieves approximately \(21\%\) of the total CPU time of density-based adaptivity with identical compliance. For a U-beam example with compliance, volume, and \(p\)-norm stress constraint, configurational-force-based adaptivity requires approximately \(40\%\) of the CPU time of density-based adaptivity for \(\Delta k=5\), while yielding equal or better enforcement of stress constraints. The paper contrasts this with density-based refinement, which refines everywhere grey, and von-Mises-based refinement, which refines stress-critical corners but can miss geometric boundaries [2507.15570].

In Kohn–Sham density functional theory, configurational forces are derived as generalized variational forces obtained from inner variations of the Kohn–Sham energy functional with respect to the position of a material point \(x\). The derivative of the ground-state energy has the form
\[
\frac{dE_0}{d\varepsilon}
=
\int_\Omega \mathbf E(x):\nabla\Theta(x)\,dx
+
\sum_I\int_{\mathbb R^3}\mathbf E'_I(x):\nabla\Theta(x)\,dx
+
F^{\rm PSP}[\Theta],
\]
where \(\mathbf E\) is the generalized Eshelby tensor, \(\mathbf E'_I\) the nuclear-self contribution, and \(F^{\rm PSP}\) the non-local pseudopotential term. By choosing \(\Theta\) as an atom-centered compact generator \(\Theta_I\), one obtains atomic forces; by choosing \(\Theta(x)=C\,x\), one obtains the stress tensor in a periodic cell [1712.05535].

The method is variationally consistent and the paper states that Pulay corrections are inherently included, so no separate Pulay-correction step is needed. It also treats pseudopotential and all-electron calculations in a single framework. In the reported higher-order FE benchmarks, forces and stresses show \(O(h^{2k-1})\) convergence for CO, CH\(_4\), N\(_2\), SiF\(_4\), Al fcc, and Li bcc, and the paper reports errors typically below \(10^{-4}\,\mathrm{Ha/Bohr}\) in forces and below \(10^{-6}\,\mathrm{Ha/Bohr}^3\) in stresses when compared with reference calculations [1712.05535].

## 6. Defects, Noether currents, and discrete line-force limits

A geometrically richer generalization appears in the mesoscopic Cosserat theory with distributed defects. The basic fields are the coframe \(e^i\) and an independent \(\mathfrak{so}(3)\)-connection \(\omega^{ij}=-\omega^{ji}\), with torsion and curvature
\[
T^i:=De^i=de^i+\omega^i{}_j\wedge e^j,
\qquad
R^{ij}:=D\omega^{ij}=d\omega^{ij}+\omega^i{}_k\wedge\omega^{kj}.
\]
The Palatini-type action leads to Euler–Lagrange equations
\[
D\,H_i+\partial_t P_i=\Sigma_i,
\qquad
D\,O_{ij}+\partial_t Q_{ij}+e_{[i}\wedge H_{j]}=M_{ij},
\]
which combine the standard force and couple balances with defect-excitation terms \(H_i\) and \(O_{ij}\) [2604.12984].

Material invariance then yields configurational currents as Noether quantities. Under infinitesimal material translations, the theory gives
\[
D S^A+\partial_t\mathscr I^A=R^A,
\]
where \(S^A\) is the configurational stress current, \(\mathscr I^A\) the configurational momentum, and \(R^A\) the configurational force density. Under material rotations, an analogous identity yields configurational moment density. The same framework also provides dynamic Bianchi transport laws,
\[
\partial_t T^i=D\,J^i+K^i{}_j\wedge e^j,
\qquad
\partial_t R^{ij}=D\,K^{ij},
\]
which connect defect transport directly to configurational forces and moments [2604.12984].

The discrete line-defect limit is especially explicit. If a defect line \(L\subset M\) carries
\[
T^i=b^i\,\delta_L,
\qquad
R^{ij}=K^{ij}\,\delta_L,
\]
then the net configurational force per unit length along \(L\) is
\[
f_A
=
\int_L R^A
=
b^i\,\Sigma_{iA}+K^{ij}\,M_{jiA}.
\]
The first term reproduces the classical Peach–Koehler force, and the second is a couple-stress correction. The worked example in the paper shows a time-dependent micro-rotation field generating torsion and curvature, which then produce defect excitations, force stress, and finally a localized configurational force that oscillates in space and decays as \(e^{-t}\) [2604.12984].

This formulation explicitly separates defect densities from classical compatibility. Torsion and curvature are treated as independent primary fields rather than constrained to vanish. The theory therefore places discrete configurational forces within a broader geometric setting in which they are not only post-processed outputs but also Noether currents linked to defect transport [2604.12984].

## 7. Statistical-mechanical interpretation and common distinctions

In jammed granular materials, “discrete configurational forces” take a different form. The Force Network Ensemble fixes particle positions and treats the contact forces \(f_{ij}\ge 0\) as the degrees of freedom, subject to force balance on every grain and fixed global stress trace. The allowed set of force configurations is a convex polytope for frictionless packings, and the phase-space volume is
\[
\Omega(P)
=
\int_{\{f_{ij}\ge 0\}}
\Biggl[\prod_{\langle i,j\rangle} df_{ij}\Biggr]
\prod_{k=1}^{dN}\delta\!\Bigl(\sum_j f_{kj}\mathbf n_{kj}\Bigr)\,
\delta\!\bigl(\mathrm{Tr}\sigma-NP\bigr).
\]
The configurational entropy is then
\[
S_f=k_B\ln\!\bigl[\Omega(P)/h_f^{N_f}\bigr],
\]
where \(h_f\) is the elementary force-cell volume, directly analogous to Planck’s constant in ordinary phase space [1912.03504].

The direct-measurement protocol proceeds by generating jammed packings, constructing the rigidity matrix, computing the left nullspace to obtain the \(\Delta Z+1\) states of self stress, intersecting the positivity cone with the unit sphere to obtain a convex spherical polytope of force solutions, and measuring its \(\Delta Z\)-dimensional volume \(V_f(P)\). The paper reports
\[
V_f(P)\approx C\,\gamma^{\Delta Z},
\qquad
\Omega(P)\approx P\,\frac{V_f(P)}{h_f^{\Delta Z}},
\]
and matches the resulting microscopic entropy with a macroscopic entropy obtained from angoricity via overlapping histograms [1912.03504].

Taken together, these results suggest that discrete configurational forces are not a single formula but a family of energy-based descriptors tied to changes in material configuration. One common misconception is to identify them with ordinary force balance; the cited works instead distinguish spatial equilibrium from material driving forces. Another misconception is that geometric admissibility alone determines evolution: the blocked-slip study shows a marked decoupling between geometric metrics and the configurational-force response, and the domain-wall study shows that lattice discreteness can create metastable pinning even when the smooth continuum force would predict monotonic expansion or shrinking [2603.24129] [2512.20453].

A further common feature is that thresholds and localization bands are problem dependent. In fracture, \(|F_{\rm CNF}/s|\) must exceed a material-specific threshold \(F_c\). In multistable metamaterials, pinning requires \(\varepsilon>\varepsilon_{\rm cr}(\Delta)\) and occurs only over a finite interval of wall sizes. In blocked slip, the paper proposes a future critical threshold \(J_c\) for slip transfer or crack nucleation. This suggests a general pattern: discrete configurational forces provide the energetic ranking and directional information, while the onset of actual evolution requires constitutive or material-specific criteria supplied by the particular application [2507.12247] [2512.20453] [2603.24129].

Source: https://www.emergentmind.com/topics/discrete-configurational-forces