src/MULTIPR/MULTIPR_GaussLoc.cxx

   1 // Project MULTIPR, IOLS WP1.2.1 - EDF/CS
   2 // Partitioning/decimation module for the SALOME v3.2 platform
   3
   4 /**
   5  * \file    MULTIPR_GaussLoc.cxx
   6  *
   7  * \brief   see MULTIPR_GaussLoc.hxx
   8  *
   9  * \author  Olivier LE ROUX - CS, Virtual Reality Dpt
  10  *
  11  * \date    01/2007
  12  */
  13
  14 //*****************************************************************************
  15 // Includes section
  16 //*****************************************************************************
  17
  18 #include "MULTIPR_GaussLoc.hxx"
  19 #include "MULTIPR_Utils.hxx"
  20 #include "MULTIPR_Exceptions.hxx"
  21
  22 #include <iostream>
  23
  24 using namespace std;
  25
  26
  27 namespace multipr
  28 {
  29
  30
  31 //*****************************************************************************
  32 // Class GaussLoc implementation
  33 //*****************************************************************************
  34
  35 GaussLoc::GaussLoc()
  36 {
  37     mRefCoo    = NULL;
  38     mGaussCoo  = NULL;
  39     mWeight    = NULL;
  40
  41     reset();
  42 }
  43
  44
  45 GaussLoc::GaussLoc(const GaussLoc& pGaussLoc)
  46 {
  47     mRefCoo   = NULL;
  48     mGaussCoo = NULL;
  49     mWeight   = NULL;
  50
  51     strcpy(mName, pGaussLoc.mName);
  52
  53     mGeom     = pGaussLoc.mGeom;
  54     mDim      = pGaussLoc.mDim;
  55     mNumNodes = pGaussLoc.mNumNodes;
  56     mNumGauss = pGaussLoc.mNumGauss;
  57
  58     if (mDim != (mGeom / 100)) throw IllegalStateException("", __FILE__, __LINE__);
  59     if (mNumNodes != (mGeom % 100)) throw IllegalStateException("", __FILE__, __LINE__);
  60
  61     if (pGaussLoc.mRefCoo != NULL)
  62     {
  63         mRefCoo = new med_float[mDim * mNumNodes];
  64         memcpy(mRefCoo, pGaussLoc.mRefCoo, sizeof(med_float) * mDim * mNumNodes);
  65     }
  66
  67     if (pGaussLoc.mGaussCoo != NULL)
  68     {
  69         mGaussCoo  = new med_float[mDim * mNumGauss];
  70         memcpy(mGaussCoo, pGaussLoc.mGaussCoo, sizeof(med_float) * mDim * mNumGauss);
  71     }
  72
  73     if (pGaussLoc.mWeight != NULL)
  74     {
  75         mWeight = new med_float[mNumGauss];
  76         memcpy(mWeight, pGaussLoc.mWeight, sizeof(med_float) * mNumGauss);
  77     }
  78 }
  79
  80
  81 GaussLoc::~GaussLoc()
  82 {
  83     reset();
  84 }
  85
  86
  87 void GaussLoc::reset()
  88 {
  89     mName[0]  = '\0';
  90     mGeom     = MED_NONE;
  91     mDim      = 0;
  92     mNumNodes = 0;
  93     mNumGauss = 0;
  94
  95     if (mRefCoo   != NULL) { delete[] mRefCoo;   mRefCoo   = NULL; }
  96     if (mGaussCoo != NULL) { delete[] mGaussCoo; mGaussCoo = NULL; }
  97     if (mWeight   != NULL) { delete[] mWeight;   mWeight   = NULL; }
  98 }
  99
 100
 101 void GaussLoc::getCoordGaussPoints(
 102     const med_float* pCooElt,
 103     med_float* pCooGaussPoints) const
 104 {
 105     // debug
 106     //printArray2D(pCooElt, mNumNodes, mDim, "Node");
 107
 108     // WARNING: assumes TETRA10 !!!
 109     // This method is not completely generic and should be extended to support all cases.
 110     if (mGeom != MED_TETRA10) throw UnsupportedOperationException("only support TETRA10 for the moment", __FILE__, __LINE__);
 111
 112     const med_float* pt1 = pCooElt;
 113     const med_float* pt2 = pt1 + mDim;
 114     const med_float* pt3 = pt2 + mDim;
 115     const med_float* pt4 = pt3 + mDim;
 116
 117     const med_float* coeff = mGaussCoo;
 118     med_float* dest        = pCooGaussPoints;
 119
 120     // for each Gauss point
 121     for (int i = 0 ; i < mNumGauss ; i++)
 122     {
 123         dest[0] = pt2[0] + (pt4[0] - pt2[0]) * coeff[0] + (pt1[0] - pt2[0]) * coeff[1] + (pt3[0] - pt2[0]) * coeff[2];
 124         dest[1] = pt2[1] + (pt4[1] - pt2[1]) * coeff[0] + (pt1[1] - pt2[1]) * coeff[1] + (pt3[1] - pt2[1]) * coeff[2];
 125         dest[2] = pt2[2] + (pt4[2] - pt2[2]) * coeff[0] + (pt1[2] - pt2[2]) * coeff[1] + (pt3[2] - pt2[2]) * coeff[2];
 126
 127         // prepare next point
 128         coeff += mDim;
 129         dest += mDim;
 130     }
 131 }
 132
 133
 134 void GaussLoc::readMED(med_idt pMEDfile, med_int pIndex)
 135 {
 136     if (pMEDfile == 0) throw IOException("pMEDfile should not be NULL", __FILE__, __LINE__);
 137     if (pIndex < 1) throw IllegalArgumentException("pIndex should be >= 1", __FILE__, __LINE__);
 138
 139     reset();
 140
 141     med_err ret = MEDgaussInfo(
 142         pMEDfile,
 143         pIndex,
 144         mName,
 145         &mGeom,
 146         &mNumGauss);
 147
 148     if (ret != 0) throw IOException("i/o error while reading Gauss localization information in MED file", __FILE__, __LINE__);
 149
 150     mDim      = mGeom / 100;
 151     mNumNodes = mGeom % 100;
 152
 153     mRefCoo   = new med_float[mDim * mNumNodes];
 154     mGaussCoo = new med_float[mDim * mNumGauss];
 155     mWeight   = new med_float[mNumGauss];
 156
 157     ret = MEDgaussLire(
 158         pMEDfile,
 159         mRefCoo,
 160         mGaussCoo,
 161         mWeight,
 162         MED_FULL_INTERLACE,
 163         mName);
 164
 165     if (ret != 0) throw IOException("i/o error while reading Gauss localization in MED file", __FILE__, __LINE__);
 166 }
 167
 168
 169 void GaussLoc::writeMED(med_idt pMEDfile)
 170 {
 171     if (pMEDfile == 0) throw IOException("", __FILE__, __LINE__);
 172     if (mNumGauss < 0) throw IllegalStateException("", __FILE__, __LINE__);
 173     if (mRefCoo == NULL) throw IllegalStateException("", __FILE__, __LINE__);
 174     if (mGaussCoo == NULL) throw IllegalStateException("", __FILE__, __LINE__);
 175     if (mWeight == NULL) throw IllegalStateException("", __FILE__, __LINE__);
 176     if (strlen(mName) > MED_TAILLE_NOM) throw IllegalStateException("", __FILE__, __LINE__);
 177
 178     med_err ret = MEDgaussEcr(
 179         pMEDfile,
 180         mGeom,
 181         mRefCoo,
 182         MED_FULL_INTERLACE,
 183         mNumGauss,
 184         mGaussCoo,
 185         mWeight,
 186         mName);
 187
 188     if (ret != 0) throw IOException("i/o error while writing Gauss localization", __FILE__, __LINE__);
 189 }
 190
 191
 192 ostream& operator<<(ostream& pOs, GaussLoc& pG)
 193 {
 194     pOs << "Gauss ref:" << endl;
 195     pOs << "    Name     =|" << pG.mName << "|" << endl;
 196     pOs << "    Geom     =" << pG.mGeom << endl;
 197     pOs << "    #Pt Gauss=" << pG.mNumGauss << endl;
 198
 199     pOs << "    Ref nodes coords.: (#nodes=" << pG.mNumNodes << " dim=" << pG.mDim << ")" << endl;
 200     for (int itNode = 0 ; itNode < pG.mNumNodes ; itNode++)
 201     {
 202         pOs << "        Node " << (itNode + 1) << ": ";
 203         for (int itDim = 0; itDim < pG.mDim ; itDim++)
 204         {
 205             pOs << pG.mRefCoo[itNode * pG.mDim + itDim] << " ";
 206         }
 207         pOs << endl;
 208     }
 209
 210     pOs << "    Gauss coords. and weight:" << endl;
 211     for (int itGauss = 0 ; itGauss < pG.mNumGauss ; itGauss++)
 212     {
 213         pOs << "        Pt " << (itGauss+1) << ": ";
 214         for (int itDim = 0; itDim < pG.mDim ; itDim++)
 215         {
 216             pOs << pG.mGaussCoo[itGauss * pG.mDim + itDim] << " ";
 217         }
 218         pOs << "weight=" << pG.mWeight[itGauss];
 219         pOs << endl;
 220     }
 221
 222     // debug
 223     //med_float res[15];
 224     //pG.getCoordGaussPoints(pG.mRefCoo, res);
 225     //printArray2D(res, pG.mNumGauss, pG.mDim, "Gauss pt");
 226
 227     return pOs;
 228 }
 229
 230
 231 } // namespace multipr
 232
 233 // EOF