Ferrite backend
The material contract knows nothing about any finite-element library, so coupling to one needs almost no adapter. What Ferrite does not provide is per-quadrature-point material state, so that bookkeeping — and an element routine — is all this extension adds.
Activated by import Ferrite alone:
import Ferrite
const EXT = Base.get_extension(
MeanFieldHomogenization, :MeanFieldHomogenizationFerriteMaterialExt
)Why a second extension
MeanFieldHomogenizationFerriteExt serves the opposite coupling and needs FerriteGmsh and Gmsh too. A structural computation that only wants a homogenized material law has no reason to pull a ~100 MB gmsh artifact into its environment — or into a documentation build.
The helpers
mfh_states(mat, dh, cv) | one initial_state per quadrature point, as the nested [cell][qp] vector Ferrite loops expect |
mfh_element!(Ke, re, cv, mat, ue, states, states_old, Δt; cache) | element stiffness and internal force, small strain, plane strain |
mfh_poro_element!(Ke, re, cvu, cvp, mat, ue, pe, …, Δt, mobility; u_range, p_range, cache) | the coupled |
annulus_grid(Ri, Ro, nr, nθ) | a structured annular sector, built by bending a rectangle — no gmsh |
cylinder_sector_grid(Ri, Ro, H, nr, nθ, nz; grading) | its 3-D twin, with geometric radial layers for a well problem |
Keep two state arrays and swap only once a step has converged; material_response never mutates its argument, so a rejected Newton iteration is undone by not swapping.
states = EXT.mfh_states(mat, dh, cv)
states_old = EXT.mfh_states(mat, dh, cv)
cache = MaterialCache()
for cell in CellIterator(dh)
reinit!(cv, cell)
fill!(Ke, 0); fill!(re, 0)
ed = celldofs(cell)
EXT.mfh_element!(Ke, re, cv, mat, u[ed],
states[cellid(cell)], states_old[cellid(cell)], Δt; cache)
assemble!(assembler, ed, Ke, re)
endA worked mechanical model is the thick-walled cylinder; the coupled element drives the ARMA 2011 well test.
The two-field element takes the mobility
for c in eachindex(states), q in eachindex(states[c])
K = transport_property(mat, states[c][q])
mob[c][q] = K === nothing ? zero(SymmetricTensor{2,3}) : to_tensors(K) / μ
endTwo traps
A DofHandler numbers dofs by cell traversal, not node order — read nodal results with evaluate_at_grid_nodes(dh, u, :u), never by indexing u with a node number. And a MaterialCache is an unlocked Dict: use one per thread, or none.
Other codes
Only Ferrite is wired today. The contract itself is backend-agnostic, so a Gridap, FEniCSx or Abaqus UMAT driver needs the same two pieces — a per-point state container and an element routine — plus, outside Julia, a way to carry ε and σ across the boundary (voigt_strain and voigt_stress exist for that). Those are on the roadmap, not in the package.