[cig-commits] r22758 - short/3D/PyLith/trunk/libsrc/pylith/faults

[brad at geodynamics.org]

Author: brad
Date: 2013-09-03 10:29:42 -0700 (Tue, 03 Sep 2013)
New Revision: 22758

Modified:
   short/3D/PyLith/trunk/libsrc/pylith/faults/FaultCohesiveDyn.cc
   short/3D/PyLith/trunk/libsrc/pylith/faults/FaultCohesiveLagrange.cc
Log:
Updated FaultCohesiveDyn and FaultCohesiveLagrange to use optimized getClosure (pass in work array as arg).

Modified: short/3D/PyLith/trunk/libsrc/pylith/faults/FaultCohesiveDyn.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/pylith/faults/FaultCohesiveDyn.cc	2013-09-03 17:22:37 UTC (rev 22757)
+++ short/3D/PyLith/trunk/libsrc/pylith/faults/FaultCohesiveDyn.cc	2013-09-03 17:29:42 UTC (rev 22758)
@@ -1823,8 +1823,10 @@
   const PetscInt numCells = cellsStratum.size();
 
   // Get sections
+  scalar_array coordsCell(numBasis*spaceDim); // :KULDGE: Update numBasis to numCorners after implementing higher order
   topology::CoordsVisitor coordsVisitor(faultDMMesh);
 
+  scalar_array dLagrangeCell(numBasis*spaceDim);
   topology::VecVisitorMesh dLagrangeVisitor(_fields->get("sensitivity dLagrange"));
 
   scalar_array residualCell(numBasis*spaceDim);
@@ -1835,16 +1837,11 @@
   // Loop over cells
   for(PetscInt c = cStart; c < cEnd; ++c) {
     // Compute geometry
-    PetscScalar *coordsCell = NULL;
-    PetscInt coordsSize = 0;
-    coordsVisitor.getClosure(&coordsCell, &coordsSize, c);assert(coordsCell);assert(numBasis*spaceDim == coordsSize);
-    _quadrature->computeGeometry(coordsCell, coordsSize, c);
-    coordsVisitor.restoreClosure(&coordsCell, &coordsSize, c);
+    coordsVisitor.getClosure(&coordsCell, c);
+    _quadrature->computeGeometry(&coordsCell[0], coordsCell.size(), c);
 
     // Restrict input fields to cell
-    PetscScalar* dLagrangeCell = NULL;
-    PetscInt dLagrangeSize = 0;
-    dLagrangeVisitor.getClosure(&dLagrangeCell, &dLagrangeSize, c);assert(dLagrangeCell);assert(numBasis*spaceDim == dLagrangeSize);
+    dLagrangeVisitor.getClosure(&dLagrangeCell, c);
 
     // Get cell geometry information that depends on cell
     const scalar_array& basis = _quadrature->basis();
@@ -1876,7 +1873,6 @@
         } // for
       } // for
     } // for
-    dLagrangeVisitor.restoreClosure(&dLagrangeCell, &dLagrangeSize, c);
 
     // Assemble cell contribution into field
     residualVisitor.setClosure(&residualCell[0], residualCell.size(), c, ADD_VALUES);

Modified: short/3D/PyLith/trunk/libsrc/pylith/faults/FaultCohesiveLagrange.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/pylith/faults/FaultCohesiveLagrange.cc	2013-09-03 17:22:37 UTC (rev 22757)
+++ short/3D/PyLith/trunk/libsrc/pylith/faults/FaultCohesiveLagrange.cc	2013-09-03 17:29:42 UTC (rev 22758)
@@ -1316,6 +1316,7 @@
   PetscScalar* orientationArray = orientationVisitor.localArray();
 
   // Get section containing coordinates of vertices
+  scalar_array coordsCell(numBasis*spaceDim); // :KULDGE: Update numBasis to numCorners after implementing higher order
   topology::CoordsVisitor coordsVisitor(faultDMMesh);
 
   // Loop over cohesive cells, computing orientation weighted by
@@ -1326,11 +1327,8 @@
     PetscInt closureSize, q = 0;
 
     // Get orientations at fault cell's vertices.
+    coordsVisitor.getClosure(&coordsCell, c);
 
-    PetscScalar *coordsCell = NULL;
-    PetscInt coordsSize = 0;
-    coordsVisitor.getClosure(&coordsCell, &coordsSize, c);assert(coordsCell);assert(numBasis*spaceDim == coordsSize);
-
     PetscErrorCode err = DMPlexGetTransitiveClosure(faultDMMesh, c, PETSC_TRUE, &closureSize, &closure);PYLITH_CHECK_ERROR(err);
 
     // Filter out non-vertices
@@ -1344,7 +1342,7 @@
     closureSize = q;
     for(PetscInt v = 0; v < closureSize; ++v) {
       // Compute Jacobian and determinant of Jacobian at vertex
-      cellGeometry.jacobian(&jacobian, &jacobianDet, coordsCell, numBasis, spaceDim, &verticesRef[v*cohesiveDim], cohesiveDim);
+      cellGeometry.jacobian(&jacobian, &jacobianDet, &coordsCell[0], numBasis, spaceDim, &verticesRef[v*cohesiveDim], cohesiveDim);
 
       // Compute orientation
       cellGeometry.orientation(&orientationVertex, jacobian, jacobianDet, up);
@@ -1355,7 +1353,6 @@
         orientationArray[ooff+d] += orientationVertex[d];
       } // for
     } // for
-    coordsVisitor.restoreClosure(&coordsCell, &coordsSize, c);
     err = DMPlexRestoreTransitiveClosure(faultDMMesh, c, PETSC_TRUE, &closureSize, &closure);PYLITH_CHECK_ERROR(err);
   } // for
   orientationVisitor.clear();
@@ -1568,6 +1565,7 @@
   topology::VecVisitorMesh areaVisitor(area);
   scalar_array areaCell(numBasis);
 
+  scalar_array coordsCell(numBasis*spaceDim); // :KULDGE: Update numBasis to numCorners after implementing higher order
   topology::CoordsVisitor coordsVisitor(faultDMMesh);
 
   // Loop over cells in fault mesh, compute area
@@ -1575,11 +1573,8 @@
     areaCell = 0.0;
 
     // Compute geometry information for current cell
-    PetscScalar* coordsCell = NULL;
-    PetscInt coordsSize = 0;
-    coordsVisitor.getClosure(&coordsCell, &coordsSize, c);assert(coordsCell);assert(numBasis*spaceDim == coordsSize);
-    _quadrature->computeGeometry(coordsCell, coordsSize, c);
-    coordsVisitor.restoreClosure(&coordsCell, &coordsSize, c);
+    coordsVisitor.getClosure(&coordsCell, c);
+    _quadrature->computeGeometry(&coordsCell[0], coordsCell.size(), c);
 
     // Get cell geometry information that depends on cell
     const scalar_array& basis = _quadrature->basis();