void pVVECT(int *key, float *V) { MULTIGRID *mg; ELEMENT *e; VERTEX *v; PreprocessingProcPtr pre; ElementVectorProcPtr foo; double *cc[MAX_CORNERS_OF_ELEM], lc[3], vv[3]; int i, k; k = (*key)-1; mg = GetCurrentMultigrid(); pre = eval[k].v->PreprocessProc; if (pre != NULL) pre(eval_name[k], mg); foo = eval[k].v->EvalProc; mg = GetCurrentMultigrid(); ClearVertexMarkers(mg); SURFACE_LOOP_BEGIN(mg, e) for (i = 0; i < CORNERS_OF_ELEM(e); i++) cc[i] = CVECT(MYVERTEX(CORNER(e, i))); for (i = 0; i < CORNERS_OF_ELEM(e); i++) { v = MYVERTEX(CORNER(e, i)); if (USED(v)) continue; SETUSED(v, 1); LocalCornerCoordinates(3, TAG(e), i, lc); foo(e, cc, lc, vv); *V++ = (float)vv[0]; *V++ = (float)vv[1]; *V++ = (float)vv[2]; } SURFACE_LOOP_END }
static void CenterOfMass (ELEMENT *e, DOUBLE *pos) { int i; V_DIM_CLEAR(pos) for(i=0; i<CORNERS_OF_ELEM(e); i++) { V_DIM_LINCOMB(1.0,pos,1.0,CVECT(MYVERTEX(CORNER(e,i))),pos) } V_DIM_SCALE(1.0/(float)CORNERS_OF_ELEM(e),pos) }
static INT ComputeSurfaceGridStats (MULTIGRID *theMG, COVISE_HEADER *covise) { NODE *theNode; INT l, n_vertices, n_elem, n_conn; /* surface grid up to current level */ n_vertices = n_elem = n_conn = 0; for (l=covise->min_level; l<=covise->max_level; l++) { ELEMENT *theElement; GRID *theGrid = GRID_ON_LEVEL(theMG,l); /* reset USED flags in all vertices to be counted */ for (theNode=FIRSTNODE(theGrid); theNode!=NULL; theNode=SUCCN(theNode)) { SETUSED(MYVERTEX(theNode),0); } /* count geometric objects */ for (theElement=FIRSTELEMENT(theGrid); theElement!=NULL; theElement=SUCCE(theElement)) { if ((EstimateHere(theElement)) || (l==covise->max_level)) { int i, coe = CORNERS_OF_ELEM(theElement); n_elem++; for (i=0; i<coe; i++) { theNode = CORNER(theElement,i); n_conn++; if (USED(MYVERTEX(theNode))) continue; SETUSED(MYVERTEX(theNode),1); if ((SONNODE(theNode)==NULL) || (l==covise->max_level)) { #ifdef ModelP if (PRIO(theNode) == PrioMaster) #endif n_vertices++; } } } } } #ifdef ModelP n_vertices = UG_GlobalSumINT(n_vertices); n_elem = UG_GlobalSumINT(n_elem); n_conn = UG_GlobalSumINT(n_conn); #endif covise->n_vertices = n_vertices; covise->n_elems = n_elem; covise->n_conns = n_conn; return(0); }
void pVSTRUC(int *knode, int *kequiv, int *kcel1, int *kcel2, int *kcel3, int *kcel4, int *knptet, int *kptet, int *knblock, int *blocks, int *kphedra, int *ksurf, int *knsurf, int *hint) { MULTIGRID *mg; ELEMENT *e; VERTEX *v; int i; *knode = *kcel1 = *kcel2 = *kcel3 = *kcel4 = *ksurf = 0; mg = GetCurrentMultigrid(); ClearVertexMarkers(mg); SURFACE_LOOP_BEGIN(mg, e) switch (TAG(e)) { case TETRAHEDRON : (*kcel1)++; break; case PYRAMID : (*kcel2)++; break; case PRISM : (*kcel3)++; break; case HEXAHEDRON : (*kcel4)++; } for (i = 0; i < CORNERS_OF_ELEM(e); i++) { v = MYVERTEX(CORNER(e, i)); if (USED(v)) continue; SETUSED(v, 1); ID(v) = *knode; /* number vertices */ (*knode)++; } /* check for domain boundary sides */ if (OBJT(e) == BEOBJ) { for (i = 0; i < SIDES_OF_ELEM(e); i++) if (SIDE_ON_BND(e, i) && !InnerBoundary(e, i)) (*ksurf)++; } SURFACE_LOOP_END *kequiv = 0; *knptet = 0; *kptet = 0; *knblock = 0; *kphedra = 0; *knsurf = 1; }
void pVGRID(float *xyz) { MULTIGRID *mg; ELEMENT *e; VERTEX *v; int i; mg = GetCurrentMultigrid(); ClearVertexMarkers(mg); SURFACE_LOOP_BEGIN(mg, e) for (i = 0; i < CORNERS_OF_ELEM(e); i++) { v = MYVERTEX(CORNER(e, i)); if (USED(v)) continue; SETUSED(v, 1); *xyz++ = XC(v); *xyz++ = YC(v); *xyz++ = ZC(v); } SURFACE_LOOP_END }
static INT SendSurfaceGrid (MULTIGRID *theMG, COVISE_HEADER *covise) { INT l, remaining, sent; TokenBuffer tb; Message* msg = new Message; msg->type = (covise_msg_type)0; printf("CoviseIF: SendSurfaceGrid start\n"); /* renumber vertex IDs */ /* TODO doesn't work in ModelP */ RenumberVertices(theMG, covise->min_level, covise->max_level); /* send surface vertices, part1: set flags */ ResetVertexFlags(theMG, covise->min_level, covise->max_level); /* send surface vertices, part2: send data */ sent = 0; /* start first buffer */ tb.reset(); remaining = MIN(covise->n_vertices,MAX_ITEMS_SENT); tb << MT_UGGRIDV; tb << remaining; for (l=covise->min_level; l<=covise->max_level; l++) { NODE *theNode; GRID *theGrid = GRID_ON_LEVEL(theMG,l); for (theNode=FIRSTNODE(theGrid); theNode!=NULL; theNode=SUCCN(theNode)) { INT vid; DOUBLE *pos; if (USED(MYVERTEX(theNode))) continue; SETUSED(MYVERTEX(theNode),1); /* extract data from vertex */ /* TODO use VXGID in ModelP */ vid = ID(MYVERTEX(theNode)); pos = CVECT(MYVERTEX(theNode)); tb << (INT32)vid; tb << (FLOAT32)pos[0]; tb << (FLOAT32)pos[1]; tb << (FLOAT32)pos[2]; remaining--; sent++; if (remaining==0) { /* send this buffer */ msg->data = (char*)tb.get_data(); msg->length = tb.get_length(); covise_connection->send_msg(msg); /* start next buffer */ tb.reset(); tb << MT_UGGRIDV; remaining = MIN(covise->n_vertices - sent, MAX_ITEMS_SENT); tb << remaining; } } } printf("CoviseIF: SendSurfaceGrid ...\n"); /* next buffer */ tb.reset(); tb << MT_UGGRIDE; remaining = MIN(covise->n_elems,MAX_ITEMS_SENT); tb << remaining; sent = 0; /* send surface elems and connectivity */ for (l=covise->min_level; l<=covise->max_level; l++) { ELEMENT *theElement; GRID *theGrid = GRID_ON_LEVEL(theMG,l); for (theElement=FIRSTELEMENT(theGrid); theElement!=NULL; theElement=SUCCE(theElement)) { if ((EstimateHere(theElement)) || (l==covise->max_level)) { int i, coe = CORNERS_OF_ELEM(theElement); tb << (INT32)coe; for (i=0; i<coe; i++) { NODE *theNode = CORNER(theElement,i); INT vid; vid = ID(MYVERTEX(theNode)); /* TODO use VXGID in ModelP */ tb << (INT32)vid; } remaining--; sent++; if (remaining==0) { /* send this buffer */ msg->data = (char*)tb.get_data(); msg->length = tb.get_length(); covise_connection->send_msg(msg); /* start next buffer */ tb.reset(); tb << MT_UGGRIDE; remaining = MIN(covise->n_elems - sent, MAX_ITEMS_SENT); tb << remaining; } } } } /* cleanup */ delete msg; printf("CoviseIF: SendSurfaceGrid stop\n"); return(0); }
void NS_DIM_PREFIX SetGhostObjectPriorities (GRID *theGrid) { ELEMENT *theElement,*theNeighbor,*SonList[MAX_SONS]; NODE *theNode; EDGE *theEdge; VECTOR *theVector; INT i,prio,*proclist,hghost,vghost; /* reset USED flag for objects of ghostelements */ for (theElement=PFIRSTELEMENT(theGrid); theElement!=NULL; theElement=SUCCE(theElement)) { SETUSED(theElement,0); SETTHEFLAG(theElement,0); for (i=0; i<EDGES_OF_ELEM(theElement); i++) { theEdge = GetEdge(CORNER(theElement,CORNER_OF_EDGE(theElement,i,0)), CORNER(theElement,CORNER_OF_EDGE(theElement,i,1))); ASSERT(theEdge != NULL); SETUSED(theEdge,0); SETTHEFLAG(theEdge,0); } if (VEC_DEF_IN_OBJ_OF_GRID(theGrid,SIDEVEC)) for (i=0; i<SIDES_OF_ELEM(theElement); i++) { theVector = SVECTOR(theElement,i); if (theVector != NULL) { SETUSED(theVector,0); SETTHEFLAG(theVector,0); } } } /* to reset also nodes which are at corners of the boundary */ /* reset of nodes need to be done through the node list */ for (theNode=PFIRSTNODE(theGrid); theNode!=NULL; theNode=SUCCN(theNode)) { SETUSED(theNode,0); SETTHEFLAG(theNode,0); SETMODIFIED(theNode,0); } /* set FLAG for objects of horizontal and vertical overlap */ for (theElement=PFIRSTELEMENT(theGrid); theElement!=NULL; theElement=SUCCE(theElement)) { if (PARTITION(theElement) == me) continue; /* check for horizontal ghost */ hghost = 0; for (i=0; i<SIDES_OF_ELEM(theElement); i++) { theNeighbor = NBELEM(theElement,i); if (theNeighbor == NULL) continue; if (PARTITION(theNeighbor) == me) { hghost = 1; break; } } /* check for vertical ghost */ vghost = 0; GetAllSons(theElement,SonList); for (i=0; SonList[i]!=NULL; i++) { if (PARTITION(SonList[i]) == me) { vghost = 1; break; } } /* one or both of vghost and hghost should be true here */ /* except for elements which will be disposed during Xfer */ if (vghost) SETTHEFLAG(theElement,1); if (hghost) SETUSED(theElement,1); for (i=0; i<CORNERS_OF_ELEM(theElement); i++) { theNode = CORNER(theElement,i); if (vghost) SETTHEFLAG(theNode,1); if (hghost) SETUSED(theNode,1); } for (i=0; i<EDGES_OF_ELEM(theElement); i++) { theEdge = GetEdge(CORNER_OF_EDGE_PTR(theElement,i,0), CORNER_OF_EDGE_PTR(theElement,i,1)); ASSERT(theEdge != NULL); if (vghost) SETTHEFLAG(theEdge,1); if (hghost) SETUSED(theEdge,1); } if (VEC_DEF_IN_OBJ_OF_GRID(theGrid,SIDEVEC)) for (i=0; i<SIDES_OF_ELEM(theElement); i++) { theVector = SVECTOR(theElement,i); if (theVector != NULL) { if (vghost) SETTHEFLAG(theVector,1); if (hghost) SETUSED(theVector,1); } } } DEBUG_TIME(0); /* set USED flag for objects of master elements */ /* reset FLAG for objects of master elements */ for (theElement=PFIRSTELEMENT(theGrid); theElement!=NULL; theElement=SUCCE(theElement)) { if (PARTITION(theElement) != me) continue; SETUSED(theElement,0); SETTHEFLAG(theElement,0); for (i=0; i<CORNERS_OF_ELEM(theElement); i++) { theNode = CORNER(theElement,i); SETUSED(theNode,0); SETTHEFLAG(theNode,0); SETMODIFIED(theNode,1); } for (i=0; i<EDGES_OF_ELEM(theElement); i++) { theEdge = GetEdge(CORNER_OF_EDGE_PTR(theElement,i,0), CORNER_OF_EDGE_PTR(theElement,i,1)); ASSERT(theEdge != NULL); SETUSED(theEdge,0); SETTHEFLAG(theEdge,0); } if (VEC_DEF_IN_OBJ_OF_GRID(theGrid,SIDEVEC)) for (i=0; i<SIDES_OF_ELEM(theElement); i++) { theVector = SVECTOR(theElement,i); if (theVector != NULL) { SETUSED(theVector,0); SETTHEFLAG(theVector,0); } } } DEBUG_TIME(0); /* set object priorities for ghostelements */ for (theElement=PFIRSTELEMENT(theGrid); theElement!=NULL; theElement=SUCCE(theElement)) { if (PARTITION(theElement) == me) continue; if (USED(theElement) || THEFLAG(theElement)) { prio = PRIO_CALC(theElement); PRINTDEBUG(gm,1,("SetGhostObjectPriorities(): e=" EID_FMTX " new prio=%d\n", EID_PRTX(theElement),prio)) SETEPRIOX(theElement,prio); if (VEC_DEF_IN_OBJ_OF_GRID(theGrid,ELEMVEC)) { theVector = EVECTOR(theElement); if (theVector != NULL) SETPRIOX(theVector,prio); } } /* set edge priorities */ for (i=0; i<EDGES_OF_ELEM(theElement); i++) { theEdge = GetEdge(CORNER(theElement,CORNER_OF_EDGE(theElement,i,0)), CORNER(theElement,CORNER_OF_EDGE(theElement,i,1))); ASSERT(theEdge != NULL); if (USED(theEdge) || THEFLAG(theEdge)) { PRINTDEBUG(dddif,3,(PFMT " dddif_SetGhostObjectPriorities():" " downgrade edge=" EDID_FMTX " from=%d to PrioHGhost\n", me,EDID_PRTX(theEdge),prio)); EDGE_PRIORITY_SET(theGrid,theEdge,PRIO_CALC(theEdge)); } else EDGE_PRIORITY_SET(theGrid,theEdge,PrioMaster); } #ifdef __THREEDIM__ /* if one(all) of the side nodes is (are) a hghost (vghost) node */ /* then its a hghost (vghost) side vector */ if (VEC_DEF_IN_OBJ_OF_GRID(theGrid,SIDEVEC)) for (i=0; i<SIDES_OF_ELEM(theElement); i++) { if (USED(theVector) || THEFLAG(theVector)) SETPRIOX(theVector,PRIO_CALC(theVector)); } #endif } /* to set also nodes which are at corners of the boundary */ /* set them through the node list */ for (theNode=PFIRSTNODE(theGrid); theNode!=NULL; theNode=SUCCN(theNode)) { /* check if its a master node */ if (USED(theNode) || THEFLAG(theNode)) { PRINTDEBUG(dddif,3,(PFMT " dddif_SetGhostObjectPriorities():" " downgrade node=" ID_FMTX " from=%d to PrioHGhost\n", me,ID_PRTX(theNode),prio)); /* set node priorities of node to ghost */ NODE_PRIORITY_SET(theGrid,theNode,PRIO_CALC(theNode)) } else if (MODIFIED(theNode) == 0)
static INT TecplotCommand (INT argc, char **argv) { INT i,j,k,v; /* counters etc. */ INT counter; /* for formatting output */ char item[1024],it[256]; /* item buffers */ INT ic=0; /* item length */ VECTOR *vc; /* a vector pointer */ ELEMENT *el; /* an element pointer */ MULTIGRID *mg; /* our multigrid */ char filename[NAMESIZE]; /* file name for output file */ PFILE *pf; /* the output file pointer */ INT nv; /* number of variables (eval functions) */ EVALUES *ev[MAXVARIABLES]; /* pointers to eval function descriptors */ char ev_name[MAXVARIABLES][NAMESIZE]; /* names for eval functions */ char s[NAMESIZE]; /* name of eval proc */ char zonename[NAMESIZE+7] = ""; /* name for zone (initialized to empty string) */ INT numNodes; /* number of data points */ INT numElements; /* number of elements */ INT gnumNodes; /* number of data points globally */ INT gnumElements; /* number of elements globallay */ PreprocessingProcPtr pre; /* pointer to prepare function */ ElementEvalProcPtr eval; /* pointer to evaluation function */ DOUBLE *CornersCoord[MAX_CORNERS_OF_ELEM]; /* pointers to coordinates */ DOUBLE LocalCoord[DIM]; /* is one of the corners local coordinates */ DOUBLE local[DIM]; /* local coordinate in DOUBLE */ DOUBLE value; /* returned by user eval proc */ INT oe,on; INT saveGeometry; /* save geometry flag */ /* get current multigrid */ mg = GetCurrentMultigrid(); if (mg==NULL) { PrintErrorMessage('W',"tecplot","no multigrid open\n"); return (OKCODE); } /* scan options */ nv = 0; saveGeometry = 0; for(i=1; i<argc; i++) { switch(argv[i][0]) { case 'e' : /* read eval proc */ if (nv>=MAXVARIABLES) { PrintErrorMessage('E',"tecplot","too many variables specified\n"); break; } sscanf(argv[i],"e %s", s); ev[nv] = GetElementValueEvalProc(s); if (ev[nv]==NULL) { PrintErrorMessageF('E',"tecplot","could not find eval proc %s\n",s); break; } if (sscanf(argv[i+1],"s %s", s) == 1) { strcpy(ev_name[nv],s); i++; } else strcpy(ev_name[nv],ev[nv]->v.name); nv++; break; case 'z' : sscanf(argv[i],"z %s", zonename+3); memcpy(zonename, "T=\"", 3); memcpy(zonename+strlen(zonename), "\", \0", 4); break; case 'g' : sscanf(argv[i],"g %d", &saveGeometry); if (saveGeometry<0) saveGeometry=0; if (saveGeometry>1) saveGeometry=1; break; } } if (nv==0) UserWrite("tecplot: no variables given, printing mesh data only\n"); /* get file name and open output file */ if (sscanf(argv[0],expandfmt(CONCAT3(" tecplot %",NAMELENSTR,"[ -~]")),filename)!=1) { PrintErrorMessage('E',"tecplot","could not read name of logfile"); return(PARAMERRORCODE); } pf = pfile_open(filename); if (pf==NULL) return(PARAMERRORCODE); /********************************/ /* TITLE */ /********************************/ ic = 0; sprintf(it,"TITLE = \"UG TECPLOT OUTPUT\"\n"); strcpy(item+ic,it); ic+=strlen(it); sprintf(it,"VARIABLES = \"X\", \"Y\""); strcpy(item+ic,it); ic+=strlen(it); if (DIM==3) { sprintf(it,", \"Z\""); strcpy(item+ic,it); ic+=strlen(it); } for (i=0; i<nv; i++) { sprintf(it,", \"%s\"",ev[i]->v.name); strcpy(item+ic,it); ic+=strlen(it); } sprintf(it,"\n"); strcpy(item+ic,it); ic+=strlen(it); pfile_master_puts(pf,item); ic=0; /********************************/ /* compute sizes */ /********************************/ /* clear VCFLAG on all levels */ for (k=0; k<=TOPLEVEL(mg); k++) for (vc=FIRSTVECTOR(GRID_ON_LEVEL(mg,k)); vc!=NULL; vc=SUCCVC(vc)) SETVCFLAG(vc,0); /* run thru all levels of elements and set index */ numNodes = numElements = 0; for (k=0; k<=TOPLEVEL(mg); k++) for (el=FIRSTELEMENT(GRID_ON_LEVEL(mg,k)); el!=NULL; el=SUCCE(el)) { if (!EstimateHere(el)) continue; /* process finest level elements only */ numElements++; /* increase element counter */ for (i=0; i<CORNERS_OF_ELEM(el); i++) { vc = NVECTOR(CORNER(el,i)); if (VCFLAG(vc)) continue; /* we have this one already */ VINDEX(vc) = ++numNodes; /* number of data points, begins with 1 ! */ SETVCFLAG(vc,1); /* tag vector as visited */ } } #ifdef ModelP gnumNodes = TPL_GlobalSumINT(numNodes); gnumElements = TPL_GlobalSumINT(numElements); on=get_offset(numNodes); oe=get_offset(numElements); /* clear VCFLAG on all levels */ for (k=0; k<=TOPLEVEL(mg); k++) for (vc=FIRSTVECTOR(GRID_ON_LEVEL(mg,k)); vc!=NULL; vc=SUCCVC(vc)) SETVCFLAG(vc,0); /* number in unique way */ for (k=0; k<=TOPLEVEL(mg); k++) for (el=FIRSTELEMENT(GRID_ON_LEVEL(mg,k)); el!=NULL; el=SUCCE(el)) { if (!EstimateHere(el)) continue; /* process finest level elements only */ for (i=0; i<CORNERS_OF_ELEM(el); i++) { vc = NVECTOR(CORNER(el,i)); if (VCFLAG(vc)) continue; /* we have this one already */ VINDEX(vc) += on; /* add offset */ SETVCFLAG(vc,1); /* tag vector as visited */ } } #else gnumNodes = numNodes; gnumElements = numElements; oe=on=0; #endif /********************************/ /* write ZONE data */ /* uses FEPOINT for data */ /* uses QUADRILATERAL in 2D */ /* and BRICK in 3D */ /********************************/ /* write zone record header */ if (DIM==2) sprintf(it,"ZONE %sN=%d, E=%d, F=FEPOINT, ET=QUADRILATERAL\n", zonename, gnumNodes,gnumElements); if (DIM==3) sprintf(it,"ZONE %sN=%d, E=%d, F=FEPOINT, ET=BRICK\n", zonename, gnumNodes,gnumElements); strcpy(item+ic,it); ic+=strlen(it); pfile_master_puts(pf,item); ic=0; /* write data in FEPOINT format, i.e. all variables of a node per line*/ for (k=0; k<=TOPLEVEL(mg); k++) for (vc=FIRSTVECTOR(GRID_ON_LEVEL(mg,k)); vc!=NULL; vc=SUCCVC(vc)) SETVCFLAG(vc,0); /* clear all flags */ counter=0; for (k=0; k<=TOPLEVEL(mg); k++) for (el=FIRSTELEMENT(GRID_ON_LEVEL(mg,k)); el!=NULL; el=SUCCE(el)) { if (!EstimateHere(el)) continue; /* process finest level elements only */ for (i=0; i<CORNERS_OF_ELEM(el); i++) CornersCoord[i] = CVECT(MYVERTEX(CORNER(el,i))); /* x,y,z of corners */ for (i=0; i<CORNERS_OF_ELEM(el); i++) { vc = NVECTOR(CORNER(el,i)); if (VCFLAG(vc)) continue; /* we have this one alre ady */ SETVCFLAG(vc,1); /* tag vector as visited */ sprintf(it,"%g",(double)XC(MYVERTEX(CORNER(el,i)))); strcpy(item+ic,it); ic+=strlen(it); sprintf(it," %g",(double)YC(MYVERTEX(CORNER(el,i)))); strcpy(item+ic,it); ic+=strlen(it); if (DIM == 3) { sprintf(it," %g",(double)ZC(MYVERTEX(CORNER(el,i)))); strcpy(item+ic,it); ic+=strlen(it); } /* now all the user variables */ /* get local coordinate of corner */ LocalCornerCoordinates(DIM,TAG(el),i,local); for (j=0; j<DIM; j++) LocalCoord[j] = local[j]; for (v=0; v<nv; v++) { pre = ev[v]->PreprocessProc; eval = ev[v]->EvalProc; /* execute prepare function */ /* This is not really equivalent to the FEBLOCK-version sinc we call "pre" more often than there. D.Werner */ if (pre!=NULL) pre(ev_name[v],mg); /* call eval function */ value = eval(el,(const DOUBLE **)CornersCoord,LocalCoord); sprintf(it," %g",value); strcpy(item+ic,it); ic+=strlen(it); } sprintf(it,"\n"); strcpy(item+ic,it); ic+=strlen(it); pfile_tagged_puts(pf,item,counter+on); ic=0; counter++; } } pfile_sync(pf); /* end of segment */ sprintf(it,"\n"); strcpy(item+ic,it); ic+=strlen(it); pfile_master_puts(pf,item); ic=0; /* finally write the connectivity list */ counter=0; for (k=0; k<=TOPLEVEL(mg); k++) for (el=FIRSTELEMENT(GRID_ON_LEVEL(mg,k)); el!=NULL; el=SUCCE(el)) { if (!EstimateHere(el)) continue; /* process finest level elements only */ switch(DIM) { case 2 : switch(TAG(el)) { case TRIANGLE : sprintf(it,"%d %d %d %d\n", VINDEX(NVECTOR(CORNER(el,0))), VINDEX(NVECTOR(CORNER(el,1))), VINDEX(NVECTOR(CORNER(el,2))), VINDEX(NVECTOR(CORNER(el,2))) ); break; case QUADRILATERAL : sprintf(it,"%d %d %d %d\n", VINDEX(NVECTOR(CORNER(el,0))), VINDEX(NVECTOR(CORNER(el,1))), VINDEX(NVECTOR(CORNER(el,2))), VINDEX(NVECTOR(CORNER(el,3))) ); break; default : UserWriteF("tecplot: unknown 2D element type with tag(el) = %d detected. Aborting further processing of command tecplot\n", TAG(el)); return CMDERRORCODE; break; } break; case 3 : switch(TAG(el)) { case HEXAHEDRON : sprintf(it,"%d %d %d %d " "%d %d %d %d\n", VINDEX(NVECTOR(CORNER(el,0))), VINDEX(NVECTOR(CORNER(el,1))), VINDEX(NVECTOR(CORNER(el,2))), VINDEX(NVECTOR(CORNER(el,3))), VINDEX(NVECTOR(CORNER(el,4))), VINDEX(NVECTOR(CORNER(el,5))), VINDEX(NVECTOR(CORNER(el,6))), VINDEX(NVECTOR(CORNER(el,7))) ); break; case TETRAHEDRON : sprintf(it,"%d %d %d %d " "%d %d %d %d\n", VINDEX(NVECTOR(CORNER(el,0))), VINDEX(NVECTOR(CORNER(el,1))), VINDEX(NVECTOR(CORNER(el,2))), VINDEX(NVECTOR(CORNER(el,2))), VINDEX(NVECTOR(CORNER(el,3))), VINDEX(NVECTOR(CORNER(el,3))), VINDEX(NVECTOR(CORNER(el,3))), VINDEX(NVECTOR(CORNER(el,3))) ); break; case PYRAMID : sprintf(it,"%d %d %d %d " "%d %d %d %d\n", VINDEX(NVECTOR(CORNER(el,0))), VINDEX(NVECTOR(CORNER(el,1))), VINDEX(NVECTOR(CORNER(el,2))), VINDEX(NVECTOR(CORNER(el,3))), VINDEX(NVECTOR(CORNER(el,4))), VINDEX(NVECTOR(CORNER(el,4))), VINDEX(NVECTOR(CORNER(el,4))), VINDEX(NVECTOR(CORNER(el,4))) ); break; case PRISM : sprintf(it,"%d %d %d %d " "%d %d %d %d\n", VINDEX(NVECTOR(CORNER(el,0))), VINDEX(NVECTOR(CORNER(el,1))), VINDEX(NVECTOR(CORNER(el,2))), VINDEX(NVECTOR(CORNER(el,2))), VINDEX(NVECTOR(CORNER(el,3))), VINDEX(NVECTOR(CORNER(el,4))), VINDEX(NVECTOR(CORNER(el,5))), VINDEX(NVECTOR(CORNER(el,5))) ); break; default : UserWriteF("tecplot: unknown 3D element type with tag(el) = %d detected. Aborting further processing of command tecplot\n", TAG(el)); return CMDERRORCODE; break; } break; } strcpy(item+ic,it); ic+=strlen(it); pfile_tagged_puts(pf,item,counter+oe); ic=0; counter++; } pfile_sync(pf); /* end of segment */ /********************************/ /* GEOMETRY */ /* we will do this later, since */ /* domain interface will change */ /********************************/ pfile_close(pf); return(OKCODE); }