ibgomark.cxx

来自「有限元学习研究用源代码(老外的),供科研人员参考」· CXX 代码 · 共 485 行 · 第 1/2 页

#include "cogenoctree.hxx"
#include "ibgoctree.hxx"
#include "cogwarnings.hxx"

wzFloat ibgOctree::LongSideIntersectionFactor = 0.45;
// wzFloat ibgOctree::LongSideIntersectionFactor = 0.25;

// 	void CogenOctree::shift(ibgoctree gg, wzFloat delta) const;
// Shifts the faces, edges and vertices which are in an absolute
// delta-environment of nodes, lines and planes to the exact position
// of these nodes, lines and planes.
// Hides nodes under faces resp. faces/nodes under edges.

//	void CogenOctree::hideIrregular(ibgoctree gg, cogeometry g) const;
// Hides nodes in boundary environment which are not regular

void CogenOctree::hideLongSideIntersections(ibgoctree gg) const
{
  ibgOctree& g = *gg;
  wzInteger nu,no,nru;	//,nro,noo,nou;
  wzIndex d,u,s,r,dr;
  wzIndex ff,ll,bu,bo; //bru,bro,bou,boo;
  wzFloat dx,dy,dx2,dy2;
  wzRangeLoop(g.faces,ff){
    ll = g.fLine[ff];	 	
    d  = g.lDirection(ll);	u  = g.down(d);
    bo = g.lNodeUp[ll];		no = g.bnBase[bo];
    bu = g.lNodeDown[ll];	nu = g.bnBase[bu];
    dx = g.LongSideIntersectionFactor*(g.nPoint[no][d] - g.nPoint[nu][d]);
    dx2 = g.chart->g_ii(g.fMetric[ff],g.fPoint[ff],d)*dx*dx;
    for(s=0;s<g.LineSides;s++){
      r = g.lsLine[d][s]; dr = g.Up(r);
      if(g.is_defined(nru=g.nNext[r][nu])){
	dy = g.nPoint[nru][dr] - g.nPoint[nu][dr];
	dy2 = g.chart->g_ii(g.fMetric[ff],g.fPoint[ff],dr)*dy*dy;
	if(dy2<dx2) goto OnLongSide;
      }
    }
    continue;
  OnLongSide:
    if(g.fPoint[ff][d]<g.nPoint[nu][d]+0.01*dx) continue;
    if(g.fPoint[ff][d]>g.nPoint[no][d]-0.01*dx) continue;
    g.bnType[bu] = ibgoHidden;
    g.bnType[bo] = ibgoHidden;
    g.fType[ ff] = ibgoHidden;
  }
}

/*

void CogenOctree::hideSkewIntersections(ibgoctree gg) const
{
  ibgOctree& g = *gg;
  wzInteger nu,no,nru,nro,noo,nou;
  wzIndex d,u,s,t,o,r,dr;
  wzIndex ff,ll,bu,bo,bru,bro,bou,boo;
  wzFloat dx,dy;
  wzRangeLoop(g.faces,ff){
    ll = g.fLine[ff];	 	
    d  = g.lDirection(ll);	u  = g.down(d);
    bo = g.lNodeUp[ll];		no = g.bnBase[bo];
    bu = g.lNodeDown[ll];	nu = g.bnBase[bu];
    dx = g.nPoint[no][d]-g.nPoint[nu][d];
    for(s=0;s<g.LineSides;s++){
      r = g.lsLine[d][s]; dr = g.Up(r);
      for(i=0,lr = g.bnNext[r][bu];i<2;i++,lr = g.bnNext[r][bu]){
	if(g.isUp(r)){
	  ff = g.lFaceDown[lr];
	}else{
	  ff = g.lFaceUp[lr];
	}
	if(ff==0) continue;
	df = g.fPoint[ff][dr]-g.nPoint[nu][dr]; if(df<0) df=-df;
	if(

      }
      
      if(g.is_defined(nru=g.nNext[r][nu])){
	dy = g.nPoint[nru][dr] - g.nPoint[nu][dr];
	if(dy<dx  && dy>-dx) goto OnLongSide;
      }
    }
    continue;
  OnLongSide:
    if(g.fPoint[ff][d]<g.nPoint[nu][d]+0.01*dx) continue;
    if(g.fPoint[ff][d]>g.nPoint[no][d]-0.01*dx) continue;
    g.bnType[bu] = ibgoHidden;
    g.bnType[bo] = ibgoHidden;
    g.fType[ ff] = ibgoHidden;
  }
}
	
*/
void CogenOctree::hideIrregular(ibgoctree gg) const
{
  ibgOctree& g = *gg;
  wzInteger nu,no,nru,nro,noo,nou;
  wzIndex d,dr,u,s,t,o,r;
  wzIndex ff,ll,bu,bo,bru,bro,bou,boo;
  wzRangeLoop(g.faces,ff){
    ll = g.fLine[ff];	 	
    d  = g.lDirection(ll);	u  = g.down(d);
    bo = g.lNodeUp[ll];		no = g.bnBase[bo];
    bu = g.lNodeDown[ll];	nu = g.bnBase[bu];
    for(s=0;s<g.LineSides;s++){
      r = g.lsLine[d][s]; dr = g.Up(r);
      if(g.is_nothing(nru=g.nNext[r][nu])){
	g.bnType[bu] = ibgoHidden;
      }else if(g.is_defined(nru)){
	bru = g.nBoundary[nru]; 
	if(!bru) bru = g.createSimpleBoundaryNode(nru);
	if(g.is_nothing(g.nNext[d][nru])){
	  g.bnType[bru] = ibgoHidden;
	}
	o = g.llOrtho[d][r];
	for(t=0;t<g.PlaneSides;t++){
	  if(g.is_nothing(nou=g.nNext[o][nru])){
	    g.bnType[bru] = ibgoHidden;
	  }else if(g.is_defined(nou)){
	    if(g.is_nothing(g.nNext[d][nou])){
	      bou = g.nBoundary[nou];
	      if(!bou) bou = g.createSimpleBoundaryNode(nou);
	      g.bnType[bou] = ibgoHidden;
	    }
	  }
	  o = g.Other(o);
	}
      }
      if(g.is_nothing(nro=g.nNext[r][no])){
	g.bnType[bo] = ibgoHidden;
      }else if(g.is_defined(nro)){
	bro = g.nBoundary[nro];
	if(!bro) bro = g.createSimpleBoundaryNode(nro);
	if(g.is_nothing(g.nNext[u][nro])){
	  g.bnType[bro] = ibgoHidden;
	}
	o = g.llOrtho[d][r];
	for(t=0;t<g.PlaneSides;t++){
	  if(g.is_nothing(noo=g.nNext[o][nro])){
	    g.bnType[bro] = ibgoHidden;
	  }else if(g.is_defined(noo)){
	    boo = g.nBoundary[noo];
	    if(g.is_nothing(g.nNext[u][noo])){
	      if(!boo) boo = g.createSimpleBoundaryNode(noo);
	      g.bnType[boo] = ibgoHidden;
	    }
	  }
	  o = g.Other(o);
	}
      }
    }
  }  
}

const int hideNothing = 0;
const int hideRefined = 1;
const int hideFew = 2;
//const int hideNearlyAll = 3;
const int hideAll = 4;

//static int	ibgBasicMode = hideRefined;
static int	ibgBasicMode = hideFew;
static float	ibgOnLongSideButClose = 0.1;
static float	ibgFaceVeryClose = 0.01;
static float	ibgFaceRelativeClose = 0.1;
static float	ibgFaceRelativeFar   = 0.9;
static float	ibgOtherFaceTooClose = 0.1;
static float	ibgEdgeTooClose = 0.3;
static float	ibgEdgeTooClosePar = 0.5;

void CogenOctreeOldmark(ibgoctree gg, cogeometry cog)
{
  ibgInteger bn,n,i,rc,ri,rb,rx,fl,ff,fr,d,dr,r,ls,ll,e,s,mark,dim,border;
  ibgInteger nu,no,nt,nr,bu,bo,ou,pe;
  const ibgIndex Default=0,FarFaceType=1,FaceType=2,CloseFaceType=3,OnFace=4,RefinedNode=5; 
  ibgFloat dx,dy,df,de,dp,de1,dp1,rt,rf,xf,yf,yr,xi,di,fi,db,fb,xu,xo;
  cogInfoOmittedFaces.report_only();
  cogInfoOmittedNodes.report_only();
  wzRangeLoop(gg->bnodes,bn){
    n = gg->bnBase[bn];
    rc = ri = rb = rx = Default; dim=border=0;
    for(i=0;i<gg->LinesUp;i++){
      if(gg->is_outside(no=gg->nUp  [i][n])){
	border++;
	if((fl=gg->bnDown[i][bn]) && (ff=gg->lFaceUp[fl])){
	  db = (xi=gg->nPoint[n][i])-gg->fPoint[ff][i];
	  wzAssert(db>-0.0001);
	  if(db<cog->Delta){
	    gg->fPoint[ff][i] = gg->nPoint[n][i];
	    gg->fPointUp[ff][i] += db;
	    gg->fPointDown[ff][i] += db;
	    rx = OnFace; continue;
	  }
	  fb = db/(xi-gg->nPoint[gg->nDown[i][n]][i]);
	  if(fb<ibgFaceRelativeClose)	{rb = CloseFaceType;}
	  else if(rb<FaceType){
	    if(fb<ibgFaceRelativeFar)	{rb = FaceType;}
	    else 			{rb = FarFaceType;}
	  }else if(rb<FarFaceType)	{rb = FarFaceType;}
	}
      }else if(gg->is_outside(nu=gg->nDown[i][n])){
	border++;
	if((fl=gg->bnUp[i][bn]) && (ff=gg->lFaceDown[fl])){
	  db = gg->fPoint[ff][i]-(xi=gg->nPoint[n][i]);
	  wzAssert(db>-0.0001);
	  if(db<cog->Delta){
	    gg->fPoint[ff][i] = gg->nPoint[n][i]; 
	    gg->fPointUp[ff][i] -= db;
	    gg->fPointDown[ff][i] -= db;
	    rx = OnFace; continue;
	  }
	  fb = db/(gg->nPoint[gg->nUp[i][n]][i]-xi);
	  if(fb<ibgFaceRelativeClose)	{rb = CloseFaceType;}
	  else if(rb<FaceType){
	    if(fb<ibgFaceRelativeFar)	{rb = FaceType;}
	    else 			{rb = FarFaceType;}
	  }else if(rb<FarFaceType)	{rb = FarFaceType;}
	}
      }else if(gg->is_undefined(no)){
	if(rx==Default) rx = RefinedNode;
      }else if(gg->is_undefined(nu)){
	if(rx==Default) rx = RefinedNode;
      }else{
	dim++;
	if((fl=gg->bnDown[i][bn]) && (ff=gg->lFaceUp[fl])){
	  di = (xi=gg->nPoint[n][i])-gg->fPoint[ff][i];
	  wzAssert(di>-0.0001);
	  if(di<cog->Delta){
	    gg->fPoint[ff][i] = gg->nPoint[n][i];
	    gg->fPointUp[ff][i] += di;
	    gg->fPointDown[ff][i] += di;
	    rx = OnFace; continue;
	  }
	  fi = di/(xi-gg->nPoint[gg->nDown[i][n]][i]);
	  if(fi<ibgFaceRelativeClose)	{ri = CloseFaceType;}
	  else if(ri<FaceType){
	    if(fi<ibgFaceRelativeFar)	{ri = FaceType;}
	    else 			{ri = FarFaceType;}
	  }else if(rb<FarFaceType)	{ri = FarFaceType;}
	}else if((fl=gg->bnUp[i][bn]) && (ff=gg->lFaceDown[fl])){
	  di = gg->fPoint[ff][i]- (xi=gg->nPoint[n][i]);
	  wzAssert(di>-0.0001);
	  if(di<cog->Delta){