MPI verified EVP 16 proc, 2562

Author: dhyan1272

src/pmpo_MPMesh.cpp

@@ -84,8 +84,8 @@ void MPMesh::calculateStress(){
 
       Vec3d strain_rate (MPsStrainRate(mp, 0), MPsStrainRate(mp, 1), MPsStrainRate(mp, 2));
       Vec3d stress(MPsStress(mp, 0), MPsStress(mp, 1), MPsStress(mp, 2));
-      //constitutive_evp(strain_rate, stress, MPsIcePressure(mp, 0), MPsRepPressure(mp, 0), MPsArea(mp, 0), elasticTimeStep, dampingTimescale);
-      constitutive_linear(strain_rate, stress);
+      constitutive_evp(strain_rate, stress, MPsIcePressure(mp, 0), MPsRepPressure(mp, 0), MPsArea(mp, 0), elasticTimeStep, dampingTimescale);
+      //constitutive_linear(strain_rate, stress);
       for (int m=0 ; m<3; m++)
         MPsStress(mp, m) = stress[m];
       /*
@@ -176,9 +176,9 @@ void MPMesh::calculateStressDivergence(){
   };
   p_MPs->parallel_for(stress_div, " stress_div_assembly");
 
-  //TO DO make to 2 fields instrad of 4 to reduce latency
   if(numProcsTot>1){
-    //Takes contribution of halo vertices and adds it in halo cells
+    
+    //Takes contribution of halo vertices and adds it in owner procs
     communicate_and_take_halo_contributions(stress_divU, numVertices, 1, 0, 0);
     communicate_and_take_halo_contributions(stress_divV, numVertices, 1, 0, 0);
     communicate_and_take_halo_contributions(divU_edge, numVertices, 1, 0, 0);
@@ -188,6 +188,7 @@ void MPMesh::calculateStressDivergence(){
     communicate_and_take_halo_contributions(stress_divV, numVertices, 1, 1, 1);
     communicate_and_take_halo_contributions(divU_edge, numVertices, 1, 1, 1);
     communicate_and_take_halo_contributions(divV_edge, numVertices, 1, 1, 1);
+    
   }
 
   auto stressDivergence = p_mesh->getMeshField<MeshF_StressDivergence>();

src/pmpo_c.cpp

@@ -620,6 +620,54 @@ void polympo_calculateMPStress_f(MPMesh_ptr p_mpmesh){
   mpMesh->calculateStress(); 
 }
 
+void polympo_setMPStress_f(MPMesh_ptr p_mpmesh, const int nComps, const int numMPs, const double* mpStressIn) {
+  Kokkos::Timer timer;
+  checkMPMeshValid(p_mpmesh);
+  auto p_MPs = ((polyMPO::MPMesh*)p_mpmesh)->p_MPs;
+  PMT_ALWAYS_ASSERT(nComps == vec3d_nEntries);
+  PMT_ALWAYS_ASSERT(numMPs >= p_MPs->getCount());
+  //PMT_ALWAYS_ASSERT(numMPs >= p_MPs->getMaxAppID());
+
+  auto mpStress = p_MPs->getData<polyMPO::MPF_Stress>();
+  auto mpAppID = p_MPs->getData<polyMPO::MPF_MP_APP_ID>();
+  kkViewHostU<const double**> mpStressIn_h(mpStressIn, nComps, numMPs);
+  Kokkos::View<double**> mpStressIn_d("mpStressDevice", vec3d_nEntries, numMPs);
+  Kokkos::deep_copy(mpStressIn_d, mpStressIn_h);
+  auto setMPStress = PS_LAMBDA(const int& elm, const int& mp, const int& mask){
+    if(mask){
+      mpStress(mp,0) = mpStressIn_d(0, mpAppID(mp));
+      mpStress(mp,1) = mpStressIn_d(1, mpAppID(mp));
+      mpStress(mp,2) = mpStressIn_d(2, mpAppID(mp));
+    }
+  };
+  p_MPs->parallel_for(setMPStress, "setMPStress");
+  pumipic::RecordTime("PolyMPO_setMPStress", timer.seconds());
+}
+
+void polympo_getMPStress_f(MPMesh_ptr p_mpmesh, const int nComps, const int numMPs, double* mpStressHost) {
+  Kokkos::Timer timer;
+  checkMPMeshValid(p_mpmesh);
+  auto p_MPs = ((polyMPO::MPMesh*)p_mpmesh)->p_MPs;
+  PMT_ALWAYS_ASSERT(nComps == vec3d_nEntries);
+  PMT_ALWAYS_ASSERT(numMPs >= p_MPs->getCount());
+  //PMT_ALWAYS_ASSERT(numMPs >= p_MPs->getMaxAppID());
+
+  auto mpStress = p_MPs->getData<polyMPO::MPF_Stress>();
+  auto mpAppID = p_MPs->getData<polyMPO::MPF_MP_APP_ID>();
+  Kokkos::View<double**> mpStressCopy("mpStressCopy", vec3d_nEntries, numMPs);
+  auto getMPStress = PS_LAMBDA(const int& elm, const int& mp, const int& mask){
+    if(mask){
+      mpStressCopy(0,mpAppID(mp)) = mpStress(mp,0);
+      mpStressCopy(1,mpAppID(mp)) = mpStress(mp,1);
+      mpStressCopy(2,mpAppID(mp)) = mpStress(mp,2);
+    }
+  };
+  p_MPs->parallel_for(getMPStress, "getMPStress");
+  kkDbl2dViewHostU arrayHost(mpStressHost, nComps, numMPs);
+  Kokkos::deep_copy(arrayHost, mpStressCopy);
+  pumipic::RecordTime("PolyMPO_getMPStress", timer.seconds());
+}
+
 void polympo_setAreaMP_f(MPMesh_ptr p_mpmesh, const int nComps, const int numMPs, double* areaMPHost){
   Kokkos::Timer timer;
   checkMPMeshValid(p_mpmesh);

src/pmpo_c.h

@@ -48,6 +48,8 @@ void polympo_setMPVel_f(MPMesh_ptr p_mpmesh, const int nComps, const int numMPs,
 void polympo_getMPVel_f(MPMesh_ptr p_mpmesh, const int nComps, const int numMPs, double* mpVelHost);
 void polympo_calculateMPStrainRate_f(MPMesh_ptr p_mpmesh);
 void polympo_calculateMPStress_f(MPMesh_ptr p_mpmesh);
+void polympo_setMPStress_f(MPMesh_ptr p_mpmesh, const int nComps, const int numMPs, const double* mpStressIn);
+void polympo_getMPStress_f(MPMesh_ptr p_mpmesh, const int nComps, const int numMPs, double* mpStressHost);
 void polympo_setAreaMP_f(MPMesh_ptr p_mpmesh, const int nComps, const int numMPs, double* areaMPHost);
 void polympo_setIcePressureMP_f(MPMesh_ptr p_mpmesh, const int nComps, const int numMPs, double* icePressureMPHost);
 void polympo_getReplacementPressureMP_f(MPMesh_ptr p_mpmesh, const int nComps, const int numMPs, double* replacementPressureMPHost);

src/pmpo_fortran.f90

@@ -377,6 +377,22 @@ subroutine polympo_calculateMPStress(mpMesh) &
     type(c_ptr), value :: mpMesh
   end subroutine
 
+  subroutine polympo_setMPStress(mpMesh, nComps, numMPs, array) &
+             bind(C, NAME='polympo_setMPStress_f')
+    use :: iso_c_binding
+    type(c_ptr), value :: mpMesh
+    integer(c_int), value :: nComps, numMPs
+    type(c_ptr), value :: array
+  end subroutine
+
+  subroutine polympo_getMPStress(mpMesh, nComps, numMPs, array) &
+             bind(C, NAME='polympo_getMPStress_f')
+    use :: iso_c_binding
+    type(c_ptr), value :: mpMesh
+    integer(c_int), value :: nComps, numMPs
+    type(c_ptr), value :: array
+  end subroutine
+
   subroutine polympo_setAreaMP(mpMesh, nComps, numMPs, array) &
              bind(C, NAME='polympo_setAreaMP_f')
     use :: iso_c_binding