DOLFIN

NonlinearProblem class, assembly problems

Asked by Knut Erik Skare

Hi.

I am trying to implement an adaptive time integrator in FEniCS and have come across some problems I hope you can help shed some light on.

Since I am changing time step in (potentially) each iteration, I would like to do as much pre-assembly as possible before I start doing the time steps. I have read the 'avoiding assembly'-chapter in the FEniCS-book, but I am looking at a nonlinear problem, rather than a linear.

My plan is to assemble all the necessary matrices and then pass THEM to the class derived from NonlinearProblem, instead of having the class assemble the system each time.

So to my question. Why doesn't the following class work?

class MyProblem(NonlinearProblem):
    def __init__(self, a, L):
        NonlinearProblem.__init__(self)
        self.L = L
        self.a = a

    def F(self, b, x):
        b = self.L

    def J(self, A, x):
        A = self.a

The remaining relevant code is as follows:
    from dolfin import *
    mesh = Rectangle(0, 0, 2.5, 2.5, 50, 50)
    V = FunctionSpace(mesh, "Lagrange", 1)

    v = TestFunction(V)
    u = TrialFunction(V)
    U = Function(V)

    l = inner(grad(U),grad(v))*dx
    a = derivative(l,U,u)

    L = Vector()
    A = Matrix()

    assemble(l,tensor=L)
    assemble(a,tensor=A)

    solver = NewtonSolver()

    problem = MySolve(A,L)
    solver.solve(problem,U.vector())

The error I get is
python: /build/buildd/dolfin-0.9.10/dolfin/la/PETScLUSolver.cpp:293: void dolfin::PETScLUSolver::set_petsc_operators(): Assertion `A->mat()' failed.

It seems like A isn't updated correctly, possibly because I am referencing it wrong. But I am somewhat lost here.

Any ideas?

Thanks in advance.

Question information

Language:
English Edit question
Status:
Solved
For:
DOLFIN Edit question
Assignee:
No assignee Edit question
Solved by:
Johan Hake
Solved:
Last query:
Last reply:
Revision history for this message
Best Johan Hake (johan-hake) said : 		#1

On Saturday November 12 2011 06:15:40 Knut Erik wrote:
> New question #178501 on DOLFIN:
> https://answers.launchpad.net/dolfin/+question/178501
>
> Hi.
>
> I am trying to implement an adaptive time integrator in FEniCS and have
> come across some problems I hope you can help shed some light on.
>
> Since I am changing time step in (potentially) each iteration, I would like
> to do as much pre-assembly as possible before I start doing the time
> steps. I have read the 'avoiding assembly'-chapter in the FEniCS-book, but
> I am looking at a nonlinear problem, rather than a linear.
>
> My plan is to assemble all the necessary matrices and then pass THEM to the
> class derived from NonlinearProblem, instead of having the class assemble
> the system each time.
>
> So to my question. Why doesn't the following class work?

Because there are no operator=() in Python.

Try:

 class MyProblem(NonlinearProblem):
     def __init__(self, a, L):
         NonlinearProblem.__init__(self)
         self.L = L
         self.a = a

     def F(self, b, x):
         b[:] = self.L

     def J(self, A, x):
         A.assign(self.a)

then you will operate directly on the tensors (A, b) passed through the
callback methods.

Johan

> The remaining relevant code is as follows:
> from dolfin import *
> mesh = Rectangle(0, 0, 2.5, 2.5, 50, 50)
> V = FunctionSpace(mesh, "Lagrange", 1)
>
> v = TestFunction(V)
> u = TrialFunction(V)
> U = Function(V)
>
> l = inner(grad(U),grad(v))*dx
> a = derivative(l,U,u)
>
> L = Vector()
> A = Matrix()
>
> assemble(l,tensor=L)
> assemble(a,tensor=A)
>
> solver = NewtonSolver()
>
> problem = MySolve(A,L)
> solver.solve(problem,U.vector())
>
> The error I get is
> python: /build/buildd/dolfin-0.9.10/dolfin/la/PETScLUSolver.cpp:293: void
> dolfin::PETScLUSolver::set_petsc_operators(): Assertion `A->mat()' failed.
>
> It seems like A isn't updated correctly, possibly because I am referencing
> it wrong. But I am somewhat lost here.
>
> Any ideas?
>
> Thanks in advance.

Revision history for this message
Knut Erik Skare (knut-erik-skare) said : 		#2

Thanks Johan Hake, that solved my question.