src/Plugins/MedReader/IO/vtkMedLocalization.cxx

   1 // Copyright (C) 2010-2011  CEA/DEN, EDF R&D
   2 //
   3 // This library is free software; you can redistribute it and/or
   4 // modify it under the terms of the GNU Lesser General Public
   5 // License as published by the Free Software Foundation; either
   6 // version 2.1 of the License.
   7 //
   8 // This library is distributed in the hope that it will be useful,
   9 // but WITHOUT ANY WARRANTY; without even the implied warranty of
  10 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  11 // Lesser General Public License for more details.
  12 //
  13 // You should have received a copy of the GNU Lesser General Public
  14 // License along with this library; if not, write to the Free Software
  15 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
  16 //
  17 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
  18 //
  19
  20 #include "vtkMedLocalization.h"
  21
  22 #include "vtkObjectFactory.h"
  23 #include "vtkDoubleArray.h"
  24 #include "vtkIntArray.h"
  25 #include "vtkFunctionParser.h"
  26
  27 #include "vtkMedUtilities.h"
  28 #include "vtkMedSetGet.h"
  29 #include "vtkMedFile.h"
  30 #include "vtkMedDriver.h"
  31 #include "vtkMedInterpolation.h"
  32 #include "vtkMedFraction.h"
  33
  34 #include "med.h"
  35
  36 // SEG2
  37 const static int SEG2_dim = 1;
  38 const static int SEG2_nnode = 2;
  39 static const int SEG2_aster2med[SEG2_nnode] =
  40 {0, 1};
  41 static const char* SEG2_varnames[SEG2_dim] = {"x"};
  42 static const char* SEG2_functions[SEG2_nnode] =
  43 {"1/2*(1-x)",
  44  "1/2*(1+x)"};
  45
  46 // SEG3
  47 const static int SEG3_dim = 1;
  48 const static int SEG3_nnode = 3;
  49 static const int SEG3_aster2med[SEG3_nnode] =
  50 {0, 1, 2};
  51 static const char* SEG3_varnames[SEG3_dim] = {"x"};
  52 static const char* SEG3_functions[SEG3_nnode] =
  53 {"-1/2*(1-x)*x",
  54   "1/2*(1+x)*x",
  55   "(1+x)*(1-x)"};
  56
  57 // SEG4
  58 const static int SEG4_dim = 1;
  59 const static int SEG4_nnode = 4;
  60 static const int SEG4_aster2med[SEG4_nnode] =
  61 {0, 1, 2, 3};
  62 static const char* SEG4_varnames[SEG4_dim] = {"x"};
  63 static const char* SEG4_functions[SEG4_nnode] =
  64 {"16/9*(1-x)*(x+1/3)*(x-1/3)",
  65 "-16/9*(1+x)*(1/3-x)*(x+1/3)",
  66 "16/27*(x-1)*(x+1)*(x-1/3)",
  67 "-16/27*(x-1)*(x+1)*(x+1/3)"};
  68
  69 // TRIA3
  70 const static int TRIA3_dim = 2;
  71 const static int TRIA3_nnode = 3;
  72 static const int TRIA3_aster2med[TRIA3_nnode] =
  73 {0, 1, 2};
  74 static const char* TRIA3_varnames[TRIA3_dim] = {"x", "y"};
  75 static const char* TRIA3_functions[TRIA3_nnode] =
  76 {"1-x-y",
  77  "x",
  78  "y"};
  79
  80 // TRIA6
  81 const static int TRIA6_dim = 2;
  82 const static int TRIA6_nnode = 6;
  83 static const int TRIA6_aster2med[TRIA6_nnode] =
  84 {0, 1, 2, 3, 4, 5};
  85 static const char* TRIA6_varnames[TRIA6_dim] = {"x", "y"};
  86 static const char* TRIA6_functions[TRIA6_nnode] =
  87 {"-(1-x-y)*(1-2*(1-x-y))",
  88  "-x*(1-2*x)",
  89  "-y*(1-2*y)",
  90  "4*x*(1-x-y)",
  91  "4*x*y",
  92  "4*y*(1-x-y)"};
  93
  94 // TRIA7
  95 const static int TRIA7_dim = 2;
  96 const static int TRIA7_nnode = 7;
  97 static const int TRIA7_aster2med[TRIA7_nnode] =
  98 {0, 1, 2, 3, 4, 5, 6};
  99 static const char* TRIA7_varnames[TRIA7_dim] = {"x", "y"};
 100 static const char* TRIA7_functions[TRIA7_nnode] =
 101 {"1-3*(x+y)+2*(x*x+y*y)+7*x*y-3*x*y*(x+y)",
 102  "x*(-1+2*x+3*y-3*y*(x+y))",
 103  "y*(-1+2*x+3*y-3*x*(x+y))",
 104  "4*x*(1-x-4*y+3*y*(x+y))",
 105  "4*x*y*(-2+3*(x+y))",
 106  "4*y*(1-4*x-y+3*x*(x+y))",
 107  "27*x*y*(1-x-y)"};
 108
 109 // QUAD4
 110 const static int QUAD4_dim = 2;
 111 const static int QUAD4_nnode = 4;
 112 static const int QUAD4_aster2med[QUAD4_nnode] =
 113 {0, 1, 2, 3};
 114 static const char* QUAD4_varnames[QUAD4_dim] = {"x", "y"};
 115 static const char* QUAD4_functions[QUAD4_nnode] =
 116 {"(1-x)*(1-y)/4",
 117  "(1+x)*(1-y)/4",
 118  "(1+x)*(1+y)/4",
 119  "(1-x)*(1+y)/4"};
 120
 121 // QUAD8
 122 const static int QUAD8_dim = 2;
 123 const static int QUAD8_nnode = 8;
 124 static const int QUAD8_aster2med[QUAD8_nnode] =
 125 {0, 1, 2, 3, 4, 5, 6, 7};
 126 static const char* QUAD8_varnames[QUAD8_dim] = {"x", "y"};
 127 static const char* QUAD8_functions[QUAD8_nnode] =
 128 {"(1-x)*(1-y)*(-1-x-y)/4",
 129  "(1+x)*(1-y)*(-1+x-y)/4",
 130  "(1+x)*(1+y)*(-1+x+y)/4",
 131  "(1-x)*(1+y)*(-1-x+y)/4",
 132  "(1-x*x)*(1-y)/2",
 133  "(1+x)*(1-y*y)/2",
 134  "(1-x*x)*(1+y)/2",
 135  "(1-x)*(1-y*y)/2"};
 136
 137 // QUAD9
 138 const static int QUAD9_dim = 2;
 139 const static int QUAD9_nnode = 9;
 140 static const int QUAD9_aster2med[QUAD9_nnode] =
 141 {0, 1, 2, 3, 4, 5, 6, 7, 8};
 142 static const char* QUAD9_varnames[QUAD9_dim] = {"x", "y"};
 143 static const char* QUAD9_functions[QUAD9_nnode] =
 144 {"x*y*(x-1)*(y-1)/4",
 145  "x*y*(x+1)*(y-1)/4",
 146  "x*y*(x+1)*(y+1)/4",
 147  "x*y*(x-1)*(y+1)/4",
 148  "(1-x*x)*y*(y-1)/2",
 149  "x*(x+1)*(1-y*y)/2",
 150  "(1-x*x)*y*(y+1)/2",
 151  "x*(x-1)*(1-y*y)/2",
 152  "(1-x*x)*(1-y*y)"};
 153
 154 // PENTA6
 155 const static int PENTA6_dim = 3;
 156 const static int PENTA6_nnode = 6;
 157 static const int PENTA6_aster2med[PENTA6_nnode] =
 158 {0, 1, 2, 3, 4, 5};
 159 static const char* PENTA6_varnames[PENTA6_dim] = {"x", "y", "z"};
 160 static const char* PENTA6_functions[PENTA6_nnode] =
 161 {"1/2*y*(1-x)",
 162  "1/2*z*(1-x)",
 163  "1/2*(1-y-z)*(1-x)",
 164  "1/2*y*(1+x)",
 165  "1/2*z*(1+x)",
 166  "1/2*(1-y-z)*(1+x)"};
 167
 168 // PENTA15
 169 const static int PENTA15_dim = 3;
 170 const static int PENTA15_nnode = 15;
 171 static const int PENTA15_aster2med[PENTA15_nnode] =
 172 {0, 1, 2, 3, 4, 5, 6, 7, 8, 12, 13, 14, 9, 10, 11};
 173 static const char* PENTA15_varnames[PENTA15_dim] = {"x", "y", "z"};
 174 static const char* PENTA15_functions[PENTA15_nnode] =
 175 {"y*(1-x)*(2*y-2-x)/2",
 176  "z*(1-x)*(2*z-2-x)/2",
 177  "(x-1)*(1-y-z)*(x+2*y+2*z)/2",
 178  "y*(1+x)*(2*y-2+x)/2",
 179  "z*(1+x)*(2*z-2+x)/2",
 180  "(-x-1)*(1-y-z)*(-x+2*y+2*z)/2",
 181  "2*y*z*(1-x)",
 182  "2*z*(1-y-z)*(1-x)",
 183  "2*y*(1-y-z)*(1-x)",
 184  "2*y*z*(1+x)",
 185  "2*z*(1-y-z)*(1+x)",
 186  "2*y*(1-y-z)*(1+x)",
 187  "y*(1-x*x)",
 188  "z*(1-x*x)",
 189  "(1-y-z)*(1-x*x)"};
 190
 191
 192 // PENTA18
 193 const static int PENTA18_dim = 3;
 194 const static int PENTA18_nnode = 18;
 195 const static int PENTA18_aster2med[PENTA18_nnode] =
 196 {0, 1, 2, 3, 4, 5, 6, 7, 8, 12, 13, 14, 9, 10, 11, 15, 16, 17};
 197 static const char* PENTA18_varnames[PENTA18_dim] = {"x", "y", "z"};
 198 static const char* PENTA18_functions[PENTA18_nnode] =
 199 {"x*y*(x−1)*(2*y−1)/2",
 200  "x*z*(x−1)*(2*z−1)/2",
 201  "x*(x−1)*(zy−1)*(2*z2*y−1)/2",
 202  "x*y*(x1)*(2*y−1)/2",
 203  "x*z*(x1)*(2*z−1)/2",
 204  "x*(x1)*(zy−1)*(2*z2*y−1)/2",
 205  "2*x*y*z*(x−1)",
 206  "−2*x*z*(x−1)*(zy−1)",
 207  "−2*x*y*(x−1)*(zy−1)",
 208  "2*x*y*z*(x1)",
 209  "−2*x*z*(x1)*(zy−1)",
 210  "−2*x*y*(x1)*(zy−1)",
 211  "y*(1−x*x)*(2*y−1)",
 212  "z*(1−x*x)*(2*z−1)",
 213  "(1−x*x)*(zy−1)*(2*z2*y−1)",
 214  "4*y*z*(1−x*x)",
 215  "4*z*(x−1)*(zy−1)",
 216  "4*y*(x−1)*(zy−1)"};
 217
 218 // HEXA8
 219 const static int HEXA8_dim = 3;
 220 const static int HEXA8_nnode = 8;
 221 static const int HEXA8_aster2med[HEXA8_nnode] =
 222 {0, 1, 2, 3, 4, 5, 6, 7};
 223 static const char* HEXA8_varnames[HEXA8_dim] = {"x", "y", "z"};
 224 static const char* HEXA8_functions[HEXA8_nnode] =
 225 {"1/8*(1-x)*(1-y)*(1-z)",
 226  "1/8*(1+x)*(1-y)*(1-z)",
 227  "1/8*(1+x)*(1+y)*(1-z)",
 228  "1/8*(1-x)*(1+y)*(1-z)",
 229  "1/8*(1-x)*(1-y)*(1+z)",
 230  "1/8*(1+x)*(1-y)*(1+z)",
 231  "1/8*(1+x)*(1+y)*(1+z)",
 232  "1/8*(1-x)*(1+y)*(1+z)"
 233    };
 234
 235 // HEXA20
 236 const static int HEXA20_dim = 3;
 237 const static int HEXA20_nnode = 20;
 238 static const int HEXA20_aster2med[HEXA20_nnode] =
 239 {0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 16, 17, 18, 19, 12, 13, 14, 15};
 240 static const char* HEXA20_varnames[HEXA20_dim] = {"x", "y", "z"};
 241 static const char* HEXA20_functions[HEXA20_nnode] =
 242 {"1/8*(1-x)*(1-y)*(1-z)*(-2-x-y-z)",
 243  "1/8*(1+x)*(1-y)*(1-z)*(-2+x-y-z)",
 244  "1/8*(1+x)*(1+y)*(1-z)*(-2+x+y-z)",
 245  "1/8*(1-x)*(1+y)*(1-z)*(-2-x+y-z)",
 246  "1/8*(1-x)*(1-y)*(1+z)*(-2-x-y+z)",
 247  "1/8*(1+x)*(1-y)*(1+z)*(-2+x-y+z)",
 248  "1/8*(1+x)*(1+y)*(1+z)*(-2+x+y+z)",
 249  "1/8*(1-x)*(1+y)*(1+z)*(-2-x+y+z)",
 250  "1/4*(1-x*x)*(1-y)*(1-z)",
 251  "1/4*(1-y*y)*(1+x)*(1-z)",
 252  "1/4*(1-x*x)*(1+y)*(1-z)",
 253  "1/4*(1-y*y)*(1-x)*(1-z)",
 254  "1/4*(1-x*x)*(1-y)*(1+z)",
 255  "1/4*(1-y*y)*(1+x)*(1+z)",
 256  "1/4*(1-x*x)*(1+y)*(1+z)",
 257  "1/4*(1-y*y)*(1-x)*(1+z)",
 258  "1/4*(1-z*z)*(1-x)*(1-y)",
 259  "1/4*(1-z*z)*(1+x)*(1-y)",
 260  "1/4*(1-z*z)*(1+x)*(1+y)",
 261  "1/4*(1-z*z)*(1-x)*(1+y)"
 262   };
 263 // HEXA27
 264 const static int HEXA27_dim = 3;
 265 const static int HEXA27_nnode = 27;
 266 static const int HEXA27_aster2med[HEXA27_nnode] =
 267 {0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 16, 17, 18, 19, 12, 13, 14, 15, 24, 22,
 268  21, 23, 20, 25, 26};
 269 static const char* HEXA27_varnames[HEXA27_dim] = {"x", "y", "z"};
 270 static const char* HEXA27_functions[HEXA27_nnode] =
 271 {"1/8*x*(x-1)*y*(y-1)*z*(z-1)",
 272  "1/8*x*(x+1)*y*(y-1)*z*(z-1)",
 273  "1/8*x*(x+1)*y*(y+1)*z*(z-1)",
 274  "1/8*x*(x-1)*y*(y+1)*z*(z-1)",
 275  "1/8*x*(x-1)*y*(y-1)*z*(z+1)",
 276  "1/8*x*(x+1)*y*(y-1)*z*(z+1)",
 277  "1/8*x*(x+1)*y*(y+1)*z*(z+1)",
 278  "1/8*x*(x-1)*y*(y+1)*z*(z+1)",
 279  "1/4*(1-x*x)*y*(y-1)*z*(z-1)",
 280  "1/4*x*(x+1)*(1-y*y)*z*(z-1)",
 281  "1/4*(1-x*x)*y*(y+1)*z*(z-1)",
 282  "1/4*x*(x-1)*(1-y*y)*z*(z-1)",
 283  "1/4*(1-x*x)*y*(y-1)*z*(z+1)",
 284  "1/4*x*(x+1)*(1-y*y)*z*(z+1)",
 285  "1/4*(1-x*x)*y*(y+1)*z*(z+1)",
 286  "1/4*x*(x-1)*(1-y*y)*z*(z+1)",
 287  "1/4*x*(x-1)*y*(y-1)*(1-z*z)",
 288  "1/4*x*(x+1)*y*(y-1)*(1-z*z)",
 289  "1/4*x*(x+1)*y*(y+1)*(1-z*z)",
 290  "1/4*x*(x-1)*y*(y+1)*(1-z*z)",
 291  "1/2*x*(x-1)*(1-y*y)*(1-z*z)",
 292  "1/2*x*(x+1)*(1-y*y)*(1-z*z)",
 293  "1/2*(1-x*x)*y*(y-1)*(1-z*z)",
 294  "1/2*(1-x*x)*y*(y+1)*(1-z*z)",
 295  "1/2*(1-x*x)*(1-y*y)*z*(z-1)",
 296  "1/2*(1-x*x)*(1-y*y)*z*(z+1)",
 297  "(1-x*x)*(1-y*y)*(1-z*z)"
 298   };
 299
 300 // TETRA4
 301 const static int TETRA4_dim = 3;
 302 const static int TETRA4_nnode = 4;
 303 static const int TETRA4_aster2med[TETRA4_nnode] =
 304 {0, 1, 2, 3};
 305 static const char* TETRA4_varnames[TETRA4_dim] = {"x", "y", "z"};
 306 static const char* TETRA4_functions[TETRA4_nnode] =
 307 {
 308   "y",
 309   "z",
 310   "1-x-y-z",
 311   "x"
 312 };
 313
 314 // TETRA10
 315 const static int TETRA10_dim = 3;
 316 const static int TETRA10_nnode = 10;
 317 static const int TETRA10_aster2med[TETRA10_nnode] =
 318 {0, 1, 2, 3, 4, 5, 6, 7, 8, 9};
 319 static const char* TETRA10_varnames[TETRA10_dim] = {"x", "y", "z"};
 320 static const char* TETRA10_functions[TETRA10_nnode] =
 321 {
 322   "y*(2*y-1)",
 323   "z*(2*z-1)",
 324   "(1-x-y-z)*(1-2*x-2*y-2*z)",
 325   "x*(2*x-1)",
 326   "4*y*z",
 327   "4*z*(1-x-y-z)",
 328   "4*y*(1-x-y-z)",
 329   "4*x*y",
 330   "4*x*z",
 331   "4*x*(1-x-y-z)"
 332 };
 333
 334 // PYRA5
 335 const static int PYRA5_dim = 3;
 336 const static int PYRA5_nnode = 5;
 337 static const int PYRA5_aster2med[PYRA5_nnode] =
 338 {0, 1, 2, 3, 4};
 339 static const char* PYRA5_varnames[PYRA5_dim] = {"x", "y", "z"};
 340 static const char* PYRA5_functions[PYRA5_nnode] =
 341 {
 342 "(-x+y+z-1)*(-x-y+z-1)/(4*(1-z))",
 343 "(-x-y+z-1)*( x-y+z-1)/(4*(1-z))",
 344 "( x-y+z-1)*( x+y+z-1)/(4*(1-z))",
 345 "( x+y+z-1)*(-x+y+z-1)/(4*(1-z))",
 346 "z"
 347 };
 348
 349 // PYRA13
 350 const static int PYRA13_dim = 3;
 351 const static int PYRA13_nnode = 13;
 352 static const int PYRA13_aster2med[PYRA13_nnode] =
 353 {0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12};
 354 static const char* PYRA13_varnames[PYRA13_dim] = {"x", "y", "z"};
 355 static const char* PYRA13_functions[PYRA13_nnode] =
 356 {
 357   "(-x+y+z-1)*(-x-y+z-1)*( x-1/2)/(2*(1-z))",
 358   "(-x-y+z-1)*( x-y+z-1)*( y-1/2)/(2*(1-z))",
 359   "( x-y+z-1)*( x+y+z-1)*(-x-1/2)/(2*(1-z))",
 360   "( x+y+z-1)*(-x+y+z-1)*(-y-1/2)/(2*(1-z))",
 361   "2*z*(z-1/2)",
 362   "-(-x+y+z-1)*(-x-y+z-1)*( x-y+z-1)/(2*(1-z))",
 363   "-(-x-y+z-1)*( x-y+z-1)*( x+y+z-1)/(2*(1-z))",
 364   "-( x-y+z-1)*( x+y+z-1)*(-x+y+z-1)/(2*(1-z))",
 365   "-( x+y+z-1)*(-x+y+z-1)*(-x-y+z-1)/(2*(1-z))",
 366   "(-x+y+z-1)*(-x-y+z-1)*z/(1-z)",
 367   "(-x-y+z-1)*( x-y+z-1)*z/(1-z)",
 368   "( x-y+z-1)*( x+y+z-1)*z/(1-z)",
 369   "( x+y+z-1)*(-x+y+z-1)*z/(1-z)"
 370 };
 371
 372 vtkCxxSetObjectMacro(vtkMedLocalization, ParentFile, vtkMedFile);
 373 vtkCxxSetObjectMacro(vtkMedLocalization, Interpolation, vtkMedInterpolation);
 374
 375 // vtkCxxRevisionMacro(vtkMedLocalization, "$Revision$")
 376 vtkStandardNewMacro(vtkMedLocalization)
 377
 378 vtkMedLocalization::vtkMedLocalization()
 379 {
 380   this->GeometryType = MED_NONE;
 381   this->NumberOfQuadraturePoint = 0;
 382   this->Weights = vtkDoubleArray::New();
 383   this->PointLocalCoordinates = vtkDoubleArray::New();
 384   this->QuadraturePointLocalCoordinates = vtkDoubleArray::New();
 385   this->ShapeFunction = vtkDoubleArray::New();
 386   this->Name = NULL;
 387   this->SectionName = NULL;
 388   this->InterpolationName = NULL;
 389   this->MedIterator = -1;
 390   this->ParentFile = NULL;
 391   this->SpaceDimension = 3;
 392   this->NumberOfCellInSection = 0;
 393   this->SectionGeometryType = MED_NONE;
 394   this->Interpolation = NULL;
 395   this->ShapeFunctionIsBuilt = 0;
 396 }
 397
 398 vtkMedLocalization::~vtkMedLocalization()
 399 {
 400   this->SetName(NULL);
 401   this->SetSectionName(NULL);
 402   this->SetInterpolationName(NULL);
 403   this->Weights->Delete();
 404   this->PointLocalCoordinates->Delete();
 405   this->QuadraturePointLocalCoordinates->Delete();
 406   this->ShapeFunction->Delete();
 407   this->SetInterpolation(NULL);
 408
 409 }
 410
 411 int vtkMedLocalization::GetSizeOfWeights()
 412 {
 413   return this->NumberOfQuadraturePoint;
 414 }
 415
 416 int vtkMedLocalization::GetSizeOfPointLocalCoordinates()
 417 {
 418   return vtkMedUtilities::GetNumberOfPoint(this->GeometryType)
 419       * vtkMedUtilities::GetDimension(this->GeometryType);
 420 }
 421
 422 int vtkMedLocalization::GetSizeOfQuadraturePointLocalCoordinates()
 423 {
 424   return this->NumberOfQuadraturePoint * vtkMedUtilities::GetDimension(
 425       this->GeometryType);
 426 }
 427
 428 int vtkMedLocalization::GetSizeOfShapeFunction()
 429 {
 430   return this->NumberOfQuadraturePoint * vtkMedUtilities::GetNumberOfPoint(
 431       this->GeometryType);
 432 }
 433
 434 void vtkMedLocalization::BuildShapeFunction()
 435 {
 436   if(this->ShapeFunctionIsBuilt)
 437     return;
 438
 439   if(this->Interpolation == NULL)
 440     {
 441     // If there is no interpolation given for this localization,
 442     // I build the default aster shape function
 443
 444     switch (this->GeometryType)
 445     {
 446       case MED_POINT1:
 447         BuildPoint1();
 448         return;
 449       case MED_SEG2:
 450         BuildAsterShapeFunction(SEG2_dim, SEG2_nnode,
 451                            (const int *) SEG2_aster2med,
 452                            (const char**)SEG2_varnames,
 453                            (const char**)SEG2_functions);
 454         break;
 455       case MED_SEG3:
 456         BuildAsterShapeFunction(SEG3_dim, SEG3_nnode,
 457                            (const int *) SEG3_aster2med,
 458                            (const char**)SEG3_varnames,
 459                            (const char**)SEG3_functions);
 460         break;
 461       case MED_SEG4:
 462         BuildAsterShapeFunction(SEG4_dim, SEG4_nnode,
 463                            (const int *) SEG4_aster2med,
 464                            (const char**)SEG4_varnames,
 465                            (const char**)SEG4_functions);
 466         break;
 467       case MED_TRIA3:
 468         BuildAsterShapeFunction(TRIA3_dim, TRIA3_nnode,
 469                            (const int *) TRIA3_aster2med,
 470                            (const char**)TRIA3_varnames,
 471                            (const char**)TRIA3_functions);
 472         break;
 473       case MED_TRIA6:
 474         BuildAsterShapeFunction(TRIA6_dim, TRIA6_nnode,
 475                            (const int *) TRIA6_aster2med,
 476                            (const char**)TRIA6_varnames,
 477                            (const char**)TRIA6_functions);
 478         break;
 479       case MED_TRIA7:
 480         BuildAsterShapeFunction(TRIA7_dim, TRIA7_nnode,
 481                            (const int *) TRIA7_aster2med,
 482                            (const char**)TRIA7_varnames,
 483                            (const char**)TRIA7_functions);
 484         break;
 485       case MED_QUAD4:
 486         BuildAsterShapeFunction(QUAD4_dim, QUAD4_nnode,
 487                            (const int *) QUAD4_aster2med,
 488                            (const char**)QUAD4_varnames,
 489                            (const char**)QUAD4_functions);
 490         break;
 491       case MED_QUAD8:
 492         BuildAsterShapeFunction(QUAD8_dim, QUAD8_nnode,
 493                            (const int *) QUAD8_aster2med,
 494                            (const char**)QUAD8_varnames,
 495                            (const char**)QUAD8_functions);
 496         break;
 497       case MED_QUAD9:
 498         BuildAsterShapeFunction(QUAD9_dim, QUAD9_nnode,
 499                            (const int *) QUAD9_aster2med,
 500                            (const char**)QUAD9_varnames,
 501                            (const char**)QUAD9_functions);
 502         break;
 503       case MED_HEXA8:
 504         BuildAsterShapeFunction(HEXA8_dim, HEXA8_nnode,
 505                            (const int *) HEXA8_aster2med,
 506                            (const char**)HEXA8_varnames,
 507                            (const char**)HEXA8_functions);
 508         break;
 509       case MED_HEXA20:
 510         BuildAsterShapeFunction(HEXA20_dim, HEXA20_nnode,
 511                            (const int *) HEXA20_aster2med,
 512                            (const char**)HEXA20_varnames,
 513                            (const char**)HEXA20_functions);
 514         break;
 515       case MED_HEXA27:
 516         BuildAsterShapeFunction(HEXA27_dim, HEXA27_nnode,
 517                            (const int *) HEXA27_aster2med,
 518                            (const char**)HEXA27_varnames,
 519                            (const char**)HEXA27_functions);
 520         break;
 521       case MED_TETRA4:
 522         BuildAsterShapeFunction(TETRA4_dim, TETRA4_nnode,
 523                            (const int *) TETRA4_aster2med,
 524                            (const char**)TETRA4_varnames,
 525                            (const char**)TETRA4_functions);
 526         break;
 527       case MED_TETRA10:
 528         BuildAsterShapeFunction(TETRA10_dim, TETRA10_nnode,
 529                            (const int *) TETRA10_aster2med,
 530                            (const char**)TETRA10_varnames,
 531                            (const char**)TETRA10_functions);
 532         break;
 533       case MED_PENTA6:
 534         BuildAsterShapeFunction(PENTA6_dim, PENTA6_nnode,
 535                            (const int *) PENTA6_aster2med,
 536                            (const char**)PENTA6_varnames,
 537                            (const char**)PENTA6_functions);
 538         break;
 539       case MED_PENTA15:
 540         BuildAsterShapeFunction(PENTA15_dim, PENTA15_nnode,
 541                            (const int *) PENTA15_aster2med,
 542                            (const char**)PENTA15_varnames,
 543                            (const char**)PENTA15_functions);
 544         break;
 545       case MED_PYRA5:
 546         BuildAsterShapeFunction(PYRA5_dim, PYRA5_nnode,
 547                            (const int *) PYRA5_aster2med,
 548                            (const char**)PYRA5_varnames,
 549                            (const char**)PYRA5_functions);
 550         break;
 551       case MED_PYRA13:
 552         BuildAsterShapeFunction(PYRA13_dim, PYRA13_nnode,
 553                            (const int *) PYRA13_aster2med,
 554                            (const char**)PYRA13_varnames,
 555                            (const char**)PYRA13_functions);
 556         break;
 557       default:
 558         vtkErrorMacro("ERROR in vtkMedLocalization::BuildShapeFunction. "
 559                       << this->GeometryType
 560                       << " : Cell geometry not supported !!! ");
 561         return;
 562       }
 563     }
 564   else
 565     {
 566     this->BuildShapeFunctionFromInterpolation();
 567     }
 568   this->ShapeFunctionIsBuilt = 1;
 569 }
 570
 571 void  vtkMedLocalization::BuildShapeFunctionFromInterpolation()
 572 {
 573   int nnodes = this->GeometryType % 100;
 574   int dim = this->GeometryType / 100;
 575   this->ShapeFunction->SetNumberOfValues(this->GetSizeOfShapeFunction());
 576
 577   int qpindex;
 578   int nodeindex;
 579
 580   vtkMedFraction* func;
 581
 582   switch(dim)
 583     {
 584     case 0 :
 585       this->ShapeFunction->SetValue(0, 1);
 586       break;
 587     default :
 588       for(qpindex=0; qpindex < this->NumberOfQuadraturePoint; qpindex++ )
 589         {
 590         double *coord = new double[dim];
 591         for(int dimid=0; dimid<dim; dimid++)
 592           {
 593           coord[dimid] = this->QuadraturePointLocalCoordinates
 594                          ->GetValue((qpindex * dim)+dimid);
 595           }
 596
 597         for(nodeindex=0; nodeindex < nnodes; nodeindex++)
 598           {
 599           func = this->Interpolation->GetBasisFunction(nodeindex);
 600           this->ShapeFunction->SetValue(
 601               qpindex*nnodes + nodeindex, func->Evaluate(coord));
 602           }
 603         }
 604     }
 605 }
 606
 607 void  vtkMedLocalization::BuildAsterShapeFunction(int dim,
 608                                  int nnodes,
 609                                  const int* aster2med,
 610                                  const char** varnames,
 611                                  const char** functions)
 612 {
 613   this->ShapeFunction->SetNumberOfValues(
 614       this->NumberOfQuadraturePoint * nnodes);
 615
 616   std::vector<vtkSmartPointer<vtkFunctionParser> > parsers;
 617   parsers.resize(nnodes);
 618   for(int nodeindex=0; nodeindex < nnodes; nodeindex++)
 619     {
 620     parsers[nodeindex] = vtkSmartPointer<vtkFunctionParser>::New();
 621     parsers[nodeindex]->SetFunction(functions[nodeindex]);
 622     }
 623
 624   for(int qpindex=0; qpindex < this->NumberOfQuadraturePoint; qpindex++ )
 625     {
 626
 627     for(int nodeindex=0; nodeindex < nnodes; nodeindex++)
 628       {
 629       int mednodeindex = aster2med[nodeindex];
 630       vtkFunctionParser* parser = parsers[mednodeindex];
 631       for(int dimid=0; dimid<dim; dimid++)
 632         {
 633         const char* varname = varnames[dimid];
 634         const double coord = this->QuadraturePointLocalCoordinates
 635                              ->GetValue((qpindex * dim)+dimid);
 636
 637         parser->SetScalarVariableValue(varname, coord);
 638         }
 639
 640       double w = parser->GetScalarResult();
 641
 642       this->ShapeFunction->SetValue(
 643           qpindex*nnodes + mednodeindex, w);
 644       }
 645     }
 646 }
 647
 648 void vtkMedLocalization::BuildPoint1()
 649 {
 650   this->Weights->SetNumberOfValues(1);
 651   this->ShapeFunction->SetNumberOfValues(1);
 652   this->Weights->SetValue(0, 1);
 653   this->ShapeFunction->SetValue(0, 1);
 654 }
 655
 656 void vtkMedLocalization::BuildCenter(med_geometry_type geometry)
 657 {
 658   this->GeometryType = geometry;
 659   this->NumberOfQuadraturePoint = 1;
 660   int npts = vtkMedUtilities::GetNumberOfPoint(this->GeometryType);
 661   this->ShapeFunction->SetNumberOfValues(npts);
 662   this->Weights->SetNumberOfValues(1);
 663   for (int i = 0; i < npts; i++)
 664     {
 665     this->ShapeFunction->SetValue(i, 1.0 / (double) npts);
 666     }
 667   this->Weights->SetValue(0, 1);
 668
 669 }
 670
 671 void vtkMedLocalization::BuildELNO(med_geometry_type geometry)
 672 {
 673   this->GeometryType = geometry;
 674   this->NumberOfQuadraturePoint = vtkMedUtilities::GetNumberOfPoint(geometry);
 675
 676   int np2 = this->NumberOfQuadraturePoint * this->NumberOfQuadraturePoint;
 677   this->ShapeFunction->SetNumberOfValues(np2);
 678   this->Weights->SetNumberOfValues(this->NumberOfQuadraturePoint);
 679
 680   for (int i = 0; i < np2; i++)
 681     {
 682     this->ShapeFunction->SetValue(i, 0);
 683     }
 684   for (int i = 0; i < this->NumberOfQuadraturePoint; i++)
 685     {
 686     this->ShapeFunction->SetValue(i + i * this->NumberOfQuadraturePoint, 1.0);
 687     }
 688   double w = 1.0 / (double) this->NumberOfQuadraturePoint;
 689   for (int i = 0; i < this->NumberOfQuadraturePoint; i++)
 690     {
 691     this->Weights->SetValue(i, w);
 692     }
 693 }
 694
 695 void vtkMedLocalization::PrintSelf(ostream& os, vtkIndent indent)
 696 {
 697   this->Superclass::PrintSelf(os, indent);
 698   PRINT_IVAR(os, indent, GeometryType);
 699   PRINT_IVAR(os, indent, NumberOfQuadraturePoint);
 700   PRINT_IVAR(os, indent, MedIterator);
 701 }