Setting up multiple subdomains for the linear elasticity example

General
1

From multiple subdomains tutorial:

kappa = Function(Q)
cells_0 = locate_entities(mesh, mesh.topology.dim, Omega_0)
cells_1 = locate_entities(mesh, mesh.topology.dim, Omega_1)
kappa.x.array[cells_0] = np.full_like(cells_0, 1, dtype=default_scalar_type)
kappa.x.array[cells_1] = np.full_like(cells_1, 0.1, dtype=default_scalar_type)

The below usage is shown:

a = inner(kappa * grad(u), grad(v)) * dx

So if I am not mistaken Omega came from the heat problem in my case I have a linear elasticity problem with multiple materials.

From linear elasticity:

def epsilon(u):
    return ufl.sym(ufl.grad(u))  # Equivalent to 0.5*(ufl.nabla_grad(u) + ufl.nabla_grad(u).T)


def sigma(u):
    return lambda_ * ufl.nabla_div(u) * ufl.Identity(len(u)) + 2 * mu * epsilon(u)

So then, in terms of expanding the epsilon and sigma functions into multiple subdomains; Is there a way to check inside epsilon function for example which cell this currently is
so the correct mu and lambda to be applied based from which domain we now ave a list of cells for as per:

kappa.x.array[cells_0] = np.full_like(cells_0, 1, dtype=default_scalar_type)
kappa.x.array[cells_1] = np.full_like(cells_1, 0.1, dtype=default_scalar_type)

If not how might two different materials be expressed for a linear elasticity example of multiple subdomains? As per:

a = ufl.inner(sigma(u), epsilon(v)) * ufl.dx

So I guess at this point one can check co-ordinates inside epsilon and make a determination based from x,y and z where one is in a model, OK. Really what things will come down to that I need to know about is how to declare epsilon and sigma for multiple materials that will be coming in from a create gmsh model with multiple subdomains how to build one like that?

I was thinking maybe a better more simple way to phrase the question might be how to enumerate the vertices of a cell so that the x,y,z co-ordinate can be checked it it appears in cell 0 for example? which I know there must be a simple way to but it doesn’t just readily come to mind right now.

Seemingly I maybe came up with a way here:

but I will admit I feel some confusion still a little bit…

This is maybe why I feel a bit confused. From slab data post inside sigma function there is a usage:

μ * ufl.sym(grad(v))

there mu, μ is a function object where a constant is set to it and μ is then ultipled by ufl.sym…

So far I see that a Function as imported from dolfinx can be seemingly used inside epsilon and sigma functions in the same fashion as slab I made at sigma function, to multiply a Function by a ufl.somefunc. OK. fine…
Is that correct and legal use of Function?

it seems there also is a way to get x,y,z inside epsilon or sigma and check if that point is inside say a hexahedron but that is related to:

Format output of get_vertex_coordinates() method of the Cell DOLFIN class - mesh - FEniCS Project

In the case of hexahedron a function to isInside(point) one has to be created as far as I know… unless someone knows if ufl for example can tell me is x,y,z inside of an element generic or otherwise?

2

Why not create lambda_ and mu_ as functions of a DG 0 space, and assign the data as done for kappa. Then you can use dolfinx.mesh.locate_entities(mesh, mesh.topology.dim, your_function)
to get the cells that you want, where your_function is a function that takes in a set of coordinates (num_points, 3), and returns an array of booleans of shape (num_points). i.e.

def left_cells(x):
   return x[1] =<0.5
marked_cells = dolfinx.mesh.locate_entities(mesh, mesh.topology.dim, left_cells)

will give you all cells whose vertices satisfies y<=0.5

1 Like