April 2015 - firedrake - imperial.ac.uk

Re: [firedrake] Mapping nodal values from parent to extruded mesh
by David Ham 04 Jun '15

04 Jun '15

Hi Tuomas, Concerning the sets of spaces which this will work for, I think that the current implementation probably works for any horizontal space but only for vertical CG. One could just try relaxing the test in that case. The reason for the restriction to vertical CG is that fs.bt_masks uses the topological association of nodes with mesh entities in order to work out which nodes are on the top or bottom of the cell. In order to allow DG in the vertical, it would be necessary to support the geometric definition (ie which basis functions do not vanish on the top/bottom). If you look at functionspace.py:78 you can see where the bottom and top masks are generated. This uses entity_closure_dofs() from FIAT. In order to use the geometric definition of dofs one would need to support using facet_support_dofs() (which is how BC maps are set up at functionspace.py:375). Currently FIAT TensorFiniteElement objects do support entity_closure_dofs() but nobody has done the legwork to get them to support facet_support_dofs. Concerning documentation, COFFEE is deplorably underdocumented (although the author is on this list so maybe this will change ;). The actual available AST nodes are only "documented" by reading the source: https://github.com/coneoproject/COFFEE/blob/master/coffee/base.py However, you are fundamentally writing a PyOP2 kernel, and the C api for that *is* documented at: http://op2.github.io/PyOP2/kernels.html so that hopefully helps somewhat. Otherwise, keep asking here or in IRC (IRC tends to get faster responses, at least between about 10am and 8pm UK time). BTW, any chance of seeing you at FEniCS '15? Cheers, David On 10 April 2015 at 02:44, Tuomas Karna <tuomas.karna(a)gmail.com> wrote: > Hi all, > > A couple of questions regarding those extruded mesh -> parent mesh copy > operations, > > A while back I got the reverse 2d->3d copy working with a hard-coded pyop2 > kernel (like in mesh extrusion). I'm not familiar with COFFEE syntax, is > there documentation/examples somewhere? > > The 3d->2d copy routine in extrusion_extraction branch checks that the > function space is CG. Does this method easily generalize to other spaces? > DGxCG prisms seem to work OK. I'm also interested in using DGxDG and > RTxCG/DG spaces. > > > Cheers, > > Tuomas > > > On 03/05/2015 10:24 AM, Tuomas Karna wrote: > > Thanks David, > > This is great, I'll try doing the reverse operation. > > - Tuomas > > On 03/05/2015 08:20 AM, David Ham wrote: > > Hi Tuomas, > > This wasn't there this morning but I've implemented one of the cases > (pulling out the top and bottom maps). The result is in the > extrusion_extraction branches of both PyOP2 and Firedrake. The 2d->3d > operation would be rather similar except that you'd have to interrogate the > fiat_element on the extruded Function in order to determine which extruded > nodes are "on top of" the 2d nodes you have. Unfortunately I won't have > time to do that one soon (huge amounts to do before SIAM CSE next week) but > feel free to have a try and complain when it doesn't work! > > Cheers, > > David > > On 4 March 2015 at 22:23, Tuomas Karna <tuomas.karna(a)gmail.com> wrote: > >> Hi all, >> >> I'd need to copy nodal values between fields on parent and extruded >> meshes. For example, copy 2d->3d (constant over vertical) or 3d->2d >> (extract surface/bottom level). The horizontal function space is the >> same. Is there an easy way to do this? Is there a map of extruded nodes >> somewhere? >> >> Thanks, >> >> Tuomas >> >> >> >> _______________________________________________ >> firedrake mailing list >> firedrake(a)imperial.ac.uk >> https://mailman.ic.ac.uk/mailman/listinfo/firedrake >> > > > > -- > Dr David Ham > Departments of Mathematics and Computing > Imperial College London > > http://www.imperial.ac.uk/people/david.ham > > > _______________________________________________ > firedrake mailing listfiredrake@imperial.ac.ukhttps://mailman.ic.ac.uk/mailman/listinfo/firedrake > > > > -- Dr David Ham Departments of Mathematics and Computing Imperial College London http://www.imperial.ac.uk/people/david.ham

4 6

unexpected keyword argument 'search'
by Colin Cotter 30 Apr '15

30 Apr '15

Dear Firedrake list, Upon executing the following: from firedrake import * n = 10 mesh = UnitSquareMesh(n,n) V1 = FunctionSpace(mesh,"BDM",2) V2 = FunctionSpace(mesh,"DG",1) W = MixedFunctionSpace((V1,V2)) F = Function(V1).project(Expression(("sin(x[0]*2*pi)*cos(x[1]*2*pi)", "-cos(x[0]*2*pi)*sin(x[1]*2*pi)"))) I get: File "/home/cjc1/firedrake/firedrake/firedrake/function.py", line 157, in project return projection.project(b, self, *args, **kwargs) File "/home/cjc1/firedrake/firedrake/firedrake/projection.py", line 69, in project f.interpolate(v) File "/home/cjc1/firedrake/firedrake/firedrake/function.py", line 221, in interpolate subset) File "/home/cjc1/firedrake/firedrake/firedrake/function.py", line 275, in _interpolate op2.par_loop(*args) File "<string>", line 2, in par_loop File "/home/cjc1/firedrake/PyOP2/pyop2/versioning.py", line 150, in modifies_arguments retval = func(*args, **kwargs) File "/home/cjc1/firedrake/PyOP2/pyop2/op2.py", line 276, in par_loop return backends._BackendSelector._backend.par_loop(kernel, iterset, *args, **kwargs) File "/home/cjc1/firedrake/PyOP2/pyop2/base.py", line 4252, in par_loop return _make_object('ParLoop', kernel, it_space, *args, **kwargs).enqueue() File "/home/cjc1/firedrake/PyOP2/pyop2/backends.py", line 72, in _make_object return _BackendSelector._backend.__dict__[obj](*args, **kwargs) File "<string>", line 2, in __init__ File "/home/cjc1/firedrake/PyOP2/pyop2/utils.py", line 122, in wrapper return f(*args, **kwargs) File "/home/cjc1/firedrake/PyOP2/pyop2/base.py", line 3906, in __init__ ast_info = ast_visit(self._kernel._ast, search=ast.FunDecl) TypeError: visit() got an unexpected keyword argument 'search' If I've done everything correctly, this is with everything on master. What do you think might be wrong? cheers --cjc

2 2

Last day for abstracts!
by David Ham 30 Apr '15

30 Apr '15

Dear all, The flood of abstracts over the last 24 hours says a lot about the power of deadlines. This is the last reminder that everyone who wishes to present a poster or a talk at FEniCS '15 needs to submit an abstract today. Abstracts should be submitted on EasyChair: https://easychair.org/conferences/?conf=fenics15 I look forward to the rest of the last minute rush. David

1 0

Re: [firedrake] Nullspace & VectorSpaceBasis
by Andrew McRae 29 Apr '15

29 Apr '15

Four comments: 1) If I understand correctly, PETSc insists that the diagonal entries of the sparse matrix are filled in. Usually, the LHS form will cause this to happen; however, that obviously isn't the case for your form! 2) It should be possible to build the nullspace you want, but I'm clueless on the actual syntax. This wouldn't get around (1) without further fixes from our end. 3) Could the equations given not be included in the weak form for the rest of the problem? 4) If not, there might be a hacky "add-a-weighted-mass-matrix-to-the-LHS" fix you could use?.... Andrew On 11 December 2014 at 12:46, Anna Kalogirou <a.kalogirou(a)leeds.ac.uk> wrote: > Dear all, > > I am trying to solve a simple problem (linear water wave equations), > where I have the Laplace equation in a square domain with Neumann > boundary conditions on the three sides and on the top boundary two > equations are satisfied. > > I have the following code: > -------------------------------------- > V = FunctionSpace(mesh,"CG",1) > > eta0 = Function(V).interpolate(Expression("cos(2*pi*x[0])")) > phi0 = Function(V) > eta1 = Function(V) > phi1 = Function(V) > eta = TrialFunction(V) > phi = TrialFunction(V) > v = TestFunction(V) > > aeta = v*eta*ds(4) > Leta = v*eta0*ds(4) + dt*inner(grad(v),grad(phi0))*dx > > eta_problem = LinearVariationalProblem(aeta,Leta,eta1) > eta_solver = LinearVariationalSolver(eta_problem) > > aphi = v*phi*ds(4) > Lphi = (v*phi0 - g*dt*v*eta1)*ds(4) > > phi_problem = LinearVariationalProblem(aphi,Lphi,phi1) > phi_solver = LinearVariationalSolver(phi_problem) > -------------------------------------- > > The error I get tells me that Object is in wrong state, and Matrix is > missing diagonal entry 39. > > I thought this had to do with the problem being singular, so I tried to > build the nullspace but I don't know how to define the VectorSpaceBasis. > Any help would be appreciated. > > Thank you very much. > > Regards, Anna. > > _______________________________________________ > firedrake mailing list > firedrake(a)imperial.ac.uk > https://mailman.ic.ac.uk/mailman/listinfo/firedrake >

3 6

Re: [firedrake] leaking memory
by David Ham 29 Apr '15

29 Apr '15

I agree about the hard limit. There is a limit to how much useful reordering we are ever going to do on the trace. Shall we set a default value for lazy_max_trace_length? Is 100 a good number? On Wed, 29 Apr 2015 at 10:00 Mitchell, Lawrence < lawrence.mitchell(a)imperial.ac.uk> wrote: > -----BEGIN PGP SIGNED MESSAGE----- > Hash: SHA1 > > On 29/04/15 08:59, Hiroe Yamazaki wrote: > > Dear all, > > > > I am experiencing a memory leak in our vertical slice model and > > wonder if there is any possibility of leaking memory associated > > with .assign(). Here I have a stand-alone code with the same > > function spaces and variable settings we use in the slice model: > > > > > https://bitbucket.org/hyamazak/test_files/src/f455e7e68722c9264845a8426ba91… > > > > What I do here is just doing > > > > new.vars += Deltavars old.vars.assign(new.vars) > > > > over and over again, and I see a constant increase of the memory > > usage both on my Linux and Macbook. I believe new.vars += Deltavars > > is not the source of leaking and wonder if there is any problem > > with doing old.vars.assign(new.vars). > > > > Do let me know if you have any ideas, > > OK, I observe the same memory leak. Here's what's going on: > > PyOP2 is running in lazy evaluation mode, which means that no > computation actually occurs until you ask to look at the resulting > data (perhaps by writing it to a file, or dropping into a solve). > Instead, we store an "execution trace" which remembers the sequence of > computations you want to do. By default this execution trace is > unbounded in length, so, since you never look at the answer in the > above code, we end up creating a trace of length 100000. This is the > cause of the memory leak. > > Arguably we should have a limit on the trace length. Can you try > setting the following parameter to see if your leak goes away? > > parameters["pyop2_options"]["lazy_max_trace_length"] = 100 > > > Cheers, > > Lawrence > > -----BEGIN PGP SIGNATURE----- > Version: GnuPG v1 > > iQEcBAEBAgAGBQJVQJ2DAAoJECOc1kQ8PEYvr8gIALFwcyFV3XIcK+S3uzKDYkZ+ > aGR5EprrM1dvn2U4nS9eGTm52a9fxSYaZq3H++JfwdfZ0J/XWCFxfHkwAoCxmWRC > 23m+sohS5om8AjhmW4hU8S2LK6RQknDiKxt14IUE5nZXcf5aIvFF6/l3/pclptvg > OBorwxDk6LqHS17/4UOhIMjkcRLEppZ3TuH/NLqo7vftgnclf77ONjWfFP8y4DIh > 6AMPmcbx1WXLUZ+UIB7llMdaNgsH98Ei7v+hlvQnrdMxFAm4AbplAoVW6PGeTJ2z > h29W6LqVddrWaE7QMp2N+r732yY/3D7Pmjk0Iw5mO2uyYkCFrBVUUxUI+e6tqs4= > =STEw > -----END PGP SIGNATURE----- > > _______________________________________________ > firedrake mailing list > firedrake(a)imperial.ac.uk > https://mailman.ic.ac.uk/mailman/listinfo/firedrake >

2 1

leaking memory
by Hiroe Yamazaki 29 Apr '15

29 Apr '15

Dear all, I am experiencing a memory leak in our vertical slice model and wonder if there is any possibility of leaking memory associated with .assign(). Here I have a stand-alone code with the same function spaces and variable settings we use in the slice model: https://bitbucket.org/hyamazak/test_files/src/f455e7e68722c9264845a8426ba91… What I do here is just doing new.vars += Deltavars old.vars.assign(new.vars) over and over again, and I see a constant increase of the memory usage both on my Linux and Macbook. I believe new.vars += Deltavars is not the source of leaking and wonder if there is any problem with doing old.vars.assign(new.vars). Do let me know if you have any ideas, all the best, Hiroe

2 1

FEniCS '15 abstracts due THURSDAY, keynote and public lectures.
by David Ham 28 Apr '15

28 Apr '15

Hello all, Thanks to all those who have submitted abstracts already. Those who are yet to submit have until Thursday to do so. I am also happy to be able to announce two of the more prominent presentations at this year's meeting. I have decided to continue the tradition of having a keynote from the wider finite element community. This year's keynote will be by Dr Jennifer Ryan of the University of East Anglia. She will talk on superconvergent filters for DG schemes. We will also hold a lecture open to people at Imperial and beyond who may be interested in the FEniCS project but who are not attending the whole workshop. The public lecture will be given by Prof Anders Logg. I look forward to the continuing flood of abstracts. Regards, David

1 0

Firedrake meeting tomorrow at the earlier time of 1500
by David Ham 27 Apr '15

27 Apr '15

Hi All, Due to Paul and my travel arrangements tomorrow evening, it would be good if we could pull tomorrow's Firedrake meeting forward to 1500. The meeting will be in the William Penney fish bowl.

1 0

CUDA support broken after update
by Hector Dearman 27 Apr '15

27 Apr '15

Hi all, I recently updated Firedrake, PyOP2, etc and now any attempt to use the CUDA backend fails with the following stack trace (the sequential backend still works): ...elided... File "/homes/hjd11/proj/Firedrake/firedrake/firedrake/functionspace.py", line 585, in __init__ super(VectorFunctionSpace, self).__init__(mesh, element, name, dim=dim, rank=1) File "/homes/hjd11/proj/Firedrake/firedrake/firedrake/functionspace.py", line 121, in __init__ with function.Function(self).dat.vec_ro as v: AttributeError: 'Dat' object has no attribute 'vec_ro' The problem seems to be this new code in firedrake/functionspace.py introduced in commit 109e128133d6622d99d80fe99144cf269e1039bf: + # Tell the DM about the layout of the global vector + with function.Function(self).dat.vec_ro as v: + self._dm.setGlobalVector(v.duplicate()) Which depends on Dat having a vec_ro which is defined on the PETSc derived Dat class but not on the CUDA one. Is the CUDA backend broken for everyone, or did I do something wrong when I updated Firedrake? Thanks! Hector Commit: https://github.com/firedrakeproject/firedrake/commit/109e128133d6622d99d80f… Code in Firedrake master: https://github.com/firedrakeproject/firedrake/blob/master/firedrake/functio… Minimal Example & Traceback: https://gist.github.com/chromy/9117a372d1b2c3b9d27d

2 1

difference between matrix-action and form-action
by Eike Mueller 27 Apr '15

27 Apr '15

Dear firedrakers, Assume I have assembled a UFL form in the HDiv space: form_uu = [a UFL form] mat_uu = assemble(form_uu,bcs=bcs) then shouldn’t the following two experssions fo exactly the same thing? version A: with v.dat.vec as y: with u.dat.vec_ro as x: mat_uu.M.handle.mult(x,y) version B: v = self._mutilde._form_action(u) as fas as I can tell, only version B gives the correct result in my case. Thanks, Eike

2 4