[cig-commits] r5106 - short/3D/PyLith/trunk/libsrc/feassemble
brad at geodynamics.org
brad at geodynamics.org
Thu Oct 26 21:52:46 PDT 2006
Author: brad
Date: 2006-10-26 21:52:45 -0700 (Thu, 26 Oct 2006)
New Revision: 5106
Added:
short/3D/PyLith/trunk/libsrc/feassemble/IntegratorInertia.cc
short/3D/PyLith/trunk/libsrc/feassemble/IntegratorInertia.hh
short/3D/PyLith/trunk/libsrc/feassemble/IntegratorInertia.icc
Removed:
short/3D/PyLith/trunk/libsrc/feassemble/IntegratorInertia3D.cc
short/3D/PyLith/trunk/libsrc/feassemble/IntegratorInertia3D.hh
short/3D/PyLith/trunk/libsrc/feassemble/IntegratorInertia3D.icc
Modified:
short/3D/PyLith/trunk/libsrc/feassemble/Makefile.am
Log:
Name change because inertia operator is dimension independent.
Copied: short/3D/PyLith/trunk/libsrc/feassemble/IntegratorInertia.cc (from rev 5105, short/3D/PyLith/trunk/libsrc/feassemble/IntegratorInertia3D.cc)
===================================================================
--- short/3D/PyLith/trunk/libsrc/feassemble/IntegratorInertia3D.cc 2006-10-27 04:49:03 UTC (rev 5105)
+++ short/3D/PyLith/trunk/libsrc/feassemble/IntegratorInertia.cc 2006-10-27 04:52:45 UTC (rev 5106)
@@ -0,0 +1,251 @@
+// -*- C++ -*-
+//
+// ======================================================================
+//
+// Brad T. Aagaard
+// U.S. Geological Survey
+//
+// {LicenseText}
+//
+// ======================================================================
+//
+
+#include <portinfo>
+
+#include "IntegratorInertia.hh" // implementation of class methods
+
+#include "Quadrature.hh" // USES Quadrature
+
+#include "petscmat.h" // USES PetscMat
+
+#include <assert.h> // USES assert()
+#include <stdexcept> // USES std::runtime_error
+
+// ----------------------------------------------------------------------
+// Constructor
+pylith::feassemble::IntegratorInertia::IntegratorInertia(void)
+{ // constructor
+} // constructor
+
+// ----------------------------------------------------------------------
+// Destructor
+pylith::feassemble::IntegratorInertia::~IntegratorInertia(void)
+{ // destructor
+} // destructor
+
+// ----------------------------------------------------------------------
+// Copy constructor.
+pylith::feassemble::IntegratorInertia::IntegratorInertia(const IntegratorInertia& i) :
+ Integrator(i)
+{ // copy constructor
+} // copy constructor
+
+// ----------------------------------------------------------------------
+// Integrate inertial term for 3-D finite elements.
+void
+pylith::feassemble::IntegratorInertia::integrateAction(
+ const ALE::Obj<real_section_type>& fieldOut,
+ const ALE::Obj<real_section_type>& fieldIn,
+ const ALE::Obj<real_section_type>& coordinates)
+{ // integrateAction
+ assert(0 != _quadrature);
+
+ // Get information about section
+ const topology_type::patch_type patch = 0;
+ const ALE::Obj<topology_type>& topology = fieldIn->getTopology();
+ const ALE::Obj<topology_type::label_sequence>& cells =
+ topology->heightStratum(patch, 0);
+ const topology_type::label_sequence::iterator cellsEnd = cells->end();
+
+ // Allocate vector for cell values (if necessary)
+ _initCellVector();
+
+ for (topology_type::label_sequence::iterator cellIter=cells->begin();
+ cellIter != cellsEnd;
+ ++cellIter) {
+ // Compute geometry information for current cell
+ _quadrature->computeGeometry(coordinates, *cellIter);
+
+ // Reset element vector to zero
+ _resetCellVector();
+
+ // Restrict input field to cell
+ const real_section_type::value_type* fieldInCell =
+ fieldIn->restrict(patch, *cellIter);
+
+ // Get cell geometry information
+ const int numQuadPts = _quadrature->numQuadPts();
+ const double* basis = _quadrature->basis();
+ const double* quadPts = _quadrature->quadPts();
+ const double* quadWts = _quadrature->quadWts();
+ const double* jacobianDet = _quadrature->jacobianDet();
+ const int numBasis = _quadrature->numCorners();
+ const int spaceDim = _quadrature->spaceDim();
+
+ // FIX THIS
+ // Hardwire mass density
+ const double density = 1.0;
+
+ // Compute action for cell
+ for (int iQuad=0; iQuad < numQuadPts; ++iQuad) {
+ const double wt = quadWts[iQuad] * jacobianDet[iQuad] * density;
+ for (int iBasis=0, iQ=iQuad*numBasis; iBasis < numBasis; ++iBasis) {
+ const int iBlock = iBasis * spaceDim;
+ double val = wt*basis[iQ+iBasis];
+ for (int jBasis=0; jBasis < numBasis; ++jBasis) {
+ val *= basis[iQ+jBasis];
+ for (int iDim=0; iDim < spaceDim; ++iDim)
+ _cellVector[iBlock+iDim] += val * fieldInCell[iBlock+iDim];
+ } // for
+ } // for
+ } // for
+ PetscErrorCode err =
+ PetscLogFlops(numQuadPts*(2+numBasis*(1+numBasis*(1+2*spaceDim))));
+ if (err)
+ throw std::runtime_error("Couldn't log PETSc flops.");
+
+ // Assemble cell contribution into field
+ fieldOut->updateAdd(patch, *cellIter, _cellVector);
+ } // for
+} // integrateAction
+
+// ----------------------------------------------------------------------
+// Compute matrix associated with operator.
+void
+pylith::feassemble::IntegratorInertia::integrate(
+ PetscMat* mat,
+ const ALE::Obj<real_section_type>& coordinates)
+{ // integrate
+ assert(0 != mat);
+ assert(0 != _quadrature);
+
+ // Get information about section
+ const topology_type::patch_type patch = 0;
+ const ALE::Obj<topology_type>& topology = coordinates->getTopology();
+ const ALE::Obj<topology_type::label_sequence>& cells =
+ topology->heightStratum(patch, 0);
+ const topology_type::label_sequence::iterator cellsEnd = cells->end();
+
+ // Allocate matrix for cell values (if necessary)
+ _initCellMatrix();
+
+ for (topology_type::label_sequence::iterator cellIter=cells->begin();
+ cellIter != cellsEnd;
+ ++cellIter) {
+ // Compute geometry information for current cell
+ _quadrature->computeGeometry(coordinates, *cellIter);
+
+ // Reset element matrix to zero
+ _resetCellMatrix();
+
+ // Get cell geometry information
+ const int numQuadPts = _quadrature->numQuadPts();
+ const double* basis = _quadrature->basis();
+ const double* quadPts = _quadrature->quadPts();
+ const double* quadWts = _quadrature->quadWts();
+ const double* jacobianDet = _quadrature->jacobianDet();
+ const int numBasis = _quadrature->numCorners();
+ const int spaceDim = _quadrature->spaceDim();
+
+ // FIX THIS
+ // Hardwire mass density
+ const double density = 1.0;
+
+ // Integrate cell
+ for (int iQuad=0; iQuad < numQuadPts; ++iQuad) {
+ const double wt = quadWts[iQuad] * jacobianDet[iQuad] * density;
+ for (int iBasis=0, iQ=iQuad*numBasis; iBasis < numBasis; ++iBasis) {
+ const int iBlock = iBasis * spaceDim;
+ double val = wt*basis[iQ+iBasis];
+ for (int jBasis=0; jBasis < numBasis; ++jBasis) {
+ const int jBlock = jBasis * spaceDim;
+ val *= basis[iQ+jBasis];
+ for (int iDim=0; iDim < spaceDim; ++iDim)
+ _cellMatrix[(iBlock+iDim)*(numBasis*spaceDim)+jBlock+iDim] += val;
+ } // for
+ } // for
+ } // for
+ PetscErrorCode err =
+ PetscLogFlops(numQuadPts*(2+numBasis*(1+numBasis*(1+spaceDim))));
+ if (err)
+ throw std::runtime_error("Couldn't log PETSc flops.");
+
+ // Assemble cell contribution into sparse matrix
+ // STUFF GOES HERE
+ } // for
+} // integrate
+
+// ----------------------------------------------------------------------
+// Compute lumped matrix associated with operator.
+void
+pylith::feassemble::IntegratorInertia::integrateLumped(
+ const ALE::Obj<real_section_type>& fieldOut,
+ const ALE::Obj<real_section_type>& coordinates)
+{ // integrateLumped
+ assert(0 != _quadrature);
+
+ // Get information about section
+ const topology_type::patch_type patch = 0;
+ const ALE::Obj<topology_type>& topology = coordinates->getTopology();
+ const ALE::Obj<topology_type::label_sequence>& cells =
+ topology->heightStratum(patch, 0);
+ const topology_type::label_sequence::iterator cellsEnd = cells->end();
+
+ // Allocate matrix for cell values (if necessary)
+ _initCellVector();
+
+ for (topology_type::label_sequence::iterator cellIter=cells->begin();
+ cellIter != cellsEnd;
+ ++cellIter) {
+ // Compute geometry information for current cell
+ _quadrature->computeGeometry(coordinates, *cellIter);
+
+ // Reset element matrix to zero
+ _resetCellVector();
+
+ // Get cell geometry information
+ const int numQuadPts = _quadrature->numQuadPts();
+ const double* basis = _quadrature->basis();
+ const double* quadPts = _quadrature->quadPts();
+ const double* quadWts = _quadrature->quadWts();
+ const double* jacobianDet = _quadrature->jacobianDet();
+ const int numBasis = _quadrature->numCorners();
+ const int spaceDim = _quadrature->spaceDim();
+
+ // FIX THIS
+ // Hardwire mass density
+ const double density = 1.0;
+
+ // Compute lumped mass matrix for cell
+ double sumdiag = 0;
+ double diagScale = 0;
+ for (int iQuad=0; iQuad < numQuadPts; ++iQuad) {
+ const double wt = quadWts[iQuad] * jacobianDet[iQuad] * density;
+ for (int iBasis=0, iQ=iQuad*numBasis; iBasis < numBasis; ++iBasis) {
+ const int iBlock = iBasis * spaceDim;
+ const double val = wt*basis[iQ+iBasis]*basis[iQ+iBasis];
+ for (int iDim=0; iDim < spaceDim; ++iDim)
+ _cellVector[iBlock+iDim] += val;
+ sumdiag += val;
+ } // for
+ diagScale += numBasis*wt;
+ } // for
+ diagScale /= sumdiag*numBasis;
+
+ for (int iBasis=0; iBasis < numBasis; ++iBasis) {
+ const int iBlock = iBasis * spaceDim;
+ for (int iDim=0; iDim < numBasis; ++iDim)
+ _cellVector[iBlock+iDim] *= diagScale;
+ } // for
+ PetscErrorCode err =
+ PetscLogFlops(numQuadPts*(2+numBasis*2) + numBasis);
+ if (err)
+ throw std::runtime_error("Couldn't log PETSc flops.");
+
+ // Assemble cell contribution into field
+ fieldOut->updateAdd(patch, *cellIter, _cellVector);
+ } // for
+} // integrateLumped
+
+
+// End of file
Copied: short/3D/PyLith/trunk/libsrc/feassemble/IntegratorInertia.hh (from rev 5105, short/3D/PyLith/trunk/libsrc/feassemble/IntegratorInertia3D.hh)
===================================================================
--- short/3D/PyLith/trunk/libsrc/feassemble/IntegratorInertia3D.hh 2006-10-27 04:49:03 UTC (rev 5105)
+++ short/3D/PyLith/trunk/libsrc/feassemble/IntegratorInertia.hh 2006-10-27 04:52:45 UTC (rev 5106)
@@ -0,0 +1,94 @@
+// -*- C++ -*-
+//
+// ======================================================================
+//
+// Brad T. Aagaard
+// U.S. Geological Survey
+//
+// {LicenseText}
+//
+// ======================================================================
+//
+
+/**
+ * @file pylith/feassemble/IntegratorInertia.hh
+ *
+ * @brief Integrate inertial term for 3-D finite elements.
+ */
+
+#if !defined(pylith_feassemble_integratorinertia_hh)
+#define pylith_feassemble_integratorinertia_hh
+
+#include "Integrator.hh"
+
+namespace pylith {
+ namespace feassemble {
+ class IntegratorInertia;
+ class TestIntegratorInertia;
+ } // feassemble
+} // pylith
+
+class pylith::feassemble::IntegratorInertia : public Integrator
+{ // Integrator1D
+ friend class TestIntegratorInertia; // unit testing
+
+// PUBLIC MEMBERS ///////////////////////////////////////////////////////
+public :
+
+ /// Constructor
+ IntegratorInertia(void);
+
+ /// Destructor
+ ~IntegratorInertia(void);
+
+ /// Create a copy of this object.
+ Integrator* clone(void) const;
+
+ /** Integrate inertial term for 3-D finite elements.
+ *
+ * @param fieldOut Output field
+ * @param fieldIn Input field
+ * @param coordinates Field of cell vertex coordinates
+ */
+ void integrateAction(const ALE::Obj<ALE::Mesh::real_section_type>& fieldOut,
+ const ALE::Obj<ALE::Mesh::real_section_type>& fieldIn,
+ const ALE::Obj<ALE::Mesh::real_section_type>& coordinates);
+
+ /** Compute matrix associated with operator.
+ *
+ * @param mat Sparse matrix
+ * @param coordinates Field of cell vertex coordinates
+ */
+ void integrate(PetscMat* mat,
+ const ALE::Obj<real_section_type>& coordinates);
+
+ /** Compute lumped mass matrix.
+ *
+ * @param fieldOut Lumped mass matrix as field
+ * @param coordinates Field of cell vertex coordinates
+ */
+ void integrateLumped(const ALE::Obj<real_section_type>& fieldOut,
+ const ALE::Obj<real_section_type>& coordinates);
+
+// PROTECTED METHODS ////////////////////////////////////////////////////
+protected :
+
+ /** Copy constructor.
+ *
+ * @param i Integrator to copy
+ */
+ IntegratorInertia(const IntegratorInertia& i);
+
+// PRIVATE METHODS //////////////////////////////////////////////////////
+private :
+
+ /// Not implemented
+ const IntegratorInertia& operator=(const IntegratorInertia&);
+
+}; // IntegratorInertia
+
+#include "IntegratorInertia.icc" // inline methods
+
+#endif // pylith_feassemble_integratorineria_hh
+
+// End of file
Copied: short/3D/PyLith/trunk/libsrc/feassemble/IntegratorInertia.icc (from rev 5104, short/3D/PyLith/trunk/libsrc/feassemble/IntegratorInertia3D.icc)
===================================================================
--- short/3D/PyLith/trunk/libsrc/feassemble/IntegratorInertia3D.icc 2006-10-26 23:37:16 UTC (rev 5104)
+++ short/3D/PyLith/trunk/libsrc/feassemble/IntegratorInertia.icc 2006-10-27 04:52:45 UTC (rev 5106)
@@ -0,0 +1,26 @@
+// -*- C++ -*-
+//
+// ======================================================================
+//
+// Brad T. Aagaard
+// U.S. Geological Survey
+//
+// {LicenseText}
+//
+// ======================================================================
+//
+
+#if !defined(pylith_feassemble_integratorinertia_hh)
+#error "IntegratorInertia.icc must be included only from IntegratorInertia.hh"
+#else
+
+// Create a copy of this object.
+inline
+pylith::feassemble::Integrator*
+pylith::feassemble::IntegratorInertia::clone(void) const {
+ return new IntegratorInertia(*this);
+} // clone
+
+#endif
+
+// End of file
Deleted: short/3D/PyLith/trunk/libsrc/feassemble/IntegratorInertia3D.cc
===================================================================
--- short/3D/PyLith/trunk/libsrc/feassemble/IntegratorInertia3D.cc 2006-10-27 04:49:03 UTC (rev 5105)
+++ short/3D/PyLith/trunk/libsrc/feassemble/IntegratorInertia3D.cc 2006-10-27 04:52:45 UTC (rev 5106)
@@ -1,251 +0,0 @@
-// -*- C++ -*-
-//
-// ======================================================================
-//
-// Brad T. Aagaard
-// U.S. Geological Survey
-//
-// {LicenseText}
-//
-// ======================================================================
-//
-
-#include <portinfo>
-
-#include "IntegratorInertia3D.hh" // implementation of class methods
-
-#include "Quadrature.hh" // USES Quadrature
-
-#include "petscmat.h" // USES PetscMat
-
-#include <assert.h> // USES assert()
-#include <stdexcept> // USES std::runtime_error
-
-// ----------------------------------------------------------------------
-// Constructor
-pylith::feassemble::IntegratorInertia3D::IntegratorInertia3D(void)
-{ // constructor
-} // constructor
-
-// ----------------------------------------------------------------------
-// Destructor
-pylith::feassemble::IntegratorInertia3D::~IntegratorInertia3D(void)
-{ // destructor
-} // destructor
-
-// ----------------------------------------------------------------------
-// Copy constructor.
-pylith::feassemble::IntegratorInertia3D::IntegratorInertia3D(const IntegratorInertia3D& i) :
- Integrator(i)
-{ // copy constructor
-} // copy constructor
-
-// ----------------------------------------------------------------------
-// Integrate inertial term for 3-D finite elements.
-void
-pylith::feassemble::IntegratorInertia3D::integrateAction(
- const ALE::Obj<real_section_type>& fieldOut,
- const ALE::Obj<real_section_type>& fieldIn,
- const ALE::Obj<real_section_type>& coordinates)
-{ // integrateAction
- assert(0 != _quadrature);
-
- // Get information about section
- const topology_type::patch_type patch = 0;
- const ALE::Obj<topology_type>& topology = fieldIn->getTopology();
- const ALE::Obj<topology_type::label_sequence>& cells =
- topology->heightStratum(patch, 0);
- const topology_type::label_sequence::iterator cellsEnd = cells->end();
-
- // Allocate vector for cell values (if necessary)
- _initCellVector();
-
- for (topology_type::label_sequence::iterator cellIter=cells->begin();
- cellIter != cellsEnd;
- ++cellIter) {
- // Compute geometry information for current cell
- _quadrature->computeGeometry(coordinates, *cellIter);
-
- // Reset element vector to zero
- _resetCellVector();
-
- // Restrict input field to cell
- const real_section_type::value_type* fieldInCell =
- fieldIn->restrict(patch, *cellIter);
-
- // Get cell geometry information
- const int numQuadPts = _quadrature->numQuadPts();
- const double* basis = _quadrature->basis();
- const double* quadPts = _quadrature->quadPts();
- const double* quadWts = _quadrature->quadWts();
- const double* jacobianDet = _quadrature->jacobianDet();
- const int numBasis = _quadrature->numCorners();
- const int spaceDim = _quadrature->spaceDim();
-
- // FIX THIS
- // Hardwire mass density
- const double density = 1.0;
-
- // Compute action for cell
- for (int iQuad=0; iQuad < numQuadPts; ++iQuad) {
- const double wt = quadWts[iQuad] * jacobianDet[iQuad] * density;
- for (int iBasis=0, iQ=iQuad*numBasis; iBasis < numBasis; ++iBasis) {
- const int iBlock = iBasis * spaceDim;
- double val = wt*basis[iQ+iBasis];
- for (int jBasis=0; jBasis < numBasis; ++jBasis) {
- val *= basis[iQ+jBasis];
- for (int iDim=0; iDim < spaceDim; ++iDim)
- _cellVector[iBlock+iDim] += val * fieldInCell[iBlock+iDim];
- } // for
- } // for
- } // for
- PetscErrorCode err =
- PetscLogFlops(numQuadPts*(2+numBasis*(1+numBasis*(1+2*spaceDim))));
- if (err)
- throw std::runtime_error("Couldn't log PETSc flops.");
-
- // Assemble cell contribution into field
- fieldOut->updateAdd(patch, *cellIter, _cellVector);
- } // for
-} // integrateAction
-
-// ----------------------------------------------------------------------
-// Compute matrix associated with operator.
-void
-pylith::feassemble::IntegratorInertia3D::integrate(
- PetscMat* mat,
- const ALE::Obj<real_section_type>& coordinates)
-{ // integrate
- assert(0 != mat);
- assert(0 != _quadrature);
-
- // Get information about section
- const topology_type::patch_type patch = 0;
- const ALE::Obj<topology_type>& topology = coordinates->getTopology();
- const ALE::Obj<topology_type::label_sequence>& cells =
- topology->heightStratum(patch, 0);
- const topology_type::label_sequence::iterator cellsEnd = cells->end();
-
- // Allocate matrix for cell values (if necessary)
- _initCellMatrix();
-
- for (topology_type::label_sequence::iterator cellIter=cells->begin();
- cellIter != cellsEnd;
- ++cellIter) {
- // Compute geometry information for current cell
- _quadrature->computeGeometry(coordinates, *cellIter);
-
- // Reset element matrix to zero
- _resetCellMatrix();
-
- // Get cell geometry information
- const int numQuadPts = _quadrature->numQuadPts();
- const double* basis = _quadrature->basis();
- const double* quadPts = _quadrature->quadPts();
- const double* quadWts = _quadrature->quadWts();
- const double* jacobianDet = _quadrature->jacobianDet();
- const int numBasis = _quadrature->numCorners();
- const int spaceDim = _quadrature->spaceDim();
-
- // FIX THIS
- // Hardwire mass density
- const double density = 1.0;
-
- // Integrate cell
- for (int iQuad=0; iQuad < numQuadPts; ++iQuad) {
- const double wt = quadWts[iQuad] * jacobianDet[iQuad] * density;
- for (int iBasis=0, iQ=iQuad*numBasis; iBasis < numBasis; ++iBasis) {
- const int iBlock = iBasis * spaceDim;
- double val = wt*basis[iQ+iBasis];
- for (int jBasis=0; jBasis < numBasis; ++jBasis) {
- const int jBlock = jBasis * spaceDim;
- val *= basis[iQ+jBasis];
- for (int iDim=0; iDim < spaceDim; ++iDim)
- _cellMatrix[(iBlock+iDim)*(numBasis*spaceDim)+jBlock+iDim] += val;
- } // for
- } // for
- } // for
- PetscErrorCode err =
- PetscLogFlops(numQuadPts*(2+numBasis*(1+numBasis*(1+spaceDim))));
- if (err)
- throw std::runtime_error("Couldn't log PETSc flops.");
-
- // Assemble cell contribution into sparse matrix
- // STUFF GOES HERE
- } // for
-} // integrate
-
-// ----------------------------------------------------------------------
-// Compute lumped matrix associated with operator.
-void
-pylith::feassemble::IntegratorInertia3D::integrateLumped(
- const ALE::Obj<real_section_type>& fieldOut,
- const ALE::Obj<real_section_type>& coordinates)
-{ // integrateLumped
- assert(0 != _quadrature);
-
- // Get information about section
- const topology_type::patch_type patch = 0;
- const ALE::Obj<topology_type>& topology = coordinates->getTopology();
- const ALE::Obj<topology_type::label_sequence>& cells =
- topology->heightStratum(patch, 0);
- const topology_type::label_sequence::iterator cellsEnd = cells->end();
-
- // Allocate matrix for cell values (if necessary)
- _initCellVector();
-
- for (topology_type::label_sequence::iterator cellIter=cells->begin();
- cellIter != cellsEnd;
- ++cellIter) {
- // Compute geometry information for current cell
- _quadrature->computeGeometry(coordinates, *cellIter);
-
- // Reset element matrix to zero
- _resetCellVector();
-
- // Get cell geometry information
- const int numQuadPts = _quadrature->numQuadPts();
- const double* basis = _quadrature->basis();
- const double* quadPts = _quadrature->quadPts();
- const double* quadWts = _quadrature->quadWts();
- const double* jacobianDet = _quadrature->jacobianDet();
- const int numBasis = _quadrature->numCorners();
- const int spaceDim = _quadrature->spaceDim();
-
- // FIX THIS
- // Hardwire mass density
- const double density = 1.0;
-
- // Compute lumped mass matrix for cell
- double sumdiag = 0;
- double diagScale = 0;
- for (int iQuad=0; iQuad < numQuadPts; ++iQuad) {
- const double wt = quadWts[iQuad] * jacobianDet[iQuad] * density;
- for (int iBasis=0, iQ=iQuad*numBasis; iBasis < numBasis; ++iBasis) {
- const int iBlock = iBasis * spaceDim;
- const double val = wt*basis[iQ+iBasis]*basis[iQ+iBasis];
- for (int iDim=0; iDim < spaceDim; ++iDim)
- _cellVector[iBlock+iDim] += val;
- sumdiag += val;
- } // for
- diagScale += numBasis*wt;
- } // for
- diagScale /= sumdiag*numBasis;
-
- for (int iBasis=0; iBasis < numBasis; ++iBasis) {
- const int iBlock = iBasis * spaceDim;
- for (int iDim=0; iDim < numBasis; ++iDim)
- _cellVector[iBlock+iDim] *= diagScale;
- } // for
- PetscErrorCode err =
- PetscLogFlops(numQuadPts*(2+numBasis*2) + numBasis);
- if (err)
- throw std::runtime_error("Couldn't log PETSc flops.");
-
- // Assemble cell contribution into field
- fieldOut->updateAdd(patch, *cellIter, _cellVector);
- } // for
-} // integrateLumped
-
-
-// End of file
Deleted: short/3D/PyLith/trunk/libsrc/feassemble/IntegratorInertia3D.hh
===================================================================
--- short/3D/PyLith/trunk/libsrc/feassemble/IntegratorInertia3D.hh 2006-10-27 04:49:03 UTC (rev 5105)
+++ short/3D/PyLith/trunk/libsrc/feassemble/IntegratorInertia3D.hh 2006-10-27 04:52:45 UTC (rev 5106)
@@ -1,94 +0,0 @@
-// -*- C++ -*-
-//
-// ======================================================================
-//
-// Brad T. Aagaard
-// U.S. Geological Survey
-//
-// {LicenseText}
-//
-// ======================================================================
-//
-
-/**
- * @file pylith/feassemble/IntegratorInertia3D.hh
- *
- * @brief Integrate inertial term for 3-D finite elements.
- */
-
-#if !defined(pylith_feassemble_integratorinertia3d_hh)
-#define pylith_feassemble_integratorinertia3d_hh
-
-#include "Integrator.hh"
-
-namespace pylith {
- namespace feassemble {
- class IntegratorInertia3D;
- class TestIntegratorInertia3D;
- } // feassemble
-} // pylith
-
-class pylith::feassemble::IntegratorInertia3D : public Integrator
-{ // Integrator1D
- friend class TestIntegratorInertia3D; // unit testing
-
-// PUBLIC MEMBERS ///////////////////////////////////////////////////////
-public :
-
- /// Constructor
- IntegratorInertia3D(void);
-
- /// Destructor
- ~IntegratorInertia3D(void);
-
- /// Create a copy of this object.
- Integrator* clone(void) const;
-
- /** Integrate inertial term for 3-D finite elements.
- *
- * @param fieldOut Output field
- * @param fieldIn Input field
- * @param coordinates Field of cell vertex coordinates
- */
- void integrateAction(const ALE::Obj<ALE::Mesh::real_section_type>& fieldOut,
- const ALE::Obj<ALE::Mesh::real_section_type>& fieldIn,
- const ALE::Obj<ALE::Mesh::real_section_type>& coordinates);
-
- /** Compute matrix associated with operator.
- *
- * @param mat Sparse matrix
- * @param coordinates Field of cell vertex coordinates
- */
- void integrate(PetscMat* mat,
- const ALE::Obj<real_section_type>& coordinates);
-
- /** Compute lumped mass matrix.
- *
- * @param fieldOut Lumped mass matrix as field
- * @param coordinates Field of cell vertex coordinates
- */
- void integrateLumped(const ALE::Obj<real_section_type>& fieldOut,
- const ALE::Obj<real_section_type>& coordinates);
-
-// PROTECTED METHODS ////////////////////////////////////////////////////
-protected :
-
- /** Copy constructor.
- *
- * @param i Integrator to copy
- */
- IntegratorInertia3D(const IntegratorInertia3D& i);
-
-// PRIVATE METHODS //////////////////////////////////////////////////////
-private :
-
- /// Not implemented
- const IntegratorInertia3D& operator=(const IntegratorInertia3D&);
-
-}; // IntegratorInertia3D
-
-#include "IntegratorInertia3D.icc" // inline methods
-
-#endif // pylith_feassemble_integratorineria3d_hh
-
-// End of file
Deleted: short/3D/PyLith/trunk/libsrc/feassemble/IntegratorInertia3D.icc
===================================================================
--- short/3D/PyLith/trunk/libsrc/feassemble/IntegratorInertia3D.icc 2006-10-27 04:49:03 UTC (rev 5105)
+++ short/3D/PyLith/trunk/libsrc/feassemble/IntegratorInertia3D.icc 2006-10-27 04:52:45 UTC (rev 5106)
@@ -1,26 +0,0 @@
-// -*- C++ -*-
-//
-// ======================================================================
-//
-// Brad T. Aagaard
-// U.S. Geological Survey
-//
-// {LicenseText}
-//
-// ======================================================================
-//
-
-#if !defined(pylith_feassemble_integratorinertia3d_hh)
-#error "IntegratorInertia3D.icc must be included only from IntegratorInertia3D.hh"
-#else
-
-// Create a copy of this object.
-inline
-pylith::feassemble::Integrator*
-pylith::feassemble::IntegratorInertia3D::clone(void) const {
- return new IntegratorInertia3D(*this);
-} // clone
-
-#endif
-
-// End of file
Modified: short/3D/PyLith/trunk/libsrc/feassemble/Makefile.am
===================================================================
--- short/3D/PyLith/trunk/libsrc/feassemble/Makefile.am 2006-10-27 04:49:03 UTC (rev 5105)
+++ short/3D/PyLith/trunk/libsrc/feassemble/Makefile.am 2006-10-27 04:52:45 UTC (rev 5106)
@@ -18,7 +18,7 @@
libpylithfeassemble_la_SOURCES = \
Integrator.cc \
IntegratorElasticity3D.cc \
- IntegratorInertia3D.cc \
+ IntegratorInertia.cc \
Quadrature.cc \
Quadrature1D.cc \
Quadrature1Din2D.cc \
@@ -31,8 +31,8 @@
Integrator.hh \
IntegratorElasticity3D.hh \
IntegratorElasticity3D.icc \
- IntegratorInertia3D.hh \
- IntegratorInertia3D.icc \
+ IntegratorInertia.hh \
+ IntegratorInertia.icc \
Quadrature.hh \
Quadrature.icc \
Quadrature1D.hh \
More information about the cig-commits
mailing list