How to define spatially varying materials

Hi. I am trying to learn to how to define the materials of a 3D model on FENICSx platform, but with a property of spatially varying elastic modulus. I have created the material file .inp by Abaqus (every element of the model has been defined by its own different material paramters), But I don’t know how to do next.

Can you provide any effective ideas? Thank you very much!

See for instance: I/O from XDMF/HDF5 files in dolfin-x - #48 by dokken
or
Dolfinx discontinous expression - #3 by dokken

Thanks a lot, Prof Dokken! Let me learn first.

Hi, Prof Dokken. Sorry to bother you.
I have test the code you provided in Dolfinx discontinous expression

import dolfinx.io
import dolfinx
from mpi4py import MPI

mesh = dolfinx.mesh.create_unit_square(MPI.COMM_WORLD, 5, 5)
V = dolfinx.fem.FunctionSpace(mesh, ("DG", 0))
v = dolfinx.fem.Function(V)
x = V.tabulate_dof_coordinates()

for i in range(x.shape[0]):
    midpoint = x[i,:]
    if midpoint[0]> 0.5 and midpoint[1]>0.25:
        v.vector.setValueLocal(i, 2)
    else:
        v.vector.setValueLocal(i, 1)

The code use the midpoint to define the function space. I’m not sure whether the list of x(V.tabulate_dof_coordinate) in DG mesh space is arranged by the element index. Or is there any API of dolfinx to get the list of element index so that the corresponding list of material parameters can be generated and attached to the each fem element?

Thanks a lot!

The dofs of a DG-space is aligned with the local cell index.
If you get an unexpected result, I would suggest you add a print statement inside your loop to see what the midpoint is.
With the snippet you have supplied above, it seems like you have gotten a sensible result?

Thank you for your quick reply.
Since the dof is aligned with the cell index, I think I know how to do next.

By the way, I have tried the second method about MeshFunction mentioned in Dolfinx discontinous expression

import dolfinx.io
import dolfinx
from mpi4py import MPI

mesh = dolfinx.mesh.create_unit_square(MPI.COMM_WORLD, 5, 5)
V = dolfinx.fem.FunctionSpace(mesh, ("DG", 0))
v2 = dolfinx.Function(V)
mf = dolfinx.MeshFunction("size_t", mesh, mesh.topology.dim, 0)
mf.mark(lambda x: numpy.logical_and(x[0] < 0.5, x[1]>0.2), 2)
for i in range(x.shape[0]):
    if mf.values[i] == 2:
        v2.vector.setValueLocal(i, 2)
    else:
        v2.vector.setValueLocal(i, 1)

But it came to an Error about the API of MeshFunction in the dolfinx library.

AttributeError: module 'dolfinx' has no attribute 'MeshFunction'

Could you please tell me where it did go wrong? Thanks!

The post is a bit out of date, as it is from 2020. MeshFunction has been replaced with dolfinx.mesh.meshtags. see for instance my demos for examples: Defining subdomains for different materials — FEniCSx tutorial

1 Like

Thank you! So kind of you.

Sorry, Prof Dokken. I’ve tried to define the material by the API of tabulate_dof_coordinates and setValueLocal. But the results are mismatched. The coordinate array given by x = V.tabulate_dof_coordinates() is not aligned with the element index defined by the xdmf file of model.

How can I resolve this problem? Thanks !

Please create a minimal reproducible example, as I cannot say what goes wrong from just looking at a screenshot of the result.

Hi, Prof.Dokken. A minimal reproducible example of a cube model containing 8 elments (Fig.1 shows the hexahedral element id ) is given as below.
In test #1, I was planning to assign the material parameters [1,2,3,4,5,6,7,8] to the corresponding element #1,#2,…,#8 . But it came to the mismatch(Fig.2).


When the element id is rearranged (Fig.3) in .msh file, the assignment of parameters is also mismatched (Fig.4).


Here is the code snippet

para= [1,2,3,4,5,6,7,8]
with dolfinx.io.XDMFFile(MPI.COMM_WORLD, 'test1.xdmf', "r") as xdmf:
    mesh = xdmf.read_mesh(name="Grid")

V = dolfinx.fem.FunctionSpace(mesh, ("DG", 0))
v = dolfinx.fem.Function(V)
x = V.tabulate_dof_coordinates()

for i in range(x.shape[0]):
    v.vector.setValueLocal(i, para[i])

It seemed that the returned list given by x=V.tabulate_dof_coordinates() is not aligned with the element of input .msh file How can I get the right parameter assignment?
Hope to get your help once more.Thanks a lot!

You can access the original cell index as done here: Wrong order of dofs after creating mesh from lists of points and cells - #5 by dokken

Prof Dokken, I’ve learned that the order of nodes of input .xdmf file are truly rearranged. By geometry.input_global_indices, the mapping relationships of rearrangement for the nodes can be obtained. And the elements index are not changed.
3

So, the problem should be due to the space v(v=dolfinx.fem.Function). It seems that the v are not arranged with the order of element index after the loop v.vetor.setValueLocal. If so, I want to know how I could write the funtion v according to the right order.
Thanks a lot!

See for instance I/O from XDMF/HDF5 files in dolfin-x - #48 by dokken
and
Wrong order of dofs after creating mesh from lists of points and cells - #2 by dokken

Yes, I have tried to learn that case by the access you provided. The reordered mapping of nodes can truly be obtained. I also checked that the element order was still unchanged, the only changes happened in the nodes order.

But the question is that I have Initial information of value list for each element arranged by the element index. I don’t know why it went wrong when I assigned the values to the function and conducted the output. I thought I have found the answer before by the access How to assign different material properties to elements (Linear Elasticity) But after I tried to perform my adapted code. The results were unsatisfactory (Fig.4). It indeed assigned the value to the corresponding element, but the element was wrong. It actually changed the element id instead of assigning the value to the correct element.

I just don’t know why the value is mismatched with the mesh after I write output. Sorry, I’m extremely confused.

The msh file of this minimal example is attched here. I just want to assign the value [1,2,3,4,5,6,7,8] respectively to the 1st,2nd,3rd… element in the msh file.

Thanks a lot!

Consider the following:


import dolfinx.io.gmshio
import meshio
from mpi4py import MPI

def create_mesh(mesh, cell_type, prune_z=False):
    cells = mesh.get_cells_type(cell_type)
    points = mesh.points[:,:2] if prune_z else mesh.points
    out_mesh = meshio.Mesh(points=points, cells={cell_type: cells})
    return out_mesh

mm = meshio.read("mesh.msh")
m_out = create_mesh(mm, "hexahedron")
original_mesh_cells = m_out.points[m_out.cells[0].data]
meshio.write("mesh.xdmf", m_out)


para= [1,2,3,4,5,6,7,8]

with dolfinx.io.XDMFFile(MPI.COMM_WORLD, 'mesh.xdmf', "r") as xdmf:
    msh = xdmf.read_mesh(name="Grid")

V = dolfinx.fem.FunctionSpace(msh, ("DG", 0))
v = dolfinx.fem.Function(V)
x = V.tabulate_dof_coordinates()

o_cell_idx =  msh.topology.original_cell_index
for i in range(x.shape[0]):
    print(i, original_mesh_cells[o_cell_idx[i]], x[i])
    v.x.array[i] = para[o_cell_idx[i]]


with dolfinx.io.XDMFFile(msh.comm, "output.xdmf", "w") as xdmf:
    xdmf.write_mesh(msh)
    xdmf.write_function(v)
2 Likes

It works. That is absolutely what I am looking for and solves my problem! :sob:
Thank you very much!