/**
 * @file
 *
 * Created by Ataollah Mesgarnejad on 3/31/15.
 * Copyright 2011 Ataollah Mesgarnejad. All rights reserved.
 *
 * This is a test distribute field
 *
 */


static char help[] = "An example of the usage of PetscSF to scatter data back and forth between a \
a dist DM and its uniformly distributed parallel counterpart.\n";

#include <petscdmplex.h>
#include <petscsnes.h>
#include <petscsf.h>
#include <exodusII.h>
#include <petsc-private/dmimpl.h>     /*I      "petscdm.h"     I*/
#include <petsc-private/dmpleximpl.h>
#include <petsc-private/vecimpl.h>
#include <petscsf.h>
#include <petscviewerhdf5.h>
#include <petscis.h>

/*void linear(const PetscReal coords[], PetscScalar *u, void *ctx)
 {
 u[0]=coords[0]*coords[1];
 }*/

#undef __FUNCT__
#define __FUNCT__ "main"
int main(int argc,char **argv)
{
    PetscErrorCode      ierr;
    PetscFE             fem;
    DM                  dm, distDM;
    char                ifilename[PETSC_MAX_PATH_LEN];
    PetscBool           flg;
    int                 CPU_word_size = 0,IO_word_size = 0,exoid;
    PetscViewer         stdoutViewer;
    float               version;
    PetscMPIInt         numproc,rank;
    PetscInt            dim;
    PetscInt            numFields = 1;
    PetscInt            numComp[1] = {1};
    PetscInt            numDof[3] = {1, 0, 0}; /*{Vertex, Edge, Cell} */
    PetscInt            bcFields[1] = {0}, numBC=0;
    IS                  bcPoints[1] = {NULL};
    PetscSection        seqSection, distSection;
    PetscSF             pointSF;
    MPI_Comm            comm;
    
    ierr = PetscInitialize(&argc, &argv, (char*)0, help);CHKERRQ(ierr);
    ierr = PetscViewerASCIIGetStdout(PETSC_COMM_SELF,&stdoutViewer);CHKERRQ(ierr);
    comm = PETSC_COMM_WORLD;
    
    ierr = PetscOptionsGetString(PETSC_NULL,"-i",ifilename,sizeof ifilename,&flg);CHKERRQ(ierr);
    
    ierr = MPI_Comm_size(PETSC_COMM_WORLD,&numproc);
    // if (numproc < 2) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_MIN_VALUE, "This example is supposed to be ran in parallel");
    ierr = MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
    if (!rank)
    {
        exoid = ex_open(ifilename,EX_READ,&CPU_word_size,&IO_word_size,&version);
        if (exoid <= 0) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_LIB,"ex_open(\"%s\",...) did not return a valid file ID",ifilename);
    }
    else
    {
        exoid = -1;                 /* Not used */
    }
    ierr = DMPlexCreateExodus(PETSC_COMM_WORLD,exoid,PETSC_FALSE,&dm);CHKERRQ(ierr);
    ierr = DMGetDimension(dm, &dim);CHKERRQ(ierr);
    ierr = DMPlexCreateSection(dm,      dim, numFields  , numComp, numDof, numBC, bcFields, bcPoints,PETSC_NULL, &seqSection);CHKERRQ(ierr);
    ierr = DMSetDefaultSection(dm, seqSection);CHKERRQ(ierr);
    

    ierr = DMPlexDistribute(dm, 1, &pointSF, &distDM);CHKERRQ(ierr);

    ierr = DMPlexCreateSection(distDM, dim, numFields  , numComp, numDof, numBC, bcFields, bcPoints,PETSC_NULL, &distSection);CHKERRQ(ierr);
    ierr = DMSetDefaultSection(distDM, distSection);CHKERRQ(ierr);

    ierr = PetscPrintf(PETSC_COMM_WORLD,"==== PointSF (from distribution)\n"); CHKERRQ(ierr);
    ierr = PetscSFView(pointSF, PETSC_VIEWER_STDOUT_WORLD); CHKERRQ(ierr);
    PetscSynchronizedFlush(PETSC_COMM_WORLD, PETSC_STDOUT);

    PetscInt            nroots,nleaves;
    const PetscInt      *localPoints;
    PetscInt            *copyoflocalPoints;
    PetscSFNode         *copyofremotePoints;
    const PetscSFNode   *remotePoints;
    
    ierr = PetscSFGetGraph(pointSF, &nroots, &nleaves, &localPoints, &remotePoints); CHKERRQ(ierr);
    PetscMalloc2(nleaves, &copyoflocalPoints,nleaves, &copyofremotePoints); CHKERRQ(ierr);
    for (int i=0;i<nleaves;++i)
    {
        copyoflocalPoints[i] = *(localPoints+i);
        copyofremotePoints[i].index = remotePoints[i].index;
        copyofremotePoints[i].rank  = remotePoints[i].rank;
    }

    ierr = DMDestroy(&dm);CHKERRQ(ierr);
    ierr = DMDestroy(&distDM);CHKERRQ(ierr);
    
    ierr = PetscFinalize();
    return 0;
}