How to restrict solutions or stresses to subdomains?

Please note that you should not use CG 1 function spaces for Von Mises stresses, as the stresses are a DG 0 function (if the displacement is CG 1. I.e.

VDG = FunctionSpace(mesh, "DG", 0)
von_Mises = project(von_Mises, V)

You can compute the stresses on a subdomain by projecting on a subdomain, as follows:

s = sigma(u) - (1./3)*tr(sigma(u))*Identity(d)
von_Mises = sqrt(3./2*inner(s, s))
u, v = TrialFunction(VDG), TestFunction(VDG)
a = inner(u, v)*dx
L = inner(von_Mises, v)*dx(1)
stress = Function(VDG)
solve(a==L, stress)

would only compute the stresses on domain 1, and be zero everywhere else.
This can be optimized in a similar fashion to How to extract values in a surface of 3D mesh to another 2D mesh? - #4 by dokken
where you would use dx(1) in the definition of a, then assemble it into A, and use ident_zeros.

Another solution would be to use MeshViews

1 Like