Exterior trace on an internal material interface using dx(skeleton=True)

General
1

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

\int_{\partial\Omega_{c}} \left(\boldsymbol{\sigma}_{\mathrm{ext}}\boldsymbol{n}_{\mathrm{mech}}\right) \cdot \boldsymbol{v}_{\mathrm{mech}}\,\mathrm{d}s,

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

2

I think yes (its a bit hard from reading code and not testing it :wink: )

I prefer this way:

interface = mesh.Materials("exterior").Boundaries() * mesh.Materials("mechanical").Boundaries()
LinearForm(BoundaryFromVolumeCF(mesh.MaterialCF({"exterior" : stress})) * specialcf.normal(mesh.Materials("mechanical")) * v * ds(interface)).Assemble()

Stress is a gridfunction of another solution or given function right?
specialcf.normal(domain) gives you outward pointing normal of domain
with BoundaryFromVolume you get (a random) volume side. With the inner MaterialCF you defiine the stress only on one side and then the chosen one is the one where it is definedon.
Best
Christopher

3

Thanks Christopher, that is very helpful!

Best,
Yash