geometrie/src/vct_sommet.cpp

#include"gestionversion.h"
//---------------------------------------------------------------------------


#pragma hdrstop

#include "vct_sommet.h"
#include "mg_sommet.h"
#include "mg_point.h"
#include "vct_point.h"
#include <iomanip>

//---------------------------------------------------------------------------

#pragma package(smart_init)


VCT_SOMMET::VCT_SOMMET(MG_SOMMET* sommet):VCT_ELEMENT_TOPOLOGIQUE(sommet)
{

int indx_premier_ptctr;
MG_POINT* point=sommet->get_point();

TPL_LISTE_ENTITE<double> nurbs_params;
point->get_param_NURBS(indx_premier_ptctr,nurbs_params) ;

OT_VECTEUR_4DD V;
           V[0]=nurbs_params.get(indx_premier_ptctr+  0);
           V[1]=nurbs_params.get(indx_premier_ptctr+  1);
           V[2]=nurbs_params.get(indx_premier_ptctr+  2);
           V[3]=nurbs_params.get(indx_premier_ptctr+  3);
lst_points.insert(lst_points.end(),V);
OT_VECTEUR_4DD VCT_NUL(0.,0.,0.,0.);
lst_vecteurs.insert(lst_vecteurs.end(), VCT_NUL);

nb_points=1;

}


VCT_SOMMET::VCT_SOMMET(VCT_SOMMET& mdd):VCT_ELEMENT_TOPOLOGIQUE(mdd.elem_topo)
{
lst_vecteurs=mdd.lst_vecteurs;
lst_points=mdd.lst_points;
nb_points=mdd.nb_points ;
}


VCT_SOMMET::~ VCT_SOMMET()
{

}


std::vector<OT_VECTEUR_4DD> &VCT_SOMMET::get_vecteurs()
{
return lst_vecteurs;
}


OT_TENSEUR VCT_SOMMET:: calcule_tenseur_metrique()
{
OT_TENSEUR tns(lst_vecteurs);
return tns;
}

OT_VECTEUR_4DD VCT_SOMMET::calcule_barycentre()
{

OT_VECTEUR_4DD barycentre(0,0,0,0);

for(int i=0;i<nb_points;i++)
 {
 barycentre+=lst_vecteurs[i];
 }

barycentre*=1./nb_points;
return barycentre;
}


int VCT_SOMMET::get_nb_points()
 {
 return nb_points  ;
 }


OT_TENSEUR  VCT_SOMMET::calcule_covariance(void)
{

 std::vector<OT_VECTEUR_4DD> Points_Centre_au_barycentre(nb_points);
 OT_VECTEUR_4DD barycentre=calcule_barycentre();

   for( int i=0;i<nb_points;i++)
    {
    Points_Centre_au_barycentre[i]=lst_points[i]-barycentre;
    }

 OT_TENSEUR COVARIANCE(4);

     for(int r=0;r<4;r++) {
       for(int s=0;s<4;s++) {
             for(int i=0;i<nb_points;i++)  {
                 OT_VECTEUR_4DD v1=Points_Centre_au_barycentre[i];
                 for(int j=0;j<nb_points;j++) {
                    OT_VECTEUR_4DD v2=Points_Centre_au_barycentre[j];
                    COVARIANCE(r,s)= COVARIANCE(r,s)+v1[i]*v2[j];
                    }
               }
               COVARIANCE(r,s)=1./nb_points*COVARIANCE(r,s);
         }
      }
 return COVARIANCE;
}


void  VCT_SOMMET::calcule_axes_dinertie(OT_VECTEUR_4DD& vp1,OT_VECTEUR_4DD& vp2,OT_VECTEUR_4DD& vp3,OT_VECTEUR_4DD& vp4)
{
 int n=4;
 int nrot;
 OT_TENSEUR v(4);
 OT_TENSEUR cov=calcule_covariance();
 OT_VECTEUR_4DD d;

 double2 zro=0.0;

 if(cov(0,0)==zro&&cov(1,0)==zro&&cov(2,0)==zro&&cov(3,0)==zro)
   {
   cov(0,0)=1.0;
   cov(1,0)=zro;
   cov(2,0)=zro;
   cov(3,0)=zro;
   }
 if(cov(0,1)==zro&&cov(1,1)==zro&&cov(2,1)==zro&&cov(3,1)==zro)
   {
   cov(0,1)=zro;
   cov(1,1)=1.0;
   cov(2,1)=zro;
   cov(3,1)=zro;
   }
  if(cov(0,2)==zro&&cov(1,2)==zro&&cov(2,2)==zro&&cov(3,2)==zro)
   {
   cov(0,2)=zro;
   cov(1,2)=zro;
   cov(2,2)=1.0;
   cov(3,2)=zro;
   }
  if(cov(0,3)==zro&&cov(1,3)==zro&&cov(2,3)==zro&&cov(3,3)==zro)
   {
   cov(0,3)=zro;
   cov(1,3)=zro;
   cov(2,3)=zro;
   cov(3,3)=1.0;
   }


cov.get_orthogonalisation(cov,d,v,n,nrot);

}


OT_TENSEUR VCT_SOMMET:: calcule_tenseur_inertie_au_barycentre()  //repere globale
{
 OT_TENSEUR Inertie(4);
 OT_TENSEUR cov=calcule_covariance();
 double2 nb_pts=nb_points;
 double2 MOINS_UN=-1.0;
 Inertie=cov*MOINS_UN;
 Inertie=Inertie*nb_pts;

 Inertie(0,0)=(cov(1,1)+cov(2,2)+cov(3,3))*nb_pts;
 Inertie(1,1)=(cov(0,0)+cov(2,2)+cov(3,3))*nb_pts;
 Inertie(2,2)=(cov(1,1)+cov(0,0)+cov(3,3))*nb_pts;
 Inertie(3,3)=(cov(1,1)+cov(2,2)+cov(0,0))*nb_pts;
 return Inertie;
}


OT_TENSEUR VCT_SOMMET:: calcule_tenseur_inertie_au_pt(OT_VECTEUR_4DD& POINT)
{

OT_VECTEUR_4DD BARYCENTRE=this->calcule_barycentre();
OT_TENSEUR TENS=this->calcule_tenseur_inertie_au_barycentre();

OT_VECTEUR_4DD ABC=POINT-BARYCENTRE;

double2 a=ABC[0] ;
double2 b=ABC[1] ;
double2 c=ABC[2] ;
double2 d=ABC[3] ;

TENS(0,0)=TENS(0,0)+b*b+c*c+d*d;
TENS(1,1)=TENS(1,1)+a*a+c*c+d*d;
TENS(2,2)=TENS(2,2)+a*a+b*b+d*d;
TENS(3,3)=TENS(3,3)+a*a+b*b+d*d;

TENS(0,1)=TENS(0,1)-nb_points*a*b;
TENS(0,2)=TENS(0,2)-nb_points*a*c;
TENS(0,3)=TENS(0,3)-nb_points*a*d;

TENS(1,0)=TENS(1,0)-nb_points*b*a;
TENS(1,2)=TENS(1,2)-nb_points*b*c;
TENS(1,3)=TENS(1,3)-nb_points*b*d;

TENS(2,0)=TENS(2,0)-nb_points*c*a;
TENS(2,1)=TENS(2,1)-nb_points*c*b;
TENS(2,3)=TENS(2,3)-nb_points*c*d;

TENS(3,0)=TENS(3,0)-nb_points*d*a;
TENS(3,1)=TENS(3,1)-nb_points*d*b;
TENS(3,2)=TENS(3,2)-nb_points*d*c;

return TENS;
}


OT_TENSEUR VCT_SOMMET::calcule_tenseur_inertie_base_locale()
{
// inertie calcul� dans le repere globale
 OT_VECTEUR_4DD v1,v2,v3,v4;
 OT_TENSEUR INERTIE;
 OT_TENSEUR I_GLOBALE;

 this->calcule_axes_dinertie(v1,v2,v3,v4);
 INERTIE=this->calcule_tenseur_inertie_au_barycentre();

 I_GLOBALE= INERTIE;

 OT_TENSEUR P(4);
 OT_TENSEUR I=I_GLOBALE;
 OT_TENSEUR P_TRSPOSE;
 OT_TENSEUR I_LOCALE;

 for(int i=0;i<4;i++) {
    for(int j=0;j<4;j++)
    {
    if (j==0) P(i,j)= v1[i];
    if (j==1) P(i,j)= v2[i];
    if (j==2) P(i,j)= v3[i];
    if (j==3) P(i,j)= v4[i];
    }
  }

 P_TRSPOSE=P.transpose();

 I_LOCALE=P_TRSPOSE*I;
 I_LOCALE=I_LOCALE*P;
  return I_LOCALE;

}


void VCT_SOMMET::calcule_tenseur_inertie_base_entite( OT_VECTEUR_4DD& g1g2,OT_TENSEUR& tens2, VCT& vct_f)
{
OT_VECTEUR_4DD v1,v2,v3,v4,w1,w2,w3,w4,G1,G2,G1G2;
OT_TENSEUR P(4);
OT_TENSEUR pt;
OT_TENSEUR IV;
OT_TENSEUR Q;
OT_TENSEUR qt;
OT_TENSEUR R;
OT_TENSEUR Rt;
OT_TENSEUR IW;

this->calcule_axes_dinertie(v1,v2,v3,v4);
vct_f.calcule_axes_dinertie(w1,w2,w3,w4);


 for(int i=0;i<4;i++) {
    for(int j=0;j<4;j++)
    {
    if (j==0) P(i,j)= v1[i];
    if (j==1) P(i,j)= v2[i];
    if (j==2) P(i,j)= v3[i];
    if (j==3) P(i,j)= v4[i];
    }
  }

 for(int i=0;i<4;i++) {
    for(int j=0;j<4;j++)
    {
    if (j==0) P(i,j)= w1[i];
    if (j==1) P(i,j)= w2[i];
    if (j==2) P(i,j)= w3[i];
    if (j==3) P(i,j)= w4[i];
    }
  }

pt=P.transpose();
qt=Q.transpose();
R=pt*Q;
Rt=qt*P;
IW=vct_f.calcule_tenseur_inertie_base_locale();

IV=R*IW;
IV=IV*Rt;


G1=this->calcule_barycentre();
G2=vct_f.calcule_barycentre();

G1G2=G2-G1;

OT_VECTEUR_4DD G1G2_LOC;  //Calcul du vecteur G1G2 dans la base locale;
  for(int i=0;i<4;i++)
  {
    for(int j=0;j<4;j++)
    {
    G1G2_LOC[i]= G1G2_LOC[i]+pt(i,j)*G1G2[j];
    }
  }


OT_TENSEUR MASSE_CONCENTRE(4);

MASSE_CONCENTRE(0,0)=nb_points*(G1G2_LOC[1]*G1G2_LOC[1]+G1G2_LOC[2]*G1G2_LOC[2]+G1G2_LOC[3]*G1G2_LOC[3]);
MASSE_CONCENTRE(0,1)=-1*nb_points*(G1G2_LOC[0]*G1G2_LOC[1]);
MASSE_CONCENTRE(0,2)=-1*nb_points*(G1G2_LOC[0]*G1G2_LOC[2]);
MASSE_CONCENTRE(0,3)=-1*nb_points*(G1G2_LOC[0]*G1G2_LOC[3]);

MASSE_CONCENTRE(1,0)=-1*nb_points*(G1G2_LOC[1]*G1G2_LOC[0]);
MASSE_CONCENTRE(1,1)=   nb_points*(G1G2_LOC[0]*G1G2_LOC[0]+G1G2_LOC[2]*G1G2_LOC[2]+G1G2_LOC[3]*G1G2_LOC[3]);
MASSE_CONCENTRE(1,2)=-1*nb_points*(G1G2_LOC[1]*G1G2_LOC[2]);
MASSE_CONCENTRE(1,3)=-1*nb_points*(G1G2_LOC[1]*G1G2_LOC[3]);

MASSE_CONCENTRE(2,0)=-1*nb_points*(G1G2_LOC[2]*G1G2_LOC[0]);
MASSE_CONCENTRE(2,1)=-1*nb_points*(G1G2_LOC[2]*G1G2_LOC[1]);
MASSE_CONCENTRE(2,2)=nb_points*(G1G2_LOC[0]*G1G2_LOC[0]+G1G2_LOC[1]*G1G2_LOC[1]+G1G2_LOC[3]*G1G2_LOC[3]);
MASSE_CONCENTRE(2,3)=-1*nb_points*(G1G2_LOC[0]*G1G2_LOC[2]);

MASSE_CONCENTRE(3,0)=-1*nb_points*(G1G2_LOC[0]*G1G2_LOC[3]);
MASSE_CONCENTRE(3,1)=-1*nb_points*(G1G2_LOC[1]*G1G2_LOC[3]);
MASSE_CONCENTRE(3,2)=-1*nb_points*(G1G2_LOC[2]*G1G2_LOC[3]);
MASSE_CONCENTRE(3,3)=nb_points*(G1G2_LOC[0]*G1G2_LOC[0]+G1G2_LOC[1]*G1G2_LOC[1]+G1G2_LOC[2]*G1G2_LOC[2]);

IV=IV+MASSE_CONCENTRE;
}


void VCT_SOMMET::enregistrer(ostream& o)
{ /*
  int compt=0;
  o<<"%%"<<((MG_SOMMET*)elem_topo)->get_id()<<"="<<"(POINTS_CONTROL_SOMMET,"<< endl;
  o<<"%%"<<((MG_ARETE*)elem_topo)->get_id()<<"="<<"(POLYGONE_DE_CONTROL_SOMMET,"<< endl;
  compt=0;
  for(unsigned int i=0;i<LISTE_TOPO_VECTEUR_SOMMET.size();i++)
    {
      o<<" %  "<<setw(15)<<LISTE_TOPO_VECTEUR_SOMMET[i];//<< " % , ";
       compt++;
       if(compt==4) {
       o<<endl;
       compt=0;}
    } 
     o<<")" <<endl;
     */
}


Revision:	69
Committed:	Thu Mar 27 13:20:26 2008 UTC (17 years, 1 month ago) by souaissa
Original Path:	*magic/lib/geometrie/geometrie/src/vct_sommet.cpp*
File size:	8521 byte(s)
Log Message:	Mise à jour des classes de la vectorisation et des calsses: ot_mathematique,ot_tenseur,ot_doubleprecision dans outil
#	User	Rev	Content
1	souaissa	66	#include"gestionversion.h"
2			//---------------------------------------------------------------------------
3
4
5			#pragma hdrstop
6
7			#include "vct_sommet.h"
8			#include "mg_sommet.h"
9			#include "mg_point.h"
10			#include "vct_point.h"
11			#include <iomanip>
12
13			//---------------------------------------------------------------------------
14
15			#pragma package(smart_init)
16
17
18
19			VCT_SOMMET::VCT_SOMMET(MG_SOMMET* sommet):VCT_ELEMENT_TOPOLOGIQUE(sommet)
20			{
21
22	souaissa	69	int indx_premier_ptctr;
23	souaissa	66	MG_POINT* point=sommet->get_point();
24	souaissa	69
25			TPL_LISTE_ENTITE<double> nurbs_params;
26	souaissa	66	point->get_param_NURBS(indx_premier_ptctr,nurbs_params) ;
27
28	souaissa	69	OT_VECTEUR_4DD V;
29			V[0]=nurbs_params.get(indx_premier_ptctr+ 0);
30			V[1]=nurbs_params.get(indx_premier_ptctr+ 1);
31			V[2]=nurbs_params.get(indx_premier_ptctr+ 2);
32			V[3]=nurbs_params.get(indx_premier_ptctr+ 3);
33			lst_points.insert(lst_points.end(),V);
34			OT_VECTEUR_4DD VCT_NUL(0.,0.,0.,0.);
35			lst_vecteurs.insert(lst_vecteurs.end(), VCT_NUL);
36	souaissa	66
37	souaissa	69	nb_points=1;
38
39	souaissa	66	}
40
41	souaissa	67
42	souaissa	66	VCT_SOMMET::VCT_SOMMET(VCT_SOMMET& mdd):VCT_ELEMENT_TOPOLOGIQUE(mdd.elem_topo)
43			{
44	souaissa	69	lst_vecteurs=mdd.lst_vecteurs;
45			lst_points=mdd.lst_points;
46			nb_points=mdd.nb_points ;
47	souaissa	66	}
48
49
50			VCT_SOMMET::~ VCT_SOMMET()
51			{
52
53			}
54
55
56
57	souaissa	69	std::vector<OT_VECTEUR_4DD> &VCT_SOMMET::get_vecteurs()
58	souaissa	66	{
59	souaissa	69	return lst_vecteurs;
60	souaissa	66	}
61
62
63	souaissa	69	OT_TENSEUR VCT_SOMMET:: calcule_tenseur_metrique()
64	souaissa	66	{
65	souaissa	69	OT_TENSEUR tns(lst_vecteurs);
66			return tns;
67	souaissa	66	}
68
69	souaissa	69	OT_VECTEUR_4DD VCT_SOMMET::calcule_barycentre()
70			{
71	souaissa	66
72	souaissa	69	OT_VECTEUR_4DD barycentre(0,0,0,0);
73	souaissa	66
74	souaissa	69	for(int i=0;i<nb_points;i++)
75			{
76			barycentre+=lst_vecteurs[i];
77			}
78	souaissa	66
79	souaissa	69	barycentre*=1./nb_points;
80			return barycentre;
81			}
82	souaissa	66
83
84
85	souaissa	69	int VCT_SOMMET::get_nb_points()
86			{
87			return nb_points ;
88			}
89	souaissa	66
90
91	souaissa	69	OT_TENSEUR VCT_SOMMET::calcule_covariance(void)
92			{
93
94			std::vector<OT_VECTEUR_4DD> Points_Centre_au_barycentre(nb_points);
95			OT_VECTEUR_4DD barycentre=calcule_barycentre();
96
97			for( int i=0;i<nb_points;i++)
98			{
99			Points_Centre_au_barycentre[i]=lst_points[i]-barycentre;
100			}
101
102			OT_TENSEUR COVARIANCE(4);
103
104			for(int r=0;r<4;r++) {
105			for(int s=0;s<4;s++) {
106			for(int i=0;i<nb_points;i++) {
107			OT_VECTEUR_4DD v1=Points_Centre_au_barycentre[i];
108			for(int j=0;j<nb_points;j++) {
109			OT_VECTEUR_4DD v2=Points_Centre_au_barycentre[j];
110			COVARIANCE(r,s)= COVARIANCE(r,s)+v1[i]*v2[j];
111			}
112			}
113			COVARIANCE(r,s)=1./nb_points*COVARIANCE(r,s);
114			}
115			}
116			return COVARIANCE;
117	souaissa	66	}
118
119
120	souaissa	69	void VCT_SOMMET::calcule_axes_dinertie(OT_VECTEUR_4DD& vp1,OT_VECTEUR_4DD& vp2,OT_VECTEUR_4DD& vp3,OT_VECTEUR_4DD& vp4)
121	souaissa	66	{
122	souaissa	69	int n=4;
123			int nrot;
124			OT_TENSEUR v(4);
125			OT_TENSEUR cov=calcule_covariance();
126			OT_VECTEUR_4DD d;
127
128			double2 zro=0.0;
129
130			if(cov(0,0)==zro&&cov(1,0)==zro&&cov(2,0)==zro&&cov(3,0)==zro)
131			{
132			cov(0,0)=1.0;
133			cov(1,0)=zro;
134			cov(2,0)=zro;
135			cov(3,0)=zro;
136			}
137			if(cov(0,1)==zro&&cov(1,1)==zro&&cov(2,1)==zro&&cov(3,1)==zro)
138			{
139			cov(0,1)=zro;
140			cov(1,1)=1.0;
141			cov(2,1)=zro;
142			cov(3,1)=zro;
143			}
144			if(cov(0,2)==zro&&cov(1,2)==zro&&cov(2,2)==zro&&cov(3,2)==zro)
145			{
146			cov(0,2)=zro;
147			cov(1,2)=zro;
148			cov(2,2)=1.0;
149			cov(3,2)=zro;
150			}
151			if(cov(0,3)==zro&&cov(1,3)==zro&&cov(2,3)==zro&&cov(3,3)==zro)
152			{
153			cov(0,3)=zro;
154			cov(1,3)=zro;
155			cov(2,3)=zro;
156			cov(3,3)=1.0;
157			}
158
159
160			cov.get_orthogonalisation(cov,d,v,n,nrot);
161
162	souaissa	66	}
163
164	souaissa	69
165			OT_TENSEUR VCT_SOMMET:: calcule_tenseur_inertie_au_barycentre() //repere globale
166	souaissa	66	{
167	souaissa	69	OT_TENSEUR Inertie(4);
168			OT_TENSEUR cov=calcule_covariance();
169			double2 nb_pts=nb_points;
170			double2 MOINS_UN=-1.0;
171			Inertie=cov*MOINS_UN;
172			Inertie=Inertie*nb_pts;
173	souaissa	66
174	souaissa	69	Inertie(0,0)=(cov(1,1)+cov(2,2)+cov(3,3))*nb_pts;
175			Inertie(1,1)=(cov(0,0)+cov(2,2)+cov(3,3))*nb_pts;
176			Inertie(2,2)=(cov(1,1)+cov(0,0)+cov(3,3))*nb_pts;
177			Inertie(3,3)=(cov(1,1)+cov(2,2)+cov(0,0))*nb_pts;
178			return Inertie;
179			}
180	souaissa	66
181	souaissa	69
182			OT_TENSEUR VCT_SOMMET:: calcule_tenseur_inertie_au_pt(OT_VECTEUR_4DD& POINT)
183	souaissa	66	{
184
185	souaissa	69	OT_VECTEUR_4DD BARYCENTRE=this->calcule_barycentre();
186			OT_TENSEUR TENS=this->calcule_tenseur_inertie_au_barycentre();
187	souaissa	66
188	souaissa	69	OT_VECTEUR_4DD ABC=POINT-BARYCENTRE;
189	souaissa	66
190	souaissa	69	double2 a=ABC[0] ;
191			double2 b=ABC[1] ;
192			double2 c=ABC[2] ;
193			double2 d=ABC[3] ;
194	souaissa	66
195	souaissa	69	TENS(0,0)=TENS(0,0)+bb+cc+d*d;
196			TENS(1,1)=TENS(1,1)+aa+cc+d*d;
197			TENS(2,2)=TENS(2,2)+aa+bb+d*d;
198			TENS(3,3)=TENS(3,3)+aa+bb+d*d;
199	souaissa	66
200	souaissa	69	TENS(0,1)=TENS(0,1)-nb_pointsab;
201			TENS(0,2)=TENS(0,2)-nb_pointsac;
202			TENS(0,3)=TENS(0,3)-nb_pointsad;
203	souaissa	67
204	souaissa	69	TENS(1,0)=TENS(1,0)-nb_pointsba;
205			TENS(1,2)=TENS(1,2)-nb_pointsbc;
206			TENS(1,3)=TENS(1,3)-nb_pointsbd;
207	souaissa	67
208	souaissa	69	TENS(2,0)=TENS(2,0)-nb_pointsca;
209			TENS(2,1)=TENS(2,1)-nb_pointscb;
210			TENS(2,3)=TENS(2,3)-nb_pointscd;
211	souaissa	66
212	souaissa	69	TENS(3,0)=TENS(3,0)-nb_pointsda;
213			TENS(3,1)=TENS(3,1)-nb_pointsdb;
214			TENS(3,2)=TENS(3,2)-nb_pointsdc;
215	souaissa	66
216	souaissa	69	return TENS;
217			}
218	souaissa	66
219	souaissa	69
220			OT_TENSEUR VCT_SOMMET::calcule_tenseur_inertie_base_locale()
221			{
222			// inertie calcul� dans le repere globale
223			OT_VECTEUR_4DD v1,v2,v3,v4;
224			OT_TENSEUR INERTIE;
225			OT_TENSEUR I_GLOBALE;
226
227			this->calcule_axes_dinertie(v1,v2,v3,v4);
228			INERTIE=this->calcule_tenseur_inertie_au_barycentre();
229
230			I_GLOBALE= INERTIE;
231
232			OT_TENSEUR P(4);
233			OT_TENSEUR I=I_GLOBALE;
234			OT_TENSEUR P_TRSPOSE;
235			OT_TENSEUR I_LOCALE;
236
237			for(int i=0;i<4;i++) {
238			for(int j=0;j<4;j++)
239			{
240			if (j==0) P(i,j)= v1[i];
241			if (j==1) P(i,j)= v2[i];
242			if (j==2) P(i,j)= v3[i];
243			if (j==3) P(i,j)= v4[i];
244			}
245			}
246
247			P_TRSPOSE=P.transpose();
248
249			I_LOCALE=P_TRSPOSE*I;
250			I_LOCALE=I_LOCALE*P;
251			return I_LOCALE;
252
253	souaissa	66	}
254
255	souaissa	69
256
257
258			void VCT_SOMMET::calcule_tenseur_inertie_base_entite( OT_VECTEUR_4DD& g1g2,OT_TENSEUR& tens2, VCT& vct_f)
259	souaissa	66	{
260	souaissa	69	OT_VECTEUR_4DD v1,v2,v3,v4,w1,w2,w3,w4,G1,G2,G1G2;
261			OT_TENSEUR P(4);
262			OT_TENSEUR pt;
263			OT_TENSEUR IV;
264			OT_TENSEUR Q;
265			OT_TENSEUR qt;
266			OT_TENSEUR R;
267			OT_TENSEUR Rt;
268			OT_TENSEUR IW;
269	souaissa	66
270	souaissa	69	this->calcule_axes_dinertie(v1,v2,v3,v4);
271			vct_f.calcule_axes_dinertie(w1,w2,w3,w4);
272
273
274			for(int i=0;i<4;i++) {
275			for(int j=0;j<4;j++)
276			{
277			if (j==0) P(i,j)= v1[i];
278			if (j==1) P(i,j)= v2[i];
279			if (j==2) P(i,j)= v3[i];
280			if (j==3) P(i,j)= v4[i];
281			}
282			}
283
284			for(int i=0;i<4;i++) {
285			for(int j=0;j<4;j++)
286			{
287			if (j==0) P(i,j)= w1[i];
288			if (j==1) P(i,j)= w2[i];
289			if (j==2) P(i,j)= w3[i];
290			if (j==3) P(i,j)= w4[i];
291			}
292			}
293
294			pt=P.transpose();
295			qt=Q.transpose();
296			R=pt*Q;
297			Rt=qt*P;
298			IW=vct_f.calcule_tenseur_inertie_base_locale();
299
300			IV=R*IW;
301			IV=IV*Rt;
302
303
304			G1=this->calcule_barycentre();
305			G2=vct_f.calcule_barycentre();
306
307			G1G2=G2-G1;
308
309			OT_VECTEUR_4DD G1G2_LOC; //Calcul du vecteur G1G2 dans la base locale;
310			for(int i=0;i<4;i++)
311			{
312			for(int j=0;j<4;j++)
313			{
314			G1G2_LOC[i]= G1G2_LOC[i]+pt(i,j)*G1G2[j];
315			}
316			}
317
318
319
320			OT_TENSEUR MASSE_CONCENTRE(4);
321
322			MASSE_CONCENTRE(0,0)=nb_points(G1G2_LOC[1]G1G2_LOC[1]+G1G2_LOC[2]G1G2_LOC[2]+G1G2_LOC[3]G1G2_LOC[3]);
323			MASSE_CONCENTRE(0,1)=-1nb_points(G1G2_LOC[0]*G1G2_LOC[1]);
324			MASSE_CONCENTRE(0,2)=-1nb_points(G1G2_LOC[0]*G1G2_LOC[2]);
325			MASSE_CONCENTRE(0,3)=-1nb_points(G1G2_LOC[0]*G1G2_LOC[3]);
326
327			MASSE_CONCENTRE(1,0)=-1nb_points(G1G2_LOC[1]*G1G2_LOC[0]);
328			MASSE_CONCENTRE(1,1)= nb_points(G1G2_LOC[0]G1G2_LOC[0]+G1G2_LOC[2]G1G2_LOC[2]+G1G2_LOC[3]G1G2_LOC[3]);
329			MASSE_CONCENTRE(1,2)=-1nb_points(G1G2_LOC[1]*G1G2_LOC[2]);
330			MASSE_CONCENTRE(1,3)=-1nb_points(G1G2_LOC[1]*G1G2_LOC[3]);
331
332			MASSE_CONCENTRE(2,0)=-1nb_points(G1G2_LOC[2]*G1G2_LOC[0]);
333			MASSE_CONCENTRE(2,1)=-1nb_points(G1G2_LOC[2]*G1G2_LOC[1]);
334			MASSE_CONCENTRE(2,2)=nb_points(G1G2_LOC[0]G1G2_LOC[0]+G1G2_LOC[1]G1G2_LOC[1]+G1G2_LOC[3]G1G2_LOC[3]);
335			MASSE_CONCENTRE(2,3)=-1nb_points(G1G2_LOC[0]*G1G2_LOC[2]);
336
337			MASSE_CONCENTRE(3,0)=-1nb_points(G1G2_LOC[0]*G1G2_LOC[3]);
338			MASSE_CONCENTRE(3,1)=-1nb_points(G1G2_LOC[1]*G1G2_LOC[3]);
339			MASSE_CONCENTRE(3,2)=-1nb_points(G1G2_LOC[2]*G1G2_LOC[3]);
340			MASSE_CONCENTRE(3,3)=nb_points(G1G2_LOC[0]G1G2_LOC[0]+G1G2_LOC[1]G1G2_LOC[1]+G1G2_LOC[2]G1G2_LOC[2]);
341
342			IV=IV+MASSE_CONCENTRE;
343	souaissa	66	}
344
345	souaissa	69
346	souaissa	66	void VCT_SOMMET::enregistrer(ostream& o)
347	souaissa	69	{ /*
348	souaissa	66	int compt=0;
349			o<<"%%"<<((MG_SOMMET*)elem_topo)->get_id()<<"="<<"(POINTS_CONTROL_SOMMET,"<< endl;
350			o<<"%%"<<((MG_ARETE*)elem_topo)->get_id()<<"="<<"(POLYGONE_DE_CONTROL_SOMMET,"<< endl;
351			compt=0;
352			for(unsigned int i=0;i<LISTE_TOPO_VECTEUR_SOMMET.size();i++)
353			{
354			o<<" % "<<setw(15)<<LISTE_TOPO_VECTEUR_SOMMET[i];//<< " % , ";
355			compt++;
356			if(compt==4) {
357			o<<endl;
358			compt=0;}
359			}
360			o<<")" <<endl;
361	souaissa	69	*/
362	souaissa	66	}
363
364
365
366