mailleur2d_structure.cpp

Aller à la documentation de ce fichier.

//####//------------------------------------------------------------
//####//------------------------------------------------------------
//####//  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/
//####//------------------------------------------------------------
//####//------------------------------------------------------------
//####//
//####//       mailleur2d_structure.cpp
//####//
//####//------------------------------------------------------------
//####//------------------------------------------------------------
//####//    COPYRIGHT 2000-2024
//####//  jeu 13 jun 2024 11:58:55 EDT
//####//------------------------------------------------------------
//####//------------------------------------------------------------
#include <stdio.h>
#include <math.h>
#include "mg_gestionnaire.h"
#include "mailleur2d_structure.h"
#include "mg_maillage.h"
#include "fem_maillage.h"
#include "fem_hexa20.h"
#include "mg_geometrie.h"
#include "mg_sommet_noeud.h"
#include "mg_arete_element.h"
#include "mg_face_element.h"
#include "mg_volume_element.h"
#include "tpl_octree.h"
#include "mg_maillage_outils.h"
#include <vector>
 
 
struct integrale
{
    double ti;
    double li;
};
 
MAILLEUR2D_STRUCTURE::MAILLEUR2D_STRUCTURE(double xmintmp, double xmaxtmp, double ymintmp, double ymaxtmp, int ntmp, int mtmp, MG_GESTIONNAIRE* gesttmp,int avecgeo,double unit):MAILLEUR(false),xmin(xmintmp),xmax(xmaxtmp),ymin(ymintmp),ymax(ymaxtmp),n(ntmp),m(mtmp),gest(gesttmp),geovirtuel(avecgeo),densitevariable(false),unite(unit),nbfonction(0)
{  
}    
MAILLEUR2D_STRUCTURE::MAILLEUR2D_STRUCTURE( double xmintmp, double xmaxtmp, double ymintmp, double ymaxtmp, double (*fx)(double),double (*fy)(double),MG_GESTIONNAIRE* gesttmp,int avecgeo,double unit):MAILLEUR(false),xmin(xmintmp),xmax(xmaxtmp),ymin(ymintmp),ymax(ymaxtmp),gest(gesttmp),geovirtuel(avecgeo),enx(fx),eny(fy),densitevariable(true),unite(unit),nbfonction(0)
{
}
 
MAILLEUR2D_STRUCTURE::MAILLEUR2D_STRUCTURE( MAILLEUR2D_STRUCTURE& mdd):MAILLEUR(mdd),
xmin(mdd.xmin),xmax(mdd.xmax),ymin(mdd.ymin),ymax(mdd.ymax),n(mdd.n),m(mdd.m),gest(mdd.gest),nbfonction(mdd.nbfonction)
{  
}  
 
MAILLEUR2D_STRUCTURE::~MAILLEUR2D_STRUCTURE()
{
}
 
 
void MAILLEUR2D_STRUCTURE::discretise(std::vector<double> &tab,double min,double max,int &nb,double (*en)(double))
{
    double t1=0;
    double t2=1;
    std::vector<integrale> tab_integrale;
    double ti,tii,t;
    double longueur_metrique=0.;
    tii=t1;
    double dt=1./pas;
    do
    {
        ti=tii;
        tii=ti+(t2-t1)*dt;
        if (tii>t2) tii=t2;
        t=0.7886751345*ti+0.2113248654*tii;
        double facteur1=(max-min)/(*en)(t);
        t=0.7886751345*tii+0.2113248654*ti;
        double facteur2=(max-min)/(*en)(t);
        longueur_metrique=longueur_metrique+0.5*(tii-ti)*(facteur1+facteur2);
        integrale pas;
        pas.ti=ti;
        pas.li=longueur_metrique;
        tab_integrale.insert(tab_integrale.end(),pas);
    }
    while (tii<t2);
    nb=(int)floor(longueur_metrique);
    if (longueur_metrique-floor(longueur_metrique)>0.5) nb++;
    if (nb<1)  nb=1;
    double valeur_cible=longueur_metrique/nb;
    tab.push_back(min);
    int nbcree=0;
    int nb_pas_integrale=tab_integrale.size();
    for (int i=0;i<nb_pas_integrale;i++)
    {
        while ((tab_integrale[i].li>(nbcree+1)*valeur_cible) && (nbcree<nb-1))
        {
            double ti=tab_integrale[i].ti;
            double tii=t2;
            double li=0.;
            if (i!=0) li=tab_integrale[i-1].li;
            double lii=tab_integrale[i].li;
            if (i!=nb_pas_integrale-1) tii=tab_integrale[i+1].ti;
            t=ti+(tii-ti)/(lii-li)*((nbcree+1)*valeur_cible-li);
            double x=min+t*(max-min);
            tab.push_back(x);
            nbcree++;
        }
    }
    tab.push_back(max);
    
}
void MAILLEUR2D_STRUCTURE::ajouter_fonction_geometrie(int nb)
{
nbfonction=nb;  
}
 
int MAILLEUR2D_STRUCTURE::maille(MG_GROUPE_TOPOLOGIQUE* mggt)
{
 affiche((char*)"Création d'un maillage MG structuré");
 double intx, inty;
 std::vector<double> tabx,taby;
 if (densitevariable==false)
    {
    intx=fabs(xmax-xmin)/n;
    inty=fabs(ymax-ymin)/m;
    }
 else
    {
    discretise(tabx,xmin,xmax,n,enx);  
    discretise(taby,ymin,ymax,m,eny);  
    }
 
 MG_GEOMETRIE* geo=NULL;
 MG_SOMMET_NOEUD *som1=NULL;
 MG_SOMMET_NOEUD *som2=NULL;
 MG_SOMMET_NOEUD *som3=NULL;
 MG_SOMMET_NOEUD *som4=NULL;
 MG_ARETE_ELEMENT* are1=NULL;
 MG_ARETE_ELEMENT* are2=NULL;
 MG_ARETE_ELEMENT* are3=NULL;
 MG_ARETE_ELEMENT* are4=NULL;
 MG_FACE_ELEMENT* face1=NULL;
 if (geovirtuel==true) 
          {     
          geo=new MG_GEOMETRIE((char*)"VIRTUEL",(char*)"VIRTUEL");
          geo->change_valeur_unite(unite);
          gest->ajouter_mg_geometrie(geo);
          som1=new MG_SOMMET_NOEUD(NULL);
          geo->ajouter_mg_sommet(som1);
          som2=new MG_SOMMET_NOEUD(NULL);
          geo->ajouter_mg_sommet(som2);
          som3=new MG_SOMMET_NOEUD(NULL);
          geo->ajouter_mg_sommet(som3);
          som4=new MG_SOMMET_NOEUD(NULL);
          geo->ajouter_mg_sommet(som4);
          are1=new MG_ARETE_ELEMENT();
          geo->ajouter_mg_arete(are1);
          are2=new MG_ARETE_ELEMENT();
          geo->ajouter_mg_arete(are2);
          are3=new MG_ARETE_ELEMENT();
          geo->ajouter_mg_arete(are3);
          are4=new MG_ARETE_ELEMENT();
          geo->ajouter_mg_arete(are4);
          face1=new MG_FACE_ELEMENT();
          geo->ajouter_mg_face(face1);
          for (int i=0;i<nbfonction;i++)
            {
            MG_GEOM_FONCTION *fonc=new MG_GEOM_FONCTION(0);
            geo->ajouter_mg_geom_fonction(fonc);
            }
          }
 MG_MAILLAGE* mai=new MG_MAILLAGE(geo);
 gest->ajouter_mg_maillage(mai); 
 std::vector<MG_NOEUD*> vecnod;
     for (int j=0; j<=m;j++)
      {
       for (int i=0; i<=n;i++)
         { 
           MG_ELEMENT_TOPOLOGIQUE *topo=NULL;
           if (geovirtuel)
            {
            if ((i==0) && (j==0)) topo=som1;
            else if ((i==0) && (j==m)) topo=som2;
            else if ((i==0) && (j==m)) topo=som3;
            else if ((i==0) && (j==0)) topo=som4;
            else if (j==0) topo=are1;
            else if (i==n) topo=are2;
            else if (j==m) topo=are3;
            else if (j==0) topo=are4;
            else topo=face1;
            }
    double x,y;
          if (densitevariable==false)
            {
            y=ymin+(j*inty);
            x=xmin+(i*intx);
            }
          else
            {
            x=tabx[i];  
            y=taby[j];  
            }
          MG_NOEUD* nou=mai->ajouter_mg_noeud(topo,x,y,0.,MAGIC::ORIGINE::MAILLEUR_AUTO);
          vecnod.push_back(nou);
          if (geovirtuel)
            {
            if ((i==0) && (j==0)) som1->change_mg_noeud(nou);
            else if ((i==n) && (j==0) ) som2->change_mg_noeud(nou);
            else if ((i==n) && (j==m) ) som3->change_mg_noeud(nou);
            else if ((i==0) && (j==m) ) som4->change_mg_noeud(nou);
            }
         }
    }
if (geovirtuel)
  {
  for(int i=0; i<n;i++)
    {
    int j=0;
    MG_NOEUD* nod1=vecnod[i+((n+1)*j)];
    MG_NOEUD* nod2=vecnod[(i+1)+((n+1)*j)]; 
    mai->ajouter_mg_segment(are1,nod1,nod2,MAGIC::ORIGINE::MAILLEUR_AUTO);
    j=m;
    MG_NOEUD* nod3=vecnod[i+((n+1)*j)];
    MG_NOEUD* nod4=vecnod[(i+1)+((n+1)*j)]; 
    mai->ajouter_mg_segment(are3,nod3,nod4,MAGIC::ORIGINE::MAILLEUR_AUTO);
    }
  for(int j=0; j<m;j++)
    {
    int i=0;
    MG_NOEUD* nod1=vecnod[i+((n+1)*j)];
    MG_NOEUD* nod2=vecnod[i+((n+1)*(j+1))]; 
    mai->ajouter_mg_segment(are4,nod1,nod2,MAGIC::ORIGINE::MAILLEUR_AUTO);
    i=n;
    MG_NOEUD* nod3=vecnod[i+((n+1)*j)];
    MG_NOEUD* nod4=vecnod[i+((n+1)*(j+1))]; 
    mai->ajouter_mg_segment(are2,nod3,nod4,MAGIC::ORIGINE::MAILLEUR_AUTO);
    }
  }
 for (int i=0;i<n;i++)
    for (int j=0;j<m;j++)
      {
      MG_NOEUD* nod1=vecnod[i+((n+1)*j)];
      MG_NOEUD* nod2=vecnod[(i+1)+((n+1)*j)]; 
      MG_NOEUD* nod3=vecnod[(i+1)+((n+1)*(j+1))]; 
      MG_NOEUD* nod4=vecnod[i+((n+1)*(j+1))];
      mai->ajouter_mg_quadrangle(face1,nod1,nod4,nod3,nod2,MAGIC::ORIGINE::MAILLEUR_AUTO);
      }
  
 BOITE_3D b(xmin,ymin,0.,xmax,ymax,0.);
 mai->change_param_structure(b,n,m,0);
 /*affiche((char*)"Création d'un maillage FEM structuré");
 FEM_MAILLAGE* fem = new FEM_MAILLAGE(geo, mai,1);
 gest->ajouter_fem_maillage(fem);
 fem->construire();*/
 return OK;
}
 
 
void MAILLEUR2D_STRUCTURE::ajuste(class MG_MAILLAGE* mai, MG_MAILLAGE* maibase, int nb)
{
    LISTE_MG_NOEUD::iterator it;
    double xmin=1e300,ymin=1e300,zmin=1e300;
    double xmax=-1e300,ymax=-1e300,zmax=-1e300;
    TPL_LISTE_ENTITE<MG_NOEUD*> lst;
    for (MG_NOEUD* no=maibase->get_premier_noeud(it);no!=NULL;no=maibase->get_suivant_noeud(it))
        {
        if (no->get_x()<xmin) xmin=no->get_x();    
        if (no->get_y()<ymin) ymin=no->get_y();    
        if (no->get_z()<zmin) zmin=no->get_z();    
        if (no->get_x()>xmax) xmax=no->get_x();    
        if (no->get_y()>ymax) ymax=no->get_y();    
        if (no->get_z()>zmax) zmax=no->get_z();    
        lst.ajouter(no);
        }
    BOITE_3D b(xmin,ymin,zmin,xmax,ymax,zmax);
    b.change_grosseur(1.05);
    TPL_OCTREE<MG_TRIANGLE*,MG_NOEUD*> oc;
    oc.initialiser(&lst,1,b.get_xmin(),b.get_xmin(),b.get_zmin(),b.get_xmax(),b.get_xmax(),b.get_zmax());
    LISTE_MG_TRIANGLE::iterator itt;
    for (MG_TRIANGLE* tri=maibase->get_premier_triangle(itt);tri!=NULL;tri=maibase->get_suivant_triangle(itt))
        oc.inserer(tri);
    
    for (MG_NOEUD* no=mai->get_premier_noeud(it);no!=NULL;no=mai->get_suivant_noeud(it))
            {
            double x=no->get_x();    
            double y=no->get_y();    
            double z=no->get_z();    
            TPL_MAP_ENTITE<MG_TRIANGLE*> lsttrouve;
            oc.rechercher(x,y,z,1e-7,lsttrouve);
            TPL_MAP_ENTITE<MG_TRIANGLE*>::ITERATEUR itf;
            bool trouve=false;
            for (MG_TRIANGLE* tri=lsttrouve.get_premier(itf);tri!=NULL;tri=lsttrouve.get_suivant(itf))
                {
                MG_MAILLAGE_OUTILS ot;
                int res=ot.estdansletriangle(tri,x,y,z);
                int val=ot.compare_etat_triangle(res,MG_MAILLAGE_OUTILS::INTERIEUR);
                if (val>0) trouve=true;
                }
            if (trouve) no->change_nouveau_numero(1); else no->change_nouveau_numero(0);
            }
    LISTE_MG_QUADRANGLE::iterator itq;
    std::vector<unsigned long> listeasupprimer;
    for (MG_QUADRANGLE* quad=mai->get_premier_quadrangle(itq);quad!=NULL;quad=mai->get_suivant_quadrangle(itq))
            {
            MG_NOEUD* no1=quad->get_noeud1();    
            MG_NOEUD* no2=quad->get_noeud2();    
            MG_NOEUD* no3=quad->get_noeud3();    
            MG_NOEUD* no4=quad->get_noeud4();    
            int val=no1->get_nouveau_numero()+no2->get_nouveau_numero()+no3->get_nouveau_numero()+no4->get_nouveau_numero();
            if (val<nb) 
                listeasupprimer.push_back(quad->get_id());
            }
    for (int i=0;i<listeasupprimer.size();i++)
        mai->supprimer_mg_quadrangleid(listeasupprimer[i]);
    
    
    
}