1 |
/* |
|
2 |
****************************************************************************** |
/* $Id$ */ |
3 |
* * |
|
4 |
* COPYRIGHT ACcESS 2003,2004,2005 - All Rights Reserved * |
/******************************************************* |
5 |
* * |
* |
6 |
* This software is the property of ACcESS. No part of this code * |
* Copyright 2003-2007 by ACceSS MNRF |
7 |
* may be copied in any form or by any means without the expressed written * |
* Copyright 2007 by University of Queensland |
8 |
* consent of ACcESS. Copying, use or modification of this software * |
* |
9 |
* by any unauthorised person is illegal unless that person has a software * |
* http://esscc.uq.edu.au |
10 |
* license agreement with ACcESS. * |
* 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 |
|
|
24 |
|
|
25 |
/**************************************************************/ |
/**************************************************************/ |
26 |
|
|
|
/* Author: gross@access.edu.au */ |
|
|
/* Version: $Id$ |
|
|
|
|
|
/**************************************************************/ |
|
|
|
|
27 |
#include "RectangularMesh.h" |
#include "RectangularMesh.h" |
28 |
|
|
29 |
/**************************************************************/ |
Finley_Mesh* Finley_RectangularMesh_Hex20(dim_t* numElements, |
30 |
|
double* Length, |
31 |
Finley_Mesh* Finley_RectangularMesh_Hex20(dim_t* numElements,double* Length,bool_t* periodic,index_t order,bool_t useElementsOnFace) { |
bool_t* periodic, |
32 |
dim_t N0,N1,N2,NE0,NE1,NE2,i0,i1,i2,k,totalNECount,faceNECount,NDOF0,NDOF1,NDOF2,NFaceElements,NUMNODES; |
index_t order, |
33 |
index_t node0; |
index_t reduced_order, |
34 |
|
bool_t useElementsOnFace, |
35 |
|
bool_t useFullElementOrder, |
36 |
|
bool_t optimize) |
37 |
|
{ |
38 |
|
#define N_PER_E 2 |
39 |
|
#define DIM 3 |
40 |
|
dim_t N0,N1,N2,NE0,NE1,NE2,i0,i1,i2,k,Nstride0,Nstride1,Nstride2, local_NE0, local_NE1, local_NE2; |
41 |
|
dim_t totalNECount,faceNECount,NDOF0,NDOF1,NDOF2,NFaceElements, local_N0, local_N1, local_N2, NN; |
42 |
|
index_t node0, myRank, e_offset0, e_offset1, e_offset2, offset0, offset1, offset2, global_i0, global_i1, global_i2; |
43 |
Finley_Mesh* out; |
Finley_Mesh* out; |
44 |
|
Paso_MPIInfo *mpi_info = NULL; |
45 |
char name[50]; |
char name[50]; |
46 |
double time0=Finley_timer(); |
double time0=Finley_timer(); |
47 |
|
|
48 |
|
/* get MPI information */ |
49 |
|
mpi_info = Paso_MPIInfo_alloc( MPI_COMM_WORLD ); |
50 |
|
if (! Finley_noError()) { |
51 |
|
return NULL; |
52 |
|
} |
53 |
|
myRank=mpi_info->rank; |
54 |
|
|
55 |
|
/* set up the global dimensions of the mesh */ |
56 |
|
|
57 |
NE0=MAX(1,numElements[0]); |
NE0=MAX(1,numElements[0]); |
58 |
NE1=MAX(1,numElements[1]); |
NE1=MAX(1,numElements[1]); |
59 |
NE2=MAX(1,numElements[2]); |
NE2=MAX(1,numElements[2]); |
60 |
N0=2*NE0+1; |
N0=N_PER_E*NE0+1; |
61 |
N1=2*NE1+1; |
N1=N_PER_E*NE1+1; |
62 |
N2=2*NE2+1; |
N2=N_PER_E*NE2+1; |
63 |
|
|
64 |
NFaceElements=0; |
/* allocate mesh: */ |
|
if (!periodic[0]) { |
|
|
NDOF0=N0; |
|
|
NFaceElements+=2*NE1*NE2; |
|
|
} else { |
|
|
NDOF0=N0-1; |
|
|
} |
|
|
if (!periodic[1]) { |
|
|
NDOF1=N1; |
|
|
NFaceElements+=2*NE0*NE2; |
|
|
} else { |
|
|
NDOF1=N1-1; |
|
|
} |
|
|
if (!periodic[2]) { |
|
|
NDOF2=N2; |
|
|
NFaceElements+=2*NE0*NE1; |
|
|
} else { |
|
|
NDOF2=N2-1; |
|
|
} |
|
|
|
|
|
/* allocate mesh: */ |
|
|
|
|
65 |
sprintf(name,"Rectangular %d x %d x %d mesh",N0,N1,N2); |
sprintf(name,"Rectangular %d x %d x %d mesh",N0,N1,N2); |
66 |
out=Finley_Mesh_alloc(name,3,order); |
out=Finley_Mesh_alloc(name,DIM,order, reduced_order, mpi_info); |
67 |
if (! Finley_noError()) return NULL; |
if (! Finley_noError()) { |
68 |
|
Paso_MPIInfo_free( mpi_info ); |
|
out->Elements=Finley_ElementFile_alloc(Hex20,out->order); |
|
|
if (useElementsOnFace) { |
|
|
out->FaceElements=Finley_ElementFile_alloc(Hex20Face,out->order); |
|
|
out->ContactElements=Finley_ElementFile_alloc(Hex20Face_Contact,out->order); |
|
|
} else { |
|
|
out->FaceElements=Finley_ElementFile_alloc(Rec8,out->order); |
|
|
out->ContactElements=Finley_ElementFile_alloc(Rec8_Contact,out->order); |
|
|
} |
|
|
out->Points=Finley_ElementFile_alloc(Point1,out->order); |
|
|
if (! Finley_noError()) { |
|
|
Finley_Mesh_dealloc(out); |
|
69 |
return NULL; |
return NULL; |
70 |
} |
} |
71 |
|
|
72 |
|
if (useFullElementOrder) { |
73 |
/* allocate tables: */ |
/* Finley_setError(SYSTEM_ERROR,"full element order for Hex elements is not supported yet."); */ |
74 |
|
Finley_Mesh_setElements(out,Finley_ElementFile_alloc(Hex27, |
75 |
Finley_NodeFile_allocTable(out->Nodes,N0*N1*N2); |
out->order, |
76 |
Finley_ElementFile_allocTable(out->Elements,NE0*NE1*NE2); |
out->reduced_order, |
77 |
Finley_ElementFile_allocTable(out->FaceElements,NFaceElements); |
mpi_info)); |
78 |
|
if (useElementsOnFace) { |
79 |
|
Finley_setError(SYSTEM_ERROR,"rich elements for Hex27 elements is not supported yet."); |
80 |
|
} else { |
81 |
|
Finley_Mesh_setFaceElements(out,Finley_ElementFile_alloc(Rec9, |
82 |
|
out->order, |
83 |
|
out->reduced_order, |
84 |
|
mpi_info)); |
85 |
|
Finley_Mesh_setContactElements(out,Finley_ElementFile_alloc(Rec9_Contact, |
86 |
|
out->order, |
87 |
|
out->reduced_order, |
88 |
|
mpi_info)); |
89 |
|
} |
90 |
|
|
91 |
|
} else { |
92 |
|
Finley_Mesh_setElements(out,Finley_ElementFile_alloc(Hex20,out->order,out->reduced_order,mpi_info)); |
93 |
|
if (useElementsOnFace) { |
94 |
|
Finley_Mesh_setFaceElements(out,Finley_ElementFile_alloc(Hex20Face, |
95 |
|
out->order, |
96 |
|
out->reduced_order, |
97 |
|
mpi_info)); |
98 |
|
Finley_Mesh_setContactElements(out,Finley_ElementFile_alloc(Hex20Face_Contact, |
99 |
|
out->order, |
100 |
|
out->reduced_order, |
101 |
|
mpi_info)); |
102 |
|
} else { |
103 |
|
Finley_Mesh_setFaceElements(out,Finley_ElementFile_alloc(Rec8, |
104 |
|
out->order, |
105 |
|
out->reduced_order, |
106 |
|
mpi_info)); |
107 |
|
Finley_Mesh_setContactElements(out,Finley_ElementFile_alloc(Rec8_Contact, |
108 |
|
out->order, |
109 |
|
out->reduced_order, |
110 |
|
mpi_info)); |
111 |
|
} |
112 |
|
} |
113 |
|
Finley_Mesh_setPoints(out,Finley_ElementFile_alloc(Point1, |
114 |
|
out->order, |
115 |
|
out->reduced_order, |
116 |
|
mpi_info)); |
117 |
if (! Finley_noError()) { |
if (! Finley_noError()) { |
118 |
Finley_Mesh_dealloc(out); |
Paso_MPIInfo_free( mpi_info ); |
119 |
|
Finley_Mesh_free(out); |
120 |
return NULL; |
return NULL; |
121 |
} |
} |
|
|
|
|
#pragma omp parallel for private(i0,i1,i2,k) |
|
|
for (i2=0;i2<N2;i2++) { |
|
|
for (i1=0;i1<N1;i1++) { |
|
|
for (i0=0;i0<N0;i0++) { |
|
|
k=i0+N0*i1+N0*N1*i2; |
|
|
out->Nodes->Coordinates[INDEX2(0,k,3)]=DBLE(i0)/DBLE(N0-1)*Length[0]; |
|
|
out->Nodes->Coordinates[INDEX2(1,k,3)]=DBLE(i1)/DBLE(N1-1)*Length[1]; |
|
|
out->Nodes->Coordinates[INDEX2(2,k,3)]=DBLE(i2)/DBLE(N2-1)*Length[2]; |
|
|
out->Nodes->Id[k]=k; |
|
|
out->Nodes->Tag[k]=0; |
|
|
out->Nodes->degreeOfFreedom[k]=(i0%NDOF0) +N0*(i1%NDOF1) +N0*N1*(i2%NDOF2); |
|
|
} |
|
|
} |
|
|
} |
|
122 |
|
|
123 |
/* tags for the faces: */ |
/* work out the largest dimension */ |
124 |
|
if (N2==MAX3(N0,N1,N2)) { |
125 |
|
Nstride0=1; |
126 |
|
Nstride1=N0; |
127 |
|
Nstride2=N0*N1; |
128 |
|
local_NE0=NE0; |
129 |
|
e_offset0=0; |
130 |
|
local_NE1=NE1; |
131 |
|
e_offset1=0; |
132 |
|
Paso_MPIInfo_Split(mpi_info,NE2,&local_NE2,&e_offset2); |
133 |
|
} else if (N1==MAX3(N0,N1,N2)) { |
134 |
|
Nstride0=N2; |
135 |
|
Nstride1=N0*N2; |
136 |
|
Nstride2=1; |
137 |
|
local_NE0=NE0; |
138 |
|
e_offset0=0; |
139 |
|
Paso_MPIInfo_Split(mpi_info,NE1,&local_NE1,&e_offset1); |
140 |
|
local_NE2=NE2; |
141 |
|
e_offset2=0; |
142 |
|
} else { |
143 |
|
Nstride0=N1*N2; |
144 |
|
Nstride1=1; |
145 |
|
Nstride2=N1; |
146 |
|
Paso_MPIInfo_Split(mpi_info,NE0,&local_NE0,&e_offset0); |
147 |
|
local_NE1=NE1; |
148 |
|
e_offset1=0; |
149 |
|
local_NE2=NE2; |
150 |
|
e_offset2=0; |
151 |
|
} |
152 |
|
offset0=e_offset0*N_PER_E; |
153 |
|
offset1=e_offset1*N_PER_E; |
154 |
|
offset2=e_offset2*N_PER_E; |
155 |
|
local_N0=local_NE0*N_PER_E+1; |
156 |
|
local_N1=local_NE1*N_PER_E+1; |
157 |
|
local_N2=local_NE2*N_PER_E+1; |
158 |
|
|
159 |
|
/* get the number of surface elements */ |
160 |
|
|
161 |
|
NFaceElements=0; |
162 |
if (!periodic[2]) { |
if (!periodic[2]) { |
163 |
for (i1=0;i1<N1;i1++) { |
NDOF2=N2; |
164 |
for (i0=0;i0<N0;i0++) { |
if (offset2==0) NFaceElements+=local_NE1*local_NE0; |
165 |
out->Nodes->Tag[i0+N0*i1+N0*N1*0]+=100; |
if (local_NE2+e_offset2 == NE2) NFaceElements+=local_NE1*local_NE0; |
166 |
out->Nodes->Tag[i0+N0*i1+N0*N1*(N2-1)]+=200; |
} else { |
167 |
} |
NDOF2=N2-1; |
|
} |
|
|
} |
|
|
if (!periodic[1]) { |
|
|
for (i2=0;i2<N2;i2++) { |
|
|
for (i0=0;i0<N0;i0++) { |
|
|
out->Nodes->Tag[i0+N0*0+N0*N1*i2]+=10; |
|
|
out->Nodes->Tag[i0+N0*(N1-1)+N0*N1*i2]+=20; |
|
|
} |
|
|
} |
|
168 |
} |
} |
169 |
|
|
170 |
if (!periodic[0]) { |
if (!periodic[0]) { |
171 |
for (i2=0;i2<N2;i2++) { |
NDOF0=N0; |
172 |
for (i1=0;i1<N1;i1++) { |
if (e_offset0 == 0) NFaceElements+=local_NE1*local_NE2; |
173 |
out->Nodes->Tag[0+N0*i1+N0*N1*i2]+=1; |
if (local_NE0+e_offset0 == NE0) NFaceElements+=local_NE1*local_NE2; |
|
out->Nodes->Tag[(N0-1)+N0*i1+N0*N1*i2]+=2; |
|
|
} |
|
|
} |
|
|
} |
|
|
|
|
|
/* set the elements: */ |
|
|
|
|
|
#pragma omp parallel for private(i0,i1,i2,k,node0) |
|
|
for (i2=0;i2<NE2;i2++) { |
|
|
for (i1=0;i1<NE1;i1++) { |
|
|
for (i0=0;i0<NE0;i0++) { |
|
|
k=i0+NE0*i1+NE0*NE1*i2; |
|
|
node0=2*i0+2*i1*N0+2*N0*N1*i2; |
|
|
|
|
|
out->Elements->Id[k]=k; |
|
|
out->Elements->Tag[k]=0; |
|
|
out->Elements->Color[k]=COLOR_MOD(i0)+3*COLOR_MOD(i1)+9*COLOR_MOD(i2); |
|
|
|
|
|
out->Elements->Nodes[INDEX2(0,k,20)]=node0; |
|
|
out->Elements->Nodes[INDEX2(1,k,20)]=node0+2; |
|
|
out->Elements->Nodes[INDEX2(2,k,20)]=node0+2*N0+2; |
|
|
out->Elements->Nodes[INDEX2(3,k,20)]=node0+2*N0; |
|
|
out->Elements->Nodes[INDEX2(4,k,20)]=node0+2*N0*N1; |
|
|
out->Elements->Nodes[INDEX2(5,k,20)]=node0+2*N0*N1+2; |
|
|
out->Elements->Nodes[INDEX2(6,k,20)]=node0+2*N0*N1+2*N0+2; |
|
|
out->Elements->Nodes[INDEX2(7,k,20)]=node0+2*N0*N1+2*N0; |
|
|
out->Elements->Nodes[INDEX2(8,k,20)]=node0+1; |
|
|
out->Elements->Nodes[INDEX2(9,k,20)]=node0+N0+2; |
|
|
out->Elements->Nodes[INDEX2(10,k,20)]=node0+2*N0+1; |
|
|
out->Elements->Nodes[INDEX2(11,k,20)]=node0+N0; |
|
|
out->Elements->Nodes[INDEX2(12,k,20)]=node0+N0*N1; |
|
|
out->Elements->Nodes[INDEX2(13,k,20)]=node0+N0*N1+2; |
|
|
out->Elements->Nodes[INDEX2(14,k,20)]=node0+N0*N1+2*N0+2; |
|
|
out->Elements->Nodes[INDEX2(15,k,20)]=node0+N0*N1+2*N0; |
|
|
out->Elements->Nodes[INDEX2(16,k,20)]=node0+2*N0*N1+1; |
|
|
out->Elements->Nodes[INDEX2(17,k,20)]=node0+2*N0*N1+N0+2; |
|
|
out->Elements->Nodes[INDEX2(18,k,20)]=node0+2*N0*N1+2*N0+1; |
|
|
out->Elements->Nodes[INDEX2(19,k,20)]=node0+2*N0*N1+N0; |
|
|
} |
|
|
} |
|
|
} |
|
|
out->Elements->minColor=0; |
|
|
out->Elements->maxColor=COLOR_MOD(0)+3*COLOR_MOD(0)+9*COLOR_MOD(0); |
|
|
|
|
|
/* face elements: */ |
|
|
|
|
|
if (useElementsOnFace) { |
|
|
NUMNODES=20; |
|
174 |
} else { |
} else { |
175 |
NUMNODES=8; |
NDOF0=N0-1; |
176 |
} |
} |
177 |
totalNECount=NE0*NE1*NE2; |
if (!periodic[1]) { |
178 |
faceNECount=0; |
NDOF1=N1; |
179 |
|
if (e_offset1 == 0) NFaceElements+=local_NE0*local_NE2; |
180 |
/* these are the quadrilateral elements on boundary 1 (x3=0): */ |
if (local_NE1+e_offset1 == NE1) NFaceElements+=local_NE0*local_NE2; |
181 |
|
} else { |
182 |
if (!periodic[2]) { |
NDOF1=N1-1; |
|
/* ** elements on boundary 100 (x3=0): */ |
|
|
#pragma omp parallel for private(i0,i1,k,node0) |
|
|
for (i1=0;i1<NE1;i1++) { |
|
|
for (i0=0;i0<NE0;i0++) { |
|
|
k=i0+NE0*i1+faceNECount; |
|
|
node0=2*i0+2*i1*N0; |
|
|
|
|
|
out->FaceElements->Id[k]=i0+NE0*i1+totalNECount; |
|
|
out->FaceElements->Tag[k]=100; |
|
|
out->FaceElements->Color[k]=(i0%2)+2*(i1%2); |
|
|
|
|
|
if (useElementsOnFace) { |
|
|
out->FaceElements->Nodes[INDEX2(0,k,NUMNODES)]=node0; |
|
|
out->FaceElements->Nodes[INDEX2(1,k,NUMNODES)]=node0+2*N0; |
|
|
out->FaceElements->Nodes[INDEX2(2,k,NUMNODES)]=node0+2*N0+2; |
|
|
out->FaceElements->Nodes[INDEX2(3,k,NUMNODES)]=node0+2; |
|
|
out->FaceElements->Nodes[INDEX2(4,k,NUMNODES)]=node0+2*N0*N1; |
|
|
out->FaceElements->Nodes[INDEX2(5,k,NUMNODES)]=node0+2*N0*N1+2*N0; |
|
|
out->FaceElements->Nodes[INDEX2(6,k,NUMNODES)]=node0+2*N0*N1+2*N0+2; |
|
|
out->FaceElements->Nodes[INDEX2(7,k,NUMNODES)]=node0+2*N0*N1+2; |
|
|
out->FaceElements->Nodes[INDEX2(8,k,NUMNODES)]=node0+N0; |
|
|
out->FaceElements->Nodes[INDEX2(9,k,NUMNODES)]=node0+2*N0+1; |
|
|
out->FaceElements->Nodes[INDEX2(10,k,NUMNODES)]=node0+N0+2; |
|
|
out->FaceElements->Nodes[INDEX2(11,k,NUMNODES)]=node0+1; |
|
|
out->FaceElements->Nodes[INDEX2(12,k,NUMNODES)]=node0+N0*N1; |
|
|
out->FaceElements->Nodes[INDEX2(13,k,NUMNODES)]=node0+N0*N1+2*N0; |
|
|
out->FaceElements->Nodes[INDEX2(14,k,NUMNODES)]=node0+N0*N1+2*N0+2; |
|
|
out->FaceElements->Nodes[INDEX2(15,k,NUMNODES)]=node0+N0*N1+2; |
|
|
out->FaceElements->Nodes[INDEX2(16,k,NUMNODES)]=node0+2*N0*N1+N0; |
|
|
out->FaceElements->Nodes[INDEX2(17,k,NUMNODES)]=node0+2*N0*N1+2*N0+1; |
|
|
out->FaceElements->Nodes[INDEX2(18,k,NUMNODES)]=node0+2*N0*N1+N0+2; |
|
|
out->FaceElements->Nodes[INDEX2(19,k,NUMNODES)]=node0+2*N0*N1+1; |
|
|
} else { |
|
|
out->FaceElements->Nodes[INDEX2(0,k,NUMNODES)]=node0; |
|
|
out->FaceElements->Nodes[INDEX2(1,k,NUMNODES)]=node0+2*N0; |
|
|
out->FaceElements->Nodes[INDEX2(2,k,NUMNODES)]=node0+2*N0+2; |
|
|
out->FaceElements->Nodes[INDEX2(3,k,NUMNODES)]=node0+2; |
|
|
out->FaceElements->Nodes[INDEX2(4,k,NUMNODES)]=node0+N0; |
|
|
out->FaceElements->Nodes[INDEX2(5,k,NUMNODES)]=node0+2*N0+1; |
|
|
out->FaceElements->Nodes[INDEX2(6,k,NUMNODES)]=node0+N0+2; |
|
|
out->FaceElements->Nodes[INDEX2(7,k,NUMNODES)]=node0+1; |
|
|
} |
|
|
} |
|
|
} |
|
|
totalNECount+=NE1*NE0; |
|
|
faceNECount+=NE1*NE0; |
|
|
|
|
|
/* ** elements on boundary 200 (x3=1) */ |
|
|
#pragma omp parallel for private(i0,i1,k,node0) |
|
|
for (i1=0;i1<NE1;i1++) { |
|
|
for (i0=0;i0<NE0;i0++) { |
|
|
k=i0+NE0*i1+faceNECount; |
|
|
node0=2*i0+2*i1*N0+2*N0*N1*(NE2-1); |
|
|
|
|
|
out->FaceElements->Id[k]=i0+NE0*i1+totalNECount; |
|
|
out->FaceElements->Tag[k]=200; |
|
|
out->FaceElements->Color[k]=(i0%2)+2*(i1%2)+4; |
|
|
|
|
|
if (useElementsOnFace) { |
|
|
out->FaceElements->Nodes[INDEX2(0,k,NUMNODES)]=node0+2*N0*N1; |
|
|
out->FaceElements->Nodes[INDEX2(1,k,NUMNODES)]=node0+2*N0*N1+2; |
|
|
out->FaceElements->Nodes[INDEX2(2,k,NUMNODES)]=node0+2*N0*N1+2*N0+2; |
|
|
out->FaceElements->Nodes[INDEX2(3,k,NUMNODES)]=node0+2*N0*N1+2*N0; |
|
|
|
|
|
out->FaceElements->Nodes[INDEX2(4,k,NUMNODES)]=node0; |
|
|
out->FaceElements->Nodes[INDEX2(5,k,NUMNODES)]=node0+2; |
|
|
out->FaceElements->Nodes[INDEX2(6,k,NUMNODES)]=node0+2*N0+2; |
|
|
out->FaceElements->Nodes[INDEX2(7,k,NUMNODES)]=node0+2*N0; |
|
|
|
|
|
out->FaceElements->Nodes[INDEX2(8,k,NUMNODES)]=node0+2*N0*N1+1; |
|
|
out->FaceElements->Nodes[INDEX2(9,k,NUMNODES)]=node0+2*N0*N1+N0+2; |
|
|
out->FaceElements->Nodes[INDEX2(10,k,NUMNODES)]=node0+2*N0*N1+2*N0+1; |
|
|
out->FaceElements->Nodes[INDEX2(11,k,NUMNODES)]=node0+2*N0*N1+N0; |
|
|
|
|
|
out->FaceElements->Nodes[INDEX2(12,k,NUMNODES)]=node0+N0*N1; |
|
|
out->FaceElements->Nodes[INDEX2(13,k,NUMNODES)]=node0+N0*N1+2; |
|
|
out->FaceElements->Nodes[INDEX2(14,k,NUMNODES)]=node0+N0*N1+2*N0+2; |
|
|
out->FaceElements->Nodes[INDEX2(15,k,NUMNODES)]=node0+N0*N1+2*N0; |
|
|
|
|
|
out->FaceElements->Nodes[INDEX2(16,k,NUMNODES)]=node0+1; |
|
|
out->FaceElements->Nodes[INDEX2(17,k,NUMNODES)]=node0+N0+2; |
|
|
out->FaceElements->Nodes[INDEX2(18,k,NUMNODES)]=node0+2*N0+1; |
|
|
out->FaceElements->Nodes[INDEX2(19,k,NUMNODES)]=node0+N0; |
|
|
|
|
|
} else { |
|
|
out->FaceElements->Nodes[INDEX2(0,k,NUMNODES)]=node0+2*N0*N1; |
|
|
out->FaceElements->Nodes[INDEX2(1,k,NUMNODES)]=node0+2*N0*N1+2; |
|
|
out->FaceElements->Nodes[INDEX2(2,k,NUMNODES)]=node0+2*N0*N1+2*N0+2; |
|
|
out->FaceElements->Nodes[INDEX2(3,k,NUMNODES)]=node0+2*N0*N1+2*N0; |
|
|
out->FaceElements->Nodes[INDEX2(4,k,NUMNODES)]=node0+2*N0*N1+1; |
|
|
out->FaceElements->Nodes[INDEX2(5,k,NUMNODES)]=node0+2*N0*N1+N0+2; |
|
|
out->FaceElements->Nodes[INDEX2(6,k,NUMNODES)]=node0+2*N0*N1+2*N0+1; |
|
|
out->FaceElements->Nodes[INDEX2(7,k,NUMNODES)]=node0+2*N0*N1+N0; |
|
|
} |
|
|
} |
|
|
} |
|
|
totalNECount+=NE1*NE0; |
|
|
faceNECount+=NE1*NE0; |
|
183 |
} |
} |
184 |
if (!periodic[0]) { |
|
185 |
/* ** elements on boundary 001 (x1=0): */ |
/* allocate tables: */ |
186 |
|
|
187 |
#pragma omp parallel for private(i1,i2,k,node0) |
Finley_NodeFile_allocTable(out->Nodes,local_N0*local_N1*local_N2); |
188 |
for (i2=0;i2<NE2;i2++) { |
Finley_ElementFile_allocTable(out->Elements,local_NE0*local_NE1*local_NE2); |
189 |
for (i1=0;i1<NE1;i1++) { |
Finley_ElementFile_allocTable(out->FaceElements,NFaceElements); |
190 |
k=i1+NE1*i2+faceNECount; |
|
191 |
node0=2*i1*N0+2*N0*N1*i2; |
if (Finley_noError()) { |
192 |
|
/* create nodes */ |
193 |
|
|
194 |
out->FaceElements->Id[k]=i1+NE1*i2+totalNECount; |
#pragma omp parallel for private(i0,i1,i2,k,global_i0,global_i1,global_i2) |
195 |
out->FaceElements->Tag[k]=1; |
for (i2=0;i2<local_N2;i2++) { |
196 |
out->FaceElements->Color[k]=(i2%2)+2*(i1%2)+8; |
for (i1=0;i1<local_N1;i1++) { |
197 |
|
for (i0=0;i0<local_N0;i0++) { |
198 |
if (useElementsOnFace) { |
k=i0+local_N0*i1+local_N0*local_N1*i2; |
199 |
out->FaceElements->Nodes[INDEX2(0,k,NUMNODES)]=node0; |
global_i0=i0+offset0; |
200 |
out->FaceElements->Nodes[INDEX2(1,k,NUMNODES)]=node0+2*N0*N1; |
global_i1=i1+offset1; |
201 |
out->FaceElements->Nodes[INDEX2(2,k,NUMNODES)]=node0+2*N0*N1+2*N0; |
global_i2=i2+offset2; |
202 |
out->FaceElements->Nodes[INDEX2(3,k,NUMNODES)]=node0+2*N0; |
out->Nodes->Coordinates[INDEX2(0,k,DIM)]=DBLE(global_i0)/DBLE(N0-1)*Length[0]; |
203 |
|
out->Nodes->Coordinates[INDEX2(1,k,DIM)]=DBLE(global_i1)/DBLE(N1-1)*Length[1]; |
204 |
out->FaceElements->Nodes[INDEX2(4,k,NUMNODES)]=node0+2; |
out->Nodes->Coordinates[INDEX2(2,k,DIM)]=DBLE(global_i2)/DBLE(N2-1)*Length[2]; |
205 |
out->FaceElements->Nodes[INDEX2(5,k,NUMNODES)]=node0+2*N0*N1+2; |
out->Nodes->Id[k]=Nstride0*global_i0+Nstride1*global_i1+Nstride2*global_i2; |
206 |
out->FaceElements->Nodes[INDEX2(6,k,NUMNODES)]=node0+2*N0*N1+2*N0+2; |
out->Nodes->Tag[k]=0; |
207 |
out->FaceElements->Nodes[INDEX2(7,k,NUMNODES)]=node0+2*N0+2; |
out->Nodes->globalDegreesOfFreedom[k]=Nstride0*(global_i0%NDOF0) |
208 |
|
+Nstride1*(global_i1%NDOF1) |
209 |
out->FaceElements->Nodes[INDEX2(8,k,NUMNODES)]=node0+N0*N1; |
+Nstride2*(global_i2%NDOF2); |
|
out->FaceElements->Nodes[INDEX2(9,k,NUMNODES)]=node0+2*N0*N1+N0; |
|
|
out->FaceElements->Nodes[INDEX2(10,k,NUMNODES)]=node0+N0*N1+2*N0; |
|
|
out->FaceElements->Nodes[INDEX2(11,k,NUMNODES)]=node0+N0; |
|
|
|
|
|
out->FaceElements->Nodes[INDEX2(12,k,NUMNODES)]=node0+1; |
|
|
out->FaceElements->Nodes[INDEX2(13,k,NUMNODES)]=node0+2*N0*N1+1; |
|
|
out->FaceElements->Nodes[INDEX2(14,k,NUMNODES)]=node0+2*N0*N1+2*N0+1; |
|
|
out->FaceElements->Nodes[INDEX2(15,k,NUMNODES)]=node0+2*N0+1; |
|
|
|
|
|
out->FaceElements->Nodes[INDEX2(16,k,NUMNODES)]=node0+N0*N1+2; |
|
|
out->FaceElements->Nodes[INDEX2(17,k,NUMNODES)]=node0+2*N0*N1+N0+2; |
|
|
out->FaceElements->Nodes[INDEX2(18,k,NUMNODES)]=node0+N0*N1+2*N0+2; |
|
|
out->FaceElements->Nodes[INDEX2(19,k,NUMNODES)]=node0+N0+2; |
|
|
|
|
|
} else { |
|
|
out->FaceElements->Nodes[INDEX2(0,k,NUMNODES)]=node0; |
|
|
out->FaceElements->Nodes[INDEX2(1,k,NUMNODES)]=node0+2*N0*N1; |
|
|
out->FaceElements->Nodes[INDEX2(2,k,NUMNODES)]=node0+2*N0*N1+2*N0; |
|
|
out->FaceElements->Nodes[INDEX2(3,k,NUMNODES)]=node0+2*N0; |
|
|
out->FaceElements->Nodes[INDEX2(4,k,NUMNODES)]=node0+N0*N1; |
|
|
out->FaceElements->Nodes[INDEX2(5,k,NUMNODES)]=node0+2*N0*N1+N0; |
|
|
out->FaceElements->Nodes[INDEX2(6,k,NUMNODES)]=node0+N0*N1+2*N0; |
|
|
out->FaceElements->Nodes[INDEX2(7,k,NUMNODES)]=node0+N0; |
|
210 |
} |
} |
211 |
} |
} |
212 |
} |
} |
213 |
totalNECount+=NE1*NE2; |
/* set the elements: */ |
214 |
faceNECount+=NE1*NE2; |
NN=out->Elements->numNodes; |
215 |
|
#pragma omp parallel for private(i0,i1,i2,k,node0) |
216 |
/* ** elements on boundary 002 (x1=1): */ |
for (i2=0;i2<local_NE2;i2++) { |
217 |
|
for (i1=0;i1<local_NE1;i1++) { |
218 |
#pragma omp parallel for private(i1,i2,k,node0) |
for (i0=0;i0<local_NE0;i0++) { |
219 |
for (i2=0;i2<NE2;i2++) { |
|
220 |
for (i1=0;i1<NE1;i1++) { |
k=i0+local_NE0*i1+local_NE0*local_NE1*i2; |
221 |
k=i1+NE1*i2+faceNECount; |
node0=Nstride0*N_PER_E*(i0+e_offset0)+Nstride1*N_PER_E*(i1+e_offset1)+Nstride2*N_PER_E*(i2+e_offset2); |
|
node0=2*(NE0-1)+2*i1*N0+2*N0*N1*i2 ; |
|
222 |
|
|
223 |
out->FaceElements->Id[k]=i1+NE1*i2+totalNECount; |
out->Elements->Id[k]=(i0+e_offset0)+NE0*(i1+e_offset1)+NE0*NE1*(i2+e_offset2); |
224 |
out->FaceElements->Tag[k]=2; |
out->Elements->Tag[k]=0; |
225 |
out->FaceElements->Color[k]=(i2%2)+2*(i1%2)+12; |
out->Elements->Owner[k]=myRank; |
226 |
|
|
227 |
if (useElementsOnFace) { |
out->Elements->Nodes[INDEX2(0,k,NN)] =node0 ; |
228 |
out->FaceElements->Nodes[INDEX2(0,k,NUMNODES)]=node0+2; |
out->Elements->Nodes[INDEX2(1,k,NN)] =node0+ 2*Nstride0; |
229 |
out->FaceElements->Nodes[INDEX2(1,k,NUMNODES)]=node0+2*N0+2; |
out->Elements->Nodes[INDEX2(2,k,NN)] =node0+ 2*Nstride1+2*Nstride0; |
230 |
out->FaceElements->Nodes[INDEX2(2,k,NUMNODES)]=node0+2*N0*N1+2*N0+2; |
out->Elements->Nodes[INDEX2(3,k,NN)] =node0+ 2*Nstride1; |
231 |
out->FaceElements->Nodes[INDEX2(3,k,NUMNODES)]=node0+2*N0*N1+2; |
out->Elements->Nodes[INDEX2(4,k,NN)] =node0+2*Nstride2 ; |
232 |
|
out->Elements->Nodes[INDEX2(5,k,NN)] =node0+2*Nstride2 +2*Nstride0; |
233 |
out->FaceElements->Nodes[INDEX2(4,k,NUMNODES)]=node0; |
out->Elements->Nodes[INDEX2(6,k,NN)] =node0+2*Nstride2+2*Nstride1+2*Nstride0; |
234 |
out->FaceElements->Nodes[INDEX2(5,k,NUMNODES)]=node0+2*N0; |
out->Elements->Nodes[INDEX2(7,k,NN)] =node0+2*Nstride2+2*Nstride1 ; |
235 |
out->FaceElements->Nodes[INDEX2(6,k,NUMNODES)]=node0+2*N0*N1+2*N0; |
out->Elements->Nodes[INDEX2(8,k,NN)] =node0+ 1*Nstride0; |
236 |
out->FaceElements->Nodes[INDEX2(7,k,NUMNODES)]=node0+2*N0*N1; |
out->Elements->Nodes[INDEX2(9,k,NN)] =node0+ 1*Nstride1+2*Nstride0; |
237 |
|
out->Elements->Nodes[INDEX2(10,k,NN)]=node0+ 2*Nstride1+1*Nstride0; |
238 |
out->FaceElements->Nodes[INDEX2(8,k,NUMNODES)]=node0+N0+2; |
out->Elements->Nodes[INDEX2(11,k,NN)]=node0+ 1*Nstride1 ; |
239 |
out->FaceElements->Nodes[INDEX2(9,k,NUMNODES)]=node0+N0*N1+2*N0+2; |
out->Elements->Nodes[INDEX2(12,k,NN)]=node0+1*Nstride2 ; |
240 |
out->FaceElements->Nodes[INDEX2(10,k,NUMNODES)]=node0+2*N0*N1+N0+2; |
out->Elements->Nodes[INDEX2(13,k,NN)]=node0+1*Nstride2 +2*Nstride0; |
241 |
out->FaceElements->Nodes[INDEX2(11,k,NUMNODES)]=node0+N0*N1+2; |
out->Elements->Nodes[INDEX2(14,k,NN)]=node0+1*Nstride2+2*Nstride1+2*Nstride0; |
242 |
|
out->Elements->Nodes[INDEX2(15,k,NN)]=node0+1*Nstride2+2*Nstride1 ; |
243 |
out->FaceElements->Nodes[INDEX2(12,k,NUMNODES)]=node0+1; |
out->Elements->Nodes[INDEX2(16,k,NN)]=node0+2*Nstride2 +1*Nstride0; |
244 |
out->FaceElements->Nodes[INDEX2(13,k,NUMNODES)]=node0+2*N0+1; |
out->Elements->Nodes[INDEX2(17,k,NN)]=node0+2*Nstride2+1*Nstride1+2*Nstride0; |
245 |
out->FaceElements->Nodes[INDEX2(14,k,NUMNODES)]=node0+2*N0*N1+2*N0+1; |
out->Elements->Nodes[INDEX2(18,k,NN)]=node0+2*Nstride2+2*Nstride1+1*Nstride0; |
246 |
out->FaceElements->Nodes[INDEX2(15,k,NUMNODES)]=node0+2*N0*N1+1; |
out->Elements->Nodes[INDEX2(19,k,NN)]=node0+2*Nstride2+1*Nstride1 ; |
247 |
|
if (useFullElementOrder) { |
248 |
out->FaceElements->Nodes[INDEX2(16,k,NUMNODES)]=node0+N0; |
out->Elements->Nodes[INDEX2(20,k,NN)]=node0+ 1*Nstride1+1*Nstride0; |
249 |
out->FaceElements->Nodes[INDEX2(17,k,NUMNODES)]=node0+N0*N1+2*N0; |
out->Elements->Nodes[INDEX2(21,k,NN)]=node0+1*Nstride2 +1*Nstride0; |
250 |
out->FaceElements->Nodes[INDEX2(18,k,NUMNODES)]=node0+2*N0*N1+N0; |
out->Elements->Nodes[INDEX2(22,k,NN)]=node0+1*Nstride2+1*Nstride1+2*Nstride0; |
251 |
out->FaceElements->Nodes[INDEX2(19,k,NUMNODES)]=node0+N0*N1; |
out->Elements->Nodes[INDEX2(23,k,NN)]=node0+1*Nstride2+2*Nstride1+1*Nstride0; |
252 |
|
out->Elements->Nodes[INDEX2(24,k,NN)]=node0+1*Nstride2+1*Nstride1 ; |
253 |
} else { |
out->Elements->Nodes[INDEX2(25,k,NN)]=node0+2*Nstride2+1*Nstride1+1*Nstride0; |
254 |
out->FaceElements->Nodes[INDEX2(0,k,NUMNODES)]=node0+2; |
out->Elements->Nodes[INDEX2(26,k,NN)]=node0+1*Nstride2+1*Nstride1+1*Nstride0; |
255 |
out->FaceElements->Nodes[INDEX2(1,k,NUMNODES)]=node0+2*N0+2; |
} |
|
out->FaceElements->Nodes[INDEX2(2,k,NUMNODES)]=node0+2*N0*N1+2*N0+2; |
|
|
out->FaceElements->Nodes[INDEX2(3,k,NUMNODES)]=node0+2*N0*N1+2; |
|
|
out->FaceElements->Nodes[INDEX2(4,k,NUMNODES)]=node0+N0+2; |
|
|
out->FaceElements->Nodes[INDEX2(5,k,NUMNODES)]=node0+N0*N1+2*N0+2; |
|
|
out->FaceElements->Nodes[INDEX2(6,k,NUMNODES)]=node0+2*N0*N1+N0+2; |
|
|
out->FaceElements->Nodes[INDEX2(7,k,NUMNODES)]=node0+N0*N1+2; |
|
256 |
} |
} |
|
|
|
257 |
} |
} |
258 |
} |
} |
259 |
totalNECount+=NE1*NE2; |
/* face elements */ |
260 |
faceNECount+=NE1*NE2; |
NN=out->FaceElements->numNodes; |
261 |
} |
totalNECount=NE0*NE1*NE2; |
262 |
if (!periodic[1]) { |
faceNECount=0; |
263 |
/* ** elements on boundary 010 (x2=0): */ |
/* these are the quadrilateral elements on boundary 1 (x3=0): */ |
264 |
|
if (!periodic[2]) { |
265 |
#pragma omp parallel for private(i0,i2,k,node0) |
/* ** elements on boundary 100 (x3=0): */ |
266 |
for (i2=0;i2<NE2;i2++) { |
if (offset2==0) { |
267 |
for (i0=0;i0<NE0;i0++) { |
#pragma omp parallel for private(i0,i1,k,node0) |
268 |
k=i0+NE0*i2+faceNECount; |
for (i1=0;i1<local_NE1;i1++) { |
269 |
node0=2*i0+2*N0*N1*i2; |
for (i0=0;i0<local_NE0;i0++) { |
270 |
|
|
271 |
out->FaceElements->Id[k]=i2+NE2*i0+totalNECount; |
k=i0+local_NE0*i1+faceNECount; |
272 |
out->FaceElements->Tag[k]=10; |
node0=Nstride0*N_PER_E*(i0+e_offset0)+Nstride1*N_PER_E*(i1+e_offset1); |
273 |
out->FaceElements->Color[k]=(i2%2)+2*(i0%2)+16; |
|
274 |
|
out->FaceElements->Id[k]=(i0+e_offset0)+NE0*(i1+e_offset1)+totalNECount; |
275 |
if (useElementsOnFace) { |
out->FaceElements->Tag[k]=100; |
276 |
out->FaceElements->Nodes[INDEX2(0,k,NUMNODES)]=node0; |
out->FaceElements->Owner[k]=myRank; |
277 |
out->FaceElements->Nodes[INDEX2(1,k,NUMNODES)]=node0+2; |
|
278 |
out->FaceElements->Nodes[INDEX2(2,k,NUMNODES)]=node0+2*N1*N0+2; |
if (useElementsOnFace) { |
279 |
out->FaceElements->Nodes[INDEX2(3,k,NUMNODES)]=node0+2*N1*N0; |
out->FaceElements->Nodes[INDEX2(0,k,NN)] =node0 ; |
280 |
|
out->FaceElements->Nodes[INDEX2(1,k,NN)] =node0 +2*Nstride1 ; |
281 |
out->FaceElements->Nodes[INDEX2(4,k,NUMNODES)]=node0+2*N0; |
out->FaceElements->Nodes[INDEX2(2,k,NN)] =node0 +2*Nstride1+2*Nstride0; |
282 |
out->FaceElements->Nodes[INDEX2(5,k,NUMNODES)]=node0+2*N0+2; |
out->FaceElements->Nodes[INDEX2(3,k,NN)] =node0+ 2*Nstride0 ; |
283 |
out->FaceElements->Nodes[INDEX2(6,k,NUMNODES)]=node0+2*N1*N0+2*N0+2; |
out->FaceElements->Nodes[INDEX2(4,k,NN)] =node0+2*Nstride2 ; |
284 |
out->FaceElements->Nodes[INDEX2(7,k,NUMNODES)]=node0+2*N1*N0+2*N0; |
out->FaceElements->Nodes[INDEX2(5,k,NN)] =node0+2*Nstride2+2*Nstride1 ; |
285 |
|
out->FaceElements->Nodes[INDEX2(6,k,NN)] =node0+2*Nstride2+2*Nstride1+2*Nstride0; |
286 |
out->FaceElements->Nodes[INDEX2(8,k,NUMNODES)]=node0+1; |
out->FaceElements->Nodes[INDEX2(7,k,NN)] =node0+2*Nstride2 +2*Nstride0; |
287 |
out->FaceElements->Nodes[INDEX2(9,k,NUMNODES)]=node0+N0*N1+2; |
out->FaceElements->Nodes[INDEX2(8,k,NN)] =node0+ 1*Nstride1 ; |
288 |
out->FaceElements->Nodes[INDEX2(10,k,NUMNODES)]=node0+2*N1*N0+1; |
out->FaceElements->Nodes[INDEX2(9,k,NN)] =node0+ 2*Nstride1+1*Nstride0; |
289 |
out->FaceElements->Nodes[INDEX2(11,k,NUMNODES)]=node0+N1*N0; |
out->FaceElements->Nodes[INDEX2(10,k,NN)]=node0+ 1*Nstride1+2*Nstride0; |
290 |
|
out->FaceElements->Nodes[INDEX2(11,k,NN)]=node0+ 1*Nstride0; |
291 |
out->FaceElements->Nodes[INDEX2(12,k,NUMNODES)]=node0+N0; |
out->FaceElements->Nodes[INDEX2(12,k,NN)]=node0+1*Nstride2 ; |
292 |
out->FaceElements->Nodes[INDEX2(13,k,NUMNODES)]=node0+N0+2; |
out->FaceElements->Nodes[INDEX2(13,k,NN)]=node0+1*Nstride2+2*Nstride1 ; |
293 |
out->FaceElements->Nodes[INDEX2(14,k,NUMNODES)]=node0+2*N1*N0+N0+2; |
out->FaceElements->Nodes[INDEX2(14,k,NN)]=node0+1*Nstride2+2*Nstride1+2*Nstride0; |
294 |
out->FaceElements->Nodes[INDEX2(15,k,NUMNODES)]=node0+2*N1*N0+N0; |
out->FaceElements->Nodes[INDEX2(15,k,NN)]=node0+1*Nstride2 +2*Nstride0; |
295 |
|
out->FaceElements->Nodes[INDEX2(16,k,NN)]=node0+2*Nstride2+1*Nstride1; |
296 |
out->FaceElements->Nodes[INDEX2(16,k,NUMNODES)]=node0+2*N0+1; |
out->FaceElements->Nodes[INDEX2(17,k,NN)]=node0+2*Nstride2+2*Nstride1+1*Nstride0; |
297 |
out->FaceElements->Nodes[INDEX2(17,k,NUMNODES)]=node0+N0*N1+2*N0+2; |
out->FaceElements->Nodes[INDEX2(18,k,NN)]=node0+2*Nstride2+1*Nstride1+2*Nstride0; |
298 |
out->FaceElements->Nodes[INDEX2(18,k,NUMNODES)]=node0+2*N1*N0+2*N0+1; |
out->FaceElements->Nodes[INDEX2(19,k,NN)]=node0+2*Nstride2 +1*Nstride0; |
299 |
out->FaceElements->Nodes[INDEX2(19,k,NUMNODES)]=node0+N1*N0+2*N0; |
} else { |
300 |
|
out->FaceElements->Nodes[INDEX2(0,k,NN)] =node0 ; |
301 |
} else { |
out->FaceElements->Nodes[INDEX2(1,k,NN)] =node0+ 2*Nstride1 ; |
302 |
out->FaceElements->Nodes[INDEX2(0,k,NUMNODES)]=node0; |
out->FaceElements->Nodes[INDEX2(2,k,NN)] =node0+ 2*Nstride1+2*Nstride0; |
303 |
out->FaceElements->Nodes[INDEX2(1,k,NUMNODES)]=node0+2; |
out->FaceElements->Nodes[INDEX2(3,k,NN)] =node0+ 2*Nstride0; |
304 |
out->FaceElements->Nodes[INDEX2(2,k,NUMNODES)]=node0+2*N1*N0+2; |
out->FaceElements->Nodes[INDEX2(4,k,NN)] =node0+ 1*Nstride1 ; |
305 |
out->FaceElements->Nodes[INDEX2(3,k,NUMNODES)]=node0+2*N1*N0; |
out->FaceElements->Nodes[INDEX2(5,k,NN)] =node0+ 2*Nstride1+1*Nstride0; |
306 |
out->FaceElements->Nodes[INDEX2(4,k,NUMNODES)]=node0+1; |
out->FaceElements->Nodes[INDEX2(6,k,NN)] =node0+ 1*Nstride1+2*Nstride0; |
307 |
out->FaceElements->Nodes[INDEX2(5,k,NUMNODES)]=node0+N0*N1+2; |
out->FaceElements->Nodes[INDEX2(7,k,NN)] =node0+ 1*Nstride0; |
308 |
out->FaceElements->Nodes[INDEX2(6,k,NUMNODES)]=node0+2*N1*N0+1; |
if (useFullElementOrder){ |
309 |
out->FaceElements->Nodes[INDEX2(7,k,NUMNODES)]=node0+N1*N0; |
out->FaceElements->Nodes[INDEX2(8,k,NN)] =node0+ 1*Nstride1+1*Nstride0; |
310 |
} |
} |
311 |
|
} |
312 |
|
} |
313 |
|
} |
314 |
|
faceNECount+=local_NE1*local_NE0; |
315 |
|
} |
316 |
|
totalNECount+=NE1*NE0; |
317 |
|
/* ** elements on boundary 200 (x3=1) */ |
318 |
|
if (local_NE2+e_offset2 == NE2) { |
319 |
|
#pragma omp parallel for private(i0,i1,k,node0) |
320 |
|
for (i1=0;i1<local_NE1;i1++) { |
321 |
|
for (i0=0;i0<local_NE0;i0++) { |
322 |
|
|
323 |
|
k=i0+local_NE0*i1+faceNECount; |
324 |
|
node0=Nstride0*N_PER_E*(i0+e_offset0)+Nstride1*N_PER_E*(i1+e_offset1)+Nstride2*N_PER_E*(NE2-1); |
325 |
|
|
326 |
|
out->FaceElements->Id[k]=(i0+e_offset0)+NE0*(i1+e_offset1)+totalNECount; |
327 |
|
out->FaceElements->Tag[k]=200; |
328 |
|
out->FaceElements->Owner[k]=myRank; |
329 |
|
if (useElementsOnFace) { |
330 |
|
out->FaceElements->Nodes[INDEX2(0,k,NN)] =node0+2*Nstride2 ; |
331 |
|
out->FaceElements->Nodes[INDEX2(1,k,NN)] =node0+2*Nstride2+ 2*Nstride0; |
332 |
|
out->FaceElements->Nodes[INDEX2(2,k,NN)] =node0+2*Nstride2+2*Nstride1+2*Nstride0; |
333 |
|
out->FaceElements->Nodes[INDEX2(3,k,NN)] =node0+2*Nstride2+2*Nstride1 ; |
334 |
|
|
335 |
|
out->FaceElements->Nodes[INDEX2(4,k,NN)] =node0 ; |
336 |
|
out->FaceElements->Nodes[INDEX2(5,k,NN)] =node0+2*Nstride0 ; |
337 |
|
out->FaceElements->Nodes[INDEX2(6,k,NN)] =node0+ 2*Nstride1+2*Nstride0; |
338 |
|
out->FaceElements->Nodes[INDEX2(7,k,NN)] =node0+ 2*Nstride1; |
339 |
|
|
340 |
|
out->FaceElements->Nodes[INDEX2(8,k,NN)] =node0+2*Nstride2+ 1*Nstride0; |
341 |
|
out->FaceElements->Nodes[INDEX2(9,k,NN)] =node0+2*Nstride2+1*Nstride1+2*Nstride0; |
342 |
|
out->FaceElements->Nodes[INDEX2(10,k,NN)]=node0+2*Nstride2+2*Nstride1+1*Nstride0; |
343 |
|
out->FaceElements->Nodes[INDEX2(11,k,NN)]=node0+2*Nstride2+1*Nstride1 ; |
344 |
|
|
345 |
|
out->FaceElements->Nodes[INDEX2(12,k,NN)]=node0+1*Nstride2; |
346 |
|
out->FaceElements->Nodes[INDEX2(13,k,NN)]=node0+1*Nstride2 +2*Nstride0; |
347 |
|
out->FaceElements->Nodes[INDEX2(14,k,NN)]=node0+1*Nstride2+2*Nstride1+2*Nstride0; |
348 |
|
out->FaceElements->Nodes[INDEX2(15,k,NN)]=node0+1*Nstride2+2*Nstride1 ; |
349 |
|
|
350 |
|
out->FaceElements->Nodes[INDEX2(16,k,NN)]=node0+ 1*Nstride0; |
351 |
|
out->FaceElements->Nodes[INDEX2(17,k,NN)]=node0+ 1*Nstride1+2*Nstride0; |
352 |
|
out->FaceElements->Nodes[INDEX2(18,k,NN)]=node0+ 2*Nstride1+1*Nstride0; |
353 |
|
out->FaceElements->Nodes[INDEX2(19,k,NN)]=node0+ 1*Nstride1 ; |
354 |
|
|
355 |
|
} else { |
356 |
|
out->FaceElements->Nodes[INDEX2(0,k,NN)] =node0+2*Nstride2 ; |
357 |
|
out->FaceElements->Nodes[INDEX2(1,k,NN)] =node0+2*Nstride2 +2*Nstride0; |
358 |
|
out->FaceElements->Nodes[INDEX2(2,k,NN)] =node0+2*Nstride2+2*Nstride1+2*Nstride0; |
359 |
|
out->FaceElements->Nodes[INDEX2(3,k,NN)] =node0+2*Nstride2+2*Nstride1 ; |
360 |
|
out->FaceElements->Nodes[INDEX2(4,k,NN)] =node0+2*Nstride2 +1*Nstride0; |
361 |
|
out->FaceElements->Nodes[INDEX2(5,k,NN)] =node0+2*Nstride2+1*Nstride1+2*Nstride0; |
362 |
|
out->FaceElements->Nodes[INDEX2(6,k,NN)] =node0+2*Nstride2+2*Nstride1+1*Nstride0; |
363 |
|
out->FaceElements->Nodes[INDEX2(7,k,NN)] =node0+2*Nstride2+1*Nstride1 ; |
364 |
|
if (useFullElementOrder){ |
365 |
|
out->FaceElements->Nodes[INDEX2(8,k,NN)] =node0+2*Nstride2+1*Nstride1+1*Nstride0; |
366 |
|
} |
367 |
|
} |
368 |
|
} |
369 |
|
} |
370 |
|
faceNECount+=local_NE1*local_NE0; |
371 |
} |
} |
372 |
|
totalNECount+=NE1*NE0; |
373 |
} |
} |
374 |
totalNECount+=NE0*NE2; |
if (!periodic[0]) { |
375 |
faceNECount+=NE0*NE2; |
/* ** elements on boundary 001 (x1=0): */ |
376 |
|
|
377 |
/* ** elements on boundary 020 (x2=1): */ |
if (e_offset0 == 0) { |
378 |
|
#pragma omp parallel for private(i1,i2,k,node0) |
379 |
|
for (i2=0;i2<local_NE2;i2++) { |
380 |
|
for (i1=0;i1<local_NE1;i1++) { |
381 |
|
|
382 |
|
k=i1+local_NE1*i2+faceNECount; |
383 |
|
node0=Nstride1*N_PER_E*(i1+e_offset1)+Nstride2*N_PER_E*(i2+e_offset2); |
384 |
|
out->FaceElements->Id[k]=(i1+e_offset1)+NE1*(i2+e_offset2)+totalNECount; |
385 |
|
out->FaceElements->Tag[k]=1; |
386 |
|
out->FaceElements->Owner[k]=myRank; |
387 |
|
|
388 |
|
if (useElementsOnFace) { |
389 |
|
out->FaceElements->Nodes[INDEX2(0,k,NN)] =node0 ; |
390 |
|
out->FaceElements->Nodes[INDEX2(1,k,NN)] =node0+2*Nstride2 ; |
391 |
|
out->FaceElements->Nodes[INDEX2(2,k,NN)] =node0+2*Nstride2+2*Nstride1 ; |
392 |
|
out->FaceElements->Nodes[INDEX2(3,k,NN)] =node0+2*Nstride1 ; |
393 |
|
|
394 |
|
out->FaceElements->Nodes[INDEX2(4,k,NN)] =node0+2*Nstride0 ; |
395 |
|
out->FaceElements->Nodes[INDEX2(5,k,NN)] =node0+2*Nstride2+2*Nstride0 ; |
396 |
|
out->FaceElements->Nodes[INDEX2(6,k,NN)] =node0+2*Nstride2+2*Nstride1+2*Nstride0; |
397 |
|
out->FaceElements->Nodes[INDEX2(7,k,NN)] =node0+2*Nstride1+2*Nstride0 ; |
398 |
|
|
399 |
|
out->FaceElements->Nodes[INDEX2(8,k,NN)] =node0+1*Nstride2 ; |
400 |
|
out->FaceElements->Nodes[INDEX2(9,k,NN)] =node0+2*Nstride2+1*Nstride1 ; |
401 |
|
out->FaceElements->Nodes[INDEX2(10,k,NN)]=node0+1*Nstride2+2*Nstride1 ; |
402 |
|
out->FaceElements->Nodes[INDEX2(11,k,NN)]=node0+ 1*Nstride1 ; |
403 |
|
|
404 |
|
out->FaceElements->Nodes[INDEX2(12,k,NN)]=node0+ 1*Nstride0; |
405 |
|
out->FaceElements->Nodes[INDEX2(13,k,NN)]=node0+2*Nstride2 +1*Nstride0; |
406 |
|
out->FaceElements->Nodes[INDEX2(14,k,NN)]=node0+2*Nstride2+2*Nstride1+1*Nstride0; |
407 |
|
out->FaceElements->Nodes[INDEX2(15,k,NN)]=node0+2*Nstride1+ 1*Nstride0; |
408 |
|
|
409 |
|
out->FaceElements->Nodes[INDEX2(16,k,NN)]=node0+1*Nstride2+ 2*Nstride0; |
410 |
|
out->FaceElements->Nodes[INDEX2(17,k,NN)]=node0+2*Nstride2+1*Nstride1+2*Nstride0; |
411 |
|
out->FaceElements->Nodes[INDEX2(18,k,NN)]=node0+1*Nstride2+2*Nstride1+2*Nstride0; |
412 |
|
out->FaceElements->Nodes[INDEX2(19,k,NN)]=node0+1*Nstride1+ 2*Nstride0; |
413 |
|
|
414 |
|
} else { |
415 |
|
out->FaceElements->Nodes[INDEX2(0,k,NN)] =node0 ; |
416 |
|
out->FaceElements->Nodes[INDEX2(1,k,NN)] =node0+2*Nstride2 ; |
417 |
|
out->FaceElements->Nodes[INDEX2(2,k,NN)] =node0+2*Nstride2+2*Nstride1 ; |
418 |
|
out->FaceElements->Nodes[INDEX2(3,k,NN)] =node0+ 2*Nstride1 ; |
419 |
|
out->FaceElements->Nodes[INDEX2(4,k,NN)] =node0+1*Nstride2 ; |
420 |
|
out->FaceElements->Nodes[INDEX2(5,k,NN)] =node0+2*Nstride2+1*Nstride1 ; |
421 |
|
out->FaceElements->Nodes[INDEX2(6,k,NN)] =node0+1*Nstride2+2*Nstride1 ; |
422 |
|
out->FaceElements->Nodes[INDEX2(7,k,NN)] =node0+ 1*Nstride1 ; |
423 |
|
if (useFullElementOrder){ |
424 |
|
out->FaceElements->Nodes[INDEX2(8,k,NN)] =node0+1*Nstride2+1*Nstride1 ; |
425 |
|
} |
426 |
|
} |
427 |
|
} |
428 |
|
} |
429 |
|
faceNECount+=local_NE1*local_NE2; |
430 |
|
} |
431 |
|
totalNECount+=NE1*NE2; |
432 |
|
|
433 |
#pragma omp parallel for private(i0,i2,k,node0) |
/* ** elements on boundary 002 (x1=1): */ |
434 |
for (i2=0;i2<NE2;i2++) { |
if (local_NE0+e_offset0 == NE0) { |
435 |
for (i0=0;i0<NE0;i0++) { |
#pragma omp parallel for private(i1,i2,k,node0) |
436 |
k=i0+NE0*i2+faceNECount; |
for (i2=0;i2<local_NE2;i2++) { |
437 |
node0=2*i0+2*(NE1-1)*N0+2*N0*N1*i2; |
for (i1=0;i1<local_NE1;i1++) { |
438 |
|
k=i1+local_NE1*i2+faceNECount; |
439 |
|
node0=Nstride0*N_PER_E*(NE0-1)+Nstride1*N_PER_E*(i1+e_offset1)+Nstride2*N_PER_E*(i2+e_offset2); |
440 |
|
out->FaceElements->Id[k]=(i1+e_offset1)+NE1*(i2+e_offset2)+totalNECount; |
441 |
|
out->FaceElements->Tag[k]=2; |
442 |
|
out->FaceElements->Owner[k]=myRank; |
443 |
|
|
444 |
out->FaceElements->Id[k]=i2+NE2*i0+totalNECount; |
if (useElementsOnFace) { |
445 |
out->FaceElements->Tag[k]=20; |
out->FaceElements->Nodes[INDEX2(0,k,NN)]=node0+ 2*Nstride0; |
446 |
out->FaceElements->Color[k]=(i2%2)+2*(i0%2)+20; |
out->FaceElements->Nodes[INDEX2(1,k,NN)]=node0+ 2*Nstride1+2*Nstride0; |
447 |
|
out->FaceElements->Nodes[INDEX2(2,k,NN)]=node0+2*Nstride2+2*Nstride1+2*Nstride0; |
448 |
if (useElementsOnFace) { |
out->FaceElements->Nodes[INDEX2(3,k,NN)]=node0+2*Nstride2+ 2*Nstride0; |
449 |
out->FaceElements->Nodes[INDEX2(0,k,NUMNODES)]=node0+2*N0; |
|
450 |
out->FaceElements->Nodes[INDEX2(1,k,NUMNODES)]=node0+2*N0*N1+2*N0; |
out->FaceElements->Nodes[INDEX2(4,k,NN)]=node0 ; |
451 |
out->FaceElements->Nodes[INDEX2(2,k,NUMNODES)]=node0+2*N0*N1+2*N0+2; |
out->FaceElements->Nodes[INDEX2(5,k,NN)]=node0+ 2*Nstride1 ; |
452 |
out->FaceElements->Nodes[INDEX2(3,k,NUMNODES)]=node0+2*N0+2; |
out->FaceElements->Nodes[INDEX2(6,k,NN)]=node0+2*Nstride2+2*Nstride1 ; |
453 |
|
out->FaceElements->Nodes[INDEX2(7,k,NN)]=node0+2*Nstride2 ; |
454 |
out->FaceElements->Nodes[INDEX2(4,k,NUMNODES)]=node0; |
|
455 |
out->FaceElements->Nodes[INDEX2(5,k,NUMNODES)]=node0+2*N0*N1; |
out->FaceElements->Nodes[INDEX2(8,k,NN)]=node0+ 1*Nstride1+2*Nstride0; |
456 |
out->FaceElements->Nodes[INDEX2(6,k,NUMNODES)]=node0+2*N0*N1+2; |
out->FaceElements->Nodes[INDEX2(9,k,NN)]=node0+1*Nstride2+2*Nstride1+2*Nstride0; |
457 |
out->FaceElements->Nodes[INDEX2(7,k,NUMNODES)]=node0+2; |
out->FaceElements->Nodes[INDEX2(10,k,NN)]=node0+2*Nstride2+1*Nstride1+2*Nstride0; |
458 |
|
out->FaceElements->Nodes[INDEX2(11,k,NN)]=node0+1*Nstride2+ 2*Nstride0; |
459 |
out->FaceElements->Nodes[INDEX2(8,k,NUMNODES)]=node0+N1*N0+2*N0; |
|
460 |
out->FaceElements->Nodes[INDEX2(9,k,NUMNODES)]=node0+2*N1*N0+2*N0+1; |
out->FaceElements->Nodes[INDEX2(12,k,NN)]=node0+ 1*Nstride0; |
461 |
out->FaceElements->Nodes[INDEX2(10,k,NUMNODES)]=node0+N0*N1+2*N0+2; |
out->FaceElements->Nodes[INDEX2(13,k,NN)]=node0+ 2*Nstride1+1*Nstride0; |
462 |
out->FaceElements->Nodes[INDEX2(11,k,NUMNODES)]=node0+2*N0+1; |
out->FaceElements->Nodes[INDEX2(14,k,NN)]=node0+2*Nstride2+2*Nstride1+1*Nstride0; |
463 |
|
out->FaceElements->Nodes[INDEX2(15,k,NN)]=node0+2*Nstride2+ 1*Nstride0; |
464 |
out->FaceElements->Nodes[INDEX2(12,k,NUMNODES)]=node0+N0; |
|
465 |
out->FaceElements->Nodes[INDEX2(13,k,NUMNODES)]=node0+2*N0*N1+N0; |
out->FaceElements->Nodes[INDEX2(16,k,NN)]=node0+ 1*Nstride1 ; |
466 |
out->FaceElements->Nodes[INDEX2(14,k,NUMNODES)]=node0+2*N0*N1+N0+2; |
out->FaceElements->Nodes[INDEX2(17,k,NN)]=node0+1*Nstride2+2*Nstride1 ; |
467 |
out->FaceElements->Nodes[INDEX2(15,k,NUMNODES)]=node0+N0+2; |
out->FaceElements->Nodes[INDEX2(18,k,NN)]=node0+2*Nstride2+1*Nstride1 ; |
468 |
|
out->FaceElements->Nodes[INDEX2(19,k,NN)]=node0+1*Nstride2 ; |
469 |
out->FaceElements->Nodes[INDEX2(16,k,NUMNODES)]=node0+N1*N0; |
|
470 |
out->FaceElements->Nodes[INDEX2(17,k,NUMNODES)]=node0+2*N1*N0+1; |
} else { |
471 |
out->FaceElements->Nodes[INDEX2(18,k,NUMNODES)]=node0+N0*N1+2; |
out->FaceElements->Nodes[INDEX2(0,k,NN)]=node0 +2*Nstride0; |
472 |
out->FaceElements->Nodes[INDEX2(19,k,NUMNODES)]=node0+1; |
out->FaceElements->Nodes[INDEX2(1,k,NN)]=node0+ 2*Nstride1+2*Nstride0; |
473 |
} else { |
out->FaceElements->Nodes[INDEX2(2,k,NN)]=node0+2*Nstride2+2*Nstride1+2*Nstride0; |
474 |
out->FaceElements->Nodes[INDEX2(0,k,NUMNODES)]=node0+2*N0; |
out->FaceElements->Nodes[INDEX2(3,k,NN)]=node0+2*Nstride2+ 2*Nstride0; |
475 |
out->FaceElements->Nodes[INDEX2(1,k,NUMNODES)]=node0+2*N0*N1+2*N0; |
out->FaceElements->Nodes[INDEX2(4,k,NN)]=node0+ 1*Nstride1+2*Nstride0; |
476 |
out->FaceElements->Nodes[INDEX2(2,k,NUMNODES)]=node0+2*N0*N1+2*N0+2; |
out->FaceElements->Nodes[INDEX2(5,k,NN)]=node0+1*Nstride2+2*Nstride1+2*Nstride0; |
477 |
out->FaceElements->Nodes[INDEX2(3,k,NUMNODES)]=node0+2*N0+2; |
out->FaceElements->Nodes[INDEX2(6,k,NN)]=node0+2*Nstride2+1*Nstride1+2*Nstride0; |
478 |
out->FaceElements->Nodes[INDEX2(4,k,NUMNODES)]=node0+N1*N0+2*N0; |
out->FaceElements->Nodes[INDEX2(7,k,NN)]=node0+1*Nstride2 +2*Nstride0; |
479 |
out->FaceElements->Nodes[INDEX2(5,k,NUMNODES)]=node0+2*N1*N0+2*N0+1; |
if (useFullElementOrder){ |
480 |
out->FaceElements->Nodes[INDEX2(6,k,NUMNODES)]=node0+N0*N1+2*N0+2; |
out->FaceElements->Nodes[INDEX2(8,k,NN)] =node0+1*Nstride2+1*Nstride1+2*Nstride0; |
481 |
out->FaceElements->Nodes[INDEX2(7,k,NUMNODES)]=node0+2*N0+1; |
} |
482 |
|
} |
483 |
|
|
484 |
|
} |
485 |
|
} |
486 |
|
faceNECount+=local_NE1*local_NE2; |
487 |
} |
} |
488 |
} |
totalNECount+=NE1*NE2; |
489 |
|
} |
490 |
|
if (!periodic[1]) { |
491 |
|
/* ** elements on boundary 010 (x2=0): */ |
492 |
|
if (e_offset1 == 0) { |
493 |
|
#pragma omp parallel for private(i0,i2,k,node0) |
494 |
|
for (i2=0;i2<local_NE2;i2++) { |
495 |
|
for (i0=0;i0<local_NE0;i0++) { |
496 |
|
k=i0+local_NE0*i2+faceNECount; |
497 |
|
node0=Nstride0*N_PER_E*(i0+e_offset0)+Nstride2*N_PER_E*(i2+e_offset2); |
498 |
|
|
499 |
|
out->FaceElements->Id[k]=(i2+e_offset2)+NE2*(e_offset0+i0)+totalNECount; |
500 |
|
out->FaceElements->Tag[k]=10; |
501 |
|
out->FaceElements->Owner[k]=myRank; |
502 |
|
if (useElementsOnFace) { |
503 |
|
out->FaceElements->Nodes[INDEX2(0,k,NN)]=node0 ; |
504 |
|
out->FaceElements->Nodes[INDEX2(1,k,NN)]=node0+ 2*Nstride0; |
505 |
|
out->FaceElements->Nodes[INDEX2(2,k,NN)]=node0+2*Nstride2 +2*Nstride0; |
506 |
|
out->FaceElements->Nodes[INDEX2(3,k,NN)]=node0+2*Nstride2 ; |
507 |
|
|
508 |
|
out->FaceElements->Nodes[INDEX2(4,k,NN)]=node0+ 2*Nstride1 ; |
509 |
|
out->FaceElements->Nodes[INDEX2(5,k,NN)]=node0+2*Nstride1+ 2*Nstride0; |
510 |
|
out->FaceElements->Nodes[INDEX2(6,k,NN)]=node0+2*Nstride2+2*Nstride1+2*Nstride0; |
511 |
|
out->FaceElements->Nodes[INDEX2(7,k,NN)]=node0+2*Nstride2+2*Nstride1 ; |
512 |
|
|
513 |
|
out->FaceElements->Nodes[INDEX2(8,k,NN)]=node0+ 1*Nstride0; |
514 |
|
out->FaceElements->Nodes[INDEX2(9,k,NN)]=node0+1*Nstride2+ 2*Nstride0; |
515 |
|
out->FaceElements->Nodes[INDEX2(10,k,NN)]=node0+2*Nstride2+ 1*Nstride0; |
516 |
|
out->FaceElements->Nodes[INDEX2(11,k,NN)]=node0+1*Nstride2 ; |
517 |
|
|
518 |
|
out->FaceElements->Nodes[INDEX2(12,k,NN)]=node0+ 1*Nstride1 ; |
519 |
|
out->FaceElements->Nodes[INDEX2(13,k,NN)]=node0+ 1*Nstride1+2*Nstride0; |
520 |
|
out->FaceElements->Nodes[INDEX2(14,k,NN)]=node0+2*Nstride2+1*Nstride1+2*Nstride0; |
521 |
|
out->FaceElements->Nodes[INDEX2(15,k,NN)]=node0+2*Nstride2+1*Nstride1 ; |
522 |
|
|
523 |
|
out->FaceElements->Nodes[INDEX2(16,k,NN)]=node0+ 2*Nstride1+1*Nstride0; |
524 |
|
out->FaceElements->Nodes[INDEX2(17,k,NN)]=node0+1*Nstride2+2*Nstride1+2*Nstride0; |
525 |
|
out->FaceElements->Nodes[INDEX2(18,k,NN)]=node0+2*Nstride2+2*Nstride1+1*Nstride0; |
526 |
|
out->FaceElements->Nodes[INDEX2(19,k,NN)]=node0+1*Nstride2+2*Nstride1 ; |
527 |
|
|
528 |
|
} else { |
529 |
|
out->FaceElements->Nodes[INDEX2(0,k,NN)]=node0 ; |
530 |
|
out->FaceElements->Nodes[INDEX2(1,k,NN)]=node0+ 2*Nstride0; |
531 |
|
out->FaceElements->Nodes[INDEX2(2,k,NN)]=node0+2*Nstride2+ 2*Nstride0; |
532 |
|
out->FaceElements->Nodes[INDEX2(3,k,NN)]=node0+2*Nstride2 ; |
533 |
|
out->FaceElements->Nodes[INDEX2(4,k,NN)]=node0+ 1*Nstride0; |
534 |
|
out->FaceElements->Nodes[INDEX2(5,k,NN)]=node0+1*Nstride2+ 2*Nstride0; |
535 |
|
out->FaceElements->Nodes[INDEX2(6,k,NN)]=node0+2*Nstride2+ 1*Nstride0; |
536 |
|
out->FaceElements->Nodes[INDEX2(7,k,NN)]=node0+1*Nstride2 ; |
537 |
|
if (useFullElementOrder){ |
538 |
|
out->FaceElements->Nodes[INDEX2(8,k,NN)] =node0+1*Nstride2+ 1*Nstride0; |
539 |
|
} |
540 |
|
} |
541 |
|
} |
542 |
|
} |
543 |
|
faceNECount+=local_NE0*local_NE2; |
544 |
|
} |
545 |
|
totalNECount+=NE0*NE2; |
546 |
|
/* ** elements on boundary 020 (x2=1): */ |
547 |
|
if (local_NE1+e_offset1 == NE1) { |
548 |
|
#pragma omp parallel for private(i0,i2,k,node0) |
549 |
|
for (i2=0;i2<local_NE2;i2++) { |
550 |
|
for (i0=0;i0<local_NE0;i0++) { |
551 |
|
k=i0+local_NE0*i2+faceNECount; |
552 |
|
node0=Nstride0*N_PER_E*(i0+e_offset0)+Nstride1*N_PER_E*(NE1-1)+Nstride2*N_PER_E*(i2+e_offset2); |
553 |
|
|
554 |
|
out->FaceElements->Id[k]=(i2+e_offset2)+NE2*(i0+e_offset0)+totalNECount; |
555 |
|
out->FaceElements->Tag[k]=20; |
556 |
|
out->FaceElements->Owner[k]=myRank; |
557 |
|
|
558 |
|
if (useElementsOnFace) { |
559 |
|
out->FaceElements->Nodes[INDEX2(0,k,NN)]=node0+ 2*Nstride1 ; |
560 |
|
out->FaceElements->Nodes[INDEX2(1,k,NN)]=node0+2*Nstride2+2*Nstride1 ; |
561 |
|
out->FaceElements->Nodes[INDEX2(2,k,NN)]=node0+2*Nstride2+2*Nstride1+2*Nstride0; |
562 |
|
out->FaceElements->Nodes[INDEX2(3,k,NN)]=node0+2*Nstride1+2*Nstride0 ; |
563 |
|
|
564 |
|
out->FaceElements->Nodes[INDEX2(4,k,NN)]=node0 ; |
565 |
|
out->FaceElements->Nodes[INDEX2(5,k,NN)]=node0+2*Nstride2 ; |
566 |
|
out->FaceElements->Nodes[INDEX2(6,k,NN)]=node0+2*Nstride2+ 2*Nstride0; |
567 |
|
out->FaceElements->Nodes[INDEX2(7,k,NN)]=node0+ 2*Nstride0; |
568 |
|
|
569 |
|
out->FaceElements->Nodes[INDEX2(8,k,NN)]=node0+1*Nstride2+2*Nstride1 ; |
570 |
|
out->FaceElements->Nodes[INDEX2(9,k,NN)]=node0+2*Nstride2+2*Nstride1+1*Nstride0; |
571 |
|
out->FaceElements->Nodes[INDEX2(10,k,NN)]=node0+1*Nstride2+2*Nstride1+2*Nstride0; |
572 |
|
out->FaceElements->Nodes[INDEX2(11,k,NN)]=node0+ 2*Nstride1+1*Nstride0; |
573 |
|
|
574 |
|
out->FaceElements->Nodes[INDEX2(12,k,NN)]=node0+ 1*Nstride1 ; |
575 |
|
out->FaceElements->Nodes[INDEX2(13,k,NN)]=node0+2*Nstride2+1*Nstride1 ; |
576 |
|
out->FaceElements->Nodes[INDEX2(14,k,NN)]=node0+2*Nstride2+1*Nstride1+2*Nstride0; |
577 |
|
out->FaceElements->Nodes[INDEX2(15,k,NN)]=node0+ 1*Nstride1+2*Nstride0; |
578 |
|
|
579 |
|
out->FaceElements->Nodes[INDEX2(16,k,NN)]=node0+1*Nstride2 ; |
580 |
|
out->FaceElements->Nodes[INDEX2(17,k,NN)]=node0+2*Nstride2 +1*Nstride0; |
581 |
|
out->FaceElements->Nodes[INDEX2(18,k,NN)]=node0+1*Nstride2 +2*Nstride0; |
582 |
|
out->FaceElements->Nodes[INDEX2(19,k,NN)]=node0+ 1*Nstride0; |
583 |
|
} else { |
584 |
|
out->FaceElements->Nodes[INDEX2(0,k,NN)]=node0+ 2*Nstride1 ; |
585 |
|
out->FaceElements->Nodes[INDEX2(1,k,NN)]=node0+2*Nstride2+2*Nstride1 ; |
586 |
|
out->FaceElements->Nodes[INDEX2(2,k,NN)]=node0+2*Nstride2+2*Nstride1+2*Nstride0; |
587 |
|
out->FaceElements->Nodes[INDEX2(3,k,NN)]=node0+ 2*Nstride1+2*Nstride0; |
588 |
|
out->FaceElements->Nodes[INDEX2(4,k,NN)]=node0+1*Nstride2+2*Nstride1 ; |
589 |
|
out->FaceElements->Nodes[INDEX2(5,k,NN)]=node0+2*Nstride2+2*Nstride1+1*Nstride0; |
590 |
|
out->FaceElements->Nodes[INDEX2(6,k,NN)]=node0+1*Nstride2+2*Nstride1+2*Nstride0; |
591 |
|
out->FaceElements->Nodes[INDEX2(7,k,NN)]=node0+ 2*Nstride1+1*Nstride0; |
592 |
|
if (useFullElementOrder){ |
593 |
|
out->FaceElements->Nodes[INDEX2(8,k,NN)] =node0+1*Nstride2+2*Nstride1+1*Nstride0; |
594 |
|
} |
595 |
|
} |
596 |
|
} |
597 |
|
} |
598 |
|
faceNECount+=local_NE0*local_NE2; |
599 |
|
} |
600 |
|
totalNECount+=NE0*NE2; |
601 |
|
} |
602 |
|
/* add tag names */ |
603 |
|
Finley_Mesh_addTagMap(out,"top", 200); |
604 |
|
Finley_Mesh_addTagMap(out,"bottom", 100); |
605 |
|
Finley_Mesh_addTagMap(out,"left", 1); |
606 |
|
Finley_Mesh_addTagMap(out,"right", 2); |
607 |
|
Finley_Mesh_addTagMap(out,"front", 10); |
608 |
|
Finley_Mesh_addTagMap(out,"back", 20); |
609 |
|
|
610 |
|
/* prepare mesh for further calculatuions:*/ |
611 |
|
if (Finley_noError()) { |
612 |
|
Finley_Mesh_resolveNodeIds(out); |
613 |
|
} |
614 |
|
if (Finley_noError()) { |
615 |
|
Finley_Mesh_prepare(out, optimize); |
616 |
} |
} |
|
totalNECount+=NE0*NE2; |
|
|
faceNECount+=NE0*NE2; |
|
617 |
} |
} |
|
out->FaceElements->minColor=0; |
|
|
out->FaceElements->maxColor=24; |
|
|
|
|
|
/* face elements done: */ |
|
|
|
|
|
/* condense the nodes: */ |
|
|
|
|
|
Finley_Mesh_resolveNodeIds(out); |
|
|
|
|
|
/* prepare mesh for further calculatuions:*/ |
|
|
Finley_Mesh_prepare(out) ; |
|
618 |
|
|
619 |
|
if (!Finley_noError()) { |
620 |
|
Finley_Mesh_free(out); |
621 |
|
} |
622 |
|
/* free up memory */ |
623 |
|
Paso_MPIInfo_free( mpi_info ); |
624 |
#ifdef Finley_TRACE |
#ifdef Finley_TRACE |
625 |
printf("timing: mesh generation: %.4e sec\n",Finley_timer()-time0); |
printf("timing: mesh generation: %.4e sec\n",Finley_timer()-time0); |
626 |
#endif |
#endif |
627 |
|
|
|
if (! Finley_noError()) { |
|
|
Finley_Mesh_dealloc(out); |
|
|
return NULL; |
|
|
} |
|
628 |
return out; |
return out; |
629 |
} |
} |
|
|
|
|
/* |
|
|
* Revision 1.3 2005/09/01 03:31:35 jgs |
|
|
* Merge of development branch dev-02 back to main trunk on 2005-09-01 |
|
|
* |
|
|
* Revision 1.2.2.2 2005/09/07 06:26:19 gross |
|
|
* the solver from finley are put into the standalone package paso now |
|
|
* |
|
|
* Revision 1.2.2.1 2005/08/24 02:02:18 gross |
|
|
* timing output switched off. solver output can be swiched through getSolution(verbose=True) now. |
|
|
* |
|
|
* Revision 1.2 2005/07/08 04:07:52 jgs |
|
|
* Merge of development branch back to main trunk on 2005-07-08 |
|
|
* |
|
|
* Revision 1.1.1.1.2.1 2005/06/29 02:34:51 gross |
|
|
* some changes towards 64 integers in finley |
|
|
* |
|
|
* Revision 1.1.1.1 2004/10/26 06:53:57 jgs |
|
|
* initial import of project esys2 |
|
|
* |
|
|
* Revision 1.1.1.1 2004/06/24 04:00:40 johng |
|
|
* Initial version of eys using boost-python. |
|
|
* |
|
|
* |
|
|
*/ |
|
|
|
|