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

Diff of /trunk/finley/src/Mesh_rec8.c

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

revision 782 by bcumming, Tue Jul 18 00:47:47 2006 UTC revision 787 by bcumming, Wed Jul 26 01:46:45 2006 UTC
# Line 208  Finley_Mesh* Finley_RectangularMesh_Rec8 Line 208  Finley_Mesh* Finley_RectangularMesh_Rec8
208        
209    /*  allocate tables: */    /*  allocate tables: */
210    Finley_NodeFile_allocTable(out->Nodes,N0*N1);    Finley_NodeFile_allocTable(out->Nodes,N0*N1);
211      Finley_ElementFile_allocTable(out->Elements,NE0*NE1);
212      Finley_ElementFile_allocTable(out->FaceElements,NFaceElements);
213  #ifdef PASO_MPI  #ifdef PASO_MPI
214    Finley_NodeDistribution_allocTable( out->Nodes->degreeOfFreedomDistribution, NDOF0*NDOF1, 0, 0 );    Finley_NodeDistribution_allocTable( out->Nodes->degreeOfFreedomDistribution, NDOF0*NDOF1, 0, 0 );
215      Finley_ElementDistribution_allocTable( out->Elements->elementDistribution, NE0*NE1, NE0*NE1);
216      Finley_ElementDistribution_allocTable( out->FaceElements->elementDistribution, NFaceElements, NFaceElements );
217  #endif  #endif
   Finley_ElementFile_allocTable(out->Elements,NE0*NE1);  
   Finley_ElementFile_allocTable(out->FaceElements,NFaceElements);  
218    if (! Finley_noError()) {    if (! Finley_noError()) {
219        Finley_Mesh_dealloc(out);        Finley_Mesh_dealloc(out);
220        return NULL;        return NULL;
# Line 433  Finley_Mesh* Finley_RectangularMesh_Rec8 Line 435  Finley_Mesh* Finley_RectangularMesh_Rec8
435    /*   condense the nodes: */    /*   condense the nodes: */
436    Finley_Mesh_resolveNodeIds(out);    Finley_Mesh_resolveNodeIds(out);
437            
438      if (! Finley_noError()) {
439          Finley_Mesh_dealloc(out);
440          return NULL;
441      }
442    /* prepare mesh for further calculatuions:*/    /* prepare mesh for further calculatuions:*/
443    Finley_Mesh_prepare(out) ;    Finley_Mesh_prepare(out) ;
444    
# Line 460  Finley_Mesh* Finley_RectangularMesh_Rec8 Line 466  Finley_Mesh* Finley_RectangularMesh_Rec8
466    double time0=Finley_timer();    double time0=Finley_timer();
467    bool_t periodicLocal[2], domLeft=FALSE, domRight=FALSE, domInternal=FALSE, boundaryLeft=FALSE, boundaryRight=FALSE;    bool_t periodicLocal[2], domLeft=FALSE, domRight=FALSE, domInternal=FALSE, boundaryLeft=FALSE, boundaryRight=FALSE;
468    Paso_MPIInfo *mpi_info = NULL;    Paso_MPIInfo *mpi_info = NULL;
469        index_t faceStart, numFaceLeft;
470    
471      /* these are used in constructing the face elements, and give the permutation      /* these are used in constructing the face elements, and give the permutation
472         of the node order of a normal element for each face */         of the node order of a normal element for each face */
# Line 546  Finley_Mesh* Finley_RectangularMesh_Rec8 Line 553  Finley_Mesh* Finley_RectangularMesh_Rec8
553        
554    /*  allocate tables: */    /*  allocate tables: */
555    Finley_NodeFile_allocTable(out->Nodes,N0t*N1);    Finley_NodeFile_allocTable(out->Nodes,N0t*N1);
   Finley_NodeDistribution_allocTable( out->Nodes->degreeOfFreedomDistribution, NDOF0t*NDOF1, NDOF1*3, 0 );  
556    Finley_ElementFile_allocTable(out->Elements,numElementsLocal*NE1);    Finley_ElementFile_allocTable(out->Elements,numElementsLocal*NE1);
557    Finley_ElementFile_allocTable(out->FaceElements,NFaceElements);    Finley_ElementFile_allocTable(out->FaceElements,NFaceElements);
558        
559      Finley_NodeDistribution_allocTable( out->Nodes->degreeOfFreedomDistribution, NDOF0t*NDOF1, NDOF1*3, 0 );
560      Finley_ElementDistribution_allocTable( out->Elements->elementDistribution, numElementsLocal*NE1, NE1*(numElementsLocal-boundaryRight*(!periodic[1])) );
561      Finley_ElementDistribution_allocTable( out->FaceElements->elementDistribution, NFaceElements, NFaceElements-2*boundaryRight*(!periodic[1]) );
562    if (! Finley_noError()) {    if (! Finley_noError()) {
563        Finley_Mesh_dealloc(out);        Finley_Mesh_dealloc(out);
564        return NULL;        return NULL;
# Line 556  Finley_Mesh* Finley_RectangularMesh_Rec8 Line 566  Finley_Mesh* Finley_RectangularMesh_Rec8
566    
567    /*  set nodes: */    /*  set nodes: */
568    #pragma omp parallel for private(i0,i1,k)    #pragma omp parallel for private(i0,i1,k)
569    for (k=0,i1=0;i1<N1;i1++) {    for (i1=0;i1<N1;i1++) {
570      for (i0=0;i0<N0t;i0++, k++) {      for (i0=0;i0<N0t;i0++, k++) {
571                k=i0+i1*N0t;
572        out->Nodes->Coordinates[INDEX2(0,k,2)]=DBLE((firstNodeConstruct + i0) % N0)/DBLE(N0-1)*Length[0];        out->Nodes->Coordinates[INDEX2(0,k,2)]=DBLE((firstNodeConstruct + i0) % N0)/DBLE(N0-1)*Length[0];
573        out->Nodes->Coordinates[INDEX2(1,k,2)]=DBLE(i1)/DBLE(N1-1)*Length[1];        out->Nodes->Coordinates[INDEX2(1,k,2)]=DBLE(i1)/DBLE(N1-1)*Length[1];
574        out->Nodes->Id[k]=i0 + i1*N0t;        out->Nodes->Id[k]=i0 + i1*N0t;
# Line 661  Finley_Mesh* Finley_RectangularMesh_Rec8 Line 672  Finley_Mesh* Finley_RectangularMesh_Rec8
672        k=i0+numElementsLocal*i1;        k=i0+numElementsLocal*i1;
673              node0 = (periodicLocal[0] && !i0) ? 2*i1*N0t :  2*i0+2*i1*N0t+periodicLocal[0] ;              node0 = (periodicLocal[0] && !i0) ? 2*i1*N0t :  2*i0+2*i1*N0t+periodicLocal[0] ;
674    
675        out->Elements->Id[k]=k;        out->Elements->Id[k]=((firstNodeConstruct/2+i0)%NE0) * NE1 + i1;
676        out->Elements->Tag[k]=0;        out->Elements->Tag[k]=0;
677        out->Elements->Color[k]=COLOR_MOD(i0)+3*COLOR_MOD(i1);        out->Elements->Color[k]=COLOR_MOD(i0)+3*COLOR_MOD(i1);
678              out->Elements->Dom[k] = ELEMENT_INTERNAL;              out->Elements->Dom[k] = ELEMENT_INTERNAL;
# Line 693  Finley_Mesh* Finley_RectangularMesh_Rec8 Line 704  Finley_Mesh* Finley_RectangularMesh_Rec8
704       NUMNODES=3;       NUMNODES=3;
705    }    }
706        
707        Finley_ElementFile_setDomainFlags( out->FaceElements );
708    totalNECount=numElementsLocal*NE1;    totalNECount=numElementsLocal*NE1;
709    faceNECount=0;    faceNECount=0;
710      faceNEBoundary = totalNECount + numElementsInternal*(!periodic[1])*2 + NE1*(!periodic[0])*(domLeft+domRight);      faceNEBoundary = totalNECount + numElementsInternal*(!periodic[1])*2 + NE1*(!periodic[0])*(domLeft+domRight);
711        
712        numFaceLeft = domLeft*(!periodic[0])*NE1;
713        faceStart = out->FaceElements->elementDistribution->vtxdist[mpi_info->rank];
714    if (!periodic[1]) {    if (!periodic[1]) {
715       /* **  elements on boundary 010 (x2=0): */       /* **  elements on boundary 010 (x1=0): */
716        
717       #pragma omp parallel for private(i0,k,node0)       #pragma omp parallel for private(i0,k,kk,facePerm)
718       for (i0=0;i0<numElementsLocal; i0++) {       for (i0=0;i0<numElementsLocal; i0++) {
719         k=i0+faceNECount;         k=i0+faceNECount;
720               kk=i0;               kk=i0;
721               facePerm = useElementsOnFace ? face0 : face0nc;               facePerm = useElementsOnFace ? face0 : face0nc;
722    
723         out->FaceElements->Id[k]=i0+totalNECount;         out->FaceElements->Id[k]=faceStart + numFaceLeft + i0*2;
724         out->FaceElements->Tag[k]=10;         out->FaceElements->Tag[k]=10;
725         out->FaceElements->Color[k]=i0%2;         out->FaceElements->Color[k]=i0%2;
726               out->FaceElements->Dom[k]=ELEMENT_INTERNAL;               out->FaceElements->Dom[k]=ELEMENT_INTERNAL;
# Line 719  Finley_Mesh* Finley_RectangularMesh_Rec8 Line 733  Finley_Mesh* Finley_RectangularMesh_Rec8
733               if( periodicLocal[0] )               if( periodicLocal[0] )
734                   out->FaceElements->Dom[faceNECount+1]=ELEMENT_BOUNDARY;                   out->FaceElements->Dom[faceNECount+1]=ELEMENT_BOUNDARY;
735           }           }
736           if( boundaryRight )           if( boundaryRight ){
737               out->FaceElements->Dom[faceNECount+numElementsLocal-1]=ELEMENT_BOUNDARY;               out->FaceElements->Dom[faceNECount+numElementsLocal-1]=ELEMENT_BOUNDARY;
738                 out->FaceElements->Id[faceNECount+numElementsLocal-1]=out->FaceElements->elementDistribution->vtxdist[mpi_info->rank+1];
739             }
740       totalNECount+=numElementsLocal;       totalNECount+=numElementsLocal;
741       faceNECount+=numElementsLocal;       faceNECount+=numElementsLocal;
742        
743       /* **  elements on boundary 020 (x2=1): */       /* **  elements on boundary 020 (x1=1): */
744        
745       #pragma omp parallel for private(i0,k,node0)       #pragma omp parallel for private(i0,k,kk,facePerm)
746       for (i0=0;i0<numElementsLocal;i0++) {       for (i0=0;i0<numElementsLocal;i0++) {
747         k=i0+faceNECount;         k=i0+faceNECount;
748               kk=i0 + numElementsLocal*(NE1-1);               kk=i0 + numElementsLocal*(NE1-1);
749               facePerm = useElementsOnFace ? face1 : face1nc;               facePerm = useElementsOnFace ? face1 : face1nc;
750    
751         out->FaceElements->Id[k]=i0+totalNECount;         out->FaceElements->Id[k]=faceStart + numFaceLeft + i0*2+1;
752         out->FaceElements->Tag[k]=20;         out->FaceElements->Tag[k]=20;
753         out->FaceElements->Color[k]=i0%2+2;         out->FaceElements->Color[k]=i0%2+2;
754               out->FaceElements->Dom[k]=ELEMENT_INTERNAL;               out->FaceElements->Dom[k]=ELEMENT_INTERNAL;
# Line 745  Finley_Mesh* Finley_RectangularMesh_Rec8 Line 761  Finley_Mesh* Finley_RectangularMesh_Rec8
761               if( periodicLocal[0] )               if( periodicLocal[0] )
762                   out->FaceElements->Dom[faceNECount+1]=ELEMENT_BOUNDARY;                   out->FaceElements->Dom[faceNECount+1]=ELEMENT_BOUNDARY;
763           }           }
764           if( boundaryRight )           if( boundaryRight ){
765               out->FaceElements->Dom[faceNECount+numElementsLocal-1]=ELEMENT_BOUNDARY;               out->FaceElements->Dom[faceNECount+numElementsLocal-1]=ELEMENT_BOUNDARY;
766                 out->FaceElements->Id[faceNECount+numElementsLocal-1]=out->FaceElements->elementDistribution->vtxdist[mpi_info->rank+1]+1;
767             }
768       totalNECount+=numElementsLocal;       totalNECount+=numElementsLocal;
769       faceNECount+=numElementsLocal;       faceNECount+=numElementsLocal;
770    }    }
771    if (domLeft && !periodicLocal[0]) {    if (domLeft && !periodicLocal[0]) {
772       /* **  elements on boundary 001 (x1=0): */       /* **  elements on boundary 001 (x0=0): */
773        
774       #pragma omp parallel for private(i1,k,node0)       #pragma omp parallel for private(i1,k,kk,facePerm)
775       for (i1=0;i1<NE1;i1++) {       for (i1=0;i1<NE1;i1++) {
776         k=i1+faceNECount;         k=i1+faceNECount;
777               kk=i1*numElementsLocal;               kk=i1*numElementsLocal;
778               facePerm = useElementsOnFace ? face2 : face2nc;               facePerm = useElementsOnFace ? face2 : face2nc;
779    
780         out->FaceElements->Id[k]=i1+totalNECount;         out->FaceElements->Id[k]=faceStart + i1;
781         out->FaceElements->Tag[k]=1;         out->FaceElements->Tag[k]=1;
782         out->FaceElements->Color[k]=(i1%2)+4;         out->FaceElements->Color[k]=(i1%2)+4;
783               out->FaceElements->Dom[k]=ELEMENT_INTERNAL;               out->FaceElements->Dom[k]=ELEMENT_INTERNAL;
# Line 770  Finley_Mesh* Finley_RectangularMesh_Rec8 Line 788  Finley_Mesh* Finley_RectangularMesh_Rec8
788       totalNECount+=NE1;       totalNECount+=NE1;
789       faceNECount+=NE1;       faceNECount+=NE1;
790        
791       /* **  elements on boundary 002 (x1=1): */       /* **  elements on boundary 002 (x0=1): */
792    }    }
793    if (domRight && !periodicLocal[1]) {    if (domRight && !periodicLocal[1]) {
794       #pragma omp parallel for private(i1,k,node0)       #pragma omp parallel for private(i1,k,kk,facePerm)
795       for (i1=0;i1<NE1;i1++) {       for (i1=0;i1<NE1;i1++) {
796         k=i1+faceNECount;         k=i1+faceNECount;
797               kk=(i1+1)*numElementsLocal-1;               kk=(i1+1)*numElementsLocal-1;
798               facePerm = useElementsOnFace ? face3 : face3nc;               facePerm = useElementsOnFace ? face3 : face3nc;
799        
800         out->FaceElements->Id[k]=i1+totalNECount;         out->FaceElements->Id[k]=faceStart + numElementsLocal*2*(!periodic[1])+i1;
801         out->FaceElements->Tag[k]=2;         out->FaceElements->Tag[k]=2;
802         out->FaceElements->Color[k]=(i1%2)+6;         out->FaceElements->Color[k]=(i1%2)+6;
803               out->FaceElements->Dom[k]=ELEMENT_INTERNAL;               out->FaceElements->Dom[k]=ELEMENT_INTERNAL;
# Line 794  Finley_Mesh* Finley_RectangularMesh_Rec8 Line 812  Finley_Mesh* Finley_RectangularMesh_Rec8
812    out->FaceElements->maxColor=7;    out->FaceElements->maxColor=7;
813      out->FaceElements->numElements = faceNECount;      out->FaceElements->numElements = faceNECount;
814    
     Finley_ElementFile_setDomainFlags( out->FaceElements );  
   
815    /* setup distribution info for other elements */    /* setup distribution info for other elements */
816      Finley_ElementFile_setDomainFlags( out->ContactElements );      Finley_ElementFile_setDomainFlags( out->ContactElements );
817      Finley_ElementFile_setDomainFlags( out->Points );      Finley_ElementFile_setDomainFlags( out->Points );
818    
819        Finley_Mesh_prepareElementDistribution( out );
820    
821      /* reorder the degrees of freedom */      /* reorder the degrees of freedom */
822      Finley_Mesh_resolveDegreeOfFreedomOrder( out, TRUE );      Finley_Mesh_resolveDegreeOfFreedomOrder( out, TRUE );
823            

Legend:
Removed from v.782  
changed lines
  Added in v.787

  ViewVC Help
Powered by ViewVC 1.1.26