ibg2simplex.c

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

/* last edit: Ilja Schmelzer -------------- 20-OCT-1994 15:49:58.82	*/
/************************************************************************/
/*                                                                      */
/*  <<< I B G >>> - Intersection - Based Grid generation package 	*/
/*                                                                      */
/*  Version 1.1 by Ilja Schmelzer   [e-mail]       */
/*                                                                      */
/*  to be distributed under IBG license conditions (see "readme.ibg")	*/
/*                                                                      */
/************************************************************************/
/*
	This program writes out the grid in the SimplexGrid 1.0 data format
 	described in:

		The Simplex Grid Data Format Version 1.0

	which is part of the IBG documentation and available as

		ibg/doc/simplex.*.

	The selection mechanism (ibgOutputType, ibgOSelectedRegion and so on)
	is declared in "ibgoutput.h". The related data will be handled in the
	following files:

		"ibgoutput.c"		general default initialization
		"ibgexamples.c"		application-specific modification
		"ibgcontrol.c"		pattern for user control of the
					selection parameter (trivial).

        You can also modify this file changing the printf calls to
	create files in your preferred file format.
*/

#define ibgNformatSimplex	"SimplexGrid "
#define ibgNversionSimplex		"1.0"
#include <string.h>
#include <stdio.h>
#include "ibgg.h"
#include "ibgdefault.h"
#include "ibgoutput.h"
#include "ibgapplication.h"
#include "ibglib.h"

#define big1 1.e20

int ibgToSimplex(ibGrid *gg, FILE *file);

ibGrid *ibgOtSimplex(FILE *file);

int ibgToSimplex(ibGrid *gg, FILE *file)
{int gdim,d,ic,cc,cl,cr,u,ul,ur,in,nn,n;
 int *cs,*cn,*nnum,*cnum,i,nmax,creg,cfac,rc;
 ibgCellType t,ts;
 ibgPoint *point;
 rc = ibgSuccess;
/* writes the header of the data file: */
 fprintf(file,ibgNformatSimplex);
 fprintf(file," ");
 fprintf(file,ibgNversionSimplex);
 fprintf(file,"\n");
 fprintf(file,"\n");
 fprintf(file,"%d dimensional grid \n",gdim=gg->gridDimension);
/* creation of the list of all cells and points which have to be written out: */
 cnum = calloc((gg->lastCell+1),sizeof(int));	/* zero initialization */
 nnum = calloc((gg->lastPoint+1),sizeof(int));
 creg = 0;	/* current region number */
 if(ibgOutputType[ibgSRegion]){
	foribgRegionCells(*gg,ic,t,u){
		if(ibgOSelectedRegion(u)){
	 		cnum[ic] = ++creg;
			cn = ibgridCPointList(*gg,ic);
			for(n=0;n<ibgcPoints(t);n++){
				nnum[cn[n]] = 1;	/* marking the points */
			}
		}
	}endibgRegionCells(*gg,ic,t,u);
 }
 cfac = 0;	/* extra counter for region and face cells */
 if(ibgOutputType[ibgSFace]){
	foribgFaceCells(*gg,ic,t,u){
		ul = ibgridCellSegment(*gg,cl = ibgridCellLeft( *gg,ic,t));
		ur = ibgridCellSegment(*gg,cr = ibgridCellRight(*gg,ic,t));
		if(ibgOSelectedFace(u,ul,ur)){
			cnum[ic] = ++cfac;
			cn = ibgridCPointList(*gg,ic);
			for(n=0;n<ibgcPoints(t);n++){
				nnum[cn[n]] = 1;	/* marking the points */
			}
		}
	}endibgFaceCells(*gg,ic,t,u);
 }
 nmax = 0;
 foribgPoints(*gg,in,t,u,point){
 	if(nnum[in]) nnum[in] = ++nmax;		/* enumeration of the points */
 }endibgPoints(*gg,in,t,u,point);
/* now we have computed the number of cells and points we want to write out: */
 fprintf(file,"%d points \n",nmax);
 fprintf(file,"%d cells \n",creg);
 fprintf(file,"%d boundary cells \n",cfac);
/* now we write out the point data: */
 foribgPoints(*gg,in,t,u,point){
       	if(nnum[in]==0) continue;
		/* the current number of the point is nnum[in] */
	for(d=0;d<gdim;d++) fprintf(file,"%f ",(float)(ibgpX(*point)[d]));
	fprintf(file,"\n");
 }endibgPoints(*gg,in,t,u,point);
/* now we write out region cells - here only for triangle and tetrahedron */
 if(creg){
	foribgRegionCells(*gg,ic,t,u){
		if(cnum[ic]==0) continue;
			/* the current reg. cell number in file is cnum[ic] */
		if(gg->gridDimension==2)
			{if(t!=ibgc2Triangle)	{rc=ibgError;continue;}}
		else if(gg->gridDimension==3)
			{if(t!=ibgc3Tetrahedron){rc=ibgError;continue;}}
	/* loop over the points of the cell: */
		cn = ibgridCPointList(*gg,ic);
		for(i=0;i<ibgcPoints(t);i++){
			nn = cn[i];		/* nn is number in grid */
			ibgassert(nn>0);
			nn = nnum[nn];		/* nn is now number in file */
			ibgassert(nn>0);
		 	fprintf(file,"%d ",nn);
		}
	 	fprintf(file,"%d ",u);
	/* loop over the neighbours of the cell: */
		cs = ibgridCSideList(*gg,ic);
		for(i=0;i<ibgcSides(t);i++){
			cc = cs[i];
	/* cc is number of neighbour in the grid. The related number in the
		data file is cnum[cc]. Remark that cnum[ic] is zero
		if the cell is not written out and a face and a region cell
		may have the same cnum[cc].
	*/
			ts = ibgridCellType(*gg,cc);
			if(ibgcCODIM(ts)==ibgSRegion){
			 	fprintf(file,"%d ",cnum[cc]);
			}else if(ibgcCODIM(ts)==ibgSFace){
			 	fprintf(file,"-%d ",cnum[cc]);
			}else{
			 	fprintf(file," 0 ");
				ibgfatal;
			}
		}
		fprintf(file,"\n");
	}endibgRegionCells(*gg,ic,t,u);
 }
/* now we write out the face cells: */
 if(cfac){
	foribgFaceCells(*gg,ic,t,u){
		if(cnum[ic]==0) continue;
			/* the current face cell number in file is cnum[ic] */
		if(gg->gridDimension==2)
			{if(t!=ibgc2Edge)	{rc=ibgError;continue;}}
		else if(gg->gridDimension==3)
			{if(t!=ibgc3Triangle)	{rc=ibgError;continue;}}
	/* loop over the points of the cell: */
		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);
		}
		if(u < ibgSMAX)	fprintf(file,"%d ",u);
		else	       	fprintf(file,"%d ",0);
		cl = ibgridCellLeft( *gg,ic,t);
		cr = ibgridCellRight(*gg,ic,t);
		if(creg){
			if(cnum[cl]) fprintf(file,"%d ",cl);
			else fprintf(file,"-%d ",ibgridCellSegment(*gg,cl));
			if(cnum[cr]) fprintf(file,"%d ",cr);
			else fprintf(file,"-%d ",ibgridCellSegment(*gg,cr));
		}else{
		 	fprintf(file,"-%d ",ibgridCellSegment(*gg,cl));
		 	fprintf(file,"-%d ",ibgridCellSegment(*gg,cr));
		}
		cs = ibgridCSideList(*gg,ic);
		for(i=0;i<ibgcSides(t);i++){
			cc = cs[i]; ts = ibgridCellType(*gg,cc);
			if(ibgcCODIM(ts)==ibgSFace){
			 	fprintf(file,"%d ",cnum[cc]);
			}else if(ibgcCODIM(ts)==ibgSLine){
			 	fprintf(file,"0 ");
			}else	ibgfatal;
		}
		fprintf(file,"\n");
	}endibgFaceCells(*gg,ic,t,u);
 }
 free(nnum);
 free(cnum);
 return rc;
}

/*
This function needs an input file in the data format SimplexGrid 1.0
*/

ibGrid *ibgOtSimplex(FILE *file)
{int ccn,ccs,cc,cb,c,cl,cr,u,ul,ur,n,nn,rc,i,j,k;
 int points,rpoints,fpoints,cells,rcells,fcells;
 int t,rtyp,cpoints,csides,num[9],gdim;
 ibgCellType ftyp;
 float x,y,z;
 ibgFloat xmin = big1,xmax = -big1;
 ibgFloat ymin = big1,ymax = -big1;
 ibgFloat zmin = big1,zmax = -big1;
 ibgPoint *point;
#define buflen 200
 char buffer[buflen];
 ibGrid *gg;
 fgets(buffer,buflen,file);
 if(!strcmp(buffer,ibgNformatSimplex)) return ibgNULL;
 fgets(buffer,buflen,file);
 fgets(buffer,buflen,file);
 rc = sscanf(buffer,"%d",&gdim);	/* the grid dimension */
 if(rc != 1) return ibgNULL;
 fgets(buffer,buflen,file);
 rc = sscanf(buffer,"%d",&points);	/* the number of points */
 if(rc != 1) return ibgNULL;
 fgets(buffer,buflen,file);
 rc = sscanf(buffer,"%d",&rcells);	/* the number of region cells */
 if(rc != 1) return ibgNULL;
 fgets(buffer,buflen,file);
 rc = sscanf(buffer,"%d",&fcells);	/* the number of boundary face cells */
 if(rc != 1) return ibgNULL;
/* memory allocation: */
 gg=malloc(sizeof(ibGrid));
 gg->spaceDimension  	= ibgDIM;	/* 3 coordinates also for 2D grid */
 gg->gridDimension      = gdim;
 gg->maxPoint 	= gg->lastPoint	= points;
 gg->Point	= malloc(IBGPOINTSIZE*(gg->maxPoint+1));
 gg->firstRegion= gg->firstCell	= gg->firstPoint	= 1;
 gg->lastRegion	= rcells;
 gg->firstFace	= rcells+1;
 gg->lastFace	= gg->lastCell	= gg->maxCell	= cells = rcells+fcells;
 gg->firstLine	= gg->firstNode	= 1;
 gg->lastLine	= gg->lastNode	= 0;
 if(ibgridGDIM(*gg)==1){
       	cpoints = 2*rcells + 1;
	csides = 2*rcells + 2*fcells + 1;
	rtyp = ibgc1Edge;	rpoints = 2;
	ftyp = ibgc1Node;	fpoints = 1;