src/simplefoamwriter.cpp

   1 //
   2 // ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
   3 // +                                                                      +
   4 // + This file is part of enGrid.                                         +
   5 // +                                                                      +
   6 // + Copyright 2008,2009 Oliver Gloth                                     +
   7 // +                                                                      +
   8 // + enGrid is free software: you can redistribute it and/or modify       +
   9 // + it under the terms of the GNU General Public License as published by +
  10 // + the Free Software Foundation, either version 3 of the License, or    +
  11 // + (at your option) any later version.                                  +
  12 // +                                                                      +
  13 // + enGrid is distributed in the hope that it will be useful,            +
  14 // + but WITHOUT ANY WARRANTY; without even the implied warranty of       +
  15 // + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the        +
  16 // + GNU General Public License for more details.                         +
  17 // +                                                                      +
  18 // + You should have received a copy of the GNU General Public License    +
  19 // + along with enGrid. If not, see <http://www.gnu.org/licenses/>.       +
  20 // +                                                                      +
  21 // ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
  22 //
  23
  24 #include "simplefoamwriter.h"
  25 #include "guimainwindow.h"
  26
  27 #include <QFileInfo>
  28 #include <QDir>
  29
  30 bool SimpleFoamWriter::face_t::operator<(const face_t &F) const
  31 {
  32   bool less = false;
  33   if (bc < F.bc) {
  34     less = true;
  35   } else if (bc == F.bc) {
  36     if (owner < F.owner) {
  37       less = true;
  38     } else if (owner == F.owner) {
  39       if (neighbour < F.neighbour) {
  40         less = true;
  41       }
  42     }
  43   }
  44   return less;
  45 }
  46
  47 vec3_t SimpleFoamWriter::face_t::normal(vtkUnstructuredGrid *m_Grid)
  48 {
  49   if (node.size() < 3) EG_BUG;
  50   vec3_t xc(0,0,0);
  51   QVector<vec3_t> x(node.size());
  52   for (int i = 0; i < node.size(); ++i) {
  53     m_Grid->GetPoint(node[i],x[i].data());
  54     xc += x[i];
  55   }
  56   xc *= 1.0/node.size();
  57   vec3_t n(0,0,0);
  58   for (int i = 0; i < node.size()-1; ++i) {
  59     vec3_t a = x[i]   - xc;
  60     vec3_t b = x[i+1] - xc;
  61     n += 0.5*(a.cross(b));
  62   }
  63   vec3_t a = x[node.size()-1] - xc;
  64   vec3_t b = x[0]             - xc;
  65   n += 0.5*(a.cross(b));
  66   return n;
  67 }
  68
  69
  70 SimpleFoamWriter::SimpleFoamWriter()
  71 {
  72   setFormat("Foam boundary files(boundary)");
  73   setExtension("");
  74 }
  75
  76 vtkIdType SimpleFoamWriter::getNeigh(int i_cells, int i_neigh)
  77 {
  78   l2g_t cells = getPartCells();
  79   l2l_t c2c   = getPartC2C();
  80   int n = c2c[i_cells][i_neigh];
  81   if (n >= 0) return cells[n];
  82   EG_BUG;
  83   return -1;
  84 }
  85
  86 void SimpleFoamWriter::addFace(face_t F)
  87 {
  88   if (isVolume(F.neighbour,m_Grid)) {
  89     if (F.neighbour > F.owner) {
  90       F.bc = 0;
  91       m_LFaces.append(F);
  92     }
  93   } else {
  94     F.bc = m_BC->GetValue(F.neighbour);
  95     F.neighbour = -1;
  96     m_LFaces.append(F);
  97   }
  98 }
  99
 100 void SimpleFoamWriter::createFaces()
 101 {
 102   l2g_t cells = getPartCells();
 103
 104   m_LFaces.clear();
 105   EG_VTKDCC(vtkIntArray, cell_code,   m_Grid, "cell_code");
 106   m_BC = cell_code;
 107   m_Eg2Of.fill(-1,cells.size());
 108   int Nvol = 0;
 109   foreach(int i_cells, cells) {
 110     vtkIdType *pts;
 111     vtkIdType  Npts;
 112     m_Grid->GetCellPoints(cells[i_cells], Npts, pts);
 113     vtkIdType type_cell = m_Grid->GetCellType(cells[i_cells]);
 114     vtkIdType id_cell = cells[i_cells];
 115
 116     // tetras
 117     //
 118     if (type_cell == VTK_TETRA) {
 119       m_Eg2Of[id_cell] = Nvol++;
 120       {
 121         face_t F(3,id_cell,getNeigh(i_cells,0));
 122         F.node[0] = pts[2]; F.node[1] = pts[1]; F.node[2] = pts[0];
 123         addFace(F);
 124       }
 125       {
 126         face_t F(3,id_cell,getNeigh(i_cells,1));
 127         F.node[0] = pts[0]; F.node[1] = pts[1]; F.node[2] = pts[3];
 128         addFace(F);
 129       }
 130       {
 131         face_t F(3,id_cell,getNeigh(i_cells,2));
 132         F.node[0] = pts[0]; F.node[1] = pts[3]; F.node[2] = pts[2];
 133         addFace(F);
 134       }
 135       {
 136         face_t F(3,id_cell,getNeigh(i_cells,3));
 137         F.node[0] = pts[1]; F.node[1] = pts[2]; F.node[2] = pts[3];
 138         addFace(F);
 139       }
 140     }
 141
 142     // prisms
 143     //
 144     if (type_cell == VTK_WEDGE) {
 145       m_Eg2Of[id_cell] = Nvol++;
 146       {
 147         face_t F(3,id_cell,getNeigh(i_cells,0));
 148         F.node[0] = pts[0]; F.node[1] = pts[1]; F.node[2] = pts[2];
 149         addFace(F);
 150       }
 151       {
 152         face_t F(3,id_cell,getNeigh(i_cells,1));
 153         F.node[0] = pts[3]; F.node[1] = pts[5]; F.node[2] = pts[4];
 154         addFace(F);
 155       }
 156       {
 157         face_t F(4,id_cell,getNeigh(i_cells,2));
 158         F.node[0] = pts[3]; F.node[1] = pts[4]; F.node[2] = pts[1]; F.node[3] = pts[0];
 159         addFace(F);
 160       }
 161       {
 162         face_t F(4,id_cell,getNeigh(i_cells,3));
 163         F.node[0] = pts[1]; F.node[1] = pts[4]; F.node[2] = pts[5]; F.node[3] = pts[2];
 164         addFace(F);
 165       }
 166       {
 167         face_t F(4,id_cell,getNeigh(i_cells,4));
 168         F.node[0] = pts[0]; F.node[1] = pts[2]; F.node[2] = pts[5]; F.node[3] = pts[3];
 169         addFace(F);
 170       }
 171     }
 172
 173     // hexes
 174     //
 175     if (type_cell == VTK_HEXAHEDRON) {
 176       m_Eg2Of[id_cell] = Nvol++;
 177       {
 178         face_t F(4,id_cell,getNeigh(i_cells,0),0);
 179         F.node[0] = pts[3]; F.node[1] = pts[2]; F.node[2] = pts[1]; F.node[3] = pts[0];
 180         addFace(F);
 181       }
 182       {
 183         face_t F(4,id_cell,getNeigh(i_cells,1),0);
 184         F.node[0] = pts[4]; F.node[1] = pts[5]; F.node[2] = pts[6]; F.node[3] = pts[7];
 185         addFace(F);
 186       }
 187       {
 188         face_t F(4,id_cell,getNeigh(i_cells,2),0);
 189         F.node[0] = pts[0]; F.node[1] = pts[1]; F.node[2] = pts[5]; F.node[3] = pts[4];
 190         addFace(F);
 191       }
 192       {
 193         face_t F(4,id_cell,getNeigh(i_cells,3),0);
 194         F.node[0] = pts[3]; F.node[1] = pts[7]; F.node[2] = pts[6]; F.node[3] = pts[2];
 195         addFace(F);
 196       }
 197       {
 198         face_t F(4,id_cell,getNeigh(i_cells,4),0);
 199         F.node[0] = pts[0]; F.node[1] = pts[4]; F.node[2] = pts[7]; F.node[3] = pts[3];
 200         addFace(F);
 201       }
 202       {
 203         face_t F(4,id_cell,getNeigh(i_cells,5),0);
 204         F.node[0] = pts[1]; F.node[1] = pts[2]; F.node[2] = pts[6]; F.node[3] = pts[5];
 205         addFace(F);
 206       }
 207     }
 208   }
 209
 210   m_Faces.resize(m_LFaces.size());
 211   qCopy(m_LFaces.begin(), m_LFaces.end(), m_Faces.begin());
 212   qSort(m_Faces);
 213
 214 }
 215
 216
 217 void SimpleFoamWriter::writePoints()
 218 {
 219   l2g_t nodes = getPartNodes();
 220
 221   QString filename = m_Path + "points";
 222   QFile file(filename);
 223   file.open(QIODevice::WriteOnly);
 224   QTextStream f(&file);
 225   f << "/*--------------------------------*- C++ -*----------------------------------*\\\n";
 226   f << "| =========                 |                                                 |\n";
 227   f << "| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |\n";
 228   f << "|  \\    /   O peration     | Version:  1.5                                   |\n";
 229   f << "|   \\  /    A nd           | Web:      http://www.OpenFOAM.org               |\n";
 230   f << "|    \\/     M anipulation  |                                                 |\n";
 231   f << "\\*---------------------------------------------------------------------------*/\n\n";
 232   f << "FoamFile\n";
 233   f << "{\n";
 234   f << "    version     2.0;\n";
 235   f << "    format      ascii;\n";
 236   f << "    class       vectorField;\n";
 237   f << "    location    \"constant/polyMesh\";\n";
 238   f << "    object      points;\n";
 239   f << "}\n\n";
 240   f << "// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //\n\n";
 241   f << nodes.size() << "\n(\n";
 242   foreach(int i_nodes, nodes) {
 243     vtkIdType id_node = nodes[i_nodes];
 244     vec3_t x;
 245     m_Grid->GetPoint(id_node,x.data());
 246     f.setRealNumberPrecision(16);
 247     f << "(" << x[0] << " " << x[1] << " " << x[2] << ")\n";
 248   }
 249   f << ")\n\n";
 250   f << "// ************************************************************************* //\n\n\n";
 251 }
 252
 253 void SimpleFoamWriter::writeFaces()
 254 {
 255   QString filename = m_Path + "faces";
 256   QFile file(filename);
 257   file.open(QIODevice::WriteOnly);
 258   QTextStream f(&file);
 259   f << "/*--------------------------------*- C++ -*----------------------------------*\\\n";
 260   f << "| =========                 |                                                 |\n";
 261   f << "| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |\n";
 262   f << "|  \\    /   O peration     | Version:  1.5                                   |\n";
 263   f << "|   \\  /    A nd           | Web:      http://www.OpenFOAM.org               |\n";
 264   f << "|    \\/     M anipulation  |                                                 |\n";
 265   f << "\\*---------------------------------------------------------------------------*/\n\n";
 266   f << "FoamFile\n";
 267   f << "{\n";
 268   f << "    version     2.0;\n";
 269   f << "    format      ascii;\n";
 270   f << "    class       faceList;\n";
 271   f << "    location    \"constant/polyMesh\";\n";
 272   f << "    object      faces;\n";
 273   f << "}\n\n";
 274   f << "// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //\n\n";
 275   f << m_Faces.size() << "\n(\n";
 276   foreach (face_t F, m_Faces) {
 277     f << F.node.size() << "(";
 278     for (int i = 0; i < F.node.size(); ++i) {
 279       f << F.node[i];
 280       if (i == F.node.size()-1) {
 281         f << ")\n";
 282       } else {
 283         f << " ";
 284       }
 285     }
 286   }
 287   f << ")\n\n";
 288   f << "// ************************************************************************* //\n\n\n";
 289 }
 290
 291 void SimpleFoamWriter::writeOwner()
 292 {
 293   QString filename = m_Path + "owner";
 294   QFile file(filename);
 295   file.open(QIODevice::WriteOnly);
 296   QTextStream f(&file);
 297   f << "/*--------------------------------*- C++ -*----------------------------------*\\\n";
 298   f << "| =========                 |                                                 |\n";
 299   f << "| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |\n";
 300   f << "|  \\    /   O peration     | Version:  1.5                                   |\n";
 301   f << "|   \\  /    A nd           | Web:      http://www.OpenFOAM.org               |\n";
 302   f << "|    \\/     M anipulation  |                                                 |\n";
 303   f << "\\*---------------------------------------------------------------------------*/\n\n";
 304   f << "FoamFile\n";
 305   f << "{\n";
 306   f << "    version     2.0;\n";
 307   f << "    format      ascii;\n";
 308   f << "    class       labelList;\n";
 309   f << "    location    \"constant/polyMesh\";\n";
 310   f << "    object      owner;\n";
 311   f << "}\n\n";
 312   f << "// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //\n\n";
 313   f << m_Faces.size() << "\n(\n";
 314   foreach (face_t F, m_Faces) {
 315     f << m_Eg2Of[F.owner] << "\n";
 316   }
 317   f << ")\n\n";
 318   f << "// ************************************************************************* //\n\n\n";
 319 }
 320
 321 void SimpleFoamWriter::writeNeighbour()
 322 {
 323   QString filename = m_Path + "neighbour";
 324   QFile file(filename);
 325   file.open(QIODevice::WriteOnly);
 326   QTextStream f(&file);
 327   f << "/*--------------------------------*- C++ -*----------------------------------*\\\n";
 328   f << "| =========                 |                                                 |\n";
 329   f << "| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |\n";
 330   f << "|  \\    /   O peration     | Version:  1.5                                   |\n";
 331   f << "|   \\  /    A nd           | Web:      http://www.OpenFOAM.org               |\n";
 332   f << "|    \\/     M anipulation  |                                                 |\n";
 333   f << "\\*---------------------------------------------------------------------------*/\n\n";
 334   f << "FoamFile\n";
 335   f << "{\n";
 336   f << "    version     2.0;\n";
 337   f << "    format      ascii;\n";
 338   f << "    class       labelList;\n";
 339   f << "    location    \"constant/polyMesh\";\n";
 340   f << "    object      neighbour;\n";
 341   f << "}\n\n";
 342   f << "// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //\n\n";
 343   f << m_Faces.size() << "\n(\n";
 344   foreach (face_t F, m_Faces) {
 345     if (F.neighbour == -1) {
 346       f << "-1\n";
 347     } else {
 348       f << m_Eg2Of[F.neighbour] << "\n";
 349     }
 350   }
 351   f << ")\n\n";
 352   f << "// ************************************************************************* //\n\n\n";
 353 }
 354
 355 void SimpleFoamWriter::writeBoundary(int faces_offset)
 356 {
 357   QString filename = m_Path + "boundary";
 358   QFile file(filename);
 359   file.open(QIODevice::WriteOnly);
 360   QTextStream f(&file);
 361   f << "/*--------------------------------*- C++ -*----------------------------------*\\\n";
 362   f << "| =========                 |                                                 |\n";
 363   f << "| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |\n";
 364   f << "|  \\    /   O peration     | Version:  1.5                                   |\n";
 365   f << "|   \\  /    A nd           | Web:      http://www.OpenFOAM.org               |\n";
 366   f << "|    \\/     M anipulation  |                                                 |\n";
 367   f << "\\*---------------------------------------------------------------------------*/\n\n";
 368   f << "FoamFile\n";
 369   f << "{\n";
 370   f << "    version     2.0;\n";
 371   f << "    format      ascii;\n";
 372   f << "    class       polyBoundaryMesh;\n";
 373   f << "    location    \"constant/polyMesh\";\n";
 374   f << "    object      boundary;\n";
 375   f << "}\n\n";
 376   f << "// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //\n\n";
 377   QSet<int> bcs;
 378   foreach (face_t F, m_Faces) {
 379     if (F.bc != 0) {
 380       bcs.insert(F.bc);
 381     }
 382   }
 383   f << bcs.size() << "\n(\n";
 384   QVector<patch_t> patch(bcs.size());
 385   int N_bc = 0;
 386   foreach (int bc, bcs) {
 387     int nFaces = 0;
 388     int startFace = -1;
 389     for (int i = 0; i < m_Faces.size(); ++i) {
 390       if (m_Faces[i].bc == bc) {
 391         ++nFaces;
 392         if (startFace == -1) {
 393           startFace = i;
 394         }
 395       }
 396     }
 397     if (startFace == -1) {
 398       EG_BUG;
 399     }
 400     patch_t P;
 401     P.bc = bc;
 402     P.startFace = startFace;
 403     P.nFaces = nFaces;
 404     patch[N_bc] = P;
 405     ++N_bc;
 406   }
 407   qSort(patch);
 408   foreach (patch_t P, patch) {
 409     BoundaryCondition BC = getBC(P.bc);
 410     QString num;
 411     num.setNum(P.bc);
 412     QString name = BC.getName();
 413     QString type = BC.getType();
 414     if (GuiMainWindow::pointer()->physicalTypeDefined(type)) {
 415       PhysicalBoundaryCondition PBC = GuiMainWindow::pointer()->getPhysicalBoundaryCondition(type);
 416       type = PBC.getFoamType();
 417     }
 418     if (name == "unknown") {
 419       name = "BC_" + num;
 420     }
 421     f << "    " << name << "\n";
 422     f << "    {\n";
 423     f << "        type        " << type << ";\n";
 424     f << "        nFaces      " << P.nFaces << ";\n";
 425     f << "        startFace   " << P.startFace + faces_offset << ";\n";
 426     f << "    }\n";
 427   }
 428   f << ")\n\n";
 429   f << "// ************************************************************************* //\n\n\n";
 430 }
 431
 432 void SimpleFoamWriter::operateOnGivenFileName()
 433 {
 434   if (isValid()) {
 435     QString p1 = getFileName();
 436     QString p2 = p1 + "/constant";
 437     QDir d1(p1);
 438     QDir d2(p2);
 439     if (!d1.exists()) {
 440       EG_BUG;
 441     }
 442     if (!d2.exists()) {
 443       d1.mkdir("constant");
 444       d2 = QDir(p2);
 445     }
 446     d1 = d2;
 447     p1 = p2;
 448     p2 = p1 + "/polyMesh";
 449     d2 = QDir(p2);
 450     if (!d2.exists()) {
 451       d1.mkdir("polyMesh");
 452     }
 453     m_Path = getFileName() + "/constant/polyMesh/";
 454     if (!QDir(m_Path).exists()) {
 455       EG_BUG;
 456     }
 457     createFaces();
 458     writePoints();
 459     writeFaces();
 460     writeOwner();
 461     writeNeighbour();
 462     writeBoundary();
 463   }
 464 }
 465
 466 void SimpleFoamWriter::operate()
 467 {
 468   try {
 469     readOutputDirectory();
 470     if (isValid()) {
 471       QString p1 = getFileName();
 472       QString p2 = p1 + "/constant";
 473       QDir d1(p1);
 474       QDir d2(p2);
 475       if (!d1.exists()) {
 476         EG_BUG;
 477       }
 478       if (!d2.exists()) {
 479         d1.mkdir("constant");
 480         d2 = QDir(p2);
 481       }
 482       d1 = d2;
 483       p1 = p2;
 484       p2 = p1 + "/polyMesh";
 485       d2 = QDir(p2);
 486       if (!d2.exists()) {
 487         d1.mkdir("polyMesh");
 488       }
 489       m_Path = getFileName() + "/constant/polyMesh/";
 490       if (!QDir(m_Path).exists()) {
 491         EG_BUG;
 492       }
 493       createFaces();
 494       writePoints();
 495       writeFaces();
 496       writeOwner();
 497       writeNeighbour();
 498       writeBoundary();
 499     }
 500   } catch (Error err) {
 501     err.display();
 502   }
 503 }