← Back to team overview

dolfin team mailing list archive

Re: dolfin <-> scipy <-> pyamg

 

On Tue, Feb 24, 2009 at 10:27:51PM +0100, kent-and@xxxxxxxxx wrote:
> 
> Did you know that dolfin can be used together with scipy and pyamg to solve
> the eg. the Poisson equation with AMG in Python. The attached example by
> Luke Olson does this! This is so cool!

Nice.

> Should we add a webpage at  www.fenics.org where we list up software
> packages that
> can be used with FEniCS software (and nice examples where possible).
> There are quite a few, eg.  VMTK and Netgen can create meshes.
> Itk/Debiosee is getting some support at the moment. Scipy and PyAMG. Of
> course all the linear algebra backends and  mesh formats etc.

Yes, that would be good. See if you can find a good location for the
page on the wiki.

-- 
Anders


> Kent

> # Luke Olson
> # lukeo@xxxxxxxx
> # http://www.pyamg.com
> 
> from dolfin import *
> 
> from scipy.sparse import csr_matrix
> from pyamg import smoothed_aggregation_solver
> from pylab import figure, show, semilogy
> 
> ############################################################
> # Dolfin
> dolfin_set("linear algebra backend","uBLAS")
> 
> mesh = UnitSquare(75, 75)
> V = FunctionSpace(mesh, "CG", 1)
> 
> # Define Dirichlet boundary (x = 0 or x = 1)
> class DirichletBoundary(SubDomain):
>     def inside(self, x, on_boundary):
>         return x[0] < DOLFIN_EPS or x[0] > 1.0 - DOLFIN_EPS
> u0 = Constant(mesh, 0.0)
> bc = DirichletBC(V, u0, DirichletBoundary())
> 
> v = TestFunction(V)
> u = TrialFunction(V)
> f = Function(V, "500.0 * exp(-(pow(x[0] - 0.5, 2) + pow(x[1] - 0.5, 2)) / 0.02)")
> a = dot(grad(v), grad(u))*dx
> L = v*f*dx
> 
> A, rhs = assemble_system(a,L,bc,mesh)
> 
> ############################################################
> # PyAMG
> (row,col,data) = A.data()   # get sparse data
> n = A.size(0)
> Asp = csr_matrix( (data,col.view('intc'),row.view('intc')), shape=(n,n))
> b = rhs.data()
> 
> residuals = []
> 
> ml = smoothed_aggregation_solver(Asp,max_coarse=10)
> x = ml.solve(b,tol=1e-10,accel='cg',residuals=residuals)
> 
> residuals = residuals/residuals[0]
> print ml
> 
> figure(2)
> semilogy(residuals)
> show()
> _______________________________________________
> DOLFIN-dev mailing list
> DOLFIN-dev@xxxxxxxxxx
> http://www.fenics.org/mailman/listinfo/dolfin-dev

Attachment: signature.asc
Description: Digital signature


Follow ups

References