Commit f84bcb85 authored by Christophe Geuzaine's avatar Christophe Geuzaine

cleanup

parent 3a08b7a2
Pipeline #1850 passed with stage
in 17 minutes 36 seconds
......@@ -3,15 +3,20 @@
// See the LICENSE.txt file for license information. Please report all
// bugs and problems to the public mailing list <gmsh@onelab.info>.
//
// GModelIO_CGNS.cpp - Copyright (C) 2008-2017 S. Guzik, B. Gorissen, K. Hillewaert, C. Geuzaine, J.-F. Remacle
// CGNSFunctions.cpp - Copyright (C) 2008-2018 S. Guzik, B. Gorissen,
// K. Hillewaert, C. Geuzaine, J.-F. Remacle
#include "GmshConfig.h"
#include "GmshMessage.h"
#include "GModel.h"
#include "CGNSOptions.h"
#if defined(HAVE_LIBCGNS)
#include <vector>
#include <sstream>
#include <fstream>
#include <iomanip>
#include "GmshMessage.h"
#include "GModel.h"
#include "CGNSOptions.h"
#include "BasisFactory.h"
#include "SVector3.h"
#include "fullMatrix.h"
......@@ -23,84 +28,65 @@ namespace CGNS {
using namespace CGNS;
#include <vector>
// -----------------------------------------------------------------------------
int parentFromCGNSType(ElementType_t cgnsType) {
switch(cgnsType) {
int parentFromCGNSType(ElementType_t cgnsType)
{
switch(cgnsType){
case NODE:
return TYPE_PNT;
break;
case BAR_2: case BAR_3: case BAR_4: case BAR_5:
return TYPE_LIN;
break;
case TRI_3: case TRI_6: case TRI_9: case TRI_10: case TRI_12: case TRI_15:
return TYPE_TRI;
break;
case QUAD_4: case QUAD_9: case QUAD_16: case QUAD_25:
case QUAD_8: case QUAD_12: case QUAD_P4_16:
return TYPE_QUA;
break;
case TETRA_4: case TETRA_10: case TETRA_20: case TETRA_35:
case TETRA_16: case TETRA_22:
case TETRA_34:
return TYPE_TET;
break;
case PYRA_5: case PYRA_14: case PYRA_30: case PYRA_55:
case PYRA_13: case PYRA_21: case PYRA_P4_29:
case PYRA_29: case PYRA_50:
return TYPE_PYR;
break;
// pentahedra
case PENTA_6: case PENTA_18: case PENTA_40: case PENTA_75:
case PENTA_15: case PENTA_24: case PENTA_33:
case PENTA_38: case PENTA_66:
return TYPE_PRI;
break;
// hexahedra
case HEXA_8: case HEXA_27: case HEXA_64: case HEXA_125:
case HEXA_20: case HEXA_32: case HEXA_44:
// case HEXA_26:
case HEXA_56: case HEXA_98:
return TYPE_HEX;
break;
case MIXED:
case NGON_n:
case NFACE_n:
case ElementTypeUserDefined:
case ElementTypeNull:
Msg::Warning("Finding parent type for unsupported CGNS element type %i",cgnsType);
Msg::Warning("Finding parent type for unsupported CGNS element type %i",
cgnsType);
return -1;
break;
}
return -1;
}
int orderFromCGNSType(ElementType_t cgnsType) {
int orderFromCGNSType(ElementType_t cgnsType)
{
switch( cgnsType ) {
switch(cgnsType){
case NODE:
return 0;
break;
case BAR_2: case TRI_3: case QUAD_4: case TETRA_4: case PYRA_5: case PENTA_6: case HEXA_8:
return 1;
break;
case BAR_3: case TRI_6: case QUAD_9: case TETRA_10: case PYRA_14: case PENTA_18: case HEXA_27:
case QUAD_8: case PYRA_13: case PENTA_15: case HEXA_20:
// case HEXA_26:
return 2;
break;
case BAR_4: case TRI_10: case QUAD_16: case TETRA_20: case PYRA_30: case PENTA_40: case HEXA_64:
case TRI_9: case QUAD_12: case TETRA_16: case PYRA_21: case PENTA_24: case HEXA_32:
case PYRA_29: case PENTA_38: case HEXA_56:
return 3;
break;
case BAR_5: case TRI_15: case QUAD_25: case TETRA_35: case PYRA_55: case PENTA_75: case HEXA_125:
case TRI_12: case QUAD_P4_16: case TETRA_22: case PYRA_P4_29: case PENTA_33: case HEXA_44:
case TETRA_34: case PYRA_50: case PENTA_66: case HEXA_98:
return 4;
break;
case MIXED:
case NGON_n:
case NFACE_n:
......@@ -108,14 +94,14 @@ int orderFromCGNSType(ElementType_t cgnsType) {
case ElementTypeNull:
Msg::Warning("Finding order for unsupported CGNS element type %i",cgnsType);
return -1;
break;
}
return -1;
}
bool completeCGNSType(ElementType_t cgnsType) {
bool completeCGNSType(ElementType_t cgnsType)
{
switch( cgnsType ) {
switch(cgnsType){
// complete elements
case NODE:
case BAR_2: case TRI_3: case QUAD_4: case TETRA_4: case PYRA_5: case PENTA_6: case HEXA_8:
......@@ -123,26 +109,22 @@ bool completeCGNSType(ElementType_t cgnsType) {
case BAR_4: case TRI_10: case QUAD_16: case TETRA_20: case PYRA_30: case PENTA_40: case HEXA_64:
case BAR_5: case TRI_15: case QUAD_25: case TETRA_35: case PYRA_55: case PENTA_75: case HEXA_125:
return true;
break;
// serendipity edge elements
case QUAD_8: case PYRA_13: case PENTA_15: case HEXA_20:
case QUAD_12: case TETRA_16: case PYRA_21: case PENTA_24: case HEXA_32:
case QUAD_8: case PYRA_13: case PENTA_15: case HEXA_20:
case QUAD_12: case TETRA_16: case PYRA_21: case PENTA_24: case HEXA_32:
case TRI_12: case QUAD_P4_16: case TETRA_22: case PYRA_P4_29: case PENTA_33: case HEXA_44:
return false;
break;
// serendipity elements
// case HEXA_26:
case PYRA_29: case PENTA_38: case HEXA_56:
case PYRA_29: case PENTA_38: case HEXA_56:
case TETRA_34: case PYRA_50: case PENTA_66: case HEXA_98:
default:
return false;
}
return false;
}
int tagFromCGNSType(ElementType_t cgnsType) {
int tagFromCGNSType(ElementType_t cgnsType)
{
return ElementType::getType(parentFromCGNSType(cgnsType),
orderFromCGNSType(cgnsType),
!completeCGNSType(cgnsType));
......@@ -153,7 +135,8 @@ std::vector<SVector3> generatePointsCGNS(int,int,bool,bool);
void addEdgePointsCGNS(const SVector3 p0,
const SVector3 p1,
int order,
std::vector<SVector3>& points) {
std::vector<SVector3>& points)
{
double ds = 1./order;
for (int i=1;i<order;i++) {
double f = ds*i;
......@@ -166,8 +149,8 @@ void addTriPointsCGNS(const SVector3 p0,
const SVector3 p2,
int order,
bool equidistant,
std::vector<SVector3>& points) {
std::vector<SVector3>& points)
{
std::vector<SVector3> triPoints = generatePointsCGNS(TYPE_TRI,
order-3,
true,true);
......@@ -186,10 +169,8 @@ void addTriPointsCGNS(const SVector3 p0,
}
}
void addQuaPointsCGNS(int order,std::vector<SVector3>& points) {
void addQuaPointsCGNS(int order,std::vector<SVector3>& points)
{
if (order > 2) {
double scale = double (order -2) / double(order);
......@@ -214,8 +195,8 @@ void addQuaPointsCGNS(const SVector3 p0,
const SVector3 p2,
const SVector3 p3,
int order,
std::vector<SVector3>& points) {
std::vector<SVector3>& points)
{
std::vector<SVector3> quaPoints;
addQuaPointsCGNS(order,quaPoints);
......@@ -232,8 +213,8 @@ void addQuaPointsCGNS(const SVector3 p0,
}
}
void print(std::vector<SVector3>& points,const char* title,int iStart,int iEnd=-1) {
void print(std::vector<SVector3>& points,const char* title,int iStart,int iEnd=-1)
{
iEnd = iEnd == -1 ? points.size() : iEnd;
std::cout << title << std::endl;
......@@ -244,12 +225,11 @@ void print(std::vector<SVector3>& points,const char* title,int iStart,int iEnd=-
}
}
std::vector<SVector3> generatePointsCGNS(int parentType,
int order,
bool complete,
bool equidistant) {
bool equidistant)
{
std::vector<SVector3> pp;
if (order == 0) pp.push_back(SVector3(0,0,0));
......@@ -479,9 +459,8 @@ std::vector<SVector3> generatePointsCGNS(int parentType,
return pp;
}
fullMatrix<double> getTransformationToGmsh(ElementType_t cgnsType) {
fullMatrix<double> getTransformationToGmsh(ElementType_t cgnsType)
{
int parent = parentFromCGNSType(cgnsType);
int order = orderFromCGNSType(cgnsType);
bool complete = completeCGNSType(cgnsType);
......@@ -508,12 +487,8 @@ fullMatrix<double> getTransformationToGmsh(ElementType_t cgnsType) {
return tfo;
}
#include <sstream>
#include <fstream>
#include <iomanip>
int* getRenumberingToGmsh(ElementType_t cgnsType) {
int* getRenumberingToGmsh(ElementType_t cgnsType)
{
int parent = parentFromCGNSType(cgnsType);
int order = orderFromCGNSType(cgnsType);
bool complete = completeCGNSType(cgnsType);
......@@ -553,8 +528,8 @@ int* getRenumberingToGmsh(ElementType_t cgnsType) {
return renum;
}
int gridLocationDimCGNS(GridLocation_t lt) {
int gridLocationDimCGNS(GridLocation_t lt)
{
switch (lt) {
case CellCenter:
return 3;
......@@ -579,9 +554,8 @@ int gridLocationDimCGNS(GridLocation_t lt) {
return -1;
}
GridLocation_t unstructuredGridLocationCGNS(int dim) {
GridLocation_t unstructuredGridLocationCGNS(int dim)
{
switch (dim) {
case 3:
return CellCenter;
......@@ -596,9 +570,7 @@ GridLocation_t unstructuredGridLocationCGNS(int dim) {
return Vertex;
break;
}
return GridLocationNull;
}
#endif
......@@ -513,13 +513,10 @@ double GEdge::parFromPoint(const SPoint3 &P) const
return t;
}
bool GEdge::refineProjection(const SVector3& Q,
double& u,
int MaxIter,
double relax,
double tol,
double& err) const {
bool GEdge::refineProjection(const SVector3 &Q, double &u,
int MaxIter, double relax, double tol,
double &err) const
{
double maxDist = tol * CTX::instance()->lc;
SVector3 P = position(u);
......@@ -551,7 +548,6 @@ bool GEdge::refineProjection(const SVector3& Q,
if (err <= maxDist) return true;
return false;
}
bool GEdge::XYZToU(const double X, const double Y, const double Z,
......
......@@ -190,13 +190,9 @@ class GEdge : public GEntity {
// return the parmater location on the edge given a point in space
// that is on the edge
virtual double parFromPoint(const SPoint3 &P) const;
virtual bool refineProjection(const SVector3& Q,
double& u,
int MaxIter,
double relax,
double tol,double& err) const;
virtual bool refineProjection(const SVector3 &Q, double &u,
int MaxIter, double relax, double tol,
double& err) const;
// compute the parameter U from a point XYZ
virtual bool XYZToU(const double X, const double Y, const double Z,
......
......@@ -217,15 +217,11 @@ void GEntity::addVerticesInSet(std::set<MVertex*>&vtcs,bool closure) const
void GEntity::updateCorrespondingVertices()
{
if (_meshMaster != this && affineTransform.size() == 16) {
correspondingVertices.clear();
closestVertexFinder cvf(_meshMaster,true);
if (cvf.getNbVtcs()) {
std::vector<double> tfo = affineTransform;
std::vector<double> inv;
invertAffineTransformation(tfo,inv);
......@@ -235,14 +231,12 @@ void GEntity::updateCorrespondingVertices()
std::set<MVertex*>::iterator vIter = vtcs.begin();
for (;vIter!=vtcs.end();++vIter) {
MVertex* tv = *vIter;
// double tgt[4] = {tv->x(),tv->y(),tv->z(),1};
// double xyz[4] = {0,0,0,0};
// int idx = 0;
// for (int i=0;i<3;i++) for (int j=0;j<4;j++) tgt[i] += inv[idx++] * ori[j];
MVertex* sv = cvf(tv->point(),inv);
correspondingVertices[tv] = sv;
......@@ -265,33 +259,24 @@ void GEntity::updateCorrespondingVertices()
void GEntity::copyMasterCoordinates()
{
if (_meshMaster != this && affineTransform.size() == 16) {
std::map<MVertex*,MVertex*>::iterator cvIter = correspondingVertices.begin();
for (;cvIter!=correspondingVertices.end();++cvIter) {
MVertex* tv = cvIter->first;
MVertex* sv = cvIter->second;
double src[4] = {sv->x(),sv->y(),sv->z(),1};
double tgt[3] = {0,0,0};
int idx = 0;
for (int i=0;i<3;i++) for (int j=0;j<4;j++) tgt[i] += affineTransform[idx++] * src[j];
tv->x() = tgt[0];
tv->y() = tgt[1];
tv->z() = tgt[2];
}
cvIter = correspondingHOPoints.begin();
for (;cvIter!=correspondingHOPoints.end();++cvIter) {
MVertex* tv = cvIter->first;
MVertex* sv = cvIter->second;
......@@ -301,13 +286,9 @@ void GEntity::copyMasterCoordinates()
int idx = 0;
for (int i=0;i<3;i++) for (int j=0;j<4;j++) tgt[i] += affineTransform[idx++] * src[j];
tv->x() = tgt[0];
tv->y() = tgt[1];
tv->z() = tgt[2];
}
}
}
......@@ -2127,15 +2127,15 @@ void GFace::alignElementsWithMaster()
MVertex* tv = face->getVertex(j);
std::map<MVertex*,MVertex*>::iterator cIter = correspondingVertices.find(tv);
if (cIter == correspondingVertices.end()) throw;
vtcs.push_back(cIter->second);
if (cIter != correspondingVertices.end())
vtcs.push_back(cIter->second);
}
MFace mf(vtcs);
std::set<MFace,Less_Face>::iterator sIter = srcFaces.find(mf);
if (sIter==srcFaces.end()) throw;
if (sIter == srcFaces.end()) continue;
MFace of = *sIter;
......@@ -2149,19 +2149,17 @@ void GFace::alignElementsWithMaster()
case 3:
{
MTriangle* tri = dynamic_cast<MTriangle*>(face);
if (!tri) throw;
tri->reorient(orientation,swap);
if (tri)
tri->reorient(orientation,swap);
break;
}
case 4:
{
MQuadrangle* qua = dynamic_cast<MQuadrangle*>(face);
if (!qua) throw;
qua->reorient(orientation,swap);
if (qua)
qua->reorient(orientation,swap);
break;
}
default:
throw;
}
}
}
......
......@@ -121,6 +121,11 @@ class GModel {
void _createFMInternals();
void _deleteFMInternals();
#if defined(HAVE_LIBCGNS)
int _readCGNSStructured(const std::string& name);
int _readCGNSUnstructured(const std::string& name);
#endif
// characteristic length (mesh size) fields
FieldManager *_fields;
......@@ -734,25 +739,6 @@ class GModel {
// GAMBIT neutral mesh file (.neu)
int writeNEU(const std::string &name, bool saveAll, double scalingFactor);
#if defined HAVE_LIBCGNS
protected:
int readCGNSBase(const std::string& name,int& nbases) const;
int readCGNSStructured (const std::string& name);
int readCGNSUnstructured(const std::string& name);
int addCGNSPoints(const std::string&,
int fileIndex,
int baseIndex,
int zoneIndex,
CGNS::cgsize_t nbPoints,
int dim,
double scale,
GEntity* ge,
int& pointIndex,
std::vector<MVertex*>& vertices);
#endif
};
#endif
......@@ -3,7 +3,8 @@
// See the LICENSE.txt file for license information. Please report all
// bugs and problems to the public mailing list <gmsh@onelab.info>.
//
// GModelIO_CGNS.cpp - Copyright (C) 2008-2012 S. Guzik, B. Gorissen, C. Geuzaine, J.-F. Remacle
// GModelIO_CGNS.cpp - Copyright (C) 2008-2012 S. Guzik, B. Gorissen,
// C. Geuzaine, J.-F. Remacle
#include "GmshConfig.h"
#include "GmshMessage.h"
......@@ -962,10 +963,10 @@ protected:
break;
case PointList:
case PointListDonor:
for (int i=0;i<size;i++) pts.push_back(ptsList[i]-1);
for (int i=0;i<size;i++) pts.push_back(ptsList[i]-1);
break;
default:
throw;
break;
}
}
......@@ -999,11 +1000,11 @@ public: // constructors
std::vector<int> tgtPts;
addPoints(tType,tSize,tPts,tgtPts);
if (tgtPts.size() != srcPts.size()) throw;
for (unsigned i=0;i<tgtPts.size();i++) {
tgtToSrcPts[tgtPts[i]] = srcPts[i];
srcToTgtPts[srcPts[i]] = tgtPts[i];
if (tgtPts.size() == srcPts.size()) {
for (unsigned i=0;i<tgtPts.size();i++) {
tgtToSrcPts[tgtPts[i]] = srcPts[i];
srcToTgtPts[srcPts[i]] = tgtPts[i];
}
}
}
......@@ -1332,17 +1333,16 @@ bool readCGNSPeriodicConnections(int fileIndex,
// -----------------------------------------------------------------------------
int GModel::addCGNSPoints(const string& fileName,
int fileIndex,
int baseIndex,
int zoneIndex,
cgsize_t nbPoints,int dim,
double scale,
GEntity* entity,
int& index,
std::vector<MVertex*>& vertices) {
int addCGNSPoints(const string& fileName,
int fileIndex,
int baseIndex,
int zoneIndex,
cgsize_t nbPoints,int dim,
double scale,
GEntity* entity,
int& index,
std::vector<MVertex*>& vertices)
{
int nbDim;
if (cg_ncoords(fileIndex,baseIndex,zoneIndex,&nbDim) != CG_OK) {
......@@ -1397,10 +1397,9 @@ int GModel::addCGNSPoints(const string& fileName,
// -----------------------------------------------------------------------------