On Tue, 14 Sep 2010 09:00:39 +0200 (CEST), Tim Kroeger 
<[email protected]> wrote:
> Well, I could make it a vector instead, but then, on the other hand, I 
> guess that in most situations the user does not have any idea of which 
> order is best to use.

That's what ordering packages are for, but defining the interface to use
std::set prohibits a smart user from providing an ordering.

> > So the LinearSolver doesn't take ownership of this set?  It seems
> > horrible to make the user responsible for that.
> 
> Well, users are supposed to do things like
> 
> std::set<unsigned int> a;
> // Add values here
> linear_solver.solve_only_on(&a);

I would expect the user to have a function that selects the subdomain
and registers it.  When the stack frame goes away, it's destructor would
be called, which is why I would normally want the solver object to
reference it.  But whatever.

> Well, what I would like to do, but don't know how, is to let subrhs be 
> a vector with the same layout as rhs and with any dof being owned by 
> that processor which owns the corresponding dof in rhs.  Or, perhaps 
> better, do an automatic load balancing for subrhs.  (Likewise for 
> subsolution of course.)

Each process provides a list of global indices that it wants to own in
the sub-problem.  The length of this list determines the local size.
You probably want

  VecCreate(comm,X);
  VecSetSizes(X,local_size,PETSC_DETERMINE);
  VecSetFromOptions(X);

which will create a serial vector in serial and a parallel vector in
parallel.  The user can use a partitioning algorithm if they want to
redistribute dofs (or you could add a runtime option to use an arbitrary
algorithm, PETSc offers MatPartitioning for this).

> Also, I am unsure whether or not I have to take care about the matrix 
> rows being owned by the same processors that own the corresponding 
> vector dofs.

It's all determined by the layout of the IS.  If you don't want rows to
"move", then the IS should only consist of global indices that are
presently owned.

> > If dofs == NULL, this function is a no-op?
> 
> Yes.  You see, this function is for all solver packages except PETSc, 
> and they should just give an error message when somebody attempts to 
> use solve_only_on().  But, of course, if a NULL pointer is used, which 
> means solve on the entire domain, there is no need to give an error 
> message, because this is what is done anyway, isn't it?  (-:

I wrote that comment first, before realizing that solve_only_on was just
restricting the domain and didn't solve anything.  Perhaps something
like set_active_domain would be a better name?

Jed

------------------------------------------------------------------------------
Start uncovering the many advantages of virtual appliances
and start using them to simplify application deployment and
accelerate your shift to cloud computing.
http://p.sf.net/sfu/novell-sfdev2dev
_______________________________________________
Libmesh-users mailing list
[email protected]
https://lists.sourceforge.net/lists/listinfo/libmesh-users

Reply via email to