sg-dump.cxx

有限元学习研究用源代码(老外的),供科研人员参考

#include <time.h>
#include "ibgold.hxx"
#include "sg-dump.hxx"

int sgWriteSmp(
	       ibGrid *gg,
	       FILE *file,
	       int cmat_max,
	       int *cmat_table,
	       wzRegion(*cmat_function)(const cogPoint&p),
	       int fmat_max,
	       int *fmat_table,
	       wzFace(*fmat_function)(const cogPoint&p)
	       )
{
  int gdim,d,ic,cc,t,u,in,nn,n,*cs,*cn,*nnum,*cnum,i,
  nnodes,ncells,nfaces,mat,rc;
  ibgPoint *point;
  cogPoint pt;
  time_t date,tt;
  int id;
  rc = ibgSuccess;
  date=time(&tt);

  
  fprintf(file,"SimplexGrid");
  fprintf(file," ");
  fprintf(file,"2.0\n");
  fprintf(file,"#created by sg-cog (c) I.Schmelzer,J.Fuhrmann\n");
  fprintf(file,"#$Id: sg-dump.cxx,v 1.2 1997/11/05 15:40:13 fuhrmann Exp $\n");
  fprintf(file,"#%s",ctime(&date));
  fprintf(file,"DIMENSION\n%d\n",gdim=gg->gridDimension);
  printf("I-Data: %d dimensional grid\n",gdim=gg->gridDimension);
  cnum = (int*)calloc((gg->lastCell+1),sizeof(int));
  nnum = (int*)calloc((gg->lastPoint+1),sizeof(int));

  ncells = 0;
  foribgRegionCells(*gg,ic,t,u)
    {
      if (!cmat_table[u]) continue;
      cnum[ic] = ++ncells;
      cn = ibgridCPointList(*gg,ic);
      for(n=0;n<ibgcPoints(t);n++)
	{
	  nnum[cn[n]] = 1;
	}
    }
  endibgRegionCells(*gg,ic,t,u);
  
  nfaces = 0;
  foribgFaceCells(*gg,ic,t,u){
    cc = ibgridCellLeft(*gg,ic,t);
    if(!cmat_table[ibgridCellSegment(*gg,cc)])
      {
	cc = ibgridCellRight(*gg,ic,t);
	if(!cmat_table[ibgridCellSegment(*gg,cc)]) continue;
      }
    cnum[ic] = ++nfaces;
    cn = ibgridCPointList(*gg,ic);
    for(n=0;n<ibgcPoints(t);n++)
      {
	nnum[cn[n]] = 1;
      }
  }
  endibgFaceCells(*gg,ic,t,u);

  nnodes = 0;
  foribgPoints(*gg,in,t,u,point)
    {
      if(nnum[in]) nnum[in] = ++nnodes;
    }
  endibgPoints(*gg,in,t,u,point);

  fprintf(file,"NODES\n%d %d\n",nnodes,gg->gridDimension);
  printf("I-Data: %d nodes\n",nnodes);


  foribgPoints(*gg,in,t,u,point)
     {
       if(nnum[in]==0) continue;
       for(d=0;d<gdim;d++) fprintf(file,"%f ",(float)(ibgpX(*point)[d]));
       fprintf(file,"\n");
     }
  endibgPoints(*gg,in,t,u,point);

  

  if(ncells)
    {
      fprintf(file,"CELLS\n%d\n",ncells);
      printf("I-Data: %d cells\n",ncells);
    
      foribgRegionCells(*gg,ic,t,u){
	if(cnum[ic]==0) continue;
	if(gg->gridDimension==2)
	  {
	    if(t!=ibgc2Triangle)	{rc=ibgError;continue;}
	  }
	else if(gg->gridDimension==3)
	  {
	    if(t!=ibgc3Tetrahedron){rc=ibgError;continue;}
	  }
	cn = ibgridCPointList(*gg,ic);
	for(i=0;i<ibgcPoints(t);i++)
	  {
	    nn = cn[i];
	    ibgassert(nn>0);
	    nn = nnum[nn];
	    ibgassert(nn>0);
	    fprintf(file,"%d ",nn);
	  }


	mat=cmat_table[u];
	if (cmat_function)
	  {
	    for (id=0;id<3;id++) pt.X[id]=0.0;
	    for(i=0;i<ibgcPoints(t);i++)
	      {
		for (id=0;id<gdim;id++)
		  pt.X[id]+=ibgpX(ibgridPoint(*gg,cn[i]))[id];
	      }
	    for (id=0;id<gdim;id++) pt.X[id]/=(double)ibgcPoints(t);
	    mat=cmat_function(pt);

	  }

        fprintf(file,"%d ",mat);

	cs = ibgridCSideList(*gg,ic);
	for(i=0;i<ibgcSides(t);i++)
	  {
	    cc = cs[i];
            t = ibgridCellType(*gg,cc);
	    if(ibgcCODIM(t)==ibgSRegion)
	      {
		fprintf(file,"%d ",cnum[cc]);
	      }
	    else if(ibgcCODIM(t)==ibgSFace)
	      {
		fprintf(file,"-%d ",cnum[cc]);
	      }
	    else	ibgfatal;
	  }
	fprintf(file,"\n");
      }
      endibgRegionCells(*gg,ic,t,u);
    }
  
  if(nfaces)
    {

      fprintf(file,"FACES\n%d\n",nfaces);
      printf("I-Data: %d faces\n",nfaces);
      
      foribgFaceCells(*gg,ic,t,u)
	{
	  if(cnum[ic]==0) continue;
	  if(gg->gridDimension==2)
	    {
	      if(t!=ibgc2Edge)	{rc=ibgError;continue;}
	    }
	  else if(gg->gridDimension==3)
	    {
	      if(t!=ibgc3Triangle)      {rc=ibgError;continue;}
	    }
	  
	  
	  cc = ibgridCellLeft(*gg,ic,t);
	  if(!cmat_table[ibgridCellSegment(*gg,cc)])
	    {
	      cc = ibgridCellRight(*gg,ic,t);
	      if(!cmat_table[ibgridCellSegment(*gg,cc)]) continue;
	    }
	  
	  
	  cn = ibgridCPointList(*gg,ic);
	  for(i=0;i<ibgcPoints(t);i++)
	    {
	      nn = cn[i];
	      ibgassert(nn>0);
	      nn = nnum[nn];
	      ibgassert(nn>0);
	      fprintf(file,"%d ",nn);
	    }
	  
	  mat=0;
	  if(u < fmat_max)
	    mat=fmat_table[u];
	  
	  if (fmat_function)
	    {
	      for (id=0;id<3;id++) pt.X[id]=0.0;
	      for(i=0;i<ibgcPoints(t);i++)
		{
		  for (id=0;id<gdim;id++)
		    pt.X[id]+=ibgpX(ibgridPoint(*gg,cn[i]))[id];
		}
	      for (id=0;id<gdim;id++) pt.X[id]/=(double)ibgcPoints(t);
	      mat=fmat_function(pt);
	    }

	  
	  fprintf(file,"%d ",mat);
	  
	  cc = ibgridCellLeft(*gg,ic,t);
	  fprintf(file,"%d ",cmat_table[ibgridCellSegment(*gg,cc)]);
	  cc = ibgridCellRight(*gg,ic,t);
	  fprintf(file,"%d ",cmat_table[ibgridCellSegment(*gg,cc)]);
	  cs = ibgridCSideList(*gg,ic);
	  for(i=0;i<ibgcSides(t);i++)
	    {
	      cc = cs[i]; t = ibgridCellType(*gg,cc);
	      if(ibgcCODIM(t)==ibgSFace)
		{
		fprintf(file,"%d ",cnum[cc]);
	      }
	    else if(ibgcCODIM(t)==ibgSLine){
	      fprintf(file,"0 ");
	    }
	    else	ibgfatal;
	  }
	fprintf(file,"\n");
      }
      endibgFaceCells(*gg,ic,t,u);
    }
  fprintf(file,"END\n");
  printf("I-Data: end\n",nfaces);
  
  free(nnum);
  free(cnum);
  return rc;
}