Setting options on matrices obtained from a DM

2009-07-22 Thread Jed Brown
If I get a matrix from a DM, I usually do something like

  PetscOptionsGetString(NULL,-q1_mat_type,mtype,sizeof(mtype),NULL);
  DMGetMatrix(dm,mtype,Jp);

so that the type can be changed on the command line (I frequently switch
between preconditioners that can use BAIJ and those that cannot).  But,
MatSetFromOptions never gets called.  If I add

  MatSetOptionsPrefix(Jp,q1_);
  MatSetFromOptions(Jp);

after the above, then preallocation will be blown if I change anything
on the command line.  This includes setting block size for AIJ formats
(I think BAIJ should implement MatSetBlockSize and just confirm that it
matches, but that's a separate issue) so MatSetValuesBlocked cannot
even assemble the correct matrix.

To rectify this, I think that DMGetMatrix either must have enough
information to call MatSetFromOptions (like the prefix) or preallocation
should be separated from setting the sizes (so the user can get the
matrix, call MatSetFromOptions, and then ask the DM to preallocate).
In the former case, we have functionality that is essentially
inaccessible from code (PetscOptionsSetValue is a workaround).  The
latter significantly complicates matters.

Note that since the user may get multiple matrices from a DM, to be used
for different purposes, it's not okay to use the DM's prefix as the
Mat's prefix.

Any preferences here?  I'm happy with any solution in which the user
doesn't have to do an inordinate amount of work to make the following
work.

  -q1_mat_type baij
  -q1_mat_type aij -mat_no_inode
  -q1_mat_type crl


Jed

-- next part --
A non-text attachment was scrubbed...
Name: signature.asc
Type: application/pgp-signature
Size: 260 bytes
Desc: OpenPGP digital signature
URL: 
http://lists.mcs.anl.gov/pipermail/petsc-dev/attachments/20090722/785e6308/attachment.pgp


Setting options on matrices obtained from a DM

2009-07-22 Thread Jed Brown
Barry Smith wrote:
 
 My eventual goal is to remove all hidden MatCreate()'s and
 instead have Mat's passed in as arguments. So MatLoad() would take a
 Mat as an argument, as would DAGetMatrix() One could then set as much
 or a as little as a Mat as they like before passing it in.  They could
 set the type, they could set the sizes, they could set the prefix.

Right, though setting sizes before DAGetMatrix() wouldn't be meaningful.

I think this is the direction you want? Will it satisfy all your needs?
 
MatCreate().
MatSetOptionsPrefix()
MatSetFromOptions()
DAGetMatrix()

This will work at present, but it's not perfect (without Mat
adjustments).  Some options would be set as a consequence of the
MatSetSizes() which occurs within DMGetMatrix() (ops-create() is saved
until the sizes are known).  If a matrix type implements
MatSetFromOptions (currently just SchurComplement), that function would
never get called in this scheme (because that operation won't be set
until the sizes are known).  This could be handled by modifying the
matrix types so that MatCreate_XXX() could be called before
MatSetSizes(), but this is nontrivial.

Jed

-- next part --
A non-text attachment was scrubbed...
Name: signature.asc
Type: application/pgp-signature
Size: 260 bytes
Desc: OpenPGP digital signature
URL: 
http://lists.mcs.anl.gov/pipermail/petsc-dev/attachments/20090722/1bfbe283/attachment.pgp


Setting options on matrices obtained from a DM

2009-07-22 Thread Jed Brown
Barry Smith wrote:

Could you call MatSetFromOptions() twice, once before the
 DAGetMatrix() call to set the type and then after the DAGetMatrix() to
 set particular options for what you set?

This is okay as long as we don't get more global matrix options (this
state appears to be stable).  It's a little clumsy to have to call that
function twice, and because preallocation comes last in the usual idiom:

  MatCreate
  MatSetSizes
  MatSetFromOptions
  MatXXXSetPreallocation

If the matrix type wants anything (programmatic) to be done before
preallocation, we don't have a way to do it with your scheme.  That is,
although we might have set the type before calling DAGetMatrix(), the
type isn't actually set until the sizes are known, so we can't call any
MatXXXSetYYY() before calling DAGetMatrix().  At the moment, I can't
think of a case were this would be a problem.  The same applies to
anything that could be set by MatSetFromOptions_XXX, but needs to be
done before preallocation.  The hypothetical

  DMGetMatrix   /* sets sizes, not preallocation or type */
  MatSetOptionsPrefix
  MatSetFromOptions
  DMSetMatPreallocation

makes it possible to set all the options in MatSetFromOptions, and
allows any MatXXXSetYYY to be called before preallocation.  I can't
think of a case where this is more limiting, but that doesn't mean it
doesn't exist.  It does mean two DM functions rather than one which is
not ideal.

We decided along time ago to cluster all the option setting in
 XXXSetFromOptions() but it does impose some limitations. Too many
 limitations. Are there other models? Having random XXGetOption() in
 any old method is dangerous because it is hard to know what and where
 options are going to be set so I think we are stuck with the
 XXXSetFromOptions() model.

I agree that keeping as much as possible in XXXSetFromOptions is the
best alternative.  I've seen the mess caused by having too many places
in which options are consulted.  The PETSc creation model is much nicer
to work with than factory objects which enforce more sequencing.

Jed

-- next part --
A non-text attachment was scrubbed...
Name: signature.asc
Type: application/pgp-signature
Size: 260 bytes
Desc: OpenPGP digital signature
URL: 
http://lists.mcs.anl.gov/pipermail/petsc-dev/attachments/20090722/5b139481/attachment.pgp


SNESSetUp() and matrix-free options

2009-07-22 Thread Lisandro Dalcin
On Wed, Jul 22, 2009 at 12:26 AM, Barry Smithbsmith at mcs.anl.gov wrote:

 Can we at least agree on the simple rules below:

 1) KSPSetFromOptions() requires a previous call to KSPSetOperators()

 2) SNESSetFromOptions() requires a previous call to:
 ?a) SNESSetFunction(), and also SNESSetJacobian() unless matrix-free is
 used.
 ?b) just SNESSetJacobian() if the problem is actually linear (a case
 that SNESComputeFunction() seems to handle)

 3) TSSetFromOptions() should have similar requirements in order to
 being able to satisfy (1) or (2) depending on the problem type.

 Barry, does all this make sense?

 ?Sometimes, so what went wrong with the broken TS example and can it fit in
 this paradgm?


Well, now the example is not broken because I changed SNES to create a
the missing snes-vec_func from the Jacobian matrix if that is
available... All this works because TSSetFromOptions() calls
SNESSetOperators() before issuing SNESSetFromOptions(). However,
nonlinear TS has an asymmetry with SNES and with TS itself:
TSSetRHSFunction() should take a Vec from the user, like
SNESSetFunction. This way, the TS will have a Vec to feed SNES.


 ? This is tough when there are nested sub KSP inside PCs and picking these
 is done at
 the command line.

OK.. You are right...  this is a nightmare...

Asking users to call SNESSetJacobian() should not be required for
matrix-free, but calling KSPSetFromOptions() without having the
operators set in KSP is broken...

Did you realize that -ksp_fischer_guess is broken if the operators are
not set in KSP? I know, using Fischer's method is a nonsense for
SNES... But perhaps they do make sense for a sub-KSP inside a PC...



-- 
Lisandro Dalc?n
---
Centro Internacional de M?todos Computacionales en Ingenier?a (CIMEC)
Instituto de Desarrollo Tecnol?gico para la Industria Qu?mica (INTEC)
Consejo Nacional de Investigaciones Cient?ficas y T?cnicas (CONICET)
PTLC - G?emes 3450, (3000) Santa Fe, Argentina
Tel/Fax: +54-(0)342-451.1594



Setting options on matrices obtained from a DM

2009-07-22 Thread Lisandro Dalcin
On Wed, Jul 22, 2009 at 12:20 AM, Barry Smithbsmith at mcs.anl.gov wrote:

 On Jul 21, 2009, at 7:47 PM, Jed Brown wrote:

 Barry Smith wrote:

 ? ?This is just a limitation of the current implementation due to the way
 the design evolved over time.
 There is nothing intrinsic to the abstract design of PETSc that prevents the
 type from being
 properly processed before the sizes are set.


I agree 100% with you.. the only issue are tons of code to review and fix :-(


 ? This is really yucky API you propose here; that no sane person would
 propose except to work
 around a shortfall in the current implementation.


We could start fixing Vec, next go for Mat ...


-- 
Lisandro Dalc?n
---
Centro Internacional de M?todos Computacionales en Ingenier?a (CIMEC)
Instituto de Desarrollo Tecnol?gico para la Industria Qu?mica (INTEC)
Consejo Nacional de Investigaciones Cient?ficas y T?cnicas (CONICET)
PTLC - G?emes 3450, (3000) Santa Fe, Argentina
Tel/Fax: +54-(0)342-451.1594



SNESSetUp() and matrix-free options

2009-07-22 Thread Barry Smith

On Jul 22, 2009, at 9:39 AM, Lisandro Dalcin wrote:

 On Wed, Jul 22, 2009 at 12:26 AM, Barry Smithbsmith at mcs.anl.gov  
 wrote:

 Can we at least agree on the simple rules below:

 1) KSPSetFromOptions() requires a previous call to KSPSetOperators()

 2) SNESSetFromOptions() requires a previous call to:
  a) SNESSetFunction(), and also SNESSetJacobian() unless matrix- 
 free is
 used.
  b) just SNESSetJacobian() if the problem is actually linear (a case
 that SNESComputeFunction() seems to handle)

 3) TSSetFromOptions() should have similar requirements in order to
 being able to satisfy (1) or (2) depending on the problem type.

 Barry, does all this make sense?

  Sometimes, so what went wrong with the broken TS example and can  
 it fit in
 this paradgm?


 Well, now the example is not broken because I changed SNES to create a
 the missing snes-vec_func from the Jacobian matrix if that is
 available... All this works because TSSetFromOptions() calls
 SNESSetOperators() before issuing SNESSetFromOptions().

This does not seem right. Having TSSetFromOptions() calling  
SNESSetJacobian() and KSPSetOperators()
seems like a bad design. Maybe some decade we should revisit this  
business.

   Barry


 However,
 nonlinear TS has an asymmetry with SNES and with TS itself:
 TSSetRHSFunction() should take a Vec from the user, like
 SNESSetFunction. This way, the TS will have a Vec to feed SNES.


   This is tough when there are nested sub KSP inside PCs and  
 picking these
 is done at
 the command line.

 OK.. You are right...  this is a nightmare...

 Asking users to call SNESSetJacobian() should not be required for
 matrix-free, but calling KSPSetFromOptions() without having the
 operators set in KSP is broken...

 Did you realize that -ksp_fischer_guess is broken if the operators are
 not set in KSP? I know, using Fischer's method is a nonsense for
 SNES... But perhaps they do make sense for a sub-KSP inside a PC...



 -- 
 Lisandro Dalc?n
 ---
 Centro Internacional de M?todos Computacionales en Ingenier?a (CIMEC)
 Instituto de Desarrollo Tecnol?gico para la Industria Qu?mica (INTEC)
 Consejo Nacional de Investigaciones Cient?ficas y T?cnicas (CONICET)
 PTLC - G?emes 3450, (3000) Santa Fe, Argentina
 Tel/Fax: +54-(0)342-451.1594




Setting options on matrices obtained from a DM

2009-07-22 Thread Barry Smith

On Jul 22, 2009, at 9:56 AM, Lisandro Dalcin wrote:

 On Wed, Jul 22, 2009 at 12:20 AM, Barry Smithbsmith at mcs.anl.gov  
 wrote:

 On Jul 21, 2009, at 7:47 PM, Jed Brown wrote:

 Barry Smith wrote:

This is just a limitation of the current implementation due to  
 the way
 the design evolved over time.
 There is nothing intrinsic to the abstract design of PETSc that  
 prevents the
 type from being
 properly processed before the sizes are set.


 I agree 100% with you.. the only issue are tons of code to review  
 and fix :-(


   This is really yucky API you propose here; that no sane person  
 would
 propose except to work
 around a shortfall in the current implementation.


 We could start fixing Vec,

Ok. VecCreate_Seq()/VecCreate_Seq_Private() and  
VecCreate_MPI_Private() do the following

1) create the Vec_xxx struct and copy over the function tables  
(setting the methods)
2) use PetscMapSetUp() to determine the local/global sizes
3) allocate space for the array if needed.

   What we want is if the sizes are not yet set then it should still  
copy over the function table (1), but delay (2) and (3) until the sizes
are set. This could be done by removing the following from VecSetType()
   if (vec-map-n  0  vec-map-N  0) {
 vec-ops-create = r;
   } else {
 ierr = (*r)(vec);CHKERRQ(ierr);
   }
and instead have each individual VecCreate_xxx() do the method setup  
and then
if sizes set) do the PetscMapSetUp() and allocate space
else) set the vec-ops-create pointer to a new function that does the  
PetscMapSetUp() and allocate space

In other words split the VecCreate_xxx() into two functions and call  
the second when the sizes are available.

   Something similar for Matrices?

Barry




Barry


 next go for Mat ...


 -- 
 Lisandro Dalc?n
 ---
 Centro Internacional de M?todos Computacionales en Ingenier?a (CIMEC)
 Instituto de Desarrollo Tecnol?gico para la Industria Qu?mica (INTEC)
 Consejo Nacional de Investigaciones Cient?ficas y T?cnicas (CONICET)
 PTLC - G?emes 3450, (3000) Santa Fe, Argentina
 Tel/Fax: +54-(0)342-451.1594




Setting options on matrices obtained from a DM

2009-07-22 Thread Lisandro Dalcin
On Wed, Jul 22, 2009 at 12:19 PM, Barry Smithbsmith at mcs.anl.gov wrote:

 On Jul 22, 2009, at 9:56 AM, Lisandro Dalcin wrote:

 On Wed, Jul 22, 2009 at 12:20 AM, Barry Smithbsmith at mcs.anl.gov wrote:

 On Jul 21, 2009, at 7:47 PM, Jed Brown wrote:

 Barry Smith wrote:

 ? This is just a limitation of the current implementation due to the way
 the design evolved over time.
 There is nothing intrinsic to the abstract design of PETSc that prevents
 the
 type from being
 properly processed before the sizes are set.


 I agree 100% with you.. the only issue are tons of code to review and fix
 :-(


 ?This is really yucky API you propose here; that no sane person would
 propose except to work
 around a shortfall in the current implementation.


 We could start fixing Vec,

 ? Ok. VecCreate_Seq()/VecCreate_Seq_Private() and VecCreate_MPI_Private() do
 the following

 1) create the Vec_xxx struct and copy over the function tables (setting the
 methods)
 2) use PetscMapSetUp() to determine the local/global sizes
 3) allocate space for the array if needed.

 ?What we want is if the sizes are not yet set then it should still copy over
 the function table (1), but delay (2) and (3) until the sizes
 are set. This could be done by removing the following from VecSetType()
 ?if (vec-map-n  0  vec-map-N  0) {
 ? ?vec-ops-create = r;
 ?} else {
 ? ?ierr = (*r)(vec);CHKERRQ(ierr);
 ?}
 and instead have each individual VecCreate_xxx() do the method setup and
 then
 if sizes set) do the PetscMapSetUp() and allocate space
 else) set the vec-ops-create pointer to a new function that does the
 PetscMapSetUp() and allocate space

 In other words split the VecCreate_xxx() into two functions and call the
 second when the sizes are available.


Something in such direction should work... However...


 ?Something similar for Matrices?


If you do MatSetType(A, MATAIJ), you end-up having a MATSEQAIJ or a
MPIMPIAIJ... Vec does not have such aliasing on type names ...


I'll try to give look at Vec to figure out what to do in a way that
can be ported to Mat taking into account my previous comment.


-- 
Lisandro Dalc?n
---
Centro Internacional de M?todos Computacionales en Ingenier?a (CIMEC)
Instituto de Desarrollo Tecnol?gico para la Industria Qu?mica (INTEC)
Consejo Nacional de Investigaciones Cient?ficas y T?cnicas (CONICET)
PTLC - G?emes 3450, (3000) Santa Fe, Argentina
Tel/Fax: +54-(0)342-451.1594



Setting options on matrices obtained from a DM

2009-07-22 Thread Jed Brown
Lisandro Dalcin wrote:

 If you do MatSetType(A, MATAIJ), you end-up having a MATSEQAIJ or a
 MPIMPIAIJ... Vec does not have such aliasing on type names ...

How is this a problem?  You aren't matching type names, you're just
setting a function to be called after the sizes are known.  Since this
happens in MatCreate_{Seq,MPI}AIJ, it doesn't matter whether MatSetType
was called with MATAIJ.

Also note that we currently have mat-ops-setsizes (only implemented by
SeqDense) so we should either use it or remove it, as long as we value
one-way-to-do-it.


Jed

-- next part --
A non-text attachment was scrubbed...
Name: signature.asc
Type: application/pgp-signature
Size: 260 bytes
Desc: OpenPGP digital signature
URL: 
http://lists.mcs.anl.gov/pipermail/petsc-dev/attachments/20090722/9135cfd4/attachment.pgp


Setting options on matrices obtained from a DM

2009-07-22 Thread Matthew Knepley
No comments? Not even This is complete shit!?

  Matt

On Tue, Jul 21, 2009 at 8:56 PM, Matthew Knepley knepley at gmail.com wrote:

 Actually, this seems like the same problem that Lisandro is having, just
 with
 different functions. I propose making data structures do the work for us
 rather than complicated organization in an imperative program.

 We could use the same mechanism we use in configure to handle issues of
 object setup. We have a setup object that takes

   a) an object to be setup

   b) a set of functions to be called for setup, and

   c) any functions which must be called prior to each given function

 This way we can flexibly add as many functions as necessary, and then they
 can be topologically sorted and executed.

 There will be some implementation issues in C, due to severe limitations
 with
 calling conventions, but I think these can all be solved.

 Comments?

Matt


 On Tue, Jul 21, 2009 at 7:47 PM, Jed Brown jed at 59a2.org wrote:

 Barry Smith wrote:

 Could you call MatSetFromOptions() twice, once before the
  DAGetMatrix() call to set the type and then after the DAGetMatrix() to
  set particular options for what you set?

 This is okay as long as we don't get more global matrix options (this
 state appears to be stable).  It's a little clumsy to have to call that
 function twice, and because preallocation comes last in the usual idiom:

  MatCreate
  MatSetSizes
  MatSetFromOptions
  MatXXXSetPreallocation

 If the matrix type wants anything (programmatic) to be done before
 preallocation, we don't have a way to do it with your scheme.  That is,
 although we might have set the type before calling DAGetMatrix(), the
 type isn't actually set until the sizes are known, so we can't call any
 MatXXXSetYYY() before calling DAGetMatrix().  At the moment, I can't
 think of a case were this would be a problem.  The same applies to
 anything that could be set by MatSetFromOptions_XXX, but needs to be
 done before preallocation.  The hypothetical

  DMGetMatrix   /* sets sizes, not preallocation or type */
  MatSetOptionsPrefix
  MatSetFromOptions
  DMSetMatPreallocation

 makes it possible to set all the options in MatSetFromOptions, and
 allows any MatXXXSetYYY to be called before preallocation.  I can't
 think of a case where this is more limiting, but that doesn't mean it
 doesn't exist.  It does mean two DM functions rather than one which is
 not ideal.

 We decided along time ago to cluster all the option setting in
  XXXSetFromOptions() but it does impose some limitations. Too many
  limitations. Are there other models? Having random XXGetOption() in
  any old method is dangerous because it is hard to know what and where
  options are going to be set so I think we are stuck with the
  XXXSetFromOptions() model.

 I agree that keeping as much as possible in XXXSetFromOptions is the
 best alternative.  I've seen the mess caused by having too many places
 in which options are consulted.  The PETSc creation model is much nicer
 to work with than factory objects which enforce more sequencing.

 Jed




 --
 What most experimenters take for granted before they begin their
 experiments is infinitely more interesting than any results to which their
 experiments lead.
 -- Norbert Wiener




-- 
What most experimenters take for granted before they begin their experiments
is infinitely more interesting than any results to which their experiments
lead.
-- Norbert Wiener
-- next part --
An HTML attachment was scrubbed...
URL: 
http://lists.mcs.anl.gov/pipermail/petsc-dev/attachments/20090722/0a1e13cc/attachment.html


Vec: types, sizes, and block sizes...

2009-07-22 Thread Lisandro Dalcin
Starting a new thread ...

Here you have simple example about creating Vec's. PETSc let me set
first the type, or first the size and block size... However...
setting FIRST the TYPE is currently broken...

Anything we are going to do to fix this will require to handle the
case above... At first, the only way I can imagine of doing that is by
using a macro VecPreallocated() (more or less like in Mat), in charge
of checking and setting-up the sizesblocksizes and the allocated
array. But then this macro have to be called on almost all Vec
interface calls...

Note: the code below have to be run in 2(two) processes to show the
problems (at second block, when creating vector y)...

#include petscvec.h

#undef __FUNCT__
#define __FUNCT__ main
int main(int argc,char **argv)
{
  PetscErrorCode ierr;
  ierr = PetscInitialize(argc,argv,0,0);CHKERRQ(ierr);
  if (1) {
Vec x;
ierr = VecCreate(PETSC_COMM_WORLD, x);CHKERRQ(ierr);
ierr = VecSetSizes(x, PETSC_DECIDE, 9);CHKERRQ(ierr);
ierr = VecSetBlockSize(x, 3);CHKERRQ(ierr);
ierr = VecSetType(x, VECMPI);CHKERRQ(ierr);
ierr = VecDestroy(x);CHKERRQ(ierr);
  }
  if (1) {
Vec y;
ierr = VecCreate(PETSC_COMM_WORLD, y);CHKERRQ(ierr);
ierr = VecSetType(y, VECMPI);CHKERRQ(ierr);
ierr = VecSetSizes(y, PETSC_DECIDE, 9);CHKERRQ(ierr);
ierr = VecSetBlockSize(y, 3);CHKERRQ(ierr);
ierr = VecDestroy(y);CHKERRQ(ierr);
  }
  ierr = PetscFinalize();CHKERRQ(ierr);
  return 0;
}



-- 
Lisandro Dalc?n
---
Centro Internacional de M?todos Computacionales en Ingenier?a (CIMEC)
Instituto de Desarrollo Tecnol?gico para la Industria Qu?mica (INTEC)
Consejo Nacional de Investigaciones Cient?ficas y T?cnicas (CONICET)
PTLC - G?emes 3450, (3000) Santa Fe, Argentina
Tel/Fax: +54-(0)342-451.1594



Setting options on matrices obtained from a DM

2009-07-22 Thread Jed Brown
Matthew Knepley wrote:
 No comments? Not even This is complete shit!?

Heh, it might be overkill and, if I understand you correctly, I think it
could be a headache.  There is a reasonable amount of code in the
various interfaces to ensure some sequencing, but delayed evaluation is
tough to reason about.  When I call a function, I expect strict
evaluation, but if you somehow store that call away until it's
dependencies are satisfied, I'm likely to see confusing behavior.

 On Tue, Jul 21, 2009 at 8:56 PM, Matthew Knepley knepley at gmail.com wrote:
 
 Actually, this seems like the same problem that Lisandro is having, just
 with
 different functions. I propose making data structures do the work for us
 rather than complicated organization in an imperative program.

 We could use the same mechanism we use in configure to handle issues of
 object setup. We have a setup object that takes

   a) an object to be setup

   b) a set of functions to be called for setup, and

   c) any functions which must be called prior to each given function

How would you specify these?  In the general case, we have a partial
ordering on setup functions.  I suppose you are suggesting that any time
you register a setup function, all dependencies would be given
explicitly (hence added to the DAG, or referenced if they were already
there).

It seems to me that we have rather few setup steps that require a
particular sequence, some of which are not necessary.  I would rather
refactor to eliminate unnecessary dependencies than deal with delayed
evaluation.  If that looks especially complex, some dependency-handling
object might help (probably just within the implementations for which it
is necessarily complex), but I don't think it should be visible to the
user.

I think it would help to enumerate the dependencies that really are
fundamental, and those that could (with work) be eliminated.  Here's a
start, I'm certainly missing some.

Vec

  * sizes before type, not fundamental

Mat

  * sizes before type, not fundamental, currently cached when type is
set first which causes other problems
  * type before preallocation, fundamental

KSP/PC

  * operators set before SetFromOptions, is this fundamental?
  * operators set before SetUp, fundamental

SNES

  * function/Jacobian before SetFromOptions, is this fundamental?
  * function/Jacobian before SetUp, fundamental

TS

  * similar to SNES


Some cases, such as PC_ASM and PC_FieldSplit, would need different
options management in order to set all the options in SetFromOptions.
They create new matrices and KSP objects from the preconditioning matrix
which doesn't exist before Solve/SetUp.  It would be possible to move
most of this into PCSetFromOptions_XXX if we could call
KSPSetFromOptions before KSPSetOperators.  It would also make it
possible for the user to control these programmatically (currently the
inner KSP cannot be created until solve, at which point it's too late
for the user to manipulate it).

Ideally we would have no dependencies between personality (types and
options) and mechanics (for lack of a better term -- sizes, operators,
preallocation).  The former can influence the latter (changing the type
changes the mechanics), but I think a mandatory ordering is bad (except
where necessary, setting an option for a specific type obviously must
happen after that type is set).


Jed

-- next part --
A non-text attachment was scrubbed...
Name: signature.asc
Type: application/pgp-signature
Size: 260 bytes
Desc: OpenPGP digital signature
URL: 
http://lists.mcs.anl.gov/pipermail/petsc-dev/attachments/20090722/f03e611c/attachment.pgp


Vec: types, sizes, and block sizes...

2009-07-22 Thread Barry Smith

 VecCreate()
 VecSetXXX
 VecSetYYY
 VecSetBlockSize()
 VecSetUp()   (called also automatically as needed by later vector  
routines, like PCSetUp()/KSPSetUp()/SNESSetUp()/TSSetUp() are called)
   use the beast

 To make this universal we need to fix Mat to have the  
MatXXXSetPreallocation() to just stash the preallocation information  
not do the allocation on the spot, but wait until the MatSetUp().


Barry


On Jul 22, 2009, at 1:47 PM, Lisandro Dalcin wrote:

 Starting a new thread ...

 Here you have simple example about creating Vec's. PETSc let me set
 first the type, or first the size and block size... However...
 setting FIRST the TYPE is currently broken...

 Anything we are going to do to fix this will require to handle the
 case above... At first, the only way I can imagine of doing that is by
 using a macro VecPreallocated() (more or less like in Mat), in charge
 of checking and setting-up the sizesblocksizes and the allocated
 array. But then this macro have to be called on almost all Vec
 interface calls...

 Note: the code below have to be run in 2(two) processes to show the
 problems (at second block, when creating vector y)...

 #include petscvec.h

 #undef __FUNCT__
 #define __FUNCT__ main
 int main(int argc,char **argv)
 {
  PetscErrorCode ierr;
  ierr = PetscInitialize(argc,argv,0,0);CHKERRQ(ierr);
  if (1) {
Vec x;
ierr = VecCreate(PETSC_COMM_WORLD, x);CHKERRQ(ierr);
ierr = VecSetSizes(x, PETSC_DECIDE, 9);CHKERRQ(ierr);
ierr = VecSetBlockSize(x, 3);CHKERRQ(ierr);
ierr = VecSetType(x, VECMPI);CHKERRQ(ierr);
ierr = VecDestroy(x);CHKERRQ(ierr);
  }
  if (1) {
Vec y;
ierr = VecCreate(PETSC_COMM_WORLD, y);CHKERRQ(ierr);
ierr = VecSetType(y, VECMPI);CHKERRQ(ierr);
ierr = VecSetSizes(y, PETSC_DECIDE, 9);CHKERRQ(ierr);
ierr = VecSetBlockSize(y, 3);CHKERRQ(ierr);
ierr = VecDestroy(y);CHKERRQ(ierr);
  }
  ierr = PetscFinalize();CHKERRQ(ierr);
  return 0;
 }



 -- 
 Lisandro Dalc?n
 ---
 Centro Internacional de M?todos Computacionales en Ingenier?a (CIMEC)
 Instituto de Desarrollo Tecnol?gico para la Industria Qu?mica (INTEC)
 Consejo Nacional de Investigaciones Cient?ficas y T?cnicas (CONICET)
 PTLC - G?emes 3450, (3000) Santa Fe, Argentina
 Tel/Fax: +54-(0)342-451.1594




Vec: types, sizes, and block sizes...

2009-07-22 Thread Lisandro Dalcin
On Wed, Jul 22, 2009 at 6:53 PM, Barry Smithbsmith at mcs.anl.gov wrote:

 ? ?VecCreate()
 ? ?VecSetXXX
 ? ?VecSetYYY
 ? ?VecSetBlockSize()
 ? ?VecSetUp() ? (called also automatically as needed by later vector
 routines, like PCSetUp()/KSPSetUp()/SNESSetUp()/TSSetUp() are called)
 ? ? ?use the beast


OK, that was what I talked about... this change requires work from our
side, but the end-user API does not change at all... Nice...

Barry, let's elaborate a bit more all this:

1) VecCreate_XXX have stop allocating the array.
2) VecOpts needs to gain a 'setup' slot. The  new slot is filled with
VecSetUp_XXX(), which basically do array allocation, but handling the
case of user-provided arrays with VecCreateXXXWithArray().
3) VecSetUp() will require at least a previous call to
VecSetSizes/VecSetBlocSize(). If the type is not set, it choses SEQ or
MPI depending on the size of Vec's MPI_Comm. After all this,
vec-ops-setup() is called.
4) Fix almost all Vec interface routines to call VecSetUp().

Am I missing something?


-- 
Lisandro Dalc?n
---
Centro Internacional de M?todos Computacionales en Ingenier?a (CIMEC)
Instituto de Desarrollo Tecnol?gico para la Industria Qu?mica (INTEC)
Consejo Nacional de Investigaciones Cient?ficas y T?cnicas (CONICET)
PTLC - G?emes 3450, (3000) Santa Fe, Argentina
Tel/Fax: +54-(0)342-451.1594