mirror of
https://github.com/JuliaFEM/JuliaFEM.jl.git
synced 2026-08-06 04:21:33 +00:00
chore(test): delete minimal THME DOF smoke test
Remove unused thermo-hydro-mechanical-electric DOF scaffolding. - Drop `test/dofs/test_thme_minimal.jl`.
This commit is contained in:
Jukka Aho
@@ -1,118 +0,0 @@
|
||||
"""
|
||||
🎯 THE SIMPLEST POSSIBLE MULTI-PHYSICS: 50 Lines of Real Code
|
||||
|
||||
This demonstrates that COMPLEX PHYSICS ≠ COMPLEX CODE in JuliaFEM.
|
||||
|
||||
Physics: Fully coupled THME (Thermo-Hydro-Mechanical-Electric):
|
||||
- 4 fields (T, u, p, φ) on 4 entity types (Vertex, Vertex, Cell, Edge)
|
||||
- 6 bidirectional couplings = 12 coupling blocks
|
||||
- Applications: Geothermal, nuclear waste, CO2 sequestration, piezoelectric sensors
|
||||
|
||||
The ENTIRE implementation: < 100 lines (including comments!)
|
||||
"""
|
||||
|
||||
using JuliaFEM, Test, Tensors, LinearAlgebra, SparseArrays
|
||||
|
||||
@testset "🎯 Minimal Multi-Physics" begin
|
||||
# 1. Create mesh (2 tetrahedra, 5 nodes)
|
||||
nodes = [Vec{3}((0.,0.,0.)), Vec{3}((1.,0.,0.)), Vec{3}((0.5,1.,0.)),
|
||||
Vec{3}((0.5,0.5,1.)), Vec{3}((1.5,0.5,0.5))]
|
||||
connectivity = [(UInt32(1),UInt32(2),UInt32(3),UInt32(4)),
|
||||
(UInt32(2),UInt32(3),UInt32(4),UInt32(5))]
|
||||
mesh = Mesh{Tetrahedron{4}}(nodes, connectivity)
|
||||
|
||||
# 2. Define multi-field element (ONE line!)
|
||||
S = @DOFSet{T::DOF{Temperature,Vertex}, u::DOF{Displacement{3},Vertex},
|
||||
p::DOF{Pressure,Cell}, φ::DOF{ElectricPotential,Edge}}
|
||||
|
||||
# 3. Create elements (TWO lines!)
|
||||
mgr = DOFManager(mesh)
|
||||
register_fields!(mgr, S)
|
||||
elements = create_elements!(mgr, Element{Tetrahedron{4}, Lagrange{1}, S})
|
||||
|
||||
# 4. Initialize global system
|
||||
K = spzeros(mgr.total_dofs, mgr.total_dofs)
|
||||
F = zeros(mgr.total_dofs)
|
||||
|
||||
# 5. Material parameters (real physics!)
|
||||
params = (κ=1.0, E=10.0, ν=0.25, k=1.0, σ_e=1.0, # Diagonal blocks
|
||||
α_T=1e-5, α_p=1e-3, β_T=1e-6, ζ=1e-7, S=1e-6, e=1e-8) # Couplings
|
||||
|
||||
# 6. Assembly loop (ONE element assembles ALL physics!)
|
||||
for (idx, elem) in enumerate(elements)
|
||||
# Get local structure
|
||||
n_loc = local_dof_count(elem)
|
||||
T_loc, u_loc, p_loc, φ_loc = field_dof_range(elem,:T), field_dof_range(elem,:u),
|
||||
field_dof_range(elem,:p), field_dof_range(elem,:φ)
|
||||
K_loc, F_loc = zeros(n_loc, n_loc), zeros(n_loc)
|
||||
|
||||
# Get geometry and basis
|
||||
X = ntuple(i -> nodes[mesh.connectivity[idx][i]], 4)
|
||||
ip = integration_points(Gauss{1}(), Tetrahedron{4}())[1]
|
||||
dN_dξ = get_basis_derivatives(Tetrahedron{4}(), Lagrange{1}(), ip.ξ)
|
||||
J = sum(X[i] ⊗ dN_dξ[i] for i in 1:4)
|
||||
∇N = ntuple(i -> transpose(inv(J)) ⋅ dN_dξ[i], 4)
|
||||
vol = det(J) * ip.weight
|
||||
|
||||
# Physics assembly (simplified but REAL!)
|
||||
λ, μ = params.E*params.ν/((1+params.ν)*(1-2params.ν)), params.E/(2*(1+params.ν))
|
||||
C = λ*one(Tensor{2,3})⊗one(Tensor{2,3}) + μ*(otimesu(one(Tensor{2,3}),one(Tensor{2,3}))+
|
||||
otimesl(one(Tensor{2,3}),one(Tensor{2,3})))
|
||||
|
||||
# Diagonal blocks (4 physics)
|
||||
for i in 1:4, j in 1:4
|
||||
K_loc[T_loc[i],T_loc[j]] += params.κ * (∇N[i]⋅∇N[j]) * vol # Thermal
|
||||
end
|
||||
for k in 1:4, l in 1:4, α in 1:3, β in 1:3
|
||||
B_kα = 0.5*(∇N[k]⊗basevec(Vec{3},α)+basevec(Vec{3},α)⊗∇N[k])
|
||||
B_lβ = 0.5*(∇N[l]⊗basevec(Vec{3},β)+basevec(Vec{3},β)⊗∇N[l])
|
||||
K_loc[u_loc[3(k-1)+α],u_loc[3(l-1)+β]] += dcontract(B_kα,dcontract(C,B_lβ))*vol # Mechanical
|
||||
end
|
||||
K_loc[p_loc[1],p_loc[1]] += params.k * vol # Hydraulic
|
||||
for i in 1:6; K_loc[φ_loc[i],φ_loc[i]] += params.σ_e*vol/6; end # Electric
|
||||
|
||||
# Coupling blocks (6 bidirectional = 12 blocks!)
|
||||
c_T = params.α_T*params.E/(1-2params.ν)
|
||||
for i in 1:4, k in 1:4, α in 1:3 # T↔u (thermal expansion)
|
||||
v = c_T*(∇N[i]⋅basevec(Vec{3},α))*norm(∇N[k])*vol
|
||||
K_loc[T_loc[i],u_loc[3(k-1)+α]] += v; K_loc[u_loc[3(k-1)+α],T_loc[i]] += v
|
||||
end
|
||||
for k in 1:4, α in 1:3 # u↔p (Biot poroelasticity)
|
||||
v = params.α_p*norm(∇N[k])*vol/3
|
||||
K_loc[u_loc[3(k-1)+α],p_loc[1]] += v; K_loc[p_loc[1],u_loc[3(k-1)+α]] += v
|
||||
end
|
||||
for i in 1:4 # T↔p (thermal pressurization)
|
||||
v = params.β_T*norm(∇N[i])*vol
|
||||
K_loc[T_loc[i],p_loc[1]] += v; K_loc[p_loc[1],T_loc[i]] += v
|
||||
end
|
||||
for j in 1:6 # p↔φ (electro-osmotic)
|
||||
v = params.ζ*vol/6
|
||||
K_loc[φ_loc[j],p_loc[1]] += v; K_loc[p_loc[1],φ_loc[j]] += v
|
||||
end
|
||||
for i in 1:4, j in 1:6 # T↔φ (Seebeck/Peltier)
|
||||
v = params.S*norm(∇N[i])*vol/6
|
||||
K_loc[T_loc[i],φ_loc[j]] += v; K_loc[φ_loc[j],T_loc[i]] += v
|
||||
end
|
||||
for k in 1:4, α in 1:3, j in 1:6 # u↔φ (piezoelectric)
|
||||
v = params.e*norm(∇N[k])*vol/6
|
||||
K_loc[u_loc[3(k-1)+α],φ_loc[j]] += v; K_loc[φ_loc[j],u_loc[3(k-1)+α]] += v
|
||||
end
|
||||
|
||||
# Scatter to global (type-safe!)
|
||||
dof_map = local_to_global_map(elem)
|
||||
for i in 1:n_loc, j in 1:n_loc; K[dof_map[i],dof_map[j]] += K_loc[i,j]; end
|
||||
end
|
||||
|
||||
# 7. Solve (standard linear algebra)
|
||||
for dof in [1, count_field_dofs(mgr,:T)+1:count_field_dofs(mgr,:T)+3...] # Fix node 1
|
||||
K[dof,:] .= 0; K[:,dof] .= 0; K[dof,dof] = 1; F[dof] = 0
|
||||
end
|
||||
F[5:end] .= 0.01 # Apply load
|
||||
sol = K \ F
|
||||
|
||||
@test all(isfinite.(sol))
|
||||
println("✅ Solved $(mgr.total_dofs)-DOF fully coupled THME system!")
|
||||
println(" 12 coupling blocks assembled in ~80 lines of code")
|
||||
println(" Applications: Geothermal, nuclear waste, smart materials")
|
||||
println("\n🎯 KEY INSIGHT: Complex physics ≠ Complex code!")
|
||||
end
|
||||
Reference in New Issue
Block a user