Changeset 06e153

Timestamp:

May 10, 2012, 2:17:09 PM (14 years ago)

Author:

Julian Iseringhausen <isering@…>

Children:

fcf7f6

Parents:

cd0fed

Message:

vmg: Improved serial performance.

git-svn-id: ​https://svn.version.fz-juelich.de/scafacos/trunk@1805 5161e1c8-67bf-11de-9fd5-51895aff932f

Location:

src/units/particle

Files:

• bspline.hpp (modified) (1 diff)
• interface_particles.cpp (modified) (1 diff)

src/units/particle/bspline.hpp

@@ -48 +48 @@
     int c = 0;
 
-    Index i;
-    Vector dir;
+    const int index_global_x = (p.Pos().X() - grid.Extent().Begin().X()) / grid.Extent().MeshWidth().X();
+    const int index_global_y = (p.Pos().Y() - grid.Extent().Begin().Y()) / grid.Extent().MeshWidth().Y();
+    const int index_global_z = (p.Pos().Z() - grid.Extent().Begin().Z()) / grid.Extent().MeshWidth().Z();
 
-    const Index index_global = static_cast<Index>((p.Pos() - grid.Extent().Begin()) / grid.Extent().MeshWidth());
-    assert(index_global.X() >= grid.Global().LocalBegin().X() && index_global.X() < grid.Global().LocalEnd().X());
-    assert(index_global.Y() >= grid.Global().LocalBegin().Y() && index_global.Y() < grid.Global().LocalEnd().Y());
-    assert(index_global.Z() >= grid.Global().LocalBegin().Z() && index_global.Z() < grid.Global().LocalEnd().Z());
+    assert(index_global_x >= grid.Global().LocalBegin().X() && index_global_x < grid.Global().LocalEnd().X());
+    assert(index_global_y >= grid.Global().LocalBegin().Y() && index_global_y < grid.Global().LocalEnd().Y());
+    assert(index_global_z >= grid.Global().LocalBegin().Z() && index_global_z < grid.Global().LocalEnd().Z());
 
-    const Index index_local = index_global - grid.Global().LocalBegin() + grid.Local().Begin();
-    assert(index_local.X() >= grid.Local().Begin().X() && index_local.X() < grid.Local().End().X());
-    assert(index_local.Y() >= grid.Local().Begin().Y() && index_local.Y() < grid.Local().End().Y());
-    assert(index_local.Z() >= grid.Local().Begin().Z() && index_local.Z() < grid.Local().End().Z());
+    const int index_local_x = index_global_x - grid.Global().LocalBegin().X() + grid.Local().Begin().X();
+    const int index_local_y = index_global_y - grid.Global().LocalBegin().Y() + grid.Local().Begin().Y();
+    const int index_local_z = index_global_z - grid.Global().LocalBegin().Z() + grid.Local().Begin().Z();
 
-    const Vector pos_beg = p.Pos() - grid.Extent().Begin() - grid.Extent().MeshWidth() * (index_global - near_field_cells);
-    const Vector& h = grid.Extent().MeshWidth();
+    assert(index_local_x >= grid.Local().Begin().X() && index_local_x < grid.Local().End().X());
+    assert(index_local_y >= grid.Local().Begin().Y() && index_local_y < grid.Local().End().Y());
+    assert(index_local_z >= grid.Local().Begin().Z() && index_local_z < grid.Local().End().Z());
+
+    const vmg_float pos_beg_x = p.Pos().X() - grid.Extent().Begin().X() - grid.Extent().MeshWidth().X() * (index_global_x - near_field_cells);
+    const vmg_float pos_beg_y = p.Pos().Y() - grid.Extent().Begin().Y() - grid.Extent().MeshWidth().Y() * (index_global_y - near_field_cells);
+    const vmg_float pos_beg_z = p.Pos().Z() - grid.Extent().Begin().Z() - grid.Extent().MeshWidth().Z() * (index_global_z - near_field_cells);
+
+    const vmg_float h_x = grid.Extent().MeshWidth().X();
+    const vmg_float h_y = grid.Extent().MeshWidth().Y();
+    const vmg_float h_z = grid.Extent().MeshWidth().Z();
 
     // Iterate over all grid points which lie in the support of the interpolating B-Spline
-    dir[0] = pos_beg[0];
-    for (i.X()=-1*near_field_cells; i.X()<=near_field_cells; ++i.X()) {
-      dir[1] = pos_beg[1];
-      for (i.Y()=-1*near_field_cells; i.Y()<=near_field_cells; ++i.Y()) {
-        dir[2] = pos_beg[2];
-        for (i.Z()=-1*near_field_cells; i.Z()<=near_field_cells; ++i.Z()) {
+
+    vmg_float dir_x, dir_y, dir_z;
+    dir_x = pos_beg_x;
+    for (int i=-1*near_field_cells; i<=near_field_cells; ++i) {
+      dir_y = pos_beg_y;
+      for (int j=-1*near_field_cells; j<=near_field_cells; ++j) {
+        dir_z = pos_beg_z;
+        for (int k=-1*near_field_cells; k<=near_field_cells; ++k) {
 
           // Compute distance from grid point to particle
-          temp_val = EvaluateSpline(dir.Length());
+          temp_val = EvaluateSpline(std::sqrt(dir_x*dir_x+dir_y*dir_y+dir_z*dir_z));
           vals[c++] = temp_val * p.Charge();
           int_val += temp_val;
 
-          dir[2] -= h[2];
+          dir_z -= h_z;
         }
-        dir[1] -= h[1];
+        dir_y -= h_y;
       }
-      dir[0] -= h[0];
+      dir_x -= h_x;
     }
 
     // Reciprocal value of the numerically integrated spline
-    int_val = 1.0 / (int_val * grid.Extent().MeshWidth().Product());
+    int_val = 1.0 / (int_val * h_x * h_y * h_z);
 
     c = 0;
-    for (i.X()=-1*near_field_cells; i.X()<=near_field_cells; ++i.X())
-      for (i.Y()=-1*near_field_cells; i.Y()<=near_field_cells; ++i.Y())
-        for (i.Z()=-1*near_field_cells; i.Z()<=near_field_cells; ++i.Z())
-          grid(index_local+i) += vals[c++] * int_val;
+    for (int i=-1*near_field_cells; i<=near_field_cells; ++i)
+      for (int j=-1*near_field_cells; j<=near_field_cells; ++j)
+        for (int k=-1*near_field_cells; k<=near_field_cells; ++k)
+          grid(index_local_x + i,
+               index_local_y + j,
+               index_local_z + k) += vals[c++] * int_val;
 
   }

src/units/particle/interface_particles.cpp

@@ -84 +84 @@
 
   // Assign charge values to the right hand side
-  for (index.X()=0; index.X()<grid.Local().Size().X(); ++index.X())
-    for (index.Y()=0; index.Y()<grid.Local().Size().Y(); ++index.Y())
-      for (index.Z()=0; index.Z()<grid.Local().Size().Z(); ++index.Z())
-        grid(index + grid.Local().Begin()) = 4.0 * Math::pi * particle_grid.GetVal(index + particle_grid.Local().Begin());
+  for (int i=0; i<grid.Local().Size().X(); ++i)
+    for (int j=0; j<grid.Local().Size().Y(); ++j)
+      for (int k=0; k<grid.Local().Size().Z(); ++k)
+        grid(grid.Local().Begin().X() + i,
+             grid.Local().Begin().Y() + j,
+             grid.Local().Begin().Z() + k) = 4.0 * Math::pi *
+          particle_grid.GetVal(particle_grid.Local().Begin().X() + i,
+                               particle_grid.Local().Begin().Y() + j,
+                               particle_grid.Local().Begin().Z() + k);
 
 #ifdef DEBUG_OUTPUT