[cig-commits] r7014 - in short/3D/PyLith/trunk: examples/twotri3 modulesrc/solver pylith/feassemble pylith/solver

[knepley at geodynamics.org]

Author: knepley
Date: 2007-05-31 10:32:08 -0700 (Thu, 31 May 2007)
New Revision: 7014

Modified:
   short/3D/PyLith/trunk/examples/twotri3/pylithapp.cfg
   short/3D/PyLith/trunk/modulesrc/solver/solver.pyxe.src
   short/3D/PyLith/trunk/pylith/feassemble/FIATCell.py
   short/3D/PyLith/trunk/pylith/solver/SolverLinear.py
Log:
Now twotri3 is working again (tensor elements incomplete)


Modified: short/3D/PyLith/trunk/examples/twotri3/pylithapp.cfg
===================================================================
--- short/3D/PyLith/trunk/examples/twotri3/pylithapp.cfg	2007-05-31 16:48:04 UTC (rev 7013)
+++ short/3D/PyLith/trunk/examples/twotri3/pylithapp.cfg	2007-05-31 17:32:08 UTC (rev 7014)
@@ -7,6 +7,7 @@
 [pylithapp.journal.info]
 eqdeformation = 1
 explicit = 1
+implicit = 1
 petsc = 1
 meshioascii = 1
 homogeneous = 1
@@ -28,7 +29,7 @@
 # problem
 # ----------------------------------------------------------------------
 [pylithapp.timedependent]
-total_time = 1.0*s
+total_time = 2.0*s
 default_dt = 1.0*s
 dimension = 2
 

Modified: short/3D/PyLith/trunk/modulesrc/solver/solver.pyxe.src
===================================================================
--- short/3D/PyLith/trunk/modulesrc/solver/solver.pyxe.src	2007-05-31 16:48:04 UTC (rev 7013)
+++ short/3D/PyLith/trunk/modulesrc/solver/solver.pyxe.src	2007-05-31 17:32:08 UTC (rev 7014)
@@ -90,7 +90,7 @@
         PetscError(__LINE__,__FUNCT__,__FILE__,__SDIR__,err,0," ");
         throw std::runtime_error("Could not create KSP object.");
       } // if
-      //err = KSPSetFromOptions(*ksp);
+      err = KSPSetFromOptions(*ksp);
       if (err) {
         PetscError(__LINE__,__FUNCT__,__FILE__,__SDIR__,err,0," ");
         throw std::runtime_error("Could not set KSP options.");

Modified: short/3D/PyLith/trunk/pylith/feassemble/FIATCell.py
===================================================================
--- short/3D/PyLith/trunk/pylith/feassemble/FIATCell.py	2007-05-31 16:48:04 UTC (rev 7013)
+++ short/3D/PyLith/trunk/pylith/feassemble/FIATCell.py	2007-05-31 17:32:08 UTC (rev 7014)
@@ -80,6 +80,46 @@
     return
 
 
+  def initializeTensor(self, dim):
+    """
+    Initialize reference finite-element cell from a tensor product of 1D Lagrange elements.
+    """
+    quadrature = self._setupQuadrature()
+    basisFns = self._setupBasisFns()
+    
+    # Evaluate basis functions at quadrature points
+    quadpts = quadrature.get_points()
+    basis = numpy.array(basisFns.tabulate(quadpts)).transpose()
+    self.basis = numpy.reshape(basis.flatten(), basis.shape)
+
+    # Evaluate derivatives of basis functions at quadrature points
+    import FIAT.shapes
+    dim = FIAT.shapes.dimension(basisFns.base.shape)
+    basisDeriv = numpy.array([basisFns.deriv_all(d).tabulate(quadpts) \
+                              for d in range(dim)]).transpose()
+    self.basisDeriv = numpy.reshape(basisDeriv.flatten(), basisDeriv.shape)
+
+    self.quadPts = numpy.array(quadrature.get_points())
+    self.quadWts = numpy.array(quadrature.get_weights())
+
+    self.cellDim = dim
+    self.numCorners = len(basisFns)
+    self.numQuadPts = len(quadrature.get_weights())
+
+    self._info.line("Basis (quad pts):")
+    self._info.line(self.basis)
+    self._info.line("Basis derivatives (quad pts):")
+    self._info.line(self.basisDeriv)
+    self._info.line("Quad pts:")
+    self._info.line(quadrature.get_points())
+    self._info.line("Quad wts:")
+    self._info.line(quadrature.get_weights())
+
+    self._info.log()
+
+    return
+
+
   # PRIVATE METHODS ////////////////////////////////////////////////////
 
   def _setupQuadrature(self):

Modified: short/3D/PyLith/trunk/pylith/solver/SolverLinear.py
===================================================================
--- short/3D/PyLith/trunk/pylith/solver/SolverLinear.py	2007-05-31 16:48:04 UTC (rev 7013)
+++ short/3D/PyLith/trunk/pylith/solver/SolverLinear.py	2007-05-31 17:32:08 UTC (rev 7014)
@@ -66,6 +66,7 @@
     """
     Solve linear system.
     """
+    self._info.log("Solving linear equations.")
     self.cppHandle.solve(fieldOut, jacobian, fieldIn)
     return