Projecting gradient

General
1

I have potential values at the centroid of a DG0 space. I found the electric field for the same. if I have n number of centroids, the electric field components are n*3( for 100 centroid, the potential values obtained are 300). A minimal code would be. I have potentials for my mesh from a different software which is stored as phi_values.

        from dolfin import *
        mesh = UnitSquareMesh(10, 10)

         V = functionSpace(mesh.'CG',1)
         proj = Function(V)
        proj.vector().set_local(phi_values)
        phi0 =  proj.vector().get_local()

        subV = FunctionSpace(mesh, "DG", 0) 
        subproj = project(proj, subV, solver_type="cg", preconditioner_type="amg")
        coords= subV.tabulate_dof_coordinates()   
        phi_cen = subproj.vector().get_local()
                         
        # Finding EF 

        V1 =  VectorFunctionSpace(self.geometry.mesh,'DG',0)
        V_phi = V
        phi_values = self.phi.vector().get_local()
        phi = Function(V_phi)
        phi.vector().set_local(phi_values)
        grad_phi = project(grad(phi), V,solver_type="cg", preconditioner_type="amg")
        E = Function(V1)
        E.vector().set_local(-grad_phi.vector().get_local())
        E_values = E.vector().get_local()
        E_values = np.column_stack(E_values)
        File("electric_field.pvd") << E

I could get the electric field, which is n3 potential values. I am confused about the values. How do I interpret it ? i mean electric field for a 3D space has 3 components, so at each centroid, it should have 3 components , so the shape should be (n,3) rather than n3.

2

Your code is with many mistakes: indentation, copied from a class with self keyword, upper and lower case typos, etc… Make sure the code is reproducible.

To answer your question: Function is a class which contains a data array of PETSc.Vec. The PETSc.Vec is a concept of vector, which only involves one-dimensional data. So it has n3 values (we say the size should be 3n) in your case but not multidimensional (n, 3) data. If you want to extract each component of it, you may want to use Function.sub(d) where d traverses from 0 to 2 in 3D space.

See for documentation:
https://fenicsproject.org/olddocs/dolfin/2019.1.0/python/_autogenerated/dolfin.cpp.function.html?highlight=functionspace#dolfin.cpp.function.Function.sub

3

so if the function space is CG with order one, it should have multidimensional data?

4

No, Function always has a vector (PETSc.Vec) data which only has a one-dimensional data concept.

5

So how do I find the electric field values corresponding to each cell in the mesh ? I was trying to get the electric field at centroids of the mesh and if I am correct I get the ef at centroids. How do I extract it for each cell from 3n?

6

I am only familar with the FEniCSx part. Within a specific partition in a parallel run, you can traverse all the cells and extract the associated DOFs cell_dofs with a poosible implementation:

num_cells = msh.topology.index_map(
        msh.topology.dim).size_local + msh.topology.index_map(
        msh.topology.dim).num_ghosts

cells = np.arange(num_cells)
cell_dofs = []

for cell in cells:
    cell_dofs.append(V.dofmap.cell_dofs(cell))

For FEniCS implementation, see the documentation: https://fenicsproject.org/olddocs/dolfin/2019.1.0/python/_autogenerated/dolfin.cpp.fem.html?highlight=cell_dof#dolfin.cpp.fem.DofMap.cell_dofs