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

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

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

revision 777 by gross, Wed Jul 12 08:54:45 2006 UTC revision 853 by gross, Wed Sep 20 05:56:36 2006 UTC
# Line 34  void Finley_Assemble_integrate(Finley_No Line 34  void Finley_Assemble_integrate(Finley_No
34      type_t datacase;      type_t datacase;
35      index_t node_offset;      index_t node_offset;
36      bool_t reducedIntegrationOrder;      bool_t reducedIntegrationOrder;
     dim_t q,e,i;  
37      if (nodes==NULL || elements==NULL) return;      if (nodes==NULL || elements==NULL) return;
38      type_t data_type=getFunctionSpaceType(data);      type_t data_type=getFunctionSpaceType(data);
39      dim_t numComps=getDataPointSize(data);      dim_t numComps=getDataPointSize(data);
# Line 68  void Finley_Assemble_integrate(Finley_No Line 67  void Finley_Assemble_integrate(Finley_No
67          /* now we can start */          /* now we can start */
68    
69          if (Finley_noError()) {          if (Finley_noError()) {
70                dim_t q,e,i;
71                double *out_local=NULL, rtmp,*data_array=NULL;
72              for (q=0;q<numComps;q++) out[q]=0;              for (q=0;q<numComps;q++) out[q]=0;
73              #pragma omp parallel private(q,i)              #pragma omp parallel private(q,i,rtmp,data_array,out_local)
74              {              {
75                  double out_local[numComps], rtmp;                  out_local=THREAD_MEMALLOC(numComps,double);
76                  double* data_array;                  if (! Finley_checkPtr(out_local) ) {
77                       /* initialize local result */
78                  /* initialize local result */  
79                   for (i=0;i<numComps;i++) out_local[i]=0;
80              for (i=0;i<numComps;i++) out_local[i]=0;  
81                       /* open the element loop */
82                  /* open the element loop */    
83                       if (isExpanded(data)) {
84                  if (isExpanded(data)) {                         #pragma omp for private(e) schedule(static)
85                      #pragma omp for private(e) schedule(static)                         for(e=0;e<elements->numElements;e++) {
86                      for(e=0;e<elements->numElements;e++) {                              data_array=getSampleData(data,e);
87                           data_array=getSampleData(data,e);                              for (q=0;q<jac->ReferenceElement->numQuadNodes;q++) {
88                           for (q=0;q<jac->ReferenceElement->numQuadNodes;q++) {                                    for (i=0;i<numComps;i++) out_local[i]+=data_array[INDEX2(i,q,numComps)]*jac->volume[INDEX2(q,e,jac->ReferenceElement->numQuadNodes)];
89                                 for (i=0;i<numComps;i++) out_local[i]+=data_array[INDEX2(i,q,numComps)]*jac->volume[INDEX2(q,e,jac->ReferenceElement->numQuadNodes)];                              }
90                           }                         }
91                      }                     } else {
92                  } else {                        #pragma omp for private(e) schedule(static)
93                     #pragma omp for private(e) schedule(static)                        for(e=0;e<elements->numElements;e++) {
94                     for(e=0;e<elements->numElements;e++) {                             data_array=getSampleData(data,e);
95                          data_array=getSampleData(data,e);                             rtmp=0.;
96                          rtmp=0.;                             for (q=0;q<jac->ReferenceElement->numQuadNodes;q++) rtmp+=jac->volume[INDEX2(q,e,jac->ReferenceElement->numQuadNodes)];
97                          for (q=0;q<jac->ReferenceElement->numQuadNodes;q++) rtmp+=jac->volume[INDEX2(q,e,jac->ReferenceElement->numQuadNodes)];                             for (i=0;i<numComps;i++) out_local[i]+=data_array[i]*rtmp;
98                          for (i=0;i<numComps;i++) out_local[i]+=data_array[i]*rtmp;                        }
99                     }                     }
100                       /* add local results to global result */
101                       #pragma omp critical
102                       for (i=0;i<numComps;i++) out[i]+=out_local[i];
103                  }                  }
104                  /* add local results to global result */                  THREAD_MEMFREE(out_local);
                 #pragma omp critical  
                 for (i=0;i<numComps;i++) out[i]+=out_local[i];  
   
105              }              }
106         }         }
107     }     }

Legend:
Removed from v.777  
changed lines
  Added in v.853

  ViewVC Help
Powered by ViewVC 1.1.26