/[escript]/trunk/finley/src/Mesh_joinFaces.c
ViewVC logotype

Contents of /trunk/finley/src/Mesh_joinFaces.c

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1388 - (show annotations)
Fri Jan 11 07:45:58 2008 UTC (11 years, 11 months ago) by trankine
File MIME type: text/plain
File size: 6003 byte(s)
And get the *(&(*&(* name right
1
2 /* $Id$ */
3
4 /*******************************************************
5 *
6 * Copyright 2003-2007 by ACceSS MNRF
7 * Copyright 2007 by University of Queensland
8 *
9 * http://esscc.uq.edu.au
10 * Primary Business: Queensland, Australia
11 * Licensed under the Open Software License version 3.0
12 * http://www.opensource.org/licenses/osl-3.0.php
13 *
14 *******************************************************/
15
16 /**************************************************************/
17
18 /* Finley: Mesh */
19
20 /* detects faces in the mesh that match and replaces it by step elements */
21
22 /**************************************************************/
23
24 #include "Mesh.h"
25
26 /**************************************************************/
27
28 void Finley_Mesh_joinFaces(Finley_Mesh* self,double safety_factor,double tolerance, bool_t optimize) {
29
30 char error_msg[LenErrorMsg_MAX];
31 index_t e0,e1,*elem1=NULL,*elem0=NULL,*elem_mask=NULL,*matching_nodes_in_elem1=NULL;
32 Finley_ElementFile *newFaceElementsFile=NULL,*newContactElementsFile=NULL;
33 dim_t e,i,numPairs, NN, NN_Contact,c, new_numFaceElements;
34
35 if (self->MPIInfo->size>1) {
36 Finley_setError(TYPE_ERROR,"Finley_Mesh_joinFaces: MPI is not supported yet.");
37 return;
38 }
39 if (self->FaceElements==NULL) return;
40
41 if (self->FaceElements->ReferenceElement->Type->numNodesOnFace<=0) {
42 sprintf(error_msg,"Finley_Mesh_joinFaces:joining faces cannot be applied to face elements of type %s",self->FaceElements->ReferenceElement->Type->Name);
43 Finley_setError(TYPE_ERROR,error_msg);
44 return;
45 }
46 if (self->ContactElements==NULL) {
47 Finley_setError(TYPE_ERROR,"Finley_Mesh_joinFaces: no contact element file present.");
48 return;
49 }
50
51 NN=self->FaceElements->ReferenceElement->Type->numNodes;
52 NN_Contact=self->ContactElements->ReferenceElement->Type->numNodes;
53
54 if (2*NN!=NN_Contact) {
55 sprintf(error_msg,"Finley_Mesh_joinFaces:contact element file for %s cannot hold elements created from face elements %s",
56 self->ContactElements->ReferenceElement->Type->Name,self->FaceElements->ReferenceElement->Type->Name);
57 Finley_setError(TYPE_ERROR,error_msg);
58 return;
59 }
60
61 /* allocate work arrays */
62 elem1=TMPMEMALLOC(self->FaceElements->numElements,index_t);
63 elem0=TMPMEMALLOC(self->FaceElements->numElements,index_t);
64 elem_mask=TMPMEMALLOC(self->FaceElements->numElements,index_t);
65 matching_nodes_in_elem1=TMPMEMALLOC(self->FaceElements->numElements*NN,index_t);
66
67 if (!(Finley_checkPtr(elem1) || Finley_checkPtr(elem0) || Finley_checkPtr(elem_mask) || Finley_checkPtr(matching_nodes_in_elem1))) {
68
69 /* find the matching face elements */
70 Finley_Mesh_findMatchingFaces(self->Nodes,self->FaceElements,safety_factor,tolerance,&numPairs,elem0,elem1,matching_nodes_in_elem1);
71 if (Finley_noError()) {
72 /* get a list of the face elements to be kept */
73 #pragma omp parallel for private(e) schedule(static)
74 for(e=0;e<self->FaceElements->numElements;e++) elem_mask[e]=1;
75 for(e=0;e<numPairs;e++) {
76 elem_mask[elem0[e]]=0;
77 elem_mask[elem1[e]]=0;
78 }
79 new_numFaceElements=0;
80 /* OMP */
81 for(e=0;e<self->FaceElements->numElements;e++) {
82 if (elem_mask[e]>0) {
83 elem_mask[new_numFaceElements]=e;
84 new_numFaceElements++;
85 }
86 }
87 /* allocate new face element and Contact element files */
88 newContactElementsFile=Finley_ElementFile_alloc(self->ContactElements->ReferenceElement->Type->TypeId,self->ContactElements->order,self->ContactElements->reduced_order, self->MPIInfo);
89 newFaceElementsFile=Finley_ElementFile_alloc(self->FaceElements->ReferenceElement->Type->TypeId,self->FaceElements->order,self->FaceElements->reduced_order, self->MPIInfo);
90 if (Finley_noError()) {
91 Finley_ElementFile_allocTable(newContactElementsFile,numPairs+self->ContactElements->numElements);
92 Finley_ElementFile_allocTable(newFaceElementsFile,new_numFaceElements);
93 }
94 /* copy the old elements over */
95 if (Finley_noError()) {
96 /* get the face elements which are still in use:*/
97 Finley_ElementFile_gather(elem_mask,self->FaceElements,newFaceElementsFile);
98 /* get the Contact elements which are still in use:*/
99 Finley_ElementFile_copyTable(0,newContactElementsFile,0,0,self->ContactElements);
100 c=self->ContactElements->numElements;
101 /* OMP */
102 for (e=0;e<numPairs;e++) {
103 e0=elem0[e];
104 e1=elem1[e];
105 newContactElementsFile->Id[c]=MIN(self->FaceElements->Id[e0],self->FaceElements->Id[e1]);
106 newContactElementsFile->Tag[c]=MIN(self->FaceElements->Tag[e0],self->FaceElements->Tag[e1]);
107 newContactElementsFile->Color[c]=e;
108 for (i=0;i<NN;i++) newContactElementsFile->Nodes[INDEX2(i,c,NN_Contact)]=self->FaceElements->Nodes[INDEX2(i,e0,NN)];
109 for (i=0;i<NN;i++) newContactElementsFile->Nodes[INDEX2(i+NN,c,NN_Contact)]=matching_nodes_in_elem1[INDEX2(i,e,NN)];
110 c++;
111 }
112 newContactElementsFile->minColor=0;
113 newContactElementsFile->maxColor=numPairs-1;
114 }
115 /* set new face and Contact elements */
116 if (Finley_noError()) {
117
118 Finley_ElementFile_free(self->FaceElements);
119 self->FaceElements=newFaceElementsFile;
120 Finley_ElementFile_free(self->ContactElements);
121 self->ContactElements=newContactElementsFile;
122 Finley_Mesh_prepare(self, optimize);
123
124 } else {
125 Finley_ElementFile_free(newFaceElementsFile);
126 Finley_ElementFile_free(newContactElementsFile);
127 }
128 }
129 }
130 TMPMEMFREE(elem1);
131 TMPMEMFREE(elem0);
132 TMPMEMFREE(matching_nodes_in_elem1);
133 TMPMEMFREE(elem_mask);
134 }

Properties

Name Value
svn:eol-style native
svn:keywords Author Date Id Revision

  ViewVC Help
Powered by ViewVC 1.1.26