How to split subdomains according to numpy array?

mesh
1

I have a 3D numpy array consisting of 0 and 1. And I want according to this array create two subdomain using mesh function and SubDomain. I think that it can be done with this, but I don’t understand how to implement this to numpy array.
Assumed cubic mesh.

2

Hey, Olga!
What is the shape of your 3D array? Do you know how to find it? I mean, analytical expressions for the shape. If yes, pls, share with me

And do you need to use SubDomain and MeshFunction? Tell more about what you want. If you’d like to set coefficient kappa (the one in your link), I think, it can be done not using these two - at least, I know how to do it in 2d case

3

Thanks for the answer.

I do not have to use SubDomain and MeshFunction

I have cube shape array (nx, ny, nz) which represent semiconductor material in capacitor. The value in each cell of array indicates the presence or absence of an electron in the trap. So, it can be any array containing 0 and 1, for example:

import numpy as np
nx = 10
ny = 10
nz = 10
nparray = np.random.choice(2, (nx,ny,nz),  p=[0.95, 0.05])

Now, I want to create a mesh with same dimensions (nx,ny,nz) and use array as source function when solving PDE.

It would be nice if you show how to do it in 2d case.

4

Consider this

import numpy as np
from dolfin import *

nx = ny = nz = 10
locs = np.random.choice(2, (nx,ny,nz),  p=[0.95, 0.05])

mesh = UnitCubeMesh(nx-1, ny-1, nz-1)
V = FunctionSpace(mesh, 'P', 1)
f = Function(V)

v2d = vertex_to_dof_map(V)
f.vector()[v2d] = locs.flatten()

# or to assign some specific values in electron locations
some_value = 100.
f.vector()[v2d[locs.flatten() > 0.]] = 100.

When I said analytic expressions, I meant the relation between your array shape = (1) and mesh dimensions; or (1) and number of elements along dimensions in your mesh, or (1) and number of vertices and so on. The counterpart in these relations depends on the approach you choose to implement the solution. I don’t know much about your intentions, so I took the simplest case, which may be not relevant to your problem. This is Lagrange finite element of 1st degree, where node values are set in the vertices: i.e. mesh_size_along_x_axis = nx - 1. In other functional spaces, or using different approaches (via MeshFunction), the relation will be different, so, this fact is worth considering

5

Oh, I need the array element to be in the central node of the cube for Lagrange finite element of 2st degree. But, I can use the grid two times denser and redefine my numpy array according to this.
Thank you very much for your help.

6

I find that in fenics array collapsed into one dimension in another order, so, you suppose to use locs.flatten('F') instead locs.flatten().

7

Interesting. Have you managed to solve your problem?

8

Almost, I have a solution for a cubic grid, now I need to change it a little and the problem will be solved.