September 2015 - firedrake

Multi-grid preconditioning
by Justin Chang 16 Nov '15

16 Nov '15

Hi everyone, Is there support for Multi-grid preconditioning? I think I saw somewhere where this was not supported yet. Because applying -pc_type gamg for my mixed forms returns me errors. Thanks, Justin

2 3

Installing Firedrake on an HPC machine
by Justin Chang 14 Nov '15

14 Nov '15

Hi everyone, I have installed firedrake on my university's HPC machine, and whenever i attempt to run any Firedrake program, I get this error: -------------------------------------------------------------------------- An MPI process has executed an operation involving a call to the "fork()" system call to create a child process. Open MPI is currently operating in a condition that could result in memory corruption or other system errors; your MPI job may hang, crash, or produce silent data corruption. The use of fork() (or system() or other calls that create child processes) is strongly discouraged. The process that invoked fork was: Local host: compute-0-0 (PID 28214) MPI_COMM_WORLD rank: 0 If you are *absolutely sure* that your application will successfully and correctly survive a call to fork(), you may disable this warning by setting the mpi_warn_on_fork MCA parameter to 0. -------------------------------------------------------------------------- Why is this happening? Program seems fine when I only use a few processes but when I use two or more compute nodes, my program hangs. I'll do my best to explain how I built firedrake: 1) Our system has openmpi-1.8.3 and gcc-4.9.2 available. so those are the only modules I have loaded, everything else was downloaded. 2) I built Python 2.7.10 from source and installed it in ${HOME}/.local 3) pip was obtained via python get-pip.py 4) PETSc was installed --with-cc=mpmicc, --with-cxx=mpicxx, and has all the necessary external packages 5) petsc4py was installed via git clone ... and set up with python setup.py build/install 6) PyOP2 was built with python setup.py build_ext --inplace 7) Firedrake was built via make and everything else was installed via 'pip install <package>'. I had to manually install swig and pcre. Any ideas what this may be? And/or how I can circumvent this issue? Thanks, Justin

2 20

Dirichlet BC on extruded mesh with discontinuous function space
by Buesing, Henrik 16 Oct '15

16 Oct '15

Dear Firedrakers, I'm constructing a DG space on an extruded mesh. meshbase = RectangleMesh(Nx, Ny, Lx, Ly, quadrilateral=True) mesh = ExtrudedMesh(meshbase, Nz, Delta_z) horiz_elt = FiniteElement("DG", quadrilateral, 0) vert_elt = FiniteElement("DG", interval, 0) elt = OuterProductElement(horiz_elt, vert_elt) DG = FunctionSpace(mesh, elt) W=DG*DG Now I would like to impose Dirchlet BC on one of the boundaries. I tried bc0 = DirichletBC(W.sub(0), Expression("1.0e7"), 2,method="geometric") and get a ValueError: Geometric boundary conditions are not yet supported on extruded meshes Using the standard topological way, seems to have no effect at all. Is there a way around this? I would be fine with fixing the value on the most outer cell layer, if that is possible... Henrik -- Dipl.-Math. Henrik Büsing Institute for Applied Geophysics and Geothermal Energy E.ON Energy Research Center RWTH Aachen University ------------------------------------------------------ Mathieustr. 10 | Tel +49 (0)241 80 49907 52074 Aachen, Germany | Fax +49 (0)241 80 49889 ------------------------------------------------------ http://www.eonerc.rwth-aachen.de/GGE hbuesing(a)eonerc.rwth-aachen.de ------------------------------------------------------

3 4

Re: [firedrake] Dirichlet BC on extruded mesh with discontinuous function space
by David Ham 09 Oct '15

09 Oct '15

I believe this functionality exists in Tuomas' branch. I'll be in Oregon next week and merging this is on my to do list. On Wed, 23 Sep 2015 at 11:34 Buesing, Henrik <HBuesing(a)eonerc.rwth-aachen.de> wrote: > > > *Von:* firedrake-bounces(a)imperial.ac.uk [mailto: > firedrake-bounces(a)imperial.ac.uk] *Im Auftrag von *Andrew McRae > *Gesendet:* 23 September 2015 12:16 > > > *An:* firedrake(a)imperial.ac.uk > *Betreff:* Re: [firedrake] Dirichlet BC on extruded mesh with > discontinuous function space > > > > Yes, but it should still be possible to set the value of the function in > those cells. Didn't Tuomas make some progress on this functionality? I > don't remember what the current status is. > > *[Buesing, Henrik] That would be nice. Then I can just prescribe the right > value at those cells. If this is enforced during the KSP, SNES solves > everything should be fine. * > > As a workaround, Miklos's suggestion may suffice. > > On 23 September 2015 at 11:09, Miklos Homolya <m.homolya14(a)imperial.ac.uk> > wrote: > > Dear Henrik, > > > > You’ve got a DG0 space, i.e. all your nodes are in the cell centres. > Neither topologically nor geometrically have you any nodes on the boundary > faces. > > > > Maybe you want to impose your BC weakly? > > > > Regards, > > Miklos > > > > On 23 Sep 2015, at 12:05, Buesing, Henrik <HBuesing(a)eonerc.rwth-aachen.de> > wrote: > > > > Dear Firedrakers, > > > > I’m constructing a DG space on an extruded mesh. > > > > > > > > *meshbase = RectangleMesh(Nx, Ny, Lx, Ly, quadrilateral=True) mesh = > ExtrudedMesh(meshbase, Nz, Delta_z) horiz_elt = FiniteElement("DG", > quadrilateral, 0) vert_elt = FiniteElement("DG", interval, 0) elt = > OuterProductElement(horiz_elt, vert_elt)* > > *DG = FunctionSpace(mesh, elt)* > > *W=DG*DG* > > > > Now I would like to impose Dirchlet BC on one of the boundaries. > > I tried > > *bc0 = DirichletBC(W.sub(0), Expression("1.0e7"), 2,method="geometric")* > > > > and get a > > > > *ValueError: Geometric boundary conditions are not yet supported on > extruded meshes* > > > > Using the standard topological way, seems to have no effect at all. > > > > Is there a way around this? I would be fine with fixing the value on the > most outer cell layer, if that is possible… > > > > Henrik > > > > > > > > > > > > -- > > Dipl.-Math. Henrik Büsing > > Institute for Applied Geophysics and Geothermal Energy > > E.ON Energy Research Center > > RWTH Aachen University > > ------------------------------------------------------ > > Mathieustr. 10 | Tel +49 (0)241 80 49907 > > 52074 Aachen, Germany | Fax +49 (0)241 80 49889 > > ------------------------------------------------------ > > http://www.eonerc.rwth-aachen.de/GGE > > hbuesing(a)eonerc.rwth-aachen.de > > ------------------------------------------------------ > > > > _______________________________________________ > firedrake mailing list > firedrake(a)imperial.ac.uk > https://mailman.ic.ac.uk/mailman/listinfo/firedrake > > >

2 2

Re: [firedrake] Solving this system of nonlinear equations
by David Ham 02 Oct '15

02 Oct '15

I'm not totally sure that's obvious. I agree that it's not the theoretically optimal solution, however it is well engineered which may make it a convenient solution, and if it's fast enough, then that's just fine. On Fri, 4 Sep 2015 at 11:22 Cotter, Colin J <colin.cotter(a)imperial.ac.uk> wrote: > PETSc shouldn't be used for this, as the problem is solved independently > at each node. > > --cjc > > On 4 September 2015 at 10:40, Ham, David A <david.ham(a)imperial.ac.uk> > wrote: > >> OK. I think I understand. In essence, assembling PETSc matrices is a core >> PyOP2 competence, so it *can* do this, although it'll be a little clunky >> because we don't have the full glory of Firedrake on top. >> >> I think it would help if you try to read some of the PyOP2 documentation >> for a bit of background: >> >> http://op2.github.io/PyOP2/concepts.html >> http://op2.github.io/PyOP2/kernels.html >> http://op2.github.io/PyOP2/linear_algebra.html >> http://op2.github.io/PyOP2/mixed.html >> >> The last one is important in this case because you're actually solving a >> mixed system over the mixed function space (A, B, C). >> >> Inconveniently, that means we need 3 kernels to assemble the residual and >> 7 kernels (one for each nonzero entry) to assemble the Jacobian. >> >> Once you have the code for assembling, I think it should be >> straightforward to provide PETSc4py with the right callbacks so that its >> variational inequality solver works. >> >> (Another possibility would be to bypass PyOP2 and write your own code >> accessing PETSc4py directly. That might actually be easier but assembly >> will be slower - which you may or may not care about). >> >> I don't have time right this instant to start writing out the Kernels (I >> have to teach on Monday and I'm not ready yet.) I'll try to get back to >> this later, or maybe someone else on the list will dive in and save me ;) >> >> Cheers, >> >> David >> > On Fri, 4 Sep 2015 at 09:37 Justin Chang <jychang48(a)gmail.com> wrote: >> > Hi everyone, >>> >>> Thanks for the input. So I was leaning more towards employing the newton >>> method for solving this. And yes, all the vectors live in the same function >>> space. If the system of equations were solved at each node i, the residual >>> algebraic vector of size 3x1 would be: >>> >>> r = [c_A(i)+c_C(i)-psi_A(i); >>> c_B(i)+c_C(i)-psi_B(i); >>> k1*c_A(i)*c_B(i)-c_C(i)] >>> >>> and the Jacobian 3x3 matrix would be: >>> >>> J = [1, 0, 1; >>> 0, 1, 1; >>> k1*c_B(i), k1*c_A(i), -1] >>> >>> Ideally I would like to have these expressed as a PETSc matrix and >>> vector because I would like to employ PETSc’s variational inequality solver >>> to enforce lower/upper bounds for c_A/B/C. >>> >>> As for the pyop2 option: >>> >>> void solve_system(double *c_A, double *c_B, double *c_C, double *k1, >>>>> double *psi_A, double *psi_B) >>>>> { >>>>> >>>>> // In here, you can dereference these variables to get the components >>>>> // of each variable at the node, e.g. >>>>> >>>>> // Set first component of c_A to sum of first two components of psi_A. >>>>> c_A[0] = psi_A[0] + psi_A[1]; >>>>> } >>>>> >>>>> I guess you can either write the newton (or fixed point) iteration out >>>>> by hand in the kernel (it shouldn't be too difficult), or else make library >>>>> calls (getting the correct link lines is possible I think). >>>>> >>>>> Now in your code you do: >>>>> >>>>> solve_kernel = op2.Kernel(kernel_code, "solve_system") >>>>> >>>>> op2.par_loop(solve_kernel, >>>>> c_A.dof_dset.node_set, >>>>> c_A.dat(op2.WRITE), >>>>> c_B.dat(op2.WRITE), >>>>> c_C.dat(op2.WRITE), >>>>> k1.dat(op2.READ), >>>>> psi_A.dat(op2.READ), >>>>> psi_B.dat(op2.READ)) >>>>> >>>>> >>>>> Does this make sense? >>>>> >>>> >>> I am guessing void solve_system(…) is a separate C file called >>> solve_system.c but how exactly do I get the components of each variable at >>> the node? So here’s how I did this with scipy: >>> >>> ================================ >>> >>> # Advection-diffusion firedrake code >>> # Compute psi_A, psi_B, and k1 ... >>> … >>> >>> from scipy.optimize import solve >>> import math >>> >>> # Solve system of nonlinear reactions >>> def equations(p,psi_A,psi_B,k_1): >>> c_A,c_B,c_C = p >>> return (c_A+c_C-psi_A, c_B+c_C-psi_B, c_C-k_1*c_A*c_B) >>> >>> # define vectors for primary species >>> c_A_vec = np.zeros(len(psi_A.vector())) >>> c_B_vec = np.zeros(len(psi_A.vector())) >>> >>> # define vector for secondary species >>> c_C_vec = np.zeros(len(psi_A.vector())) >>> >>> # Iterate through each node >>> for i in range(len(psi_A.vector())): >>> c_A_i,c_B_i,c_C_i = >>> fsolve(equations,(0.0,0.0,0.0),args=(psi_A.vector()[i],psi_B.vector()[i],k1) >>> c_A_vec[i] = c_A_i >>> c_B_vec[i] = c_B_i >>> c_C_vec[i] = c_C_i >>> >>> # Output >>> c_A = Function(Q) >>> c_B = Function(Q) >>> c_C = Function(Q) >>> c_A.vector()[:] = c_A_vec >>> c_B.vector()[:] = c_B_vec >>> c_C.vector()[:] = c_C_vec >>> File(“Solution_c_A.pvd”) << c_A >>> File(“Solution_c_B.pvd”) << c_B >>> File(“Solution_c_C.pvd”) << c_C >>> >>> =================================== >>> >>> The PyOP2 equivalent of the above would be a good starting point I >>> suppose. >>> >>> Thanks all, >>> Justin >>> >>> > -- > http://www.imperial.ac.uk/people/colin.cotter > > www.cambridge.org/9781107663916 > > >

3 12

Re: [firedrake] The wonderful world of firedrake-install
by David Ham 02 Oct '15

02 Oct '15

Hi Justin, Sorry for the delay. I'm rushing backwards and forwards between the University of Reading and Imperial this week so lots of balls are getting dropped. The hdf5 problem is yet another "pip doesn't really look after dependencies" issue. I think we need to forcibly move the PETSc install earlier in the install process (because it's PETSc which installs the hdf5). I'll give that a try when I'm on a real internet connection. Your problem on your cluster is that the just in time compilation of one of your forms is failing. If you read the error log closely, it tells you which files failed to compile (in this case it's /scratch/306103/pyop2-cache-uid3003/c47dfce75aa72266c548fad137321af9.c) It also tells you where the related files containing the compile log and compile errors are. Could you post those files please? Cheers, David On Mon, 7 Sep 2015 at 00:01 Justin Chang <jychang48(a)gmail.com> wrote: > David, > > The error is still there even with the firedrake-install > from install-cython-workaround. For some strange reason I cannot link > the mip library to Cython using pip install ... but I compile Cython 0.23 > from source (./configure --prefix=/home/jchang23/.local && make && make > install) it links properly. I suspect this may have more to do with how my > local cluster's path dependencies behave. > > But now, i get this other error as shown in the attached > firedrake-install.log. I am guessing this has to do with not having an HDF5 > library available within my custom installed Python? > > Also, I note from inspection that apt-get normally installs the > blas/lapack libraries for you, so I added the tag '--download-fblaslapack' > to the petsc configure options. Is that the right thing to do or should I > somehow install blas/lapack as a whole so that all packages depend on it? > > Finally, I finally got firedrake to install properly when I loaded the > intel MPI compiled python libraries. They happened to have all of the most > up-to-date dependencies, and installed with 'python firedrake-install > --user --no_package_manager --disable_ssh'. When I run my programs on the > head node (and avoid getting yelled at), the RT0 code I have attached runs > fine, but when I submit a batch script as shown in runCompare, i get the > outputs as shown in output_error/test. Can you guys dissect what's going on > or would this be more of an issue that needs to be discussed with the > system admins? > > Thanks, > Justin > > > On Fri, Sep 4, 2015 at 4:14 AM, David Ham <David.Ham(a)imperial.ac.uk> > wrote: > >> Hi Justin, >> >> The prospective fix is in install-cython-workaround . Can you try that >> and tell me if it works for you? It will be merged to master as soon as >> Travis gives the green light. >> >> Cheers, >> >> David >> >> On Fri, 4 Sep 2015 at 09:42 Ham, David A <david.ham(a)imperial.ac.uk> >> wrote: >> >>> Hi Justin, >>> >>> Ah yes. I saw something like this when attempting to set up travis. The >>> problem appears to be that the dependency analysis in pip is broken so it >>> doesn't realise it has to build and install Cython before it builds h5py. I >>> think the workaround is to force Cython installation in the install script >>> before looking at the requirements file. I will push a branch presently and >>> merge it as soon as it passes tests. >>> >>> Cheers, >>> >>> David >>> >>> On Fri, 4 Sep 2015 at 08:17 Justin Chang <jychang48(a)gmail.com> wrote: >>> >>>> Hi David >>>> >>>> Thanks for the updated script. So I am using the system provided >>>> openmpi/1.8.3 module but built my own Python 2.7.10 and pip from source. I >>>> am getting an error as shown in the attached firedrake-install.log >>>> >>>> I also got a similar issue when I built my own openmpi library, >>>> although in that case I didn’t have the cython related errors and only had >>>> the ones related to hdf5/h5py >>>> >>>> Any idea what this is? >>>> >>>> Thanks, >>>> Justin >>>> >>>> >>>> On Sep 3, 2015, at 8:54 AM, David Ham <David.Ham(a)imperial.ac.uk> wrote: >>>> >>>> The branch is merged now so I believe this is now fixed. >>>> >>>> David >>>> >>>> On Thu, 3 Sep 2015 at 12:09 David Ham <David.Ham(a)imperial.ac.uk> wrote: >>>> >>>> It turns out that this is for the usual reason: I am an idiot. A fix >>>> for this is in the further-install-fixes branch and I will merge the branch >>>> as soon as tests pass. >>>> >>>> (What's actually going on is that I am abusing CalledProcessError when >>>> I should really raise a special exception for this case. I have now created >>>> a new exception and raised that instead). >>>> >>>> Cheers, >>>> >>>> David >>>> >>>> On Thu, 3 Sep 2015 at 01:56 Justin Chang <jychang48(a)gmail.com> wrote: >>>> >>>> Thanks David, >>>> >>>> I actually forgot about the --no_package_manager option. But when I >>>> tried installing with it, i get this error: >>>> >>>> $ python firedrake-install --no_package_manager >>>> Traceback (most recent call last): >>>> File "firedrake-install", line 265, in <module> >>>> raise subprocess.CalledProcessError >>>> TypeError: __init__() takes at least 3 arguments (1 given) >>>> >>>> Even if I add the --disable_ssh and --developer tags I get the above >>>> error. Am I suppose to do something else with this? >>>> >>>> Thanks, >>>> Justin >>>> >>>> On Tue, Sep 1, 2015 at 3:57 AM, David Ham <David.Ham(a)imperial.ac.uk> >>>> wrote: >>>> >>>> Thanks, >>>> >>>> On Sun, 30 Aug 2015 at 08:15 Justin Chang <jychang48(a)gmail.com> wrote: >>>> >>>> I agree that hype and a sparse direct solver like mumps is necessary. >>>> The latter needs --download-metis --download-parmetis >>>> --download-scalapack. >>>> >>>> >>>> Added and currently testing. >>>> >>>> >>>> More comments: >>>> >>>> 4) If mumps is to be included, the install script has to ensure that >>>> one has cmake >=2.5 when configuring metis, otherwise PETSc configure >>>> will return an error. Perhaps include "brew install cmake" and "sudo >>>> apt-get install cmake" in the script? >>>> >>>> >>>> Also added. Thanks for catching this. >>>> >>>> >>>> 5) Should --download-exodusii (and it's dependencies --download-netcdf >>>> --download-hdf5) also be included in the PETSC_CONFIGURE_OPTIONS? I >>>> imagine that if one wants to solve problems on large-scale >>>> unstructured grids, .exo files would be more feasible to use over .msh >>>> files. Or is exodusii already accounted for in firedrake/PyOP2/etc? >>>> >>>> >>>> And this. >>>> >>>> >>>> 6) If I run the script on a pristine MacOSX, I get an error saying >>>> "OSError: Unable to find SWIG installation. Please install SWIG >>>> version 2.0.0 or higher." I am guessing "brew install swig" was >>>> somehow missed in the firedrake-install script? >>>> >>>> Ah yes. Fecking swig. I've added the dependency to the new branch. >>>> Hopefully the swig dependency is going away very soon. >>>> >>>> >>>> >>>> 7) Perhaps not as important, but do y'all think it's possible to make >>>> this script more non-Ubuntu HPC system friendly? For isntance, either >>>> a) require the user to install his or her own local OpenMPI and Python >>>> libraries from source, b) let script download the tarballs/git >>>> repositories remotely and install it for you (kind of like how PETSc >>>> handles external packages through --download-<package>), or c) enable >>>> the user to point to the system provided OpenMPI and Python libraries. >>>> Same would need to be done for CMake, SWIG, and PCRE as these packages >>>> cannot simply be obtained from pip. >>>> >>>> >>>> I think (c) is the case right now. If you run the script on a >>>> non-Ubuntu system or pass the --no_package_manager option then it will just >>>> assume that you have installed the dependencies somewhere visible and will >>>> get on with pip. >>>> >>>> I think (b) is [a] hard and [b] in the case of MPI probably a bad idea. >>>> Supercomputers are such a diverse bunch of hacked up distributions that I >>>> think that installing the compiled dependencies automatically in a way >>>> which is portable would require a lot of install script work. In the case >>>> of MPI, you really need to use the MPI build which talks to the fast >>>> interconnect on your supercomputer, so automatically building vanilla >>>> OpenMPI is a Bad Thing. >>>> >>>> Cheers, >>>> >>>> David >>>> >>>> >>>> >>>> >>>> Thanks, >>>> Justin >>>> >>>> On Fri, Aug 28, 2015 at 7:20 AM, Lawrence Mitchell >>>> <lawrence.mitchell(a)imperial.ac.uk> wrote: >>>> > >>>> >> On 28 Aug 2015, at 14:01, David Ham <David.Ham(a)imperial.ac.uk> >>>> wrote: >>>> >> >>>> >> Hi All, >>>> >> >>>> >> I have a branch going through testing now which will handle >>>> PETSC_CONFIGURE_OPTIONS in a smarter way. Basically we'll just make sure >>>> that the things which are required are there, and we'll also honour >>>> anything the user has already set. I've also updated docs to say that >>>> that's what happens. >>>> >> >>>> >> I am also open to the suggestion that we could add more >>>> configuration options to the default set. Would anyone like to suggest what >>>> they should be? >>>> > >>>> > I think --download-hypre=1 at least. Maybe also a sparse direct >>>> solver (mumps?) which I think needs --download-metis --download-parmetis >>>> --download-mumps (maybe some others?) >>>> > >>>> > Lawrence >>>> > >>>> > >>>> > _______________________________________________ >>>> > 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 >>>> >>>> >>>> _______________________________________________ >>>> 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 2

Re: [firedrake] Expession
by Andrew McRae 01 Oct '15

01 Oct '15

In the scalar approach, are you also redefining Q at each time step? ...just checking On 30 September 2015 at 15:46, Anna Kalogirou <a.kalogirou(a)leeds.ac.uk> wrote: > Yes I update the time component of the expession at every time step. > > > On 30/09/15 15:38, Andrew McRae wrote: > > Are you performing the interpolation every time step? > > On 30 September 2015 at 15:36, Anna Kalogirou <a.kalogirou(a)leeds.ac.uk> > wrote: > >> The solution behaves as expected when I define Q as a function, but in >> the other case there seems to be some kind of instability and I get a >> very strange result. >> >> Anna. >> >> >> On 30/09/15 15:28, Lawrence Mitchell wrote: >> > -----BEGIN PGP SIGNED MESSAGE----- >> > Hash: SHA1 >> > >> > On 30/09/15 15:19, Anna Kalogirou wrote: >> >> Dear all, >> >> >> >> I have a very simple, yet not so obvious - to me at least - >> >> question: >> >> >> >> If I define an expression which does not depend on space variables, >> >> for example >> >> >> >> Q_expr = Expression("A*sin(2*pi*F*t)", t=t, A=A, F=F) >> >> Q.interpolate(Q_expr) >> >> >> >> shouldn't this be equivalent to just writing Q as a scalar, >> >> >> >> Q = A*sin(2*pi*F*t) ? >> >> >> >> This quantity Q is then used in a weak form to multiply the test >> >> function v on a boundary, e.g. Q*v*ds_v(1). >> >> >> >> I get different results for the two cases. Clearly, in the 1st case >> >> Q is a function while in the 2nd case it is a scalar, but I don't >> >> see why this would give different results. I might be missing >> >> something obvious here... >> > These ought to be equivalent. But can you describe how they are >> > different? >> > >> > >> > Lawrence. >> > >> > -----BEGIN PGP SIGNATURE----- >> > Version: GnuPG v1 >> > >> > iQEcBAEBAgAGBQJWC/GPAAoJECOc1kQ8PEYvGncIAM7WI2JkrZ5Y1wZMyJAHNDJG >> > Ca476ZQcNCzKsZ/ntPL9SoXsbXu4Uu2mIVuT4r30VBc/NAEuBdkKppVmiwS3DTsP >> > qzLVP2nRL1ENxYWOUk8pWryE/jffyPP9ny5gHoRLPP/09vL3e4y8jKplw0bmKMgZ >> > UVvcp+kL1hqFujvCddxpym8gg4Y3UrMTdD2MaYJ8WHaFeqB9ZPNWomv24ZKKWhHU >> > HAms5sgsNedvdrgIY32AH2AARRysk4AsLsfjYymsiXPrOyPrlG5V2Muf4gKdqID1 >> > MRWiL7p6ZVxLnOytF4NtFRUItYLKcOrv487IIIHVfl+TVvmKO+IJpXDyx5wyZcY= >> > =EkO8 >> > -----END PGP SIGNATURE----- >> > >> > _______________________________________________ >> > firedrake mailing list >> > firedrake(a)imperial.ac.uk >> > https://mailman.ic.ac.uk/mailman/listinfo/firedrake >> >> -- >> >> Dr Anna Kalogirou >> Research Fellow >> School of Mathematics >> University of Leeds >> >> http://www1.maths.leeds.ac.uk/~matak/ >> >> >> _______________________________________________ >> firedrake mailing list >> firedrake(a)imperial.ac.uk >> https://mailman.ic.ac.uk/mailman/listinfo/firedrake >> > > > > _______________________________________________ > firedrake mailing listfiredrake@imperial.ac.ukhttps://mailman.ic.ac.uk/mailman/listinfo/firedrake > > > -- > > Dr Anna Kalogirou > Research Fellow > School of Mathematics > University of Leeds > > http://www1.maths.leeds.ac.uk/~matak/ > >

3 3

Re: [firedrake] Firedrake issue
by David Ham 30 Sep '15

30 Sep '15

OK. What I think is going on here is that mathematically your test functions are only defined on the surface. So having the test function inside the z integral does not really make sense. What we would need to do is integrate over z first and then feed that in as a term in the boundary integral. That's not something Firedrake can do yet. We should give some thought to what it will take and maybe that's something you might be able to contribute. Regards, David On Mon, 28 Sep 2015 at 07:04 Tomasz Salwa [RPG] <mmtjs(a)leeds.ac.uk> wrote: > Hi David, > > > This is the question, I don't know how to formulate weak form properly, > provided I have the discrete one. It was my guess so far. I used the test > function defined for the 3D domain (should it be different?). > > Let's stick to the discrete equations. For index \alpha traversing the > free surface nodes, the red coupling terms contribute only to the corner > nodes which are at the free surface and coupling surface. Then, in the sum > over 'm' (coupling surface), for CG1 elements, there will be a contribution > from corner nodes and the neighbouring nodes below only. The question is > how to input it so that FD gets it right. > > > Does it clarify a bit? > > > Regards, > > Tomasz > > > ------------------------------ > *From:* David Ham <David.Ham(a)imperial.ac.uk> > *Sent:* 28 September 2015 14:15 > > *To:* Tomasz Salwa [RPG]; firedrake > *Subject:* Re: Firedrake issue > Hi Tomasz, > > I've gone over your working, and I can't understand how your weak form > works. The test function \delta phi seems to simultaneously be defined over > the xy free surface and the yz coupling surface. Since each discrete > equation corresponds to one matrix row, I don't understand to which matrix > rows the contributions are supposed to contribute. > > Can you enlighten me? > > Regards, > > David > > On Thu, 24 Sep 2015 at 08:06 Tomasz Salwa [RPG] <mmtjs(a)leeds.ac.uk> wrote: > >> Hi David, >> >> >> Concerning the FD update, it will take a while, as I encounter problems >> installing petsc4py and I need to wait for IT guy. So I'm working on the >> version of 15/06/15 for the moment. >> >> >> I'm new to github, but this should work: >> https://github.com/mmtjs/lin_coupled >> <https://github.com/mmtjs/lin_coupled> >> mmtjs/lin_coupled >> Contribute to lin_coupled development by creating an account on GitHub. >> Read more... <https://github.com/mmtjs/lin_coupled> >> I had an error in the previous formulation, so I redirived it and at the >> moment the solver is partitioned i.e. it doesn't require one matrix for the >> whole system, but just mappings of corresponding boundary nodes from one >> domain to the other (these work). >> >> >> The issue is in lines 175-179, where the coupling terms are commented out >> ... >> a_eta = trial * v * ds_t >> L_eta = ( eta + dt*dot( n_w, grad(phi) ) ) * v * ds_t #+ dt/rho * phi * >> v * ds_v(2) - dt * Pw * v * ds_v(2) >> ... >> More on it in the previous pdf file. >> >> Cheers, >> Tomasz >> >> >> ------------------------------ >> *From:* David Ham <David.Ham(a)imperial.ac.uk> >> *Sent:* 23 September 2015 10:27 >> *To:* Tomasz Salwa [RPG]; firedrake >> *Subject:* Re: Firedrake issue >> >> Hi Tomasz, >> >> Welcome back. >> >> Let's get through a few preliminaries. >> >> First, if you haven't used firedrake for 2-3 months then it will be quite >> out of date. You need to update to the current version. You can do this >> using your existing install and update all the packages, or you can use the >> new installer to get a new installation which will also have an attached >> update script. See http://firedrakeproject.org/download.html >> Obtaining Firedrake — Firedrake 0.11.0-658-ge79e125-dirty ... >> Running firedrake-install with no arguments will install firedrake in a >> python virtualenv created in a firedrake subdirectory of the current >> directory. >> Read more... <http://firedrakeproject.org/download.html> >> >> Next, it's much better to have this sort of conversation on the Firedrake >> list. This keeps everyone in the loop and enables other people to intervene >> when I say something stupid ;). I have taken the liberty of adding you to >> the list and I'm replying to the list now so that this conversation ends up >> there. >> >> Now, we need to be talking about code, not only maths, so please push >> what you're working on to a git repository on github or bitbucket so we can >> see what's going on. >> >> Finally we can start to think about the actual problem. I think that last >> time we talked about this, we discussed setting up a mixed mass matrix >> across the two domains and checking that projection works. Does that still >> work (and can you post some code to a repository which does this)? >> >> Once we get past that then we need to talk about how to assemble forms on >> the boundary facets. This will involve doing some mucking around with maps, >> because at the moment the facets belong to two different meshes. We need to >> think a little about how to do that, but we can do so when the first steps >> work. >> >> Cheers, >> >> David >> >> >> On Tue, 22 Sep 2015 at 12:54 Tomasz Salwa [RPG] <mmtjs(a)leeds.ac.uk> >> wrote: >> >>> Dear David, >>> >>> >>> I was busy in the recent months and could not work on the >>> fluid-structure interaction problem. Now I returned to it and reformulated >>> it. >>> >>> >>> I am getting the following error: >>> >>> >>> RuntimeError: Nonlinear solve failed to converge after 0 nonlinear >>> iterations. >>> Reason: >>> Inner linear solve failed to converge after 10000 iterations with >>> reason: DIVERGED_MAX_IT >>> >>> >>> but the problem is with the proper formulation in Firedrake. >>> >>> I include the pdf with it. Would you be able to help me with this? >>> >>> >>> Regards, >>> >>> Tomasz >>> >>> >>> >>> ----- >>> >>> Tomasz Salwa >>> >>> PhD Research Student >>> >>> University of Leeds >>> >>> http://www1.maths.leeds.ac.uk/~mmtjs >>> >>

1 0

Re: [firedrake] Expession
by Andrew McRae 30 Sep '15

30 Sep '15

Are you performing the interpolation every time step? On 30 September 2015 at 15:36, Anna Kalogirou <a.kalogirou(a)leeds.ac.uk> wrote: > The solution behaves as expected when I define Q as a function, but in > the other case there seems to be some kind of instability and I get a > very strange result. > > Anna. > > > On 30/09/15 15:28, Lawrence Mitchell wrote: > > -----BEGIN PGP SIGNED MESSAGE----- > > Hash: SHA1 > > > > On 30/09/15 15:19, Anna Kalogirou wrote: > >> Dear all, > >> > >> I have a very simple, yet not so obvious - to me at least - > >> question: > >> > >> If I define an expression which does not depend on space variables, > >> for example > >> > >> Q_expr = Expression("A*sin(2*pi*F*t)", t=t, A=A, F=F) > >> Q.interpolate(Q_expr) > >> > >> shouldn't this be equivalent to just writing Q as a scalar, > >> > >> Q = A*sin(2*pi*F*t) ? > >> > >> This quantity Q is then used in a weak form to multiply the test > >> function v on a boundary, e.g. Q*v*ds_v(1). > >> > >> I get different results for the two cases. Clearly, in the 1st case > >> Q is a function while in the 2nd case it is a scalar, but I don't > >> see why this would give different results. I might be missing > >> something obvious here... > > These ought to be equivalent. But can you describe how they are > > different? > > > > > > Lawrence. > > > > -----BEGIN PGP SIGNATURE----- > > Version: GnuPG v1 > > > > iQEcBAEBAgAGBQJWC/GPAAoJECOc1kQ8PEYvGncIAM7WI2JkrZ5Y1wZMyJAHNDJG > > Ca476ZQcNCzKsZ/ntPL9SoXsbXu4Uu2mIVuT4r30VBc/NAEuBdkKppVmiwS3DTsP > > qzLVP2nRL1ENxYWOUk8pWryE/jffyPP9ny5gHoRLPP/09vL3e4y8jKplw0bmKMgZ > > UVvcp+kL1hqFujvCddxpym8gg4Y3UrMTdD2MaYJ8WHaFeqB9ZPNWomv24ZKKWhHU > > HAms5sgsNedvdrgIY32AH2AARRysk4AsLsfjYymsiXPrOyPrlG5V2Muf4gKdqID1 > > MRWiL7p6ZVxLnOytF4NtFRUItYLKcOrv487IIIHVfl+TVvmKO+IJpXDyx5wyZcY= > > =EkO8 > > -----END PGP SIGNATURE----- > > > > _______________________________________________ > > firedrake mailing list > > firedrake(a)imperial.ac.uk > > https://mailman.ic.ac.uk/mailman/listinfo/firedrake > > -- > > Dr Anna Kalogirou > Research Fellow > School of Mathematics > University of Leeds > > http://www1.maths.leeds.ac.uk/~matak/ > > > _______________________________________________ > firedrake mailing list > firedrake(a)imperial.ac.uk > https://mailman.ic.ac.uk/mailman/listinfo/firedrake >

2 1

Re: [firedrake] Solve a Variational problem in a part of the domain
by Colin Cotter 30 Sep '15

30 Sep '15

Well, can you sketch out how you would solve this system, then we can see how to build it. all the best --cjc On 17 August 2015 at 13:33, Anna Kalogirou <a.kalogirou(a)leeds.ac.uk> wrote: > Yes, obviously I wouldn't be able to write > > assemble(u*dx) > > for the trial function u, but I just wanted to show you what I have. > > Is it too difficult to design such a non-conventional solver? > > Thanks, > > Anna. > > > > On 17/08/15 13:04, Colin Cotter wrote: > > Hi Anna, > Not as written, no. I think that we would need to design a bespoke > solver for that. > > all the best > --cjc > > On 17 August 2015 at 12:18, Anna Kalogirou <a.kalogirou(a)leeds.ac.uk> > wrote: > >> Hi all, >> >> I have to solve a problem with the bilinear form similar to the following: >> >> a = ( v*u + dt*inner(grad(u),grad(v)) )*dx + ( v*dt*assemble(u*dx) )*dx >> >> The last term is essentially the product of the integral of u and the >> integral of test function v. >> Is it even possible to solve this in Firedrake? >> >> Best, >> >> Anna. >> >> >> On 14/08/15 15:43, Lawrence Mitchell wrote: >> >> On 14 Aug 2015, at 11:55, Anna Kalogirou <a.kalogirou(a)leeds.ac.uk> <a.kalogirou(a)leeds.ac.uk> wrote: >> >> Hi, >> >> I was planning to define the heavyside function as a DG0 function. >> >> Any ideas about solving equation (1b), which contains a time-update for both a function and a scalar? >> >> I haven't been fully following along, but it looks like the formulation in (1f) introduces, effectively, a global coupling. Is that right? Assuming you didn't have any implementation constraints, how would you go about solving this? >> >> Lawrence >> >> >> >> _______________________________________________ >> firedrake mailing listfiredrake@imperial.ac.ukhttps://mailman.ic.ac.uk/mailman/listinfo/firedrake >> >> >> -- >> >> Dr Anna Kalogirou >> Research Fellow >> School of Mathematics >> University of Leeds >> >> http://www1.maths.leeds.ac.uk/~matak/ >> >> > > > -- > http://www.imperial.ac.uk/people/colin.cotter > > www.cambridge.org/9781107663916 > > > > > _______________________________________________ > firedrake mailing listfiredrake@imperial.ac.ukhttps://mailman.ic.ac.uk/mailman/listinfo/firedrake > > > -- > > Dr Anna Kalogirou > Research Fellow > School of Mathematics > University of Leeds > > http://www1.maths.leeds.ac.uk/~matak/ > > -- http://www.imperial.ac.uk/people/colin.cotter www.cambridge.org/9781107663916

3 3