src/INTERPOLATION/MEDMEM_InterpolationTools.hxx

   1 // Copyright (C) 2005  OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN,
   2 // CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS
   3 //
   4 // This library is free software; you can redistribute it and/or
   5 // modify it under the terms of the GNU Lesser General Public
   6 // License as published by the Free Software Foundation; either
   7 // version 2.1 of the License.
   8 //
   9 // This library is distributed in the hope that it will be useful
  10 // but WITHOUT ANY WARRANTY; without even the implied warranty of
  11 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  12 // Lesser General Public License for more details.
  13 //
  14 // You should have received a copy of the GNU Lesser General Public
  15 // License along with this library; if not, write to the Free Software
  16 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
  17 //
  18 // See http://www.salome-platform.org/
  19 //
  20 #ifndef MEDMEM_INTERPOLATION_TOOLS_HXX
  21 #define MEDMEM_INTERPOLATION_TOOLS_HXX
  22
  23 #include "MEDMEM_define.hxx"
  24
  25 #define _TEMPLATE_ template <class CHAMP,class VALEURCHAMP,class NUAGENOEUD, class NOEUD, class NUAGEMAILLE>
  26 #define _PARAM_ CHAMP,VALEURCHAMP,NUAGENOEUD,NOEUD,NUAGEMAILLE
  27
  28 // ces macros définissent pour une face carrée plane la fonction projection sur cette face, non normalisée, la numérotation est dans le cas d'un hexaedre
  29 #define face2367(x,y,z) ((x6*(-y2 + y3) + x3*(y2 - y6) + x2*(-y3 + y6))*z - x6*y3*z2 + x3*y6*z2 + x6*y2*z3 - x2*y6*z3 - x3*y2*z6 + x2*y3*z6 + y*(x6*(z2 - z3) + x2*(z3 - z6) + x3*(-z2 + z6)) + x*(y6*(-z2 + z3) + y3*(z2 - z6) + y2*(-z3 + z6)))
  30 #define face4567(x,y,z) ((x6*(-y4 + y5) + x5*(y4 - y6) + x4*(-y5 + y6))*z - x6*y5*z4 + x5*y6*z4 + x6*y4*z5 - x4*y6*z5 - x5*y4*z6 + x4*y5*z6 + y*(x6*(z4 - z5) + x4*(z5 - z6) + x5*(-z4 + z6)) + x*(y6*(-z4 + z5) + y5*(z4 - z6) + y4*(-z5 + z6)))
  31 #define face1256(x,y,z) ((x5*(-y1 + y2) + x2*(y1 - y5) + x1*(-y2 + y5))*z - x5*y2*z1 + x2*y5*z1 + x5*y1*z2 - x1*y5*z2 - x2*y1*z5 + x1*y2*z5 + y*(x5*(z1 - z2) + x1*(z2 - z5) + x2*(-z1 + z5)) + x*(y5*(-z1 + z2) + y2*(z1 - z5) + y1*(-z2 + z5)))
  32 #define face0347(x,y,z) ((x4*(-y0 + y3) + x3*(y0 - y4) + x0*(-y3 + y4))*z - x4*y3*z0 + x3*y4*z0 + x4*y0*z3 - x0*y4*z3 - x3*y0*z4 + x0*y3*z4 + y*(x4*(z0 - z3) + x0*(z3 - z4) + x3*(-z0 + z4)) + x*(y4*(-z0 + z3) + y3*(z0 - z4) + y0*(-z3 + z4)))
  33 #define face0145(x,y,z) ((x4*(-y0 + y1) + x1*(y0 - y4) + x0*(-y1 + y4))*z - x4*y1*z0 + x1*y4*z0 + x4*y0*z1 - x0*y4*z1 - x1*y0*z4 + x0*y1*z4 + y*(x4*(z0 - z1) + x0*(z1 - z4) + x1*(-z0 + z4)) + x*(y4*(-z0 + z1) + y1*(z0 - z4) + y0*(-z1 + z4)))
  34 #define face0123(x,y,z) ((x2*(-y0 + y1) + x1*(y0 - y2) + x0*(-y1 + y2))*z - x2*y1*z0 + x1*y2*z0 + x2*y0*z1 - x0*y2*z1 - x1*y0*z2 + x0*y1*z2 + y*(x2*(z0 - z1) + x0*(z1 - z2) + x1*(-z0 + z2)) + x*(y2*(-z0 + z1) + y1*(z0 - z2) + y0*(-z1 + z2)))
  35 // des macros définissent pour une face triangulaire orientée vers l'extérieur de la maille la fonction de projection, non normalisée ( =(12^13).1M )
  36 #define face(x1,y1,z1,x2,y2,z2,x3,y3,z3,x,y,z)  ( ((y2-y1)*(z3-z1)-(z2-z1)*(y3-y1))*(x-x1)+((z2-z1)*(x3-x1)-(x2-x1)*(z3-z1))*(y-y1)+((x2-x1)*(y3-y1)-(y2-y1)*(x3-x1))*(z-z1) )
  37 #define projection(x1,y1,z1,x2,y2,z2,x3,y3,z3,x0,y0,z0,x,y,z) (face(x1,y1,z1,x2,y2,z2,x3,y3,z3,x,y,z)/face(x1,y1,z1,x2,y2,z2,x3,y3,z3,x0,y0,z0))
  38
  39 // DECLARATION
  40
  41 _TEMPLATE_ class Calcul_Interpolation
  42 {
  43 protected :
  44         NUAGENOEUD * noeuds;
  45         NUAGEMAILLE * mailles;
  46         CHAMP * champ;
  47 public :
  48         Calcul_Interpolation():noeuds(NULL),mailles(NULL),champ(NULL) {}
  49         Calcul_Interpolation(NUAGENOEUD * nn,NUAGEMAILLE * nm,CHAMP * c):noeuds(nn),mailles(nm),champ(c) {}
  50         ~Calcul_Interpolation() {}
  51         virtual VALEURCHAMP operator() (const NOEUD & n, int num_maille) {cerr<<"APPEL OPERATOR() DE LA CLASSE MERE CALCUL_INTERPOLATION => EXIT(-1)"<<endl;exit(-1);}
  52 };
  53
  54 _TEMPLATE_ class Calcul_Interpolation_P0;
  55
  56 _TEMPLATE_ class Calcul_Interpolation_Tria3;
  57 _TEMPLATE_ class Calcul_Interpolation_Tria6;
  58 _TEMPLATE_ class Calcul_Interpolation_Quad4;
  59 _TEMPLATE_ class Calcul_Interpolation_Quad8;
  60 _TEMPLATE_ class Calcul_Interpolation_Tetra4;
  61 _TEMPLATE_ class Calcul_Interpolation_Tetra10;
  62 _TEMPLATE_ class Calcul_Interpolation_Hexa8;
  63 _TEMPLATE_ class Calcul_Interpolation_Hexa20;
  64 _TEMPLATE_ class Calcul_Interpolation_Penta6;
  65 _TEMPLATE_ class Calcul_Interpolation_Penta15;
  66 _TEMPLATE_ class Calcul_Interpolation_Pyra5;
  67 _TEMPLATE_ class Calcul_Interpolation_Pyra13;
  68
  69 _TEMPLATE_ class Calcul_Hybride
  70 {
  71 protected :
  72         NUAGEMAILLE * mailles;
  73         map<int,Calcul_Interpolation<_PARAM_> *> fonctions;
  74 public :
  75         Calcul_Hybride():mailles(NULL) {}
  76         Calcul_Hybride(NUAGENOEUD * nn,NUAGEMAILLE * nm,CHAMP * c);
  77         ~Calcul_Hybride() {}
  78         VALEURCHAMP operator() (const NOEUD & n, int num_maille);
  79 };
  80
  81 //CODE
  82
  83 _TEMPLATE_ Calcul_Hybride<_PARAM_>::Calcul_Hybride(NUAGENOEUD * nn,NUAGEMAILLE * nm,CHAMP * c):mailles(nm)
  84         {
  85         fonctions[ MED_EN::MED_TRIA3  ]=new Calcul_Interpolation_Tria3  <_PARAM_>(nn,nm,c);
  86         fonctions[ MED_EN::MED_TRIA6  ]=new Calcul_Interpolation_Tria6  <_PARAM_>(nn,nm,c);
  87         fonctions[ MED_EN::MED_QUAD4  ]=new Calcul_Interpolation_Quad4  <_PARAM_>(nn,nm,c);
  88         fonctions[ MED_EN::MED_QUAD8  ]=new Calcul_Interpolation_Quad8  <_PARAM_>(nn,nm,c);
  89         fonctions[ MED_EN::MED_TETRA4 ]=new Calcul_Interpolation_Tetra4 <_PARAM_>(nn,nm,c);
  90         fonctions[ MED_EN::MED_TETRA10]=new Calcul_Interpolation_Tetra10<_PARAM_>(nn,nm,c);
  91         fonctions[ MED_EN::MED_HEXA8  ]=new Calcul_Interpolation_Hexa8  <_PARAM_>(nn,nm,c);
  92         fonctions[ MED_EN::MED_HEXA20 ]=new Calcul_Interpolation_Hexa20 <_PARAM_>(nn,nm,c);
  93         fonctions[ MED_EN::MED_PENTA6 ]=new Calcul_Interpolation_Penta6 <_PARAM_>(nn,nm,c);
  94         fonctions[ MED_EN::MED_PENTA15]=new Calcul_Interpolation_Penta15<_PARAM_>(nn,nm,c);
  95         fonctions[ MED_EN::MED_PYRA5  ]=new Calcul_Interpolation_Pyra5  <_PARAM_>(nn,nm,c);
  96         fonctions[ MED_EN::MED_PYRA13 ]=new Calcul_Interpolation_Pyra13 <_PARAM_>(nn,nm,c);
  97         }
  98
  99 _TEMPLATE_ VALEURCHAMP Calcul_Hybride<_PARAM_>::operator() (const NOEUD & n, int num_maille)
 100         {
 101         return fonctions[(*mailles)[num_maille].DONNE_TYPE_MED_MAILLE()]->operator()(n,num_maille);
 102         }
 103
 104 _TEMPLATE_ class Calcul_Interpolation_P0      : public Calcul_Interpolation<_PARAM_>
 105 {
 106 public : Calcul_Interpolation_P0(NUAGENOEUD * nn,NUAGEMAILLE * nm,CHAMP * c):Calcul_Interpolation<_PARAM_>(nn,nm,c) {}
 107 public : VALEURCHAMP operator() (const NOEUD & n, int num_maille)
 108         {
 109         return (*Calcul_Interpolation<_PARAM_>::champ)[num_maille];
 110         }
 111 };
 112 _TEMPLATE_ class Calcul_Interpolation_Tria3   : public Calcul_Interpolation<_PARAM_>
 113 {
 114 public : Calcul_Interpolation_Tria3(NUAGENOEUD * nn,NUAGEMAILLE * nm,CHAMP * c):Calcul_Interpolation<_PARAM_>(nn,nm,c) {}
 115 public : VALEURCHAMP operator() (const NOEUD & n, int num_maille)
 116         {
 117         int num0=(*Calcul_Interpolation<_PARAM_>::mailles)[num_maille][0];
 118         int num1=(*Calcul_Interpolation<_PARAM_>::mailles)[num_maille][1];
 119         int num2=(*Calcul_Interpolation<_PARAM_>::mailles)[num_maille][2];
 120
 121         double x0=(*Calcul_Interpolation<_PARAM_>::noeuds)[num0][0];
 122         double y0=(*Calcul_Interpolation<_PARAM_>::noeuds)[num0][1];
 123         double x1=(*Calcul_Interpolation<_PARAM_>::noeuds)[num1][0];
 124         double y1=(*Calcul_Interpolation<_PARAM_>::noeuds)[num1][1];
 125         double x2=(*Calcul_Interpolation<_PARAM_>::noeuds)[num2][0];
 126         double y2=(*Calcul_Interpolation<_PARAM_>::noeuds)[num2][1];
 127
 128         VALEURCHAMP v0=(*Calcul_Interpolation<_PARAM_>::champ)[num0];
 129         VALEURCHAMP v1=(*Calcul_Interpolation<_PARAM_>::champ)[num1];
 130         VALEURCHAMP v2=(*Calcul_Interpolation<_PARAM_>::champ)[num2];
 131
 132         double x=n[0];
 133         double y=n[1];
 134
 135         double lambda0=(y1-y2)*x+(x2-x1)*y+(x1*y2-x2*y1);
 136         double lambda1=(y2-y0)*x+(x0-x2)*y+(x2*y0-x0*y2);
 137         double lambda2=(y0-y1)*x+(x1-x0)*y+(x0*y1-x1*y0);
 138
 139         double delta = (x2-x1)*y0+(x0-x2)*y1+(x1-x0)*y2;
 140
 141         VALEURCHAMP retour(v0.SIZE());
 142
 143         retour=(1/delta)*(lambda0*v0+lambda1*v1+lambda2*v2);
 144
 145         return retour; //
 146         }
 147 };
 148 _TEMPLATE_ class Calcul_Interpolation_Tria6   : public Calcul_Interpolation<_PARAM_>
 149 {
 150 public : Calcul_Interpolation_Tria6(NUAGENOEUD * nn,NUAGEMAILLE * nm,CHAMP * c):Calcul_Interpolation<_PARAM_>(nn,nm,c) {}
 151 public : VALEURCHAMP operator() (const NOEUD & n, int num_maille)
 152         {
 153         // ON SUPPOSE IMPLICITEMENT QUE LES NOEUDS SUPPLEMENTAIRES SONT BIEN DES NOEUDS MILIEUX
 154         int num0 =(*Calcul_Interpolation<_PARAM_>::mailles)[num_maille][0];
 155         int num1 =(*Calcul_Interpolation<_PARAM_>::mailles)[num_maille][1];
 156         int num2 =(*Calcul_Interpolation<_PARAM_>::mailles)[num_maille][2];
 157         int num01=(*Calcul_Interpolation<_PARAM_>::mailles)[num_maille][3];
 158         int num12=(*Calcul_Interpolation<_PARAM_>::mailles)[num_maille][4];
 159         int num20=(*Calcul_Interpolation<_PARAM_>::mailles)[num_maille][5];
 160
 161         double x0=(*Calcul_Interpolation<_PARAM_>::noeuds)[num0][0];
 162         double y0=(*Calcul_Interpolation<_PARAM_>::noeuds)[num0][1];
 163         double x1=(*Calcul_Interpolation<_PARAM_>::noeuds)[num1][0];
 164         double y1=(*Calcul_Interpolation<_PARAM_>::noeuds)[num1][1];
 165         double x2=(*Calcul_Interpolation<_PARAM_>::noeuds)[num2][0];
 166         double y2=(*Calcul_Interpolation<_PARAM_>::noeuds)[num2][1];
 167
 168         VALEURCHAMP v0=(*Calcul_Interpolation<_PARAM_>::champ)[num0];
 169         VALEURCHAMP v1=(*Calcul_Interpolation<_PARAM_>::champ)[num1];
 170         VALEURCHAMP v2=(*Calcul_Interpolation<_PARAM_>::champ)[num2];
 171         VALEURCHAMP v01=(*Calcul_Interpolation<_PARAM_>::champ)[num01];
 172         VALEURCHAMP v12=(*Calcul_Interpolation<_PARAM_>::champ)[num12];
 173         VALEURCHAMP v20=(*Calcul_Interpolation<_PARAM_>::champ)[num20];
 174
 175         double x=n[0];
 176         double y=n[1];
 177
 178         double lambda0=(y1-y2)*x+(x2-x1)*y+(x1*y2-x2*y1);
 179         double lambda1=(y2-y0)*x+(x0-x2)*y+(x2*y0-x0*y2);
 180         double lambda2=(y0-y1)*x+(x1-x0)*y+(x0*y1-x1*y0);
 181
 182         double delta = (x2-x1)*y0+(x0-x2)*y1+(x1-x0)*y2;
 183
 184
 185         // VALEURCHAMP retour(v0.SIZE()); //
 186
 187         return  2*(lambda0*lambda0*v0+
 188                    lambda1*lambda1*v1+
 189                    lambda2*lambda2*v2+
 190                    2*(lambda0*lambda1*v01+
 191                       lambda1*lambda2*v12+
 192                       lambda2*lambda0*v20))/(delta*delta)+
 193                 (lambda0*v0+lambda1*v1+lambda2*v2)/(-delta);
 194
 195         // return retour; //
 196         }
 197 };
 198 _TEMPLATE_ class Calcul_Interpolation_Quad4   : public Calcul_Interpolation<_PARAM_>
 199 {
 200 public : Calcul_Interpolation_Quad4(NUAGENOEUD * nn,NUAGEMAILLE * nm,CHAMP * c):Calcul_Interpolation<_PARAM_>(nn,nm,c) {}
 201 public : VALEURCHAMP operator() (const NOEUD & n, int num_maille)
 202         {
 203         int num0=(*Calcul_Interpolation<_PARAM_>::mailles)[num_maille][0];
 204         int num1=(*Calcul_Interpolation<_PARAM_>::mailles)[num_maille][1];
 205         int num2=(*Calcul_Interpolation<_PARAM_>::mailles)[num_maille][2];
 206         int num3=(*Calcul_Interpolation<_PARAM_>::mailles)[num_maille][3];
 207
 208         double x0=(*Calcul_Interpolation<_PARAM_>::noeuds)[num0][0];
 209         double y0=(*Calcul_Interpolation<_PARAM_>::noeuds)[num0][1];
 210         double x1=(*Calcul_Interpolation<_PARAM_>::noeuds)[num1][0];
 211         double y1=(*Calcul_Interpolation<_PARAM_>::noeuds)[num1][1];
 212         double x2=(*Calcul_Interpolation<_PARAM_>::noeuds)[num2][0];
 213         double y2=(*Calcul_Interpolation<_PARAM_>::noeuds)[num2][1];
 214         double x3=(*Calcul_Interpolation<_PARAM_>::noeuds)[num3][0];
 215         double y3=(*Calcul_Interpolation<_PARAM_>::noeuds)[num3][1];
 216
 217         VALEURCHAMP v0=(*Calcul_Interpolation<_PARAM_>::champ)[num0];
 218         VALEURCHAMP v1=(*Calcul_Interpolation<_PARAM_>::champ)[num1];
 219         VALEURCHAMP v2=(*Calcul_Interpolation<_PARAM_>::champ)[num2];
 220         VALEURCHAMP v3=(*Calcul_Interpolation<_PARAM_>::champ)[num3];
 221
 222         double x=n[0];
 223         double y=n[1];
 224
 225
 226         double mu0=-((x3*(y1-y2)+x1*(y2-y3)+x2*(y3-y1))*x*y+(x2*y2*(y1-y3)+x3*(-y1+y2)*y3+x1*y1*(y3-y2))*x+(x2*x3*(y2-y3)+x1*(x2*(y1-y2)+x3*(y3-y1)))*y+(x2*x3*y1*(y3-y2)+x1*(x3*y2*(y1-y3)+x2*(y2-y1)*y3)));
 227         double mu1=(x0*(y2-y3)+x2*(y3-y0)+x3*(y0-y2))*x*y+(x3*y3*(y2-y0)+x0*(-y2+y3)*y0+x2*y2*(y0-y3))*x+(x3*x0*(y3-y0)+x2*(x3*(y2-y3)+x0*(y0-y2)))*y+(x2*x0*y2*(y0-y2)+x2*(x0*y2*(y2-y0)+x3*(y3-y2)*y0));
 228         double mu2=-((x1*(y3-y0)+x3*(y0-y1)+x0*(y1-y3))*x*y+(x0*y0*(y3-y1)+x1*(-y3+y0)*y1+x3*y3*(y1-y0))*x+(x0*x1*(y0-y1)+x3*(x0*(y3-y0)+x1*(y1-y3)))*y+(x3*x1*y3*(y1-y2)+x3*(x1*y2*(y3-y1)+x0*(y0-y3)*y1)));
 229         double mu3=(x2*(y0-y1)+x0*(y1-y2)+x1*(y2-y0))*x*y+(x1*y1*(y0-y2)+x2*(-y0+y1)*y2+x0*y0*(y2-y1))*x+(x1*x2*(y1-y2)+x0*(x1*(y0-y1)+x2*(y2-y0)))*y+(x0*x2*y0*(y2-y2)+x0*(x2*y2*(y0-y2)+x1*(y1-y0)*y2));
 230
 231         double delta=(y0-y2)*(y1-y3)*(x0*x2+x1*x3)-(x1*x2+x0*x3)*(y1-y2)*(y0-y3)-(x0*x1+x2*x3)*(y0-y1)*(y2-y3);
 232
 233         /*
 234         cout<<"  ### Pour ( "<<x<<" , "<<y<<" )"<<endl;
 235         cout<<"  ### delta = "<<delta<<endl;
 236         cout<<"  ### Mu0 = "<<mu0<<endl;
 237         cout<<"  ### Mu1 = "<<mu1<<endl;
 238         cout<<"  ### Mu2 = "<<mu2<<endl;
 239         cout<<"  ### Mu3 = "<<mu3<<endl;
 240         */
 241         VALEURCHAMP retour(v0.SIZE()); //
 242
 243         retour=(mu0*v0+mu1*v1+mu2*v2+mu3*v3)/delta;
 244
 245         return retour; //
 246         }
 247 };
 248 _TEMPLATE_ class Calcul_Interpolation_Quad8   : public Calcul_Interpolation<_PARAM_>
 249 {
 250 public : Calcul_Interpolation_Quad8(NUAGENOEUD * nn,NUAGEMAILLE * nm,CHAMP * c):Calcul_Interpolation<_PARAM_>(nn,nm,c) {}
 251 public : VALEURCHAMP operator() (const NOEUD & n, int num_maille)
 252         {
 253         cerr<<"Interpolation Q2 pas encore implémentée"<<endl;
 254         exit(-1);
 255         }
 256 };
 257 _TEMPLATE_ class Calcul_Interpolation_Tetra4  : public Calcul_Interpolation<_PARAM_>
 258 {
 259 public : Calcul_Interpolation_Tetra4(NUAGENOEUD * nn,NUAGEMAILLE * nm,CHAMP * c):Calcul_Interpolation<_PARAM_>(nn,nm,c) {}
 260 public : VALEURCHAMP operator() (const NOEUD & n, int num_maille)
 261         {
 262         int num0=(*Calcul_Interpolation<_PARAM_>::mailles)[num_maille][0];
 263         int num1=(*Calcul_Interpolation<_PARAM_>::mailles)[num_maille][1];
 264         int num2=(*Calcul_Interpolation<_PARAM_>::mailles)[num_maille][2];
 265         int num3=(*Calcul_Interpolation<_PARAM_>::mailles)[num_maille][3];
 266
 267         double x0=(*Calcul_Interpolation<_PARAM_>::noeuds)[num0][0];
 268         double y0=(*Calcul_Interpolation<_PARAM_>::noeuds)[num0][1];
 269         double z0=(*Calcul_Interpolation<_PARAM_>::noeuds)[num0][2];
 270         double x1=(*Calcul_Interpolation<_PARAM_>::noeuds)[num1][0];
 271         double y1=(*Calcul_Interpolation<_PARAM_>::noeuds)[num1][1];
 272         double z1=(*Calcul_Interpolation<_PARAM_>::noeuds)[num1][2];
 273         double x2=(*Calcul_Interpolation<_PARAM_>::noeuds)[num2][0];
 274         double y2=(*Calcul_Interpolation<_PARAM_>::noeuds)[num2][1];
 275         double z2=(*Calcul_Interpolation<_PARAM_>::noeuds)[num2][2];
 276         double x3=(*Calcul_Interpolation<_PARAM_>::noeuds)[num3][0];
 277         double y3=(*Calcul_Interpolation<_PARAM_>::noeuds)[num3][1];
 278         double z3=(*Calcul_Interpolation<_PARAM_>::noeuds)[num3][2];
 279
 280         VALEURCHAMP v0=(*Calcul_Interpolation<_PARAM_>::champ)[num0];
 281         VALEURCHAMP v1=(*Calcul_Interpolation<_PARAM_>::champ)[num1];
 282         VALEURCHAMP v2=(*Calcul_Interpolation<_PARAM_>::champ)[num2];
 283         VALEURCHAMP v3=(*Calcul_Interpolation<_PARAM_>::champ)[num3];
 284
 285         double x=n[0];
 286         double y=n[1];
 287         double z=n[2];
 288
 289         double lambda0=face(x1,y1,z1,x2,y2,z2,x3,y3,z3,x,y,z)/face(x1,y1,z1,x2,y2,z2,x3,y3,z3,x0,y0,z0);
 290         double lambda1=face(x0,y0,z0,x2,y2,z2,x3,y3,z3,x,y,z)/face(x0,y0,z0,x2,y2,z2,x3,y3,z3,x1,y1,z1);
 291         double lambda2=face(x1,y1,z1,x0,y0,z0,x3,y3,z3,x,y,z)/face(x1,y1,z1,x0,y0,z0,x3,y3,z3,x2,y2,z2);
 292         double lambda3=face(x1,y1,z1,x2,y2,z2,x0,y0,z0,x,y,z)/face(x1,y1,z1,x2,y2,z2,x0,y0,z0,x3,y3,z3);
 293
 294         VALEURCHAMP retour(v0.SIZE()); //
 295
 296         retour=(lambda0*v0+lambda1*v1+lambda2*v2+lambda3*v3);
 297
 298         return retour; //
 299         }
 300 };
 301 _TEMPLATE_ class Calcul_Interpolation_Tetra10 : public Calcul_Interpolation<_PARAM_>
 302 {
 303 public : Calcul_Interpolation_Tetra10(NUAGENOEUD * nn,NUAGEMAILLE * nm,CHAMP * c):Calcul_Interpolation<_PARAM_>(nn,nm,c) {}
 304 public : VALEURCHAMP operator() (const NOEUD & n, int num_maille)
 305         {
 306         // ON SUPPOSE IMPLICITEMENT QUE LES NOEUDS SUPPLEMENTAIRES SONT BIEN DES NOEUDS MILIEUX
 307         int num0 =(*Calcul_Interpolation<_PARAM_>::mailles)[num_maille][0];
 308         int num1 =(*Calcul_Interpolation<_PARAM_>::mailles)[num_maille][1];
 309         int num2 =(*Calcul_Interpolation<_PARAM_>::mailles)[num_maille][2];
 310         int num3 =(*Calcul_Interpolation<_PARAM_>::mailles)[num_maille][3];
 311         int num01=(*Calcul_Interpolation<_PARAM_>::mailles)[num_maille][4];
 312         int num02=(*Calcul_Interpolation<_PARAM_>::mailles)[num_maille][6];
 313         int num03=(*Calcul_Interpolation<_PARAM_>::mailles)[num_maille][7];
 314         int num12=(*Calcul_Interpolation<_PARAM_>::mailles)[num_maille][5];
 315         int num13=(*Calcul_Interpolation<_PARAM_>::mailles)[num_maille][8];
 316         int num23=(*Calcul_Interpolation<_PARAM_>::mailles)[num_maille][9];
 317
 318         double x0=(*Calcul_Interpolation<_PARAM_>::noeuds)[num0][0];double y0=(*Calcul_Interpolation<_PARAM_>::noeuds)[num0][1];double z0=(*Calcul_Interpolation<_PARAM_>::noeuds)[num0][2];
 319         double x1=(*Calcul_Interpolation<_PARAM_>::noeuds)[num1][0];double y1=(*Calcul_Interpolation<_PARAM_>::noeuds)[num1][1];double z1=(*Calcul_Interpolation<_PARAM_>::noeuds)[num1][2];
 320         double x2=(*Calcul_Interpolation<_PARAM_>::noeuds)[num2][0];double y2=(*Calcul_Interpolation<_PARAM_>::noeuds)[num2][1];double z2=(*Calcul_Interpolation<_PARAM_>::noeuds)[num2][2];
 321         double x3=(*Calcul_Interpolation<_PARAM_>::noeuds)[num3][0];double y3=(*Calcul_Interpolation<_PARAM_>::noeuds)[num3][1];double z3=(*Calcul_Interpolation<_PARAM_>::noeuds)[num3][2];
 322
 323         VALEURCHAMP v0=(*Calcul_Interpolation<_PARAM_>::champ)[num0];
 324         VALEURCHAMP v1=(*Calcul_Interpolation<_PARAM_>::champ)[num1];
 325         VALEURCHAMP v2=(*Calcul_Interpolation<_PARAM_>::champ)[num2];
 326         VALEURCHAMP v3=(*Calcul_Interpolation<_PARAM_>::champ)[num3];
 327         VALEURCHAMP v01=(*Calcul_Interpolation<_PARAM_>::champ)[num01];
 328         VALEURCHAMP v02=(*Calcul_Interpolation<_PARAM_>::champ)[num02];
 329         VALEURCHAMP v03=(*Calcul_Interpolation<_PARAM_>::champ)[num03];
 330         VALEURCHAMP v12=(*Calcul_Interpolation<_PARAM_>::champ)[num12];
 331         VALEURCHAMP v13=(*Calcul_Interpolation<_PARAM_>::champ)[num13];
 332         VALEURCHAMP v23=(*Calcul_Interpolation<_PARAM_>::champ)[num23];
 333
 334         double x=n[0];
 335         double y=n[1];
 336         double z=n[2];
 337
 338         double lambda0=face(x1,y1,z1,x2,y2,z2,x3,y3,z3,x,y,z)/face(x1,y1,z1,x2,y2,z2,x3,y3,z3,x0,y0,z0);
 339         double lambda1=face(x0,y0,z0,x2,y2,z2,x3,y3,z3,x,y,z)/face(x0,y0,z0,x2,y2,z2,x3,y3,z3,x1,y1,z1);
 340         double lambda2=face(x1,y1,z1,x0,y0,z0,x3,y3,z3,x,y,z)/face(x1,y1,z1,x0,y0,z0,x3,y3,z3,x2,y2,z2);
 341         double lambda3=face(x1,y1,z1,x2,y2,z2,x0,y0,z0,x,y,z)/face(x1,y1,z1,x2,y2,z2,x0,y0,z0,x3,y3,z3);
 342
 343         /*
 344         cout<<"  ### Pour ( "<<x<<" , "<<y<<" , "<<z<<" )"<<endl;
 345         cout<<"  ### Lambda0 = "<<lambda0<<endl;
 346         cout<<"  ### Lambda1 = "<<lambda1<<endl;
 347         cout<<"  ### Lambda2 = "<<lambda2<<endl;
 348         cout<<"  ### Lambda3 = "<<lambda3<<endl;
 349         cout<<"  ### Controle = "<<(lambda0+lambda1+lambda2+lambda3)<<endl;
 350         */
 351
 352         VALEURCHAMP retour(v0.SIZE()); //
 353
 354         retour=2*(lambda0*lambda0*v0+
 355                   lambda1*lambda1*v1+
 356                   lambda2*lambda2*v2+
 357                   lambda3*lambda3*v3+
 358                   2*(lambda0*lambda1*v01+
 359                      lambda0*lambda2*v02+
 360                      lambda0*lambda3*v03+
 361                      lambda1*lambda2*v12+
 362                      lambda1*lambda3*v13+
 363                      lambda2*lambda3*v23))+
 364                (-1.0)*(lambda0*v0+lambda1*v1+lambda2*v2+lambda3*v3);
 365
 366         return retour; //
 367         }
 368 };
 369 _TEMPLATE_ class Calcul_Interpolation_Hexa8   : public Calcul_Interpolation<_PARAM_>
 370 {
 371 public : Calcul_Interpolation_Hexa8(NUAGENOEUD * nn,NUAGEMAILLE * nm,CHAMP * c):Calcul_Interpolation<_PARAM_>(nn,nm,c) {}
 372 public : VALEURCHAMP operator() (const NOEUD & n, int num_maille)
 373         {
 374
 375         // Tres probablement numériquement mauvais, à revoir
 376
 377         int num0=(*Calcul_Interpolation<_PARAM_>::mailles)[num_maille][0];
 378         int num1=(*Calcul_Interpolation<_PARAM_>::mailles)[num_maille][1];
 379         int num2=(*Calcul_Interpolation<_PARAM_>::mailles)[num_maille][2];
 380         int num3=(*Calcul_Interpolation<_PARAM_>::mailles)[num_maille][3];
 381         int num4=(*Calcul_Interpolation<_PARAM_>::mailles)[num_maille][4];
 382         int num5=(*Calcul_Interpolation<_PARAM_>::mailles)[num_maille][5];
 383         int num6=(*Calcul_Interpolation<_PARAM_>::mailles)[num_maille][6];
 384         int num7=(*Calcul_Interpolation<_PARAM_>::mailles)[num_maille][7];
 385
 386
 387         double x0=(*Calcul_Interpolation<_PARAM_>::noeuds)[num0][0];double y0=(*Calcul_Interpolation<_PARAM_>::noeuds)[num0][1];double z0=(*Calcul_Interpolation<_PARAM_>::noeuds)[num0][2];
 388         double x1=(*Calcul_Interpolation<_PARAM_>::noeuds)[num1][0];double y1=(*Calcul_Interpolation<_PARAM_>::noeuds)[num1][1];double z1=(*Calcul_Interpolation<_PARAM_>::noeuds)[num1][2];
 389         double x2=(*Calcul_Interpolation<_PARAM_>::noeuds)[num2][0];double y2=(*Calcul_Interpolation<_PARAM_>::noeuds)[num2][1];double z2=(*Calcul_Interpolation<_PARAM_>::noeuds)[num2][2];
 390         double x3=(*Calcul_Interpolation<_PARAM_>::noeuds)[num3][0];double y3=(*Calcul_Interpolation<_PARAM_>::noeuds)[num3][1];double z3=(*Calcul_Interpolation<_PARAM_>::noeuds)[num3][2];
 391         double x4=(*Calcul_Interpolation<_PARAM_>::noeuds)[num4][0];double y4=(*Calcul_Interpolation<_PARAM_>::noeuds)[num4][1];double z4=(*Calcul_Interpolation<_PARAM_>::noeuds)[num4][2];
 392         double x5=(*Calcul_Interpolation<_PARAM_>::noeuds)[num5][0];double y5=(*Calcul_Interpolation<_PARAM_>::noeuds)[num5][1];double z5=(*Calcul_Interpolation<_PARAM_>::noeuds)[num5][2];
 393         double x6=(*Calcul_Interpolation<_PARAM_>::noeuds)[num6][0];double y6=(*Calcul_Interpolation<_PARAM_>::noeuds)[num6][1];double z6=(*Calcul_Interpolation<_PARAM_>::noeuds)[num6][2];
 394         double x7=(*Calcul_Interpolation<_PARAM_>::noeuds)[num7][0];double y7=(*Calcul_Interpolation<_PARAM_>::noeuds)[num7][1];double z7=(*Calcul_Interpolation<_PARAM_>::noeuds)[num7][2];
 395
 396
 397         VALEURCHAMP v0=(*Calcul_Interpolation<_PARAM_>::champ)[num0];
 398         VALEURCHAMP v1=(*Calcul_Interpolation<_PARAM_>::champ)[num1];
 399         VALEURCHAMP v2=(*Calcul_Interpolation<_PARAM_>::champ)[num2];
 400         VALEURCHAMP v3=(*Calcul_Interpolation<_PARAM_>::champ)[num3];
 401         VALEURCHAMP v4=(*Calcul_Interpolation<_PARAM_>::champ)[num4];
 402         VALEURCHAMP v5=(*Calcul_Interpolation<_PARAM_>::champ)[num5];
 403         VALEURCHAMP v6=(*Calcul_Interpolation<_PARAM_>::champ)[num6];
 404         VALEURCHAMP v7=(*Calcul_Interpolation<_PARAM_>::champ)[num7];
 405
 406         double x=n[0];
 407         double y=n[1];
 408         double z=n[2];
 409
 410         // double mu0=face1256(x,y,z)*face2367(x,y,z)*face4567(x,y,z)/(face1256(x0,y0,z0)*face2367(x0,y0,z0)*face4567(x0,y0,z0));
 411         // double mu1=face0347(x,y,z)*face2367(x,y,z)*face4567(x,y,z)/(face0347(x1,y1,z1)*face2367(x1,y1,z1)*face4567(x1,y1,z1));
 412         // double mu2=face0347(x,y,z)*face0145(x,y,z)*face4567(x,y,z)/(face0347(x2,y2,z2)*face0145(x2,y2,z2)*face4567(x2,y2,z2));
 413         // double mu3=face1256(x,y,z)*face0145(x,y,z)*face4567(x,y,z)/(face1256(x3,y3,z3)*face0145(x3,y3,z3)*face4567(x3,y3,z3));
 414         // double mu4=face1256(x,y,z)*face2367(x,y,z)*face0123(x,y,z)/(face1256(x4,y4,z4)*face2367(x4,y4,z4)*face0123(x4,y4,z4));
 415         // double mu5=face0347(x,y,z)*face2367(x,y,z)*face0123(x,y,z)/(face0347(x5,y5,z5)*face2367(x5,y5,z5)*face0123(x5,y5,z5));
 416         // double mu6=face0347(x,y,z)*face0145(x,y,z)*face0123(x,y,z)/(face0347(x6,y6,z6)*face0145(x6,y6,z6)*face0123(x6,y6,z6));
 417         // double mu7=face1256(x,y,z)*face0145(x,y,z)*face0123(x,y,z)/(face1256(x7,y7,z7)*face0145(x7,y7,z7)*face0123(x7,y7,z7));
 418
 419         double mu0=projection(x1,y1,z1,x2,y2,z2,x6,y6,z6,x0,y0,z0,x,y,z)*projection(x2,y2,z2,x3,y3,z3,x7,y7,z7,x0,y0,z0,x,y,z)*projection(x4,y4,z4,x5,y5,z5,x6,y6,z6,x0,y0,z0,x,y,z);
 420         double mu1=projection(x0,y0,z0,x3,y3,z3,x7,y7,z7,x1,y1,z1,x,y,z)*projection(x2,y2,z2,x3,y3,z3,x7,y7,z7,x1,y1,z1,x,y,z)*projection(x4,y4,z4,x5,y5,z5,x6,y6,z6,x1,y1,z1,x,y,z);
 421         double mu2=projection(x0,y0,z0,x3,y3,z3,x7,y7,z7,x2,y2,z2,x,y,z)*projection(x0,y0,z0,x1,y1,z1,x5,y5,z5,x2,y2,z2,x,y,z)*projection(x4,y4,z4,x5,y5,z5,x6,y6,z6,x2,y2,z2,x,y,z);
 422         double mu3=projection(x1,y1,z1,x2,y2,z2,x6,y6,z6,x3,y3,z3,x,y,z)*projection(x0,y0,z0,x1,y1,z1,x5,y5,z5,x3,y3,z3,x,y,z)*projection(x4,y4,z4,x6,y6,z6,x7,y7,z7,x3,y3,z3,x,y,z);
 423         double mu4=projection(x2,y2,z2,x6,y6,z6,x5,y5,z5,x4,y4,z4,x,y,z)*projection(x2,y2,z2,x3,y3,z3,x7,y7,z7,x4,y4,z4,x,y,z)*projection(x0,y0,z0,x1,y1,z1,x2,y2,z2,x4,y4,z4,x,y,z);
 424         double mu5=projection(x3,y3,z3,x7,y7,z7,x4,y4,z4,x5,y5,z5,x,y,z)*projection(x2,y2,z2,x3,y3,z3,x6,y6,z6,x5,y5,z5,x,y,z)*projection(x0,y0,z0,x1,y1,z1,x2,y2,z2,x5,y5,z5,x,y,z);
 425         double mu6=projection(x3,y3,z3,x7,y7,z7,x4,y4,z4,x6,y6,z6,x,y,z)*projection(x1,y1,z1,x5,y5,z5,x4,y4,z4,x6,y6,z6,x,y,z)*projection(x0,y0,z0,x1,y1,z1,x2,y2,z2,x6,y6,z6,x,y,z);
 426         double mu7=projection(x2,y2,z2,x6,y6,z6,x5,y5,z5,x7,y7,z7,x,y,z)*projection(x1,y1,z1,x5,y5,z5,x4,y4,z4,x7,y7,z7,x,y,z)*projection(x0,y0,z0,x1,y1,z1,x3,y3,z3,x7,y7,z7,x,y,z);
 427
 428         VALEURCHAMP retour(v0.SIZE()); //
 429
 430         retour=(mu0*v0+mu1*v1+mu2*v2+mu3*v3+mu4*v4+mu5*v5+mu6*v6+mu7*v7);
 431
 432         return retour; //
 433         }
 434 };
 435 _TEMPLATE_ class Calcul_Interpolation_Hexa20  : public Calcul_Interpolation<_PARAM_>
 436 {
 437 public : Calcul_Interpolation_Hexa20(NUAGENOEUD * nn,NUAGEMAILLE * nm,CHAMP * c):Calcul_Interpolation<_PARAM_>(nn,nm,c) {}
 438 public : VALEURCHAMP operator() (const NOEUD & n, int num_maille)
 439         {
 440         cerr<<"Interpolation H2 pasencore implémentée"<<endl;
 441         exit(-1);
 442         }
 443 };
 444 _TEMPLATE_ class Calcul_Interpolation_Penta6  : public Calcul_Interpolation<_PARAM_>
 445 {
 446 public : Calcul_Interpolation_Penta6(NUAGENOEUD * nn,NUAGEMAILLE * nm,CHAMP * c):Calcul_Interpolation<_PARAM_>(nn,nm,c) {}
 447 public : VALEURCHAMP operator() (const NOEUD & n, int num_maille)
 448         {
 449         int num0=(*Calcul_Interpolation<_PARAM_>::mailles)[num_maille][0];
 450         int num1=(*Calcul_Interpolation<_PARAM_>::mailles)[num_maille][1];
 451         int num2=(*Calcul_Interpolation<_PARAM_>::mailles)[num_maille][2];
 452         int num3=(*Calcul_Interpolation<_PARAM_>::mailles)[num_maille][3];
 453         int num4=(*Calcul_Interpolation<_PARAM_>::mailles)[num_maille][4];
 454         int num5=(*Calcul_Interpolation<_PARAM_>::mailles)[num_maille][5];
 455
 456         double x0=(*Calcul_Interpolation<_PARAM_>::noeuds)[num0][0];double y0=(*Calcul_Interpolation<_PARAM_>::noeuds)[num0][1];double z0=(*Calcul_Interpolation<_PARAM_>::noeuds)[num0][2];
 457         double x1=(*Calcul_Interpolation<_PARAM_>::noeuds)[num1][0];double y1=(*Calcul_Interpolation<_PARAM_>::noeuds)[num1][1];double z1=(*Calcul_Interpolation<_PARAM_>::noeuds)[num1][2];
 458         double x2=(*Calcul_Interpolation<_PARAM_>::noeuds)[num2][0];double y2=(*Calcul_Interpolation<_PARAM_>::noeuds)[num2][1];double z2=(*Calcul_Interpolation<_PARAM_>::noeuds)[num2][2];
 459         double x3=(*Calcul_Interpolation<_PARAM_>::noeuds)[num3][0];double y3=(*Calcul_Interpolation<_PARAM_>::noeuds)[num3][1];double z3=(*Calcul_Interpolation<_PARAM_>::noeuds)[num3][2];
 460         double x4=(*Calcul_Interpolation<_PARAM_>::noeuds)[num4][0];double y4=(*Calcul_Interpolation<_PARAM_>::noeuds)[num4][1];double z4=(*Calcul_Interpolation<_PARAM_>::noeuds)[num4][2];
 461         double x5=(*Calcul_Interpolation<_PARAM_>::noeuds)[num5][0];double y5=(*Calcul_Interpolation<_PARAM_>::noeuds)[num5][1];double z5=(*Calcul_Interpolation<_PARAM_>::noeuds)[num5][2];
 462
 463         VALEURCHAMP v0=(*Calcul_Interpolation<_PARAM_>::champ)[num0];
 464         VALEURCHAMP v1=(*Calcul_Interpolation<_PARAM_>::champ)[num1];
 465         VALEURCHAMP v2=(*Calcul_Interpolation<_PARAM_>::champ)[num2];
 466         VALEURCHAMP v3=(*Calcul_Interpolation<_PARAM_>::champ)[num3];
 467         VALEURCHAMP v4=(*Calcul_Interpolation<_PARAM_>::champ)[num4];
 468         VALEURCHAMP v5=(*Calcul_Interpolation<_PARAM_>::champ)[num5];
 469
 470         double x=n[0];
 471         double y=n[1];
 472         double z=n[2];
 473
 474         double mu0=face(x1,y1,z1,x2,y2,z2,x5,y5,z5,x,y,z)*face(x3,y3,z3,x4,y4,z4,x5,y5,z5,x,y,z)/(face(x1,y1,z1,x2,y2,z2,x5,y5,z5,x0,y0,z0)*face(x3,y3,z3,x4,y4,z4,x5,y5,z5,x0,y0,z0));
 475         double mu1=face(x0,y0,z0,x2,y2,z2,x5,y5,z5,x,y,z)*face(x3,y3,z3,x4,y4,z4,x5,y5,z5,x,y,z)/(face(x0,y0,z0,x2,y2,z2,x5,y5,z5,x1,y1,z1)*face(x3,y3,z3,x4,y4,z4,x5,y5,z5,x1,y1,z1));
 476         double mu2=face(x0,y0,z0,x1,y1,z1,x4,y4,z4,x,y,z)*face(x3,y3,z3,x4,y4,z4,x5,y5,z5,x,y,z)/(face(x0,y0,z0,x1,y1,z1,x4,y4,z4,x2,y2,z2)*face(x3,y3,z3,x4,y4,z4,x5,y5,z5,x2,y2,z2));
 477         double mu3=face(x1,y1,z1,x2,y2,z2,x5,y5,z5,x,y,z)*face(x0,y0,z0,x1,y1,z1,x2,y2,z2,x,y,z)/(face(x1,y1,z1,x2,y2,z2,x5,y5,z5,x3,y3,z3)*face(x0,y0,z0,x1,y1,z1,x2,y2,z2,x3,y3,z3));
 478         double mu4=face(x0,y0,z0,x2,y2,z2,x5,y5,z5,x,y,z)*face(x0,y0,z0,x1,y1,z1,x2,y2,z2,x,y,z)/(face(x0,y0,z0,x2,y2,z2,x5,y5,z5,x4,y4,z4)*face(x0,y0,z0,x1,y1,z1,x2,y2,z2,x4,y4,z4));
 479         double mu5=face(x0,y0,z0,x1,y1,z1,x4,y4,z4,x,y,z)*face(x0,y0,z0,x1,y1,z1,x2,y2,z2,x,y,z)/(face(x0,y0,z0,x1,y1,z1,x4,y4,z4,x5,y5,z5)*face(x0,y0,z0,x1,y1,z1,x2,y2,z2,x5,y5,z5));
 480
 481         VALEURCHAMP retour(v0.SIZE()); //
 482
 483         retour=(mu0*v0+mu1*v1+mu2*v2+mu3*v3+mu4*v4+mu5*v5);
 484
 485         return retour; //
 486         }
 487 };
 488 _TEMPLATE_ class Calcul_Interpolation_Penta15 : public Calcul_Interpolation<_PARAM_>
 489 {
 490 public : Calcul_Interpolation_Penta15(NUAGENOEUD * nn,NUAGEMAILLE * nm,CHAMP * c):Calcul_Interpolation<_PARAM_>(nn,nm,c) {}
 491 public : VALEURCHAMP operator() (const NOEUD & n, int num_maille)
 492         {
 493         cerr<<"Interpolation Pe2 pasencore implémentée"<<endl;
 494         exit(-1);
 495         }
 496 };
 497 _TEMPLATE_ class Calcul_Interpolation_Pyra5   : public Calcul_Interpolation<_PARAM_>
 498 {
 499 public : Calcul_Interpolation_Pyra5(NUAGENOEUD * nn,NUAGEMAILLE * nm,CHAMP * c):Calcul_Interpolation<_PARAM_>(nn,nm,c) {}
 500 public : VALEURCHAMP operator() (const NOEUD & n, int num_maille)
 501         {
 502         // NON TESTE
 503         int num0=(*Calcul_Interpolation<_PARAM_>::mailles)[num_maille][0];
 504         int num1=(*Calcul_Interpolation<_PARAM_>::mailles)[num_maille][1];
 505         int num2=(*Calcul_Interpolation<_PARAM_>::mailles)[num_maille][2];
 506         int num3=(*Calcul_Interpolation<_PARAM_>::mailles)[num_maille][3];
 507         int num4=(*Calcul_Interpolation<_PARAM_>::mailles)[num_maille][4];
 508
 509         double x0=(*Calcul_Interpolation<_PARAM_>::noeuds)[num0][0];double y0=(*Calcul_Interpolation<_PARAM_>::noeuds)[num0][1];double z0=(*Calcul_Interpolation<_PARAM_>::noeuds)[num0][2];
 510         double x1=(*Calcul_Interpolation<_PARAM_>::noeuds)[num1][0];double y1=(*Calcul_Interpolation<_PARAM_>::noeuds)[num1][1];double z1=(*Calcul_Interpolation<_PARAM_>::noeuds)[num1][2];
 511         double x2=(*Calcul_Interpolation<_PARAM_>::noeuds)[num2][0];double y2=(*Calcul_Interpolation<_PARAM_>::noeuds)[num2][1];double z2=(*Calcul_Interpolation<_PARAM_>::noeuds)[num2][2];
 512         double x3=(*Calcul_Interpolation<_PARAM_>::noeuds)[num3][0];double y3=(*Calcul_Interpolation<_PARAM_>::noeuds)[num3][1];double z3=(*Calcul_Interpolation<_PARAM_>::noeuds)[num3][2];
 513         double x4=(*Calcul_Interpolation<_PARAM_>::noeuds)[num4][0];double y4=(*Calcul_Interpolation<_PARAM_>::noeuds)[num4][1];double z4=(*Calcul_Interpolation<_PARAM_>::noeuds)[num4][2];
 514
 515         VALEURCHAMP v0=(*Calcul_Interpolation<_PARAM_>::champ)[num0];
 516         VALEURCHAMP v1=(*Calcul_Interpolation<_PARAM_>::champ)[num1];
 517         VALEURCHAMP v2=(*Calcul_Interpolation<_PARAM_>::champ)[num2];
 518         VALEURCHAMP v3=(*Calcul_Interpolation<_PARAM_>::champ)[num3];
 519         VALEURCHAMP v4=(*Calcul_Interpolation<_PARAM_>::champ)[num4];
 520
 521         double x=n[0];
 522         double y=n[1];
 523         double z=n[2];
 524
 525         double mu0=face(x1,y1,z1,x2,y2,z2,x4,y4,z4,x,y,z)*face(x2,y2,z2,x3,y3,z3,x4,y4,z4,x,y,z)/(face(x1,y1,z1,x2,y2,z2,x4,y4,z4,x0,y0,z0)*face(x2,y2,z2,x3,y3,z3,x4,y4,z4,x0,y0,z0));
 526         double mu1=face(x0,y0,z0,x3,y3,z3,x4,y4,z4,x,y,z)*face(x2,y2,z2,x3,y3,z3,x4,y4,z4,x,y,z)/(face(x0,y0,z0,x3,y3,z3,x4,y4,z4,x1,y1,z1)*face(x2,y2,z2,x3,y3,z3,x4,y4,z4,x1,y1,z1));
 527         double mu2=face(x0,y0,z0,x1,y1,z1,x4,y4,z4,x,y,z)*face(x0,y0,z0,x3,y3,z3,x4,y4,z4,x,y,z)/(face(x0,y0,z0,x1,y1,z1,x4,y4,z4,x2,y2,z2)*face(x0,y0,z0,x3,y3,z3,x4,y4,z4,x2,y2,z2));
 528         double mu3=face(x0,y0,z0,x1,y1,z1,x4,y4,z4,x,y,z)*face(x1,y1,z1,x2,y2,z2,x4,y4,z4,x,y,z)/(face(x0,y0,z0,x1,y1,z1,x4,y4,z4,x3,y3,z3)*face(x1,y1,z1,x2,y2,z2,x4,y4,z4,x3,y3,z3));
 529         double mu4=face(x0,y0,z0,x1,y1,z1,x2,y2,z2,x,y,z)/face(x0,y0,z0,x1,y1,z1,x2,y2,z2,x4,y4,z4);
 530
 531         VALEURCHAMP retour(v0.SIZE()); //
 532
 533         retour=(mu0*v0+mu1*v1+mu2*v2+mu3*v3+mu4*v4);
 534
 535         return retour; //
 536         }
 537 };
 538 _TEMPLATE_ class Calcul_Interpolation_Pyra13  : public Calcul_Interpolation<_PARAM_>
 539 {
 540 public : Calcul_Interpolation_Pyra13(NUAGENOEUD * nn,NUAGEMAILLE * nm,CHAMP * c):Calcul_Interpolation<_PARAM_>(nn,nm,c) {}
 541 public : VALEURCHAMP operator() (const NOEUD & n, int num_maille)
 542         {
 543         cerr<<"Interpolation Py2 pasencore implémentée"<<endl;
 544         exit(-1);
 545         }
 546 };
 547
 548 #undef _TEMPLATE_
 549 #undef _PARAM_
 550
 551 #endif