qsimpcomp.cpp

来自「算断裂的」· C++ 代码 · 共 1,059 行 · 第 1/3 页

// ------------------------------------------------------------------
// qsimpcomp.cpp
//
// This file contains member functions of classes SimpComp and 
// SimpComp_Under_Construction.  See simpcomp.h for more information.
// ------------------------------------------------------------------
// Author: Stephen A. Vavasis
// Copyright (c) 1999 by Cornell University.  All rights reserved.
// 
// See the accompanying file 'Copyright' for authorship information,
// the terms of the license governing this software, and disclaimers
// concerning this software.
// ------------------------------------------------------------------
// This file is part of the QMG software.  
// Version 2.0 of QMG, release date September 3, 1999.
// ------------------------------------------------------------------

#include "qsimpcomp.h"

extern "C" {
#include <stdio.h>
#include <stdlib.h>
  void cptc_free_mesh(Cptc_MshDesc* mesh);
}



namespace {
  using namespace QMG;
  void* malloc1(size_t sz) {
    char tmp[100];
    void* p = malloc(sz);
    if (p == 0 && sz > 0) {
      sprintf(tmp, "Out of memory on malloc of size %u", sz);
      throw_error(tmp);
    }
    return p;
  }
}



#ifdef _MSC_VER
static void unusedfunc() {
  using namespace QMG;
  map<int, int> q;
  q.insert(pair<int,int>(2,3));
  if (q.begin() != q.end())
    throw_error("");
}
#endif


int 
QMG::SimpComplex::brep_levelsize(int level) const {
#ifdef RANGECHECK
  if (level < 0 || level > mesh_ -> embedded_dimension)
    throw_error("levelsize arg out of range");
#endif
  if (level == 0)
    return mesh_ -> num_top_vertices;
  else if (level == 1)
    return mesh_ -> num_top_edges;
  else if (level == 2)
    return mesh_ -> num_top_surfaces;
  else
    return mesh_ -> num_top_regions;
}

QMG::string
QMG::SimpComplex::lookup_prop_val(const string& prop) const {
  for (int i = 0; i < mesh_ -> num_prop_val; ++i) {
    if (compare_nocase(prop, string(mesh_ -> props[i])) == 0)
      return string(mesh_ -> vals[i]);
  }
  return string("");
}


QMG::Real 
QMG::SimpComplex::real_coord_o(VertexOrdinalIndex vnumo, int jj) const {
#ifdef RANGECHECK
  if (jj < 0 || jj >= mesh_ -> embedded_dimension)
    throw_error("Range check");
  return mesh_ -> nodes[process_ordinal_vnum_(vnumo)].coord[jj];
#else
  return mesh_ -> nodes[vnumo].coord[jj];
#endif
}

void 
QMG::SimpComplex::check_fspec_(const Brep::Face_Spec& fspec) const {
  int fdim = fspec.fdim();
  Brep::FaceIndex faceind = fspec.faceind();
  if (fdim < 0 || fdim > mesh_ -> intrinsic_dimension || 
    fdim >= mesh_ -> embedded_dimension)
    throw_error("Fspec dim out of range");
  if (faceind < 0)
    throw_error("Face ind out of range");
  if (fdim == 0) {
    if (faceind >= mesh_ -> num_top_vertices)
      throw_error("Face ind out of range");
  }
  else if (fdim == 1) {
    if (faceind >= mesh_ -> num_top_edges)
      throw_error("Face ind out of range");
  }
  else {
    if (faceind >= mesh_ -> num_top_surfaces)
      throw_error("Face ind out of range");
  }
}

void
QMG::SimpComplex::check_seqno_(const Brep::Face_Spec& fspec, int seqno) const {
  if (seqno < 0)
    throw_error("seqno out of range");
  int fdim = fspec.fdim();
  Brep::FaceIndex faceind = fspec.faceind();
  if (fdim == 0) {
    if (seqno >= static_cast<int>(mesh_ -> vtxinfo[faceind].vtx >= 0))
      throw_error("seqno out of range");
  } 
  else if (fdim == 1) {
    if (seqno >= mesh_ -> edgeinfo[faceind].num_mnode)
      throw_error("seqno out of range");
  } 
  else { 
    if (seqno >= mesh_ -> surfinfo[faceind].num_mnode)
      throw_error("seqno out of range");
  } 
}


int 
QMG::SimpComplex::num_node_on_face(const Brep::Face_Spec& fspec) const {
  int fdim = fspec.fdim();
  Brep::FaceIndex faceind = fspec.faceind();
#ifdef RANGECHECK
  check_fspec2_(fspec);
#endif
  if (fdim == mesh_ -> embedded_dimension)
    return 0;
  else if (fdim == 0)
    return static_cast<int>(mesh_ -> vtxinfo[faceind].vtx >= 0);
  else if (fdim == 1)
    return mesh_ -> edgeinfo[faceind].num_mnode;
  else if (fdim == 2)
    return mesh_ -> surfinfo[faceind].num_mnode;
  return 0;
}

QMG::SimpComplex::VertexGlobalIndex 
QMG::SimpComplex::node_on_face(const Brep::Face_Spec& fspec, int seqno) const {
  int fdim = fspec.fdim();
  Brep::FaceIndex faceind = fspec.faceind();
#ifdef RANGECHECK
  check_fspec_(fspec);
  check_seqno_(fspec, seqno);
#endif
  if (fdim == 0) {
    return mesh_ -> vtxinfo[faceind].vtx;
  }
  else if (fdim == 1) {
    return mesh_ -> edgeinfo[faceind].incs[seqno].id;
  }
  else {
    return mesh_ -> surfinfo[faceind].incs[seqno].id;
  }
}


QMG::Brep::PatchIndex 
QMG::SimpComplex::patchind_on_face(const Brep::Face_Spec& fspec, int seqno) const {
  int fdim = fspec.fdim();
  Brep::FaceIndex faceind = fspec.faceind();
#ifdef RANGECHECK
  check_fspec_(fspec);
  check_seqno_(fspec, seqno);
#endif
  if (fdim == 0) {
    return 0;
  }
  else if (fdim == 1) {
    return mesh_ -> edgeinfo[faceind].incs[seqno].curve_index;
  }
  else {
    return mesh_ -> surfinfo[faceind].incs[seqno].patch_index;
  }
}

QMG::Real 
QMG::SimpComplex::param_coord_on_face(const Brep::Face_Spec& fspec, int seqno, int jj) const {
  int fdim = fspec.fdim();
  Brep::FaceIndex faceind = fspec.faceind();
#ifdef RANGECHECK
  check_fspec_(fspec);
  check_seqno_(fspec, seqno);
  if (jj < 0 || jj >= fdim)
    throw_error("param coord index out of range");
#endif
  if (fdim == 1) {
    return mesh_ -> edgeinfo[faceind].incs[seqno].param;
  }
  else {
    return mesh_ -> surfinfo[faceind].incs[seqno].param[jj];
  }
}

void 
QMG::SimpComplex::check_fspec2_(const Brep::Face_Spec& fspec) const {
  int fdim = fspec.fdim();
  Brep::FaceIndex faceind = fspec.faceind();
  if (fdim < 0 || fdim > mesh_ -> intrinsic_dimension)
    throw_error("Fspec dim out of range");
  if (faceind < 0)
    throw_error("Face ind out of range");
  if (fdim == 0) {
    if (faceind >= mesh_ -> num_top_vertices)
      throw_error("Face ind out of range");
  }
  else if (fdim == 1) {
    if (faceind >= mesh_ -> num_top_edges)
      throw_error("Face ind out of range");
  }
  else if (fdim == 2) {
    if (faceind >= mesh_ -> num_top_surfaces)
      throw_error("Face ind out of range");
  }
  else {
    if (faceind >= mesh_ -> num_top_regions)
      throw_error("Face ind out of range");
  }
}



int 
QMG::SimpComplex::num_meshface_on_face(const Brep::Face_Spec& fspec) const {
  int fdim = fspec.fdim();
  Brep::FaceIndex faceind = fspec.faceind();
#ifdef RANGECHECK
  check_fspec2_(fspec);
#endif
  if (fdim == 0)    
    return static_cast<int>(mesh_ -> vtxinfo[faceind].vtx >= 0);
  else if (fdim == 1)
    return mesh_ -> edgeinfo[faceind].num_medge;
  else if (fdim == 2)
    return mesh_ -> surfinfo[faceind].num_mtrifacet;
  else
    return mesh_ -> regioninfo[faceind].num_tet;
}

QMG::SimpComplex::VertexGlobalIndex 
QMG::SimpComplex::node_of_meshface_on_face(const Brep::Face_Spec& fspec, 
    int seqno, 
    int jj) const {
  int fdim = fspec.fdim();
  Brep::FaceIndex faceind = fspec.faceind();
#ifdef RANGECHECK
  int nm1 = num_meshface_on_face(fspec);
  if (seqno < 0 || seqno >= nm1)
    throw_error("seqno out of range");
  if (jj < 0 || jj > fspec.fdim())
    throw_error("Vertex selector out of range");
#endif
  if (fdim == 0)
    return mesh_ -> vtxinfo[faceind].vtx;
  else if (fdim == 1) 
    return mesh_ -> edgeinfo[faceind].edgelist[seqno][jj];
  else if (fdim == 2)
    return mesh_ -> surfinfo[faceind].trifacetlist[seqno][jj];
  else {
    return mesh_ -> regioninfo[faceind].tets[seqno][jj];
   // map<int,int> q;
    //q.insert(pair<int,int>(2,3));
    //return q.begin() == q.end();
  }
}




void
QMG::SimpComplex_From_FrontEnd::make_global_to_ordinal_() const {
  if (global_to_ordinal_) throw_error("Make global to ordinal called twice");
  global_to_ordinal_ = new GlobalToOrdinalMapType;
  for (VertexOrdinalIndex vnumo = 0;
  vnumo < mesh_ -> num_nodes;
  ++vnumo) {
    VertexGlobalIndex vnum = mesh_ -> nodes[vnumo].id;
    if (vnum < 0)
      throw_error("global vertex id out of range");
    pair<GlobalToOrdinalMapType::iterator, bool> rval =
      global_to_ordinal_ -> 
      insert(pair<VertexGlobalIndex, VertexOrdinalIndex>(vnum,vnumo));
    if (!rval.second)
      throw_error("global vertex id occurs twice in input");
  }
}


QMG::SimpComplex_From_FrontEnd::~SimpComplex_From_FrontEnd() {
  delete global_to_ordinal_;
}



QMG::SimpComplex_From_FrontEnd::
SimpComplex_From_FrontEnd(SimpComplex_From_FrontEnd& other) :
SimpComplex(other.mesh_, other.global_to_ordinal_) {
  other.global_to_ordinal_ = 0;
}


QMG::SimpComplex_From_FrontEnd::
SimpComplex_From_FrontEnd(const SimpComplex_From_FrontEnd_Returnval& sucr) :
SimpComplex(sucr.sref_.mesh_, sucr.sref_.global_to_ordinal_) {
  sucr.sref_.global_to_ordinal_ = 0;
}




bool 
QMG::SimpComplex_Under_Construction::Face_Bdry_Inc_::
operator<(const Face_Bdry_Inc_& other) const {
  if (meshv[0] < other.meshv[0]) return true;
  if (meshv[0] > other.meshv[0]) return false;
  if (meshv[1] < other.meshv[1]) return true;
  if (meshv[1] > other.meshv[1]) return false;
  return meshv[2] < other.meshv[2];
}

bool
QMG::SimpComplex_Under_Construction::Vertex_Bdry_Inc_::
operator<(const Vertex_Bdry_Inc_& other) const {
  if (vnum < other.vnum) return true;
  if (vnum > other.vnum) return false;
  return fspec < other.fspec;
}



void 
QMG::SimpComplex_Under_Construction::clobber_() {
  mesh_ = 0;
  global_to_ordinal_ = 0;
  vertinclookup_ = 0;  
  faceinclookup_ = 0;
  max_node_edge_.resize(0);
  max_edge_edge_.resize(0);