PETSc Segmentation Violation
Dear all, I am trying to run a code which returns the following error, before even entering the time loop: [0]PETSC ERROR: ------------------------------------------------------------------------ [0]PETSC ERROR: Caught signal number 11 SEGV: Segmentation Violation, probably memory access out of range [0]PETSC ERROR: Try option -start_in_debugger or -on_error_attach_debugger [0]PETSC ERROR: or see http://www.mcs.anl.gov/petsc/documentation/faq.html#valgrind [0]PETSC ERROR: or try http://valgrind.org on GNU/linux and Apple Mac OS X to find memory corruption errors [0]PETSC ERROR: configure using --with-debugging=yes, recompile, link, and run [0]PETSC ERROR: to get more information on the crash. -------------------------------------------------------------------------- MPI_ABORT was invoked on rank 0 in communicator MPI_COMM_WORLD with errorcode 59. NOTE: invoking MPI_ABORT causes Open MPI to kill all MPI processes. You may or may not see output from other processes, depending on exactly when Open MPI kills them. -------------------------------------------------------------------------- I tracked down the problem and I found that it is caused by the following: Matrix = assemble(u*v*dx).M.handle where u is a trial and v is a test function. Can any of you help with what might be the issue here? Thank you, Anna. -- Dr Anna Kalogirou Research Fellow School of Mathematics University of Leeds http://www1.maths.leeds.ac.uk/~matak/
Anna, Can you send us this code so that we can reproduce this problem? Thanks, Justin On Thu, Jun 2, 2016 at 11:56 AM, Anna Kalogirou <a.kalogirou@leeds.ac.uk> wrote:
Dear all,
I am trying to run a code which returns the following error, before even entering the time loop:
[0]PETSC ERROR: ------------------------------------------------------------------------ [0]PETSC ERROR: Caught signal number 11 SEGV: Segmentation Violation, probably memory access out of range [0]PETSC ERROR: Try option -start_in_debugger or -on_error_attach_debugger [0]PETSC ERROR: or see http://www.mcs.anl.gov/petsc/documentation/faq.html#valgrind [0]PETSC ERROR: or try http://valgrind.org on GNU/linux and Apple Mac OS X to find memory corruption errors [0]PETSC ERROR: configure using --with-debugging=yes, recompile, link, and run [0]PETSC ERROR: to get more information on the crash. -------------------------------------------------------------------------- MPI_ABORT was invoked on rank 0 in communicator MPI_COMM_WORLD with errorcode 59.
NOTE: invoking MPI_ABORT causes Open MPI to kill all MPI processes. You may or may not see output from other processes, depending on exactly when Open MPI kills them. --------------------------------------------------------------------------
I tracked down the problem and I found that it is caused by the following: Matrix = assemble(u*v*dx).M.handle
where u is a trial and v is a test function.
Can any of you help with what might be the issue here?
Thank you,
Anna.
--
Dr Anna Kalogirou Research Fellow School of Mathematics University of Leeds
http://www1.maths.leeds.ac.uk/~matak/
_______________________________________________ firedrake mailing list firedrake@imperial.ac.uk https://mailman.ic.ac.uk/mailman/listinfo/firedrake
Good morning, The code can be found here: https://bitbucket.org/annakalog/buoy2d/src/d05fb55ab1b3b6d75486e5c4b9a268933... Thanks, Anna. On 02/06/16 20:01, Justin Chang wrote:
Anna,
Can you send us this code so that we can reproduce this problem?
Thanks, Justin
On Thu, Jun 2, 2016 at 11:56 AM, Anna Kalogirou <a.kalogirou@leeds.ac.uk <mailto:a.kalogirou@leeds.ac.uk>> wrote:
Dear all,
I am trying to run a code which returns the following error, before even entering the time loop:
[0]PETSC ERROR: ------------------------------------------------------------------------ [0]PETSC ERROR: Caught signal number 11 SEGV: Segmentation Violation, probably memory access out of range [0]PETSC ERROR: Try option -start_in_debugger or -on_error_attach_debugger [0]PETSC ERROR: or see http://www.mcs.anl.gov/petsc/documentation/faq.html#valgrind [0]PETSC ERROR: or try http://valgrind.org on GNU/linux and Apple Mac OS X to find memory corruption errors [0]PETSC ERROR: configure using --with-debugging=yes, recompile, link, and run [0]PETSC ERROR: to get more information on the crash. -------------------------------------------------------------------------- MPI_ABORT was invoked on rank 0 in communicator MPI_COMM_WORLD with errorcode 59.
NOTE: invoking MPI_ABORT causes Open MPI to kill all MPI processes. You may or may not see output from other processes, depending on exactly when Open MPI kills them. --------------------------------------------------------------------------
I tracked down the problem and I found that it is caused by the following: Matrix = assemble(u*v*dx).M.handle
where u is a trial and v is a test function.
Can any of you help with what might be the issue here?
Thank you,
Anna.
--
Dr Anna Kalogirou Research Fellow School of Mathematics University of Leeds
http://www1.maths.leeds.ac.uk/~matak/ <http://www1.maths.leeds.ac.uk/%7Ematak/>
_______________________________________________ firedrake mailing list firedrake@imperial.ac.uk <mailto:firedrake@imperial.ac.uk> https://mailman.ic.ac.uk/mailman/listinfo/firedrake
_______________________________________________ firedrake mailing list firedrake@imperial.ac.uk https://mailman.ic.ac.uk/mailman/listinfo/firedrake
Hi Anna,
On 3 Jun 2016, at 09:59, Anna Kalogirou <a.kalogirou@leeds.ac.uk> wrote:
Good morning,
The code can be found here: https://bitbucket.org/annakalog/buoy2d/src/d05fb55ab1b3b6d75486e5c4b9a268933...
Running the buoy-swe.py in this directory fails for me with compilation errors in your expression for mu_bar: diff --git a/Inequality constraint/parameters.py b/Inequality constraint/parameters.py index b38c0f6..f35b94b 100755 --- a/Inequality constraint/parameters.py +++ b/Inequality constraint/parameters.py @@ -85,7 +85,7 @@ def def_Hb(Hb, x, L, Lp, H0, rho, Mass, a, Nx): def def_lam_mu_bar(lambda_bar, mu_bar, x, L, Lp, H0, g, rho, Mass, a): lambda_bar.interpolate(Expression("-0.5*(1.0+copysign(1.0,x[0]-Lp))*g*(sqrt(2*M*tan(a)/rho)+tan(a)*(x[0]-L))", Lp=Lp, g=g, M=Mass, a=a, rho=rho, L=L)) - mu_bar.interpolate(Expression("0.5*(1.0-copysign(1.0,x[0]-Lp))*(sqrt(tan(a)*(Lp-x[0]))", Lp=Lp, a=a)) + mu_bar.interpolate(Expression("0.5*(1.0-copysign(1.0,x[0]-Lp))*(sqrt(tan(a)*(Lp-x[0])))", Lp=Lp, a=a)) return (lambda_bar, mu_bar); Then I get a problem because making the mixed matrix in solver_mu doesn't work right. This seems to fix things for me: diff --git a/Inequality constraint/solvers.py b/Inequality constraint/solvers.py index 789be96..8987bca 100755 --- a/Inequality constraint/solvers.py +++ b/Inequality constraint/solvers.py @@ -58,22 +58,24 @@ def solver_mu(eta0, phi0, Z0, W0, LM_0_5, mu, v, v1, dt, Hb, lam_bar, mu_bar, g, # Build the matrix free operator B class MatrixFreeB(object): - def __init__(self, A, Q1, Q2): + def __init__(self, A, Q1, Q2, scale): self.A = A self.Q1 = Q1 self.Q2 = Q2 + self.scale = scale # Compute y = alpha Q1 + y = Q2^T*Q1 + A*x def mult(self, mat, x, y): self.A.mult(x, y) alpha = self.Q2.dot(x) y.axpy(alpha, self.Q1) + y.scale(self.scale) # Define MatrixFreeB class operator B B = PETSc.Mat().create() B.setSizes(*A.getSizes()) B.setType(B.Type.PYTHON) - B.setPythonContext(MatrixFreeB(A, rho/Mass*Q, Q)) + B.setPythonContext(MatrixFreeB(A, rho/Mass*Q, Q, (dt/2))) B.setUp() ''' # Print LHS 'matrix' B to inspect its values @@ -84,7 +86,7 @@ def solver_mu(eta0, phi0, Z0, W0, LM_0_5, mu, v, v1, dt, Hb, lam_bar, mu_bar, g, print 'LHS matrix A + Q*Q^T \n------------\n' + str( dense_B ) + ' \n' ''' # Define combined matrix - mixed_matrix = PETSc.Mat().createNest([[(dt/2)*B, (2/dt)*B_mu], [B_mu, C_lam]]) + mixed_matrix = PETSc.Mat().createNest([[B, (2/dt)*B_mu], [B_mu, C_lam]]) # For small problems, we don't need a preconditioner at all mu_solver = PETSc.KSP().create() Then one of your forms in buoy-swe.py is wrong, since one of the terms contains a test function but the other doesn't (the RHS just before the solvers_SV call). I don't know what this should be though. After fixing that you'll probably find that the mu_solver now doesn't work, because the RHS that you feed in comes from V but the solver works on the space W. But I didn't get any of the segfaults you described! Assuming you have run firedrake-update successfully, can you run firedrake-clean and see if that removes those problems? Cheers, Lawrence
Hi Lawrence, Thanks for the suggestions. Indeed, it doesn’t give any segmentation faults when running it from my mac laptop. On my linux machine, I ran firedrake-update, firedrake-clean and make alltest without any problems, but it still complains about those segmentation faults. Do you have any idea what the problem on linux might be? Now regarding the code itself, I fixed a few mistakes and updated it on bitbucket. It is now running but it’s blowing up. I suspect is due to bad preconditioning (or lack of any). How can I set up the preconditioner for the mixed system via setting the correct PETSc options? Thanks, Anna. On 3 Jun 2016, at 10:22, Lawrence Mitchell <lawrence.mitchell@imperial.ac.uk> wrote:
Hi Anna,
On 3 Jun 2016, at 09:59, Anna Kalogirou <a.kalogirou@leeds.ac.uk> wrote:
Good morning,
The code can be found here: https://bitbucket.org/annakalog/buoy2d/src/d05fb55ab1b3b6d75486e5c4b9a268933...
Running the buoy-swe.py in this directory fails for me with compilation errors in your expression for mu_bar:
diff --git a/Inequality constraint/parameters.py b/Inequality constraint/parameters.py index b38c0f6..f35b94b 100755 --- a/Inequality constraint/parameters.py +++ b/Inequality constraint/parameters.py @@ -85,7 +85,7 @@ def def_Hb(Hb, x, L, Lp, H0, rho, Mass, a, Nx): def def_lam_mu_bar(lambda_bar, mu_bar, x, L, Lp, H0, g, rho, Mass, a):
lambda_bar.interpolate(Expression("-0.5*(1.0+copysign(1.0,x[0]-Lp))*g*(sqrt(2*M*tan(a)/rho)+tan(a)*(x[0]-L))", Lp=Lp, g=g, M=Mass, a=a, rho=rho, L=L)) - mu_bar.interpolate(Expression("0.5*(1.0-copysign(1.0,x[0]-Lp))*(sqrt(tan(a)*(Lp-x[0]))", Lp=Lp, a=a)) + mu_bar.interpolate(Expression("0.5*(1.0-copysign(1.0,x[0]-Lp))*(sqrt(tan(a)*(Lp-x[0])))", Lp=Lp, a=a))
return (lambda_bar, mu_bar);
Then I get a problem because making the mixed matrix in solver_mu doesn't work right. This seems to fix things for me:
diff --git a/Inequality constraint/solvers.py b/Inequality constraint/solvers.py index 789be96..8987bca 100755 --- a/Inequality constraint/solvers.py +++ b/Inequality constraint/solvers.py @@ -58,22 +58,24 @@ def solver_mu(eta0, phi0, Z0, W0, LM_0_5, mu, v, v1, dt, Hb, lam_bar, mu_bar, g,
# Build the matrix free operator B class MatrixFreeB(object): - def __init__(self, A, Q1, Q2): + def __init__(self, A, Q1, Q2, scale): self.A = A self.Q1 = Q1 self.Q2 = Q2 + self.scale = scale
# Compute y = alpha Q1 + y = Q2^T*Q1 + A*x def mult(self, mat, x, y): self.A.mult(x, y) alpha = self.Q2.dot(x) y.axpy(alpha, self.Q1) + y.scale(self.scale)
# Define MatrixFreeB class operator B B = PETSc.Mat().create() B.setSizes(*A.getSizes()) B.setType(B.Type.PYTHON) - B.setPythonContext(MatrixFreeB(A, rho/Mass*Q, Q)) + B.setPythonContext(MatrixFreeB(A, rho/Mass*Q, Q, (dt/2))) B.setUp() ''' # Print LHS 'matrix' B to inspect its values @@ -84,7 +86,7 @@ def solver_mu(eta0, phi0, Z0, W0, LM_0_5, mu, v, v1, dt, Hb, lam_bar, mu_bar, g, print 'LHS matrix A + Q*Q^T \n------------\n' + str( dense_B ) + ' \n' ''' # Define combined matrix - mixed_matrix = PETSc.Mat().createNest([[(dt/2)*B, (2/dt)*B_mu], [B_mu, C_lam]]) + mixed_matrix = PETSc.Mat().createNest([[B, (2/dt)*B_mu], [B_mu, C_lam]])
# For small problems, we don't need a preconditioner at all mu_solver = PETSc.KSP().create()
Then one of your forms in buoy-swe.py is wrong, since one of the terms contains a test function but the other doesn't (the RHS just before the solvers_SV call). I don't know what this should be though.
After fixing that you'll probably find that the mu_solver now doesn't work, because the RHS that you feed in comes from V but the solver works on the space W.
But I didn't get any of the segfaults you described! Assuming you have run firedrake-update successfully, can you run firedrake-clean and see if that removes those problems?
Cheers,
Lawrence _______________________________________________ firedrake mailing list firedrake@imperial.ac.uk https://mailman.ic.ac.uk/mailman/listinfo/firedrake
On 3 Jun 2016, at 16:21, Anna Kalogirou <A.Kalogirou@leeds.ac.uk> wrote:
Hi Lawrence,
Thanks for the suggestions. Indeed, it doesn’t give any segmentation faults when running it from my mac laptop. On my linux machine, I ran firedrake-update, firedrake-clean and make alltest without any problems, but it still complains about those segmentation faults. Do you have any idea what the problem on linux might be?
Not sure. I would run under gdb and see where the segfault occurs: gdb --args python buoy-swe.py (gdb) run ... wait for segfault (gdb) bt
Now regarding the code itself, I fixed a few mistakes and updated it on bitbucket. It is now running but it’s blowing up. I suspect is due to bad preconditioning (or lack of any). How can I set up the preconditioner for the mixed system via setting the correct PETSc options?
If you're blowing up, but the solves are converging, then preconditioning is not likely to be the problem. Instead the numerics might be screwy. Or, by blowing up, do you mean that the solvers don't converge? Lawrence
I attach the gdb outcome. It seems to complain about a missing file /tmp/pip-HUEqVU-build/src/vec/is/utils/isltog.c. Can someone understand what this file is doing? Thanks, Anna. On 03/06/16 16:31, Lawrence Mitchell wrote:
On 3 Jun 2016, at 16:21, Anna Kalogirou <A.Kalogirou@leeds.ac.uk> wrote:
Hi Lawrence,
Thanks for the suggestions. Indeed, it doesn’t give any segmentation faults when running it from my mac laptop. On my linux machine, I ran firedrake-update, firedrake-clean and make alltest without any problems, but it still complains about those segmentation faults. Do you have any idea what the problem on linux might be? Not sure. I would run under gdb and see where the segfault occurs:
gdb --args python buoy-swe.py (gdb) run ... wait for segfault (gdb) bt
Now regarding the code itself, I fixed a few mistakes and updated it on bitbucket. It is now running but it’s blowing up. I suspect is due to bad preconditioning (or lack of any). How can I set up the preconditioner for the mixed system via setting the correct PETSc options? If you're blowing up, but the solves are converging, then preconditioning is not likely to be the problem. Instead the numerics might be screwy. Or, by blowing up, do you mean that the solvers don't converge?
Lawrence
_______________________________________________ firedrake mailing list firedrake@imperial.ac.uk https://mailman.ic.ac.uk/mailman/listinfo/firedrake
On 3 Jun 2016, at 16:55, Anna Kalogirou <a.kalogirou@leeds.ac.uk> wrote:
I attach the gdb outcome. It seems to complain about a missing file /tmp/pip-HUEqVU-build/src/vec/is/utils/isltog.c. Can someone understand what this file is doing? #0 ISLocalToGlobalMappingApply (mapping=0xbe00a3e606158d48, N=2, ^^^^^^^^^^^^^^^^^^^
This looks like a very suspicious address. But I suspect that the problem is that somehow something has gone wrong well before this. This is happening when assembling a matrix. Can you tell which one it is? Can you narrow down the set of operations you need to do in your buoy setup such that you get the error? What has probably happened is that somewhere some memory has been overwritten, but it's difficult to guess where that might have occurred. Lawrence
Anna, I would be interested seeing /tmp/pyop2-cache-uid375330/bd937acb61afd8e72052253aad28bcd6*.c Can you send that file to us? Thanks, Miklós
On 3 Jun 2016, at 16:55, Anna Kalogirou <a.kalogirou@leeds.ac.uk> wrote:
I attach the gdb outcome. It seems to complain about a missing file /tmp/pip-HUEqVU-build/src/vec/is/utils/isltog.c. Can someone understand what this file is doing?
Thanks, Anna. On 03/06/16 16:31, Lawrence Mitchell wrote:
On 3 Jun 2016, at 16:21, Anna Kalogirou <A.Kalogirou@leeds.ac.uk> <mailto:A.Kalogirou@leeds.ac.uk> wrote:
Hi Lawrence,
Thanks for the suggestions. Indeed, it doesn’t give any segmentation faults when running it from my mac laptop. On my linux machine, I ran firedrake-update, firedrake-clean and make alltest without any problems, but it still complains about those segmentation faults. Do you have any idea what the problem on linux might be? Not sure. I would run under gdb and see where the segfault occurs:
gdb --args python buoy-swe.py (gdb) run ... wait for segfault (gdb) bt
Now regarding the code itself, I fixed a few mistakes and updated it on bitbucket. It is now running but it’s blowing up. I suspect is due to bad preconditioning (or lack of any). How can I set up the preconditioner for the mixed system via setting the correct PETSc options? If you're blowing up, but the solves are converging, then preconditioning is not likely to be the problem. Instead the numerics might be screwy. Or, by blowing up, do you mean that the solvers don't converge?
Lawrence
_______________________________________________ firedrake mailing list firedrake@imperial.ac.uk <mailto:firedrake@imperial.ac.uk> https://mailman.ic.ac.uk/mailman/listinfo/firedrake <https://mailman.ic.ac.uk/mailman/listinfo/firedrake> <Python segmentation fault.txt>_______________________________________________ firedrake mailing list firedrake@imperial.ac.uk https://mailman.ic.ac.uk/mailman/listinfo/firedrake
Hi Miklos, You can find the file attached. Thanks, Anna. On 04/06/16 09:08, Miklós Homolya wrote:
Anna,
I would be interested seeing /tmp/pyop2-cache-uid375330/bd937acb61afd8e72052253aad28bcd6*.c Can you send that file to us?
Thanks, Miklós
On 3 Jun 2016, at 16:55, Anna Kalogirou <a.kalogirou@leeds.ac.uk <mailto:a.kalogirou@leeds.ac.uk>> wrote:
I attach the gdb outcome. It seems to complain about a missing file /tmp/pip-HUEqVU-build/src/vec/is/utils/isltog.c. Can someone understand what this file is doing?
Thanks, Anna. On 03/06/16 16:31, Lawrence Mitchell wrote:
On 3 Jun 2016, at 16:21, Anna Kalogirou<A.Kalogirou@leeds.ac.uk> wrote:
Hi Lawrence,
Thanks for the suggestions. Indeed, it doesn’t give any segmentation faults when running it from my mac laptop. On my linux machine, I ran firedrake-update, firedrake-clean and make alltest without any problems, but it still complains about those segmentation faults. Do you have any idea what the problem on linux might be? Not sure. I would run under gdb and see where the segfault occurs:
gdb --args python buoy-swe.py (gdb) run ... wait for segfault (gdb) bt
Now regarding the code itself, I fixed a few mistakes and updated it on bitbucket. It is now running but it’s blowing up. I suspect is due to bad preconditioning (or lack of any). How can I set up the preconditioner for the mixed system via setting the correct PETSc options? If you're blowing up, but the solves are converging, then preconditioning is not likely to be the problem. Instead the numerics might be screwy. Or, by blowing up, do you mean that the solvers don't converge?
Lawrence
_______________________________________________ firedrake mailing list firedrake@imperial.ac.uk https://mailman.ic.ac.uk/mailman/listinfo/firedrake
<Python segmentation fault.txt>_______________________________________________ firedrake mailing list firedrake@imperial.ac.uk <mailto:firedrake@imperial.ac.uk> https://mailman.ic.ac.uk/mailman/listinfo/firedrake
_______________________________________________ firedrake mailing list firedrake@imperial.ac.uk https://mailman.ic.ac.uk/mailman/listinfo/firedrake
participants (5)
-
Anna Kalogirou
-
Anna Kalogirou
-
Justin Chang
-
Lawrence Mitchell
-
Miklós Homolya