Some refactoring: Moved instanciation of 2D patches to static create method in ASM2D-scope and Discretization enum to ASM-scope

git-svn-id: http://svn.sintef.no/trondheim/IFEM/trunk@1280 e10b68d5-8a6e-419e-a041-bce267b0401d

Apps/HDF5toVTx/HDF5toVTx.C

@@ -14,56 +14,47 @@
 #include "HDF5Writer.h"
 #include "XMLWriter.h"
 #include "StringUtils.h"
-#include <sstream>
-#include <stdlib.h>
-#include <string.h>
 #include "ASMs1D.h"
 #include "ASMs2D.h"
 #include "ASMs3D.h"
-#if HAS_LRSPLINE == 1
-#include "LR/ASMu2D.h"
-#endif
 #include "ElementBlock.h"
 #include "VTF.h"
 #include "VTU.h"
-
+#include <sstream>
+#include <stdlib.h>
+#include <string.h>
 #include <assert.h>
 
 typedef std::map< std::string,std::vector<XMLWriter::Entry> > ProcessList;
-typedef std::map< std::string, std::vector<int> > VTFList;
+typedef std::map< std::string,std::vector<int> > VTFList;
 
-std::vector<ASMbase*> readBasis(const std::string& name, 
+
-                                int patches, HDF5Writer& hdf, int dim, int level)
+std::vector<ASMbase*> readBasis (const std::string& name, 
+				 int patches, HDF5Writer& hdf,
+				 int dim, int level)
 {
+  unsigned char nf = 1;
+  ASM::Discretization ptype;
   std::vector<ASMbase*> result;
   for (int i=0;i<patches;++i) {
-    std::stringstream geom;
+    std::stringstream geom, basis;
     geom << '/' << level << "/basis/";
     geom << name;
     geom << "/";
     geom << i+1;
     std::string out;
     hdf.readString(geom.str(),out);
-    std::stringstream basis;
+    ptype = out.substr(0,10) == "# LRSPLINE" ? ASM::LRSpline : ASM::Spline;
     basis << out;
-    if (out.substr(0,10) == "# LRSPLINE") {
+    if (dim == 1)
-      switch (dim) {
+      result.push_back(new ASMs1D(basis,1,1));
-#if HAS_LRSPLINE == 1
+    else if (dim == 2) {
-        case 2:
+      result.push_back(ASM2D::create(ptype,&nf));
-          result.push_back(new ASMu2D(basis,2,1));
+      assert(result.back());
-          break;
+      result.back()->read(basis);
-#endif
-        default:
-          assert(0);
-      }
-    } else {
-      if (dim == 1)
-        result.push_back(new ASMs1D(basis,1,1));
-      if (dim == 2)
-        result.push_back(new ASMs2D(basis,2,1));
-      if (dim == 3)
-        result.push_back(new ASMs3D(basis,false,1));
     }
+    else if (dim == 3)
+      result.push_back(new ASMs3D(basis,false,1));
     result.back()->generateFEMTopology();
   }
 
@@ -157,26 +148,21 @@ void writePatchGeometry(ASMbase* patch, int id, VTF& myVtf, int* nViz, int block
 }
 
 
-#define TRY(x,y) { x* t = dynamic_cast<x*>(y); \
-                   if (t) \
-                     t->getGridParameters(gpar[k],k,n[k]-1); \
-                 }
 std::vector<RealArray*> generateFEModel(std::vector<ASMbase*> patches,
                                         int dims, int* n)
 {
   std::vector<RealArray*> result;
+  result.reserve(patches.size());
   for (size_t i=0;i<patches.size();++i) {
     RealArray* gpar = new RealArray[dims];
     for (int k=0;k<dims;++k) {
       if (dims == 2) {
-#if HAS_LRSPLINE == 1
+        ASM2D* patch = dynamic_cast<ASM2D*>(patches[i]);
-        TRY(ASMu2D,patches[i])
+        if (patch) patch->getGridParameters(gpar[k],k,n[k]-1);
-#endif
-        TRY(ASMs2D,patches[i])
       }
-      if (dims == 3) {
+      else if (dims == 3) {
-        ASMs3D* patch = (ASMs3D*)patches[i];
+        ASMs3D* patch = dynamic_cast<ASMs3D*>(patches[i]);
-        patch->getGridParameters(gpar[k],k,n[k]-1);
+        if (patch) patch->getGridParameters(gpar[k],k,n[k]-1);
       }
     }
     result.push_back(gpar);

src/SIM/SIM1D.C

@@ -40,10 +40,10 @@ bool SIM1D::parse (char* keyWord, std::istream& is)
     for (int i = 0; i < npatch && (cline = utl::readLine(is)); i++)
     {
       switch (discretization) {
-      case Lagrange:
+      case ASM::Lagrange:
 	pch = new ASMs1DLag(strtok(cline," "),1,nf);
 	break;
-      case Spectral:
+      case ASM::Spectral:
 	pch = new ASMs1DSpec(strtok(cline," "),1,nf);
 	break;
       default:
@@ -380,16 +380,17 @@ void SIM1D::readPatches (std::istream& isp)
   for (int patchNo = 1; isp.good(); patchNo++)
   {
     std::cout <<"Reading patch "<< patchNo << std::endl;
-    switch (discretization) {
+    switch (discretization)
-      case Lagrange:
+      {
+      case ASM::Lagrange:
         pch = new ASMs1DLag(isp,1,nf);
         break;
-      case Spectral:
+      case ASM::Spectral:
         pch = new ASMs1DSpec(isp,1,nf);
         break;
       default:
         pch = new ASMs1D(isp,1,nf);
-    }
+      }
     if (pch->empty())
       delete pch;
     else

src/SIM/SIM2D.C

@@ -12,12 +12,7 @@
 //==============================================================================
 
 #include "SIM2D.h"
-#include "ASMs2Dmx.h"
+#include "ASMs2D.h"
-#include "ASMs2DmxLag.h"
-#include "ASMs2DSpec.h"
-#ifdef HAS_LRSPLINE
-#include "LR/ASMu2D.h"
-#endif
 #include "Functions.h"
 #include "Utilities.h"
 #include <fstream>
@@ -115,14 +110,12 @@ bool SIM2D::parse (char* keyWord, std::istream& is)
       }
 
       if (isn.good() && pid > 0)
-      {
         if (!static_cast<ASMs2D*>(myModel[pid-1])->assignNodeNumbers(n))
         {
           std::cerr <<" *** SIM2D::parse: Failed to assign node numbers"
                     <<" for patch "<< patch << std::endl;
           return false;
         }
-      }
     }
   }
 
@@ -386,6 +379,7 @@ bool SIM2D::parse (char* keyWord, std::istream& is)
     if (ignoreDirichlet) return true; // Ignore all boundary conditions
     if (!this->createFEMmodel()) return false;
 
+    ASM2D* pch = 0;
     int nfix = atoi(keyWord+9);
     std::cout <<"\nNumber of fixed points: "<< nfix << std::endl;
     for (int i = 0; i < nfix && (cline = utl::readLine(is)); i++)
@@ -393,15 +387,16 @@ bool SIM2D::parse (char* keyWord, std::istream& is)
       int patch = atoi(strtok(cline," "));
       double rx = atof(strtok(NULL," "));
       double ry = atof(strtok(NULL," "));
-      int bcode = (cline = strtok(NULL," ")) ? atoi(cline) : 123;
+      int bcode = (cline = strtok(NULL," ")) ? atoi(cline) : 12;
 
       int pid = this->getLocalPatchIndex(patch);
-      if (pid < 1) continue;
+      if (pid > 0 && (pch = dynamic_cast<ASM2D*>(myModel[pid-1])))
-
+      {
-      std::cout <<"\tConstraining P"<< patch
+	std::cout <<"\tConstraining P"<< patch
-		<<" point at "<< rx <<" "<< ry
+		  <<" point at "<< rx <<" "<< ry
-		<<" with code "<< bcode << std::endl;
+		  <<" with code "<< bcode << std::endl;
-      static_cast<ASMs2D*>(myModel[pid-1])->constrainNode(rx,ry,bcode);
+	pch->constrainNode(rx,ry,bcode);
+      }
     }
   }
 
@@ -479,69 +474,54 @@ void SIM2D::setQuadratureRule (size_t ng)
 }
 
 
-void SIM2D::readPatch (const char* patchFile, int pchInd)
+bool SIM2D::readPatch (const char* patchFile, int pchInd)
 {
-  ASMbase* pch = 0;
+  std::ifstream is(patchFile);
-  switch (discretization) {
+  if (!is.good())
-  case Lagrange:
-    if (nf[1] > 0)
-      pch = new ASMs2DmxLag(patchFile,2,nf[0],nf[1]);
-    else
-      pch = new ASMs2DLag(patchFile,2,nf[0]);
-    break;
-  case Spectral:
-    pch = new ASMs2DSpec(patchFile,2,nf[0]);
-    break;
-#ifdef HAS_LRSPLINE
-  case LRSpline:
-    pch = new ASMu2D(patchFile,2,nf[0]);
-    break;
-#endif
-  default:
-    if (nf[1] > 0)
-      pch = new ASMs2Dmx(patchFile,2,nf[0],nf[1]);
-    else
-      pch = new ASMs2D(patchFile,2,nf[0]);
-  }
-  if (pch->empty() || this->getLocalPatchIndex(pchInd+1) < 1)
-    delete pch;
-  else
-    myModel.push_back(pch);
-}
-
-
-void SIM2D::readPatches (std::istream& isp)
-{
-  ASMbase* pch = 0;
-  for (int patchNo = 1; isp.good(); patchNo++)
   {
-    std::cout <<"Reading patch "<< patchNo << std::endl;
+    std::cerr <<" *** SIM2D: Failure opening patch file"
-    switch (discretization) {
+	      << patchFile << std::endl;
-    case Lagrange:
+    return false;
-      if (nf[1] > 0)
+  }
-	pch = new ASMs2DmxLag(isp,2,nf[0],nf[1]);
+
-      else
+  ASMbase* pch = ASM2D::create(discretization,nf);
-	pch = new ASMs2DLag(isp,2,nf[0]);
+  if (pch)
-      break;
+  {
-    case Spectral:
+    std::cout <<"\nReading patch file "<< patchFile << std::endl;
-      pch = new ASMs2DSpec(isp,2,nf[0]);
+    if (!pch->read(is))
-      break;
+    {
-    case LRSpline:
+      delete pch;
-#ifdef HAS_LRSPLINE
+      return false;
-      pch = new ASMu2D(isp,2,nf[0]);
-      break;
-#endif
-    default:
-      if (nf[1] > 0)
-	pch = new ASMs2Dmx(isp,2,nf[0],nf[1]);
-      else
-	pch = new ASMs2D(isp,2,nf[0]);
     }
-    if (pch->empty() || this->getLocalPatchIndex(patchNo) < 1)
+    else if (pch->empty() || this->getLocalPatchIndex(pchInd+1) < 1)
       delete pch;
     else
       myModel.push_back(pch);
   }
+
+  return true;
+}
+
+
+bool SIM2D::readPatches (std::istream& isp)
+{
+  ASMbase* pch = 0;
+  for (int pchInd = 1; isp.good(); pchInd++)
+    if ((pch = ASM2D::create(discretization,nf)))
+    {
+      std::cout <<"Reading patch "<< pchInd << std::endl;
+      if (!pch->read(isp))
+      {
+	delete pch;
+	return false;
+      }
+      else if (pch->empty() || this->getLocalPatchIndex(pchInd) < 1)
+        delete pch;
+      else
+        myModel.push_back(pch);
+    }
+
+  return true;
 }
 
 
@@ -560,11 +540,10 @@ void SIM2D::clonePatches (const FEModelVec& patches,
 bool SIM2D::refine (const std::vector<int>& elements,
 		    const std::vector<int>& options, const char* fName)
 {
+  ASM2D* pch = 0;
   for (size_t i = 0; i < myModel.size(); i++)
-    if (!myModel.empty())
+    if (!myModel[i]->empty() && (pch = dynamic_cast<ASM2D*>(myModel[i])))
-#ifdef HAS_LRSPLINE
+      if (!pch->refine(elements,options,fName))
-      if (!static_cast<ASMu2D*>(myModel[i])->refine(elements,options,fName))
-#endif
 	return false;
 
   isRefined = true;

src/SIM/SIM2D.h

@@ -51,10 +51,10 @@ protected:
   //! \brief Reads a patch from given input file.
   //! \param[in] patchFile Name of file to read from
   //! \param[in] pchInd 0-based index of the patch to read
-  void readPatch(const char* patchFile, int pchInd);
+  bool readPatch(const char* patchFile, int pchInd);
   //! \brief Reads patches from given input stream.
   //! \param[in] isp The file stream to read from
-  void readPatches(std::istream& isp);
+  bool readPatches(std::istream& isp);
 
   //! \brief Refines a list of elements.
   //! \param[in] elements 1-based indices of the elements to refine

src/SIM/SIM3D.C

@@ -43,13 +43,13 @@ bool SIM3D::parse (char* keyWord, std::istream& is)
     {
       cline = strtok(cline," ");
       switch (discretization) {
-      case Lagrange:
+      case ASM::Lagrange:
 	if (nf[1] > 0)
 	  pch = new ASMs3DmxLag(cline,checkRHSys,nf[0],nf[1]);
 	else
 	  pch = new ASMs3DLag(cline,checkRHSys,nf[0]);
 	break;
-      case Spectral:
+      case ASM::Spectral:
 	pch = new ASMs3DSpec(cline,checkRHSys,nf[0]);
 	break;
       default:
@@ -598,14 +598,15 @@ void SIM3D::readPatches (std::istream& isp)
   for (int patchNo = 1; isp.good(); patchNo++)
   {
     std::cout <<"Reading patch "<< patchNo << std::endl;
-    switch (discretization) {
+    switch (discretization)
-      case Lagrange:
+      {
+      case ASM::Lagrange:
         if (nf[1] > 0)
           pch = new ASMs3DmxLag(isp,checkRHSys,nf[0],nf[1]);
         else
           pch = new ASMs3DLag(isp,checkRHSys,nf[0]);
         break;
-      case Spectral:
+      case ASM::Spectral:
         pch = new ASMs3DSpec(isp,checkRHSys,nf[0]);
         break;
       default:
@@ -613,7 +614,7 @@ void SIM3D::readPatches (std::istream& isp)
           pch = new ASMs3Dmx(isp,checkRHSys,nf[0],nf[1]);
         else
           pch = new ASMs3D(isp,checkRHSys,nf[0]);
-    }
+      }
     if (pch->empty() || this->getLocalPatchIndex(patchNo) < 1)
       delete pch;
     else

src/SIM/SIMbase.C

@@ -37,10 +37,10 @@
 #include <stdio.h>
 
 
-SIMbase::Discretization SIMbase::discretization  = SIMbase::Spline;
+ASM::Discretization SIMbase::discretization  = ASM::Spline;
-bool                    SIMbase::preserveNOrder  = false;
+bool                SIMbase::preserveNOrder  = false;
-bool                    SIMbase::ignoreDirichlet = false;
+bool                SIMbase::ignoreDirichlet = false;
-int                     SIMbase::num_threads_SLU = 1;
+int                 SIMbase::num_threads_SLU = 1;
 
 
 SIMbase::SIMbase () : g2l(&myGlb2Loc)
@@ -283,7 +283,7 @@ bool SIMbase::preprocess (const std::vector<int>& ignored, bool fixDup)
   // will map the global node numbers to local node numbers on the current
   // processor. In serial simulations, the global-to-local mapping will be unity
   // unless the original global node number sequence had "holes" due to
-  // duplicated nodes and/or erast patches.
+  // duplicated nodes and/or erased patches.
   int ngnod = 0;
   int renum = 0;
   if (preserveNOrder)
@@ -1225,7 +1225,7 @@ bool SIMbase::writeGlvS (const Vector& psol, const int* nViz,
       else
 	sID[k++].push_back(nBlock);
 
-    if (discretization == Spline)
+    if (discretization == ASM::Spline || discretization == ASM::SplineC1)
     {
       // 3. Projection of secondary solution variables (tensorial splines only)
 
@@ -1305,7 +1305,7 @@ bool SIMbase::writeGlvS (const Vector& psol, const int* nViz,
     if (!myVtf->writeSblk(sID[j],myProblem->getField2Name(i,haveAsol?"FE":0),
 			  ++idBlock,iStep)) return false;
 
-  if (discretization == Spline)
+  if (discretization == ASM::Spline || discretization == ASM::SplineC1)
     for (i = 0; i < nf && !sID[j].empty(); i++, j++)
       if (!myVtf->writeSblk(sID[j],myProblem->getField2Name(i,"Projected"),
 			    ++idBlock,iStep)) return false;
@@ -1610,7 +1610,7 @@ bool SIMbase::dumpResults (const Vector& psol, double time, std::ostream& os,
     // Find all evaluation points within this patch, if any
     for (j = 0, p = myPoints.begin(); p != myPoints.end(); j++, p++)
       if (this->getLocalPatchIndex(p->patch) == (int)(i+1))
-	if (discretization == Spline || discretization == LRSpline)
+	if (discretization >= ASM::Spline)
 	{
 	  points.push_back(p->inod > 0 ? p->inod : -(j+1));
 	  for (k = 0; k < myModel[i]->getNoParamDim(); k++)
@@ -1622,7 +1622,7 @@ bool SIMbase::dumpResults (const Vector& psol, double time, std::ostream& os,
     if (points.empty()) continue; // no points in this patch
 
     myModel[i]->extractNodeVec(psol,myProblem->getSolution());
-    if (discretization == Spline || discretization == LRSpline)
+    if (discretization >= ASM::Spline)
     {
       // Evaluate the primary solution variables
       if (!myModel[i]->evalSolution(sol1,myProblem->getSolution(),params,false))
@@ -1658,7 +1658,7 @@ bool SIMbase::dumpResults (const Vector& psol, double time, std::ostream& os,
       for (k = 1; k <= sol1.rows(); k++)
 	os << std::setw(flWidth) << utl::trunc(sol1(k,j+1));
 
-      if (discretization == Spline || discretization == LRSpline)
+      if (discretization >= ASM::Spline)
       {
         if (formatted && sol2.rows() > 0)
           os <<"\n\t\tsol2 =";

src/SIM/SIMbase.h

@@ -17,6 +17,7 @@
 #include "SIMinput.h"
 #include "SystemMatrix.h"
 #include "TimeDomain.h"
+#include "ASMenums.h"
 #include "Property.h"
 #include "Function.h"
 #include "Vec3.h"
@@ -543,10 +544,7 @@ protected:
   virtual double externalEnergy(const Vectors& psol) const;
 
 public:
-  //! \brief Enum defining the available discretization methods.
+  static ASM::Discretization discretization; //!< Spatial discretization option
-  enum Discretization { Spline, Lagrange, Spectral, LRSpline };
-
-  static Discretization discretization; //!< Spatial discretization option
 
   static bool ignoreDirichlet; //!< Set to \e true for free vibration analysis
   static bool preserveNOrder;  //!< Set to \e true to preserve node ordering