/[escript]/branches/domexper/dudley/src/Mesh_hex20.c
ViewVC logotype

Diff of /branches/domexper/dudley/src/Mesh_hex20.c

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

temp/finley/src/Mesh_hex20.c revision 1387 by trankine, Fri Jan 11 07:45:26 2008 UTC trunk/finley/src/Mesh_hex20.c revision 2722 by gross, Fri Oct 16 06:45:01 2009 UTC
# Line 1  Line 1 
1    
 /* $Id$ */  
   
2  /*******************************************************  /*******************************************************
3   *  *
4   *           Copyright 2003-2007 by ACceSS MNRF  * Copyright (c) 2003-2009 by University of Queensland
5   *       Copyright 2007 by University of Queensland  * Earth Systems Science Computational Center (ESSCC)
6   *  * http://www.uq.edu.au/esscc
7   *                http://esscc.uq.edu.au  *
8   *        Primary Business: Queensland, Australia  * Primary Business: Queensland, Australia
9   *  Licensed under the Open Software License version 3.0  * Licensed under the Open Software License version 3.0
10   *     http://www.opensource.org/licenses/osl-3.0.php  * http://www.opensource.org/licenses/osl-3.0.php
11   *  *
12   *******************************************************/  *******************************************************/
13    
14    
15  /**************************************************************/  /**************************************************************/
16    
# Line 33  Finley_Mesh* Finley_RectangularMesh_Hex2 Line 32  Finley_Mesh* Finley_RectangularMesh_Hex2
32                                            index_t reduced_order,                                            index_t reduced_order,
33                                            bool_t useElementsOnFace,                                            bool_t useElementsOnFace,
34                                            bool_t useFullElementOrder,                                            bool_t useFullElementOrder,
35                                              bool_t useMacroElements,
36                                            bool_t optimize)                                            bool_t optimize)
37  {  {
38    #define N_PER_E 2    #define N_PER_E 2
# Line 43  Finley_Mesh* Finley_RectangularMesh_Hex2 Line 43  Finley_Mesh* Finley_RectangularMesh_Hex2
43    Finley_Mesh* out;    Finley_Mesh* out;
44    Paso_MPIInfo *mpi_info = NULL;    Paso_MPIInfo *mpi_info = NULL;
45    char name[50];    char name[50];
46      bool_t generateAllNodes= useFullElementOrder || useMacroElements;
47      #ifdef Finley_TRACE
48    double time0=Finley_timer();    double time0=Finley_timer();
49      #endif
50    
51    /* get MPI information */    /* get MPI information */
52    mpi_info = Paso_MPIInfo_alloc( MPI_COMM_WORLD );    mpi_info = Paso_MPIInfo_alloc( MPI_COMM_WORLD );
# Line 62  Finley_Mesh* Finley_RectangularMesh_Hex2 Line 65  Finley_Mesh* Finley_RectangularMesh_Hex2
65    N2=N_PER_E*NE2+1;    N2=N_PER_E*NE2+1;
66    
67    /*  allocate mesh: */      /*  allocate mesh: */  
68    sprintf(name,"Rectangular %d x %d x %d mesh",N0,N1,N2);    sprintf(name,"Brick %d x %d x %d mesh",N0,N1,N2);
69    out=Finley_Mesh_alloc(name,DIM,order, reduced_order, mpi_info);    out=Finley_Mesh_alloc(name,DIM,order, reduced_order, mpi_info);
70    if (! Finley_noError()) {    if (! Finley_noError()) {
71        Paso_MPIInfo_free( mpi_info );        Paso_MPIInfo_free( mpi_info );
72        return NULL;        return NULL;
73    }    }
74    
75    if (useFullElementOrder) {    if (generateAllNodes) {
76       Finley_setError(SYSTEM_ERROR,"full element order for Hex elements is not supported yet.");       /* Finley_setError(SYSTEM_ERROR,"full element order for Hex elements is not supported yet."); */
77  /*       if (useMacroElements) {
78       Finley_Mesh_setElements(out,Finley_ElementFile_alloc(Hex27,            Finley_Mesh_setElements(out,Finley_ElementFile_alloc(Hex27Macro,
79                                              out->order,                                                   out->order,
80                                              out->reduced_order,                                                   out->reduced_order,
81                                              mpi_info))`;                                                   mpi_info));
82         } else {
83              Finley_Mesh_setElements(out,Finley_ElementFile_alloc(Hex27,
84                                                     out->order,
85                                                     out->reduced_order,
86                                                     mpi_info));
87         }
88       if (useElementsOnFace) {       if (useElementsOnFace) {
89           Finley_setError(SYSTEM_ERROR,"rich elements for Hex27 elements is not supported yet.");           Finley_setError(SYSTEM_ERROR,"rich elements for Hex27 elements is not supported yet.");
90       } else {       } else {
91           Finley_Mesh_setFaceElements(out,Finley_ElementFile_alloc(Rec9,           if (useMacroElements) {
92                                                      out->order,               Finley_Mesh_setFaceElements(out,Finley_ElementFile_alloc(Rec9Macro,
93                                                      out->reduced_order,                                                          out->order,
94                                                      mpi_info));                                                          out->reduced_order,
95                                                            mpi_info));
96             } else {
97                 Finley_Mesh_setFaceElements(out,Finley_ElementFile_alloc(Rec9,
98                                                            out->order,
99                                                            out->reduced_order,
100                                                            mpi_info));
101             }
102           Finley_Mesh_setContactElements(out,Finley_ElementFile_alloc(Rec9_Contact,           Finley_Mesh_setContactElements(out,Finley_ElementFile_alloc(Rec9_Contact,
103                                                         out->order,                                                         out->order,
104                                                         out->reduced_order,                                                         out->reduced_order,
105                                                         mpi_info));                                                         mpi_info));
106       }       }
 */  
107    
108    } else  {    } else  {
109       Finley_Mesh_setElements(out,Finley_ElementFile_alloc(Hex20,out->order,out->reduced_order,mpi_info));       Finley_Mesh_setElements(out,Finley_ElementFile_alloc(Hex20,out->order,out->reduced_order,mpi_info));
# Line 154  Finley_Mesh* Finley_RectangularMesh_Hex2 Line 169  Finley_Mesh* Finley_RectangularMesh_Hex2
169    offset0=e_offset0*N_PER_E;    offset0=e_offset0*N_PER_E;
170    offset1=e_offset1*N_PER_E;    offset1=e_offset1*N_PER_E;
171    offset2=e_offset2*N_PER_E;    offset2=e_offset2*N_PER_E;
172    local_N0=local_NE0*N_PER_E+1;    local_N0=local_NE0>0 ? local_NE0*N_PER_E+1 : 0;
173    local_N1=local_NE1*N_PER_E+1;    local_N1=local_NE1>0 ? local_NE1*N_PER_E+1 : 0;
174    local_N2=local_NE2*N_PER_E+1;    local_N2=local_NE2>0 ? local_NE2*N_PER_E+1 : 0;
175    
176    /* get the number of surface elements */    /* get the number of surface elements */
177    
178    NFaceElements=0;    NFaceElements=0;
179    if (!periodic[2]) {    if (!periodic[2] && (local_NE2>0) ) {
180      NDOF2=N2;      NDOF2=N2;
181      if (offset2==0) NFaceElements+=local_NE1*local_NE0;      if (offset2==0) NFaceElements+=local_NE1*local_NE0;
182      if (local_NE2+e_offset2 == NE2) NFaceElements+=local_NE1*local_NE0;      if (local_NE2+e_offset2 == NE2) NFaceElements+=local_NE1*local_NE0;
# Line 169  Finley_Mesh* Finley_RectangularMesh_Hex2 Line 184  Finley_Mesh* Finley_RectangularMesh_Hex2
184        NDOF2=N2-1;        NDOF2=N2-1;
185    }    }
186    
187    if (!periodic[0]) {    if (!periodic[0] && (local_NE0>0) ) {
188       NDOF0=N0;       NDOF0=N0;
189       if (e_offset0 == 0) NFaceElements+=local_NE1*local_NE2;       if (e_offset0 == 0) NFaceElements+=local_NE1*local_NE2;
190       if (local_NE0+e_offset0 == NE0) NFaceElements+=local_NE1*local_NE2;       if (local_NE0+e_offset0 == NE0) NFaceElements+=local_NE1*local_NE2;
191    } else {    } else {
192        NDOF0=N0-1;        NDOF0=N0-1;
193    }    }
194    if (!periodic[1]) {    if (!periodic[1] && (local_NE1>0) ) {
195       NDOF1=N1;       NDOF1=N1;
196       if (e_offset1 == 0) NFaceElements+=local_NE0*local_NE2;       if (e_offset1 == 0) NFaceElements+=local_NE0*local_NE2;
197       if (local_NE1+e_offset1 == NE1) NFaceElements+=local_NE0*local_NE2;       if (local_NE1+e_offset1 == NE1) NFaceElements+=local_NE0*local_NE2;
198    } else {    } else {
199        NDOF1=N1-1;        NDOF1=N1-1;
200    }    }
   
201    /*  allocate tables: */    /*  allocate tables: */
202    
203    Finley_NodeFile_allocTable(out->Nodes,local_N0*local_N1*local_N2);    Finley_NodeFile_allocTable(out->Nodes,local_N0*local_N1*local_N2);
# Line 246  Finley_Mesh* Finley_RectangularMesh_Hex2 Line 260  Finley_Mesh* Finley_RectangularMesh_Hex2
260             out->Elements->Nodes[INDEX2(17,k,NN)]=node0+2*Nstride2+1*Nstride1+2*Nstride0;             out->Elements->Nodes[INDEX2(17,k,NN)]=node0+2*Nstride2+1*Nstride1+2*Nstride0;
261             out->Elements->Nodes[INDEX2(18,k,NN)]=node0+2*Nstride2+2*Nstride1+1*Nstride0;             out->Elements->Nodes[INDEX2(18,k,NN)]=node0+2*Nstride2+2*Nstride1+1*Nstride0;
262             out->Elements->Nodes[INDEX2(19,k,NN)]=node0+2*Nstride2+1*Nstride1           ;             out->Elements->Nodes[INDEX2(19,k,NN)]=node0+2*Nstride2+1*Nstride1           ;
263             if (useFullElementOrder) {             if (generateAllNodes) {
264                out->Elements->Nodes[INDEX2(20,k,NN)]=node0+           1*Nstride1+1*Nstride0;                out->Elements->Nodes[INDEX2(20,k,NN)]=node0+           1*Nstride1+1*Nstride0;
265                out->Elements->Nodes[INDEX2(21,k,NN)]=node0+1*Nstride2           +1*Nstride0;                out->Elements->Nodes[INDEX2(21,k,NN)]=node0+1*Nstride2           +1*Nstride0;
266                out->Elements->Nodes[INDEX2(22,k,NN)]=node0+1*Nstride2+1*Nstride1+2*Nstride0;                out->Elements->Nodes[INDEX2(22,k,NN)]=node0+1*Nstride2+1*Nstride1+2*Nstride0;
# Line 263  Finley_Mesh* Finley_RectangularMesh_Hex2 Line 277  Finley_Mesh* Finley_RectangularMesh_Hex2
277       totalNECount=NE0*NE1*NE2;       totalNECount=NE0*NE1*NE2;
278       faceNECount=0;       faceNECount=0;
279       /*   these are the quadrilateral elements on boundary 1 (x3=0): */       /*   these are the quadrilateral elements on boundary 1 (x3=0): */
280       if (!periodic[2]) {       if (!periodic[2] && (local_NE2>0)) {
281         /* **  elements on boundary 100 (x3=0): */         /* **  elements on boundary 100 (x3=0): */
282         if (offset2==0) {         if (offset2==0) {
283            #pragma omp parallel for private(i0,i1,k,node0)            #pragma omp parallel for private(i0,i1,k,node0)
# Line 307  Finley_Mesh* Finley_RectangularMesh_Hex2 Line 321  Finley_Mesh* Finley_RectangularMesh_Hex2
321                   out->FaceElements->Nodes[INDEX2(5,k,NN)] =node0+           2*Nstride1+1*Nstride0;                   out->FaceElements->Nodes[INDEX2(5,k,NN)] =node0+           2*Nstride1+1*Nstride0;
322                   out->FaceElements->Nodes[INDEX2(6,k,NN)] =node0+           1*Nstride1+2*Nstride0;                   out->FaceElements->Nodes[INDEX2(6,k,NN)] =node0+           1*Nstride1+2*Nstride0;
323                   out->FaceElements->Nodes[INDEX2(7,k,NN)] =node0+                      1*Nstride0;                   out->FaceElements->Nodes[INDEX2(7,k,NN)] =node0+                      1*Nstride0;
324                   if (useFullElementOrder){                   if (generateAllNodes){
325                      out->FaceElements->Nodes[INDEX2(8,k,NN)] =node0+           1*Nstride1+1*Nstride0;                      out->FaceElements->Nodes[INDEX2(8,k,NN)] =node0+           1*Nstride1+1*Nstride0;
326                   }                   }
327                }                }
# Line 363  Finley_Mesh* Finley_RectangularMesh_Hex2 Line 377  Finley_Mesh* Finley_RectangularMesh_Hex2
377                   out->FaceElements->Nodes[INDEX2(5,k,NN)] =node0+2*Nstride2+1*Nstride1+2*Nstride0;                   out->FaceElements->Nodes[INDEX2(5,k,NN)] =node0+2*Nstride2+1*Nstride1+2*Nstride0;
378                   out->FaceElements->Nodes[INDEX2(6,k,NN)] =node0+2*Nstride2+2*Nstride1+1*Nstride0;                   out->FaceElements->Nodes[INDEX2(6,k,NN)] =node0+2*Nstride2+2*Nstride1+1*Nstride0;
379                   out->FaceElements->Nodes[INDEX2(7,k,NN)] =node0+2*Nstride2+1*Nstride1           ;                   out->FaceElements->Nodes[INDEX2(7,k,NN)] =node0+2*Nstride2+1*Nstride1           ;
380                   if (useFullElementOrder){                   if (generateAllNodes){
381                   out->FaceElements->Nodes[INDEX2(8,k,NN)] =node0+2*Nstride2+1*Nstride1+1*Nstride0;                   out->FaceElements->Nodes[INDEX2(8,k,NN)] =node0+2*Nstride2+1*Nstride1+1*Nstride0;
382                   }                   }
383                }                }
# Line 373  Finley_Mesh* Finley_RectangularMesh_Hex2 Line 387  Finley_Mesh* Finley_RectangularMesh_Hex2
387         }         }
388         totalNECount+=NE1*NE0;         totalNECount+=NE1*NE0;
389       }       }
390       if (!periodic[0]) {       if (!periodic[0] && (local_NE0>0)) {
391          /* **  elements on boundary 001 (x1=0): */          /* **  elements on boundary 001 (x1=0): */
392            
393          if (e_offset0 == 0) {          if (e_offset0 == 0) {
# Line 422  Finley_Mesh* Finley_RectangularMesh_Hex2 Line 436  Finley_Mesh* Finley_RectangularMesh_Hex2
436                    out->FaceElements->Nodes[INDEX2(5,k,NN)] =node0+2*Nstride2+1*Nstride1           ;                    out->FaceElements->Nodes[INDEX2(5,k,NN)] =node0+2*Nstride2+1*Nstride1           ;
437                    out->FaceElements->Nodes[INDEX2(6,k,NN)] =node0+1*Nstride2+2*Nstride1           ;                    out->FaceElements->Nodes[INDEX2(6,k,NN)] =node0+1*Nstride2+2*Nstride1           ;
438                    out->FaceElements->Nodes[INDEX2(7,k,NN)] =node0+           1*Nstride1           ;                    out->FaceElements->Nodes[INDEX2(7,k,NN)] =node0+           1*Nstride1           ;
439                   if (useFullElementOrder){                   if (generateAllNodes){
440                      out->FaceElements->Nodes[INDEX2(8,k,NN)] =node0+1*Nstride2+1*Nstride1           ;                      out->FaceElements->Nodes[INDEX2(8,k,NN)] =node0+1*Nstride2+1*Nstride1           ;
441                   }                   }
442                 }                 }
# Line 478  Finley_Mesh* Finley_RectangularMesh_Hex2 Line 492  Finley_Mesh* Finley_RectangularMesh_Hex2
492                    out->FaceElements->Nodes[INDEX2(5,k,NN)]=node0+1*Nstride2+2*Nstride1+2*Nstride0;                    out->FaceElements->Nodes[INDEX2(5,k,NN)]=node0+1*Nstride2+2*Nstride1+2*Nstride0;
493                    out->FaceElements->Nodes[INDEX2(6,k,NN)]=node0+2*Nstride2+1*Nstride1+2*Nstride0;                    out->FaceElements->Nodes[INDEX2(6,k,NN)]=node0+2*Nstride2+1*Nstride1+2*Nstride0;
494                    out->FaceElements->Nodes[INDEX2(7,k,NN)]=node0+1*Nstride2           +2*Nstride0;                    out->FaceElements->Nodes[INDEX2(7,k,NN)]=node0+1*Nstride2           +2*Nstride0;
495                   if (useFullElementOrder){                   if (generateAllNodes){
496                      out->FaceElements->Nodes[INDEX2(8,k,NN)] =node0+1*Nstride2+1*Nstride1+2*Nstride0;                      out->FaceElements->Nodes[INDEX2(8,k,NN)] =node0+1*Nstride2+1*Nstride1+2*Nstride0;
497                   }                   }
498                 }                 }
# Line 489  Finley_Mesh* Finley_RectangularMesh_Hex2 Line 503  Finley_Mesh* Finley_RectangularMesh_Hex2
503           }           }
504           totalNECount+=NE1*NE2;           totalNECount+=NE1*NE2;
505       }       }
506       if (!periodic[1]) {       if (!periodic[1] && (local_NE1>0)) {
507          /* **  elements on boundary 010 (x2=0): */          /* **  elements on boundary 010 (x2=0): */
508          if (e_offset1 == 0) {          if (e_offset1 == 0) {
509             #pragma omp parallel for private(i0,i2,k,node0)             #pragma omp parallel for private(i0,i2,k,node0)
# Line 536  Finley_Mesh* Finley_RectangularMesh_Hex2 Line 550  Finley_Mesh* Finley_RectangularMesh_Hex2
550                    out->FaceElements->Nodes[INDEX2(5,k,NN)]=node0+1*Nstride2+           2*Nstride0;                    out->FaceElements->Nodes[INDEX2(5,k,NN)]=node0+1*Nstride2+           2*Nstride0;
551                    out->FaceElements->Nodes[INDEX2(6,k,NN)]=node0+2*Nstride2+           1*Nstride0;                    out->FaceElements->Nodes[INDEX2(6,k,NN)]=node0+2*Nstride2+           1*Nstride0;
552                    out->FaceElements->Nodes[INDEX2(7,k,NN)]=node0+1*Nstride2                      ;                    out->FaceElements->Nodes[INDEX2(7,k,NN)]=node0+1*Nstride2                      ;
553                   if (useFullElementOrder){                   if (generateAllNodes){
554                      out->FaceElements->Nodes[INDEX2(8,k,NN)] =node0+1*Nstride2+         1*Nstride0;                      out->FaceElements->Nodes[INDEX2(8,k,NN)] =node0+1*Nstride2+         1*Nstride0;
555                   }                   }
556                 }                 }
# Line 591  Finley_Mesh* Finley_RectangularMesh_Hex2 Line 605  Finley_Mesh* Finley_RectangularMesh_Hex2
605                    out->FaceElements->Nodes[INDEX2(5,k,NN)]=node0+2*Nstride2+2*Nstride1+1*Nstride0;                    out->FaceElements->Nodes[INDEX2(5,k,NN)]=node0+2*Nstride2+2*Nstride1+1*Nstride0;
606                    out->FaceElements->Nodes[INDEX2(6,k,NN)]=node0+1*Nstride2+2*Nstride1+2*Nstride0;                    out->FaceElements->Nodes[INDEX2(6,k,NN)]=node0+1*Nstride2+2*Nstride1+2*Nstride0;
607                    out->FaceElements->Nodes[INDEX2(7,k,NN)]=node0+           2*Nstride1+1*Nstride0;                    out->FaceElements->Nodes[INDEX2(7,k,NN)]=node0+           2*Nstride1+1*Nstride0;
608                   if (useFullElementOrder){                   if (generateAllNodes){
609                      out->FaceElements->Nodes[INDEX2(8,k,NN)] =node0+1*Nstride2+2*Nstride1+1*Nstride0;                      out->FaceElements->Nodes[INDEX2(8,k,NN)] =node0+1*Nstride2+2*Nstride1+1*Nstride0;
610                   }                   }
611                 }                 }

Legend:
Removed from v.1387  
changed lines
  Added in v.2722

  ViewVC Help
Powered by ViewVC 1.1.26