mtu/src/fem_element3.cpp

//####//------------------------------------------------------------
//####//------------------------------------------------------------
//####//  MAGiC
//####//  Jean Christophe Cuilliere et Vincent FRANCOIS
//####//  Departement de Genie Mecanique - UQTR
//####//------------------------------------------------------------
//####//  MAGIC est un projet de recherche de l equipe ERICCA
//####//  du departement de genie mecanique de l Universite du Quebec a Trois Rivieres
//####//  http://www.uqtr.ca/ericca
//####//  http://www.uqtr.ca/
//####//------------------------------------------------------------
//####//------------------------------------------------------------
//####//
//####//       fem_element3.cpp
//####//
//####//------------------------------------------------------------
//####//------------------------------------------------------------
//####//    COPYRIGHT 2000-2024
//####//  jeu 13 jun 2024 11:58:53 EDT
//####//------------------------------------------------------------
//####//------------------------------------------------------------


#include "gestionversion.h"
#include <math.h>
#include "fem_element3.h"
#include "xfem_element3.h"
#include "fem_noeud.h"
#include "fem_maillage.h"
#include "mg_element_maillage.h"
#include "ot_boite_3d.h"
#include "ot_mathematique.h"

FEM_ELEMENT3::FEM_ELEMENT3(unsigned long num,class MG_ELEMENT_MAILLAGE* mai):FEM_ELEMENT_MAILLAGE(num,mai)
{
}
FEM_ELEMENT3::FEM_ELEMENT3(unsigned long num,class MG_ELEMENT_TOPOLOGIQUE* topo):FEM_ELEMENT_MAILLAGE(num,topo)
{
}
FEM_ELEMENT3::FEM_ELEMENT3(unsigned long num,class MG_ELEMENT_TOPOLOGIQUE* topo,class MG_ELEMENT_MAILLAGE* mai):FEM_ELEMENT_MAILLAGE(num,topo,mai)
{
}

FEM_ELEMENT3::FEM_ELEMENT3(class MG_ELEMENT_MAILLAGE* mai):FEM_ELEMENT_MAILLAGE(mai)
{
}

FEM_ELEMENT3::FEM_ELEMENT3(class MG_ELEMENT_TOPOLOGIQUE* topo):FEM_ELEMENT_MAILLAGE(topo)
{
}
FEM_ELEMENT3::FEM_ELEMENT3(class MG_ELEMENT_TOPOLOGIQUE* topo,class MG_ELEMENT_MAILLAGE* mai):FEM_ELEMENT_MAILLAGE(topo,mai)
{
}


FEM_ELEMENT3::FEM_ELEMENT3(FEM_ELEMENT3& mdd):FEM_ELEMENT_MAILLAGE(mdd)
{
}


FEM_ELEMENT3::~FEM_ELEMENT3()
{
}

bool FEM_ELEMENT3::get_param_element_fini(double *xyz,double *uvw)
{
OT_MATRICE_3D mat;
OT_VECTEUR_3D vec;
uvw[0]=0.;
uvw[1]=0.;
uvw[2]=0.;

int compteur=0;
int ok=0;
while (ok!=1)
{
ok=0;
int cas=compteur/100-1;
int base=4;
if (compteur!=0)
    {
    int i=cas%base;  
    int j=(cas/base)%base;  
    int k=(cas/base/base)%base; 
    uvw[0]=-1.+2.*i/(base-1);
    uvw[1]=-1.+2.*j/(base-1);
    uvw[2]=-1.+2.*k/(base-1);
    }
if (cas>base*base*base) return false;
while (ok==0)
{
mat(0,0)=0.;mat(0,1)=0.;mat(0,2)=0.;
mat(1,0)=0.;mat(1,1)=0.;mat(1,2)=0.;
mat(2,0)=0.;mat(2,1)=0.;mat(2,2)=0.;
vec(0)=xyz[0];vec(1)=xyz[1];vec(2)=xyz[2];

for (int i=0;i<get_nb_fem_noeud();i++)
  {
  mat(0,0)=mat(0,0)+get_fonction_derive_interpolation(i+1,1,uvw)*get_fem_noeud(i)->get_x();
  mat(0,1)=mat(0,1)+get_fonction_derive_interpolation(i+1,2,uvw)*get_fem_noeud(i)->get_x();
  mat(0,2)=mat(0,2)+get_fonction_derive_interpolation(i+1,3,uvw)*get_fem_noeud(i)->get_x();
  mat(1,0)=mat(1,0)+get_fonction_derive_interpolation(i+1,1,uvw)*get_fem_noeud(i)->get_y();
  mat(1,1)=mat(1,1)+get_fonction_derive_interpolation(i+1,2,uvw)*get_fem_noeud(i)->get_y();
  mat(1,2)=mat(1,2)+get_fonction_derive_interpolation(i+1,3,uvw)*get_fem_noeud(i)->get_y();
  mat(2,0)=mat(2,0)+get_fonction_derive_interpolation(i+1,1,uvw)*get_fem_noeud(i)->get_z();
  mat(2,1)=mat(2,1)+get_fonction_derive_interpolation(i+1,2,uvw)*get_fem_noeud(i)->get_z();
  mat(2,2)=mat(2,2)+get_fonction_derive_interpolation(i+1,3,uvw)*get_fem_noeud(i)->get_z();
  vec(0)=vec(0)-get_fonction_interpolation(i+1,uvw)*get_fem_noeud(i)->get_x();
  vec(1)=vec(1)-get_fonction_interpolation(i+1,uvw)*get_fem_noeud(i)->get_y();
  vec(2)=vec(2)-get_fonction_interpolation(i+1,uvw)*get_fem_noeud(i)->get_z();
  }
double det=mat.get_determinant();
if (fabs(det)<1e-12) ok=3;
OT_MATRICE_3D mat1(vec,mat.get_vecteur2(),mat.get_vecteur3());
OT_MATRICE_3D mat2(mat.get_vecteur1(),vec,mat.get_vecteur3());
OT_MATRICE_3D mat3(mat.get_vecteur1(),mat.get_vecteur2(),vec);
double d1=mat1.get_determinant()/det;
double d2=mat2.get_determinant()/det;
double d3=mat3.get_determinant()/det;
compteur++;
if (compteur%100==0) ok=2;
if ((fabs(d1)<1e-8)&&(fabs(d2)<1e-8)&&(fabs(d3)<1e-8))  ok=1;
uvw[0]=uvw[0]+d1;
uvw[1]=uvw[1]+d2;
uvw[2]=uvw[2]+d3;
}
}
return true; 
}


void FEM_ELEMENT3::get_interpolation_xyz(double* uvw, double* xyz)
{
xyz[0]=0.;xyz[1]=0.;xyz[2]=0;
for (int i=0;i<get_nb_fem_noeud();i++)
  {
  xyz[0]=xyz[0]+get_fonction_interpolation(i+1,uvw)*get_fem_noeud(i)->get_x();
  xyz[1]=xyz[1]+get_fonction_interpolation(i+1,uvw)*get_fem_noeud(i)->get_y();
  xyz[2]=xyz[2]+get_fonction_interpolation(i+1,uvw)*get_fem_noeud(i)->get_z();
  }
}

double FEM_ELEMENT3::get_jacobien(double* jac,double *uv,double unite)
{   
    int nb=get_nb_fem_noeud();
      
    OT_MATRICE_3D jacobien;
    for (int i=0;i<3;i++)
    {
     jacobien(i,0)=0.;jacobien(i,1)=0.;jacobien(i,2)=0.;
      for (int k=0;k<nb;k++)
          {
      double valderiv=get_fonction_derive_interpolation(k+1,i+1,uv);
      double *xyz=get_fem_noeud(k)->get_coord();
          jacobien(i,0)=jacobien(i,0)+valderiv*xyz[0]*unite;
          jacobien(i,1)=jacobien(i,1)+valderiv*xyz[1]*unite;
          jacobien(i,2)=jacobien(i,2)+valderiv*xyz[2]*unite;
          }
    }
    jac[0]=jacobien(0,0);
    jac[1]=jacobien(0,1);
    jac[2]=jacobien(0,2);
    
    jac[3]=jacobien(1,0);
    jac[4]=jacobien(1,1);
    jac[5]=jacobien(1,2);
    
    jac[6]=jacobien(2,0);
    jac[7]=jacobien(2,1);
    jac[8]=jacobien(2,2);
    
    double det=jacobien.get_determinant();
    return det;
}

void FEM_ELEMENT3::get_inverse_jacob(double* j,double *uv,double unite)
{
    double jac[9];
    get_jacobien(jac,uv,unite);
    OT_MATRICE_3D J;
    J(0,0)=jac[0];
    J(0,1)=jac[1];
    J(0,2)=jac[2];
    J(1,0)=jac[3];
    J(1,1)=jac[4];
    J(1,2)=jac[5];
    J(2,0)=jac[6];
    J(2,1)=jac[7];
    J(2,2)=jac[8];
    OT_MATRICE_3D j_i=J.inverse();
  
    
    j[0]=j_i(0,0);
    j[1]=j_i(0,1);
    j[2]=j_i(0,2);
    
    j[3]=j_i(1,0);
    j[4]=j_i(1,1);
    j[5]=j_i(1,2);
    
    j[6]=j_i(2,0);
    j[7]=j_i(2,1);
    j[8]=j_i(2,2);
}
Revision:	1158
Committed:	Thu Jun 13 22:18:49 2024 UTC (11 months, 1 week ago) by francois
File size:	6209 byte(s)
Log Message:	compatibilité Ubuntu 22.04 Suppression des refeences à Windows Ajout d'une banière
#	User	Rev	Content
1	francois	1158	//####//------------------------------------------------------------
2			//####//------------------------------------------------------------
3			//####// MAGiC
4			//####// Jean Christophe Cuilliere et Vincent FRANCOIS
5			//####// Departement de Genie Mecanique - UQTR
6			//####//------------------------------------------------------------
7			//####// MAGIC est un projet de recherche de l equipe ERICCA
8			//####// du departement de genie mecanique de l Universite du Quebec a Trois Rivieres
9			//####// http://www.uqtr.ca/ericca
10			//####// http://www.uqtr.ca/
11			//####//------------------------------------------------------------
12			//####//------------------------------------------------------------
13			//####//
14			//####// fem_element3.cpp
15			//####//
16			//####//------------------------------------------------------------
17			//####//------------------------------------------------------------
18			//####// COPYRIGHT 2000-2024
19			//####// jeu 13 jun 2024 11:58:53 EDT
20			//####//------------------------------------------------------------
21			//####//------------------------------------------------------------
22	francois	283
23
24			#include "gestionversion.h"
25			#include <math.h>
26	francois	309	#include "fem_element3.h"
27	francois	339	#include "xfem_element3.h"
28	francois	283	#include "fem_noeud.h"
29	francois	339	#include "fem_maillage.h"
30	francois	283	#include "mg_element_maillage.h"
31	francois	481	#include "ot_boite_3d.h"
32	francois	628	#include "ot_mathematique.h"
33	francois	283
34	francois	399	FEM_ELEMENT3::FEM_ELEMENT3(unsigned long num,class MG_ELEMENT_MAILLAGE* mai):FEM_ELEMENT_MAILLAGE(num,mai)
35	francois	283	{
36			}
37	francois	399	FEM_ELEMENT3::FEM_ELEMENT3(unsigned long num,class MG_ELEMENT_TOPOLOGIQUE* topo):FEM_ELEMENT_MAILLAGE(num,topo)
38	francois	378	{
39			}
40	francois	399	FEM_ELEMENT3::FEM_ELEMENT3(unsigned long num,class MG_ELEMENT_TOPOLOGIQUE* topo,class MG_ELEMENT_MAILLAGE* mai):FEM_ELEMENT_MAILLAGE(num,topo,mai)
41	francois	378	{
42			}
43	francois	283
44	francois	399	FEM_ELEMENT3::FEM_ELEMENT3(class MG_ELEMENT_MAILLAGE* mai):FEM_ELEMENT_MAILLAGE(mai)
45	francois	283	{
46			}
47
48	francois	399	FEM_ELEMENT3::FEM_ELEMENT3(class MG_ELEMENT_TOPOLOGIQUE* topo):FEM_ELEMENT_MAILLAGE(topo)
49	francois	378	{
50			}
51	francois	399	FEM_ELEMENT3::FEM_ELEMENT3(class MG_ELEMENT_TOPOLOGIQUE* topo,class MG_ELEMENT_MAILLAGE* mai):FEM_ELEMENT_MAILLAGE(topo,mai)
52	francois	378	{
53			}
54	francois	283
55	francois	378
56	francois	399	FEM_ELEMENT3::FEM_ELEMENT3(FEM_ELEMENT3& mdd):FEM_ELEMENT_MAILLAGE(mdd)
57	francois	283	{
58			}
59
60
61
62	francois	309	FEM_ELEMENT3::~FEM_ELEMENT3()
63	francois	283	{
64			}
65
66	francois	635	bool FEM_ELEMENT3::get_param_element_fini(double xyz,double uvw)
67	francois	628	{
68			OT_MATRICE_3D mat;
69			OT_VECTEUR_3D vec;
70			uvw[0]=0.;
71			uvw[1]=0.;
72			uvw[2]=0.;
73	francois	283
74	francois	684	int compteur=0;
75	francois	628	int ok=0;
76	francois	684	while (ok!=1)
77			{
78			ok=0;
79			int cas=compteur/100-1;
80			int base=4;
81			if (compteur!=0)
82			{
83			int i=cas%base;
84			int j=(cas/base)%base;
85			int k=(cas/base/base)%base;
86			uvw[0]=-1.+2.*i/(base-1);
87			uvw[1]=-1.+2.*j/(base-1);
88			uvw[2]=-1.+2.*k/(base-1);
89			}
90			if (cas>basebasebase) return false;
91	francois	628	while (ok==0)
92			{
93			mat(0,0)=0.;mat(0,1)=0.;mat(0,2)=0.;
94			mat(1,0)=0.;mat(1,1)=0.;mat(1,2)=0.;
95			mat(2,0)=0.;mat(2,1)=0.;mat(2,2)=0.;
96			vec(0)=xyz[0];vec(1)=xyz[1];vec(2)=xyz[2];
97
98			for (int i=0;i<get_nb_fem_noeud();i++)
99			{
100			mat(0,0)=mat(0,0)+get_fonction_derive_interpolation(i+1,1,uvw)*get_fem_noeud(i)->get_x();
101			mat(0,1)=mat(0,1)+get_fonction_derive_interpolation(i+1,2,uvw)*get_fem_noeud(i)->get_x();
102			mat(0,2)=mat(0,2)+get_fonction_derive_interpolation(i+1,3,uvw)*get_fem_noeud(i)->get_x();
103			mat(1,0)=mat(1,0)+get_fonction_derive_interpolation(i+1,1,uvw)*get_fem_noeud(i)->get_y();
104			mat(1,1)=mat(1,1)+get_fonction_derive_interpolation(i+1,2,uvw)*get_fem_noeud(i)->get_y();
105			mat(1,2)=mat(1,2)+get_fonction_derive_interpolation(i+1,3,uvw)*get_fem_noeud(i)->get_y();
106			mat(2,0)=mat(2,0)+get_fonction_derive_interpolation(i+1,1,uvw)*get_fem_noeud(i)->get_z();
107			mat(2,1)=mat(2,1)+get_fonction_derive_interpolation(i+1,2,uvw)*get_fem_noeud(i)->get_z();
108			mat(2,2)=mat(2,2)+get_fonction_derive_interpolation(i+1,3,uvw)*get_fem_noeud(i)->get_z();
109			vec(0)=vec(0)-get_fonction_interpolation(i+1,uvw)*get_fem_noeud(i)->get_x();
110			vec(1)=vec(1)-get_fonction_interpolation(i+1,uvw)*get_fem_noeud(i)->get_y();
111			vec(2)=vec(2)-get_fonction_interpolation(i+1,uvw)*get_fem_noeud(i)->get_z();
112			}
113			double det=mat.get_determinant();
114	francois	684	if (fabs(det)<1e-12) ok=3;
115	francois	628	OT_MATRICE_3D mat1(vec,mat.get_vecteur2(),mat.get_vecteur3());
116			OT_MATRICE_3D mat2(mat.get_vecteur1(),vec,mat.get_vecteur3());
117			OT_MATRICE_3D mat3(mat.get_vecteur1(),mat.get_vecteur2(),vec);
118			double d1=mat1.get_determinant()/det;
119			double d2=mat2.get_determinant()/det;
120			double d3=mat3.get_determinant()/det;
121	francois	635	compteur++;
122	francois	684	if (compteur%100==0) ok=2;
123	francois	628	if ((fabs(d1)<1e-8)&&(fabs(d2)<1e-8)&&(fabs(d3)<1e-8)) ok=1;
124			uvw[0]=uvw[0]+d1;
125			uvw[1]=uvw[1]+d2;
126			uvw[2]=uvw[2]+d3;
127			}
128	francois	684	}
129	francois	674	return true;
130	francois	628	}
131
132
133			void FEM_ELEMENT3::get_interpolation_xyz(double* uvw, double* xyz)
134			{
135			xyz[0]=0.;xyz[1]=0.;xyz[2]=0;
136			for (int i=0;i<get_nb_fem_noeud();i++)
137			{
138			xyz[0]=xyz[0]+get_fonction_interpolation(i+1,uvw)*get_fem_noeud(i)->get_x();
139			xyz[1]=xyz[1]+get_fonction_interpolation(i+1,uvw)*get_fem_noeud(i)->get_y();
140			xyz[2]=xyz[2]+get_fonction_interpolation(i+1,uvw)*get_fem_noeud(i)->get_z();
141			}
142			}
143
144	francois	1104	double FEM_ELEMENT3::get_jacobien(double* jac,double *uv,double unite)
145			{
146			int nb=get_nb_fem_noeud();
147
148			OT_MATRICE_3D jacobien;
149			for (int i=0;i<3;i++)
150			{
151			jacobien(i,0)=0.;jacobien(i,1)=0.;jacobien(i,2)=0.;
152			for (int k=0;k<nb;k++)
153			{
154			double valderiv=get_fonction_derive_interpolation(k+1,i+1,uv);
155			double *xyz=get_fem_noeud(k)->get_coord();
156			jacobien(i,0)=jacobien(i,0)+valderivxyz[0]unite;
157			jacobien(i,1)=jacobien(i,1)+valderivxyz[1]unite;
158			jacobien(i,2)=jacobien(i,2)+valderivxyz[2]unite;
159			}
160			}
161			jac[0]=jacobien(0,0);
162			jac[1]=jacobien(0,1);
163			jac[2]=jacobien(0,2);
164
165			jac[3]=jacobien(1,0);
166			jac[4]=jacobien(1,1);
167			jac[5]=jacobien(1,2);
168
169			jac[6]=jacobien(2,0);
170			jac[7]=jacobien(2,1);
171			jac[8]=jacobien(2,2);
172
173			double det=jacobien.get_determinant();
174			return det;
175			}
176
177			void FEM_ELEMENT3::get_inverse_jacob(double* j,double *uv,double unite)
178			{
179			double jac[9];
180			get_jacobien(jac,uv,unite);
181			OT_MATRICE_3D J;
182			J(0,0)=jac[0];
183			J(0,1)=jac[1];
184			J(0,2)=jac[2];
185			J(1,0)=jac[3];
186			J(1,1)=jac[4];
187			J(1,2)=jac[5];
188			J(2,0)=jac[6];
189			J(2,1)=jac[7];
190			J(2,2)=jac[8];
191			OT_MATRICE_3D j_i=J.inverse();
192
193
194			j[0]=j_i(0,0);
195			j[1]=j_i(0,1);
196			j[2]=j_i(0,2);
197
198			j[3]=j_i(1,0);
199			j[4]=j_i(1,1);
200			j[5]=j_i(1,2);
201
202			j[6]=j_i(2,0);
203			j[7]=j_i(2,1);
204			j[8]=j_i(2,2);
205			}