Hello, I want to build a global vector (and matrix) from local data. The local data has a global index, which can be non-contiguous i.e. global index (0, 5, 2) is on rank 0, (1, 4, 3) is on rank 1. To keep all local data on the local part of vector I need a mapping:
Application -> PETSc (0, 5, 2) -> (0, 1, 2) (1, 4, 3) -> (3, 4, 5) It can happen that rank 1 also need to touch vector[5] (application ordering) i.e. vector[1] (petsc ordering) on rank 0. I can produce a mapping using index sets: ierr = ISCreateGeneral(PETSC_COMM_WORLD, indizes.size(), indizes.data(), PETSC_COPY_VALUES, &ISlocal); CHKERRV(ierr); ierr = ISAllGather(ISlocal, &ISglobal); CHKERRV(ierr); // Gather the IS from all processors ierr = ISLocalToGlobalMappingCreateIS(ISglobal, &ISmapping); CHKERRV(ierr); // Make it a mapping I construct std::vector indizes by walking though all local indizes. This way I have a mapping ISmapping that looks like that: 0 0 1 5 2 2 3 1 4 4 5 3 Now I can use that mapping: VecSetLocalToGlobalMapping(vec, mapping); but this does not work like it should (for me). If I use VecSetValueLocal to set element 5 it gets mapped to 3, element 1 gets mapped to 5. I want to have it reversed, accessing element 3 gets mapped to 5, element 1 to 3. This loop: for (size_t i = 0; i < indizes.size(); i++) { ierr = VecSetValueLocal(vec, indizes[i], i*10, INSERT_VALUES); CHKERRQ(ierr); } should produce a vector = [0, 10, 20, 30, 40, 50] and not [0, 50, 20, 10, 40, 30]. Are AO (application orderings) more suited for that? But afaik I can not produce a ISLocalToGlobalMapping from an AO that can be used with VecSetValueLocal, can I? Is there a way to inverse a mapping? To produce: 0 0 1 3 2 2 3 5 4 4 5 1 Or am I down the wrong path? Thanks, Florian Example code: #include <iostream> #include <vector> #include <mpi.h> #include "petscmat.h" #include "petscviewer.h" #include "petscksp.h" using namespace std; int main(int argc, char **argv) { PetscInitialize(&argc, &argv, "", NULL); PetscErrorCode ierr; Vec vec; VecCreate(PETSC_COMM_SELF, &vec); std::vector<int> indizes = {0, 5, 2, 1, 4, 3}; IS is; ISLocalToGlobalMapping mapping; ierr = ISCreateGeneral(PETSC_COMM_WORLD, indizes.size(),indizes.data(), PETSC_COPY_VALUES, &is); CHKERRQ(ierr); ierr = ISLocalToGlobalMappingCreateIS(is, &mapping); CHKERRQ(ierr); ierr = VecCreate(PETSC_COMM_WORLD, &vec); CHKERRQ(ierr); ierr = VecSetLocalToGlobalMapping(vec, mapping); CHKERRQ(ierr); ierr = ISLocalToGlobalMappingView(mapping, PETSC_VIEWER_STDOUT_SELF); CHKERRQ(ierr); ierr = VecSetSizes(vec, PETSC_DECIDE, indizes.size()); CHKERRQ(ierr); ierr = VecSetFromOptions(vec); CHKERRQ(ierr); for (size_t i = 0; i < indizes.size(); i++) { ierr = VecSetValueLocal(vec, indizes[i], i*10, INSERT_VALUES); CHKERRQ(ierr); } ierr = VecView(vec, PETSC_VIEWER_STDOUT_WORLD); CHKERRQ(ierr); PetscFinalize(); return 0; }