cogversionfull.cxx

Delaunay三角形的网格剖分程序

#define  COG_VERSION COG_VERSION_FULL

#include "cogenchart.hxx"
#include "cogencomposite.hxx"
#include "cogwarnings.hxx"

wzpoints CogenChart::generatePoints(cogeometry g,wzmetric m,wzchart c) const
{
  wzIndex i;
  cogeometry g1 = new CogeometryChart(g,coordinates,(void*) this);
  wzpoints pp = original->generatePoints(g1,m,c);
  wzRangeLoop(*pp,i){
    coordinates->callDirect(pp->Point[i].X);
  }
  return pp;
}

wzpoints CogenComposite::generatePoints(cogeometry g,wzmetric r,wzchart c) const
{
  wzpoints p = new wzPointList;
  wzpoints p1;
  wzFloat fi,fj;
  wzIndex i,j,k,cc,c0=0;
  wzRangeLoop(*this,i){
    cogInfoGenerate("\n  pointlist <> ")<<i;
    p1 = generator[i]->generatePoints(g,
				      generator[i]->refinement(),
				      generator[i]->chart());
    if(p1->length()<=0){
      cogInfoGenerate("empty");
      continue;
    }
    if(p->spaceDim<p1->spaceDim) p->spaceDim=p1->spaceDim;
    if(p->boxDim<p1->boxDim) p->boxDim=p1->boxDim;
    wzArray<wzIndex> newnode(*p1); wzIndex n=0;
    newnode[0] = 0; cc=0;
    wzRangeLoop(*p1,k){
      fi = generator[i]->chi(p1->Point[k]);
      wzRangeLoop(*this,j){
	if(i==j)fj = 0;
     	else	fj = generator[j]->chi(p1->Point[k]);
	if(fj>fi)	goto omit;
	else if(fj==fi){
	  if(j<i) 	goto omit;
	}
      }
      cc++;
      n = p->create();
      p->Point[n] = p1->Point[k];
      p->Near[n] = newnode[p1->Near[k]];
      newnode[k] = n;
      continue;
    omit:
      newnode[k] = 0;
    }
    cogInfoGenerate
      (" used <>")<<cc;
    c0 += cc;
  }
  cogInfoGenerate
      ("\n summary: <> points for Delaunay algorithm\n")<<c0;
  return p;
}