C++ and MPI how to write part of code as parallel?
===

I've been writing some code using PETSc library and now I'm going to change a part of it to be run as parallel. Most of the things what I want to parallelize is matrix initializings and the parts where I generate and calculate a large amount of values. Anyway my problem is following if I run the code with more than 1 core for some reason all parts of the code will be run as many times as how many cores I use.

This is just simple sample code where I tested PETSc and MPI

    int main(int argc, char** argv)
    {
     	time_t rawtime;
    	time ( &rawtime );
    	string sta = ctime (&rawtime);
    	cout << "Solving began..." << endl;
    	
    PetscInitialize(&argc, &argv, 0, 0);
    
      Mat            A;        /* linear system matrix */
      PetscInt       i,j,Ii,J,Istart,Iend,m = 120000,n = 3,its;
      PetscErrorCode ierr;
      PetscBool      flg = PETSC_FALSE;
      PetscScalar    v;
    #if defined(PETSC_USE_LOG)
      PetscLogStage  stage;
    #endif
    
      /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 
             Compute the matrix and right-hand-side vector that define
             the linear system, Ax = b.
         - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
      /* 
         Create parallel matrix, specifying only its global dimensions.
         When using MatCreate(), the matrix format can be specified at
         runtime. Also, the parallel partitioning of the matrix is
         determined by PETSc at runtime.
    
         Performance tuning note:  For problems of substantial size,
         preallocation of matrix memory is crucial for attaining good 
         performance. See the matrix chapter of the users manual for details.
      */
      ierr = MatCreate(PETSC_COMM_WORLD,&A);CHKERRQ(ierr);
      ierr = MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,m,n);CHKERRQ(ierr);
      ierr = MatSetFromOptions(A);CHKERRQ(ierr);
      ierr = MatMPIAIJSetPreallocation(A,5,PETSC_NULL,5,PETSC_NULL);CHKERRQ(ierr);
      ierr = MatSeqAIJSetPreallocation(A,5,PETSC_NULL);CHKERRQ(ierr);
      ierr = MatSetUp(A);CHKERRQ(ierr);
    
      /* 
         Currently, all PETSc parallel matrix formats are partitioned by
         contiguous chunks of rows across the processors.  Determine which
         rows of the matrix are locally owned. 
      */
      ierr = MatGetOwnershipRange(A,&Istart,&Iend);CHKERRQ(ierr);
    
      /* 
         Set matrix elements for the 2-D, five-point stencil in parallel.
          - Each processor needs to insert only elements that it owns
            locally (but any non-local elements will be sent to the
            appropriate processor during matrix assembly). 
          - Always specify global rows and columns of matrix entries.
    
         Note: this uses the less common natural ordering that orders first
         all the unknowns for x = h then for x = 2h etc; Hence you see J = Ii +- n
         instead of J = I +- m as you might expect. The more standard ordering
         would first do all variables for y = h, then y = 2h etc.
    
       */
    PetscMPIInt    rank;        // processor rank
    PetscMPIInt    size;        // size of communicator
    MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
    MPI_Comm_size(PETSC_COMM_WORLD,&size);
    
    cout << "Rank = " << rank << endl;
    cout << "Size = " << size << endl;
    
    cout << "Generating 2D-Array" << endl;
    
    double temp2D[120000][3];
     for (Ii=Istart; Ii<Iend; Ii++) { 
        for(J=0; J<n;J++){
          temp2D[Ii][J] = 1;
        }
      }
      cout << "Processor " << rank << " set values : " << Istart << " - " << Iend << " into 2D-Array" << endl;
      
      v = -1.0;
      for (Ii=Istart; Ii<Iend; Ii++) { 
        for(J=0; J<n;J++){
           MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);CHKERRQ(ierr);
       }
      }
      cout << "Ii = " << Ii << " processor " << rank << " and it owns: " << Istart << " - " << Iend << endl;
    
      /* 
         Assemble matrix, using the 2-step process:
           MatAssemblyBegin(), MatAssemblyEnd()
         Computations can be done while messages are in transition
         by placing code between these two statements.
      */
      ierr = MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);
      ierr = MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);
    
    	MPI_Finalize();
    cout << "No more MPI" << endl;
    return 0;
    
    }

And my real program has a couple different .cpp files. I initialize MPI in the main program what calls a function in another .cpp file where I did implement same kind of matrix filling but all the cout's what the program does before filling the matrices will be printed as many times as the number of my cores.

I can run my test program as mpiexec -n 4 test and it runs successfully but for some reason I have to run my real program as mpiexec -n 4 ./myprog

Output of my test program is following

    Solving began...
    Solving began...
    Solving began...
    Solving began...
    Rank = 0
    Size = 4
    Generating 2D-Array
    Processor 0 set values : 0 - 30000 into 2D-Array
    Rank = 2
    Size = 4
    Generating 2D-Array
    Processor 2 set values : 60000 - 90000 into 2D-Array
    Rank = 3
    Size = 4
    Generating 2D-Array
    Processor 3 set values : 90000 - 120000 into 2D-Array
    Rank = 1
    Size = 4
    Generating 2D-Array
    Processor 1 set values : 30000 - 60000 into 2D-Array
    Ii = 30000 processor 0 and it owns: 0 - 30000
    Ii = 90000 processor 2 and it owns: 60000 - 90000
    Ii = 120000 processor 3 and it owns: 90000 - 120000
    Ii = 60000 processor 1 and it owns: 30000 - 60000
    no more MPI
    no more MPI
    no more MPI
    no more MPI

