geometrie/src/mg_volume.cpp

//------------------------------------------------------------
//------------------------------------------------------------
//  MAGiC
//  Jean Christophe Cuilli?re et Vincent FRANCOIS
//  D?partement de G?nie M?canique - UQTR
//------------------------------------------------------------
//  Le projet MAGIC est un projet de recherche du d?partement
//  de g?nie m?canique de l'Universit? du Qu?bec ?
//  Trois Rivi?res
//  Les librairies ne peuvent ?tre utilis?es sans l'accord
//  des auteurs (contact : francois@uqtr.ca)
//------------------------------------------------------------
//------------------------------------------------------------
//
//       mg_volume.cpp
//
//------------------------------------------------------------
//------------------------------------------------------------
//  COPYRIGHT 2000
//  Version du 02/03/2006 ? 11H22
//------------------------------------------------------------
//------------------------------------------------------------


#include "gestionversion.h"
#include "mg_volume.h"
#include "ot_mathematique.h"
#include "mg_maillage.h"
#include "vct_volume.h"
//#include "message.h"

MG_VOLUME::MG_VOLUME(std::string idori,unsigned long num):MG_ELEMENT_TOPOLOGIQUE(num,idori),num_materiau(-1),vect(NULL),mince(false)
{
}

MG_VOLUME::MG_VOLUME(std::string idori):MG_ELEMENT_TOPOLOGIQUE(idori),num_materiau(-1),vect(NULL),mince(false)
{
}

MG_VOLUME::MG_VOLUME(MG_VOLUME& mdd):MG_ELEMENT_TOPOLOGIQUE(mdd),lst_coquille(mdd.lst_coquille),vect(mdd.vect),mince(mdd.mince)
{
}

MG_VOLUME::~MG_VOLUME()
{
}

void MG_VOLUME::ajouter_mg_coquille(class MG_COQUILLE* mgcoq)
{
    lst_coquille.insert(lst_coquille.end(),mgcoq);
}

void MG_VOLUME::supprimer_mg_coquille(class MG_COQUILLE* mgcoq)
{
    std::vector<MG_COQUILLE*>::iterator i;
    for (i=lst_coquille.begin();i!=lst_coquille.end();i++)
    {
        if ((*i)==mgcoq)
        {
            lst_coquille.erase(i);
            return;
        }
    }
}


int MG_VOLUME::get_nb_mg_coquille(void)
{
    return lst_coquille.size();
}

MG_COQUILLE* MG_VOLUME::get_mg_coquille(int num)
{
    return lst_coquille[num];
}

int MG_VOLUME::get_dimension(void)
{
    return 3;
}
VCT& MG_VOLUME::get_vectorisation(void)
{
//VCT* p=NULL;
//return *p;
    if (vect==NULL) vect=new VCT_VOLUME(this);
    return *vect;
}
void MG_VOLUME::get_topologie_sousjacente(TPL_MAP_ENTITE<MG_ELEMENT_TOPOLOGIQUE*> *lst)
{
    int nbcoq=lst_coquille.size();
    for (int i=0;i<nbcoq;i++)
    {
        MG_COQUILLE* coq=lst_coquille[i];
        int nbface=coq->get_nb_mg_coface();
        for (int j=0;j<nbface;j++)
        {
            MG_FACE* face=coq->get_mg_coface(j)->get_face();
            lst->ajouter(face);
            face->get_topologie_sousjacente(lst);
        }
    }
}


BOITE_3D MG_VOLUME::get_boite_3D(int pas_echantillon)
{
BOITE_3D boitevol;
int nb_coquille=get_nb_mg_coquille();
for (int i=0;i<nb_coquille;i++)
    {
    MG_COQUILLE* coq=get_mg_coquille(i);
    int nb_face=coq->get_nb_mg_coface();
    for (int j=0;j<nb_face;j++)
          {
          MG_FACE* face=coq->get_mg_coface(j)->get_face();
          BOITE_3D tmp=face->get_boite_3D(pas_echantillon);
          if ((i==0) && (j==0))  boitevol=tmp;
          else boitevol=tmp+boitevol;
          }
    }
return boitevol;
}

void MG_VOLUME::enregistrer(std::ostream& o,double version)
{
    o << "%" << get_id() << "=VOLUME(" << get_idoriginal() << ",(";
    for (unsigned int i=0;i<lst_coquille.size();i++)
    {
        o << "$" << lst_coquille[i]->get_id();
        if (i!=lst_coquille.size()-1) o << ",";
        else o << ")";
    }
    int nb=get_nb_ccf();
    if (version<2) 
    {
    o << "," << num_materiau << "," << nb;
    if (nb!=0)
    {
        o << ",(";
        for (int i=0;i<nb;i++)
        {
            char nom[3];
            get_type_ccf(i,nom);
            o << "(" <<  nom << "," << get_valeur_ccf(i) << ")";
            if (i!=nb-1) o << "," ;
        }
        o << ")";
    }
    o << ");" << std::endl;
    }
    else
    {
     if (version>2.2)
        {
          if (mince) o << ",1" ;
          else o << ",0" ;
          if (mince)
            {
            o << ",(" ;  
            for (int i=0;i<get_nb_face_correspondante()-1;i++)
            o <<  "$" << lst_base[i]->get_id() << ","  << "$" <<  lst_extrude[i]->get_id() << "," ; 
            o << "$" << lst_base[get_nb_face_correspondante()-1]->get_id() << ","  << "$" << lst_extrude[get_nb_face_correspondante()-1]->get_id(); 
            o << ")" ;  
            }
        }
    o << "," ;
    enregistrer_ccf(o,version);  
    o << ");" << std::endl;
    }
}

void MG_VOLUME::change_num_materiau(int num)
{
    num_materiau=num;
}

int MG_VOLUME::get_num_materiau(void)
{
    return num_materiau;
}


bool MG_VOLUME::est_un_volume_element(void)
{
return false;
}


void MG_VOLUME::ajouter_face_correspondante(MG_FACE* face1,MG_FACE* face2)
{
lst_base.push_back(face1);  
lst_extrude.push_back(face2);  
mince=true;
}

void MG_VOLUME::get_face_correspondante(int num,MG_FACE** face1,MG_FACE** face2)
{
(*face1)=lst_base[num];  
(*face2)=lst_extrude[num];  
}

int MG_VOLUME::get_nb_face_correspondante(void)
{
return lst_base.size();  
}

bool MG_VOLUME::est_mince(void)
{
return mince;  
}


void MG_VOLUME::get_propriete_massique(class MG_MAILLAGE* mai,double& volume,class OT_VECTEUR_3D& cdm,class OT_MATRICE_3D& inertieglobale,class OT_MATRICE_3D& inertiecdm,double dens)
{
    double vol=0.;
    double xg=0.;
    double yg=0.;
    double zg=0.;
    double a=0.;
    double b=0.;
    double c=0.;
    double d=0.;
    double e=0.;
    double f=0.;
    int nbcoquille=get_nb_mg_coquille();
    for (int i=0;i<nbcoquille;i++)
    {
        MG_COQUILLE* coq=get_mg_coquille(i);
        int nbface=coq->get_nb_mg_coface();
        for (int j=0;j<nbface;j++)
        {
            MG_COFACE* coface=coq->get_mg_coface(j);
            MG_FACE* face=coface->get_face();
            TPL_SET<MG_ELEMENT_MAILLAGE*>::ITERATEUR it;
            for (MG_TRIANGLE* tri=(MG_TRIANGLE*)face->get_lien_maillage()->get_premier(it);tri!=NULL;tri=(MG_TRIANGLE*)face->get_lien_maillage()->get_suivant(it))
            {
                if (mai->get_mg_triangleid(tri->get_id())==NULL) continue;
                MG_NOEUD* noeud1=tri->get_noeud1();
                MG_NOEUD* noeud2=tri->get_noeud2();
                MG_NOEUD* noeud3=tri->get_noeud3();
                double *xyz1=noeud1->get_coord();
                double *xyz2=noeud2->get_coord();
                double *xyz3=noeud3->get_coord();
                OT_VECTEUR_3D vec1(xyz1,xyz2);
                OT_VECTEUR_3D vec2(xyz1,xyz3);
                OT_VECTEUR_3D n=vec1&vec2;
                double detj=n.get_longueur();
                n.norme();
                n=-1.*coface->get_orientation()*n;
                for (int k=0;k<3;k++)
                {
                    double xsi=1./6.;
                    double eta=1./6.;
                    double wi=1./6.;
                    if (k==1) xsi=2./3.;
                    if (k==2) eta=2./3.;
                    double x=(1-xsi-eta)*xyz1[0]+xsi*xyz2[0]+eta*xyz3[0];
                    double y=(1-xsi-eta)*xyz1[1]+xsi*xyz2[1]+eta*xyz3[1];
                    double z=(1-xsi-eta)*xyz1[2]+xsi*xyz2[2]+eta*xyz3[2];
                    OT_VECTEUR_3D psi1(x,0,0);
                    OT_VECTEUR_3D psi2(0.5*x*x,0,0);
                    OT_VECTEUR_3D psi3(0,0.5*y*y,0);
                    OT_VECTEUR_3D psi4(0,0,0.5*z*z);
                    OT_VECTEUR_3D psi5(0,y*y*y/3.,z*z*z/3.);
                    OT_VECTEUR_3D psi6(x*x*x/3.,0,z*z*z/3.);
                    OT_VECTEUR_3D psi7(x*x*x/3.,y*y*y/3.,0);
                    OT_VECTEUR_3D psi8(0.,y*y*z/2.,0);
                    OT_VECTEUR_3D psi9(x*x*z/2.,0.,0);
                    OT_VECTEUR_3D psi10(x*x*y/2.,0,0);
                    vol=vol+detj*wi*(psi1*n);
                    xg=xg+wi*detj*(psi2*n);
                    yg=yg+wi*detj*(psi3*n);
                    zg=zg+wi*detj*(psi4*n);
                    a=a+wi*detj*(psi5*n);
                    b=b+wi*detj*(psi6*n);
                    c=c+wi*detj*(psi7*n);
                    d=d+wi*detj*(psi8*n);
                    e=e+wi*detj*(psi9*n);
                    f=f+wi*detj*(psi10*n);
                }
            }
        }
    }
    volume=vol;
    cdm.change_x(xg/volume);
    cdm.change_y(yg/volume);
    cdm.change_z(zg/volume);
    inertieglobale(0,0)=a*dens;
    inertieglobale(1,1)=b*dens;
    inertieglobale(2,2)=c*dens;
    inertieglobale(1,0)=f*dens;
    inertieglobale(0,1)=f*dens;
    inertieglobale(2,0)=e*dens;
    inertieglobale(0,2)=e*dens;
    inertieglobale(2,1)=d*dens;
    inertieglobale(1,2)=d*dens;
    OT_VECTEUR_3D m1(volume*(cdm.get_y()*cdm.get_y()+cdm.get_z()*cdm.get_z()),volume*cdm.get_x()*cdm.get_y(),volume*cdm.get_x()*cdm.get_z());
    OT_VECTEUR_3D m2(volume*cdm.get_x()*cdm.get_y(),volume*(cdm.get_x()*cdm.get_x()+cdm.get_z()*cdm.get_z()),volume*cdm.get_y()*cdm.get_z());
    OT_VECTEUR_3D m3(volume*cdm.get_x()*cdm.get_z(),volume*cdm.get_y()*cdm.get_z(),volume*(cdm.get_x()*cdm.get_x()+cdm.get_y()*cdm.get_y()));
    m1=(-1.)*m1;
    m2=(-1.)*m2;
    m3=(-1.)*m3;
    OT_MATRICE_3D change(m1,m2,m3);
    inertiecdm=inertieglobale+change;
}

Revision:	881
Committed:	Wed Apr 12 23:14:07 2017 UTC (8 years, 1 month ago) by francois
File size:	9129 byte(s)
Log Message:	Changement de strategie pour le mailleur couche. Il fonctionne maintenant avec des volumes minces et non des coquilles mince est une propriète de volume et ne change pas son nom. cela necessite un changement de version de fichier -> maintenant 2.3 compatible en arrière valider sur ancien exemple plus sphere creuse
#	User	Rev	Content
1	francois	283	//------------------------------------------------------------
2			//------------------------------------------------------------
3			// MAGiC
4			// Jean Christophe Cuilli?re et Vincent FRANCOIS
5			// D?partement de G?nie M?canique - UQTR
6			//------------------------------------------------------------
7			// Le projet MAGIC est un projet de recherche du d?partement
8			// de g?nie m?canique de l'Universit? du Qu?bec ?
9			// Trois Rivi?res
10			// Les librairies ne peuvent ?tre utilis?es sans l'accord
11			// des auteurs (contact : francois@uqtr.ca)
12			//------------------------------------------------------------
13			//------------------------------------------------------------
14			//
15			// mg_volume.cpp
16			//
17			//------------------------------------------------------------
18			//------------------------------------------------------------
19			// COPYRIGHT 2000
20			// Version du 02/03/2006 ? 11H22
21			//------------------------------------------------------------
22			//------------------------------------------------------------
23
24
25			#include "gestionversion.h"
26			#include "mg_volume.h"
27			#include "ot_mathematique.h"
28			#include "mg_maillage.h"
29			#include "vct_volume.h"
30			//#include "message.h"
31
32	francois	881	MG_VOLUME::MG_VOLUME(std::string idori,unsigned long num):MG_ELEMENT_TOPOLOGIQUE(num,idori),num_materiau(-1),vect(NULL),mince(false)
33	francois	283	{
34			}
35
36	francois	881	MG_VOLUME::MG_VOLUME(std::string idori):MG_ELEMENT_TOPOLOGIQUE(idori),num_materiau(-1),vect(NULL),mince(false)
37	francois	283	{
38			}
39
40	francois	881	MG_VOLUME::MG_VOLUME(MG_VOLUME& mdd):MG_ELEMENT_TOPOLOGIQUE(mdd),lst_coquille(mdd.lst_coquille),vect(mdd.vect),mince(mdd.mince)
41	francois	283	{
42			}
43
44			MG_VOLUME::~MG_VOLUME()
45			{
46			}
47
48			void MG_VOLUME::ajouter_mg_coquille(class MG_COQUILLE* mgcoq)
49			{
50			lst_coquille.insert(lst_coquille.end(),mgcoq);
51			}
52
53			void MG_VOLUME::supprimer_mg_coquille(class MG_COQUILLE* mgcoq)
54			{
55			std::vector<MG_COQUILLE*>::iterator i;
56			for (i=lst_coquille.begin();i!=lst_coquille.end();i++)
57			{
58			if ((*i)==mgcoq)
59			{
60			lst_coquille.erase(i);
61			return;
62			}
63			}
64			}
65
66
67			int MG_VOLUME::get_nb_mg_coquille(void)
68			{
69			return lst_coquille.size();
70			}
71
72			MG_COQUILLE* MG_VOLUME::get_mg_coquille(int num)
73			{
74			return lst_coquille[num];
75			}
76
77			int MG_VOLUME::get_dimension(void)
78			{
79			return 3;
80			}
81			VCT& MG_VOLUME::get_vectorisation(void)
82			{
83			//VCT* p=NULL;
84			//return *p;
85			if (vect==NULL) vect=new VCT_VOLUME(this);
86			return *vect;
87			}
88			void MG_VOLUME::get_topologie_sousjacente(TPL_MAP_ENTITE<MG_ELEMENT_TOPOLOGIQUE> lst)
89			{
90			int nbcoq=lst_coquille.size();
91			for (int i=0;i<nbcoq;i++)
92			{
93			MG_COQUILLE* coq=lst_coquille[i];
94			int nbface=coq->get_nb_mg_coface();
95			for (int j=0;j<nbface;j++)
96			{
97			MG_FACE* face=coq->get_mg_coface(j)->get_face();
98			lst->ajouter(face);
99			face->get_topologie_sousjacente(lst);
100			}
101			}
102			}
103
104
105	francois	632
106
107			BOITE_3D MG_VOLUME::get_boite_3D(int pas_echantillon)
108			{
109			BOITE_3D boitevol;
110			int nb_coquille=get_nb_mg_coquille();
111			for (int i=0;i<nb_coquille;i++)
112			{
113			MG_COQUILLE* coq=get_mg_coquille(i);
114			int nb_face=coq->get_nb_mg_coface();
115			for (int j=0;j<nb_face;j++)
116			{
117			MG_FACE* face=coq->get_mg_coface(j)->get_face();
118			BOITE_3D tmp=face->get_boite_3D(pas_echantillon);
119			if ((i==0) && (j==0)) boitevol=tmp;
120			else boitevol=tmp+boitevol;
121			}
122			}
123			return boitevol;
124			}
125
126	francois	763	void MG_VOLUME::enregistrer(std::ostream& o,double version)
127	francois	283	{
128			o << "%" << get_id() << "=VOLUME(" << get_idoriginal() << ",(";
129			for (unsigned int i=0;i<lst_coquille.size();i++)
130			{
131			o << "$" << lst_coquille[i]->get_id();
132			if (i!=lst_coquille.size()-1) o << ",";
133			else o << ")";
134			}
135			int nb=get_nb_ccf();
136	francois	763	if (version<2)
137			{
138	francois	283	o << "," << num_materiau << "," << nb;
139			if (nb!=0)
140			{
141			o << ",(";
142			for (int i=0;i<nb;i++)
143			{
144			char nom[3];
145			get_type_ccf(i,nom);
146			o << "(" << nom << "," << get_valeur_ccf(i) << ")";
147			if (i!=nb-1) o << "," ;
148			}
149			o << ")";
150			}
151			o << ");" << std::endl;
152	francois	763	}
153			else
154			{
155	francois	881	if (version>2.2)
156			{
157			if (mince) o << ",1" ;
158			else o << ",0" ;
159			if (mince)
160			{
161			o << ",(" ;
162			for (int i=0;i<get_nb_face_correspondante()-1;i++)
163			o << "$" << lst_base[i]->get_id() << "," << "$" << lst_extrude[i]->get_id() << "," ;
164			o << "$" << lst_base[get_nb_face_correspondante()-1]->get_id() << "," << "$" << lst_extrude[get_nb_face_correspondante()-1]->get_id();
165			o << ")" ;
166			}
167			}
168	francois	763	o << "," ;
169			enregistrer_ccf(o,version);
170			o << ");" << std::endl;
171			}
172	francois	283	}
173
174			void MG_VOLUME::change_num_materiau(int num)
175			{
176			num_materiau=num;
177			}
178
179			int MG_VOLUME::get_num_materiau(void)
180			{
181			return num_materiau;
182			}
183
184	francois	576
185			bool MG_VOLUME::est_un_volume_element(void)
186			{
187			return false;
188			}
189	francois	881
190
191
192
193			void MG_VOLUME::ajouter_face_correspondante(MG_FACE* face1,MG_FACE* face2)
194			{
195			lst_base.push_back(face1);
196			lst_extrude.push_back(face2);
197			mince=true;
198			}
199
200			void MG_VOLUME::get_face_correspondante(int num,MG_FACE face1,MG_FACE face2)
201			{
202			(*face1)=lst_base[num];
203			(*face2)=lst_extrude[num];
204			}
205
206			int MG_VOLUME::get_nb_face_correspondante(void)
207			{
208			return lst_base.size();
209			}
210
211			bool MG_VOLUME::est_mince(void)
212			{
213			return mince;
214			}
215
216
217
218	francois	283	void MG_VOLUME::get_propriete_massique(class MG_MAILLAGE* mai,double& volume,class OT_VECTEUR_3D& cdm,class OT_MATRICE_3D& inertieglobale,class OT_MATRICE_3D& inertiecdm,double dens)
219			{
220			double vol=0.;
221			double xg=0.;
222			double yg=0.;
223			double zg=0.;
224			double a=0.;
225			double b=0.;
226			double c=0.;
227			double d=0.;
228			double e=0.;
229			double f=0.;
230			int nbcoquille=get_nb_mg_coquille();
231			for (int i=0;i<nbcoquille;i++)
232			{
233			MG_COQUILLE* coq=get_mg_coquille(i);
234			int nbface=coq->get_nb_mg_coface();
235			for (int j=0;j<nbface;j++)
236			{
237			MG_COFACE* coface=coq->get_mg_coface(j);
238			MG_FACE* face=coface->get_face();
239			TPL_SET<MG_ELEMENT_MAILLAGE*>::ITERATEUR it;
240			for (MG_TRIANGLE* tri=(MG_TRIANGLE)face->get_lien_maillage()->get_premier(it);tri!=NULL;tri=(MG_TRIANGLE)face->get_lien_maillage()->get_suivant(it))
241			{
242			if (mai->get_mg_triangleid(tri->get_id())==NULL) continue;
243			MG_NOEUD* noeud1=tri->get_noeud1();
244			MG_NOEUD* noeud2=tri->get_noeud2();
245			MG_NOEUD* noeud3=tri->get_noeud3();
246			double *xyz1=noeud1->get_coord();
247			double *xyz2=noeud2->get_coord();
248			double *xyz3=noeud3->get_coord();
249			OT_VECTEUR_3D vec1(xyz1,xyz2);
250			OT_VECTEUR_3D vec2(xyz1,xyz3);
251			OT_VECTEUR_3D n=vec1&vec2;
252			double detj=n.get_longueur();
253			n.norme();
254			n=-1.coface->get_orientation()n;
255			for (int k=0;k<3;k++)
256			{
257			double xsi=1./6.;
258			double eta=1./6.;
259			double wi=1./6.;
260			if (k==1) xsi=2./3.;
261			if (k==2) eta=2./3.;
262			double x=(1-xsi-eta)xyz1[0]+xsixyz2[0]+eta*xyz3[0];
263			double y=(1-xsi-eta)xyz1[1]+xsixyz2[1]+eta*xyz3[1];
264			double z=(1-xsi-eta)xyz1[2]+xsixyz2[2]+eta*xyz3[2];
265			OT_VECTEUR_3D psi1(x,0,0);
266			OT_VECTEUR_3D psi2(0.5xx,0,0);
267			OT_VECTEUR_3D psi3(0,0.5yy,0);
268			OT_VECTEUR_3D psi4(0,0,0.5zz);
269			OT_VECTEUR_3D psi5(0,yyy/3.,zzz/3.);
270			OT_VECTEUR_3D psi6(xxx/3.,0,zzz/3.);
271			OT_VECTEUR_3D psi7(xxx/3.,yyy/3.,0);
272			OT_VECTEUR_3D psi8(0.,yyz/2.,0);
273			OT_VECTEUR_3D psi9(xxz/2.,0.,0);
274			OT_VECTEUR_3D psi10(xxy/2.,0,0);
275			vol=vol+detjwi(psi1*n);
276			xg=xg+widetj(psi2*n);
277			yg=yg+widetj(psi3*n);
278			zg=zg+widetj(psi4*n);
279			a=a+widetj(psi5*n);
280			b=b+widetj(psi6*n);
281			c=c+widetj(psi7*n);
282			d=d+widetj(psi8*n);
283			e=e+widetj(psi9*n);
284			f=f+widetj(psi10*n);
285			}
286			}
287			}
288			}
289			volume=vol;
290			cdm.change_x(xg/volume);
291			cdm.change_y(yg/volume);
292			cdm.change_z(zg/volume);
293			inertieglobale(0,0)=a*dens;
294			inertieglobale(1,1)=b*dens;
295			inertieglobale(2,2)=c*dens;
296			inertieglobale(1,0)=f*dens;
297			inertieglobale(0,1)=f*dens;
298			inertieglobale(2,0)=e*dens;
299			inertieglobale(0,2)=e*dens;
300			inertieglobale(2,1)=d*dens;
301			inertieglobale(1,2)=d*dens;
302			OT_VECTEUR_3D m1(volume(cdm.get_y()cdm.get_y()+cdm.get_z()cdm.get_z()),volumecdm.get_x()cdm.get_y(),volumecdm.get_x()*cdm.get_z());
303			OT_VECTEUR_3D m2(volumecdm.get_x()cdm.get_y(),volume(cdm.get_x()cdm.get_x()+cdm.get_z()cdm.get_z()),volumecdm.get_y()*cdm.get_z());
304			OT_VECTEUR_3D m3(volumecdm.get_x()cdm.get_z(),volumecdm.get_y()cdm.get_z(),volume(cdm.get_x()cdm.get_x()+cdm.get_y()*cdm.get_y()));
305			m1=(-1.)*m1;
306			m2=(-1.)*m2;
307			m3=(-1.)*m3;
308			OT_MATRICE_3D change(m1,m2,m3);
309			inertiecdm=inertieglobale+change;
310			}
311