Error to run a script converted from FEniCS
Hi all, My name is Simone Puel, a graduate student at the University of Texas at Austin. I think I have installed Firedrake correctly since I can run the example in the online tutorial. Thank you so much. I am trying to convert a script I made in FEniCS in Firedrake. I modified some syntax, but I have the following error: Traceback (most recent call last): File "/home/simone/firedrake/src/PyOP2/pyop2/caching.py", line 197, in __new__ return cls._cache_lookup(key) File "/home/simone/firedrake/src/firedrake/firedrake/tsfc_interface.py", line 53, in _cache_lookup return cls._cache.get(key) or cls._read_from_disk(key, comm) File "/home/simone/firedrake/src/firedrake/firedrake/tsfc_interface.py", line 74, in _read_from_disk raise KeyError("Object with key %s not found" % key) KeyError: 'Object with key aec71f5f7a30a3d0b4db90beac07de85 not found' During handling of the above exception, another exception occurred: Traceback (most recent call last): File "Elasticity2D_Firedrake.py", line 161, in <module> solve(lhs_varf == rhs_varf, sol, bcs=[bc1,bc2]) File "/home/simone/firedrake/src/firedrake/firedrake/solving.py", line 125, in solve _solve_varproblem(*args, **kwargs) File "/home/simone/firedrake/src/firedrake/firedrake/solving.py", line 152, in _solve_varproblem appctx=appctx) File "/home/simone/firedrake/src/firedrake/firedrake/variational_solver.py", line 335, in __init__ super(LinearVariationalSolver, self).__init__(*args, **kwargs) File "/home/simone/firedrake/src/firedrake/firedrake/variational_solver.py", line 182, in __init__ options_prefix=self.options_prefix) File "/home/simone/firedrake/src/firedrake/firedrake/solving_utils.py", line 131, in __init__ form_compiler_parameters=self.fcp) File "/home/simone/firedrake/src/firedrake/firedrake/assemble.py", line 158, in create_assembly_callable loops = tuple(loops) File "/home/simone/firedrake/src/firedrake/firedrake/assemble.py", line 228, in _assemble kernels = tsfc_interface.compile_form(f, "form", parameters=form_compiler_parameters, inverse=inverse) File "/home/simone/firedrake/src/firedrake/firedrake/tsfc_interface.py", line 212, in compile_form number_map, interface, coffee).kernels File "/home/simone/firedrake/src/PyOP2/pyop2/caching.py", line 199, in __new__ obj = make_obj() File "/home/simone/firedrake/src/PyOP2/pyop2/caching.py", line 189, in make_obj obj.__init__(*args, **kwargs) File "/home/simone/firedrake/src/firedrake/firedrake/tsfc_interface.py", line 121, in __init__ tree = tsfc_compile_form(form, prefix=name, parameters=parameters, interface=interface, coffee=coffee) File "/home/simone/firedrake/src/tsfc/tsfc/driver.py", line 57, in compile_form fd = ufl_utils.compute_form_data(form, complex_mode=complex_mode) File "/home/simone/firedrake/src/tsfc/tsfc/ufl_utils.py", line 59, in compute_form_data complex_mode=complex_mode File "/home/simone/firedrake/src/ufl/ufl/algorithms/compute_form_data.py", line 298, in compute_form_data form = apply_function_pullbacks(form) File "/home/simone/firedrake/src/ufl/ufl/algorithms/apply_function_pullbacks.py", line 267, in apply_function_pullbacks return map_integrand_dags(FunctionPullbackApplier(), expr) File "/home/simone/firedrake/src/ufl/ufl/algorithms/map_integrands.py", line 58, in map_integrand_dags form, only_integral_type) File "/home/simone/firedrake/src/ufl/ufl/algorithms/map_integrands.py", line 39, in map_integrands for itg in form.integrals()] File "/home/simone/firedrake/src/ufl/ufl/algorithms/map_integrands.py", line 39, in <listcomp> for itg in form.integrals()] File "/home/simone/firedrake/src/ufl/ufl/algorithms/map_integrands.py", line 46, in map_integrands return itg.reconstruct(function(itg.integrand())) File "/home/simone/firedrake/src/ufl/ufl/algorithms/map_integrands.py", line 57, in <lambda> return map_integrands(lambda expr: map_expr_dag(function, expr, compress), File "/home/simone/firedrake/src/ufl/ufl/corealg/map_dag.py", line 37, in map_expr_dag result, = map_expr_dags(function, [expression], compress=compress) File "/home/simone/firedrake/src/ufl/ufl/corealg/map_dag.py", line 84, in map_expr_dags r = handlers[v._ufl_typecode_](v) File "/home/simone/firedrake/src/ufl/ufl/corealg/multifunction.py", line 42, in _memoized_handler r = handler(self, o) File "/home/simone/firedrake/src/ufl/ufl/algorithms/apply_function_pullbacks.py", line 256, in form_argument return apply_single_function_pullbacks(o) File "/home/simone/firedrake/src/ufl/ufl/algorithms/apply_function_pullbacks.py", line 151, in apply_single_function_pullbacks assert len(gsh) == 1 AssertionError I made a very simple mesh in Gmsh (attached) and the script is the linear elasticity using a mixed formulation (Arnold et al., 2007). The problem is a 2D rectangle with a vertical fault where I prescribe 1 m of slip. I attached the mesh (.msh file) and the script. If anyone is able to help, I would appreciate it. Thank you so much in advance. Have a wonderful day! Best regards, Simone -- *Simone Puel* Ph.D. student in Geological Sciences at the University of Texas at Austin Jackson School of Geosciences & Institute for Geophysics (UTIG) JGB 5.332 & ROC 2.116L spuel@utexas.edu +1 (737)-230-6811
Dear Simone On Tue, 25 Jun 2019 at 00:34, Simone Puel <spuel@utexas.edu> wrote:
Hi all,
My name is Simone Puel, a graduate student at the University of Texas at Austin. I think I have installed Firedrake correctly since I can run the example in the online tutorial. Thank you so much. I am trying to convert a script I made in FEniCS in Firedrake. I modified some syntax, but I have the following error:
[...]
assert len(gsh) == 1 AssertionError
I made a very simple mesh in Gmsh (attached) and the script is the linear elasticity using a mixed formulation (Arnold et al., 2007). The problem is a 2D rectangle with a vertical fault where I prescribe 1 m of slip. I attached the mesh (.msh file) and the script. If anyone is able to help, I would appreciate it. Thank you so much in advance. Have a wonderful day!
This is, I think, I missing feature in the way UFL applies pullback to certain elements. (I'm slightly surprised the fenics version works!). In any case, I think I have a fix. Can you try with these changes in UFL: https://github.com/firedrakeproject/ufl/pull/18 And in Firedrake: https://github.com/firedrakeproject/firedrake/pull/1466 Thanks, Lawrence
On Tue, Jun 25, 2019 at 4:52 AM Lawrence Mitchell <wence@gmx.li> wrote:
Dear Simone
On Tue, 25 Jun 2019 at 00:34, Simone Puel <spuel@utexas.edu> wrote:
Hi all,
My name is Simone Puel, a graduate student at the University of Texas at
Austin. I think I have installed Firedrake correctly since I can run the example in the online tutorial. Thank you so much.
I am trying to convert a script I made in FEniCS in Firedrake. I modified some syntax, but I have the following error:
[...]
assert len(gsh) == 1 AssertionError
I made a very simple mesh in Gmsh (attached) and the script is the
linear elasticity using a mixed formulation (Arnold et al., 2007). The problem is a 2D rectangle with a vertical fault where I prescribe 1 m of slip.
I attached the mesh (.msh file) and the script. If anyone is able to help, I would appreciate it. Thank you so much in advance. Have a wonderful day!
This is, I think, I missing feature in the way UFL applies pullback to certain elements. (I'm slightly surprised the fenics version works!).
This is my chance to ask a notational question. I recently put this in PETSc so I could see how everything worked, and it looks like the operation here is to put tabulated reference values of the basis or test functions into real space for a cell integral say. Isn't that operation the pushforward (adjoint of pullback). I use the pullback when I am mapping back user function evals to evaluate the action of dual basis vectors. Am I misunderstanding the terminology? Thanks, Matt
In any case, I think I have a fix. Can you try with these changes in UFL:
https://github.com/firedrakeproject/ufl/pull/18
And in Firedrake:
https://github.com/firedrakeproject/firedrake/pull/1466
Thanks,
Lawrence
_______________________________________________ firedrake mailing list firedrake@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/>
I *think* that the concept is that the physical cell is mapped to the reference cell, so you are pulling back reference values through that mapping to get physical values. From: <firedrake-bounces@imperial.ac.uk> on behalf of Matthew Knepley <knepley@gmail.com> Date: Tuesday, 25 June 2019 at 11:37 To: Lawrence Mitchell <wence@gmx.li> Cc: firedrake <firedrake@imperial.ac.uk> Subject: Re: [firedrake] Error to run a script converted from FEniCS On Tue, Jun 25, 2019 at 4:52 AM Lawrence Mitchell <wence@gmx.li<mailto:wence@gmx.li>> wrote: Dear Simone On Tue, 25 Jun 2019 at 00:34, Simone Puel <spuel@utexas.edu<mailto:spuel@utexas.edu>> wrote:
Hi all,
My name is Simone Puel, a graduate student at the University of Texas at Austin. I think I have installed Firedrake correctly since I can run the example in the online tutorial. Thank you so much. I am trying to convert a script I made in FEniCS in Firedrake. I modified some syntax, but I have the following error:
[...]
assert len(gsh) == 1 AssertionError
I made a very simple mesh in Gmsh (attached) and the script is the linear elasticity using a mixed formulation (Arnold et al., 2007). The problem is a 2D rectangle with a vertical fault where I prescribe 1 m of slip. I attached the mesh (.msh file) and the script. If anyone is able to help, I would appreciate it. Thank you so much in advance. Have a wonderful day!
This is, I think, I missing feature in the way UFL applies pullback to certain elements. (I'm slightly surprised the fenics version works!). This is my chance to ask a notational question. I recently put this in PETSc so I could see how everything worked, and it looks like the operation here is to put tabulated reference values of the basis or test functions into real space for a cell integral say. Isn't that operation the pushforward (adjoint of pullback). I use the pullback when I am mapping back user function evals to evaluate the action of dual basis vectors. Am I misunderstanding the terminology? Thanks, Matt In any case, I think I have a fix. Can you try with these changes in UFL: https://github.com/firedrakeproject/ufl/pull/18 And in Firedrake: https://github.com/firedrakeproject/firedrake/pull/1466 Thanks, Lawrence _______________________________________________ firedrake mailing list firedrake@imperial.ac.uk<mailto:firedrake@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/>
On Tue, Jun 25, 2019 at 6:44 AM Ham, David A <david.ham@imperial.ac.uk> wrote:
I **think** that the concept is that the physical cell is mapped to the reference cell, so you are pulling back reference values through that mapping to get physical values.
Ah. I named them the opposite way :) This is a good window into how my mind works. I thought since the Jacobian maps reference to real (and phi for that matter) that the default direction of any mapping should be ref-to-real. It seems like this is the overwhelming convention in large deformations. Thanks, Matt
*From: *<firedrake-bounces@imperial.ac.uk> on behalf of Matthew Knepley < knepley@gmail.com> *Date: *Tuesday, 25 June 2019 at 11:37 *To: *Lawrence Mitchell <wence@gmx.li> *Cc: *firedrake <firedrake@imperial.ac.uk> *Subject: *Re: [firedrake] Error to run a script converted from FEniCS
On Tue, Jun 25, 2019 at 4:52 AM Lawrence Mitchell <wence@gmx.li> wrote:
Dear Simone
On Tue, 25 Jun 2019 at 00:34, Simone Puel <spuel@utexas.edu> wrote:
Hi all,
My name is Simone Puel, a graduate student at the University of Texas at
Austin. I think I have installed Firedrake correctly since I can run the example in the online tutorial. Thank you so much.
I am trying to convert a script I made in FEniCS in Firedrake. I modified some syntax, but I have the following error:
[...]
assert len(gsh) == 1 AssertionError
I made a very simple mesh in Gmsh (attached) and the script is the
linear elasticity using a mixed formulation (Arnold et al., 2007). The problem is a 2D rectangle with a vertical fault where I prescribe 1 m of slip.
I attached the mesh (.msh file) and the script. If anyone is able to help, I would appreciate it. Thank you so much in advance. Have a wonderful day!
This is, I think, I missing feature in the way UFL applies pullback to certain elements. (I'm slightly surprised the fenics version works!).
This is my chance to ask a notational question. I recently put this in PETSc so I could see how everything worked,
and it looks like the operation here is to put tabulated reference values of the basis or test functions into real space
for a cell integral say. Isn't that operation the pushforward (adjoint of pullback). I use the pullback when I am mapping
back user function evals to evaluate the action of dual basis vectors. Am I misunderstanding the terminology?
Thanks,
Matt
In any case, I think I have a fix. Can you try with these changes in UFL:
https://github.com/firedrakeproject/ufl/pull/18
And in Firedrake:
https://github.com/firedrakeproject/firedrake/pull/1466
Thanks,
Lawrence
_______________________________________________ firedrake mailing list firedrake@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/>
On 25 Jun 2019, at 11:49, Matthew Knepley <knepley@gmail.com> wrote:
On Tue, Jun 25, 2019 at 6:44 AM Ham, David A <david.ham@imperial.ac.uk> wrote: I *think* that the concept is that the physical cell is mapped to the reference cell, so you are pulling back reference values through that mapping to get physical values.
Ah. I named them the opposite way :) This is a good window into how my mind works. I thought since the Jacobian maps reference to real (and phi for that matter) that the default direction of any mapping should be ref-to-real. It seems like this is the overwhelming convention in large deformations.
If you have a mapping \phi : \hat{K} \to K Then \phi^* : Alt^k(K) \to Alt^k(\hat{K}) is the pullback defined by \phi^* w(v_1, ..., v_k) = w(\phi v_1, ... \phi v_k), v_i \in \hat{K}, w \in Alt^k(K) So at least differential geometry says the pullback goes in the opposite direction to the mapping. Lawrence
On Tue, Jun 25, 2019 at 7:22 AM Lawrence Mitchell <wence@gmx.li> wrote:
On 25 Jun 2019, at 11:49, Matthew Knepley <knepley@gmail.com> wrote:
On Tue, Jun 25, 2019 at 6:44 AM Ham, David A <david.ham@imperial.ac.uk> wrote: I *think* that the concept is that the physical cell is mapped to the reference cell, so you are pulling back reference values through that mapping to get physical values.
Ah. I named them the opposite way :) This is a good window into how my mind works. I thought since the Jacobian maps reference to real (and phi for that matter) that the default direction of any mapping should be ref-to-real. It seems like this is the overwhelming convention in large deformations.
If you have a mapping
\phi : \hat{K} \to K
Then
\phi^* : Alt^k(K) \to Alt^k(\hat{K})
is the pullback defined by
\phi^* w(v_1, ..., v_k) = w(\phi v_1, ... \phi v_k), v_i \in \hat{K}, w \in Alt^k(K)
So at least differential geometry says the pullback goes in the opposite direction to the mapping.
Yes, exactly. So if the mapping is ref-to-real, then then pulllback should be real-to-ref, which is what I needed for dual basis application. I needed ref-to-real for cell integration (since tabulations are in ref), which I called pushforward. Thanks, Matt
Lawrence
-- 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/>
Dear Lawrence and Matt, Thank you so much for your emails. I am not sure I understood correctly. In my case BDM is used for the stress (tensor), DGv for displacement (vector) and DGs for rotation (vector). If I write something like this: BDM = TensorFunctionSpace(mesh, "BDM", 1) DGv = VectorFunctionSpace(mesh, "DG", 0) DGs = VectorFunctionSpace(mesh, "DG", 0) # Mixed function spaces Vh_element = MixedElement([TensorElement(BDM), VectorElement(DGv), VectorElement(DGs)]) Vh = FunctionSpace(mesh, Vh_element) I have the following error: Traceback (most recent call last): File "Elasticity2D_Firedrake.py", line 116, in <module> BDM = TensorFunctionSpace(mesh, "BDM", 1) File "/home/simone/firedrake/src/firedrake/firedrake/functionspace.py", line 213, in TensorFunctionSpace assert sub_element.value_shape() == () AssertionError Thank you so much, probably I am a little bit confused how to use these terms. Thanks in advance, have a wonderful day. Best regards, Simone On Tue, Jun 25, 2019 at 6:33 AM Matthew Knepley <knepley@gmail.com> wrote:
On Tue, Jun 25, 2019 at 7:22 AM Lawrence Mitchell <wence@gmx.li> wrote:
On 25 Jun 2019, at 11:49, Matthew Knepley <knepley@gmail.com> wrote:
On Tue, Jun 25, 2019 at 6:44 AM Ham, David A <david.ham@imperial.ac.uk> wrote: I *think* that the concept is that the physical cell is mapped to the reference cell, so you are pulling back reference values through that mapping to get physical values.
Ah. I named them the opposite way :) This is a good window into how my mind works. I thought since the Jacobian maps reference to real (and phi for that matter) that the default direction of any mapping should be ref-to-real. It seems like this is the overwhelming convention in large deformations.
If you have a mapping
\phi : \hat{K} \to K
Then
\phi^* : Alt^k(K) \to Alt^k(\hat{K})
is the pullback defined by
\phi^* w(v_1, ..., v_k) = w(\phi v_1, ... \phi v_k), v_i \in \hat{K}, w \in Alt^k(K)
So at least differential geometry says the pullback goes in the opposite direction to the mapping.
Yes, exactly. So if the mapping is ref-to-real, then then pulllback should be real-to-ref, which is what I needed for dual basis application. I needed ref-to-real for cell integration (since tabulations are in ref), which I called pushforward.
Thanks,
Matt
Lawrence
-- 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/> _______________________________________________ firedrake mailing list firedrake@imperial.ac.uk https://mailman.ic.ac.uk/mailman/listinfo/firedrake
-- *Simone Puel* Ph.D. student in Geological Sciences at the University of Texas at Austin Jackson School of Geosciences & Institute for Geophysics (UTIG) JGB 5.332 & ROC 2.116L spuel@utexas.edu +1 (737)-230-6811
Even though I use the approach similar to the mixed formulation in the online tutorial (Mixed formulation for the Poisson equation): BDM = VectorFunctionSpace(mesh, "BDM", 1) DGv = VectorFunctionSpace(mesh, "DG", 0) DGs = FunctionSpace(mesh, "DG", 0) Vh = BDM * DGv * DGs I have the same error as before: Traceback (most recent call last): File "/home/simone/firedrake/src/PyOP2/pyop2/caching.py", line 197, in __new__ return cls._cache_lookup(key) File "/home/simone/firedrake/src/firedrake/firedrake/tsfc_interface.py", line 53, in _cache_lookup return cls._cache.get(key) or cls._read_from_disk(key, comm) File "/home/simone/firedrake/src/firedrake/firedrake/tsfc_interface.py", line 74, in _read_from_disk raise KeyError("Object with key %s not found" % key) KeyError: 'Object with key aec71f5f7a30a3d0b4db90beac07de85 not found' During handling of the above exception, another exception occurred: Traceback (most recent call last): File "Elasticity2D_Firedrake.py", line 163, in <module> solve(lhs_varf == rhs_varf, sol, bcs=[bc1,bc2]) File "/home/simone/firedrake/src/firedrake/firedrake/solving.py", line 125, in solve _solve_varproblem(*args, **kwargs) File "/home/simone/firedrake/src/firedrake/firedrake/solving.py", line 152, in _solve_varproblem appctx=appctx) File "/home/simone/firedrake/src/firedrake/firedrake/variational_solver.py", line 335, in __init__ super(LinearVariationalSolver, self).__init__(*args, **kwargs) File "/home/simone/firedrake/src/firedrake/firedrake/variational_solver.py", line 182, in __init__ options_prefix=self.options_prefix) File "/home/simone/firedrake/src/firedrake/firedrake/solving_utils.py", line 131, in __init__ form_compiler_parameters=self.fcp) File "/home/simone/firedrake/src/firedrake/firedrake/assemble.py", line 158, in create_assembly_callable loops = tuple(loops) File "/home/simone/firedrake/src/firedrake/firedrake/assemble.py", line 228, in _assemble kernels = tsfc_interface.compile_form(f, "form", parameters=form_compiler_parameters, inverse=inverse) File "/home/simone/firedrake/src/firedrake/firedrake/tsfc_interface.py", line 212, in compile_form number_map, interface, coffee).kernels File "/home/simone/firedrake/src/PyOP2/pyop2/caching.py", line 199, in __new__ obj = make_obj() File "/home/simone/firedrake/src/PyOP2/pyop2/caching.py", line 189, in make_obj obj.__init__(*args, **kwargs) File "/home/simone/firedrake/src/firedrake/firedrake/tsfc_interface.py", line 121, in __init__ tree = tsfc_compile_form(form, prefix=name, parameters=parameters, interface=interface, coffee=coffee) File "/home/simone/firedrake/src/tsfc/tsfc/driver.py", line 57, in compile_form fd = ufl_utils.compute_form_data(form, complex_mode=complex_mode) File "/home/simone/firedrake/src/tsfc/tsfc/ufl_utils.py", line 59, in compute_form_data complex_mode=complex_mode File "/home/simone/firedrake/src/ufl/ufl/algorithms/compute_form_data.py", line 298, in compute_form_data form = apply_function_pullbacks(form) File "/home/simone/firedrake/src/ufl/ufl/algorithms/apply_function_pullbacks.py", line 267, in apply_function_pullbacks return map_integrand_dags(FunctionPullbackApplier(), expr) File "/home/simone/firedrake/src/ufl/ufl/algorithms/map_integrands.py", line 58, in map_integrand_dags form, only_integral_type) File "/home/simone/firedrake/src/ufl/ufl/algorithms/map_integrands.py", line 39, in map_integrands for itg in form.integrals()] File "/home/simone/firedrake/src/ufl/ufl/algorithms/map_integrands.py", line 39, in <listcomp> for itg in form.integrals()] File "/home/simone/firedrake/src/ufl/ufl/algorithms/map_integrands.py", line 46, in map_integrands return itg.reconstruct(function(itg.integrand())) File "/home/simone/firedrake/src/ufl/ufl/algorithms/map_integrands.py", line 57, in <lambda> return map_integrands(lambda expr: map_expr_dag(function, expr, compress), File "/home/simone/firedrake/src/ufl/ufl/corealg/map_dag.py", line 37, in map_expr_dag result, = map_expr_dags(function, [expression], compress=compress) File "/home/simone/firedrake/src/ufl/ufl/corealg/map_dag.py", line 84, in map_expr_dags r = handlers[v._ufl_typecode_](v) File "/home/simone/firedrake/src/ufl/ufl/corealg/multifunction.py", line 42, in _memoized_handler r = handler(self, o) File "/home/simone/firedrake/src/ufl/ufl/algorithms/apply_function_pullbacks.py", line 256, in form_argument return apply_single_function_pullbacks(o) File "/home/simone/firedrake/src/ufl/ufl/algorithms/apply_function_pullbacks.py", line 151, in apply_single_function_pullbacks assert len(gsh) == 1 AssertionError Thank you so much, have a wonderful day. Best regards, Simone On Tue, Jun 25, 2019 at 8:03 AM Simone Puel <spuel@utexas.edu> wrote:
Dear Lawrence and Matt,
Thank you so much for your emails. I am not sure I understood correctly. In my case BDM is used for the stress (tensor), DGv for displacement (vector) and DGs for rotation (vector). If I write something like this:
BDM = TensorFunctionSpace(mesh, "BDM", 1) DGv = VectorFunctionSpace(mesh, "DG", 0) DGs = VectorFunctionSpace(mesh, "DG", 0)
# Mixed function spaces Vh_element = MixedElement([TensorElement(BDM), VectorElement(DGv), VectorElement(DGs)]) Vh = FunctionSpace(mesh, Vh_element)
I have the following error:
Traceback (most recent call last): File "Elasticity2D_Firedrake.py", line 116, in <module> BDM = TensorFunctionSpace(mesh, "BDM", 1) File "/home/simone/firedrake/src/firedrake/firedrake/functionspace.py", line 213, in TensorFunctionSpace assert sub_element.value_shape() == () AssertionError
Thank you so much, probably I am a little bit confused how to use these terms. Thanks in advance, have a wonderful day.
Best regards, Simone
On Tue, Jun 25, 2019 at 6:33 AM Matthew Knepley <knepley@gmail.com> wrote:
On Tue, Jun 25, 2019 at 7:22 AM Lawrence Mitchell <wence@gmx.li> wrote:
On 25 Jun 2019, at 11:49, Matthew Knepley <knepley@gmail.com> wrote:
On Tue, Jun 25, 2019 at 6:44 AM Ham, David A <david.ham@imperial.ac.uk> wrote: I *think* that the concept is that the physical cell is mapped to the reference cell, so you are pulling back reference values through that mapping to get physical values.
Ah. I named them the opposite way :) This is a good window into how my mind works. I thought since the Jacobian maps reference to real (and phi for that matter) that the default direction of any mapping should be ref-to-real. It seems like this is the overwhelming convention in large deformations.
If you have a mapping
\phi : \hat{K} \to K
Then
\phi^* : Alt^k(K) \to Alt^k(\hat{K})
is the pullback defined by
\phi^* w(v_1, ..., v_k) = w(\phi v_1, ... \phi v_k), v_i \in \hat{K}, w \in Alt^k(K)
So at least differential geometry says the pullback goes in the opposite direction to the mapping.
Yes, exactly. So if the mapping is ref-to-real, then then pulllback should be real-to-ref, which is what I needed for dual basis application. I needed ref-to-real for cell integration (since tabulations are in ref), which I called pushforward.
Thanks,
Matt
Lawrence
-- 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/> _______________________________________________ firedrake mailing list firedrake@imperial.ac.uk https://mailman.ic.ac.uk/mailman/listinfo/firedrake
-- *Simone Puel* Ph.D. student in Geological Sciences at the University of Texas at Austin Jackson School of Geosciences & Institute for Geophysics (UTIG) JGB 5.332 & ROC 2.116L spuel@utexas.edu +1 (737)-230-6811
-- *Simone Puel* Ph.D. student in Geological Sciences at the University of Texas at Austin Jackson School of Geosciences & Institute for Geophysics (UTIG) JGB 5.332 & ROC 2.116L spuel@utexas.edu +1 (737)-230-6811
Dear Simone,
On 25 Jun 2019, at 14:12, Simone Puel <spuel@utexas.edu> wrote:
Even though I use the approach similar to the mixed formulation in the online tutorial (Mixed formulation for the Poisson equation):
BDM = VectorFunctionSpace(mesh, "BDM", 1) DGv = VectorFunctionSpace(mesh, "DG", 0) DGs = FunctionSpace(mesh, "DG", 0)
Vh = BDM * DGv * DGs
[...] (Matt's question was a bit of a distraction). The problem here is that UFL applies a mapping from physical space to reference space before we perform numerical integration on the cell. Presently, UFL doesn't have the correct rules to do this for the VectorFunctionSpace(mesh, "BDM", 1) part. I added these in a branch (you need changes in UFL and in Firedrake) here: https://github.com/firedrakeproject/ufl/pull/18 and here: https://github.com/firedrakeproject/firedrake/pull/1466 I think if you use these branches (hopefully merged to master soon), then your model should work.
On Tue, Jun 25, 2019 at 8:03 AM Simone Puel <spuel@utexas.edu> wrote: Dear Lawrence and Matt,
Thank you so much for your emails. I am not sure I understood correctly. In my case BDM is used for the stress (tensor), DGv for displacement (vector) and DGs for rotation (vector). If I write something like this:
BDM = TensorFunctionSpace(mesh, "BDM", 1)
Just as an aside, you probably (not sure) don't want to take a TensorFunctionSpace of a BDM element for mixed elasticity. Since BDM is intrinsically vector-valued, this gives you a rank-3 tensor for the stress, which I guess you don't want? Thanks, Lawrence
Dear Lawrence, Thank you so much for the quick response. So, the idea is to leave my code as it is and modify the following scripts according with the updates you made in the links: ufl/algorithms/apply_function_pullbacks.py <https://github.com/firedrakeproject/ufl/pull/18/commits/d1472db0cec2c18acd0dd60ca1dde74ddf1df329#diff-82f3f0009b45b48c9f4b7c301585b525> ufl/finiteelement/mixedelement.py <https://github.com/firedrakeproject/ufl/pull/18/commits/2da7945b1af897810f5a215d9e1987c735582319#diff-5b1178ddf91101e9f03ce2e0e808f0cd> firedrake/functionspace.py <https://github.com/firedrakeproject/firedrake/pull/1466/commits/58807112ef773f99a4f1a03d7b3c553d6d752153#diff-e27ab14973abcfe3b4d5f97066c3c058> Is it right? Thanks a lot for your help. Have a wonderful day. Best regards, Simone On Tue, Jun 25, 2019 at 8:28 AM Lawrence Mitchell <wence@gmx.li> wrote:
Dear Simone,
On 25 Jun 2019, at 14:12, Simone Puel <spuel@utexas.edu> wrote:
Even though I use the approach similar to the mixed formulation in the online tutorial (Mixed formulation for the Poisson equation):
BDM = VectorFunctionSpace(mesh, "BDM", 1) DGv = VectorFunctionSpace(mesh, "DG", 0) DGs = FunctionSpace(mesh, "DG", 0)
Vh = BDM * DGv * DGs
[...]
(Matt's question was a bit of a distraction).
The problem here is that UFL applies a mapping from physical space to reference space before we perform numerical integration on the cell. Presently, UFL doesn't have the correct rules to do this for the
VectorFunctionSpace(mesh, "BDM", 1)
part.
I added these in a branch (you need changes in UFL and in Firedrake)
here:
https://github.com/firedrakeproject/ufl/pull/18
and here:
https://github.com/firedrakeproject/firedrake/pull/1466
I think if you use these branches (hopefully merged to master soon), then your model should work.
On Tue, Jun 25, 2019 at 8:03 AM Simone Puel <spuel@utexas.edu> wrote: Dear Lawrence and Matt,
Thank you so much for your emails. I am not sure I understood correctly.
In my case BDM is used for the stress (tensor), DGv for displacement (vector) and DGs for rotation (vector).
If I write something like this:
BDM = TensorFunctionSpace(mesh, "BDM", 1)
Just as an aside, you probably (not sure) don't want to take a TensorFunctionSpace of a BDM element for mixed elasticity. Since BDM is intrinsically vector-valued, this gives you a rank-3 tensor for the stress, which I guess you don't want?
Thanks,
Lawrence
-- *Simone Puel* Ph.D. student in Geological Sciences at the University of Texas at Austin Jackson School of Geosciences & Institute for Geophysics (UTIG) JGB 5.332 & ROC 2.116L spuel@utexas.edu +1 (737)-230-6811
On 25 Jun 2019, at 14:39, Simone Puel <spuel@utexas.edu> wrote:
Dear Lawrence,
Thank you so much for the quick response. So, the idea is to leave my code as it is and modify the following scripts according with the updates you made in the links:
ufl/algorithms/apply_function_pullbacks.py ufl/finiteelement/mixedelement.py firedrake/functionspace.py
Is it right? Thanks a lot for your help. Have a wonderful day.
Yes. Although the easier thing to do is to just change branches. If you do: cd src/ufl git fetch git checkout wence/pullback-tensor-vector-valued cd ../firedrake git fetch git checkout wence/generalise-tensor-element And confirm that the script then runs. Since we will incorporate these changes in the master branches soon, it is best to switch back (so that further updates will continue to work). You can do that by changing into the relevant directories and doing: git checkout master If you would rather just wait, I can let you know when the changes are merged (at which point a firedrake-update will get them). Thanks, Lawrence
Dear Lawrence, Thanks so much! I have done as you suggested and now I have a different error: Traceback (most recent call last): File "Elasticity2D_Firedrake.py", line 163, in <module> solve(lhs_varf == rhs_varf, sol, bcs=[bc1,bc2]) File "/home/simone/firedrake/src/firedrake/firedrake/solving.py", line 125, in solve _solve_varproblem(*args, **kwargs) File "/home/simone/firedrake/src/firedrake/firedrake/solving.py", line 153, in _solve_varproblem solver.solve() File "/home/simone/firedrake/src/firedrake/firedrake/variational_solver.py", line 249, in solve dbc.apply(self._problem.u) File "</home/simone/firedrake/lib/python3.6/site-packages/decorator.py:decorator-gen-48>", line 2, in apply File "/home/simone/firedrake/src/PyOP2/pyop2/profiling.py", line 60, in wrapper return f(*args, **kwargs) File "/home/simone/firedrake/src/firedrake/firedrake/bcs.py", line 432, in apply r.assign(self.function_arg, subset=self.node_set) File "</home/simone/firedrake/lib/python3.6/site-packages/decorator.py:decorator-gen-40>", line 2, in assign File "/home/simone/firedrake/src/firedrake/firedrake/utils.py", line 61, in wrapper return f(*args, **kwargs) File "/home/simone/firedrake/src/firedrake/firedrake/function.py", line 377, in assign assemble_expressions.Assign(self, expr), subset) File "</home/simone/firedrake/lib/python3.6/site-packages/decorator.py:decorator-gen-45>", line 2, in evaluate_expression File "/home/simone/firedrake/src/firedrake/firedrake/utils.py", line 58, in wrapper return f(*args, **kwargs) File "/home/simone/firedrake/src/firedrake/firedrake/assemble_expressions.py", line 507, in evaluate_expression e, args, _ = ExpressionWalker().walk(tree) File "/home/simone/firedrake/src/firedrake/firedrake/assemble_expressions.py", line 346, in walk return (self.visit(expr), File "/home/simone/firedrake/src/ufl/ufl/algorithms/transformer.py", line 114, in visit r = h(o) File "/home/simone/firedrake/src/firedrake/firedrake/assemble_expressions.py", line 413, in operator operands = [operands[0], self.visit(operands[1])] File "/home/simone/firedrake/src/ufl/ufl/algorithms/transformer.py", line 114, in visit r = h(o) File "/home/simone/firedrake/src/firedrake/firedrake/assemble_expressions.py", line 417, in operator operands = list(map(self.visit, o.ufl_operands)) File "/home/simone/firedrake/src/ufl/ufl/algorithms/transformer.py", line 114, in visit r = h(o) File "/home/simone/firedrake/src/firedrake/firedrake/assemble_expressions.py", line 383, in coefficient raise ValueError("Constant has mismatched shape for expression function space") ValueError: Constant has mismatched shape for expression function space I kept the function spaces as the original: BDM = VectorFunctionSpace(mesh, "BDM", 1) DGv = VectorFunctionSpace(mesh, "DG", 0) DGs = FunctionSpace(mesh, "DG", 0) Vh = BDM * DGv * DGs Is this error due to the boundary conditions that I prescribed? u0 = Constant((0.0, 0.0)) # zero displacement left, right and bottom zero_stress = as_vector([Constant((0., 0.)), Constant((0., 0.))]) # zero stress at top Thank you so much for your help, I really appreciate it. Best regards, Simone On Tue, Jun 25, 2019 at 9:11 AM Lawrence Mitchell <wence@gmx.li> wrote:
On 25 Jun 2019, at 14:39, Simone Puel <spuel@utexas.edu> wrote:
Dear Lawrence,
Thank you so much for the quick response. So, the idea is to leave my code as it is and modify the following scripts according with the updates you made in the links:
ufl/algorithms/apply_function_pullbacks.py ufl/finiteelement/mixedelement.py firedrake/functionspace.py
Is it right? Thanks a lot for your help. Have a wonderful day.
Yes. Although the easier thing to do is to just change branches.
If you do:
cd src/ufl git fetch git checkout wence/pullback-tensor-vector-valued cd ../firedrake git fetch git checkout wence/generalise-tensor-element
And confirm that the script then runs.
Since we will incorporate these changes in the master branches soon, it is best to switch back (so that further updates will continue to work).
You can do that by changing into the relevant directories and doing:
git checkout master
If you would rather just wait, I can let you know when the changes are merged (at which point a firedrake-update will get them).
Thanks,
Lawrence
-- *Simone Puel* Ph.D. student in Geological Sciences at the University of Texas at Austin Jackson School of Geosciences & Institute for Geophysics (UTIG) JGB 5.332 & ROC 2.116L spuel@utexas.edu +1 (737)-230-6811
On 25 Jun 2019, at 15:25, Simone Puel <spuel@utexas.edu> wrote:
Dear Lawrence,
Thanks so much! I have done as you suggested and now I have a different error:
Traceback (most recent call last): File "Elasticity2D_Firedrake.py", line 163, in <module> solve(lhs_varf == rhs_varf, sol, bcs=[bc1,bc2]) File "/home/simone/firedrake/src/firedrake/firedrake/solving.py", line 125, in solve _solve_varproblem(*args, **kwargs) File "/home/simone/firedrake/src/firedrake/firedrake/solving.py", line 153, in _solve_varproblem solver.solve() File "/home/simone/firedrake/src/firedrake/firedrake/variational_solver.py", line 249, in solve dbc.apply(self._problem.u) File "</home/simone/firedrake/lib/python3.6/site-packages/decorator.py:decorator-gen-48>", line 2, in apply File "/home/simone/firedrake/src/PyOP2/pyop2/profiling.py", line 60, in wrapper return f(*args, **kwargs) File "/home/simone/firedrake/src/firedrake/firedrake/bcs.py", line 432, in apply r.assign(self.function_arg, subset=self.node_set) File "</home/simone/firedrake/lib/python3.6/site-packages/decorator.py:decorator-gen-40>", line 2, in assign File "/home/simone/firedrake/src/firedrake/firedrake/utils.py", line 61, in wrapper return f(*args, **kwargs) File "/home/simone/firedrake/src/firedrake/firedrake/function.py", line 377, in assign assemble_expressions.Assign(self, expr), subset) File "</home/simone/firedrake/lib/python3.6/site-packages/decorator.py:decorator-gen-45>", line 2, in evaluate_expression File "/home/simone/firedrake/src/firedrake/firedrake/utils.py", line 58, in wrapper return f(*args, **kwargs) File "/home/simone/firedrake/src/firedrake/firedrake/assemble_expressions.py", line 507, in evaluate_expression e, args, _ = ExpressionWalker().walk(tree) File "/home/simone/firedrake/src/firedrake/firedrake/assemble_expressions.py", line 346, in walk return (self.visit(expr), File "/home/simone/firedrake/src/ufl/ufl/algorithms/transformer.py", line 114, in visit r = h(o) File "/home/simone/firedrake/src/firedrake/firedrake/assemble_expressions.py", line 413, in operator operands = [operands[0], self.visit(operands[1])] File "/home/simone/firedrake/src/ufl/ufl/algorithms/transformer.py", line 114, in visit r = h(o) File "/home/simone/firedrake/src/firedrake/firedrake/assemble_expressions.py", line 417, in operator operands = list(map(self.visit, o.ufl_operands)) File "/home/simone/firedrake/src/ufl/ufl/algorithms/transformer.py", line 114, in visit r = h(o) File "/home/simone/firedrake/src/firedrake/firedrake/assemble_expressions.py", line 383, in coefficient raise ValueError("Constant has mismatched shape for expression function space") ValueError: Constant has mismatched shape for expression function space
I kept the function spaces as the original: BDM = VectorFunctionSpace(mesh, "BDM", 1) DGv = VectorFunctionSpace(mesh, "DG", 0) DGs = FunctionSpace(mesh, "DG", 0)
Vh = BDM * DGv * DGs
Is this error due to the boundary conditions that I prescribed? u0 = Constant((0.0, 0.0)) # zero displacement left, right and bottom zero_stress = as_vector([Constant((0., 0.)), Constant((0., 0.))]) # zero stress at top
Thank you so much for your help, I really appreciate it.
I think so. You probably want zero_stress = as_tensor([[0, 0], [0, 0]]) or: zero_stress = Constant([[0, 0], [0, 0]]) Thanks, Lawrence
Dear Lawrence, Thanks, I really appreciate it. Hope it is the last one. The solver doesn't converge. If I run same mesh and same script in FEniCS (with the changes you suggested, that is the boundary conditions), the linear solver converges immediately (0.13 s). In Firedrake I have the following error: Traceback (most recent call last): File "Elasticity2D_Firedrake.py", line 163, in <module> solve(lhs_varf == rhs_varf, sol, bcs=[bc1,bc2]) File "/home/simone/firedrake/src/firedrake/firedrake/solving.py", line 125, in solve _solve_varproblem(*args, **kwargs) File "/home/simone/firedrake/src/firedrake/firedrake/solving.py", line 153, in _solve_varproblem solver.solve() File "/home/simone/firedrake/src/firedrake/firedrake/variational_solver.py", line 267, in solve solving_utils.check_snes_convergence(self.snes) File "/home/simone/firedrake/src/firedrake/firedrake/solving_utils.py", line 44, in check_snes_convergence %s""" % (snes.getIterationNumber(), msg)) firedrake.exceptions.ConvergenceError: Nonlinear solve failed to converge after 0 nonlinear iterations. Reason: Inner linear solve failed to converge after 10000 iterations with reason: DIVERGED_MAX_IT Thank you so much! Best regards, Simone On Tue, Jun 25, 2019 at 9:28 AM Lawrence Mitchell <wence@gmx.li> wrote:
On 25 Jun 2019, at 15:25, Simone Puel <spuel@utexas.edu> wrote:
Dear Lawrence,
Thanks so much! I have done as you suggested and now I have a different error:
Traceback (most recent call last): File "Elasticity2D_Firedrake.py", line 163, in <module> solve(lhs_varf == rhs_varf, sol, bcs=[bc1,bc2]) File "/home/simone/firedrake/src/firedrake/firedrake/solving.py", line 125, in solve _solve_varproblem(*args, **kwargs) File "/home/simone/firedrake/src/firedrake/firedrake/solving.py", line 153, in _solve_varproblem solver.solve() File "/home/simone/firedrake/src/firedrake/firedrake/variational_solver.py", line 249, in solve dbc.apply(self._problem.u) File "</home/simone/firedrake/lib/python3.6/site-packages/decorator.py:decorator-gen-48>", line 2, in apply File "/home/simone/firedrake/src/PyOP2/pyop2/profiling.py", line 60, in wrapper return f(*args, **kwargs) File "/home/simone/firedrake/src/firedrake/firedrake/bcs.py", line 432, in apply r.assign(self.function_arg, subset=self.node_set) File "</home/simone/firedrake/lib/python3.6/site-packages/decorator.py:decorator-gen-40>", line 2, in assign File "/home/simone/firedrake/src/firedrake/firedrake/utils.py", line 61, in wrapper return f(*args, **kwargs) File "/home/simone/firedrake/src/firedrake/firedrake/function.py", line 377, in assign assemble_expressions.Assign(self, expr), subset) File "</home/simone/firedrake/lib/python3.6/site-packages/decorator.py:decorator-gen-45>", line 2, in evaluate_expression File "/home/simone/firedrake/src/firedrake/firedrake/utils.py", line 58, in wrapper return f(*args, **kwargs) File "/home/simone/firedrake/src/firedrake/firedrake/assemble_expressions.py", line 507, in evaluate_expression e, args, _ = ExpressionWalker().walk(tree) File "/home/simone/firedrake/src/firedrake/firedrake/assemble_expressions.py", line 346, in walk return (self.visit(expr), File "/home/simone/firedrake/src/ufl/ufl/algorithms/transformer.py", line 114, in visit r = h(o) File "/home/simone/firedrake/src/firedrake/firedrake/assemble_expressions.py", line 413, in operator operands = [operands[0], self.visit(operands[1])] File "/home/simone/firedrake/src/ufl/ufl/algorithms/transformer.py", line 114, in visit r = h(o) File "/home/simone/firedrake/src/firedrake/firedrake/assemble_expressions.py", line 417, in operator operands = list(map(self.visit, o.ufl_operands)) File "/home/simone/firedrake/src/ufl/ufl/algorithms/transformer.py", line 114, in visit r = h(o) File "/home/simone/firedrake/src/firedrake/firedrake/assemble_expressions.py", line 383, in coefficient raise ValueError("Constant has mismatched shape for expression function space") ValueError: Constant has mismatched shape for expression function space
I kept the function spaces as the original: BDM = VectorFunctionSpace(mesh, "BDM", 1) DGv = VectorFunctionSpace(mesh, "DG", 0) DGs = FunctionSpace(mesh, "DG", 0)
Vh = BDM * DGv * DGs
Is this error due to the boundary conditions that I prescribed? u0 = Constant((0.0, 0.0)) # zero displacement left, right and bottom zero_stress = as_vector([Constant((0., 0.)), Constant((0., 0.))]) # zero stress at top
Thank you so much for your help, I really appreciate it.
I think so. You probably want
zero_stress = as_tensor([[0, 0], [0, 0]])
or:
zero_stress = Constant([[0, 0], [0, 0]])
Thanks,
Lawrence
-- *Simone Puel* Ph.D. student in Geological Sciences at the University of Texas at Austin Jackson School of Geosciences & Institute for Geophysics (UTIG) JGB 5.332 & ROC 2.116L spuel@utexas.edu +1 (737)-230-6811
Dear Simone, aha, this is an easier one.
On 25 Jun 2019, at 15:44, Simone Puel <spuel@utexas.edu> wrote:
Dear Lawrence,
Thanks, I really appreciate it. Hope it is the last one. The solver doesn't converge. If I run same mesh and same script in FEniCS (with the changes you suggested, that is the boundary conditions), the linear solver converges immediately (0.13 s). In Firedrake I have the following error:
Traceback (most recent call last): File "Elasticity2D_Firedrake.py", line 163, in <module> solve(lhs_varf == rhs_varf, sol, bcs=[bc1,bc2]) File "/home/simone/firedrake/src/firedrake/firedrake/solving.py", line 125, in solve _solve_varproblem(*args, **kwargs) File "/home/simone/firedrake/src/firedrake/firedrake/solving.py", line 153, in _solve_varproblem solver.solve() File "/home/simone/firedrake/src/firedrake/firedrake/variational_solver.py", line 267, in solve solving_utils.check_snes_convergence(self.snes) File "/home/simone/firedrake/src/firedrake/firedrake/solving_utils.py", line 44, in check_snes_convergence %s""" % (snes.getIterationNumber(), msg)) firedrake.exceptions.ConvergenceError: Nonlinear solve failed to converge after 0 nonlinear iterations. Reason: Inner linear solve failed to converge after 10000 iterations with reason: DIVERGED_MAX_IT
By default, the FEniCS solver is a direct factorisation, whereas ours is not (and this mixed elliptic problem is rather hard for the default preconditioner). You can select a sparse LU factorisation with the following solver options: solve(lhs_varf == rhs_varf, sol, bcs=[bc1, bc2], solver_parameters={"ksp_type": "preonly", "pc_type": "lu", "pc_factor_mat_solver_type": "mumps"}) Aside (to the list), perhaps it is worth reconsidering our choice of defaults for these solver options. Thanks, Lawrence
Dear Lawrence, Thank you so much. Modifying as you suggest, I have the following error: Traceback (most recent call last): File "Elasticity2D_Firedrake.py", line 166, in <module> "pc_factor_mat_solver_type": "mumps"}) File "/home/simone/firedrake/src/firedrake/firedrake/solving.py", line 125, in solve _solve_varproblem(*args, **kwargs) File "/home/simone/firedrake/src/firedrake/firedrake/solving.py", line 153, in _solve_varproblem solver.solve() File "/home/simone/firedrake/src/firedrake/firedrake/variational_solver.py", line 264, in solve self.snes.solve(None, work) File "PETSc/SNES.pyx", line 555, in petsc4py.PETSc.SNES.solve petsc4py.PETSc.Error: error code 56 [0] SNESSolve() line 4560 in /tmp/pip-req-build-eqqddqg6/src/snes/interface/snes.c [0] SNESSolve_KSPONLY() line 48 in /tmp/pip-req-build-eqqddqg6/src/snes/impls/ksponly/ksponly.c [0] KSPSolve() line 725 in /tmp/pip-req-build-eqqddqg6/src/ksp/ksp/interface/itfunc.c [0] KSPSetUp() line 391 in /tmp/pip-req-build-eqqddqg6/src/ksp/ksp/interface/itfunc.c [0] PCSetUp() line 932 in /tmp/pip-req-build-eqqddqg6/src/ksp/pc/interface/precon.c [0] PCSetUp_LU() line 88 in /tmp/pip-req-build-eqqddqg6/src/ksp/pc/impls/factor/lu/lu.c [0] MatGetOrdering() line 250 in /tmp/pip-req-build-eqqddqg6/src/mat/order/sorder.c [0] MatGetOrdering_ND() line 17 in /tmp/pip-req-build-eqqddqg6/src/mat/order/spnd.c [0] No support for this operation for this object type [0] Cannot get rows for matrix type nest Could you explain me in few words which solver it is better to use for my kind of problem (elasticity in mixed formulation), please? I will really appreciate it. Thanks again for your help. Best regards, Simone On Tue, Jun 25, 2019 at 10:11 AM Lawrence Mitchell <wence@gmx.li> wrote:
Dear Simone,
aha, this is an easier one.
On 25 Jun 2019, at 15:44, Simone Puel <spuel@utexas.edu> wrote:
Dear Lawrence,
Thanks, I really appreciate it. Hope it is the last one. The solver doesn't converge. If I run same mesh and same script in FEniCS (with the changes you suggested, that is the boundary conditions), the linear solver converges immediately (0.13 s). In Firedrake I have the following error:
Traceback (most recent call last): File "Elasticity2D_Firedrake.py", line 163, in <module> solve(lhs_varf == rhs_varf, sol, bcs=[bc1,bc2]) File "/home/simone/firedrake/src/firedrake/firedrake/solving.py", line 125, in solve _solve_varproblem(*args, **kwargs) File "/home/simone/firedrake/src/firedrake/firedrake/solving.py", line 153, in _solve_varproblem solver.solve() File "/home/simone/firedrake/src/firedrake/firedrake/variational_solver.py", line 267, in solve solving_utils.check_snes_convergence(self.snes) File "/home/simone/firedrake/src/firedrake/firedrake/solving_utils.py", line 44, in check_snes_convergence %s""" % (snes.getIterationNumber(), msg)) firedrake.exceptions.ConvergenceError: Nonlinear solve failed to converge after 0 nonlinear iterations. Reason: Inner linear solve failed to converge after 10000 iterations with reason: DIVERGED_MAX_IT
By default, the FEniCS solver is a direct factorisation, whereas ours is not (and this mixed elliptic problem is rather hard for the default preconditioner).
You can select a sparse LU factorisation with the following solver options:
solve(lhs_varf == rhs_varf, sol, bcs=[bc1, bc2], solver_parameters={"ksp_type": "preonly", "pc_type": "lu", "pc_factor_mat_solver_type": "mumps"})
Aside (to the list), perhaps it is worth reconsidering our choice of defaults for these solver options.
Thanks,
Lawrence
-- *Simone Puel* Ph.D. student in Geological Sciences at the University of Texas at Austin Jackson School of Geosciences & Institute for Geophysics (UTIG) JGB 5.332 & ROC 2.116L spuel@utexas.edu +1 (737)-230-6811
Ah, sorry,
On 25 Jun 2019, at 16:22, Simone Puel <spuel@utexas.edu> wrote:
Dear Lawrence,
Thank you so much. Modifying as you suggest, I have the following error:
Traceback (most recent call last): File "Elasticity2D_Firedrake.py", line 166, in <module> "pc_factor_mat_solver_type": "mumps"})
Please also add "mat_type": "aij" to the solver_parameters dictionary.
File "/home/simone/firedrake/src/firedrake/firedrake/solving.py", line 125, in solve _solve_varproblem(*args, **kwargs) File "/home/simone/firedrake/src/firedrake/firedrake/solving.py", line 153, in _solve_varproblem solver.solve() File "/home/simone/firedrake/src/firedrake/firedrake/variational_solver.py", line 264, in solve self.snes.solve(None, work) File "PETSc/SNES.pyx", line 555, in petsc4py.PETSc.SNES.solve petsc4py.PETSc.Error: error code 56 [0] SNESSolve() line 4560 in /tmp/pip-req-build-eqqddqg6/src/snes/interface/snes.c [0] SNESSolve_KSPONLY() line 48 in /tmp/pip-req-build-eqqddqg6/src/snes/impls/ksponly/ksponly.c [0] KSPSolve() line 725 in /tmp/pip-req-build-eqqddqg6/src/ksp/ksp/interface/itfunc.c [0] KSPSetUp() line 391 in /tmp/pip-req-build-eqqddqg6/src/ksp/ksp/interface/itfunc.c [0] PCSetUp() line 932 in /tmp/pip-req-build-eqqddqg6/src/ksp/pc/interface/precon.c [0] PCSetUp_LU() line 88 in /tmp/pip-req-build-eqqddqg6/src/ksp/pc/impls/factor/lu/lu.c [0] MatGetOrdering() line 250 in /tmp/pip-req-build-eqqddqg6/src/mat/order/sorder.c [0] MatGetOrdering_ND() line 17 in /tmp/pip-req-build-eqqddqg6/src/mat/order/spnd.c [0] No support for this operation for this object type [0] Cannot get rows for matrix type nest
Could you explain me in few words which solver it is better to use for my kind of problem (elasticity in mixed formulation), please? I will really appreciate it. Thanks again for your help.
For small (up to a few million dofs in 2D) sparse factorisations are hard to beat. Let's get these working before moving on to the preconditioning. I suspect you will want a block-diagonal preconditioner, for example as developed in this recent paper of Chen and collaborators https://www.ams.org/journals/mcom/2018-87-312/S0025-5718-2017-03285-3/S0025-... If your target domain has a simple geometric structure, then I think we can do the first thing they propose (eq 1.2) using geometric multigrid and overlapping Schwarz smoothers. If you need complicated geometric domains, the proposed auxiliary space approach that uses only algebraic multigrid solvers will likely be better. But I don't know a lot about these problems, so it is likely that a more detailed literature search may provide other ideas. In any case, if you have particular preconditioners in mind, we can certainly help attempting to implement them. Thanks, Lawrence
Hi Simone, The missing solver option is “mat_type”: “aij”. The default Firedrake behaviour is to build a “nest” meaning that a matrix is built for each block of the mixed system. This is better for Schur-compliment type solvers. The direct solvers PETSc interfaces to can’t deal with nests so you have to tell Firedrake to build a monolithic matrix instead. PETSc calls such matrices “aij” For 2D systems of moderate size, it’s very hard to beat direct solvers. For 3D systems and very large 2D systems you want some sort of iterative solver with the right preconditioner. In your case this is likely to be some sort of Schur compliment with multigrid on the reduced system. However that’s a complication which you can and should avoid until your system gets big. Regards, David From: <firedrake-bounces@imperial.ac.uk> on behalf of Simone Puel <spuel@utexas.edu> Date: Tuesday, 25 June 2019 at 16:22 To: Lawrence Mitchell <wence@gmx.li> Cc: firedrake <firedrake@imperial.ac.uk> Subject: Re: [firedrake] Error to run a script converted from FEniCS Dear Lawrence, Thank you so much. Modifying as you suggest, I have the following error: Traceback (most recent call last): File "Elasticity2D_Firedrake.py", line 166, in <module> "pc_factor_mat_solver_type": "mumps"}) File "/home/simone/firedrake/src/firedrake/firedrake/solving.py", line 125, in solve _solve_varproblem(*args, **kwargs) File "/home/simone/firedrake/src/firedrake/firedrake/solving.py", line 153, in _solve_varproblem solver.solve() File "/home/simone/firedrake/src/firedrake/firedrake/variational_solver.py", line 264, in solve self.snes.solve(None, work) File "PETSc/SNES.pyx", line 555, in petsc4py.PETSc.SNES.solve petsc4py.PETSc.Error: error code 56 [0] SNESSolve() line 4560 in /tmp/pip-req-build-eqqddqg6/src/snes/interface/snes.c [0] SNESSolve_KSPONLY() line 48 in /tmp/pip-req-build-eqqddqg6/src/snes/impls/ksponly/ksponly.c [0] KSPSolve() line 725 in /tmp/pip-req-build-eqqddqg6/src/ksp/ksp/interface/itfunc.c [0] KSPSetUp() line 391 in /tmp/pip-req-build-eqqddqg6/src/ksp/ksp/interface/itfunc.c [0] PCSetUp() line 932 in /tmp/pip-req-build-eqqddqg6/src/ksp/pc/interface/precon.c [0] PCSetUp_LU() line 88 in /tmp/pip-req-build-eqqddqg6/src/ksp/pc/impls/factor/lu/lu.c [0] MatGetOrdering() line 250 in /tmp/pip-req-build-eqqddqg6/src/mat/order/sorder.c [0] MatGetOrdering_ND() line 17 in /tmp/pip-req-build-eqqddqg6/src/mat/order/spnd.c [0] No support for this operation for this object type [0] Cannot get rows for matrix type nest Could you explain me in few words which solver it is better to use for my kind of problem (elasticity in mixed formulation), please? I will really appreciate it. Thanks again for your help. Best regards, Simone On Tue, Jun 25, 2019 at 10:11 AM Lawrence Mitchell <wence@gmx.li<mailto:wence@gmx.li>> wrote: Dear Simone, aha, this is an easier one.
On 25 Jun 2019, at 15:44, Simone Puel <spuel@utexas.edu<mailto:spuel@utexas.edu>> wrote:
Dear Lawrence,
Thanks, I really appreciate it. Hope it is the last one. The solver doesn't converge. If I run same mesh and same script in FEniCS (with the changes you suggested, that is the boundary conditions), the linear solver converges immediately (0.13 s). In Firedrake I have the following error:
Traceback (most recent call last): File "Elasticity2D_Firedrake.py", line 163, in <module> solve(lhs_varf == rhs_varf, sol, bcs=[bc1,bc2]) File "/home/simone/firedrake/src/firedrake/firedrake/solving.py", line 125, in solve _solve_varproblem(*args, **kwargs) File "/home/simone/firedrake/src/firedrake/firedrake/solving.py", line 153, in _solve_varproblem solver.solve() File "/home/simone/firedrake/src/firedrake/firedrake/variational_solver.py", line 267, in solve solving_utils.check_snes_convergence(self.snes) File "/home/simone/firedrake/src/firedrake/firedrake/solving_utils.py", line 44, in check_snes_convergence %s""" % (snes.getIterationNumber(), msg)) firedrake.exceptions.ConvergenceError: Nonlinear solve failed to converge after 0 nonlinear iterations. Reason: Inner linear solve failed to converge after 10000 iterations with reason: DIVERGED_MAX_IT
By default, the FEniCS solver is a direct factorisation, whereas ours is not (and this mixed elliptic problem is rather hard for the default preconditioner). You can select a sparse LU factorisation with the following solver options: solve(lhs_varf == rhs_varf, sol, bcs=[bc1, bc2], solver_parameters={"ksp_type": "preonly", "pc_type": "lu", "pc_factor_mat_solver_type": "mumps"}) Aside (to the list), perhaps it is worth reconsidering our choice of defaults for these solver options. Thanks, Lawrence -- Simone Puel Ph.D. student in Geological Sciences at the University of Texas at Austin Jackson School of Geosciences & Institute for Geophysics (UTIG) JGB 5.332 & ROC 2.116L spuel@utexas.edu<mailto:spuel@utexas.edu> +1 (737)-230-6811
Dear Lawrence and David, Thank you so much for your explanation and helps, I really appreciate it. I would like to switch from FEniCS to Firedrake and I am trying to rewrite my scripts according to the Firedrake syntax and logic. This should be a very simple problem (2D with less than 1000 dofs). However, adding the line you suggested in the solver_parameters as: solve(lhs_varf == rhs_varf, sol, bcs=[bc1, bc2], solver_parameters={"ksp_type": "preonly", "pc_type": "lu", "pc_factor_mat_solver_type": "mumps", "mat_type": "aij"}) The solver does not still converge giving the following error: Traceback (most recent call last): File "Elasticity2D_Firedrake.py", line 167, in <module> "mat_type": "aij"}) File "/home/simone/firedrake/src/firedrake/firedrake/solving.py", line 125, in solve _solve_varproblem(*args, **kwargs) File "/home/simone/firedrake/src/firedrake/firedrake/solving.py", line 153, in _solve_varproblem solver.solve() File "/home/simone/firedrake/src/firedrake/firedrake/variational_solver.py", line 267, in solve solving_utils.check_snes_convergence(self.snes) File "/home/simone/firedrake/src/firedrake/firedrake/solving_utils.py", line 44, in check_snes_convergence %s""" % (snes.getIterationNumber(), msg)) firedrake.exceptions.ConvergenceError: Nonlinear solve failed to converge after 0 nonlinear iterations. Reason: Inner linear solve failed to converge after 0 iterations with reason: DIVERGED_PCSETUP_FAILED Thank you so much! Best regards, Simone On Tue, Jun 25, 2019 at 10:36 AM Ham, David A <david.ham@imperial.ac.uk> wrote:
Hi Simone,
The missing solver option is “mat_type”: “aij”. The default Firedrake behaviour is to build a “nest” meaning that a matrix is built for each block of the mixed system. This is better for Schur-compliment type solvers. The direct solvers PETSc interfaces to can’t deal with nests so you have to tell Firedrake to build a monolithic matrix instead. PETSc calls such matrices “aij”
For 2D systems of moderate size, it’s very hard to beat direct solvers. For 3D systems and very large 2D systems you want some sort of iterative solver with the right preconditioner. In your case this is likely to be some sort of Schur compliment with multigrid on the reduced system. However that’s a complication which you can and should avoid until your system gets big.
Regards,
David
*From: *<firedrake-bounces@imperial.ac.uk> on behalf of Simone Puel < spuel@utexas.edu> *Date: *Tuesday, 25 June 2019 at 16:22 *To: *Lawrence Mitchell <wence@gmx.li> *Cc: *firedrake <firedrake@imperial.ac.uk> *Subject: *Re: [firedrake] Error to run a script converted from FEniCS
Dear Lawrence,
Thank you so much. Modifying as you suggest, I have the following error:
Traceback (most recent call last): File "Elasticity2D_Firedrake.py", line 166, in <module> "pc_factor_mat_solver_type": "mumps"}) File "/home/simone/firedrake/src/firedrake/firedrake/solving.py", line 125, in solve _solve_varproblem(*args, **kwargs) File "/home/simone/firedrake/src/firedrake/firedrake/solving.py", line 153, in _solve_varproblem solver.solve() File "/home/simone/firedrake/src/firedrake/firedrake/variational_solver.py", line 264, in solve self.snes.solve(None, work) File "PETSc/SNES.pyx", line 555, in petsc4py.PETSc.SNES.solve petsc4py.PETSc.Error: error code 56 [0] SNESSolve() line 4560 in /tmp/pip-req-build-eqqddqg6/src/snes/interface/snes.c [0] SNESSolve_KSPONLY() line 48 in /tmp/pip-req-build-eqqddqg6/src/snes/impls/ksponly/ksponly.c [0] KSPSolve() line 725 in /tmp/pip-req-build-eqqddqg6/src/ksp/ksp/interface/itfunc.c [0] KSPSetUp() line 391 in /tmp/pip-req-build-eqqddqg6/src/ksp/ksp/interface/itfunc.c [0] PCSetUp() line 932 in /tmp/pip-req-build-eqqddqg6/src/ksp/pc/interface/precon.c [0] PCSetUp_LU() line 88 in /tmp/pip-req-build-eqqddqg6/src/ksp/pc/impls/factor/lu/lu.c [0] MatGetOrdering() line 250 in /tmp/pip-req-build-eqqddqg6/src/mat/order/sorder.c [0] MatGetOrdering_ND() line 17 in /tmp/pip-req-build-eqqddqg6/src/mat/order/spnd.c [0] No support for this operation for this object type [0] Cannot get rows for matrix type nest
Could you explain me in few words which solver it is better to use for my kind of problem (elasticity in mixed formulation), please? I will really appreciate it.
Thanks again for your help.
Best regards,
Simone
On Tue, Jun 25, 2019 at 10:11 AM Lawrence Mitchell <wence@gmx.li> wrote:
Dear Simone,
aha, this is an easier one.
On 25 Jun 2019, at 15:44, Simone Puel <spuel@utexas.edu> wrote:
Dear Lawrence,
Thanks, I really appreciate it. Hope it is the last one. The solver doesn't converge. If I run same mesh and same script in FEniCS (with the changes you suggested, that is the boundary conditions), the linear solver converges immediately (0.13 s). In Firedrake I have the following error:
Traceback (most recent call last): File "Elasticity2D_Firedrake.py", line 163, in <module> solve(lhs_varf == rhs_varf, sol, bcs=[bc1,bc2]) File "/home/simone/firedrake/src/firedrake/firedrake/solving.py", line 125, in solve _solve_varproblem(*args, **kwargs) File "/home/simone/firedrake/src/firedrake/firedrake/solving.py", line 153, in _solve_varproblem solver.solve() File "/home/simone/firedrake/src/firedrake/firedrake/variational_solver.py", line 267, in solve solving_utils.check_snes_convergence(self.snes) File "/home/simone/firedrake/src/firedrake/firedrake/solving_utils.py", line 44, in check_snes_convergence %s""" % (snes.getIterationNumber(), msg)) firedrake.exceptions.ConvergenceError: Nonlinear solve failed to converge after 0 nonlinear iterations. Reason: Inner linear solve failed to converge after 10000 iterations with reason: DIVERGED_MAX_IT
By default, the FEniCS solver is a direct factorisation, whereas ours is not (and this mixed elliptic problem is rather hard for the default preconditioner).
You can select a sparse LU factorisation with the following solver options:
solve(lhs_varf == rhs_varf, sol, bcs=[bc1, bc2], solver_parameters={"ksp_type": "preonly", "pc_type": "lu", "pc_factor_mat_solver_type": "mumps"})
Aside (to the list), perhaps it is worth reconsidering our choice of defaults for these solver options.
Thanks,
Lawrence
--
*Simone Puel*
Ph.D. student in Geological Sciences at the University of Texas at Austin
Jackson School of Geosciences & Institute for Geophysics (UTIG)
JGB 5.332 & ROC 2.116L
spuel@utexas.edu
+1 (737)-230-6811
-- *Simone Puel* Ph.D. student in Geological Sciences at the University of Texas at Austin Jackson School of Geosciences & Institute for Geophysics (UTIG) JGB 5.332 & ROC 2.116L spuel@utexas.edu +1 (737)-230-6811
On Tue, Jun 25, 2019 at 11:48 AM Simone Puel <spuel@utexas.edu> wrote:
Dear Lawrence and David,
Thank you so much for your explanation and helps, I really appreciate it. I would like to switch from FEniCS to Firedrake and I am trying to rewrite my scripts according to the Firedrake syntax and logic. This should be a very simple problem (2D with less than 1000 dofs). However, adding the line you suggested in the solver_parameters as: solve(lhs_varf == rhs_varf, sol, bcs=[bc1, bc2], solver_parameters={"ksp_type": "preonly", "pc_type": "lu", "pc_factor_mat_solver_type": "mumps", "mat_type": "aij"})
The solver does not still converge giving the following error:
Traceback (most recent call last): File "Elasticity2D_Firedrake.py", line 167, in <module> "mat_type": "aij"}) File "/home/simone/firedrake/src/firedrake/firedrake/solving.py", line 125, in solve _solve_varproblem(*args, **kwargs) File "/home/simone/firedrake/src/firedrake/firedrake/solving.py", line 153, in _solve_varproblem solver.solve() File "/home/simone/firedrake/src/firedrake/firedrake/variational_solver.py", line 267, in solve solving_utils.check_snes_convergence(self.snes) File "/home/simone/firedrake/src/firedrake/firedrake/solving_utils.py", line 44, in check_snes_convergence %s""" % (snes.getIterationNumber(), msg)) firedrake.exceptions.ConvergenceError: Nonlinear solve failed to converge after 0 nonlinear iterations. Reason: Inner linear solve failed to converge after 0 iterations with reason: DIVERGED_PCSETUP_FAILED
Hi Simone, This means the factorization failed If this is elasticity, you should not have a nullspace. Are you sure the boundary conditions are correct? The other possibility is that the package you were using in FEniCS, maybe UMFPACK, is doing better pivoting than what you are using now, maybe PETSc. You can get UMFPACK by configuring PETSc with --download-suitesparse. Thanks, Matt
Thank you so much!
Best regards, Simone
On Tue, Jun 25, 2019 at 10:36 AM Ham, David A <david.ham@imperial.ac.uk> wrote:
Hi Simone,
The missing solver option is “mat_type”: “aij”. The default Firedrake behaviour is to build a “nest” meaning that a matrix is built for each block of the mixed system. This is better for Schur-compliment type solvers. The direct solvers PETSc interfaces to can’t deal with nests so you have to tell Firedrake to build a monolithic matrix instead. PETSc calls such matrices “aij”
For 2D systems of moderate size, it’s very hard to beat direct solvers. For 3D systems and very large 2D systems you want some sort of iterative solver with the right preconditioner. In your case this is likely to be some sort of Schur compliment with multigrid on the reduced system. However that’s a complication which you can and should avoid until your system gets big.
Regards,
David
*From: *<firedrake-bounces@imperial.ac.uk> on behalf of Simone Puel < spuel@utexas.edu> *Date: *Tuesday, 25 June 2019 at 16:22 *To: *Lawrence Mitchell <wence@gmx.li> *Cc: *firedrake <firedrake@imperial.ac.uk> *Subject: *Re: [firedrake] Error to run a script converted from FEniCS
Dear Lawrence,
Thank you so much. Modifying as you suggest, I have the following error:
Traceback (most recent call last): File "Elasticity2D_Firedrake.py", line 166, in <module> "pc_factor_mat_solver_type": "mumps"}) File "/home/simone/firedrake/src/firedrake/firedrake/solving.py", line 125, in solve _solve_varproblem(*args, **kwargs) File "/home/simone/firedrake/src/firedrake/firedrake/solving.py", line 153, in _solve_varproblem solver.solve() File "/home/simone/firedrake/src/firedrake/firedrake/variational_solver.py", line 264, in solve self.snes.solve(None, work) File "PETSc/SNES.pyx", line 555, in petsc4py.PETSc.SNES.solve petsc4py.PETSc.Error: error code 56 [0] SNESSolve() line 4560 in /tmp/pip-req-build-eqqddqg6/src/snes/interface/snes.c [0] SNESSolve_KSPONLY() line 48 in /tmp/pip-req-build-eqqddqg6/src/snes/impls/ksponly/ksponly.c [0] KSPSolve() line 725 in /tmp/pip-req-build-eqqddqg6/src/ksp/ksp/interface/itfunc.c [0] KSPSetUp() line 391 in /tmp/pip-req-build-eqqddqg6/src/ksp/ksp/interface/itfunc.c [0] PCSetUp() line 932 in /tmp/pip-req-build-eqqddqg6/src/ksp/pc/interface/precon.c [0] PCSetUp_LU() line 88 in /tmp/pip-req-build-eqqddqg6/src/ksp/pc/impls/factor/lu/lu.c [0] MatGetOrdering() line 250 in /tmp/pip-req-build-eqqddqg6/src/mat/order/sorder.c [0] MatGetOrdering_ND() line 17 in /tmp/pip-req-build-eqqddqg6/src/mat/order/spnd.c [0] No support for this operation for this object type [0] Cannot get rows for matrix type nest
Could you explain me in few words which solver it is better to use for my kind of problem (elasticity in mixed formulation), please? I will really appreciate it.
Thanks again for your help.
Best regards,
Simone
On Tue, Jun 25, 2019 at 10:11 AM Lawrence Mitchell <wence@gmx.li> wrote:
Dear Simone,
aha, this is an easier one.
On 25 Jun 2019, at 15:44, Simone Puel <spuel@utexas.edu> wrote:
Dear Lawrence,
Thanks, I really appreciate it. Hope it is the last one. The solver doesn't converge. If I run same mesh and same script in FEniCS (with the changes you suggested, that is the boundary conditions), the linear solver converges immediately (0.13 s). In Firedrake I have the following error:
Traceback (most recent call last): File "Elasticity2D_Firedrake.py", line 163, in <module> solve(lhs_varf == rhs_varf, sol, bcs=[bc1,bc2]) File "/home/simone/firedrake/src/firedrake/firedrake/solving.py", line 125, in solve _solve_varproblem(*args, **kwargs) File "/home/simone/firedrake/src/firedrake/firedrake/solving.py", line 153, in _solve_varproblem solver.solve() File "/home/simone/firedrake/src/firedrake/firedrake/variational_solver.py", line 267, in solve solving_utils.check_snes_convergence(self.snes) File "/home/simone/firedrake/src/firedrake/firedrake/solving_utils.py", line 44, in check_snes_convergence %s""" % (snes.getIterationNumber(), msg)) firedrake.exceptions.ConvergenceError: Nonlinear solve failed to converge after 0 nonlinear iterations. Reason: Inner linear solve failed to converge after 10000 iterations with reason: DIVERGED_MAX_IT
By default, the FEniCS solver is a direct factorisation, whereas ours is not (and this mixed elliptic problem is rather hard for the default preconditioner).
You can select a sparse LU factorisation with the following solver options:
solve(lhs_varf == rhs_varf, sol, bcs=[bc1, bc2], solver_parameters={"ksp_type": "preonly", "pc_type": "lu", "pc_factor_mat_solver_type": "mumps"})
Aside (to the list), perhaps it is worth reconsidering our choice of defaults for these solver options.
Thanks,
Lawrence
--
*Simone Puel*
Ph.D. student in Geological Sciences at the University of Texas at Austin
Jackson School of Geosciences & Institute for Geophysics (UTIG)
JGB 5.332 & ROC 2.116L
spuel@utexas.edu
+1 (737)-230-6811
-- *Simone Puel* Ph.D. student in Geological Sciences at the University of Texas at Austin Jackson School of Geosciences & Institute for Geophysics (UTIG) JGB 5.332 & ROC 2.116L spuel@utexas.edu +1 (737)-230-6811 _______________________________________________ firedrake mailing list firedrake@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/>
Dear Matt, The only differences between my scripts in FEniCS and Firedrake are: FEniCS: mesh = Mesh('Mesh_NormalFault.xml') boundaries = MeshFunction("size_t", mesh, 'Mesh_NormalFault_facet_region.xml') subdomains = MeshFunction("size_t", mesh, 'Mesh_NormalFault_physical_region.xml') ds = Measure('ds')(domain=mesh, subdomain_data=boundaries) dS = Measure('dS')(domain=mesh, subdomain_data=boundaries) BDM = VectorFunctionSpace(mesh, "BDM", 1) DGv = VectorFunctionSpace(mesh, "DG", 0) DGs = FunctionSpace(mesh, "DG", 0) Vh_element = MixedElement([BDM.ufl_element(), DGv.ufl_element(), DGs.ufl_element()]) Vh = FunctionSpace(mesh, Vh_element) u0 = Expression(("0.0", "0.0"), degree=1) # Neumann BCs zero_stress = Expression((("0.0", "0.0"), ("0.0", "0.0")), degree=5) # Dirichlet BCs A, b = assemble_system(lhs_varf, rhs_varf, bcs=[bc1,bc2]) solve(A, sol.vector(), b) Firedrake: mesh = Mesh("Mesh_NormalFault.msh") BDM = VectorFunctionSpace(mesh, "BDM", 1) DGv = VectorFunctionSpace(mesh, "DG", 0) DGs = FunctionSpace(mesh, "DG", 0) Vh = BDM * DGv * DGs u0 = Constant([0, 0]) # Neumann BCs zero_stress = Constant([[0, 0], [0, 0]]) # Dirichlet BCs solve(lhs_varf == rhs_varf, sol, bcs=[bc1, bc2], solver_parameters={"ksp_type": "preonly", "pc_type": "lu", "pc_factor_mat_solver_type": "mumps", "mat_type": "aij"}) All the rest is the same. The mesh is the same, but FEniCS requires a .xml file instead of .msh file. The variational form is the same. Probably I am doing something wrong, but the FEniCS version converges in 0.13 s. Thank you so much in advance. Have a wonderful day. Best regards, Simone On Tue, Jun 25, 2019 at 1:21 PM Matthew Knepley <knepley@gmail.com> wrote:
On Tue, Jun 25, 2019 at 11:48 AM Simone Puel <spuel@utexas.edu> wrote:
Dear Lawrence and David,
Thank you so much for your explanation and helps, I really appreciate it. I would like to switch from FEniCS to Firedrake and I am trying to rewrite my scripts according to the Firedrake syntax and logic. This should be a very simple problem (2D with less than 1000 dofs). However, adding the line you suggested in the solver_parameters as: solve(lhs_varf == rhs_varf, sol, bcs=[bc1, bc2], solver_parameters={"ksp_type": "preonly", "pc_type": "lu", "pc_factor_mat_solver_type": "mumps", "mat_type": "aij"})
The solver does not still converge giving the following error:
Traceback (most recent call last): File "Elasticity2D_Firedrake.py", line 167, in <module> "mat_type": "aij"}) File "/home/simone/firedrake/src/firedrake/firedrake/solving.py", line 125, in solve _solve_varproblem(*args, **kwargs) File "/home/simone/firedrake/src/firedrake/firedrake/solving.py", line 153, in _solve_varproblem solver.solve() File "/home/simone/firedrake/src/firedrake/firedrake/variational_solver.py", line 267, in solve solving_utils.check_snes_convergence(self.snes) File "/home/simone/firedrake/src/firedrake/firedrake/solving_utils.py", line 44, in check_snes_convergence %s""" % (snes.getIterationNumber(), msg)) firedrake.exceptions.ConvergenceError: Nonlinear solve failed to converge after 0 nonlinear iterations. Reason: Inner linear solve failed to converge after 0 iterations with reason: DIVERGED_PCSETUP_FAILED
Hi Simone,
This means the factorization failed If this is elasticity, you should not have a nullspace. Are you sure the boundary conditions are correct? The other possibility is that the package you were using in FEniCS, maybe UMFPACK, is doing better pivoting than what you are using now, maybe PETSc. You can get UMFPACK by configuring PETSc with --download-suitesparse.
Thanks,
Matt
Thank you so much!
Best regards, Simone
On Tue, Jun 25, 2019 at 10:36 AM Ham, David A <david.ham@imperial.ac.uk> wrote:
Hi Simone,
The missing solver option is “mat_type”: “aij”. The default Firedrake behaviour is to build a “nest” meaning that a matrix is built for each block of the mixed system. This is better for Schur-compliment type solvers. The direct solvers PETSc interfaces to can’t deal with nests so you have to tell Firedrake to build a monolithic matrix instead. PETSc calls such matrices “aij”
For 2D systems of moderate size, it’s very hard to beat direct solvers. For 3D systems and very large 2D systems you want some sort of iterative solver with the right preconditioner. In your case this is likely to be some sort of Schur compliment with multigrid on the reduced system. However that’s a complication which you can and should avoid until your system gets big.
Regards,
David
*From: *<firedrake-bounces@imperial.ac.uk> on behalf of Simone Puel < spuel@utexas.edu> *Date: *Tuesday, 25 June 2019 at 16:22 *To: *Lawrence Mitchell <wence@gmx.li> *Cc: *firedrake <firedrake@imperial.ac.uk> *Subject: *Re: [firedrake] Error to run a script converted from FEniCS
Dear Lawrence,
Thank you so much. Modifying as you suggest, I have the following error:
Traceback (most recent call last): File "Elasticity2D_Firedrake.py", line 166, in <module> "pc_factor_mat_solver_type": "mumps"}) File "/home/simone/firedrake/src/firedrake/firedrake/solving.py", line 125, in solve _solve_varproblem(*args, **kwargs) File "/home/simone/firedrake/src/firedrake/firedrake/solving.py", line 153, in _solve_varproblem solver.solve() File "/home/simone/firedrake/src/firedrake/firedrake/variational_solver.py", line 264, in solve self.snes.solve(None, work) File "PETSc/SNES.pyx", line 555, in petsc4py.PETSc.SNES.solve petsc4py.PETSc.Error: error code 56 [0] SNESSolve() line 4560 in /tmp/pip-req-build-eqqddqg6/src/snes/interface/snes.c [0] SNESSolve_KSPONLY() line 48 in /tmp/pip-req-build-eqqddqg6/src/snes/impls/ksponly/ksponly.c [0] KSPSolve() line 725 in /tmp/pip-req-build-eqqddqg6/src/ksp/ksp/interface/itfunc.c [0] KSPSetUp() line 391 in /tmp/pip-req-build-eqqddqg6/src/ksp/ksp/interface/itfunc.c [0] PCSetUp() line 932 in /tmp/pip-req-build-eqqddqg6/src/ksp/pc/interface/precon.c [0] PCSetUp_LU() line 88 in /tmp/pip-req-build-eqqddqg6/src/ksp/pc/impls/factor/lu/lu.c [0] MatGetOrdering() line 250 in /tmp/pip-req-build-eqqddqg6/src/mat/order/sorder.c [0] MatGetOrdering_ND() line 17 in /tmp/pip-req-build-eqqddqg6/src/mat/order/spnd.c [0] No support for this operation for this object type [0] Cannot get rows for matrix type nest
Could you explain me in few words which solver it is better to use for my kind of problem (elasticity in mixed formulation), please? I will really appreciate it.
Thanks again for your help.
Best regards,
Simone
On Tue, Jun 25, 2019 at 10:11 AM Lawrence Mitchell <wence@gmx.li> wrote:
Dear Simone,
aha, this is an easier one.
On 25 Jun 2019, at 15:44, Simone Puel <spuel@utexas.edu> wrote:
Dear Lawrence,
Thanks, I really appreciate it. Hope it is the last one. The solver doesn't converge. If I run same mesh and same script in FEniCS (with the changes you suggested, that is the boundary conditions), the linear solver converges immediately (0.13 s). In Firedrake I have the following error:
Traceback (most recent call last): File "Elasticity2D_Firedrake.py", line 163, in <module> solve(lhs_varf == rhs_varf, sol, bcs=[bc1,bc2]) File "/home/simone/firedrake/src/firedrake/firedrake/solving.py", line 125, in solve _solve_varproblem(*args, **kwargs) File "/home/simone/firedrake/src/firedrake/firedrake/solving.py", line 153, in _solve_varproblem solver.solve() File "/home/simone/firedrake/src/firedrake/firedrake/variational_solver.py", line 267, in solve solving_utils.check_snes_convergence(self.snes) File "/home/simone/firedrake/src/firedrake/firedrake/solving_utils.py", line 44, in check_snes_convergence %s""" % (snes.getIterationNumber(), msg)) firedrake.exceptions.ConvergenceError: Nonlinear solve failed to converge after 0 nonlinear iterations. Reason: Inner linear solve failed to converge after 10000 iterations with reason: DIVERGED_MAX_IT
By default, the FEniCS solver is a direct factorisation, whereas ours is not (and this mixed elliptic problem is rather hard for the default preconditioner).
You can select a sparse LU factorisation with the following solver options:
solve(lhs_varf == rhs_varf, sol, bcs=[bc1, bc2], solver_parameters={"ksp_type": "preonly", "pc_type": "lu", "pc_factor_mat_solver_type": "mumps"})
Aside (to the list), perhaps it is worth reconsidering our choice of defaults for these solver options.
Thanks,
Lawrence
--
*Simone Puel*
Ph.D. student in Geological Sciences at the University of Texas at Austin
Jackson School of Geosciences & Institute for Geophysics (UTIG)
JGB 5.332 & ROC 2.116L
spuel@utexas.edu
+1 (737)-230-6811
-- *Simone Puel* Ph.D. student in Geological Sciences at the University of Texas at Austin Jackson School of Geosciences & Institute for Geophysics (UTIG) JGB 5.332 & ROC 2.116L spuel@utexas.edu +1 (737)-230-6811 _______________________________________________ firedrake mailing list firedrake@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/>
-- *Simone Puel* Ph.D. student in Geological Sciences at the University of Texas at Austin Jackson School of Geosciences & Institute for Geophysics (UTIG) JGB 5.332 & ROC 2.116L spuel@utexas.edu +1 (737)-230-6811
Dear Simone,
On 25 Jun 2019, at 19:20, Matthew Knepley <knepley@gmail.com> wrote:
On Tue, Jun 25, 2019 at 11:48 AM Simone Puel <spuel@utexas.edu> wrote: Dear Lawrence and David,
Thank you so much for your explanation and helps, I really appreciate it. I would like to switch from FEniCS to Firedrake and I am trying to rewrite my scripts according to the Firedrake syntax and logic. This should be a very simple problem (2D with less than 1000 dofs). However, adding the line you suggested in the solver_parameters as: solve(lhs_varf == rhs_varf, sol, bcs=[bc1, bc2], solver_parameters={"ksp_type": "preonly", "pc_type": "lu", "pc_factor_mat_solver_type": "mumps", "mat_type": "aij"})
The solver does not still converge giving the following error:
Traceback (most recent call last): File "Elasticity2D_Firedrake.py", line 167, in <module> "mat_type": "aij"}) File "/home/simone/firedrake/src/firedrake/firedrake/solving.py", line 125, in solve _solve_varproblem(*args, **kwargs) File "/home/simone/firedrake/src/firedrake/firedrake/solving.py", line 153, in _solve_varproblem solver.solve() File "/home/simone/firedrake/src/firedrake/firedrake/variational_solver.py", line 267, in solve solving_utils.check_snes_convergence(self.snes) File "/home/simone/firedrake/src/firedrake/firedrake/solving_utils.py", line 44, in check_snes_convergence %s""" % (snes.getIterationNumber(), msg)) firedrake.exceptions.ConvergenceError: Nonlinear solve failed to converge after 0 nonlinear iterations. Reason: Inner linear solve failed to converge after 0 iterations with reason: DIVERGED_PCSETUP_FAILED
Hi Simone,
This means the factorization failed If this is elasticity, you should not have a nullspace. Are you sure the boundary conditions are correct? The other possibility is that the package you were using in FEniCS, maybe UMFPACK, is doing better pivoting than what you are using now, maybe PETSc. You can get UMFPACK by configuring PETSc with --download-suitesparse.
Matt is correct. We are using MUMPS to do the factorisation, we can get a bit more information by adding: "ksp_converged_reason": None to the solver parameters dict. When I do that, I get the same error, but I also see: Linear firedrake_0_ solve did not converge due to DIVERGED_PC_FAILED iterations 0 PC_FAILED due to FACTOR_OUTMEMORY Aha! MUMPS has an "esoteric" approach to memory allocation, it guesses how much space it will need, and then never attempts to reallocate. In this case, its guess is too low, and we get an error. We can help it out by upping the safety factor it uses. The relevant option is: -mat_mumps_icntl_14 "ICNTL(14): percentage increase in the estimated working space" So, let's have a guess and add "mat_mumps_icntl_14": 200 to the solver parameters dictionary. Now it works! $ python Elasticity2D_Firedrake.py Linear firedrake_0_ solve converged due to CONVERGED_ITS iterations 1 Time elapsed 0.42232704162597656 seconds Finished! So my final solver options were: solver_parameters={"ksp_type": "preonly", "pc_type": "lu", "ksp_converged_reason": None, "pc_factor_mat_solver_type": "mumps", "mat_mumps_icntl_14": 200, "mat_type": "aij"} Other sparse direct solvers are also available. The default firedrake-install does not build petsc with them, but as Matt mentions, UMFPACK is often a good option, as well as superlu_dist. I need to look up the right magic incantations to install those though. Thanks, Lawrence
Dear Lawrence, Yes now it works! Thank you so much, you all are the best! Thanks for your help, I really appreciate it. Have a wonderful day! Best regards, Simone On Wed, Jun 26, 2019 at 2:48 AM Lawrence Mitchell <wencel@gmail.com> wrote:
Dear Simone,
On 25 Jun 2019, at 19:20, Matthew Knepley <knepley@gmail.com> wrote:
On Tue, Jun 25, 2019 at 11:48 AM Simone Puel <spuel@utexas.edu> wrote: Dear Lawrence and David,
Thank you so much for your explanation and helps, I really appreciate it. I would like to switch from FEniCS to Firedrake and I am trying to rewrite my scripts according to the Firedrake syntax and logic. This should be a very simple problem (2D with less than 1000 dofs). However, adding the line you suggested in the solver_parameters as: solve(lhs_varf == rhs_varf, sol, bcs=[bc1, bc2], solver_parameters={"ksp_type": "preonly", "pc_type": "lu", "pc_factor_mat_solver_type": "mumps", "mat_type": "aij"})
The solver does not still converge giving the following error:
Traceback (most recent call last): File "Elasticity2D_Firedrake.py", line 167, in <module> "mat_type": "aij"}) File "/home/simone/firedrake/src/firedrake/firedrake/solving.py", line 125, in solve _solve_varproblem(*args, **kwargs) File "/home/simone/firedrake/src/firedrake/firedrake/solving.py", line 153, in _solve_varproblem solver.solve() File "/home/simone/firedrake/src/firedrake/firedrake/variational_solver.py", line 267, in solve solving_utils.check_snes_convergence(self.snes) File "/home/simone/firedrake/src/firedrake/firedrake/solving_utils.py", line 44, in check_snes_convergence %s""" % (snes.getIterationNumber(), msg)) firedrake.exceptions.ConvergenceError: Nonlinear solve failed to converge after 0 nonlinear iterations. Reason: Inner linear solve failed to converge after 0 iterations with reason: DIVERGED_PCSETUP_FAILED
Hi Simone,
This means the factorization failed If this is elasticity, you should not have a nullspace. Are you sure the boundary conditions are correct? The other possibility is that the package you were using in FEniCS, maybe UMFPACK, is doing better pivoting than what you are using now, maybe PETSc. You can get UMFPACK by configuring PETSc with --download-suitesparse.
Matt is correct. We are using MUMPS to do the factorisation, we can get a bit more information by adding:
"ksp_converged_reason": None
to the solver parameters dict.
When I do that, I get the same error, but I also see:
Linear firedrake_0_ solve did not converge due to DIVERGED_PC_FAILED iterations 0 PC_FAILED due to FACTOR_OUTMEMORY
Aha! MUMPS has an "esoteric" approach to memory allocation, it guesses how much space it will need, and then never attempts to reallocate. In this case, its guess is too low, and we get an error. We can help it out by upping the safety factor it uses.
The relevant option is:
-mat_mumps_icntl_14 "ICNTL(14): percentage increase in the estimated working space"
So, let's have a guess and add
"mat_mumps_icntl_14": 200
to the solver parameters dictionary.
Now it works!
$ python Elasticity2D_Firedrake.py Linear firedrake_0_ solve converged due to CONVERGED_ITS iterations 1 Time elapsed 0.42232704162597656 seconds Finished!
So my final solver options were:
solver_parameters={"ksp_type": "preonly", "pc_type": "lu", "ksp_converged_reason": None, "pc_factor_mat_solver_type": "mumps", "mat_mumps_icntl_14": 200, "mat_type": "aij"}
Other sparse direct solvers are also available. The default firedrake-install does not build petsc with them, but as Matt mentions, UMFPACK is often a good option, as well as superlu_dist. I need to look up the right magic incantations to install those though.
Thanks,
Lawrence
-- *Simone Puel* Ph.D. student in Geological Sciences at the University of Texas at Austin Jackson School of Geosciences & Institute for Geophysics (UTIG) JGB 5.332 & ROC 2.116L spuel@utexas.edu +1 (737)-230-6811
participants (5)
-
Ham, David A
-
Lawrence Mitchell
-
Lawrence Mitchell
-
Matthew Knepley
-
Simone Puel