August 2016 - firedrake - imperial.ac.uk

Re: [firedrake] Mixed system
by David Ham 12 Sep '16

12 Sep '16

Unfortunately, It isn't as simple as saying "no mpi". MPI pervades Firedrake. Unless real function spaces play nicely with communicators, Firedrake will simply crash. I'm sorry that this feature is not available right now, but it actually does require significant work. On Wed, 10 Aug 2016 at 16:51 Anna Kalogirou <A.Kalogirou(a)leeds.ac.uk> wrote: > Yes indeed, but I believe David would want to make it work for all cases. > The merge is necessary in order to proceed with my simulations, so an > update as soon as possible would be greatly appreciated. > > Best, Anna. > > > On 10 Aug 2016, at 16:29, Onno Bokhove <O.Bokhove(a)leeds.ac.uk> wrote: > > We would already be helped for the case without MPI, I think that is > correct, Anna? > > > *From:* firedrake-bounces(a)imperial.ac.uk <firedrake-bounces(a)imperial.ac.uk> > on behalf of David Ham <David.Ham(a)imperial.ac.uk> > > *Sent:* Wednesday, August 10, 2016 10:01 AM > > > *To:* firedrake > *Subject:* Re: [firedrake] Mixed system > > The short version is it depends on when I find some time to work on it. It > is not simply the case that we have new features lying around waiting to be > merged. The current proximate issue is to make the real spaces work > properly with MPI communicators. I can't guarantee that will be the last > issue. > > > On Tue, 9 Aug 2016 at 17:37 Onno Bokhove <O.Bokhove(a)leeds.ac.uk> wrote: > >> It is possible to merge the branch for the mixed system relatively soon? >> Anna's results are relevant for the grant application Colin and I aim to >> submit in Sept. on further wave-ship modelling. >> Sorry for asking! >> ------------------------------ >> *From:* firedrake-bounces(a)imperial.ac.uk < >> firedrake-bounces(a)imperial.ac.uk> on behalf of David Ham < >> David.Ham(a)imperial.ac.uk> >> *Sent:* Friday, August 5, 2016 11:37:23 AM >> *To:* firedrake >> *Subject:* Re: [firedrake] Mixed system >> >> >> >> On Fri, 5 Aug 2016 at 11:31 Anna Kalogirou <a.kalogirou(a)leeds.ac.uk> >> wrote: >> >>> Hi Kramer, >>> >>> Thank you for the reply. >>> >>> The only problem I see with the solution you suggest, is that both I and >>> lambda will be trial functions, so both terms should be on the bilinear >>> part of the variational form. >>> >> >> Not quite. I and lambda together form a Function in a mixed finite >> element space. You solve for both of them together, just like pressure and >> velocity in a fluid simulation. >> >> >>> >>> I could, of course, define a residual and solve everything using >>> nonlinear solvers, but I don't know how efficient that would be. >>> >> >> I think you are confused. None of the relevant issues change in the >> nonlinear case. After all, the nonlinear solver will differentiate the >> residual to create a bilinear form anyway. >> >> Incidentally, there is no particular reason to believe that writing the >> problem in nonlinear form would be inefficient. >> >> Regards, >> >> David >> >> >>> >>> Best, >>> >>> Anna. >>> >>> Dr Anna Kalogirou >>> Research Fellow >>> School of Mathematics >>> University of Leeds >>> >>> http://www1.maths.leeds.ac.uk/~matak/ >>> >>> On 04/08/16 17:29, Stephan Kramer wrote: >>> > On 27/07/16 11:01, Anna Kalogirou wrote: >>> >> Dear all, >>> >> >>> >> I am trying to solve the system found in the attached pdf as a mixed >>> >> system. I want to solve the same >>> >> problem as the one found here >>> >> < >>> https://bitbucket.org/annakalog/buoy2d/src/14334d3c20b9f10ed7c1246cde9e3cb6… >>> >, >>> >> >>> >> but with the use of Schur complements and not linear algebra. >>> >> >>> >> In a previous discussion about this problem, Lawrence mentioned that >>> the >>> >> test function for integral(lambda*Theta(x-Lp)dx) needs to be >>> considered >>> >> as coming from the real space of constant functions. Could you please >>> >> elaborate on this? >>> > >>> > My guess at what Lawrence means is (it's how I would probably approach >>> > it), is that you add an additional equation that says >>> > >>> > (1) I = integral(lambda*Theta(x-Lp)dx) >>> > >>> > and then substitute I in your third equation. Equation (1) can >>> > be turned into a finite element weak formulation by considering the >>> > function space of functions that are constant over the entire domain >>> > This function space is just 1-dimensional and therefore equivalent to >>> > the space of reals R and hence we denote this function space of >>> > constants simply as R. Then we can consider I to be a trial function >>> > in R, and with a test function v4 in R, we can write: >>> > >>> > v4*lambda*Theta(x-Lp)*dx == v4*I*dx/area(domain) >>> > >>> > This is a proper finite element equation that you could solve in >>> > conjuction with the other equations. I believe however that this >>> > function space R is not currently implemented in Firedrake. It is >>> > available in fenics, and is on the wishlist. >>> > >>> > Cheers >>> > Stephan >>> > >>> > _______________________________________________ >>> > firedrake mailing list >>> > firedrake(a)imperial.ac.uk >>> > https://mailman.ic.ac.uk/mailman/listinfo/firedrake >>> >>> >>> _______________________________________________ >>> firedrake mailing list >>> firedrake(a)imperial.ac.uk >>> https://mailman.ic.ac.uk/mailman/listinfo/firedrake >>> >> _______________________________________________ > firedrake mailing list > firedrake(a)imperial.ac.uk > https://mailman.ic.ac.uk/mailman/listinfo/firedrake > >

2 3

flood drake query
by Onno Bokhove 08 Sep '16

08 Sep '16

Hi, Using/adapting flood drake: - has there been an update with a gmsh example in which solid walls and inflow/outflow/specified flow conditions have been imposed at the boundaries - has the time step been made adaptable based on flood/signal speeds. - I tried daybreak over uniform beach and lowered the time step (adaptable one preferred) but am not sure whether the slope limiter causes some spurious oscillations. Does the Ern slope limiter approach add an extra time step restriction above and beyond the normal CFL* min_dx/|lambda| with CFL<1 and lambda = |u|+ sort (g h)? - Also the speed seems a bit low although I aim to run it against a 2D c-code I once made. Thank you, Onno See code snippet below. """ demo file for simple 2d shallow water equations flooding on a slope """ from __future__ import division from firedrake import * from flooddrake import * # Meshsize n = 5 mesh = UnitSquareMesh(10*n, n) # mixed function space v_h = FunctionSpace(mesh, "DG", 1) v_mu = FunctionSpace(mesh, "DG", 1) v_mv = FunctionSpace(mesh, "DG", 1) V = v_h*v_mu*v_mv # for slope limiter v_hcg = FunctionSpace(mesh, "CG", 1) v_mucg = FunctionSpace(mesh, "CG", 1) v_mvcg = FunctionSpace(mesh, "CG", 1) VCG = v_hcg*v_mucg*v_mvcg # setup free surface depth g = Function(V) x = SpatialCoordinate(V.mesh()) g.sub(0).interpolate(conditional(x[0] < 0.2, 3.0/9.8, 0.0)) # setup bed bed = Function(V) # setup actual depth w = g.assign(g - bed) x = SpatialCoordinate(V.mesh()) bed.sub(0).interpolate(0.25*x[0]) # setup source (is only a depth function) source = Function(v_h) # timestep solution = Timestepper(V, VCG, bed, source, 0.0065) solution.stepper(0, 2, w, 0.025)

3 8

Re: [firedrake] Saving txt files while running in parallel
by Andrew McRae 30 Aug '16

30 Aug '16

Yes, I guess that in Floriane's original code, each process is computing its own max_eta (since eta.dat.data contains the values for degrees of freedom on that MPI process). Then all process are trying to write (or append?) to the same file(?). You could get the separate max values by making the filename depend on op2.MPI.comm.rank, so that each process writes to a separate file. Otherwise, computing the max over all subdomains indeed requires a parallel operation. On 8 June 2016 at 09:53, Shipton, Jemma <j.shipton(a)imperial.ac.uk> wrote: > Hi Floriane, > > > We have defined a max method that works in parallel for our code... try > something like: > > > def max(f): fmax = op2.Global(1, [-1000], dtype=float) op2.par_loop(op2.Kernel("""void maxify(double *a, double *b) { a[0] = a[0] < fabs(b[0]) ? fabs(b[0]) : a[0]; }""", "maxify"), f.dof_dset.set, fmax(op2.MAX), f.dat(op2.READ)) return fmax.data[0] > > > But I agree, it would be nice if that weren't necessary. > > > Hope that helps, > > > Jemma > > ------------------------------ > *From:* firedrake-bounces(a)imperial.ac.uk <firedrake-bounces(a)imperial.ac.uk> > on behalf of Floriane Gidel [RPG] <mmfg(a)leeds.ac.uk> > *Sent:* 08 June 2016 09:33:46 > *To:* firedrake > *Subject:* [firedrake] Saving txt files while running in parallel > > > Dear all, > > > I am running a Firedrake code in parallel with 4 cores, in which I save > the maximal value of the amplitude of the wave at each time, with the > command: > > > max_eta = max(eta.dat.data) > > Eta_file.write('%-10s %-10s %-10s\n' % (t,'', eta_max)) > > > > But when opening the .txt file, I notice two issues: > > - the data are not saved continuously: time goes from 0 to 100 (let's say) > and then starts from 90 again; > > - the values of max_eta are not the maximal values of eta, except at the > beginning. So it looks like the function max takes the maximum value of eta > only in one subdomain, so that after the wave has crossed the subdomain, > the maximal value of eta goes back to the depth of rest value. > > > How can I force this command to be applied on the full domain, even if the > code is run in parallel? > > > Thanks, > > Floriane > > > >

5 18

Accessing values of a block matrix system
by William Booker 25 Aug '16

25 Aug '16

Dear Firedrakers, I'm trying to compare the matrices I get from matlab to the ones generated in Firedrake, but I've run into some issues: A = assemble(a, nest=True) print A.getValues(1,1) gives me File "compressible_stratified.py", line 137, in <module> print A.getValues(1,1) AttributeError: 'Matrix' object has no attribute 'getValues' How would I fix this ? Regards Will

3 3

Accessing values of a block matrix system
by William Booker 24 Aug '16

24 Aug '16

Dear Firedrakers, I'm trying to compare the matrices I get from matlab to the ones generated in Firedrake, but I've run into some issues: A = assemble(a, nest=True) print A.getValues(1,1) gives me File "compressible_stratified.py", line 137, in <module> print A.getValues(1,1) AttributeError: 'Matrix' object has no attribute 'getValues' How would I fix this ? Regards Will

1 0

firedrake docker image
by moritz braun 22 Aug '16

22 Aug '16

Dear All I am looking for the author of the firedrake docker image, that I am currently using. There is a strange bug in the firedrake installation on it in connection with the evaluation of a finite element function: ................. File "/home/braunm/firedrakeHOMDFEMMol/sphericalFemDomain.py", line 136, in f return waveIP.at(r) File "/usr/local/lib/python2.7/dist-packages/firedrake/function.py", line 780, in at single_eval(points[i:i+1], p_result) File "/usr/local/lib/python2.7/dist-packages/firedrake/function.py", line 758, in single_eval err = self._c_evaluate(self._ctypes, File "/usr/local/lib/python2.7/dist-packages/pyop2/utils.py", line 64, in __get__ obj.__dict__[self.__name__] = result = self.fget(obj) File "/usr/local/lib/python2.7/dist-packages/firedrake/function.py", line 705, in _c_evaluate result = make_c_evaluate(self) File "/usr/local/lib/python2.7/dist-packages/firedrake/function.py", line 842, in make_c_evaluate src = pq_utils.src_locate_cell(function_space.mesh()) File "/usr/local/lib/python2.7/dist-packages/firedrake/pointquery_utils.py", line 32, in src_locate_cell with open(path.join(path.dirname(__file__), "locate.cpp")) as f: IOError: [Errno 2] No such file or directory: '/usr/local/lib/python2.7/dist-packages/firedrake/locate.cpp' ^C Regards Moritz -- Prof M Braun Tel.:27-12-4298006/8027 Physics Department Fax.: 27-12-4293643 University of South Africa (UNISA) moritz.braun(a)gmail.com P.O. Box 392 0003 UNISA South Africa http://moritz-braun.blogspot.com

2 1

Re: [firedrake] preconditioner for time-dependent Stokes
by Colin Cotter 22 Aug '16

22 Aug '16

Hi Justin, No, it's unsteady Stokes which is a different PDE. cheers --cjc On 22 August 2016 at 04:20, Justin Chang <jychang48(a)gmail.com> wrote: > Colin, > > So I didn't read the paper and therefore don't know exactly what you're > trying to implement, but the PETSc folks have some examples on how to > precondition the Stokes equation using command-line options starting around > slide 129: > > http://www.mcs.anl.gov/petsc/documentation/tutorials/BuffaloTutorial.pdf > > Perhaps one of their examples closely resembles what you're looking for? > > Justin > > On Sat, Aug 20, 2016 at 11:12 AM, Colin Cotter < > colin.cotter(a)imperial.ac.uk> wrote: > >> Dear all, >> I'd very much like to try the preconditioner in >> Mardal, Kent-Andre, and Ragnar Winther. "Uniform preconditioners for the >> time dependent Stokes problem." *Numerische Mathematik* 98.2 (2004): >> 305-327. >> >> This is an additive preconditioner that involves summing the result of >> two different preconditioning systems. What is the best path to trying this >> right now? Is it the new matrix-free stuff? >> >> cheers >> --cjc >> >> -- >> http://www.imperial.ac.uk/people/colin.cotter >> >> www.cambridge.org/9781107663916 >> >> >> >> _______________________________________________ >> firedrake mailing list >> firedrake(a)imperial.ac.uk >> https://mailman.ic.ac.uk/mailman/listinfo/firedrake >> >> > -- http://www.imperial.ac.uk/people/colin.cotter www.cambridge.org/9781107663916

1 0

Firedrake not strong-scaling as well as a plain PETSc program?
by Justin Chang 22 Aug '16

22 Aug '16

Hi all, I am not sure if this is a firedrake issue/question per se, but I noticed that when I run PETSc's SNES ex12 (3D FEM poisson with DMPlex) with 531,441 dofs, i get the following speedup for KSPSolve(...) on up to 8 cores: 1: 4.4628 s 2: 2.0991 s 4: 1.0930 s 8: 0.6591 s That's roughly 85 % parallel efficiency. Now, when I run a similar firedrake 3D version: 1: 7.1377 s 2: 3.6969 s 4: 2.0406 s 8: 1.2939 s That's now almost 69% parallel efficiency. I used the same solver and preconditioner (GMRE with ML). The PETSC_DIR I am using for the SNES ex12 is /path/to/firedrake/lib/python2.7/site-packages/petsc so I am basically using the same petsc and compilers for both cases. Using OpenMPI-1.6.5 with the binding options "--bind-to-core --bysocket" on an Intel Xeon E5 2670 node. The efficiency for firedrake gets worse when I scale up the problem and use up to 64 cores (8 cores and 8 nodes). I get as bad as 40 % efficiency when I still maintain roughly 75% for the PETSc case. This is rather strange. Shouldn't I expect roughly the same scaling performance for both implementations? Or is this normal? Note that I disregarded the assembly of the Jacobian or Function for either case because firedrake is much faster that Matt's PetscFE. I understand that the more "serially efficient" a code is, the less parallel efficient it may be. Might this scaling issue have to do with python and/or the implementation of firedrake? FWIW, I installed my own python library because the HPC machine does not have a compatible python 2.7 Thanks, Justin

2 4

preconditioner for time-dependent Stokes
by Colin Cotter 22 Aug '16

22 Aug '16

Dear all, I'd very much like to try the preconditioner in Mardal, Kent-Andre, and Ragnar Winther. "Uniform preconditioners for the time dependent Stokes problem." *Numerische Mathematik* 98.2 (2004): 305-327. This is an additive preconditioner that involves summing the result of two different preconditioning systems. What is the best path to trying this right now? Is it the new matrix-free stuff? cheers --cjc -- http://www.imperial.ac.uk/people/colin.cotter www.cambridge.org/9781107663916

2 1

Call-back functions not working (I think)
by Justin Chang 15 Aug '16

15 Aug '16

Hi all, So I am implementing PETSc/TAO's semi smooth solver to ensure non-negative solutions for anisotropic diffusion. Attached is the code. Basically, the call-back functions require a Jacobian and a gradient/constraint evaluation routine. It begins around line 160. Run the code as: python SS_ex2.py 30 <guess> where <guess> is either use zeros as initial guess (0) or use the standard formulation solution (i.e., no bounds) as the initial guess (1). The implementation somewhat works if <guess> is 1 or an initial non-zero solution is provided. If no guess is provided, I get a funky solution where everywhere is zero except for the non-homogeneous DirichletBCs. Also, even if I have an initial guess (1), the values of the BC aren't exactly accurate either (I should have a max value of 1 at the very center of the top face, but I get 0.983941) Something tells me I am not implementing the BCs correctly. Can anyone tell me what I am possibly doing wrong here? Thanks! Justin

2 3