dolfin team mailing list archive

[Anders Logg <logg@xxxxxxxxx>]

I haven't looked much at the design of the nonlinear solver before,
but I've looked at it some now. I would suggest some minor changes.
Here are some questions/suggestions:

1. The main class should be NewtonSolver and it should solve
   a nonlinear system F(x) = 0. The class should be in src/kernel/nls.

2. The main function should be

   // Solve F(x) = 0
   uint NewtonSolver::solve(NonlinearFunction& F, Vector& x);

3. The class NewtonSolver should not know anything about PDE, bilinear
   forms, variational problems etc, it just solves the discrete
   system F(x) = 0.

4. All classes in src/kernel/nls only handle the system F(x) = 0,
   no PDEs etc. Looks like we only need NewtonSolver and
   NonlinerFunction right now.

5. Maybe NonlinearFunction::form() is not needed? NonlinearFunction
   should only need F() and J().

6. Why does NewtonSolver inherit from KrylovSolver?

7. Move all PDE stuff to src/kernel/pde. It's enough to have form
   stuff in src/kernel/form.

8. NonlinearPDE knows about forms etc and it uses NewtonSolver to
   solve the nonlinear PDE. I don't think newton_solver needs to
   be a member variable. It's probably enough to create a NewtonSolver
   in the NonlinearPDE::solve() function, same as a GMRES object is
   created in LinearPDE::solve().

9. NonlinearPDE needs to create an object of type NonlinearFunction
   that it can pass to NewtonSolver::solve(). This can be done by
   creating a sub class of NonlinearFunction that represents the
   nonlinear function for a given BilinearForm and LinearForm. The sub
   class can be called for example NonlinearPDEFunction and be put in
   a separate file in src/kernel/pde.

   The constructor of NonlinearPDEFunction takes a NonlinearPDE as
   an argument, so in NonlinearPDE::solve(), we would do something
   like

   NonlinearPDEFunction F(this);
   NewtonSolver newton;
   newton.solve(F, x);

So all in all, we need the following classes:

   NewtonSolver         - solves F(x) = 0
   NonlinearFunction    - represents the function F

   NonlinearPDE         - represents a nonlinear PDE
   NonlinearPDEFunction - represents F for a nonlinear PDE

/Anders


On Fri, Feb 24, 2006 at 02:29:36PM +0100, Garth N. Wells wrote:
> On Fri, 2006-02-24 at 14:06 +0100, Johan Jansson wrote:
> > On Fri, Feb 24, 2006 at 01:45:33PM +0100, Garth N. Wells wrote:
> > > Hi,
> > > 
> > > I'm working on integrating NonlinearPDE into PDE, and it's coming along
> > > well. I do have a question about the correntness/style of some code.
> > > I've sketched out the code below. The class NonlinearPDE has a member
> > > object of type NewtonSolver (called newton_solver). When asking
> > > newton_solver to solve the problem, I pass nonlinear_pde to it (because
> > > nonlinear_pde knows how to form it's Jacobian and RHS vector). Is this a
> > > problem? I could tidy it up by making a another class - is there a
> > > problem with embedded classes (I recall something to so with swig)?
> > > 
> > > Garth
> > > 
> > > class NonlinearPDE : public GenericPDE
> > > {
> > > public:
> > >     .
> > >     .
> > >     void solve(NonlinearPDE
> > > 
> > >     // Form Jacobian and RHS (x is current solution)
> > >     void form(Matrix& A, Vector& b, const Vector& x)
> > > 
> > >     void solve(Function &u) 
> > >     {     
> > >       newton_solver.solve(*this, u);
> > >     }
> > > 
> > > private:
> > > 
> > >   NewtonSolver newton_solver;
> > >  
> > > }
> > > 
> > > 
> > > 
> > > 
> > 
> > Hi!
> > 
> > What you are doing is perfectly all right. NonlinearPDE should also
> > inherit from NonlinearProblem if I understand correctly:
> > 
> > class NonlinearPDE : public GenericPDE, public NonlinearProblem
> > 
> > 
> 
> Exactly. I'm using the name NonlinearFunction at the moment, but
> NonlinearProblem is a better name. An even better name would not refer
> to "nonlinear". Any ideas for names for a base class for objects which
> simply know how to form a matrix and a vector for a given problem?
> 
> Garth 
> 
> 
> > 
> > SWIG has problems parsing nested class definitions:
> > 
> > class A
> > {
> >   ...
> > 
> >   class B
> >   {
> >     ...
> >   }
> > 
> >   ...
> > }
> > 
> > This is usually not a problem though, SWIG cannot see the interface of
> > B, but it can see the interface of A, and that's usually enough. We
> > shouldn't worry about this, it's a deficiency of SWIG which will
> > eventually be fixed.
> > 
> >   Johan
> 
> 

-- 
Anders Logg
Research Assistant Professor
Toyota Technological Institute at Chicago
http://www.tti-c.org/logg/