From Echoes to MeanFieldHom

MeanFieldHom (MFH) is a Julia port of Echoes [57], the C++ mean-field homogenization library (with a Python interface) developed at Cerema. If you already have Echoes scripts, this page is a direct translation guide: the same RVE/scheme/homogenize concepts, a line-by-line API table, and a worked example that calls both implementations — Echoes live from Julia via PyCall.jl — to check that they agree.

The @example blocks run at build time; the Python and PyCall blocks are shown for reference, and the Echoes numbers below were captured by running them against Echoes 1.0.

API translation table

Both libraries share the same mental model: an RVE holding a matrix phase and inclusion/crack phases, each with a geometry and property dictionary, homogenized by a named scheme. Only the surface syntax differs.

Two conventions worth flagging for a smooth transition:

• Both stiff_kmu/iso_stiffness take physical (k, μ). Where they differ is in what the resulting tensor stores: MFH's raw TensISO{3}(a, b) constructor takes the pair (3k, 2μ), not (k, μ) — see the first tutorial. Building with iso_stiffness(k, μ) (as in the table above) sidesteps this entirely.
• Symbol-string vs. type-instance schemes. Echoes selects a scheme via a module-level constant (MT, SC, …) passed to homogenize; MFH uses a scheme type instance (MoriTanaka(), SelfConsistent(; kwargs...)), which is also how solver options (abstol, maxiters, select_best) attach directly to the scheme rather than as loose homogenize keywords. A Symbol shortcut (:mt, :sc, …) is also accepted — see the schemes manual for the full alias table.

Same problem, both sides

The classic porous benchmark — a solid matrix with spherical pores, porosity $\varphi \in [0, 1]$ — makes the translation concrete. In Echoes (the classic porous benchmark from the Echoes book):

from echoes import *

ks, mus = 72., 32.
kp, mup = 1.e-6, 1.e-6

myrve = rve(matrix="SOLID")
myrve["SOLID"] = ellipsoid(shape=spheroidal(1.), symmetrize=[ISO],
                            prop={"C": stiff_kmu(ks, mus)})
myrve["PORE"] = ellipsoid(shape=spheroidal(1.), symmetrize=[ISO],
                           prop={"C": stiff_kmu(kp, mup)})

def Chom_porous(phi, scheme):
    myrve["PORE"].fraction = phi
    myrve["SOLID"].fraction = 1. - phi
    C = homogenize(prop="C", rve=myrve, scheme=scheme,
                   epsrel=1.e-6, maxnb=300, select_best=True)
    return max(C.k, 0.), max(C.mu, 0.)

The same problem, live, in MFH:

using MeanFieldHom
using TensND
using LinearAlgebra
using Plots
gr()  # headless backend; GKSwstype is set to "100" in make.jl

const ks, mus = 72.0, 32.0
const kp, mup = 1.0e-6, 1.0e-6
C_s = iso_stiffness(ks, mus)
C_p = iso_stiffness(kp, mup)

function build_rve(φ)
    r = RVE(:SOLID)
    add_matrix!(r, Spheroid(1.0), Dict(:C => C_s); symmetrize = IsoSymmetrize())
    add_phase!(r, :PORE, Spheroid(1.0), Dict(:C => C_p); fraction = φ, symmetrize = IsoSymmetrize())
    return r
end

φs = collect(0.0:0.1:1.0)
k_mt = [k_mu(homogenize(build_rve(φ), MoriTanaka(), :C))[1] for φ in φs]
k_mt[1:4]

4-element Vector{Float64}:
 72.0
 55.44385079613372
 43.065421368241765
 33.460581867908616

Line for line, the two snippets build the same RVE and call the same scheme — this is the translation to internalize for every other Echoes script.

Calling Echoes from Julia with PyCall

To cross-check MFH against a live Echoes instead of hand-copied numbers, PyCall.jl can call Echoes directly from a Julia session. This requires Echoes installed in some Python environment — any interpreter where import echoes succeeds will do, conda or otherwise:

# One-time setup: point PyCall at the interpreter where `import echoes`
# succeeds, then rebuild the binding.
ENV["PYTHON"] = raw"/path/to/that/python"
using Pkg
Pkg.build("PyCall")

using PyCall

py"""
from echoes import *

def run_echoes(phi, scheme_name):
    _SCHEMES = {"MT": MT, "SC": SC, "DIFF": DIFF}
    scheme = _SCHEMES[scheme_name]
    myrve = rve(matrix="SOLID")
    myrve["SOLID"] = ellipsoid(shape=spheroidal(1.), symmetrize=[ISO],
                                prop={"C": stiff_kmu(72., 32.)})
    myrve["PORE"] = ellipsoid(shape=spheroidal(1.), symmetrize=[ISO],
                               prop={"C": stiff_kmu(1e-6, 1e-6)})
    myrve["PORE"].fraction = phi
    myrve["SOLID"].fraction = 1. - phi
    C = homogenize(prop="C", rve=myrve, scheme=scheme,
                   epsrel=1e-6, maxnb=300, select_best=True)
    return max(C.k, 0.), max(C.mu, 0.)
"""

k_e, μ_e = py"run_echoes"(0.3, "MT")

Running this snippet against Echoes 1.0 returns (33.46058186790858, 17.626741619570563) at φ = 0.3, scheme = "MT", matching the reference table below to solver tolerance. This is the same pattern used throughout scripts/bench_echoes/ (see Production cross-checks) — a py"""...""" block defining a small helper, called from Julia like an ordinary function.

The benchmark result

The table below reproduces the captured Echoes 1.0 output for the porous benchmark above (ks, μs = 72, 32; φ = 0.0, 0.1, …, 1.0), obtained exactly as shown in the previous section. It is embedded here as literal data — reading it back requires no Echoes installation — so this page renders identically whether or not Echoes is available:

# Captured from a run of Echoes 1.0 (see the run_echoes snippet above),
# once per (φ, scheme) pair.
echoes_ref = Dict(
    "MT" => (
        k = [72.0, 55.443851, 43.065421, 33.460582, 25.791046, 19.525425, 14.31056, 9.90258, 6.127661, 2.858562, 1.0e-6],
        μ = [32.0, 26.415585, 21.685158, 17.626742, 14.106633, 11.024391, 8.303101, 5.882864, 3.716342, 1.765626, 1.0e-6],
    ),
    "DIFF" => (
        k = [72.0, 54.60602, 40.616795, 29.414749, 20.528772, 13.593911, 8.323876, 4.491934, 1.917566, 0.457355, 1.0e-6],
        μ = [32.0, 26.168105, 20.864777, 16.11255, 11.933231, 8.347851, 5.376596, 3.038716, 1.352384, 0.334464, 1.0e-6],
    ),
    "SC" => (
        k = [72.0, 53.608383, 37.156454, 22.689248, 10.27169, 0.076007, 6.0e-6, 3.0e-6, 2.0e-6, 1.0e-6, 1.0e-6],
        μ = [32.0, 25.866249, 19.629147, 13.264365, 6.737918, 0.056937, 5.0e-6, 2.0e-6, 2.0e-6, 1.0e-6, 1.0e-6],
    ),
)

MFH, computed live, against those captured references:

mfh_schemes = Dict(
    "MT" => MoriTanaka(),
    "DIFF" => DifferentialScheme(; nsteps = 300),
    "SC" => SelfConsistent(; abstol = 1.0e-6, maxiters = 300, select_best = true),
)

_relerr(a, b) = abs(b) < 1.0e-9 ? abs(a - b) : abs(a - b) / abs(b)

plt = plot(;
    xlabel = "φ (porosity)", ylabel = "k_hom",
    legend = :topright, framestyle = :box, size = (760, 480),
)
for (name, color) in (("MT", :black), ("DIFF", :orange), ("SC", :red))
    k_jl = [max(k_mu(best_fit_iso(homogenize(build_rve(φ), mfh_schemes[name], :C)))[1], 0.0) for φ in φs]
    plot!(plt, φs, k_jl; label = "MFH $name", color = color, lw = 2)
    scatter!(plt, φs, echoes_ref[name].k; label = "Echoes $name", color = color, marker = :circle, markersize = 4)
end
plt

for name in ("MT", "DIFF", "SC")
    scheme = mfh_schemes[name]
    println("== $name ==")
    for (i, φ) in enumerate(φs)
        k_jl, _ = k_mu(best_fit_iso(homogenize(build_rve(φ), scheme, :C)))
        k_jl = max(k_jl, 0.0)
        k_e = echoes_ref[name].k[i]
        println("  φ=", round(φ, digits = 2), "  k_MFH=", round(k_jl, digits = 6),
                "  k_Echoes=", k_e, "  relerr=", round(_relerr(k_jl, k_e), sigdigits = 3))
    end
end

== MT ==
  φ=0.0  k_MFH=72.0  k_Echoes=72.0  relerr=1.97e-16
  φ=0.1  k_MFH=55.443851  k_Echoes=55.443851  relerr=3.68e-9
  φ=0.2  k_MFH=43.065421  k_Echoes=43.065421  relerr=8.55e-9
  φ=0.3  k_MFH=33.460582  k_Echoes=33.460582  relerr=3.95e-9
  φ=0.4  k_MFH=25.791046  k_Echoes=25.791046  relerr=7.53e-9
  φ=0.5  k_MFH=19.525425  k_Echoes=19.525425  relerr=1.07e-8
  φ=0.6  k_MFH=14.31056  k_Echoes=14.31056  relerr=5.32e-9
  φ=0.7  k_MFH=9.90258  k_Echoes=9.90258  relerr=1.42e-8
  φ=0.8  k_MFH=6.127661  k_Echoes=6.127661  relerr=6.19e-8
  φ=0.9  k_MFH=2.858562  k_Echoes=2.858562  relerr=6.34e-8
  φ=1.0  k_MFH=1.0e-6  k_Echoes=1.0e-6  relerr=1.64e-8
== DIFF ==
  φ=0.0  k_MFH=72.0  k_Echoes=72.0  relerr=1.97e-16
  φ=0.1  k_MFH=54.611188  k_Echoes=54.60602  relerr=9.46e-5
  φ=0.2  k_MFH=40.631978  k_Echoes=40.616795  relerr=0.000374
  φ=0.3  k_MFH=29.43971  k_Echoes=29.414749  relerr=0.000849
  φ=0.4  k_MFH=20.56085  k_Echoes=20.528772  relerr=0.00156
  φ=0.5  k_MFH=13.62948  k_Echoes=13.593911  relerr=0.00262
  φ=0.6  k_MFH=8.359112  k_Echoes=8.323876  relerr=0.00423
  φ=0.7  k_MFH=4.523194  k_Echoes=4.491934  relerr=0.00696
  φ=0.8  k_MFH=1.94151  k_Echoes=1.917566  relerr=0.0125
  φ=0.9  k_MFH=0.470894  k_Echoes=0.457355  relerr=0.0296
  φ=1.0  k_MFH=1.0e-6  k_Echoes=1.0e-6  relerr=1.66e-10
== SC ==
  φ=0.0  k_MFH=72.0  k_Echoes=72.0  relerr=1.97e-16
  φ=0.1  k_MFH=53.60838  k_Echoes=53.608383  relerr=6.35e-8
  φ=0.2  k_MFH=37.156441  k_Echoes=37.156454  relerr=3.49e-7
  φ=0.3  k_MFH=22.689235  k_Echoes=22.689248  relerr=5.7e-7
  φ=0.4  k_MFH=10.27167  k_Echoes=10.27169  relerr=1.99e-6
  φ=0.5  k_MFH=0.076512  k_Echoes=0.076007  relerr=0.00665
  φ=0.6  k_MFH=6.0e-6  k_Echoes=6.0e-6  relerr=0.0476
  φ=0.7  k_MFH=3.0e-6  k_Echoes=3.0e-6  relerr=0.0222
  φ=0.8  k_MFH=2.0e-6  k_Echoes=2.0e-6  relerr=0.0685
  φ=0.9  k_MFH=1.0e-6  k_Echoes=1.0e-6  relerr=0.35
  φ=1.0  k_MFH=1.0e-6  k_Echoes=1.0e-6  relerr=0.0

Production cross-checks

scripts/bench_echoes/ runs the same PyCall pattern systematically; tolerances and measured agreement are in Cross-validation, timings in Performance vs Echoes.

Without Echoes or PyCall, each cross-check has a PyCall-free counterpart: committed *_python.json dumps, and pure-Julia transcriptions such as scripts/55_ageing_creep_dirichlet_chains.jl.