static char help[] = "Playing with block sizes.\n\n";
#include <petscvec.h>

#undef __FUNCT__
#define __FUNCT__ "main"
int main(int argc,char **argv)
{
  Vec            x1,x2;
  PetscInt       bs1,bs2,n = 10;
  PetscErrorCode ierr;
  
  ierr = PetscInitialize(&argc,&argv,(char *) 0,help);CHKERRQ(ierr); 
  ierr = PetscOptionsGetInt(PETSC_NULL,"-n",&n,PETSC_NULL);CHKERRQ(ierr);

  /*
    Everything is fine here
  */
  ierr = VecCreate(PETSC_COMM_WORLD,&x1);CHKERRQ(ierr);
  ierr = VecSetSizes(x1,PETSC_DECIDE,n);CHKERRQ(ierr);
  ierr = VecSetFromOptions(x1);CHKERRQ(ierr);

  ierr = VecGetBlockSize(x1,&bs1);CHKERRQ(ierr);
  ierr = PetscPrintf(PETSC_COMM_WORLD,"x1 block size is %i\n",bs1);CHKERRQ(ierr);

  //ierr = VecDuplicate(x1,&x2);CHKERRQ(ierr);

  ierr = VecCreate(PETSC_COMM_WORLD,&x2);CHKERRQ(ierr);
  ierr = VecSetSizes(x2,PETSC_DECIDE,n);CHKERRQ(ierr);
  ierr = VecSetFromOptions(x2);CHKERRQ(ierr);
  ierr = VecSetBlockSize(x2,2);CHKERRQ(ierr);

  ierr = VecGetBlockSize(x2,&bs2);CHKERRQ(ierr);
  ierr = PetscPrintf(PETSC_COMM_WORLD,"x2 block size is %i\n",bs2);CHKERRQ(ierr);

  ierr = VecGetBlockSize(x1,&bs1);CHKERRQ(ierr);
  ierr = PetscPrintf(PETSC_COMM_WORLD,"x1 block size is %i\n",bs1);CHKERRQ(ierr);
  
  ierr = VecDestroy(&x1);CHKERRQ(ierr);
  ierr = VecDestroy(&x2);CHKERRQ(ierr);

  /*
    Same thing where x2 is obtained from VecDuplicate
  */
  ierr = VecCreate(PETSC_COMM_WORLD,&x1);CHKERRQ(ierr);
  ierr = VecSetSizes(x1,PETSC_DECIDE,n);CHKERRQ(ierr);
  ierr = VecSetFromOptions(x1);CHKERRQ(ierr);

  ierr = VecGetBlockSize(x1,&bs1);CHKERRQ(ierr);
  ierr = PetscPrintf(PETSC_COMM_WORLD,"x1 block size is %i\n",bs1);CHKERRQ(ierr);

  ierr = VecDuplicate(x1,&x2);CHKERRQ(ierr);
  ierr = VecSetBlockSize(x2,2);CHKERRQ(ierr);

  ierr = VecGetBlockSize(x2,&bs2);CHKERRQ(ierr);
  ierr = PetscPrintf(PETSC_COMM_WORLD,"x2 block size is %i\n",bs2);CHKERRQ(ierr);

  ierr = VecGetBlockSize(x1,&bs1);CHKERRQ(ierr);
  ierr = PetscPrintf(PETSC_COMM_WORLD,"x1 block size is %i\n",bs1);CHKERRQ(ierr);
  /*
    WTF, 2???
  */
  
  ierr = VecDestroy(&x1);CHKERRQ(ierr);
  ierr = VecDestroy(&x2);CHKERRQ(ierr);
  ierr = PetscFinalize();
  PetscFunctionReturn(0);
}