1*11689aebSJed Brown #include <petsc-private/dmimpl.h> /*I "petscdm.h" I*/ 2*11689aebSJed Brown 3*11689aebSJed Brown extern PetscErrorCode VecView_Seq(Vec, PetscViewer); 4*11689aebSJed Brown extern PetscErrorCode VecView_MPI(Vec, PetscViewer); 5*11689aebSJed Brown 6*11689aebSJed Brown #undef __FUNCT__ 7*11689aebSJed Brown #define __FUNCT__ "DMCreateGlobalVector_Section_Private" 8*11689aebSJed Brown PetscErrorCode DMCreateGlobalVector_Section_Private(DM dm,Vec *vec) 9*11689aebSJed Brown { 10*11689aebSJed Brown PetscErrorCode ierr; 11*11689aebSJed Brown PetscSection gSection; 12*11689aebSJed Brown PetscInt localSize, blockSize = -1, pStart, pEnd, p; 13*11689aebSJed Brown 14*11689aebSJed Brown PetscFunctionBegin; 15*11689aebSJed Brown ierr = DMGetDefaultGlobalSection(dm, &gSection);CHKERRQ(ierr); 16*11689aebSJed Brown ierr = PetscSectionGetChart(dm->defaultSection, &pStart, &pEnd);CHKERRQ(ierr); 17*11689aebSJed Brown for(p = pStart; p < pEnd; ++p) { 18*11689aebSJed Brown PetscInt dof, cdof; 19*11689aebSJed Brown 20*11689aebSJed Brown ierr = PetscSectionGetDof(dm->defaultSection, p, &dof);CHKERRQ(ierr); 21*11689aebSJed Brown ierr = PetscSectionGetConstraintDof(dm->defaultSection, p, &cdof);CHKERRQ(ierr); 22*11689aebSJed Brown if ((blockSize < 0) && (dof > 0)) blockSize = dof-cdof; 23*11689aebSJed Brown if ((dof > 0) && (dof-cdof != blockSize)) { 24*11689aebSJed Brown blockSize = 1; 25*11689aebSJed Brown break; 26*11689aebSJed Brown } 27*11689aebSJed Brown } 28*11689aebSJed Brown ierr = PetscSectionGetConstrainedStorageSize(dm->defaultGlobalSection, &localSize);CHKERRQ(ierr); 29*11689aebSJed Brown if (localSize%blockSize) SETERRQ2(((PetscObject) dm)->comm, PETSC_ERR_ARG_WRONG, "Mismatch between blocksize %d and local storage size %d", blockSize, localSize); 30*11689aebSJed Brown ierr = VecCreate(((PetscObject) dm)->comm, vec);CHKERRQ(ierr); 31*11689aebSJed Brown ierr = VecSetSizes(*vec, localSize, PETSC_DETERMINE);CHKERRQ(ierr); 32*11689aebSJed Brown ierr = VecSetBlockSize(*vec, blockSize);CHKERRQ(ierr); 33*11689aebSJed Brown /* ierr = VecSetType(*vec, dm->vectype);CHKERRQ(ierr); */ 34*11689aebSJed Brown ierr = VecSetFromOptions(*vec);CHKERRQ(ierr); 35*11689aebSJed Brown ierr = VecSetDM(*vec, dm);CHKERRQ(ierr); 36*11689aebSJed Brown /* ierr = VecSetLocalToGlobalMapping(*vec, dm->ltogmap);CHKERRQ(ierr); */ 37*11689aebSJed Brown /* ierr = VecSetLocalToGlobalMappingBlock(*vec, dm->ltogmapb);CHKERRQ(ierr); */ 38*11689aebSJed Brown PetscFunctionReturn(0); 39*11689aebSJed Brown } 40*11689aebSJed Brown 41*11689aebSJed Brown #undef __FUNCT__ 42*11689aebSJed Brown #define __FUNCT__ "DMCreateLocalVector_Section_Private" 43*11689aebSJed Brown PetscErrorCode DMCreateLocalVector_Section_Private(DM dm,Vec *vec) 44*11689aebSJed Brown { 45*11689aebSJed Brown PetscErrorCode ierr; 46*11689aebSJed Brown PetscInt localSize, blockSize = -1, pStart, pEnd, p; 47*11689aebSJed Brown 48*11689aebSJed Brown PetscFunctionBegin; 49*11689aebSJed Brown ierr = PetscSectionGetChart(dm->defaultSection, &pStart, &pEnd);CHKERRQ(ierr); 50*11689aebSJed Brown for(p = pStart; p < pEnd; ++p) { 51*11689aebSJed Brown PetscInt dof; 52*11689aebSJed Brown 53*11689aebSJed Brown ierr = PetscSectionGetDof(dm->defaultSection, p, &dof);CHKERRQ(ierr); 54*11689aebSJed Brown if ((blockSize < 0) && (dof > 0)) blockSize = dof; 55*11689aebSJed Brown if ((dof > 0) && (dof != blockSize)) { 56*11689aebSJed Brown blockSize = 1; 57*11689aebSJed Brown break; 58*11689aebSJed Brown } 59*11689aebSJed Brown } 60*11689aebSJed Brown ierr = PetscSectionGetStorageSize(dm->defaultSection, &localSize);CHKERRQ(ierr); 61*11689aebSJed Brown ierr = VecCreate(PETSC_COMM_SELF, vec);CHKERRQ(ierr); 62*11689aebSJed Brown ierr = VecSetSizes(*vec, localSize, localSize);CHKERRQ(ierr); 63*11689aebSJed Brown ierr = VecSetBlockSize(*vec, blockSize);CHKERRQ(ierr); 64*11689aebSJed Brown ierr = VecSetFromOptions(*vec);CHKERRQ(ierr); 65*11689aebSJed Brown ierr = VecSetDM(*vec, dm);CHKERRQ(ierr); 66*11689aebSJed Brown PetscFunctionReturn(0); 67*11689aebSJed Brown } 68