Convert numpy array to Generic Matrix

DOLFIN and its linear algebra backends use sparse matrices, and
numpy.ndarray is a dense one.

I would suggest you to do the scaling either directly in the form (if
that is possible) or by toying around with getrow/setrow in the generic
matrix interface.

Have in mind that this will most probably break if you run your
application in parallel using mpi.

Johan

On 09/13/2012 05:15 PM, Karel Van Bockstal wrote:
> Question #208481 on DOLFIN changed:
> https://answers.launchpad.net/dolfin/+question/208481
>
> Description changed to:
> Hello,
>
> I have a problem with applying the boundary condition to a matrix. This
> is a simplified code with the same error:
>
> from dolfin import *
> import numpy as np
>
> #Build mesh
>
> nx = 5 #Number of cells in x direction
> ny = 5 #Number of cells in y direction
> nz = 5 #Number of cells in z direction
>
> domain = UnitCube(nx,ny,nz)
>
> numver = (nx + 1)*(ny + 1)*(nz + 1)
> # Define function spaces
> V = VectorFunctionSpace(domain, "Lagrange", 1)
>
> # Dirichlet boundary
> class DirichletBoundary(SubDomain):
> def inside(self, x, on_boundary):
> return on_boundary
>
> # Define Dirichlet boundary conditions at boundaries
> gD = Expression(('1.0','1.0','1.0'))
>
> bc = DirichletBC(V, gD, DirichletBoundary())
>
> # Define test and trial functions
> H = TrialFunction(V)
> v = TestFunction(V)
>
> F = Expression(('0.0','0.0','0.0'))
>
> # Define variational problem
> a_S = inner(nabla_grad(H), nabla_grad(v))*dx
> a_M = inner(H, v)*dx
>
> M = assemble(a_M, mesh = domain) # mass matrix
> S = assemble(a_S, mesh = domain) # stiffness matrix
>
>
> #pre-built a matrix C
>
> C = np.zeros((3*numver,3*numver))
>
> #for instance
> for j in range(0,numver):
> for i in range(numver):
> C[i,j] = 1
> C[numver + i,numver + j] = 2
> C[2*numver + i, 2*numver+j] = 3
>
> A = M.array()*C + S.array()
>
> f = interpolate(F,V)
>
> b = M*f.vector()
>
> bc.apply(A,b)
>
> #solve problem
> H = Function(V)
>
> solve(A, H.vector(),b)
>
> ###################
>
> I get the following (expected) error due to the matrix A who is a numpy matrix (C is a numpy array)
> bc.apply(A,b)
> TypeError: in method 'DirichletBC_apply', argument 2 of type 'dolfin::GenericVector &'
>
> How can I solve this? In short, I want to multiply a pre-built matrix C
> with the assembled matrix M. Afterwards, I want to apply the boundary
> conditions on A (not by hand and then not using linalg.solve(A,b)
> instead of solve). Thank you in advance.
>
> Karel Van Bockstal
>