Actual source code: dmi.c

petsc-3.4.2 2013-07-02
  1: #include <petsc-private/dmimpl.h>     /*I      "petscdm.h"     I*/

  5: PetscErrorCode DMCreateGlobalVector_Section_Private(DM dm,Vec *vec)
  6: {
  7:   PetscSection   gSection;
  8:   PetscInt       localSize, bs, blockSize = -1, pStart, pEnd, p;

 12:   DMGetDefaultGlobalSection(dm, &gSection);
 13:   PetscSectionGetChart(gSection, &pStart, &pEnd);
 14:   for (p = pStart; p < pEnd; ++p) {
 15:     PetscInt dof, cdof;

 17:     PetscSectionGetDof(gSection, p, &dof);
 18:     PetscSectionGetConstraintDof(gSection, p, &cdof);
 19:     if ((blockSize < 0) && (dof > 0) && (dof-cdof > 0)) blockSize = dof-cdof;
 20:     if ((dof > 0) && (dof-cdof != blockSize)) {
 21:       blockSize = 1;
 22:       break;
 23:     }
 24:   }
 25:   if (blockSize < 0) blockSize = 1;
 26:   MPI_Allreduce(&blockSize, &bs, 1, MPIU_INT, MPI_MIN, PetscObjectComm((PetscObject)dm));
 27:   PetscSectionGetConstrainedStorageSize(gSection, &localSize);
 28:   if (localSize%blockSize) SETERRQ2(PetscObjectComm((PetscObject)dm), PETSC_ERR_ARG_WRONG, "Mismatch between blocksize %d and local storage size %d", blockSize, localSize);
 29:   VecCreate(PetscObjectComm((PetscObject)dm), vec);
 30:   VecSetSizes(*vec, localSize, PETSC_DETERMINE);
 31:   VecSetBlockSize(*vec, bs);
 32:   /* VecSetType(*vec, dm->vectype); */
 33:   VecSetFromOptions(*vec);
 34:   VecSetDM(*vec, dm);
 35:   /* VecSetLocalToGlobalMapping(*vec, dm->ltogmap); */
 36:   /* VecSetLocalToGlobalMappingBlock(*vec, dm->ltogmapb); */
 37:   return(0);
 38: }

 42: PetscErrorCode DMCreateLocalVector_Section_Private(DM dm,Vec *vec)
 43: {
 44:   PetscSection   section;
 45:   PetscInt       localSize, blockSize = -1, pStart, pEnd, p;

 49:   DMGetDefaultSection(dm, &section);
 50:   PetscSectionGetChart(section, &pStart, &pEnd);
 51:   for (p = pStart; p < pEnd; ++p) {
 52:     PetscInt dof;

 54:     PetscSectionGetDof(section, p, &dof);
 55:     if ((blockSize < 0) && (dof > 0)) blockSize = dof;
 56:     if ((dof > 0) && (dof != blockSize)) {
 57:       blockSize = 1;
 58:       break;
 59:     }
 60:   }
 61:   PetscSectionGetStorageSize(section, &localSize);
 62:   VecCreate(PETSC_COMM_SELF, vec);
 63:   VecSetSizes(*vec, localSize, localSize);
 64:   VecSetBlockSize(*vec, blockSize);
 65:   VecSetFromOptions(*vec);
 66:   VecSetDM(*vec, dm);
 67:   return(0);
 68: }