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

[knepley at geodynamics.org]

Author: knepley
Date: 2007-05-09 12:05:31 -0700 (Wed, 09 May 2007)
New Revision: 6808

Modified:
   short/3D/PyLith/trunk/modulesrc/solver/solver.pyxe.src
   short/3D/PyLith/trunk/pylith/feassemble/FIATLagrange.py
Log:
Fixed VecScatter calls, Add FIATLagrange stuff (need tests)


Modified: short/3D/PyLith/trunk/modulesrc/solver/solver.pyxe.src
===================================================================
--- short/3D/PyLith/trunk/modulesrc/solver/solver.pyxe.src	2007-05-09 18:26:50 UTC (rev 6807)
+++ short/3D/PyLith/trunk/modulesrc/solver/solver.pyxe.src	2007-05-09 19:05:31 UTC (rev 6808)
@@ -209,10 +209,10 @@
       err = VecCreateSeqWithArray(PETSC_COMM_SELF, (*fieldIn)->size(),
                                   (*fieldIn)->restrict(), &localVec);CHKERRQ(err);
       err = VecDuplicate(localVec, &vecIn); CHKERRQ(err);
-      err = VecScatterBegin(localVec, vecIn, INSERT_VALUES, SCATTER_FORWARD,
-                            scatter);CHKERRQ(err);
-      err = VecScatterEnd(localVec, vecIn, INSERT_VALUES, SCATTER_FORWARD,
-                          scatter); CHKERRQ(err);
+      err = VecScatterBegin(scatter, localVec, vecIn, INSERT_VALUES, SCATTER_FORWARD
+                            );CHKERRQ(err);
+      err = VecScatterEnd(scatter, localVec, vecIn, INSERT_VALUES, SCATTER_FORWARD
+                          ); CHKERRQ(err);
       err = VecDestroy(localVec); CHKERRQ(err);
       err = KSPSetOperators(*ksp, *jacobian, *jacobian,
                             DIFFERENT_NONZERO_PATTERN); CHKERRQ(err);
@@ -220,10 +220,10 @@
       err = KSPSolve(*ksp, vecIn, vecOut); CHKERRQ(err);
       err = VecCreateSeqWithArray(PETSC_COMM_SELF, (*fieldOut)->size(),
                                 (*fieldOut)->restrict(), &localVec);CHKERRQ(err);
-      err = VecScatterBegin(vecOut, localVec, INSERT_VALUES, SCATTER_REVERSE,
-                            scatter); CHKERRQ(err);
-      err = VecScatterEnd(vecOut, localVec, INSERT_VALUES, SCATTER_REVERSE,
-                          scatter); CHKERRQ(err);
+      err = VecScatterBegin(scatter, vecOut, localVec, INSERT_VALUES, SCATTER_REVERSE
+                            ); CHKERRQ(err);
+      err = VecScatterEnd(scatter, vecOut, localVec, INSERT_VALUES, SCATTER_REVERSE
+                          ); CHKERRQ(err);
       err = VecDestroy(localVec); CHKERRQ(err);
     } catch (const std::exception& err) {
       PyErr_SetString(PyExc_RuntimeError,

Modified: short/3D/PyLith/trunk/pylith/feassemble/FIATLagrange.py
===================================================================
--- short/3D/PyLith/trunk/pylith/feassemble/FIATLagrange.py	2007-05-09 18:26:50 UTC (rev 6807)
+++ short/3D/PyLith/trunk/pylith/feassemble/FIATLagrange.py	2007-05-09 19:05:31 UTC (rev 6808)
@@ -60,7 +60,7 @@
     degree = pyre.inventory.int("degree", default=1)
     degree.meta['tip'] = "Degree of finite-element cell."
 
-    quadOrder = pyre.inventory.int("quad_order", default=1)
+    quadOrder = pyre.inventory.int("quad_order", default=-1)
     quadOrder.meta['tip'] = "Order of quadrature rule."
     
 
@@ -80,10 +80,21 @@
     """
     Set members based using inventory.
     """
+    import FIAT.shapes
+
     FIATCell._configure(self)
-    self.dimension = self.inventory.dimension
+    self.cellDim = self.inventory.dimension
     self.degree = self.inventory.degree
-    self.order = self.inventory.order
+    self.quadOrder = self.inventory.quadOrder
+    if self.cellDim == 1:
+      self.shape = FIAT.shapes.LINE
+    elif self.cellDim == 2:
+      self.shape = FIAT.shapes.TRIANGLE
+    elif self.cellDim == 3:
+      self.shape = FIAT.shapes.TETRAHEDRON
+    if self.quadOrder == -1:
+      self.quadOrder = 2*self.degree+1
+    self.numCorners = self.cellDim+1
     return
 
 
@@ -91,10 +102,11 @@
     """
     Setup quadrature rule for reference cell.
     """
-
-    # ADD STUFF HERE
-    raise NotImplementedError()
-
+    import FIAT.quadrature
+    self.quadrature = FIAT.quadrature.make_quadrature_by_degree(shape, self.quadOrder)
+    self.numQuadPts = len(quadrature.get_points())
+    self.quadPts = quadrature.get_points()
+    self.quadWts = quadrature.get_weights()
     return
 
 
@@ -104,9 +116,12 @@
     """
     from FIAT.Lagrange import Lagrange
 
-    # ADD STUFF HERE
-    raise NotImplementedError()
-
+    self.element = Lagrange(self.shape, self.degree)
+    points = self.quadrature.get_points()
+    basis = self.element.function_space()
+    self.numBasisFuncs = len(basis)
+    self.basis = numpy.transpose(basis.tabulate(points))
+    self.basisDeriv = numpy.transpose([basis.deriv_all(d).tabulate(points) for d in range(self.cellDim)])
     return
 
 # FACTORIES ////////////////////////////////////////////////////////////