1

I am writing a heat equation for a cube by fenics. The material will be changed during the time which is handled by a function(let say function_k) to update the property. The code is OK and now I want to parallel it. Everything is fine but the problem is the function_k which makes some problem. As I found till now the mesh will be partitioned based on mesh vertex coordinate. However, when I want to change function value I need to use nodal value which makes problem. The mesh, function space and partition parameter:

mesh = BoxMesh(Point(x0, y0, z0), Point(x1, y1, z1), nx, ny, nz)
parameters["mesh_partitioner"] = "SCOTCH"#"ParMETIS"s
V = FunctionSpace(mesh, 'CG', 1)

The material property"K" is initialised as:

tol = 1e-6
muz = 0.1
k  = Expression('x[2]>zf+tol ? k_a  : k_p' ,degree=0, zf=muz, k_a=20 , 
k_p=1, tol=tol)​

The function_k is(which will be called after solving the eq in each time step):

def function_k(mesh,f,k,lw,u,V):
    k_p  = 20       
    v2d = vertex_to_dof_map(V)
    coordinates = mesh.coordinates()
    k = interpolate(k,V)
    nodal_values_k  = k.vector()
    for i, x in enumerate(coordinates):
        if( (x[2] <= f.muz) and (x[2] > f.muz-lw)):
            nodal_values_k[v2d[i]]  = k_p
    k.vector()[:]  = nodal_values_k
    return(k)

It gives me in line : "nodal_values_k[v2d[i]] = k_p" following error: IndexError: expected indices to be in [0..106]

which is from v2d... means"vertex_to_dof_map" is not partitioned as desire....Is there any way to correct it?

Can we partitioned mesh based on degree of freedom(dof)?

Can we have "vertex_to_dof_map(V)" based on partitioned mesh?

In addition If I could write the function_k with vertex value as follow is there any way the vertex value_k to k without using dof to vertex map?

def Kfunction(mesh,f,k,lw,u,V):
    k_p  = 20       
    coordinates = mesh.coordinates()
    k = interpolate(k,V)
    vertex_values_k = k.compute_vertex_values(mesh)
    for i, x in enumerate(coordinates):
        if( (x[2] <= f.muz) and (x[2] > f.muz-lw)):
            vertex_values_k[i] = k_p

    #??? how should I assign the vertex value_k to k witouth using dof 
    # to vertex map
    return(k)
  • I don't know the answer to your question, but I recommend asking on https://www.allanswered.com/community/s/fenics-project/ for FEniCS-related questions. Many of the developers are active there, and I've found them to be helpful. – sigvaldm Mar 07 '18 at 11:11

0 Answers0