| 1 |
jgs |
82 |
/* $Id$ */ |
| 2 |
|
|
|
| 3 |
|
|
/**************************************************************/ |
| 4 |
|
|
|
| 5 |
|
|
/* writes data and mesh in an opendx file */ |
| 6 |
|
|
|
| 7 |
|
|
/**************************************************************/ |
| 8 |
|
|
|
| 9 |
|
|
/* Copyrights by ACcESS, Australia 2004 */ |
| 10 |
|
|
/* Author: Lutz Gross, gross@access.edu.au */ |
| 11 |
|
|
|
| 12 |
|
|
/**************************************************************/ |
| 13 |
|
|
|
| 14 |
|
|
#include "Finley.h" |
| 15 |
|
|
#include "Common.h" |
| 16 |
|
|
#include "Mesh.h" |
| 17 |
|
|
#include "escript/Data/DataC.h" |
| 18 |
|
|
|
| 19 |
|
|
void Finley_Mesh_saveDX(const char * filename_p, Finley_Mesh *mesh_p, escriptDataC* data_p) { |
| 20 |
|
|
/* if there is no mesh we just return */ |
| 21 |
|
|
if (mesh_p==NULL) return; |
| 22 |
|
|
/* some tables needed for reordering */ |
| 23 |
|
|
int resort[6][9]={ |
| 24 |
|
|
{0,1}, /* line */ |
| 25 |
|
|
{0,1,2}, /* triagle */ |
| 26 |
|
|
{0,1,2,3}, /* tetrahedron */ |
| 27 |
|
|
{0,3,1,2}, /* quadrilateral */ |
| 28 |
|
|
{3,0,7,4,2,1,6,5}, /* hexahedron */ |
| 29 |
|
|
}; |
| 30 |
|
|
Finley_ElementFile* elements=NULL; |
| 31 |
|
|
char elemTypeStr[32]; |
| 32 |
|
|
int i,j,k,numDXNodesPerElement,*resortIndex,isCellCentered,nodetype; |
| 33 |
|
|
double* values,rtmp; |
| 34 |
|
|
int nDim = mesh_p->Nodes->numDim; |
| 35 |
|
|
/* open the file and check handel */ |
| 36 |
|
|
FILE * fileHandle_p = fopen(filename_p, "w"); |
| 37 |
|
|
if (fileHandle_p==NULL) { |
| 38 |
|
|
Finley_ErrorCode=IO_ERROR; |
| 39 |
|
|
sprintf(Finley_ErrorMsg,"File %s could not be opened for writing.",filename_p); |
| 40 |
|
|
return; |
| 41 |
|
|
} |
| 42 |
|
|
/* positions */ |
| 43 |
|
|
fprintf(fileHandle_p, "object 1 class array type float rank 1 shape %d items %d data follows\n", |
| 44 |
|
|
nDim, mesh_p->Nodes->reducedNumNodes); |
| 45 |
|
|
for (i = 0; i < mesh_p->Nodes->numNodes; i++) { |
| 46 |
|
|
if (mesh_p->Nodes->toReduced[i]>=0) { |
| 47 |
|
|
fprintf(fileHandle_p, "%e", mesh_p->Nodes->Coordinates[INDEX2(0, i, nDim)]); |
| 48 |
|
|
for (j = 1; j < nDim; j++) fprintf(fileHandle_p, " %f",mesh_p->Nodes->Coordinates[INDEX2(j, i, nDim)]); |
| 49 |
|
|
fprintf(fileHandle_p, "\n"); |
| 50 |
|
|
} |
| 51 |
|
|
} |
| 52 |
|
|
/* connections */ |
| 53 |
|
|
/* get a pointer to the relevant mesh component */ |
| 54 |
|
|
if (isEmpty(data_p)) { |
| 55 |
|
|
elements=mesh_p->Elements; |
| 56 |
|
|
} else { |
| 57 |
|
|
switch(getFunctionSpaceType(data_p)) { |
| 58 |
|
|
case(FINLEY_DEGREES_OF_FREEDOM): |
| 59 |
|
|
nodetype=FINLEY_DEGREES_OF_FREEDOM; |
| 60 |
|
|
case(FINLEY_REDUCED_DEGREES_OF_FREEDOM): |
| 61 |
|
|
nodetype=FINLEY_REDUCED_DEGREES_OF_FREEDOM; |
| 62 |
|
|
case(FINLEY_NODES): |
| 63 |
|
|
nodetype=FINLEY_NODES; |
| 64 |
|
|
isCellCentered=FALSE; |
| 65 |
|
|
elements=mesh_p->Elements; |
| 66 |
|
|
break; |
| 67 |
|
|
case(FINLEY_ELEMENTS): |
| 68 |
|
|
isCellCentered=TRUE; |
| 69 |
|
|
elements=mesh_p->Elements; |
| 70 |
|
|
break; |
| 71 |
|
|
case(FINLEY_FACE_ELEMENTS): |
| 72 |
|
|
isCellCentered=TRUE; |
| 73 |
|
|
elements=mesh_p->FaceElements; |
| 74 |
|
|
break; |
| 75 |
|
|
case(FINLEY_POINTS): |
| 76 |
|
|
isCellCentered=TRUE; |
| 77 |
|
|
elements=mesh_p->Points; |
| 78 |
|
|
break; |
| 79 |
|
|
case(FINLEY_CONTACT_ELEMENTS_1): |
| 80 |
|
|
case(FINLEY_CONTACT_ELEMENTS_2): |
| 81 |
|
|
isCellCentered=TRUE; |
| 82 |
|
|
elements=mesh_p->ContactElements; |
| 83 |
|
|
break; |
| 84 |
|
|
default: |
| 85 |
|
|
Finley_ErrorCode=TYPE_ERROR; |
| 86 |
|
|
sprintf(Finley_ErrorMsg,"Finley does not know anything about function space type %d",getFunctionSpaceType(data_p)); |
| 87 |
|
|
break; |
| 88 |
|
|
} |
| 89 |
|
|
} |
| 90 |
|
|
/* if no element table is present jump over the connection table */ |
| 91 |
|
|
if (elements!=NULL) { |
| 92 |
|
|
ElementTypeId TypeId = elements->ReferenceElement->Type->TypeId; |
| 93 |
|
|
if (TypeId==Line2 || TypeId==Line3 || TypeId==Line4 ) { |
| 94 |
|
|
numDXNodesPerElement=2; |
| 95 |
|
|
resortIndex=resort[0]; |
| 96 |
|
|
strcpy(elemTypeStr, "lines"); |
| 97 |
|
|
} else if (TypeId==Tri3 || TypeId==Tri6 || TypeId==Tri9 || TypeId==Tri10 ) { |
| 98 |
|
|
numDXNodesPerElement = 3; |
| 99 |
|
|
resortIndex=resort[1]; |
| 100 |
|
|
strcpy(elemTypeStr, "triangles"); |
| 101 |
|
|
} else if (TypeId==Rec4 || TypeId==Rec8 || TypeId==Rec9 || TypeId==Rec12 || TypeId==Rec16 ) { |
| 102 |
|
|
numDXNodesPerElement = 4; |
| 103 |
|
|
resortIndex=resort[3]; |
| 104 |
|
|
strcpy(elemTypeStr, "quads"); |
| 105 |
|
|
} else if (TypeId==Tet4 || TypeId==Tet10 || TypeId==Tet16 ) { |
| 106 |
|
|
numDXNodesPerElement = 4; |
| 107 |
|
|
resortIndex=resort[2]; |
| 108 |
|
|
strcpy(elemTypeStr, "tetrahedra"); |
| 109 |
|
|
} else if (TypeId==Hex8 || TypeId==Hex20 || TypeId==Hex32 ) { |
| 110 |
|
|
numDXNodesPerElement = 8; |
| 111 |
|
|
resortIndex=resort[4]; |
| 112 |
|
|
strcpy(elemTypeStr, "cubes"); |
| 113 |
|
|
} else { |
| 114 |
|
|
Finley_ErrorCode=VALUE_ERROR; |
| 115 |
|
|
sprintf(Finley_ErrorMsg, "Element type %s is not supported by DX",elements->ReferenceElement->Type->Name); |
| 116 |
|
|
return; |
| 117 |
|
|
} |
| 118 |
|
|
int NN=elements->ReferenceElement->Type->numNodes; |
| 119 |
|
|
fprintf(fileHandle_p, "object 2 class array type int rank 1 shape %d items %d data follows\n",numDXNodesPerElement, elements->numElements); |
| 120 |
|
|
for (i = 0; i < elements->numElements; i++) { |
| 121 |
|
|
fprintf(fileHandle_p,"%d",mesh_p->Nodes->toReduced[mesh_p->Elements->Nodes[INDEX2(resortIndex[0], i, NN)]]); |
| 122 |
|
|
for (j = 1; j < numDXNodesPerElement; j++) { |
| 123 |
|
|
fprintf(fileHandle_p," %d",mesh_p->Nodes->toReduced[mesh_p->Elements->Nodes[INDEX2(resortIndex[j], i, NN)]]); |
| 124 |
|
|
} |
| 125 |
|
|
fprintf(fileHandle_p, "\n"); |
| 126 |
|
|
} |
| 127 |
|
|
fprintf(fileHandle_p, "attribute \"element type\" string \"%s\"\n",elemTypeStr); |
| 128 |
|
|
fprintf(fileHandle_p, "attribute \"ref\" string \"positions\"\n"); |
| 129 |
|
|
|
| 130 |
|
|
} |
| 131 |
|
|
/* data */ |
| 132 |
|
|
if (!isEmpty(data_p)) { |
| 133 |
|
|
int rank=getDataPointRank(data_p); |
| 134 |
|
|
int* shape=getDataPointShape(data_p); |
| 135 |
|
|
int nComp=getDataPointSize(data_p); |
| 136 |
|
|
fprintf(fileHandle_p, "object 3 class array type float rank %d ", rank); |
| 137 |
|
|
if (0 < rank) { |
| 138 |
|
|
fprintf(fileHandle_p, "shape "); |
| 139 |
|
|
for (i = 0; i < rank; i++) fprintf(fileHandle_p, "%d ", shape[i]); |
| 140 |
|
|
} |
| 141 |
|
|
if (isCellCentered) { |
| 142 |
|
|
int numPointsPerSample=elements->ReferenceElement->numQuadNodes; |
| 143 |
|
|
if (numPointsPerSample) { |
| 144 |
|
|
fprintf(fileHandle_p, "items %d data follows\n", elements->numElements); |
| 145 |
|
|
for (i=0;i<elements->numElements;i++) { |
| 146 |
|
|
values=getSampleData(data_p,i); |
| 147 |
|
|
for (k=0;k<nComp;k++) { |
| 148 |
|
|
rtmp=0.; |
| 149 |
|
|
for (j=0;j<numPointsPerSample;j++) rtmp+=values[INDEX2(k,j,nComp)]; |
| 150 |
|
|
fprintf(fileHandle_p, " %f", rtmp/numPointsPerSample); |
| 151 |
|
|
} |
| 152 |
|
|
fprintf(fileHandle_p, "\n"); |
| 153 |
|
|
} |
| 154 |
|
|
fprintf(fileHandle_p, "attribute \"dep\" string \"connections\"\n"); |
| 155 |
|
|
} |
| 156 |
|
|
} else { |
| 157 |
|
|
fprintf(fileHandle_p, "items %d data follows\n", mesh_p->Nodes->reducedNumNodes); |
| 158 |
|
|
for (i=0;i<mesh_p->Nodes->numNodes;i++) { |
| 159 |
|
|
if (mesh_p->Nodes->toReduced[i]>=0) { |
| 160 |
|
|
switch (nodetype) { |
| 161 |
|
|
case(FINLEY_DEGREES_OF_FREEDOM): |
| 162 |
|
|
values=getSampleData(data_p,mesh_p->Nodes->degreeOfFreedom[i]); |
| 163 |
|
|
break; |
| 164 |
|
|
case(FINLEY_REDUCED_DEGREES_OF_FREEDOM): |
| 165 |
|
|
values=getSampleData(data_p,mesh_p->Nodes->reducedDegreeOfFreedom[i]); |
| 166 |
|
|
break; |
| 167 |
|
|
case(FINLEY_NODES): |
| 168 |
|
|
values=getSampleData(data_p,i); |
| 169 |
|
|
break; |
| 170 |
|
|
} |
| 171 |
|
|
} |
| 172 |
|
|
for (k=0;k<nComp;k++) fprintf(fileHandle_p, " %f", values[k]); |
| 173 |
|
|
fprintf(fileHandle_p, "\n"); |
| 174 |
|
|
} |
| 175 |
|
|
fprintf(fileHandle_p, "attribute \"dep\" string \"positions\"\n"); |
| 176 |
|
|
} |
| 177 |
|
|
} |
| 178 |
|
|
|
| 179 |
|
|
/* and finish it up */ |
| 180 |
|
|
fprintf(fileHandle_p, "object 4 class field\n"); |
| 181 |
|
|
fprintf(fileHandle_p, "component \"positions\" value 1\n"); |
| 182 |
|
|
if (elements!=NULL) fprintf(fileHandle_p, "component \"connections\" value 2\n"); |
| 183 |
|
|
if (!isEmpty(data_p)) fprintf(fileHandle_p, "component \"data\" value 3\n"); |
| 184 |
|
|
fprintf(fileHandle_p, "end\n"); |
| 185 |
|
|
/* close the file */ |
| 186 |
|
|
fclose(fileHandle_p); |
| 187 |
|
|
return; |
| 188 |
|
|
} |
| 189 |
|
|
|
| 190 |
|
|
/* |
| 191 |
|
|
* $Log$ |
| 192 |
|
|
* Revision 1.1 2004/10/26 06:53:57 jgs |
| 193 |
|
|
* Initial revision |
| 194 |
|
|
* |
| 195 |
|
|
* Revision 1.1 2004/07/27 08:27:11 gross |
| 196 |
|
|
* Finley: saveDX added: now it is possible to write data on boundary and contact elements |
| 197 |
|
|
* |
| 198 |
|
|
*/ |
Parent Directory
| 
Revision Log