J'ai un problème concernant une implémentation d'un programme sur la méthode des éléments finis en hydrogéologie en C++, au fait, le programme doit travailler sur une matrice 595x595, et il "crashe" à la compilation ("peut être du à la taille de la matrice", je ne sais pas exactement où est le problème), je veux réécrire le code sous matlab.
Le hic c'est que j'ai 5 fichiers ".txt" que le programme doit lire chaque élément et récupérer le nombre de ligne de chaque fichier.
je cherche un code pour reprogrammer sous un script matlab, merci!
Mes requêtes:
-ouvrir les 5 fichiers .txt dans un script.
-lire chaque élément pour chaque fichier.
-récupérer le nombre de lignes et de colonnes pour chaque fichier .txt
les 5 fichiers .txt dans le dossier zippé sont :
-connectivites (ligne i: valeur(i,1) valeur(i,2) valeur(i,3))
-conditions_aux_limites (ligne i: valeur(i,1) valeur(i,1))
-transmissivites (ligne i: valeur(i,1))
-coordonnees (ligne i: valeur(i,1) valeur(i,2)
-terme_source (ligne i: valeur(i,1))
Merci de votre coopération et bonne journée
PS: voici le code projet utilisé en C++
fichier "main.cpp
#include <cstdlib>
#include
#include"outils.h"
#include <fstream>
#include <conio.h>
#include <math.h>
using namespace std;
//fonction pour le calcul du nombre de noeuds et de colonnes à partir des fichiers de données
int nbrlignes(char chemin[])
{
string ligne;
int count=0;
ifstream fi(chemin);
if(fi==0) return(0);
while(!fi.eof())
{
getline(fi, ligne);
count++;
}
fi.close();
return(count);
}
int main(int argc, char *argv[])
{
ifstream f1("coordonnees.txt");
ifstream f2("connectivites.txt");
ifstream f3("terme_source.txt");
ifstream f4("transmissivites.txt");
ifstream f5("conditions_aux_limites.txt");
ofstream g1("resultat1.txt"); //fichier de sortie des résultats de résolution par éléments finis
//Recuperer le nombre des noeuds et des triangles a partir des fichiers de donnees
int NT=nbrlignes("connectivites.txt");
int NN=nbrlignes("coordonnees.txt");
int NCL=nbrlignes("conditions_aux_limites.txt");
cout<<" Informations pratiques sur le maillage"<<endl<<endl;
cout<<" Le nombre de triangle de votre maillage est : "<<NT<<" "<<endl;
cout<<" Le nombre de noeud de votre maillage : "<<NN<<" "<<endl;
cout<<" Le nombre de noeud de conditions aux limites : "<<NCL<<" "<<endl<<endl;
cout<<" ============================================================================="<<endl;
system("PAUSE");
string reponse;
//déclarations des différentes variables
tab k(NN,NN),G(NN,NN);
vect SM(NN),SMG(NN),H(NN),det(NN),TS(NN),X(NN),Y(NN);
float T[NN];
int sommet[NT][3],CL[NCL][2];
//initialisation
for(int i=0;i<NN;i++)
{
for(int j=0;j<NN;j++)
{
G(i,j)=0;
}
}
for(int i=0;i<NN;i++)
{
SMG[i]=0;
SM[i]=0;
X[i]=0;
Y[i]=0;
}
//lecture du fichier des coordonnées
if ( !f1 )
{
cout << "le fichier des coordonnees est inexistant";
}
else
{
for(int i=0;i<NN;i++)
f1>>X[i]>>Y[i];
}
f1.close();
//lecture du fichier des connectivités
if ( !f2 )
{
cout << "le fichier des sommets est inexistant";
}
else
{
for(int i=0;i<NT;i++)
{
for(int j=0;j<3;j++)
f2>>sommet[i][j];
}
}
f2.close();
//lecture du fichier des termes source
if ( !f3 )
{
cout << "le fichier des termes source est inexistant";
}
else
{
for(int i=0;i<NN;i++)
{
f3>>TS[i];
SM[i]=TS[i];
}
}
f3.close();
//lecture du fichier des transmissivités
if ( !f4 )
{
cout << "le fichier des transmissivites est inexistant";
}
else
{
for(int i=0;i<NT;i++)
{
f4>>T[i];
}
}
f4.close();
// lecture du fichier des cnditions aux limites
if ( !f5 )
{
cout << "le fichier des conditions aux limites est inexistant";
}
else
{
for(int i=0;i<NCL;i++)
{
for(int j=0;j<2;j++)
f5>>CL[i][j];
}
}
f5.close();
//Remplissage de la matrice globale elementaire et de la matrice Globale générale
for(int i=0;i<NT;i++)
{
int i1=sommet[i][0]-1, i2=sommet[i][1]-1, i3=sommet[i][2]-1;
det[i]=(2*(((X[i2]-X[i1])*(Y[i3]-Y[i1]))-((X[i3]-X[i1])*(Y[i2]-Y[i1]))));
//élément i1,i1
k(i1,i1)=((Y[i3]-Y[i2])*(Y[i3]-Y[i2]))+((X[i3]-X[i2])*(X[i3]-X[i2]));
k(i1,i1)=(k(i1,i1))*(T[i]/det[i]);
//élément i1,i2
k(i1,i2)=((Y[i3]-Y[i2])*(Y[i1]-Y[i3]))+((X[i3]-X[i2])*(X[i1]-X[i3]));
k(i1,i2)=(k(i1,i2))*(T[i]/det[i]);
//élément i1,i3
k(i1,i3)=((Y[i2]-Y[i1])*(Y[i3]-Y[i2]))+((X[i2]-X[i1])*(X[i3]-X[i2]));
k(i1,i3)=(k(i1,i3))*(T[i]/det[i]);
//élément i2,i2
k(i2,i2)=((Y[i3]-Y[i1])*(Y[i3]-Y[i1]))+((X[i3]-X[i1])*(X[i3]-X[i1]));
k(i2,i2)=k(i2,i2)*(T[i]/det[i]);
//élemnt i2,i3
k(i2,i3)=((Y[i1]-Y[i3])*(Y[i2]-Y[i1]))+((X[i1]-X[i3])*(X[i2]-X[i1]));
k(i2,i3)=k(i2,i3)*(T[i]/det[i]);
//élément i3,i3
k(i3,i3)=((Y[i2]-Y[i1])*(Y[i2]-Y[i1]))+((X[i2]-X[i1])*(X[i2]-X[i1]));
k(i3,i3)=k(i3,i3)*(T[i]/det[i]);
//la transposé
k(i2,i1)=k(i1,i2);
k(i3,i1)=k(i1,i3);
k(i2,i1)=k(i1,i2);
k(i3,i2)=k(i2,i3);
k(i3,i2)=k(i2,i3);
//Assemblage
for(int i=0;i<NN;i++)
{
for(int j=0;j<NN;j++)
{
G(i,j)=G(i,j)+k(i,j);
k(i,j)=0;
}
}
//second membre
for(int i=0;i<NN;i++)
{
SMG[i]= SMG[i]+SM[i];
}
for(int i=0;i<NN;i++)
{
SM[i]=0;
}
}
int x;double y;
//introduction des conditions aux limites
for(int j=0;j<NCL;j++)
{
x=CL[j][0]-1;
y=CL[j][1];
for(int i=0;i<NN;i++)
{
SMG[i]=SMG[i]-((G(i,x))*y);
G(i,x)=0;
G(x,i)=0;
}
SMG[x]=y;
G(x,x)=1;
}
//resolution du système linéaire
gauss_seidel( G ,SMG ,50 ,H,NN);
//Edition du résultat dans un fichier de sortie
for(int i=0;i<NN;i++)
{
g1<<X[i]<<" "<<Y[i]<<" "<<H[i]<<endl;
}
ofstream gg("Matrice_Globale.txt");
for(int i=0;i<NN;i++)
{
for(int j=0;j<NN;j++)
{
gg<<" "<<G(i,j);
}
gg<<endl;
}
ofstream hh("Second_Membre.txt");
for(int i=0;i<NN;i++)
{
hh<<" "<<SMG[i];
hh<<"\n";
}
ofstream zz("Piezo_calcule.txt");
for(int i=0;i<NN;i++)
{
zz<<" "<<H[i];
zz<<"\n";
}
cout<<"Le second membre est"<<endl<<endl;
for(int i=0;i<NN;i++)cout<<"q"< > (nl)
{
nlig=nl;
ncol=nc;
vector<double> v(nc);
for(int i=0;i<nl;i++)
(*this)[i]=v;
}
double & tab::operator()(int i,int j)
{
if((i<0)||(i>nlig))
i=0 ;
if((j<0)||(j>ncol))
j=0 ;
return (*this)[i][j];
}
ostream & operator<<(ostream & os, tab & t)
{
for(int i=0;i<t.nlig;i++)
{
for(int j=0;j<t.ncol;j++)
{
os<<t(i,j)<<" ";
}
os<<"\n";
}
return os ;
}
ostream & operator<<(ostream & os,vect & v)
{
vector<double>::iterator iv;
for (iv=v.begin();iv!=v.end();iv++)
os<<*iv<<" ";
return os ;
}
void gauss_seidel(tab& A,vect& b,int itmax,vect& C,int nv)
{
int k,i,j;
float s;
for(k=0;k
#include <vector>
using namespace std;
class tab:vector<vector<double> >
{
public :
int nlig,ncol;
tab(int nl,int nc);
double & operator()(int,int);
friend ostream & operator<<(ostream &, tab &);
};
class vect:public vector<double>
{
public :
vect(int n):vector<double>(n){}
friend ostream & operator<<(ostream &,vect &);
};
void gauss_seidel( tab& ,vect& ,int ,vect&,int );
#endif