openMP modif

srcs/BEM/mainBEM.cpp

@@ -11,7 +11,10 @@ int main(int argc, char **argv)
         return 0;
     }
 
+    double start = omp_get_wtime();
+    int nthreads = omp_get_max_threads();
     gmsh::initialize();
+    omp_set_num_threads(nthreads);
 
     gmsh::open(argv[1]);
     gmsh::model::mesh::generate(2);

srcs/BEM/postProcessing.cpp

@@ -377,6 +377,7 @@ void computeElementData(std::vector<std::vector<std::vector<std::vector<double>>
             phi[i][j].resize(domainElementTags[i][j].size());
             electricField[i][j].resize(domainElementTags[i][j].size());
             gmsh::model::mesh::getBarycenters(domainElementTypes[i][j], domainEntityTags[i], false, false, barycenters);
+            #pragma omp parallel for schedule(static)
             for(size_t k = 0; k < domainElementTags[i][j].size(); k++)
             {
                 double x = barycenters[3*k];
@@ -386,55 +387,44 @@ void computeElementData(std::vector<std::vector<std::vector<std::vector<double>>
                 double gradPhiX = 0;
                 double gradPhiY = 0;
 
-                #pragma omp parallel shared(currentPhi, gradPhiX, gradPhiY)
+                double localCurrentPhi = 0;
+                double localGradPhiX = 0;
+                double localGradPhiY = 0;
+
+                double h, hX, hY;
+                double g, gX, gY;
+                
+                for(std::size_t el = 0; el < elementVector.size(); el++)
                 {
-                    double localCurrentPhi = 0;
-                    double localGradPhiX = 0;
-                    double localGradPhiY = 0;
-
-                    double h, hX, hY;
-                    double g, gX, gY;
-                    
-                    #pragma omp for schedule(guided) 
-                    for(std::size_t el = 0; el < elementVector.size(); el++)
+                    if(computePhi)
                     {
-                        if(computePhi)
-                        {
-                            h = computeAnalyticalH(x, y, elementVector[el]);
-                            if(!solutionType)
-                                g = computeNumericalG(x, y, elementVector[el], integrationType);
-                            else
-                                g = computeAnalyticalG(x, y, elementVector[el]);
+                        h = computeAnalyticalH(x, y, elementVector[el]);
+                        if(!solutionType)
+                            g = computeNumericalG(x, y, elementVector[el], integrationType);
+                        else
+                            g = computeAnalyticalG(x, y, elementVector[el]);
 
-                            localCurrentPhi += h*elementVector[el].bcValue[0] - g*elementVector[el].bcValue[1];
+                        localCurrentPhi += h*elementVector[el].bcValue[0] - g*elementVector[el].bcValue[1];
+                    }
+                    if(computeElectricField)
+                    {
+                        if(!solutionType)
+                        {
+                            hX = computeNumericalGradHx(x, y, elementVector[el], integrationType);
+                            gX = computeNumericalGradGx(x, y, elementVector[el], integrationType);
+                            hY = computeNumericalGradHy(x, y, elementVector[el], integrationType);
+                            gY = computeNumericalGradGy(x, y, elementVector[el], integrationType);
                         }
-                        if(computeElectricField)
+                        else
                         {
-                            if(!solutionType)
-                            {
-                                hX = computeNumericalGradHx(x, y, elementVector[el], integrationType);
-                                gX = computeNumericalGradGx(x, y, elementVector[el], integrationType);
-                                hY = computeNumericalGradHy(x, y, elementVector[el], integrationType);
-                                gY = computeNumericalGradGy(x, y, elementVector[el], integrationType);
-                            }
-                            else
-                            {
-                                hX = computeAnalyticalGradHx(x, y, elementVector[el]);
-                                gX = computeAnalyticalGradGx(x, y, elementVector[el]);
-                                hY = computeAnalyticalGradHy(x, y, elementVector[el]);
-                                gY = computeAnalyticalGradGy(x, y, elementVector[el]);
-                            }
-
-                            localGradPhiX += hX*elementVector[el].bcValue[0] - gX*elementVector[el].bcValue[1];
-                            localGradPhiY += hY*elementVector[el].bcValue[0] - gY*elementVector[el].bcValue[1];
+                            hX = computeAnalyticalGradHx(x, y, elementVector[el]);
+                            gX = computeAnalyticalGradGx(x, y, elementVector[el]);
+                            hY = computeAnalyticalGradHy(x, y, elementVector[el]);
+                            gY = computeAnalyticalGradGy(x, y, elementVector[el]);
                         }
-                    }
 
-                    #pragma omp critical
-                    {
-                        currentPhi += localCurrentPhi;
-                        gradPhiX += localGradPhiX;
-                        gradPhiY += localGradPhiY;
+                        localGradPhiX += hX*elementVector[el].bcValue[0] - gX*elementVector[el].bcValue[1];
+                        localGradPhiY += hY*elementVector[el].bcValue[0] - gY*elementVector[el].bcValue[1];
                     }
                 }