# simple example on how to use dolfin::assemble_local

842

views

3

I would like to call dolfin FEM assemble routine for each element by providing 2d-elements and triangles (extracted from GMSH).

Does anyone have tried it and if yes can you provide a simple example code which can help me get started with it?

Thanks,

- AJ

Does anyone have tried it and if yes can you provide a simple example code which can help me get started with it?

Thanks,

- AJ

### 2 Answers

6

Here is a little code snippet which assembles a global tensors from local cell tensors using the assemble_local function.

It follows the logic of Alg. 2 in Ch. 6 of "Automated Solution of Differential Equations by the Finite Element Method".

The result is verified against the default fenics implementation by direct comparision.

```
# Define model problem
import numpy as np
import dolfin as df
mesh = df.RectangleMesh(df.Point(0,0), df.Point(1,1), 4, 4)
V = df.FunctionSpace(mesh, 'Nedelec 1st kind H(curl)', 2)
u = df.TrialFunction(V)
v = df.TestFunction(V)
form = df.inner(df.grad(u), df.grad(v)) * df.dx
# Assemble global matrix from local matrices
A = np.zeros([V.dim(), V.dim()])
for cell in df.cells(V.mesh()):
A_local = df.assemble_local(form, cell)
local_to_global_map = V.dofmap().cell_dofs(cell.index())
# Copy local matrix into global matrix
for (idx_loc_1, idx_glob_1) in enumerate(local_to_global_map):
for (idx_loc_2, idx_glob_2) in enumerate(local_to_global_map):
A[idx_glob_1, idx_glob_2] += A_local[idx_loc_1, idx_loc_2]
# Verification against fenics:
A_df = df.assemble(form)
A_dense = df.as_backend_type(A_df).mat().convert('dense')
assert np.allclose(A_dense.getDenseArray() - A, 0.), "Both matrices should be the same"
```

Hope this helps and leads you on the right track!

Oh, very nice indeed!

written
13 months ago by
pf4d

Thanks, Wagnandr. That's is working great!

Following that, how can we apply boundary conditions to the Global matrix 'A' after assembling the way you explained?

I tried "bc.apply(A)" but get the following error:

It works if directly applied on "A_direct = df.assemble(form)".

Appreciate your help.

-A

Following that, how can we apply boundary conditions to the Global matrix 'A' after assembling the way you explained?

I tried "bc.apply(A)" but get the following error:

*TypeError: in method 'DirichletBC_apply', argument 2 of type 'dolfin::GenericVector &'*It works if directly applied on "A_direct = df.assemble(form)".

Appreciate your help.

-A

written
11 months ago by
AJ

Hi wagnandr,

Any suggestion on how to

Thanks,

AJ

Any suggestion on how to

**apply Boundary Conditions**on**global assembly matrix (A)**assembled using above method?Thanks,

AJ

written
9 months ago by
AJ

0

With help from old docs :

https://fenicsproject.org/qa/927/identify-and-assemble-on-boundary-elements

https://fenicsproject.org/qa/927/identify-and-assemble-on-boundary-elements

```
from fenics import *
mesh = UnitSquareMesh(1,1)
P1 = FiniteElement('CG', triangle, 1)
Q = FunctionSpace(mesh, P1)
# create a cell function for each cell :
cf = CellFunction('size_t', mesh, 0)
# for example, the first cell :
one_cell = cells(mesh).next();
cf[one_cell] = 1
# re-define the measure :
dx = Measure('dx', subdomain_data=cf)
# variational form :
u = TestFunction(Q)
v = TrialFunction(Q)
a = inner(grad(u), grad(v))*dx(1)
# assemble over one cell :
A = assemble(a)
print A.array()
```

This gives a global matrix with non-zero values for the component of the local matrix only.

Thanks! That works. But in your example, the global mesh have just one cell.

What if I have mesh with N cells and I want to assemble local matrix on each cell to form global matrix?

Any suggestion?

What if I have mesh with N cells and I want to assemble local matrix on each cell to form global matrix?

Any suggestion?

written
13 months ago by
AJ

There are two cells in this example, one upper and one lower triangle partition of the square.

I am uncertain of the best way to extract the individual local tensors; however, you can construct the global tensor by addition over all single-element assembled tensors, if that helps :

Thanks : https://fenicsproject.org/qa/3150/how-to-perform-a-matrix-multiply

I am uncertain of the best way to extract the individual local tensors; however, you can construct the global tensor by addition over all single-element assembled tensors, if that helps :

```
from fenics import *
import numpy as np
mesh = UnitSquareMesh(1,1)
P1 = FiniteElement('CG', triangle, 1)
Q = FunctionSpace(mesh, P1)
# create a cell function for each cell :
cf = CellFunction('size_t', mesh, 0)
# iterate over all cells :
for i,cell in enumerate(cells(mesh)):
cf[cell] = i;
# re-define the measure :
dx = Measure('dx', subdomain_data=cf)
# variational form :
u = TestFunction(Q)
v = TrialFunction(Q)
# local tensor list :
A_a = []
# assemble global contribution of local tensors :
print "LOCAL STIFFNESS:\n"
for i in range(mesh.num_cells()):
A_i = assemble(inner(grad(u), grad(v))*dx(i))
print A_i.array(), "\n"
A_m = as_backend_type(A_i).mat() # add local petsc matrix
A_a.append(A_m)
A = np.sum(np.array(A_a))
# convert back to dolfin matrix :
C = Matrix(PETScMatrix(A))
print "GLOBAL STIFFNESS:\n", C.array()
```

Thanks : https://fenicsproject.org/qa/3150/how-to-perform-a-matrix-multiply

written
13 months ago by
pf4d

Please login to add an answer/comment or follow this question.

I want to see if it is possible to call FEniCS assembly routines for each element and then build the global matrix from those local assembly calls.

Thanks