Changed: Can now use dumpMatlabGrid without SAM for single-patch

src/ASM/Test/TestMatlabPatch.C

@@ -12,7 +12,6 @@
 
 #include "SIM2D.h"
 #include "ASMs2D.h"
-#include "IntegrandBase.h"
 #include "ModelGenerator.h"
 #include "tinyxml.h"
 #include <fstream>
@@ -20,20 +19,17 @@
 #include "gtest/gtest.h"
 
 
-class Dummy : public IntegrandBase {};
-
-
 class SIM2D_default : public SIM2D
 {
 public:
-  SIM2D_default(int nu, int nv) : SIM2D(new Dummy())
+  SIM2D_default(int nu, int nv)
   {
     // Create a bi-unit square grid with (nu+1)x(nv+1) linear elements
     this->createDefaultModel();
     ASMs2D* pch1 = static_cast<ASMs2D*>(myModel.front());
     pch1->uniformRefine(0,nu);
     pch1->uniformRefine(1,nv);
-    EXPECT_TRUE(this->preprocess());
+    EXPECT_TRUE(this->createFEMmodel());
 
     // Create topological boundary entities (vertices and edges)
     TiXmlDocument doc;
@@ -48,7 +44,7 @@ public:
 class SIM2D_matlab : public SIM2D
 {
 public:
-  SIM2D_matlab(const char* inputFile) : SIM2D(new Dummy())
+  SIM2D_matlab(const char* inputFile)
   {
     std::string xml("<geometry><patchfile type=\"matlab\">");
     xml.append(inputFile);
@@ -60,7 +56,7 @@ public:
     const TiXmlElement* tag = doc.RootElement();
     EXPECT_TRUE(tag != nullptr);
     EXPECT_TRUE(this->parse(tag));
-    EXPECT_TRUE(this->preprocess());
+    EXPECT_TRUE(this->createFEMmodel());
   }
   virtual ~SIM2D_matlab() {}
 };

src/SIM/SIMbase.C

@@ -1047,7 +1047,7 @@ bool SIMbase::solveSystem (Vector& solution, int printSol, double* rCond,
 
   // Expand solution vector from equation ordering to DOF-ordering
   if (status)
-    status = mySam->expandSolution(*b,solution);
+    status = mySam ? mySam->expandSolution(*b,solution) : false;
 
 #if SP_DEBUG > 2
   if (printSol < 1000) printSol = 1000;
@@ -1557,7 +1557,7 @@ bool SIMbase::systemModes (std::vector<Mode>& solution,
 			   int nev, int ncv, int iop, double shift,
 			   size_t iA, size_t iB)
 {
-  if (nev < 1 || ncv <= nev || !myEqSys) return false;
+  if (nev < 1 || ncv <= nev || !myEqSys || !mySam) return false;
 
   PROFILE1("Eigenvalue analysis");
 
@@ -1818,7 +1818,7 @@ bool SIMbase::extractPatchSolution (IntegrandBase* problem,
                                     const Vectors& sol, size_t pindx) const
 {
   ASMbase* pch = this->getPatch(pindx+1);
-  if (!pch) return false;
+  if (!pch || !mySam) return false;
 
   problem->initNodeMap(pch->getGlobalNodeNums());
   for (size_t i = 0; i < sol.size() && i < problem->getNoSolutions(); i++)
@@ -1853,21 +1853,21 @@ bool SIMbase::projectAnaSol (Vector& ssol,
   FunctionBase* f = mySol->getScalarSecSol();
   if (f)
   {
-    ssol.resize(f->dim()*mySam->getNoNodes());
+    ssol.resize(f->dim()*this->getNoNodes());
     return this->project(ssol,f,0,0,0,pMethod);
   }
 
   f = mySol->getVectorSecSol();
   if (f)
   {
-    ssol.resize(f->dim()*mySam->getNoNodes());
+    ssol.resize(f->dim()*this->getNoNodes());
     return this->project(ssol,f,0,0,0,pMethod);
   }
 
   f = mySol->getStressSol();
   if (f)
   {
-    ssol.resize(f->dim()*mySam->getNoNodes());
+    ssol.resize(f->dim()*this->getNoNodes());
     return this->project(ssol,f,0,0,0,pMethod);
   }
 
@@ -1943,14 +1943,13 @@ bool SIMbase::setPatchMaterial (size_t patch)
 
 bool SIMbase::addMADOF (unsigned char basis, unsigned char nndof, bool other)
 {
-  if (!mySam) return false;
-
   int key = basis << 16 + nndof;
   if (mixedMADOFs.find(key) != mixedMADOFs.end())
     return false; // This MADOF already calculated
 
   IntVec& madof = mixedMADOFs[key];
-  madof.resize(mySam->getNoNodes()+1,0);
+  madof.resize(this->getNoNodes()+1,0);
+  if (madof.size() < 2) return false;
 
   char nType = basis <= 1 ? 'D' : 'P' + basis-2;
   for (size_t i = 0; i < myModel.size(); i++)

src/SIM/SIMoutput.C

@@ -1154,7 +1154,7 @@ bool SIMoutput::dumpMatlabGrid (std::ostream& os, const std::string& name,
                                 const std::vector<std::string>& sets,
                                 double scale) const
 {
-  if (!mySam || this->getNoParamDim() != 2)
+  if (this->getNoParamDim() != 2)
   {
     std::cerr <<" *** SIMoutput::dumpMatlabGrid: For 2D only."<< std::endl;
     return false;
@@ -1195,7 +1195,7 @@ bool SIMoutput::dumpMatlabGrid (std::ostream& os, const std::string& name,
 
   // Write out all nodes
   os <<"\n  Node=[1";
-  size_t nnod = mySam->getNoNodes();
+  size_t nnod = this->getNoNodes();
   for (size_t inod = 1; inod <= nnod; inod++)
   {
     Vec3 X = this->getNodeCoord(inod);
@@ -1207,11 +1207,11 @@ bool SIMoutput::dumpMatlabGrid (std::ostream& os, const std::string& name,
 
   // Write out all elements
   os <<"\n  Element=[1";
-  size_t nen = 0, nel = mySam->getNoElms();
+  size_t nen = 0, nel = this->getNoElms(true);
   for (size_t iel = 1; iel <= nel; iel++)
   {
     IntVec nodes;
-    if (!mySam->getElmNodes(nodes,iel))
+    if (!this->getElmNodes(nodes,iel))
       return false;
     else if (nodes.size() == 4)
       std::swap(nodes[2],nodes[3]);
@@ -1246,32 +1246,34 @@ bool SIMoutput::dumpGeometry (std::ostream& os) const
 void SIMoutput::dumpPrimSol (const Vector& psol, utl::LogStream& os,
                              bool withID) const
 {
-  if (psol.empty()) return;
+  if (psol.empty() || !mySam)
+    return;
 
-  size_t i, j, ip;
+  for (const ASMbase* pch : myModel)
-  unsigned char k, n;
-  for (i = 0; i < myModel.size(); i++)
   {
-    if (myModel[i]->empty()) continue; // skip empty patches
+    if (pch->empty()) continue; // skip empty patches
 
     Vector patchSol;
-    myModel[i]->extractNodalVec(psol,patchSol,mySam->getMADOF());
+    pch->extractNodalVec(psol,patchSol,mySam->getMADOF());
 
     if (withID)
     {
       if (myModel.size() > 1)
-        os <<"\n# Patch: "<< i+1;
+        os <<"\n# Patch: "<< pch->idx+1;
       os <<"\n# inod/gnod\tNodal Coordinates\tSolution\n";
     }
-    for (ip = 0, j = 1; j <= myModel[i]->getNoNodes(); j++)
+
+    size_t ip = 0, nnod = pch->getNoNodes();
+    for (size_t j = 1; j <= nnod; j++)
     {
-      if ((n = myModel[i]->getNodalDOFs(j)) == 0)
+      unsigned char nndof = pch->getNodalDOFs(j);
+      if (nndof == 0)
         continue;
       else if (withID)
-        os << j <<' '<< myModel[i]->getNodeID(j)
+        os << j <<" "<< pch->getNodeID(j) <<"\t\t"<< pch->getCoord(j) <<"\t\t";
-           <<"\t\t"<< myModel[i]->getCoord(j) <<"\t\t";
+
       os << utl::trunc(patchSol[ip++]);
-      for (k = 1; k < n; k++)
+      for (unsigned char k = 1; k < nndof; k++)
         os <<' '<< utl::trunc(patchSol[ip++]);
       os <<'\n';
     }
@@ -1285,6 +1287,8 @@ bool SIMoutput::dumpSolution (const Vector& psol, utl::LogStream& os) const
 {
   if (psol.empty())
     return true;
+  else if (!mySam)
+    return false;
 
   Matrix field;
   size_t i, j, k;
@@ -1488,7 +1492,7 @@ bool SIMoutput::dumpVector (const Vector& vsol, const char* fname,
   if (vsol.empty() || myPoints.empty())
     return true;
 
-  size_t ngNodes = mySam->getNoNodes();
+  size_t ngNodes = this->getNoNodes();
   size_t nComp = vsol.size() / ngNodes;
   if (nComp*ngNodes != vsol.size() || nComp == this->getNoFields())
     nComp = 0; // Using the number of primary field components
@@ -1669,7 +1673,7 @@ void SIMoutput::printNorms (const Vectors& norms, size_t w) const
                <<"\nExact relative error (%) : "<< 100.0*n(4)/n(3);
 
   size_t j = 0;
-  for (const auto& prj : opt.project)
+  for (const SIMoptions::ProjectionMap::value_type& prj : opt.project)
     if (++j < norms.size())
       this->printNormGroup(norms[j],n,prj.second);