October 2019 - firedrake - imperial.ac.uk

Preconditioner for the Stokes-Brinkmann problem
by Salazar De Troya, Miguel 24 Oct '19

24 Oct '19

Hello, I am restating this question here for better context. Last email was a bit sloppy with all the stuff coming from the petsc mail list. I am trying to solve the Stokes equation with the Brinkman term to simulate a solid material. My intention is to implement the preconditioner in this paper: https://onlinelibrary.wiley.com/doi/epdf/10.1002/fld.426 (section 2.6) where they solve for the velocity and substitute that expression in the pressure equation. They end up solving a system of the type B K^-1 B^T, i.e. the Schur complement of the problem. For this system of equations, they argue that the preconditioner in page 11 is perfect for a given constant Brinkman penalty term. Because I am solving for velocity and pressure without doing any substitution, I thought I could use a PC fieldsplit type Schur (full factorization) and provide the preconditioner in the paper to solve the Schur complement. In the paper, the action of the PC is defined by: M^{-1} q = \phi_q + \mu q \mu is some number, and \phi_q solves the BVP -div \alpha^{-1} grad \phi_q = q grad \phi_q \dot n = 0 \int \phi_q = 0 Therefore, the PC for the fieldsplit_1 needs to: 1. Discretise and solve this BVP to compute \phi_q 2. Add \mu q. I.e. the action of PCApply is: y <- \mu x + Laplace^{-1} x I followed the example in https://www.firedrakeproject.org/_modules/firedrake/preconditioners/pcd.html to try to implement this preconditioner . The code is attached below. The fieldsplit_1_ solve takes many iterations to converge whenever a Gaussian bump or a jump function models the \alpha parameter. It behaves better for \alpha constant (an assumption in the paper for the preconditioner to be perfect). I am solving the pure Neumann problem by passing the nullspace to the KSP. Are there mistakes in the way I am implementing the preconditioner? Or maybe the lack of convergence is just purely because of the method? Thanks Miguel Miguel A. Salazar de Troya Postdoctoral Researcher, Lawrence Livermore National Laboratory B141 Rm: 1085-5 Ph: 1(925) 422-6411 from firedrake import * from firedrake.petsc import PETSc from firedrake.preconditioners.base import PCBase class BrinkmannPC(PCBase): r""" Preconditioner for the Stokes-Brinkmann problem Fill details """ def initialize(self, pc): from firedrake import TrialFunction, TestFunction, dx, \ assemble, inner, grad, split, Constant, parameters from firedrake.assemble import allocate_matrix, create_assembly_callable if pc.getType() != "python": raise ValueError("Expecting PC type python") prefix = pc.getOptionsPrefix() + "brink_" # we assume P has things stuffed inside of it _, P = pc.getOperators() context = P.getPythonContext() test, trial = context.a.arguments() if test.function_space() != trial.function_space(): raise ValueError("Pressure space test and trial space differ") Q = test.function_space() p = TrialFunction(Q) q = TestFunction(Q) # Regularisation to avoid having to think about nullspaces. alpha = context.appctx.get("alpha", Constant(1.0)) stiffness = inner(1.0 / alpha * grad(p), grad(q))*dx opts = PETSc.Options() # we're inverting Kp, so default them to assembled. # Mp only needs its action, so default it to mat-free. # These can of course be overridden. default = parameters["default_matrix_type"] #default = "aij" Kp_mat_type = opts.getString(prefix+"Kp_mat_type", default) self.Mp_mat_type = opts.getString(prefix+"Mp_mat_type", "matfree") Kp = assemble(stiffness, form_compiler_parameters=context.fc_params, mat_type=Kp_mat_type, options_prefix=prefix + "Kp_") Kp_petsc = Kp.M.handle L = q * dx self.rhs = assemble(L) from firedrake import VectorSpaceBasis, Function x_axis = Function(Q).interpolate(Constant(1.0)) nullspace = VectorSpaceBasis(vecs=[x_axis]) nullspace.orthonormalize() nullspace_petsc = nullspace.nullspace() Kp_petsc.setNearNullSpace(nullspace_petsc) self.nullspace = nullspace Kksp = PETSc.KSP().create(comm=pc.comm) Kksp.incrementTabLevel(1, parent=pc) Kksp.setOptionsPrefix(prefix + "Kp_") Kksp.setOperators(Kp_petsc) Kksp.setFromOptions() Kksp.setUp() self.Kksp = Kksp mu = context.appctx.get("mu", Constant(1.0)) mass = mu*p*q*dx self.Mp = allocate_matrix(mass, form_compiler_parameters=context.fc_params, mat_type=self.Mp_mat_type, options_prefix=prefix + "Mp_") self._assemble_Mp = create_assembly_callable(mass, tensor=self.Mp, form_compiler_parameters=context.fc_params, mat_type=self.Mp_mat_type) self._assemble_Mp() Mpmat = self.Mp.petscmat self.workspace = [Mpmat.createVecLeft() for i in (0, 1)] def update(self, pc): pass def apply(self, pc, x, y): a, b = self.workspace with self.rhs.dat.vec_wo as rhs_dat: x.copy(rhs_dat) self.nullspace.orthogonalize(self.rhs) with self.rhs.dat.vec_wo as rhs_dat: self.Kksp.solve(rhs_dat, a) self.Mp.petscmat.mult(x, y) y.axpy(1.0, a) def applyTranspose(self, pc, x, y): a, b = self.workspace with self.rhs.dat.vec_wo as rhs_dat: x.copy(rhs_dat) self.nullspace.orthogonalize(self.rhs) with self.rhs.dat.vec_wo as rhs_dat: self.Kksp.solveTranspose(rhs_dat, a) self.Mp.petscmat.mult(x, y) y.axpy(1.0, a) def view(self, pc, viewer=None): super(BrinkmannPC, self).view(pc, viewer) viewer.printfASCII("KSP solver for K^-1:\n") self.Kksp.view(viewer) N = 100 mesh = UnitSquareMesh(N, N) V = VectorFunctionSpace(mesh, "CG", 2) W = FunctionSpace(mesh, "CG", 1) Z = V * W print("# DOFS {}".format(Z.dim())) w = Function(Z) u, p = split(w) v, q = TestFunctions(Z) mu = Constant(1.0) x, y = SpatialCoordinate(mesh) from ufl import And alpha = conditional(And(y > 0.4, And(x < 0.2, y < 0.8)), Constant(1e2), Constant(1e-1)) #alpha = Constant(1e1) #alpha = 1e2*exp((0.2 - (x - 0.5)*(x - 0.5) - (y - 0.5)*(y - 0.5)) / 0.05) alpha_interp = Function(W) alpha_interp.interpolate(alpha) File("alpha.pvd").write(alpha_interp) F = mu*inner(grad(u), grad(v))*dx - p*div(v)*dx - div(u)*q*dx + alpha*inner(u, v)*dx bc_value = as_vector([0, x * (x - 1)]) bcs = [DirichletBC(Z.sub(0), bc_value, (4,)), DirichletBC(Z.sub(0), zero(mesh.geometric_dimension()), (1, 2)), DirichletBC(Z.sub(0), bc_value, (3,))] def create_solver(solver_parameters, appctx=None, nullspace=nullspace): p = {} if solver_parameters is not None: p.update(solver_parameters) # Default to linear SNES p.setdefault("snes_type", "ksponly") p.setdefault("ksp_atol", 1e-6) problem = NonlinearVariationalProblem(F, w, bcs=bcs) solver = NonlinearVariationalSolver(problem, solver_parameters=p, appctx=appctx, nullspace=nullspace) return solver fieldsplit_0_lu = { "ksp_type": "preonly", "ksp_max_it": 1, "pc_type": "lu", "pc_factor_mat_solver_type": "mumps", } fieldsplit_1_lu = fieldsplit_0_lu fieldsplit_0_python = { "ksp_type": "preonly", "pc_type": "python", "pc_python_type": "firedrake.AssembledPC", "assembled_pc_type": "hypre", } fieldsplit_1_python = { "ksp_type" : "cg", "ksp_atol" : 1e-2, "pc_type": "python", "pc_python_type": "__main__.BrinkmannPC", "ksp_converged_reason" : None, "brink_Kp_ksp_type": "preonly", "brink_Kp_pc_type": "lu", "brink_Kp_pc_factor_mat_solver_type" : "mumps", } iterative = True params = { "mat_type" : "matfree" if iterative else "aij", "ksp_monitor_true_residual": None, "ksp_converged_reason": None, "ksp_max_it" : 100, "ksp_atol" : 1e-6, "ksp_type" : "fgmres", "pc_type" : "fieldsplit", "pc_fieldsplit_type" : "schur", "pc_fieldsplit_schur_fact_type": "full", "pc_fieldsplit_schur_precondition": "selfp" if not iterative else "a11", "fieldsplit_0": fieldsplit_0_python if iterative else fieldsplit_0_lu, "fieldsplit_1" : fieldsplit_1_python if iterative else fieldsplit_1_lu, } w.assign(0) appctx = {"alpha" : alpha, "mu" : mu} nullspace = MixedVectorSpaceBasis( Z, [Z.sub(0), VectorSpaceBasis(constant=True)]) solver = create_solver(params, appctx=appctx, nullspace=nullspace) solver.solve() u, p = w.split() u.rename("velocity") p.rename("pressure") File("brinkman_problem.pvd").write(u, p)

1 0

Point evaluation
by Rodrigo Piccinini 22 Oct '19

22 Oct '19

Hello, I'm trying to evaluate a Function at a point, but I get the error bellow. Could anybody help me on that? The log files are attached. My code can be found here: https://github.com/rbpiccinini/codes/blob/master/eigendrake/spe10_transient… CompilationError: Command "['mpicc', '-std=c99', '-fPIC', '-Wall', > '-march=native', '-O3', '-ffast-math', > '-I/home/piccinini/apps/firedrake/src/firedrake/firedrake', > '-I/home/piccinini/apps/firedrake/include', > '-I/home/piccinini/apps/firedrake/lib/python3.7/site-packages/petsc/include', > '-o', > '/home/piccinini/apps/firedrake/.cache/pyop2/0a/2ff88afc9cdab5dd593941b370e84c_p6855.so.tmp', > '/home/piccinini/apps/firedrake/.cache/pyop2/0a/2ff88afc9cdab5dd593941b370e84c_p6855.c', > '-shared', '-L/home/piccinini/apps/firedrake/lib', '-lspatialindex_c', > '-Wl,-rpath,/home/piccinini/apps/firedrake/lib']" return error status 1. > Unable to compile code > Thanks -- Rodrigo

2 4

Installation on Ubuntu: 18.04.3 LTS
by Bhavesh Shrimali 18 Oct '19

18 Oct '19

Hi all, I tried installing firedrake on my laptop with `Ubuntu: 18.04` but have been unsuccessful so far. These are the warnings I see in the `log` file =============================================================================== ***** WARNING : F90 (set to mpif90) found in environment variables - ignoring use ./configure F90=$F90 if you really want to use that value ****** =============================================================================== =============================================================================== ***** WARNING : Using default optimization C flags -g -O You might consider manually setting optimal optimization flags for your system with COPTFLAGS="optimization flags" see config/examples/arch-*-opt.py for examples =============================================================================== =============================================================================== ***** WARNING : Using default C++ optimization flags -g -O You might consider manually setting optimal optimization flags for your system with CXXOPTFLAGS="optimization flags" see config/examples/arch-*-opt.py for examples =============================================================================== =============================================================================== ***** WARNING : Using default FORTRAN optimization flags -g -O You might consider manually setting optimal optimization flags for your system with FOPTFLAGS="optimization flags" see config/examples/arch-*-opt.py for examples =============================================================================== and it fails to build PETSc beyond this point. ******************************************************************************* UNABLE to CONFIGURE with GIVEN OPTIONS (see configure.log for details): ------------------------------------------------------------------------------- Your libraries are from OpenMPI but it appears your mpiexec is from MPICH ******************************************************************************* I am attaching the entire `log` file for reference. Any leads into what could be possibly going wrong would be helpful. Thanks Bhavesh -- *Bhavesh Shrimali* *Ph.D. Student, **Civil and Environmental Engineering,* *University of Illinois Urbana-Champaign* *Champaign, IL, United States* *Email : bshrima2(a)illinois.edu <bshrima2(a)illinois.edu>* *Web : **https://bhaveshshrimali.github.io/ <https://bhaveshshrimali.github.io/>*

2 1

Talks from Firedrake 19
by Lawrence Mitchell 11 Oct '19

11 Oct '19

Dear all, the Firedrake '19 page (https://firedrakeproject.org/firedrake_19.html) has been updated to add talk slides from some of the talks from attendees. Unfortunately I was not able to make the in-room recording work, and so there is no audio. For those of you who presented and would like your slides to also appear, please send them to me (or provide a link to a website). Thanks, Lawrence

1 0

using different petsc version
by Florian Bruckner 09 Oct '19

09 Oct '19

Dear Firedrake-Team, this may be a question that is very easy to answer, but I didn't get it working by myself. I would like to use a different petsc version which contains some small changes. I am using firedrake inside a Docker container, where i run the firedrake-install script. I tried using PETSC_DIR with --honour-petsc-dir, but this fails with some missing include files. Is PETSC_DIR only the source directory of my own petsc version, or do I have to compile it before? Is it somehow possible to just specify an individual url / git-repository to use a different petsc version? thanks for any advice Florian

2 1

Re: [firedrake] [petsc-users] Stokes-Brinkmann equation preconditioner
by Salazar De Troya, Miguel 07 Oct '19

07 Oct '19

From: Matthew Knepley <knepley(a)gmail.com> Date: Friday, October 4, 2019 at 3:19 AM To: Lawrence Mitchell <wence(a)gmx.li> Cc: "Salazar De Troya, Miguel" <salazardetro1(a)llnl.gov>, "petsc-users(a)mcs.anl.gov" <petsc-users(a)mcs.anl.gov> Subject: Re: [petsc-users] Stokes-Brinkmann equation preconditioner On Fri, Oct 4, 2019 at 6:04 AM Lawrence Mitchell <wence(a)gmx.li<mailto:wence@gmx.li>> wrote: > On 4 Oct 2019, at 10:46, Matthew Knepley via petsc-users <petsc-users(a)mcs.anl.gov<mailto:petsc-users@mcs.anl.gov>> wrote: > > On Thu, Oct 3, 2019 at 6:34 PM Salazar De Troya, Miguel via petsc-users <petsc-users(a)mcs.anl.gov<mailto:petsc-users@mcs.anl.gov>> wrote: > I am trying to solve the Stokes equation with the Brinkman term to simulate a solid material. My intention is to implement the preconditioner in this paper:https://onlinelibrary.wiley.com/doi/epdf/10.1002/fld.426 (section 2.6) > > > Link does not work for me. Try https://onlinelibrary.wiley.com/doi/epdf/10.1002/fld.426 (mail clients are awful) > > where they solve for the velocity and substitute that expression in the pressure equation. They end up solving a system of the type B K^-1 B^T, i.e. the Schur complement of the problem. For this system of equations, they argue that the preconditioner in page 11 is perfect for a given constant Brinkman penalty term. > > > > Because I am solving for velocity and pressure without doing any substitution, I thought I could use a PC fieldsplit type Schur (full factorization) > > > Yes, this will form the exact factorization and a matrix-free form of the Schur complement. > > and provide the preconditioner in the paper to solve the Schur complement. > > > Yes, you can provide a user-defined PC for the Schur complement. > > My question is, should I provide this preconditioner through PCFieldSplitSetSchurPre or through fieldsplit_1_pc_type (probably through the Firedrake interface as inhttps://www.firedrakeproject.org/demos/stokes.py.html<http://www.firedrakeproject.org/demos/stokes.py.html>) ? > > > The name PCFieldSplitSetSchurPre seems to be very misleading. You do not use it to provide a _preconditioner_. You use it to determine > the _preconditioning matrix_ from which the actual preconditioner is built. The preconditioner itself is defined using -fieldpslit_1_pc_type. > Since I do not know what the preconditioner looks like, I cannot say what preconditioner matrix you would want. Since Firedrake can construct > any operator for you, you might not care about the matrix we pass to you. The action of the PC is defined by: M^{-1} q = \phi_q + \mu q \mu is some number, and \phi_q solves the BVP -div \alpha^{-1} grad \phi_q = q grad \phi_q \dot n = 0 \int \phi_q = 0 So your PC for the fieldsplit_1 needs to: 1. Discretise and solve this BVP to compute \phi_q 2. Add \mu q. I.e. the action of PCApply is: y <- \mu x + Laplace^{-1} x (Aside, surely a primal-dual error has been made in the analysis here) You could, I think, do this by providing the discretisation of this laplacian and using an additive PCComposite (although I don't know what PC you would use to just scale the input, it's easy to write one though). Richardson R(x, y) is y += mu (b - Ax) so if b = 0 and A = I, you get that. It seems easier just to make a PCSHELL with VecAXPY I guess. Thanks, Matt Cheers, Lawrence Why is there a need to use a PCComposite instead of applying “y <- \mu x + Laplace^{-1} x” directly in the same preconditioner? I have tried to do the latter without success, maybe that is why. I have used Firedrake and followed the example in https://www.firedrakeproject.org/_modules/firedrake/preconditioners/pcd.html . I hope it is ok to include it here given that it uses a lot of the petsc4py routines. The code is attached below. The implementation of the preconditioner for the Schur complement must be wrong somewhere because the fieldsplit_1_ solve takes many iterations to converge whenever a Gaussian bump or a jump function models the \alpha parameter. It behaves better for \alpha constant (an assumption in the paper for the preconditioner to be perfect). I am solving the pure Neumann problem by passing the nullspace to the KSP. I hope it is okay to ask here if there is something wrong in how I am applying “y <- \mu x + Laplace^{-1} x” Thanks Miguel from firedrake import * from firedrake.petsc import PETSc from firedrake.preconditioners.base import PCBase class BrinkmannPC(PCBase): r""" Preconditioner for the Stokes-Brinkmann problem Fill details """ def initialize(self, pc): from firedrake import TrialFunction, TestFunction, dx, \ assemble, inner, grad, split, Constant, parameters from firedrake.assemble import allocate_matrix, create_assembly_callable if pc.getType() != "python": raise ValueError("Expecting PC type python") prefix = pc.getOptionsPrefix() + "brink_" # we assume P has things stuffed inside of it _, P = pc.getOperators() context = P.getPythonContext() test, trial = context.a.arguments() if test.function_space() != trial.function_space(): raise ValueError("Pressure space test and trial space differ") Q = test.function_space() p = TrialFunction(Q) q = TestFunction(Q) alpha = context.appctx.get("alpha", Constant(1.0)) stiffness = inner(1.0 / alpha * grad(p), grad(q))*dx opts = PETSc.Options() # we're inverting Kp, so default them to assembled. # Mp only needs its action, so default it to mat-free. # These can of course be overridden. default = parameters["default_matrix_type"] Kp_mat_type = opts.getString(prefix+"Kp_mat_type", default) self.Mp_mat_type = opts.getString(prefix+"Mp_mat_type", "matfree") Kp = assemble(stiffness, form_compiler_parameters=context.fc_params, mat_type=Kp_mat_type, options_prefix=prefix + "Kp_") Kp_petsc = Kp.M.handle L = q * dx self.rhs = assemble(L) from firedrake import VectorSpaceBasis, Function ctsp = Function(Q).interpolate(Constant(1.0)) nullspace = VectorSpaceBasis(vecs=[ctsp]) nullspace.orthonormalize() nullspace_petsc = nullspace.nullspace() Kp_petsc.setNearNullSpace(nullspace_petsc) self.nullspace = nullspace Kksp = PETSc.KSP().create(comm=pc.comm) Kksp.incrementTabLevel(1, parent=pc) Kksp.setOptionsPrefix(prefix + "Kp_") Kksp.setOperators(Kp_petsc) Kksp.setFromOptions() Kksp.setUp() self.Kksp = Kksp mu = context.appctx.get("mu", Constant(1.0)) mass = mu*p*q*dx self.Mp = allocate_matrix(mass, form_compiler_parameters=context.fc_params, mat_type=self.Mp_mat_type, options_prefix=prefix + "Mp_") self._assemble_Mp = create_assembly_callable(mass, tensor=self.Mp, form_compiler_parameters=context.fc_params, mat_type=self.Mp_mat_type) self._assemble_Mp() Mpmat = self.Mp.petscmat self.workspace = [Mpmat.createVecLeft() for i in (0, 1)] def update(self, pc): pass def apply(self, pc, x, y): a, b = self.workspace with self.rhs.dat.vec_wo as rhs_dat: x.copy(rhs_dat) self.nullspace.orthogonalize(self.rhs) with self.rhs.dat.vec_wo as rhs_dat: self.Kksp.solve(rhs_dat, a) self.Mp.petscmat.mult(x, y) y.axpy(1.0, a) def applyTranspose(self, pc, x, y): a, b = self.workspace with self.rhs.dat.vec_wo as rhs_dat: x.copy(rhs_dat) self.nullspace.orthogonalize(self.rhs) with self.rhs.dat.vec_wo as rhs_dat: self.Kksp.solveTranspose(rhs_dat, a) self.Mp.petscmat.mult(x, y) y.axpy(1.0, a) def view(self, pc, viewer=None): super(BrinkmannPC, self).view(pc, viewer) viewer.printfASCII("KSP solver for K^-1:\n") self.Kksp.view(viewer) N = 100 mesh = UnitSquareMesh(N, N) V = VectorFunctionSpace(mesh, "CG", 2) W = FunctionSpace(mesh, "CG", 1) Z = V * W print("# DOFS {}".format(Z.dim())) w = Function(Z) u, p = split(w) v, q = TestFunctions(Z) mu = Constant(1.0) x, y = SpatialCoordinate(mesh) from ufl import And alpha = conditional(And(y > 0.4, And(x < 0.2, y < 0.8)), Constant(1e5), Constant(1e5)) alpha = Constant(1e3) alpha = 1e2*exp((0.2 - (x - 0.5)*(x - 0.5) - (y - 0.5)*(y - 0.5)) / 0.05) alpha_interp = Function(W) alpha_interp.interpolate(alpha) File("alpha.pvd").write(alpha_interp) F = mu*inner(grad(u), grad(v))*dx - p*div(v)*dx - div(u)*q*dx + alpha*inner(u, v)*dx bc_value = as_vector([0, x * (x - 1)]) bcs = [DirichletBC(Z.sub(0), bc_value, (4,)), DirichletBC(Z.sub(0), zero(mesh.geometric_dimension()), (1, 2)), DirichletBC(Z.sub(0), bc_value, (3,))] def create_solver(solver_parameters, appctx=None): p = {} if solver_parameters is not None: p.update(solver_parameters) # Default to linear SNES p.setdefault("snes_type", "ksponly") p.setdefault("ksp_atol", 1e-6) problem = NonlinearVariationalProblem(F, w, bcs=bcs) solver = NonlinearVariationalSolver(problem, solver_parameters=p, appctx=appctx) return solver fieldsplit_0_lu = { "ksp_type": "preonly", "ksp_max_it": 1, "pc_type": "lu", "pc_factor_mat_solver_type": "mumps", } fieldsplit_1_lu = fieldsplit_0_lu fieldsplit_0_python = { "ksp_type": "preonly", "pc_type": "python", "pc_python_type": "firedrake.AssembledPC", "assembled_pc_type": "hypre", } fieldsplit_1_python = { "ksp_type" : "cg", "ksp_atol" : 1e-2, "pc_type": "python", "pc_python_type": "__main__.BrinkmannPC", "ksp_converged_reason" : None, "brink_Kp_ksp_type": "preonly", "brink_Kp_pc_type": "lu", "brink_Kp_pc_factor_mat_solver_type" : "mumps", } iterative = True params = { "mat_type" : "matfree" if iterative else "aij", "ksp_monitor_true_residual": None, "ksp_converged_reason": None, "ksp_max_it" : 100, "ksp_atol" : 1e-6, "ksp_type" : "fgmres", "pc_type" : "fieldsplit", "pc_fieldsplit_type" : "schur", "pc_fieldsplit_schur_fact_type": "full", "pc_fieldsplit_schur_precondition": "selfp" if not iterative else "a11", "fieldsplit_0": fieldsplit_0_python if iterative else fieldsplit_0_lu, "fieldsplit_1" : fieldsplit_1_python if iterative else fieldsplit_1_lu, } w.assign(0) appctx = {"alpha" : alpha, "mu" : mu} solver = create_solver(params, appctx=appctx) solver.solve() u, p = w.split() u.rename("velocity") p.rename("pressure") File("out.pvd").write(u, p)

1 0

Installation fail inquiry
by Qiwei Zhan, Ph.D. 07 Oct '19

07 Oct '19

Hi colleagues, I have spent several weeks to install firedrake, but still failed. I have totally uninstall anaconda to make sure "conda deactivate", but my installation still fails, even earlier in the initial steps, mainly related to "pip is configured with locations that require TLS/SSL, however the ssl module in Python is not available." My installation command is "python3 firedrake-install" , I am using LinuxMint17. The log file is attached as "firedrake-install.log" On the opposite side, if I use Anaconda ("conda activate"), I can go to the very last steps, but fails with " firedrake/supermeshimpl.c:1606:3: error: unknown type name ‘MPI_Message’ MPI_Message ob_mpi" The log file is attached as the second file "firedrake-install-Anaconda.log". I would appreicate it, if anyone could help me with this, thank you! Thank you! Qiwei *Qiwei Zhan，Ph.D.* Peter O’Donnell, Jr. Postdoctoral Fellow Oden Institute for Computational Engineering and Sciences University of Texas at Austin Email: qiwei.zhan(a)utexas.edu; qiwei.zhan(a)alumni.duke.edu <qiwei.zhan(a)duke.edu>

3 4

Re: [firedrake] Gathering D.O.F.s onto one process
by Matthew Knepley 03 Oct '19

03 Oct '19

On Mon, Sep 30, 2019 at 10:46 PM Sentz, Peter <sentz2(a)illinois.edu> wrote: > So is there a way to return how the global degrees of freedom are > distributed across processors? > I would recommend changing your workflow so that you do not need this, since it complicates everything. However, if you still really want to do this, you can 1) Tell the DM to constructs mappings back to the original ordering https://www.mcs.anl.gov/petsc/petsc-master/docs/manualpages/DM/DMSetUseNatu… 2) After distribution, map a vector back to that ordering https://www.mcs.anl.gov/petsc/petsc-master/docs/manualpages/DMPLEX/DMPlexGl… 3) Output that vector using view() Note that DMDA outputs in the natural ordering by default since it is not expensive in the structured case. Thanks, Matt > -Peter > ------------------------------ > *From:* Matthew Knepley <knepley(a)gmail.com> > *Sent:* Monday, September 30, 2019 6:33 PM > *To:* Sentz, Peter <sentz2(a)illinois.edu> > *Cc:* Firedrake Project <firedrake(a)imperial.ac.uk> > *Subject:* Re: [firedrake] Gathering D.O.F.s onto one process > > On Mon, Sep 30, 2019 at 3:45 PM Sentz, Peter <sentz2(a)illinois.edu> wrote: > > Hi, > > u.vector().view() does not work because a Firedrake.vector.Vector has no > such method. > > However, I did try the following: > > > Q = PETSc.Viewer().createBinary('name.dat','w') > > with u.dat.vec_ro as vu: > vu.view(Q) > > Then in separate .py file loaded this .dat file and looked at the > resulting vector. Unfortunately it does not maintain the serial ordering. > This is perplexing because according to > https://www.firedrakeproject.org/demos/parprint.py.html, shouldn't the > .view() maintain serial ordering? > > > I think we may have a misunderstanding about what maintaining an ordering > means. If you output > using the PETSc parallel viewer, then it will definitely output in the > same order in which the vector > is distributed in parallel, meaning the dofs from proc0, then proc1, etc. > However, it now sounds like > this order does not match the order you have in serial, which is > completely understandable. When > you distribute meshes, it is normal to partition them to reduce > communication, which necessitates > reordering. > > Thanks, > > Matt > > > Peter > ------------------------------ > *From:* Matthew Knepley <knepley(a)gmail.com> > *Sent:* Saturday, September 28, 2019 12:14 AM > *To:* Sentz, Peter <sentz2(a)illinois.edu> > *Cc:* firedrake(a)imperial.ac.uk <firedrake(a)imperial.ac.uk> > *Subject:* Re: [firedrake] Gathering D.O.F.s onto one process > > On Thu, Sep 26, 2019 at 8:15 PM Sentz, Peter <sentz2(a)illinois.edu> wrote: > > Hello, > > I am working on a project where I need to save the Numpy arrays of degrees > of freedom of FE functions and use them in another script. I am wondering > how to get consistent results when using MPI. > > Here's a minimal example of what I would like to do: > > ---------------------------------------- > from firedrake import * > import numpy as np > > mesh = UnitSquareMesh(3, 3) > > V = FunctionSpace(mesh, "CG", 1) > u = TrialFunction(V) > v = TestFunction(V) > f = Function(V) > x, y = SpatialCoordinate(mesh) > f.interpolate((1+8*pi*pi)*cos(x*pi*2)*cos(y*pi*2)) > a = (dot(grad(v), grad(u)) + v * u) * dx > L = f * v * dx > > u = Function(V) > solve(a == L, u, solver_parameters={'ksp_type': 'cg'}) > > U = u.vector().gather() > > > Here I think you can just use > > u.vector().view(viewer) > > with a petsc4py Viewer of the appropriate type (maybe binary?). > > Thanks, > > Matt > > > if u.comm.rank == 0: > print(U) > name = 'U_global_' + str(u.comm.size) + '.npy' > np.save(name,U) > --------------------------------------------------------------------- > > If I run this in serial, and then with two processes in MPI, I get two > files 'U_global_1.npy' and 'U_global_2.npy' but they do not agree. The > first (in serial) has values: > > [-0.29551776 -0.67467152 -0.67467152 -0.15155862 0.06606523 0.06606523 > 0.36508954 0.36508954 1.33298501 1.33298501 0.06606523 -0.15155862 > 0.06606522 -0.67467153 -0.67467153 -0.29551776] > > > while the second has values: > > [-0.67467152 -0.67467152 -0.29551776 0.06606523 -0.15155862 0.06606523 > 0.06606523 -0.67467152 -0.15155862 -0.29551776 -0.67467152 0.06606523 > 1.33298501 0.36508954 0.36508954 1.33298502] > > > How do I get the parallel implementation to match the ordering of the > d.o.f.s in the serial case? And avoiding this .gather() approach would be > good, because I really only need the array gathered on a single process. > > Thanks, > > Peter > _______________________________________________ > firedrake mailing list > firedrake(a)imperial.ac.uk > https://mailman.ic.ac.uk/mailman/listinfo/firedrake > > > > -- > What most experimenters take for granted before they begin their > experiments is infinitely more interesting than any results to which their > experiments lead. > -- Norbert Wiener > > https://www.cse.buffalo.edu/~knepley/ > <http://www.cse.buffalo.edu/~knepley/> > > > > -- > What most experimenters take for granted before they begin their > experiments is infinitely more interesting than any results to which their > experiments lead. > -- Norbert Wiener > > https://www.cse.buffalo.edu/~knepley/ > <http://www.cse.buffalo.edu/~knepley/> > -- What most experimenters take for granted before they begin their experiments is infinitely more interesting than any results to which their experiments lead. -- Norbert Wiener https://www.cse.buffalo.edu/~knepley/ <http://www.cse.buffalo.edu/~knepley/>

4 4