[cig-commits] r8040 - short/3D/PyLith/trunk/libsrc/materials

willic3 at geodynamics.org willic3 at geodynamics.org
Thu Sep 27 14:18:02 PDT 2007


Author: willic3
Date: 2007-09-27 14:18:02 -0700 (Thu, 27 Sep 2007)
New Revision: 8040

Modified:
   short/3D/PyLith/trunk/libsrc/materials/ElasticIsotropic3D.cc
   short/3D/PyLith/trunk/libsrc/materials/ElasticPlaneStrain.cc
   short/3D/PyLith/trunk/libsrc/materials/ElasticPlaneStress.cc
   short/3D/PyLith/trunk/libsrc/materials/ElasticStrain1D.cc
   short/3D/PyLith/trunk/libsrc/materials/ElasticStress1D.cc
   short/3D/PyLith/trunk/libsrc/materials/MaxwellIsotropic3D.cc
Log:
Added PETSc flop logging to material routines and eliminated some
unnecessary flops.
Not quite finished with MaxwellIsotropic3D.



Modified: short/3D/PyLith/trunk/libsrc/materials/ElasticIsotropic3D.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/materials/ElasticIsotropic3D.cc	2007-09-27 19:51:03 UTC (rev 8039)
+++ short/3D/PyLith/trunk/libsrc/materials/ElasticIsotropic3D.cc	2007-09-27 21:18:02 UTC (rev 8040)
@@ -16,6 +16,8 @@
 
 #include "pylith/utils/array.hh" // USES double_array
 
+#include "petsc.h" // USES PetscLogFlopsNoCheck
+
 #include <assert.h> // USES assert()
 #include <sstream> // USES std::ostringstream
 #include <stdexcept> // USES std::runtime_error
@@ -126,6 +128,8 @@
   (*paramVals)[_ElasticIsotropic3D::pidDensity][0] = density;
   (*paramVals)[_ElasticIsotropic3D::pidMu][0] = mu;
   (*paramVals)[_ElasticIsotropic3D::pidLambda][0] = lambda;
+
+  PetscLogFlopsNoCheck(6);
 } // _dbToParameters
 
 // ----------------------------------------------------------------------
@@ -179,8 +183,8 @@
   const double mu = parameters[_ElasticIsotropic3D::pidMu][0];
   const double lambda = parameters[_ElasticIsotropic3D::pidLambda][0];
 
-  const double lambda2mu = lambda + 2.0 * mu;
-  
+  const double mu2 = 2.0*mu;
+
   const double e11 = totalStrain[0];
   const double e22 = totalStrain[1];
   const double e33 = totalStrain[2];
@@ -190,12 +194,13 @@
   
   const double s123 = lambda * (e11 + e22 + e33);
 
-  (*stress)[0] = s123 + 2.0*mu*e11;
-  (*stress)[1] = s123 + 2.0*mu*e22;
-  (*stress)[2] = s123 + 2.0*mu*e33;
-  (*stress)[3] = 2.0 * mu * e12;
-  (*stress)[4] = 2.0 * mu * e23;
-  (*stress)[5] = 2.0 * mu * e13;
+  (*stress)[0] = s123 + mu2*e11;
+  (*stress)[1] = s123 + mu2*e22;
+  (*stress)[2] = s123 + mu2*e33;
+  (*stress)[3] = mu2 * e12;
+  (*stress)[4] = mu2 * e23;
+  (*stress)[5] = mu2 * e13;
+  PetscLogFlopsNoCheck(13);
 } // _calcStress
 
 // ----------------------------------------------------------------------
@@ -218,7 +223,8 @@
   const double mu = parameters[_ElasticIsotropic3D::pidMu][0];
   const double lambda = parameters[_ElasticIsotropic3D::pidLambda][0];
 
-  const double lambda2mu = lambda + 2.0 * mu;
+  const double mu2 = 2.0 * mu;
+  const double lambda2mu = lambda + mu2;
    
   (*elasticConsts)[ 0] = lambda2mu; // C1111
   (*elasticConsts)[ 1] = lambda; // C1122
@@ -235,12 +241,13 @@
   (*elasticConsts)[12] = 0; // C3312
   (*elasticConsts)[13] = 0; // C3323
   (*elasticConsts)[14] = 0; // C3313
-  (*elasticConsts)[15] = 2.0 * mu; // C1212
+  (*elasticConsts)[15] = mu2; // C1212
   (*elasticConsts)[16] = 0; // C1223
   (*elasticConsts)[17] = 0; // C1213
-  (*elasticConsts)[18] = 2.0 * mu; // C2323
+  (*elasticConsts)[18] = mu2; // C2323
   (*elasticConsts)[19] = 0; // C2313
-  (*elasticConsts)[20] = 2.0 * mu; // C1313
+  (*elasticConsts)[20] = mu2; // C1313
+  PetscLogFlopsNoCheck(2);
 } // _calcElasticConsts
 
 

Modified: short/3D/PyLith/trunk/libsrc/materials/ElasticPlaneStrain.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/materials/ElasticPlaneStrain.cc	2007-09-27 19:51:03 UTC (rev 8039)
+++ short/3D/PyLith/trunk/libsrc/materials/ElasticPlaneStrain.cc	2007-09-27 21:18:02 UTC (rev 8040)
@@ -16,6 +16,8 @@
 
 #include "pylith/utils/array.hh" // USES double_array
 
+#include "petsc.h" // USES PetscLogFlopsNoCheck
+
 #include <assert.h> // USES assert()
 #include <sstream> // USES std::ostringstream
 #include <stdexcept> // USES std::runtime_error
@@ -127,6 +129,8 @@
   (*paramVals)[_ElasticPlaneStrain::pidDensity][0] = density;
   (*paramVals)[_ElasticPlaneStrain::pidMu][0] = mu;
   (*paramVals)[_ElasticPlaneStrain::pidLambda][0] = lambda;
+
+  PetscLogFlopsNoCheck(6);
 } // computeParameters
 
 // ----------------------------------------------------------------------
@@ -180,7 +184,7 @@
   const double mu = parameters[_ElasticPlaneStrain::pidMu][0];
   const double lambda = parameters[_ElasticPlaneStrain::pidLambda][0];
   
-  const double lambda2mu = lambda + 2.0 * mu;
+  const double mu2 = 2.0 * mu;
   
   const double e11 = totalStrain[0];
   const double e22 = totalStrain[1];
@@ -188,9 +192,11 @@
 
   const double s12 = lambda * (e11 + e22);
 
-  (*stress)[0] = s12 + 2.0*mu*e11;
-  (*stress)[1] = s12 + 2.0*mu*e22;
-  (*stress)[2] = 2.0 * mu * e12;
+  (*stress)[0] = s12 + mu2*e11;
+  (*stress)[1] = s12 + mu2*e22;
+  (*stress)[2] = mu2 * e12;
+
+  PetscLogFlopsNoCheck(8);
 } // _calcStress
 
 // ----------------------------------------------------------------------
@@ -213,14 +219,17 @@
   const double mu = parameters[_ElasticPlaneStrain::pidMu][0];
   const double lambda = parameters[_ElasticPlaneStrain::pidLambda][0];
 
-  const double lambda2mu = lambda + 2.0*mu;
+  const double mu2 = 2.0 * mu;
+  const double lambda2mu = lambda + mu2;
   
   (*elasticConsts)[0] = lambda2mu; // C1111
   (*elasticConsts)[1] = lambda; // C1122
   (*elasticConsts)[2] = 0; // C1112
   (*elasticConsts)[3] = lambda2mu; // C2222
   (*elasticConsts)[4] = 0; // C2212
-  (*elasticConsts)[5] = 2.0*mu; // C1212
+  (*elasticConsts)[5] = mu2; // C1212
+
+  PetscLogFlopsNoCheck(2);
 } // calcElasticConsts
 
 

Modified: short/3D/PyLith/trunk/libsrc/materials/ElasticPlaneStress.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/materials/ElasticPlaneStress.cc	2007-09-27 19:51:03 UTC (rev 8039)
+++ short/3D/PyLith/trunk/libsrc/materials/ElasticPlaneStress.cc	2007-09-27 21:18:02 UTC (rev 8040)
@@ -16,6 +16,8 @@
 
 #include "pylith/utils/array.hh" // USES double_array
 
+#include "petsc.h" // USES PetscLogFlopsNoCheck
+
 #include <assert.h> // USES assert()
 #include <sstream> // USES std::ostringstream
 #include <stdexcept> // USES std::runtime_error
@@ -126,6 +128,8 @@
   (*paramVals)[_ElasticPlaneStress::pidDensity][0] = density;
   (*paramVals)[_ElasticPlaneStress::pidMu][0] = mu;
   (*paramVals)[_ElasticPlaneStress::pidLambda][0] = lambda;
+
+  PetscLogFlopsNoCheck(6);
 } // _dbToParameters
 
 // ----------------------------------------------------------------------
@@ -178,16 +182,19 @@
   const double mu = parameters[_ElasticPlaneStress::pidMu][0];
   const double lambda = parameters[_ElasticPlaneStress::pidLambda][0];
 
-  const double lambda2mu = lambda + 2.0 * mu;
+  const double mu2 = 2.0 * mu;
+  const double lambda2mu = lambda + mu2;
   const double lambdamu = lambda + mu;
   
   const double e11 = totalStrain[0];
   const double e22 = totalStrain[1];
   const double e12 = totalStrain[2];
 
-  (*stress)[0] = (4.0*mu*lambdamu * e11 + 2.0*mu*lambda * e22)/lambda2mu;
-  (*stress)[1] = (2.0*mu*lambda * e11 + 4.0*mu*lambdamu * e22)/lambda2mu;
-  (*stress)[2] = 2.0 * mu * e12;
+  (*stress)[0] = (2.0*mu2*lambdamu * e11 + mu2*lambda * e22)/lambda2mu;
+  (*stress)[1] = (mu2*lambda * e11 + 2.0*mu2*lambdamu * e22)/lambda2mu;
+  (*stress)[2] = mu2 * e12;
+
+  PetscLogFlopsNoCheck(18);
 } // _calcStress
 
 // ----------------------------------------------------------------------
@@ -210,15 +217,18 @@
   const double mu = parameters[_ElasticPlaneStress::pidMu][0];
   const double lambda = parameters[_ElasticPlaneStress::pidLambda][0];
   
-  const double lambda2mu = lambda + 2.0 * mu;
-  const double c11 = 4.0 * mu * (lambda + mu) / lambda2mu;
+  const double mu2 = 2.0 * mu;
+  const double lambda2mu = lambda + mu2;
+  const double c11 = 2.0 * mu2 * (lambda + mu) / lambda2mu;
 
   (*elasticConsts)[0] = c11; // C1111
-  (*elasticConsts)[1] = 2.0 * mu * lambda / lambda2mu; // C1122
+  (*elasticConsts)[1] = mu2 * lambda / lambda2mu; // C1122
   (*elasticConsts)[2] = 0; // C1112
   (*elasticConsts)[3] = c11; // C2222
   (*elasticConsts)[4] = 0; // C2212
-  (*elasticConsts)[5] = 2.0 * mu; // C1212
+  (*elasticConsts)[5] = mu2; // C1212
+
+  PetscLogFlopsNoCheck(8);
 } // calcElasticConsts
 
 

Modified: short/3D/PyLith/trunk/libsrc/materials/ElasticStrain1D.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/materials/ElasticStrain1D.cc	2007-09-27 19:51:03 UTC (rev 8039)
+++ short/3D/PyLith/trunk/libsrc/materials/ElasticStrain1D.cc	2007-09-27 21:18:02 UTC (rev 8040)
@@ -16,6 +16,8 @@
 
 #include "pylith/utils/array.hh" // USES double_array
 
+#include "petsc.h" // USES PetscLogFlopsNoCheck
+
 #include <assert.h> // USES assert()
 
 // ----------------------------------------------------------------------
@@ -109,6 +111,8 @@
 
   (*paramVals)[_ElasticStrain1D::pidDensity][0] = density;
   (*paramVals)[_ElasticStrain1D::pidLambda2mu][0] = lambda2mu;
+
+  PetscLogFlopsNoCheck(2);
 } // _dbToParameters
 
 // ----------------------------------------------------------------------
@@ -161,6 +165,8 @@
 
   const double e11 = totalStrain[0];
   (*stress)[0] = lambda2mu * e11;
+
+  PetscLogFlopsNoCheck(1);
 } // _calcStress
 
 // ----------------------------------------------------------------------

Modified: short/3D/PyLith/trunk/libsrc/materials/ElasticStress1D.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/materials/ElasticStress1D.cc	2007-09-27 19:51:03 UTC (rev 8039)
+++ short/3D/PyLith/trunk/libsrc/materials/ElasticStress1D.cc	2007-09-27 21:18:02 UTC (rev 8040)
@@ -16,6 +16,8 @@
 
 #include "pylith/utils/array.hh" // USES double_array
 
+#include "petsc.h" // USES PetscLogFlopsNoCheck
+
 #include <assert.h> // USES assert()
 
 // ----------------------------------------------------------------------
@@ -116,6 +118,8 @@
   (*paramVals)[_ElasticStress1D::pidDensity][0] = density;
   (*paramVals)[_ElasticStress1D::pidMu][0] = mu;
   (*paramVals)[_ElasticStress1D::pidLambda][0] = lambda;
+
+  PetscLogFlopsNoCheck(6);
 } // computeParameters
 
 // ----------------------------------------------------------------------
@@ -166,6 +170,8 @@
 
   const double e11 = totalStrain[0];
   (*stress)[0] = mu * (3.0*lambda+2.0*mu) / (lambda + mu) * e11;
+
+  PetscLogFlopsNoCheck(7);
 } // _calcStress
 
 // ----------------------------------------------------------------------
@@ -186,6 +192,8 @@
   const double lambda = parameters[_ElasticStress1D::pidLambda][0];
 
   (*elasticConsts)[0] = mu * (3.0*lambda+2.0*mu) / (lambda + mu);
+
+  PetscLogFlopsNoCheck(6);
 } // _calcElasticConsts
 
 

Modified: short/3D/PyLith/trunk/libsrc/materials/MaxwellIsotropic3D.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/materials/MaxwellIsotropic3D.cc	2007-09-27 19:51:03 UTC (rev 8039)
+++ short/3D/PyLith/trunk/libsrc/materials/MaxwellIsotropic3D.cc	2007-09-27 21:18:02 UTC (rev 8040)
@@ -16,6 +16,8 @@
 
 #include "pylith/utils/array.hh" // USES double_array
 
+#include "petsc.h" // USES PetscLogFlopsNoCheck
+
 #include <assert.h> // USES assert()
 #include <sstream> // USES std::ostringstream
 #include <stdexcept> // USES std::runtime_error
@@ -137,6 +139,8 @@
   (*paramVals)[_MaxwellIsotropic3D::pidMu][0] = mu;
   (*paramVals)[_MaxwellIsotropic3D::pidLambda][0] = lambda;
   (*paramVals)[_MaxwellIsotropic3D::pidMaxwellTime][0] = maxwelltime;
+
+  PetscLogFlopsNoCheck(7);
 } // _dbToParameters
 
 // ----------------------------------------------------------------------
@@ -194,8 +198,7 @@
   const double lambda = parameters[_MaxwellIsotropic3D::pidLambda][0];
   const double maxwelltime = parameters[_MaxwellIsotropic3D::pidMaxwellTime][0];
 
-  const double lambda2mu = lambda + 2.0 * mu;
-  const double bulkmodulus = lambda + 2.0 * mu/3.0;
+  const double mu2 = 2.0 * mu;
 
   const double e11 = totalStrain[0];
   const double e22 = totalStrain[1];
@@ -205,15 +208,14 @@
   const double e13 = totalStrain[5];
   
   const double traceStrainTpdt = e11 + e22 + e33;
-  const double meanStrainTpdt = traceStrainTpdt/3.0;
   const double s123 = lambda * traceStrainTpdt;
 
-  (*stress)[0] = s123 + 2.0*mu*e11;
-  (*stress)[1] = s123 + 2.0*mu*e22;
-  (*stress)[2] = s123 + 2.0*mu*e33;
-  (*stress)[3] = 2.0 * mu * e12;
-  (*stress)[4] = 2.0 * mu * e23;
-  (*stress)[5] = 2.0 * mu * e13;
+  (*stress)[0] = s123 + mu2*e11;
+  (*stress)[1] = s123 + mu2*e22;
+  (*stress)[2] = s123 + mu2*e33;
+  (*stress)[3] = mu2 * e12;
+  (*stress)[4] = mu2 * e23;
+  (*stress)[5] = mu2 * e13;
   // std::cout << " _calcStressElastic: " << std::endl;
   // std::cout << " totalStrain: " << std::endl;
   // for (int iComp=0; iComp < _MaxwellIsotropic3D::tensorSize; ++iComp)
@@ -222,6 +224,8 @@
   // for (int iComp=0; iComp < _MaxwellIsotropic3D::tensorSize; ++iComp)
     // std::cout << "  " << (*stress)[iComp];
   // std::cout << std::endl;
+
+  PetscLogFlopsNoCheck(13);
 } // _calcStressElastic
 
 // ----------------------------------------------------------------------
@@ -248,8 +252,8 @@
   const double lambda = parameters[_MaxwellIsotropic3D::pidLambda][0];
   const double maxwelltime = parameters[_MaxwellIsotropic3D::pidMaxwellTime][0];
 
-  const double lambda2mu = lambda + 2.0 * mu;
-  const double bulkmodulus = lambda + 2.0 * mu/3.0;
+  const double mu2 = 2.0 * mu;
+  const double bulkmodulus = lambda + mu2/3.0;
 
   const double e11 = totalStrain[0];
   const double e22 = totalStrain[1];
@@ -275,6 +279,7 @@
 			      parameters[_MaxwellIsotropic3D::pidStrainT][1] +
 			      parameters[_MaxwellIsotropic3D::pidStrainT][2])/3.0;
   
+  PetscLogFlopsNoCheck(11);
   // The code below should probably be in a separate function since it
   // is used more than once.  I should also probably cover the possibility
   // that Maxwell time is zero (although this should never happen).
@@ -292,8 +297,11 @@
       fraction *= _dt/maxwelltime;
       dq += fSign*fraction/factorial;
     } // for
-  } else
+    PetscLogFlopsNoCheck(1+7*numTerms);
+  } else {
     dq = maxwelltime*(1.0-exp(-_dt/maxwelltime))/_dt;
+    PetscLogFlopsNoCheck(7);
+  } // else
 
   const double expFac = exp(-_dt/maxwelltime);
   const double elasFac = 2.0*mu;
@@ -301,6 +309,7 @@
   double devStrainT = 0.0;
   double devStressTpdt = 0.0;
   double visStrain = 0.0;
+  PetscLogFlopsNoCheck(4);
   // std::cout << " _calcStressViscoelastic: " << std::endl;
   // std::cout << " stress  totalStrain  visStrain: " << std::endl;
   for (int iComp=0; iComp < _MaxwellIsotropic3D::tensorSize; ++iComp) {
@@ -316,6 +325,7 @@
     // Temporary to get stresses and strains.
     // std::cout << "  " << (*stress)[iComp] << "  " << totalStrain[iComp] << "  " << visStrain << std:: endl;
   } // for
+  PetscLogFlopsNoCheck(11 * _MaxwellIsotropic3D::tensorSize);
 } // _calcStress
 
 // ----------------------------------------------------------------------
@@ -336,8 +346,9 @@
   const double lambda = parameters[_MaxwellIsotropic3D::pidLambda][0];
   const double maxwelltime = parameters[_MaxwellIsotropic3D::pidMaxwellTime][0];
 
-  const double lambda2mu = lambda + 2.0 * mu;
-  const double bulkmodulus = lambda + 2.0 * mu/3.0;
+  const double mu2 = 2.0 * mu;
+  const double lambda2mu = lambda + mu2;
+  const double bulkmodulus = lambda + mu2/3.0;
 
   (*elasticConsts)[ 0] = lambda2mu; // C1111
   (*elasticConsts)[ 1] = lambda; // C1122
@@ -354,12 +365,14 @@
   (*elasticConsts)[12] = 0; // C3312
   (*elasticConsts)[13] = 0; // C3323
   (*elasticConsts)[14] = 0; // C3313
-  (*elasticConsts)[15] = 2.0 * mu; // C1212
+  (*elasticConsts)[15] = mu2; // C1212
   (*elasticConsts)[16] = 0; // C1223
   (*elasticConsts)[17] = 0; // C1213
-  (*elasticConsts)[18] = 2.0 * mu; // C2323
+  (*elasticConsts)[18] = mu2; // C2323
   (*elasticConsts)[19] = 0; // C2313
-  (*elasticConsts)[20] = 2.0 * mu; // C1313
+  (*elasticConsts)[20] = mu2; // C1313
+
+  PetscLogFlopsNoCheck(4);
 } // _calcElasticConstsElastic
 
 // ----------------------------------------------------------------------
@@ -381,12 +394,14 @@
   const double lambda = parameters[_MaxwellIsotropic3D::pidLambda][0];
   const double maxwelltime = parameters[_MaxwellIsotropic3D::pidMaxwellTime][0];
 
-  const double lambda2mu = lambda + 2.0 * mu;
-  const double bulkmodulus = lambda + 2.0 * mu/3.0;
+  const double mu2 = 2.0 * mu;
+  const double bulkmodulus = lambda + mu2/3.0;
 
   const double timeFrac = 1.0e-5;
   const int numTerms = 5;
   double dq = 0.0;
+
+  PetscLogFlopsNoCheck(3);
   if(maxwelltime < timeFrac*_dt) {
     double fSign = 1.0;
     double factorial = 1.0;
@@ -398,8 +413,11 @@
       fraction *= _dt/maxwelltime;
       dq += fSign*fraction/factorial;
     } // for
-  } else
+    PetscLogFlopsNoCheck(1+7*numTerms);
+  } else {
     dq = maxwelltime*(1.0-exp(-_dt/maxwelltime))/_dt;
+    PetscLogFlopsNoCheck(7);
+  } // else
 
   const double visFac = mu*dq/3.0;
   (*elasticConsts)[ 0] = bulkmodulus + 4.0*visFac; // C1111
@@ -423,6 +441,8 @@
   (*elasticConsts)[18] = (*elasticConsts)[15]; // C2323
   (*elasticConsts)[19] = 0; // C2313
   (*elasticConsts)[20] = (*elasticConsts)[15]; // C1313
+
+  PetscLogFlopsNoCheck(7);
 } // _calcElasticConstsViscoelastic
 
 // ----------------------------------------------------------------------
@@ -454,14 +474,15 @@
   const double diag[] = { 1.0, 1.0, 1.0, 0.0, 0.0, 0.0 };
 
   // Temporary to get stresses.
-  double_array stress(6);
-  _calcStressElastic(&stress, (*parameters), totalStrain);
+  // double_array stress(6);
+  // _calcStressElastic(&stress, (*parameters), totalStrain);
 
   for (int iComp=0; iComp < _MaxwellIsotropic3D::tensorSize; ++iComp) {
     (*parameters)[_MaxwellIsotropic3D::pidStrainT][iComp] = totalStrain[iComp];
     (*parameters)[_MaxwellIsotropic3D::pidVisStrain][iComp] =
       totalStrain[iComp] - diag[iComp]*meanStrainTpdt;
   } // for
+  PetscLogFlopsNoCheck(5 * _MaxwellIsotropic3D::tensorSize);
 //   std::cout << std::endl;
 //   std::cout << " updateStateElastic: "<< std::endl;
 //   std::cout << " StrainT  VisStrain  Stress: " << std::endl;
@@ -473,6 +494,8 @@
   _needNewJacobian = true;
 } // _updateStateElastic
 
+// **************  Finish adding PETSc logging from here *******************
+
 // ----------------------------------------------------------------------
 // Update state variables.
 void



More information about the cig-commits mailing list