Restricted dof (part 2)

Asked by dbeacham on 2010-03-10

I've finally got round to calculating the dofs on a boundary from a mesh function, so that I can go about trying to solve a pde purely on the boundary. However, I can't quite work out how best to play about with the assembled system:

If I have:

a = ( .. )*ds(i)
L = ( .. )*ds(i)

A = assemble(a, exterior_facet_domains=meshfn)
b = assemble(L, exterior_facet_domains=meshfn)

then, as far as I can tell I have two basic options, but I'm not sure what is best, and can't get (2) to work:

1. Do something along the lines of A_sub = A.vector()[rows, cols] (where rows = cols, is the dofs on which I wish to solve), ditto for b, and then solve using numpy (maybe I can put it back into a GenericMatrix and use one of the dolfin solvers, but I've certainly changed the structure and don't know how to do it anyway).

2. Adjust A,b appropriately as:
A.ident(unused_dofs)
b[unused_dofs] = 0.0

But this method throws up an error from petsc about the object not being in the correct state.

Matrix of size 10201 x 10201 has 0 nonzero entries.
[0]PETSC ERROR: --------------------- Error Message ------------------------------------
[0]PETSC ERROR: Object is in wrong state!
[0]PETSC ERROR: Matrix is missing diagonal entry in row 10200!
[0]PETSC ERROR: ------------------------------------------------------------------------
[0]PETSC ERROR: Petsc Release Version 3.0.0, Patch 11, Mon Feb 1 11:01:51 CST 2010
[0]PETSC ERROR: See docs/changes/index.html for recent updates.
[0]PETSC ERROR: See docs/faq.html for hints about trouble shooting.
[0]PETSC ERROR: See docs/index.html for manual pages.
[0]PETSC ERROR: ------------------------------------------------------------------------
[0]PETSC ERROR: Unknown Name on a linux-gnu named lappytop-dos by dbeacham Wed Mar 10 17:49:18 2010
[0]PETSC ERROR: Libraries linked from /home/dbeacham/FEM/stable/build/lib
[0]PETSC ERROR: Configure run at Thu Feb 11 11:37:57 2010
[0]PETSC ERROR: Configure options COPTFLAGS=-O3 --with-debugging=0 --with-shared=1 --with-clanguage=cxx --download-umfpack=1 --prefix=/home/dbeacham/FEM/stable/build
[0]PETSC ERROR: ------------------------------------------------------------------------
[0]PETSC ERROR: MatZeroRows_SeqAIJ() line 1367 in src/mat/impls/aij/seq/aij.c
[0]PETSC ERROR: MatZeroRows() line 4856 in src/mat/interface/matrix.c
[0]PETSC ERROR: MatZeroRowsIS() line 4913 in src/mat/interface/matrix.c

Thanks for any help,
David

Question information

Language:
English Edit question
Status:
Answered
For:
DOLFIN Edit question
Assignee:
No assignee Edit question
Last query:
2010-03-10
Last reply:
2010-04-21
Anders Logg (logg) said : #1

You could try using the new function A.ident_zeros() which should fix
all the rows which are zero (inserting one on the diagonal).

--
Anders

On Wed, Mar 10, 2010 at 06:13:39PM -0000, dbeacham wrote:
> New question #103905 on DOLFIN:
> https://answers.launchpad.net/dolfin/+question/103905
>
> I've finally got round to calculating the dofs on a boundary from a mesh function, so that I can go about trying to solve a pde purely on the boundary. However, I can't quite work out how best to play about with the assembled system:
>
> If I have:
>
> a = ( .. )*ds(i)
> L = ( .. )*ds(i)
>
> A = assemble(a, exterior_facet_domains=meshfn)
> b = assemble(L, exterior_facet_domains=meshfn)
>
> then, as far as I can tell I have two basic options, but I'm not sure what is best, and can't get (2) to work:
>
> 1. Do something along the lines of A_sub = A.vector()[rows, cols] (where rows = cols, is the dofs on which I wish to solve), ditto for b, and then solve using numpy (maybe I can put it back into a GenericMatrix and use one of the dolfin solvers, but I've certainly changed the structure and don't know how to do it anyway).
>
> 2. Adjust A,b appropriately as:
> A.ident(unused_dofs)
> b[unused_dofs] = 0.0
>
> But this method throws up an error from petsc about the object not being in the correct state.
>
> Matrix of size 10201 x 10201 has 0 nonzero entries.
> [0]PETSC ERROR: --------------------- Error Message ------------------------------------
> [0]PETSC ERROR: Object is in wrong state!
> [0]PETSC ERROR: Matrix is missing diagonal entry in row 10200!
> [0]PETSC ERROR: ------------------------------------------------------------------------
> [0]PETSC ERROR: Petsc Release Version 3.0.0, Patch 11, Mon Feb 1 11:01:51 CST 2010
> [0]PETSC ERROR: See docs/changes/index.html for recent updates.
> [0]PETSC ERROR: See docs/faq.html for hints about trouble shooting.
> [0]PETSC ERROR: See docs/index.html for manual pages.
> [0]PETSC ERROR: ------------------------------------------------------------------------
> [0]PETSC ERROR: Unknown Name on a linux-gnu named lappytop-dos by dbeacham Wed Mar 10 17:49:18 2010
> [0]PETSC ERROR: Libraries linked from /home/dbeacham/FEM/stable/build/lib
> [0]PETSC ERROR: Configure run at Thu Feb 11 11:37:57 2010
> [0]PETSC ERROR: Configure options COPTFLAGS=-O3 --with-debugging=0 --with-shared=1 --with-clanguage=cxx --download-umfpack=1 --prefix=/home/dbeacham/FEM/stable/build
> [0]PETSC ERROR: ------------------------------------------------------------------------
> [0]PETSC ERROR: MatZeroRows_SeqAIJ() line 1367 in src/mat/impls/aij/seq/aij.c
> [0]PETSC ERROR: MatZeroRows() line 4856 in src/mat/interface/matrix.c
> [0]PETSC ERROR: MatZeroRowsIS() line 4913 in src/mat/interface/matrix.c
>
> Thanks for any help,
> David
>
> You received this question notification because you are a member of
> DOLFIN Team, which is an answer contact for DOLFIN.
>
> _______________________________________________
> Mailing list: https://launchpad.net/~dolfin
> Post to : <email address hidden>
> Unsubscribe : https://launchpad.net/~dolfin
> More help : https://help.launchpad.net/ListHelp

Anders Logg (logg) said : #2

Try again now (still using ident_zeros). I have just pushed a fix for ident_zeros.

Can you help with this problem?

Provide an answer of your own, or ask dbeacham for more information if necessary.

To post a message you must log in.