Hi,
I have a mesh containing two material subdomains, which I will call mechanical and exterior. The displacement test function is only defined on the mechanical subdomain.
I want to assemble
where:
- \partial\Omega_{c} is the internal interface between the two materials
- \boldsymbol{\sigma}_{\mathrm{ext}} is evaluated on the exterior side
- \boldsymbol{n}_{\mathrm{mech}} points outwards from the mechanical subdomain
- \boldsymbol{v}_{\mathrm{mech}} is the mechanical test-function
The stress is a material-wise discontinuous coefficient function. The relevant part of my current implementation is:
mechanical_other = is_mechanical.Other()
normal = (is_mechanical - mechanical_other) * specialcf.normal(mesh.dim)
transmission_pair = (
is_mechanical * is_exterior.Other()
+ mechanical_other * is_exterior
)
exterior_stress = (
is_mechanical * stress.Other()
+ mechanical_other * stress
)
test_on_mechanical = (
is_mechanical * v
+ mechanical_other * v.Other()
)
rhs = LinearForm(fes)
rhs += (
transmission_pair
* InnerProduct(exterior_stress * normal, test_on_mechanical)
* dx(skeleton=True, bonus_intorder=6)
)
rhs.Assemble()
Here, is_mechanical and is_exterior are material-wise cutoff coefficient functions, and fes is defined only on the mechanical subdomain.
As a check, I created a separate mesh containing only the mechanical domain. The interface then became a genuine external boundary, so I could impose the known traction directly using an ordinary boundary integral. The displacement obtained from that reference problem was very close to the result from the skeleton expression above.
I also compared the direct dx(skeleton=True) expression with:
rhs += SymbolicLFI(
transmission_pair
* InnerProduct(exterior_stress * normal, test_on_mechanical),
VOL,
skeleton=True,
bonus_intorder=6,
)
The two expressions produced identical assembled right-hand-side vectors.
Could someone confirm whether the first expression is the recommended way to combine the exterior stress trace, outward mechanical normal and mechanical test-function trace?
In particular, does including both facet orientations in transmission_pair correctly select the required traces, or does it count the interface twice?
Thanks,
Yash
)