<commit_before>#include "apfMPAS.h"
#include <apfMesh2.h>
#include <apfNumbering.h>
#include <diffMC/src/parma_ghostOwner.h>
#include <netcdf>

namespace apf {

using namespace netCDF;

struct MpasFile
{
  ~MpasFile()
  {
    delete [] x;
    delete [] y;
    delete [] z;
    delete [] conn;
  }
  int elements;
  int nodes;
  int nodesPerElement;
  double* x;
  double* y;
  double* z;
/* conn[i * vertexDegree + j] - 1
   is the index of the j'th node
   of element i */
  int* conn;
};

static int readDim(NcFile& f, const char* name)
{
  return f.getDim(name).getSize();
}

static double* readDoubles(NcFile& f, const char* name, int n)
{
  double* a = new double[n];
  NcVar v = f.getVar(name);
  v.getVar(a);
  return a;
}

static int* readInts(NcFile& f, const char* name, int n)
{
  int* a = new int[n];
  NcVar v = f.getVar(name);
  v.getVar(a);
  return a;
}

void readMpasFile(MpasFile& out, const char* filename)
{
  NcFile in(filename, NcFile::read);
  /* this is the dual of a hexagonal mesh,
     hence the reversing of terms */
  out.nodes = readDim(in, "nCells");
  printf("%d nodes\n", out.nodes);
  out.elements = readDim(in, "nVertices");
  printf("%d elements\n", out.elements);
  out.nodesPerElement = readDim(in, "vertexDegree");
  out.x = readDoubles(in, "xCell", out.nodes);
  out.y = readDoubles(in, "yCell", out.nodes);
  out.z = readDoubles(in, "zCell", out.nodes);
  out.conn = readInts(in, "cellsOnVertex", out.elements * out.nodesPerElement);
}

static int getType(int nodesPerElement)
{
  if (nodesPerElement == 3)
    return apf::Mesh::TRIANGLE;
  if (nodesPerElement == 4)
    return apf::Mesh::QUAD;
  abort();
}

void numberInitialNodes(MpasFile& in, apf::Mesh2* out,
    apf::MeshEntity** v)
{
  apf::Numbering* n = apf::createNumbering(out, "mpas_id", out->getShape(), 1);
  for (int i = 0; i < in.nodes; ++i)
    apf::number(n, v[i], 0, 0, i);
}

void buildInitialMesh(MpasFile& in, apf::Mesh2* out)
{
  apf::MeshEntity** v = new apf::MeshEntity* [in.nodes];
/* fake classification for now, deriveMdsModel will fix this up */
  apf::ModelEntity* c = out->findModelEntity(2, 0);
  int t = getType(in.nodesPerElement);
  for (int i = 0; i < in.nodes; ++i)
    v[i] = out->createVert(c);
  for (int i = 0; i < in.nodes; ++i) {
    apf::Vector3 x(in.x[i], in.y[i], in.z[i]);
    out->setPoint(v[i], 0, x);
  }
  for (int i = 0; i < in.elements; ++i) {
    apf::Downward ev;
    for (int j = 0; j < in.nodesPerElement; ++j) {
      int vi = in.conn[i * in.nodesPerElement + j] - 1;
/* zero indices (-1 after conversion) indicate that an MPAS
   vertex does not have an adjacent cell in that direction
   (it is a boundary vertex). In our dual mesh, we simply
   omit that dual element. This is correct and should not
   cause problems. */
      if (vi < 0)
        goto bad_vi;
      assert(vi < in.nodes);
      ev[j] = v[vi];
    }
    apf::buildElement(out, c, t, ev);
bad_vi:
    continue;
  }
  numberInitialNodes(in, out, v);
  delete [] v;
}

void removeIsolatedNodes(apf::Mesh2* m)
{
/* another problem with the dual: MPAS cells without adjacent
   cells become vertices with no adjacent elements in the
   dual mesh, which our codes don't handle.
   We have to give up on these cells, but tell the user about it.
   We also need to be aware of these omitted cells going forward */
  apf::MeshIterator* it = m->begin(0);
  int n = 0;
  apf::MeshEntity* e;
  while ((e = m->iterate(it)))
    if ( ! m->countUpward(e)) {
      m->destroy(e);
      ++n;
    }
  m->end(it);
  if (n)
    fprintf(stderr,
        "warning: removed %d isolated nodes from MPAS dual mesh\n",
        n);
}

void loadMpasMesh(apf::Mesh2* m, const char* filename)
{
  MpasFile f;
  readMpasFile(f, filename);
  buildInitialMesh(f, m);
  removeIsolatedNodes(m);
}


void writeMpasAssignments(apf::Mesh2* m, const char* filename) {
  NcFile in(filename, NcFile::read);
  /* this is the dual of a hexagonal mesh,
     hence the reversing of terms */
  int numMpasVtx = readDim(in, "nCells");
  
  apf::Numbering* n = apf::createNumbering(m, "mpas_id", m->getShape(), 1);

  // collect N/#parts contiguous vertex assignments on each process (and deal with any remainders)
  
  // assign missing vertices to a random part id
  
  // use MPI IO to write the contiguous blocks to a single graph.info.part.<#parts> file
  // see https://gist.github.com/cwsmith/166d5beb400f3a8136f7 and the comments
}

}<commit_msg>vertex passing and missing vertex finding with PCU error<commit_after>#include "apfMPAS.h"
#include <apfMesh2.h>
#include <apfNumbering.h>
#include <diffMC/src/parma_ghostOwner.h>
#include <netcdf>

namespace apf {

using namespace netCDF;

struct MpasFile
{
  ~MpasFile()
  {
    delete [] x;
    delete [] y;
    delete [] z;
    delete [] conn;
  }
  int elements;
  int nodes;
  int nodesPerElement;
  double* x;
  double* y;
  double* z;
/* conn[i * vertexDegree + j] - 1
   is the index of the j'th node
   of element i */
  int* conn;
};

static int readDim(NcFile& f, const char* name)
{
  return f.getDim(name).getSize();
}

static double* readDoubles(NcFile& f, const char* name, int n)
{
  double* a = new double[n];
  NcVar v = f.getVar(name);
  v.getVar(a);
  return a;
}

static int* readInts(NcFile& f, const char* name, int n)
{
  int* a = new int[n];
  NcVar v = f.getVar(name);
  v.getVar(a);
  return a;
}

void readMpasFile(MpasFile& out, const char* filename)
{
  NcFile in(filename, NcFile::read);
  /* this is the dual of a hexagonal mesh,
     hence the reversing of terms */
  out.nodes = readDim(in, "nCells");
  printf("%d nodes\n", out.nodes);
  out.elements = readDim(in, "nVertices");
  printf("%d elements\n", out.elements);
  out.nodesPerElement = readDim(in, "vertexDegree");
  out.x = readDoubles(in, "xCell", out.nodes);
  out.y = readDoubles(in, "yCell", out.nodes);
  out.z = readDoubles(in, "zCell", out.nodes);
  out.conn = readInts(in, "cellsOnVertex", out.elements * out.nodesPerElement);
}

static int getType(int nodesPerElement)
{
  if (nodesPerElement == 3)
    return apf::Mesh::TRIANGLE;
  if (nodesPerElement == 4)
    return apf::Mesh::QUAD;
  abort();
}

void numberInitialNodes(MpasFile& in, apf::Mesh2* out,
    apf::MeshEntity** v)
{
  apf::Numbering* n = apf::createNumbering(out, "mpas_id", out->getShape(), 1);
  for (int i = 0; i < in.nodes; ++i)
    apf::number(n, v[i], 0, 0, i);
}

void buildInitialMesh(MpasFile& in, apf::Mesh2* out)
{
  apf::MeshEntity** v = new apf::MeshEntity* [in.nodes];
/* fake classification for now, deriveMdsModel will fix this up */
  apf::ModelEntity* c = out->findModelEntity(2, 0);
  int t = getType(in.nodesPerElement);
  for (int i = 0; i < in.nodes; ++i)
    v[i] = out->createVert(c);
  for (int i = 0; i < in.nodes; ++i) {
    apf::Vector3 x(in.x[i], in.y[i], in.z[i]);
    out->setPoint(v[i], 0, x);
  }
  for (int i = 0; i < in.elements; ++i) {
    apf::Downward ev;
    for (int j = 0; j < in.nodesPerElement; ++j) {
      int vi = in.conn[i * in.nodesPerElement + j] - 1;
/* zero indices (-1 after conversion) indicate that an MPAS
   vertex does not have an adjacent cell in that direction
   (it is a boundary vertex). In our dual mesh, we simply
   omit that dual element. This is correct and should not
   cause problems. */
      if (vi < 0)
        goto bad_vi;
      assert(vi < in.nodes);
      ev[j] = v[vi];
    }
    apf::buildElement(out, c, t, ev);
bad_vi:
    continue;
  }
  numberInitialNodes(in, out, v);
  delete [] v;
}

void removeIsolatedNodes(apf::Mesh2* m)
{
/* another problem with the dual: MPAS cells without adjacent
   cells become vertices with no adjacent elements in the
   dual mesh, which our codes don't handle.
   We have to give up on these cells, but tell the user about it.
   We also need to be aware of these omitted cells going forward */
  apf::MeshIterator* it = m->begin(0);
  int n = 0;
  apf::MeshEntity* e;
  while ((e = m->iterate(it)))
    if ( ! m->countUpward(e)) {
      m->destroy(e);
      ++n;
    }
  m->end(it);
  if (n)
    fprintf(stderr,
        "warning: removed %d isolated nodes from MPAS dual mesh\n",
        n);
}

void loadMpasMesh(apf::Mesh2* m, const char* filename)
{
  MpasFile f;
  readMpasFile(f, filename);
  buildInitialMesh(f, m);
  removeIsolatedNodes(m);
}


void writeMpasAssignments(apf::Mesh2* m, const char* filename) {
  NcFile in(filename, NcFile::read);
  /* this is the dual of a hexagonal mesh,
     hence the reversing of terms */
  int numMpasVtx = readDim(in, "nCells");
  
  apf::Numbering* n = apf::createNumbering(m, "mpas_id", m->getShape(), 1);

  // collect N/#parts contiguous vertex assignments on each process (and deal with any remainders)
  int numPerPart=numMpasVtx/PCU_Comm_Peers()+1;
  apf::MeshIterator* itr = m->begin(0);
  apf::MeshEntity* e;
  std::map<int,int> vtxs;
  PCU_Comm_Begin();
  while (e=m->iterate(itr)) {
    if (!parma::isOwned(m,e))
      continue;
    int num=getNumber(n,e,0,0);
    int target = num/numPerPart;
    if (target==PCU_Comm_Self())
      vtxs[num]=parma::getOwner(m,e);
    else 
      PCU_COMM_PACK(target,num);
  }

  PCU_Comm_Send();
  while (PCU_Comm_Listen()) {
    int num;
    PCU_COMM_UNPACK(num);
    int owner =PCU_Comm_Sender();
    
    vtxs[num]=owner;
  }
  // assign missing vertices to a random part id
  int count=0;
  for (int i=numPerPart*PCU_Comm_Self();i<numPerPart*(PCU_Comm_Self()+1)&&i<numMpasVtx;i++)
    if (vtxs.find(i)==vtxs.end())  {
      vtxs[i]=0; //to be random
      count++;
    }
  fprintf(stdout,"missing vertices found %d\n",count);
  // use MPI IO to write the contiguous blocks to a single graph.info.part.<#parts> file
  // see https://gist.github.com/cwsmith/166d5beb400f3a8136f7 and the comments
}

}
{'fraction_non_alphanumeric': 0.14899536321483772, 'fraction_numerical': 0.008964451313755796, 'mean_word_length': 3.4058102587380845, 'pattern_counts': {'":': 0, '<': 30, '<?xml version=': 0, '>': 36, 'https://': 2, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'dir': 'git-commits-cleaned', 'id': '1149521', 'n_tokens_mistral': 3745, 'n_tokens_neox': 3660, 'n_words': 1049}