finley/src/Mesh_rec4.c


/* $Id$ */

/*******************************************************
 *
 *           Copyright 2003-2007 by ACceSS MNRF
 *       Copyright 2007 by University of Queensland
 *
 *                http://esscc.uq.edu.au
 *        Primary Business: Queensland, Australia
 *  Licensed under the Open Software License version 3.0
 *     http://www.opensource.org/licenses/osl-3.0.php
 *
 *******************************************************/

/**************************************************************/

/*   Finley: generates rectangular meshes  */

/*   Generates a numElements[0] x numElements[1] mesh with first order elements (Rec4) in the rectangle */
/*   [0,Length[0]] x [0,Length[1]]. order is the desired accuracy of the integration scheme. */


/**************************************************************/

#include "RectangularMesh.h"


Finley_Mesh* Finley_RectangularMesh_Rec4(dim_t* numElements,
                                          double* Length,
                                          bool_t* periodic,
                                          index_t order, 
                                          index_t reduced_order, 
                                          bool_t useElementsOnFace,
                                          bool_t useFullElementOrder,
                                          bool_t optimize)
{
  #define N_PER_E 1
  #define DIM 2
  dim_t N0,N1,NE0,NE1,i0,i1,k,Nstride0,Nstride1, local_NE0, local_NE1, local_N0, local_N1, global_i0, global_i1;
  index_t offset0, offset1, e_offset0, e_offset1;
  dim_t totalNECount,faceNECount,NDOF0,NDOF1,NFaceElements,NN;
  index_t node0, myRank;
  Finley_Mesh* out;
  Paso_MPIInfo *mpi_info = NULL;
  char name[50];
  double time0=Finley_timer();

  /* get MPI information */
  mpi_info = Paso_MPIInfo_alloc( MPI_COMM_WORLD );
  if (! Finley_noError()) {
        return NULL;
  }
  myRank=mpi_info->rank;

  /* set up the global dimensions of the mesh */

  NE0=MAX(1,numElements[0]);
  NE1=MAX(1,numElements[1]);
  N0=N_PER_E*NE0+1;
  N1=N_PER_E*NE1+1;

  /*  allocate mesh: */  
  sprintf(name,"Rectangular %d x %d mesh",N0,N1);
  out=Finley_Mesh_alloc(name,DIM,order, reduced_order, mpi_info);
  if (! Finley_noError()) { 
      Paso_MPIInfo_free( mpi_info );
      return NULL;
  }

  Finley_Mesh_setElements(out,Finley_ElementFile_alloc(Rec4,out->order,out->reduced_order,mpi_info));
  if (useElementsOnFace) {
         Finley_Mesh_setFaceElements(out,Finley_ElementFile_alloc(Rec4Face,
                                                                  out->order,
                                                                  out->reduced_order,
                                                                  mpi_info));
         Finley_Mesh_setContactElements(out,Finley_ElementFile_alloc(Rec4Face_Contact,
                                                                    out->order,
                                                                    out->reduced_order,
                                                                    mpi_info));
  } else {
         Finley_Mesh_setFaceElements(out,Finley_ElementFile_alloc(Line2,
                                                                  out->order,
                                                                  out->reduced_order,
                                                                  mpi_info));
         Finley_Mesh_setContactElements(out,Finley_ElementFile_alloc(Line2_Contact,
                                                                     out->order,
                                                                     out->reduced_order,
                                                                     mpi_info));
  }
  Finley_Mesh_setPoints(out,Finley_ElementFile_alloc(Point1,
                                                 out->order,
                                                 out->reduced_order,
                                                 mpi_info));
  if (! Finley_noError()) {
      Paso_MPIInfo_free( mpi_info );
      Finley_Mesh_free(out);
      return NULL;
  }

  /* work out the largest dimension */
  if (N1==MAX(N0,N1)) {
     Nstride0=1;
     Nstride1=N0;
     local_NE0=NE0;
     e_offset0=0;
     Paso_MPIInfo_Split(mpi_info,NE1,&local_NE1,&e_offset1);
  } else {
     Nstride0=N1;
     Nstride1=1;
     Paso_MPIInfo_Split(mpi_info,NE0,&local_NE0,&e_offset0);
     local_NE1=NE1;
     e_offset1=0;
  }
  offset0=e_offset0*N_PER_E;
  offset1=e_offset1*N_PER_E;
  local_N0=local_NE0>0 ? local_NE0*N_PER_E+1 : 0;
  local_N1=local_NE1>0 ? local_NE1*N_PER_E+1 : 0;

  /* get the number of surface elements */

  NFaceElements=0;
  if (!periodic[0] && (local_NE0>0)) {
     NDOF0=N0;
     if (e_offset0 == 0) NFaceElements+=local_NE1;
     if (local_NE0+e_offset0 == NE0) NFaceElements+=local_NE1;
  } else {
      NDOF0=N0-1;
  }
  if (!periodic[1]  && (local_NE1>0)) {
     NDOF1=N1;
     if (e_offset1 == 0) NFaceElements+=local_NE0;
     if (local_NE1+e_offset1 == NE1) NFaceElements+=local_NE0;
  } else {
      NDOF1=N1-1;
  }

  /*  allocate tables: */

  Finley_NodeFile_allocTable(out->Nodes,local_N0*local_N1);
  Finley_ElementFile_allocTable(out->Elements,local_NE0*local_NE1);
  Finley_ElementFile_allocTable(out->FaceElements,NFaceElements);

  if (Finley_noError()) {
     /* create nodes */
     #pragma omp parallel for private(i0,i1,k,global_i0,global_i1)
     for (i1=0;i1<local_N1;i1++) {
       for (i0=0;i0<local_N0;i0++) {
           k=i0+local_N0*i1;
           global_i0=i0+offset0;
           global_i1=i1+offset1;
           out->Nodes->Coordinates[INDEX2(0,k,DIM)]=DBLE(global_i0)/DBLE(N0-1)*Length[0];
           out->Nodes->Coordinates[INDEX2(1,k,DIM)]=DBLE(global_i1)/DBLE(N1-1)*Length[1];
           out->Nodes->Id[k]=Nstride0*global_i0+Nstride1*global_i1;
           out->Nodes->Tag[k]=0;
           out->Nodes->globalDegreesOfFreedom[k]=Nstride0*(global_i0%NDOF0) 
                                               +Nstride1*(global_i1%NDOF1);
       }
     }
     /*   set the elements: */
     NN=out->Elements->numNodes;
     #pragma omp parallel for private(i0,i1,k,node0) 
     for (i1=0;i1<local_NE1;i1++) {
         for (i0=0;i0<local_NE0;i0++) {
           
           k=i0+local_NE0*i1;        
           node0=Nstride0*N_PER_E*(i0+e_offset0)+Nstride1*N_PER_E*(i1+e_offset1);

           out->Elements->Id[k]=(i0+e_offset0)+NE0*(i1+e_offset1);
           out->Elements->Tag[k]=0;
           out->Elements->Owner[k]=myRank;

           out->Elements->Nodes[INDEX2(0,k,NN)]=node0;
           out->Elements->Nodes[INDEX2(1,k,NN)]=node0+Nstride0;
           out->Elements->Nodes[INDEX2(2,k,NN)]=node0+Nstride1+Nstride0;
           out->Elements->Nodes[INDEX2(3,k,NN)]=node0+Nstride1;
         }
     }
     /* face elements */
     NN=out->FaceElements->numNodes;
     totalNECount=NE0*NE1;
     faceNECount=0;
     if (!periodic[0]  && (local_NE0>0)) {
        /* **  elements on boundary 001 (x1=0): */
     
        if (e_offset0 == 0) {
           #pragma omp parallel for private(i1,k,node0) 
           for (i1=0;i1<local_NE1;i1++) {
      
               k=i1+faceNECount;
               node0=Nstride1*N_PER_E*(i1+e_offset1);

               out->FaceElements->Id[k]=i1+e_offset1+totalNECount;
               out->FaceElements->Tag[k]=1;
               out->FaceElements->Owner[k]=myRank;
               if (useElementsOnFace) {
                  out->FaceElements->Nodes[INDEX2(0,k,NN)]=node0+Nstride1;
                  out->FaceElements->Nodes[INDEX2(1,k,NN)]=node0;
                  out->FaceElements->Nodes[INDEX2(2,k,NN)]=node0+Nstride0;
                  out->FaceElements->Nodes[INDEX2(3,k,NN)]=node0+Nstride1+Nstride0;
               } else {
                  out->FaceElements->Nodes[INDEX2(0,k,NN)]=node0+Nstride1;
                  out->FaceElements->Nodes[INDEX2(1,k,NN)]=node0;
               }
           }
           faceNECount+=local_NE1;
        }
        totalNECount+=NE1;
        /* **  elements on boundary 002 (x1=1): */
        if (local_NE0+e_offset0 == NE0) {
           #pragma omp parallel for private(i1,k,node0) 
           for (i1=0;i1<local_NE1;i1++) {
               k=i1+faceNECount;
               node0=Nstride0*N_PER_E*(NE0-1)+Nstride1*N_PER_E*(i1+e_offset1);

               out->FaceElements->Id[k]=(i1+e_offset1)+totalNECount;
               out->FaceElements->Tag[k]=2;
               out->FaceElements->Owner[k]=myRank;

               if (useElementsOnFace) {
                  out->FaceElements->Nodes[INDEX2(0,k,NN)]=node0+Nstride0;
                  out->FaceElements->Nodes[INDEX2(1,k,NN)]=node0+Nstride1+Nstride0;
                  out->FaceElements->Nodes[INDEX2(2,k,NN)]=node0+Nstride1;
                  out->FaceElements->Nodes[INDEX2(3,k,NN)]=node0;
               } else {
                  out->FaceElements->Nodes[INDEX2(0,k,NN)]=node0+Nstride0;
                  out->FaceElements->Nodes[INDEX2(1,k,NN)]=node0+Nstride1+Nstride0;
               }
           }
           faceNECount+=local_NE1;
         }
         totalNECount+=NE1;
     }
     if (!periodic[1]  && (local_NE1>0)) {
        /* **  elements on boundary 010 (x2=0): */
        if (e_offset1 == 0) {
           #pragma omp parallel for private(i0,k,node0) 
           for (i0=0;i0<local_NE0;i0++) {
               k=i0+faceNECount;
               node0=Nstride0*N_PER_E*(i0+e_offset0);
         
               out->FaceElements->Id[k]=e_offset0+i0+totalNECount;
               out->FaceElements->Tag[k]=10;
               out->FaceElements->Owner[k]=myRank;
       
               if (useElementsOnFace) {
                   out->FaceElements->Nodes[INDEX2(0,k,NN)]=node0;
                   out->FaceElements->Nodes[INDEX2(1,k,NN)]=node0+Nstride0;
                   out->FaceElements->Nodes[INDEX2(2,k,NN)]=node0+Nstride1+Nstride0;
                   out->FaceElements->Nodes[INDEX2(3,k,NN)]=node0+Nstride1;
               } else {
                   out->FaceElements->Nodes[INDEX2(0,k,NN)]=node0;
                   out->FaceElements->Nodes[INDEX2(1,k,NN)]=node0+Nstride0;
               }
           }
           faceNECount+=local_NE0;
        }
        totalNECount+=NE0;
        /* **  elements on boundary 020 (x2=1): */
        if (local_NE1+e_offset1 == NE1) {
           #pragma omp parallel for private(i0,k,node0) 
           for (i0=0;i0<local_NE0;i0++) {
               k=i0+faceNECount;
               node0=Nstride0*N_PER_E*(i0+e_offset0)+Nstride1*N_PER_E*(NE1-1);
   
               out->FaceElements->Id[k]=i0+e_offset0+totalNECount;
               out->FaceElements->Tag[k]=20;
               out->FaceElements->Owner[k]=myRank;
               if (useElementsOnFace) {
                    out->FaceElements->Nodes[INDEX2(0,k,NN)]=node0+Nstride1+Nstride0;
                    out->FaceElements->Nodes[INDEX2(1,k,NN)]=node0+Nstride1;
                    out->FaceElements->Nodes[INDEX2(2,k,NN)]=node0;
                    out->FaceElements->Nodes[INDEX2(3,k,NN)]=node0+Nstride0;
               } else {
                    out->FaceElements->Nodes[INDEX2(0,k,NN)]=node0+Nstride1+Nstride0;
                    out->FaceElements->Nodes[INDEX2(1,k,NN)]=node0+Nstride1;
               }
           }
           faceNECount+=local_NE0;
        }
        totalNECount+=NE0;
     }
     /* add tag names */
     Finley_Mesh_addTagMap(out,"top", 20);
     Finley_Mesh_addTagMap(out,"bottom", 10);
     Finley_Mesh_addTagMap(out,"left", 1);
     Finley_Mesh_addTagMap(out,"right", 2);
   
     /* prepare mesh for further calculatuions:*/
     if (Finley_noError()) {
         Finley_Mesh_resolveNodeIds(out);
     }
     if (Finley_noError()) {
         Finley_Mesh_prepare(out, optimize);
     }
  }

  if (!Finley_noError()) {
      Finley_Mesh_free(out);
  }
  /* free up memory */
  Paso_MPIInfo_free( mpi_info );  
  #ifdef Finley_TRACE
  printf("timing: mesh generation: %.4e sec\n",Finley_timer()-time0);
  #endif

  return out;
}
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: generates rectangular meshes */
19
20	/* Generates a numElements[0] x numElements[1] mesh with first order elements (Rec4) in the rectangle */
21	/* [0,Length[0]] x [0,Length[1]]. order is the desired accuracy of the integration scheme. */
22
23
24	/**************************************************************/
25
26	#include "RectangularMesh.h"
27
28
29	Finley_Mesh* Finley_RectangularMesh_Rec4(dim_t* numElements,
30	double* Length,
31	bool_t* periodic,
32	index_t order,
33	index_t reduced_order,
34	bool_t useElementsOnFace,
35	bool_t useFullElementOrder,
36	bool_t optimize)
37	{
38	#define N_PER_E 1
39	#define DIM 2
40	dim_t N0,N1,NE0,NE1,i0,i1,k,Nstride0,Nstride1, local_NE0, local_NE1, local_N0, local_N1, global_i0, global_i1;
41	index_t offset0, offset1, e_offset0, e_offset1;
42	dim_t totalNECount,faceNECount,NDOF0,NDOF1,NFaceElements,NN;
43	index_t node0, myRank;
44	Finley_Mesh* out;
45	Paso_MPIInfo *mpi_info = NULL;
46	char name[50];
47	double time0=Finley_timer();
48
49	/* get MPI information */
50	mpi_info = Paso_MPIInfo_alloc( MPI_COMM_WORLD );
51	if (! Finley_noError()) {
52	return NULL;
53	}
54	myRank=mpi_info->rank;
55
56	/* set up the global dimensions of the mesh */
57
58	NE0=MAX(1,numElements[0]);
59	NE1=MAX(1,numElements[1]);
60	N0=N_PER_E*NE0+1;
61	N1=N_PER_E*NE1+1;
62
63	/* allocate mesh: */
64	sprintf(name,"Rectangular %d x %d mesh",N0,N1);
65	out=Finley_Mesh_alloc(name,DIM,order, reduced_order, mpi_info);
66	if (! Finley_noError()) {
67	Paso_MPIInfo_free( mpi_info );
68	return NULL;
69	}
70
71	Finley_Mesh_setElements(out,Finley_ElementFile_alloc(Rec4,out->order,out->reduced_order,mpi_info));
72	if (useElementsOnFace) {
73	Finley_Mesh_setFaceElements(out,Finley_ElementFile_alloc(Rec4Face,
74	out->order,
75	out->reduced_order,
76	mpi_info));
77	Finley_Mesh_setContactElements(out,Finley_ElementFile_alloc(Rec4Face_Contact,
78	out->order,
79	out->reduced_order,
80	mpi_info));
81	} else {
82	Finley_Mesh_setFaceElements(out,Finley_ElementFile_alloc(Line2,
83	out->order,
84	out->reduced_order,
85	mpi_info));
86	Finley_Mesh_setContactElements(out,Finley_ElementFile_alloc(Line2_Contact,
87	out->order,
88	out->reduced_order,
89	mpi_info));
90	}
91	Finley_Mesh_setPoints(out,Finley_ElementFile_alloc(Point1,
92	out->order,
93	out->reduced_order,
94	mpi_info));
95	if (! Finley_noError()) {
96	Paso_MPIInfo_free( mpi_info );
97	Finley_Mesh_free(out);
98	return NULL;
99	}
100
101	/* work out the largest dimension */
102	if (N1==MAX(N0,N1)) {
103	Nstride0=1;
104	Nstride1=N0;
105	local_NE0=NE0;
106	e_offset0=0;
107	Paso_MPIInfo_Split(mpi_info,NE1,&local_NE1,&e_offset1);
108	} else {
109	Nstride0=N1;
110	Nstride1=1;
111	Paso_MPIInfo_Split(mpi_info,NE0,&local_NE0,&e_offset0);
112	local_NE1=NE1;
113	e_offset1=0;
114	}
115	offset0=e_offset0*N_PER_E;
116	offset1=e_offset1*N_PER_E;
117	local_N0=local_NE0>0 ? local_NE0*N_PER_E+1 : 0;
118	local_N1=local_NE1>0 ? local_NE1*N_PER_E+1 : 0;
119
120	/* get the number of surface elements */
121
122	NFaceElements=0;
123	if (!periodic[0] && (local_NE0>0)) {
124	NDOF0=N0;
125	if (e_offset0 == 0) NFaceElements+=local_NE1;
126	if (local_NE0+e_offset0 == NE0) NFaceElements+=local_NE1;
127	} else {
128	NDOF0=N0-1;
129	}
130	if (!periodic[1] && (local_NE1>0)) {
131	NDOF1=N1;
132	if (e_offset1 == 0) NFaceElements+=local_NE0;
133	if (local_NE1+e_offset1 == NE1) NFaceElements+=local_NE0;
134	} else {
135	NDOF1=N1-1;
136	}
137
138	/* allocate tables: */
139
140	Finley_NodeFile_allocTable(out->Nodes,local_N0*local_N1);
141	Finley_ElementFile_allocTable(out->Elements,local_NE0*local_NE1);
142	Finley_ElementFile_allocTable(out->FaceElements,NFaceElements);
143
144	if (Finley_noError()) {
145	/* create nodes */
146	#pragma omp parallel for private(i0,i1,k,global_i0,global_i1)
147	for (i1=0;i1<local_N1;i1++) {
148	for (i0=0;i0<local_N0;i0++) {
149	k=i0+local_N0*i1;
150	global_i0=i0+offset0;
151	global_i1=i1+offset1;
152	out->Nodes->Coordinates[INDEX2(0,k,DIM)]=DBLE(global_i0)/DBLE(N0-1)*Length[0];
153	out->Nodes->Coordinates[INDEX2(1,k,DIM)]=DBLE(global_i1)/DBLE(N1-1)*Length[1];
154	out->Nodes->Id[k]=Nstride0global_i0+Nstride1global_i1;
155	out->Nodes->Tag[k]=0;
156	out->Nodes->globalDegreesOfFreedom[k]=Nstride0*(global_i0%NDOF0)
157	+Nstride1*(global_i1%NDOF1);
158	}
159	}
160	/* set the elements: */
161	NN=out->Elements->numNodes;
162	#pragma omp parallel for private(i0,i1,k,node0)
163	for (i1=0;i1<local_NE1;i1++) {
164	for (i0=0;i0<local_NE0;i0++) {
165
166	k=i0+local_NE0*i1;
167	node0=Nstride0N_PER_E(i0+e_offset0)+Nstride1N_PER_E(i1+e_offset1);
168
169	out->Elements->Id[k]=(i0+e_offset0)+NE0*(i1+e_offset1);
170	out->Elements->Tag[k]=0;
171	out->Elements->Owner[k]=myRank;
172
173	out->Elements->Nodes[INDEX2(0,k,NN)]=node0;
174	out->Elements->Nodes[INDEX2(1,k,NN)]=node0+Nstride0;
175	out->Elements->Nodes[INDEX2(2,k,NN)]=node0+Nstride1+Nstride0;
176	out->Elements->Nodes[INDEX2(3,k,NN)]=node0+Nstride1;
177	}
178	}
179	/* face elements */
180	NN=out->FaceElements->numNodes;
181	totalNECount=NE0*NE1;
182	faceNECount=0;
183	if (!periodic[0] && (local_NE0>0)) {
184	/* ** elements on boundary 001 (x1=0): */
185
186	if (e_offset0 == 0) {
187	#pragma omp parallel for private(i1,k,node0)
188	for (i1=0;i1<local_NE1;i1++) {
189
190	k=i1+faceNECount;
191	node0=Nstride1N_PER_E(i1+e_offset1);
192
193	out->FaceElements->Id[k]=i1+e_offset1+totalNECount;
194	out->FaceElements->Tag[k]=1;
195	out->FaceElements->Owner[k]=myRank;
196	if (useElementsOnFace) {
197	out->FaceElements->Nodes[INDEX2(0,k,NN)]=node0+Nstride1;
198	out->FaceElements->Nodes[INDEX2(1,k,NN)]=node0;
199	out->FaceElements->Nodes[INDEX2(2,k,NN)]=node0+Nstride0;
200	out->FaceElements->Nodes[INDEX2(3,k,NN)]=node0+Nstride1+Nstride0;
201	} else {
202	out->FaceElements->Nodes[INDEX2(0,k,NN)]=node0+Nstride1;
203	out->FaceElements->Nodes[INDEX2(1,k,NN)]=node0;
204	}
205	}
206	faceNECount+=local_NE1;
207	}
208	totalNECount+=NE1;
209	/* ** elements on boundary 002 (x1=1): */
210	if (local_NE0+e_offset0 == NE0) {
211	#pragma omp parallel for private(i1,k,node0)
212	for (i1=0;i1<local_NE1;i1++) {
213	k=i1+faceNECount;
214	node0=Nstride0N_PER_E(NE0-1)+Nstride1N_PER_E(i1+e_offset1);
215
216	out->FaceElements->Id[k]=(i1+e_offset1)+totalNECount;
217	out->FaceElements->Tag[k]=2;
218	out->FaceElements->Owner[k]=myRank;
219
220	if (useElementsOnFace) {
221	out->FaceElements->Nodes[INDEX2(0,k,NN)]=node0+Nstride0;
222	out->FaceElements->Nodes[INDEX2(1,k,NN)]=node0+Nstride1+Nstride0;
223	out->FaceElements->Nodes[INDEX2(2,k,NN)]=node0+Nstride1;
224	out->FaceElements->Nodes[INDEX2(3,k,NN)]=node0;
225	} else {
226	out->FaceElements->Nodes[INDEX2(0,k,NN)]=node0+Nstride0;
227	out->FaceElements->Nodes[INDEX2(1,k,NN)]=node0+Nstride1+Nstride0;
228	}
229	}
230	faceNECount+=local_NE1;
231	}
232	totalNECount+=NE1;
233	}
234	if (!periodic[1] && (local_NE1>0)) {
235	/* ** elements on boundary 010 (x2=0): */
236	if (e_offset1 == 0) {
237	#pragma omp parallel for private(i0,k,node0)
238	for (i0=0;i0<local_NE0;i0++) {
239	k=i0+faceNECount;
240	node0=Nstride0N_PER_E(i0+e_offset0);
241
242	out->FaceElements->Id[k]=e_offset0+i0+totalNECount;
243	out->FaceElements->Tag[k]=10;
244	out->FaceElements->Owner[k]=myRank;
245
246	if (useElementsOnFace) {
247	out->FaceElements->Nodes[INDEX2(0,k,NN)]=node0;
248	out->FaceElements->Nodes[INDEX2(1,k,NN)]=node0+Nstride0;
249	out->FaceElements->Nodes[INDEX2(2,k,NN)]=node0+Nstride1+Nstride0;
250	out->FaceElements->Nodes[INDEX2(3,k,NN)]=node0+Nstride1;
251	} else {
252	out->FaceElements->Nodes[INDEX2(0,k,NN)]=node0;
253	out->FaceElements->Nodes[INDEX2(1,k,NN)]=node0+Nstride0;
254	}
255	}
256	faceNECount+=local_NE0;
257	}
258	totalNECount+=NE0;
259	/* ** elements on boundary 020 (x2=1): */
260	if (local_NE1+e_offset1 == NE1) {
261	#pragma omp parallel for private(i0,k,node0)
262	for (i0=0;i0<local_NE0;i0++) {
263	k=i0+faceNECount;
264	node0=Nstride0N_PER_E(i0+e_offset0)+Nstride1N_PER_E(NE1-1);
265
266	out->FaceElements->Id[k]=i0+e_offset0+totalNECount;
267	out->FaceElements->Tag[k]=20;
268	out->FaceElements->Owner[k]=myRank;
269	if (useElementsOnFace) {
270	out->FaceElements->Nodes[INDEX2(0,k,NN)]=node0+Nstride1+Nstride0;
271	out->FaceElements->Nodes[INDEX2(1,k,NN)]=node0+Nstride1;
272	out->FaceElements->Nodes[INDEX2(2,k,NN)]=node0;
273	out->FaceElements->Nodes[INDEX2(3,k,NN)]=node0+Nstride0;
274	} else {
275	out->FaceElements->Nodes[INDEX2(0,k,NN)]=node0+Nstride1+Nstride0;
276	out->FaceElements->Nodes[INDEX2(1,k,NN)]=node0+Nstride1;
277	}
278	}
279	faceNECount+=local_NE0;
280	}
281	totalNECount+=NE0;
282	}
283	/* add tag names */
284	Finley_Mesh_addTagMap(out,"top", 20);
285	Finley_Mesh_addTagMap(out,"bottom", 10);
286	Finley_Mesh_addTagMap(out,"left", 1);
287	Finley_Mesh_addTagMap(out,"right", 2);
288
289	/* prepare mesh for further calculatuions:*/
290	if (Finley_noError()) {
291	Finley_Mesh_resolveNodeIds(out);
292	}
293	if (Finley_noError()) {
294	Finley_Mesh_prepare(out, optimize);
295	}
296	}
297
298	if (!Finley_noError()) {
299	Finley_Mesh_free(out);
300	}
301	/* free up memory */
302	Paso_MPIInfo_free( mpi_info );
303	#ifdef Finley_TRACE
304	printf("timing: mesh generation: %.4e sec\n",Finley_timer()-time0);
305	#endif
306
307	return out;
308	}
Name	Value
svn:eol-style	native
svn:keywords	Author Date Id Revision