programmation en C++ de léquation de la chaleur

programmation en C++ de léquation de la chaleur - C++ - Programmation

Marsh Posté le 24-02-2007 à 19:00:29    

voila je ne suis qu'un petit étudiant en physique qui doit programmer en C++, je connai un pe le C mais la je rame vraiment. En gros je dois modéliser la diffusion de la chaleur sur une plaque a 2 dimension. cette plaque est laissé à température constant que l'on va nommer T0. L'aspect physique du truc je le maitrise et je sais comment faire.
alors il faut que je crée une matrice de taille J*J dont l'algorythme est:
     for(i=0;i<J;i=i++)    { M[i][i]=(1+4*mu);}
     for(i=0;i<J-1;i=i++)  { M[i][i+1]=-mu;}  
     for(i=0;i<J-1;i=i++)  { M[i+1][i]=-mu;}
     for(i=0;i<J-I;i=i++)  { M[i][i+I]=-mu;}
     for(i=0;i<J-I;i=i++)  { M[I+i][i]=-mu;}
et les autres élément sont tous nulle et mu est un double on a aussi I*I=J (genre I=2, J=4)
 
bon kan la matrice est construite il fau que je l'inverse, je pense faire du gauss jordan ou un carrémen prendre la définition de l'inversion de matrice.
puis kan cette matrice est inversée je doi la multiplier avec un vecteur et bien sur récupérer la solution qui est un vecteur J ligne.
ce vecteur J ligne je le remultipli avec la matrice qui a ete inversé et je recupere un vecteur que je remultipli et ainsi de suite... je trouve donc la chaleur à chaque instant. de facon plus calculatoire ca donne: M ma matrice invM son inverse, T1 & T2 la température tel que T2=(T1+T0)M avec T0 ki dépen de la température des bords. et aussi information utile, je connai la premiere température T1 ki représente en fait la température de la plaque a l'instan 0. De plus j'ai quand mem essayé de commencé et j'ai fai un petit programme (ki ne marche pa) pour construire M et multiplier M avec un vecteur ce programme est le suivan:
#include <iostream>  
#include <fstream>  
#include <math.h>  
 
using namespace std;  
const int I=2;
const int J=4;
 
class vecteur
{
double vect[J];
 public  :
 vecteur(double T0) {int i,j;
for(i=0;i<J;i++) vect[i]=T0;}
friend vecteur prod(Matrice,vecteur);
};  
               
class Matrice  
{  
 
double Mat[J][J];
public  :
Matrice (double mu) {int i,j;
for(i=0;i<J;i=i++)  {for(j=0;j<J;j=j++)   { Mat[i][j]=0;}}  
for(i=0;i<J;i=i++)  {for(i=0;i<J;i=i++)   { Mat[i][i]=(1+4*mu);}}  
for(i=0;i<J-1;i=i++){for(i=0;i<J-1;i=i++) { Mat[i][i+1]=-mu;}}  
for(i=0;i<J-1;i=i++){for(i=0;i<J-1;i=i++) { Mat[i+1][i]=-mu;}}  
for(i=0;i<J-I;i=i++){for(i=0;i<J-I;i=i++) { Mat[i][i+I]=-mu;}}  
for(i=0;i<J-I;i=i++){for(i=0;i<J-I;i=i++) { Mat[I+i][i]=-mu;}}  
                    }
friend vecteur prod(Matrice,vecteur);
};                                    
 
 
vecteur prod(matrice m,vecteur x)
{int i,j;
double som;
vecteur res;  
for(i=0;i<J;i=i++)  { for(j=0;som=0;j<J;j=j++) som+=m.Mat[i][j]*x.vect[j];
                      res.vect[i]=som;
                      return(res);
                    }
 
   
 
 main ()  
{  
int i,j;
double mu=1;  
double T0=0;  
vecteur T;
Matrice M;
res=prod(M,T);
for(i=0;i<J;i++) {cout<< res[i] <<"  ";} cout<< endl;}
 
}  
bon ba jaimerai ke kelkun maide vou pouvez me toper sur msn g pa bocou de tps pour finir ca et c capital.  

Reply

Marsh Posté le 24-02-2007 à 19:00:29   

Reply

Marsh Posté le 24-02-2007 à 19:13:27    

ce programme est peu etre incohérent mais jai fai comme jle sentai si ya kelkun ki pe mdir pourkoi ca marche pa ce serai vraimen cool.

Reply

Marsh Posté le 24-02-2007 à 19:17:16    

et les valeur de T0 et mu sont arbitraire on peut mettre otre choz

Reply

Marsh Posté le 24-02-2007 à 23:23:41    

euh... c'est quoi qui marche pas ? ca compile ?


Message édité par nrth le 24-02-2007 à 23:23:57

---------------
BarrePoint.org : News pour techno-geeks.
Reply

Marsh Posté le 25-02-2007 à 07:57:09    

Putain, je vois que l'apprentissage de la prog aux étudiants en Physique, c'est tjrs aussi catastrophique.  :pfff: Enfin, maintenant, ils enseignent le C++, c'est déjà ça. A mon époque, c'était encore du Fortran 77. [:petrus75]

 

Un code de classe matrice, -qui marche, elle -, t'est gracieusement offert ici. Il te reste juste à implémenter la multiplication (inspire-toi du code de l'opérateur = pour cette dernière), la transposée puis Gauss-Jordan, ce qui devrait être très facile avec cette base. Merci de poster le code ici quand tu as fini, hisoitre d'en faire bénéficier d'autres personnes. Après tout, on t'a aidé.
Lis tout le topic, il donne d'autres directions intéressantes orientées calcul scientifique (Blitz++, uBlas, Boost ou lib de NT² de Joel F qui a l'air très simpel à utiliser (voir sa signature)). Mais comme ton truc est un exercice, je ne suis pa sûr que ton prof accepte cette solution.


Message édité par el muchacho le 25-02-2007 à 08:35:36
Reply

Marsh Posté le 25-02-2007 à 18:29:27    

non ca compile pa ca me di kil conai pa le type de retour de la fonction friend alors ke je lui donne vect.... compren po
 

Reply

Marsh Posté le 25-02-2007 à 19:17:40    

pierre581 a écrit :

non ca compile pa ca me di kil conai pa le type de retour de la fonction friend alors ke je lui donne vect.... compren po


merci d'éviter le style SMS stp, par respect pour ceux qui te lisent et qui te répondent.


---------------
J'ai un string dans l'array (Paris Hilton)
Reply

Marsh Posté le 25-02-2007 à 19:19:18    

res est pas défini dans ton main déjà...(et j n'a pas de valeur)


---------------
BarrePoint.org : News pour techno-geeks.
Reply

Marsh Posté le 25-02-2007 à 21:17:52    

Salut,
il y a clairement une accolade fermante en trop au dernier for()...

Reply

Marsh Posté le 25-02-2007 à 23:22:48    

désolé pour le style sms

Reply

Marsh Posté le 25-02-2007 à 23:22:48   

Reply

Marsh Posté le 26-02-2007 à 08:30:52    

for(i=0;i<J-1;i=i++){for(i=0;i<J-1;i=i++) { Mat[i][i+1]=-mu;}}
 
tu fais comment la différence entre les i ?

Reply

Marsh Posté le 26-02-2007 à 11:33:37    

En fait au début je ne faisais qu'une boucle vu que c'est juste pour remplir des termes appartenant à une meme diagonale mais ça ne marchais pas (ce problème m'a fait d'ailleurs perdre énormément de temps) puis qu'en j'ai mis  une autre boucle ma matrice se remplissait correctement, donc je pense que ce n'est pas un probleme.
L'erreur se situe au niveau des fonctions amies il ne veut pas reconnaitre l'objet vecteur.

Reply

Marsh Posté le 26-02-2007 à 11:39:55    

non ...
 
PUTAIN JE LES HAIE CES PROFS A LA CON QUI ONT RIEN D4AUTRES A FOUTRE QUE D4ENSEIGNER FRIEND ALORS QU4EN DES ANNEES DE C++ J4AI JAMAIS AU A UTILISER.

Reply

Marsh Posté le 26-02-2007 à 11:41:19    

j'utilise quoi alors pour utiliser deux classes avec une meme fonction

Reply

Marsh Posté le 26-02-2007 à 11:45:57    

tu fais marcher ta boîte à cerveau en terme de visibilité.
Et fais un joli

Code :
  1. vecteur& vecteur::operator*=(const Matrice &m)
  2. {
  3.   // des trucs
  4. }
  5. vecteur operator*(const Matrice &m, const vecteur &v)
  6. {
  7.   vecteur res(v);
  8.   res *= m;
  9.   return res;
  10. }


 
toutes façons, ton code compile même pas et en est loin.

Reply

Marsh Posté le 26-02-2007 à 11:54:25    

Just my two cents...
 
Tu n'as pas écris exactement ton problème, mais si tu considères l'équation de la chaleur utilisée dans 99.99% des cas, ta matrice M est fausse, c'est des 1-2*mu que tu as sur la diagonale. De surcroît Taz a raison, tes boucles en i ne sont pas nécessaires, tu devrais pouvoir tout faire en une seule boucle du type:

Code :
  1. for(i=0; i<J-1; ++i)
  2. {
  3.    for (j=0; j<J-1; ++j)
  4.    {
  5.       if(j==i) Mat[i][j]=1-2*mu;
  6.       else if(j==i+1 || j==i-1) Mat[i][j]=mu;
  7.       else Mat[i][j] = 0.;
  8.    }
  9. }


 
Enfin faire un schéma explicite pour ce type de problème, avec absolument aucun contrôle sur le pas de temps, c'est l'échec assuré le schéma n'étant pas inconditionellement stable. Tout ca pour dire, que même si ton code compile, ou marche, la résolution mathématique est assez loin d'être correcte (m'enfin je dis ca je dis rien).  Bon courage!

Reply

Marsh Posté le 26-02-2007 à 12:03:36    

ok mais je te promet que sur la diagonale c'est 1+4mu cela vient du fait que j'utilise une méthode implicite à 2 dimensions sinon je suis d'accord au niveau des boucles en i. En fait j'avais fait ca parceque ça marchais et ça m'affichais la bonne matrice. Mais je prend note merci

Reply

Marsh Posté le 26-02-2007 à 12:09:42    

Tu veux sans doute dire méthode explicite à 2 dimensions? Au temps pour moi je n'ai que survoler ton code (quoiqu'il en soit en général le pb se formule avec des 1-qqchose plutôt que l'inverse). De totu les manières tu vas avoir des pb de stabilité et des oscillations, si tu utilises un CN ou un schémas implicite, tu auras qqchose d'inconditionellement stable. Cheers!    

Reply

Marsh Posté le 26-02-2007 à 12:10:35    

c'est d'ailleurs parce que c'est la méthode implicite que je dois inverser la matrice  

Reply

Marsh Posté le 26-02-2007 à 12:13:51    

l autre méthode je l'ai fait en C et tu a raison elle est instable à partir d'une certaine valeur de mu la température diverge.
 

Reply

Marsh Posté le 26-02-2007 à 20:52:49    

Utilise mon code et les exemples de Taz, je te dis, et commence par écrire une vraie lib matricielle en complétant ce que je t'ai donné avant de te lancer dans la résolution de l'équation.
Bref, apprends à programmer.
Là, ce que tu fais, c'est n'importe quoi, c'est du code dégueu, tellement dégueu qu'on n'a pas envie de le lire et que tu fais des erreurs dans tous les sens.


Message édité par el muchacho le 26-02-2007 à 20:54:49
Reply

Marsh Posté le 26-02-2007 à 22:56:29    

Taz a écrit :

PUTAIN JE LES HAIE CES PROFS A LA CON QUI ONT RIEN D4AUTRES A FOUTRE QUE D4ENSEIGNER FRIEND ALORS QU4EN DES ANNEES DE C++ J4AI JAMAIS AU A UTILISER.


Alors on a pas du tout la même expérience ...

Reply

Marsh Posté le 26-02-2007 à 23:03:34    

el muchacho c'est justement ce que j'essaie de faire mais c'est pas intuitif,dans le prog que tu m'a filé, avant de l'utilisé, il faut que je saisisse toute les chses qu'il y a dedans.Mais je vais m'en sortir avec un peu de percéverence, laisse moi juste un peu de temps et soit un peu plus cool et moins méchant.

Reply

Marsh Posté le 26-02-2007 à 23:12:14    

mais merci quand meme

Reply

Marsh Posté le 27-02-2007 à 16:00:32    

Si tu te limites a des schemas d'ordre deux alors tu auras toujours une matrice tridiagonale, en general (et c'est le cas pour ton pb) a diagonale dominante. Tout ca pour dire qu'un structure efficace consiste a creer une classe TriMatrix (avec trois vecteurs pour les trois diagonales). Ensuite il te suffit d'écrire les différentes opérations/opérateurs pour cette classe et surtout une méthode d'inversion basée sur cette structure en particulier (en gros descente/remontée) qui t'inverse ta matrice en 2*N et tu auras un code a peu près efficace.


Message édité par ElDesdichado le 27-02-2007 à 16:01:17
Reply

Marsh Posté le 27-02-2007 à 16:39:40    

je suis entierement d'accord avec toi mais ce que tu me di c'est pour une résolution impicite unidimensionnelle, la mon but c'est de faire cette méthode en deux dimensions. J'ai donc trouvé une manière simple qui est de développer ma matrice représentant ma température sur ma plaque selon les colonnes. C'est pour cela que ma plaque qui est de dimension I*I me donne un vecteur de J=I*I lignes et donc ma matrice (de dimension J*J) qui me permet de passer d'une température T(t+dt) à T(t) n'est plus tridiagonale. elle a bien sur des élément tridiagonaux mais il faut lui rajouter deux autres diagonales. C'est pour cela que j'avais fait un algo tout moche pour vous faire une idée de la tete de ma matrice. Je pense que mon approche physique est bonne.La j'essai d'apprendre à programmer. Je posterai mon programme à la fin quand il sera fini(de facon il faut le finir mon passage en M2 en dépend !!).
Et tu as raison pour résoudre l'inversion d'une matrice tridiagonale il suffit de prendre un truc tel que l'algorithme de thomas ou meme la méthode de résolution de cholewski (je crois!!!). Mais la je pense qu'il me faut un pivot de gauss.

Reply

Marsh Posté le 27-02-2007 à 17:35:01    

ce que je te dis en deux dim peut s'appliquer pour n'importequelle matrice de Toeplitz i.e. a bande avec pleins de 0 ou tu ne stockes que les diagonales, un peu d'algèbre et ca fonctionne nettement plus rapidement. Sinon Choleski va être plus lent, il vaut mieux une méthode Gauss directe avec un peu de réflexion tu arrives a la solution la plus rapide.

Reply

Marsh Posté le 28-02-2007 à 21:08:18    

pierre581 a écrit :

el muchacho c'est justement ce que j'essaie de faire mais c'est pas intuitif,dans le prog que tu m'a filé, avant de l'utilisé, il faut que je saisisse toute les chses qu'il y a dedans.Mais je vais m'en sortir avec un peu de percéverence, laisse moi juste un peu de temps et soit un peu plus cool et moins méchant.


OK désolé d'avoir été aggressif, on est un peu trop réactifs, par ici.  :jap: Tu devrais pouvoir y arriver, ce code est très simple, sans la moindre ruse particulière, et pas franchement optimisé, mais c'est une base de départ. Par contre, il faudra sans doute la modifier pas mal si tu veux implémenter les méthodes de ElDesdichado (que je ne connais pas).


Message édité par el muchacho le 28-02-2007 à 21:27:57
Reply

Marsh Posté le 28-02-2007 à 23:46:58    

Bon j'ai commencé à programmer, et c'est pas une mince affaire que de rester des heures sur le PC pour en fin de compte devoir poster une fois de plus car je suis bloqué. Ce qui est positif c'est que si ce programme est débogé j'ai déja tous les algo pour finir en temps et en heure.
Bon mon programme ne réinvente rien il essai seulement de multiplier une ligne d'une matrice par un réel. mon but est bien sur par la suite de finir la méthode des pivots de gauss pour inverser la matrice (puis résoudre cette équation de la chaleur).
Bon je pense que ca commence à être mieu mais bon ca marche pas et je sais qu'il y a un problème au niveau de la fonction dilatation mais des heures durant jai cherché et je trouve pas.
 
 
Voici le matrix.hh

Code :
  1. using namespace std;
  2. #include <iostream>
  3. #include <math.h>
  4. const int I=4;
  5. template <class T>
  6. class matrix
  7. {
  8. T mat[I*I][I*I];
  9. public:
  10. matrix(T t[I*I][I*I]); 
  11. matrix();
  12. void affiche_mat();
  13. matrix<T> dilatation(int,double);
  14. };


voici le matrix.cpp

Code :
  1. using namespace std;
  2. #include <iostream>
  3. #include <math.h>
  4. #include "Matrix.hh"
  5. template <class T>
  6. matrix<T>::matrix<T>(T t[N][N])
  7. { int i; int j;
  8. for(i=0;i<I*I;i++)
  9.    {
  10.      for(j=0;j<I*I;j++) mat[i][j]=t[i][j];
  11.    }
  12. }
  13. template <class T>
  14. matrix<T>::matrix<T>()
  15. { int i; int j;
  16. for(i=0;i<I*I;i++)
  17.    {
  18.      for(j=0;j<I*I;j++) mat[i][j]=0;
  19.    }
  20. }
  21. template <class T>
  22. void matrix<T>::affiche_mat()
  23. {
  24. int i,j;
  25. for(i=0;i<I*I;i++)
  26.    {
  27.      for(j=0;j<I*I;j++) cout<<mat[i][j]<<" ";
  28.      cout<<endl;
  29.    }
  30.  
  31.  
  32. }
  33. //dilatation d'une ligne
  34. template <class T>
  35. matrix<T> matrix<T>::dilatation(int i,double k)
  36. {
  37.        matrix<T> dil_i_k();
  38.        int m; int n;
  39.        for(n=0;n<I*I;n++)
  40.        {for (m=0;m<I*I;m++)
  41.             {
  42.             if (m!=i) {dil_i_k.mat[m][n]=mat[m][n];}
  43.             else      {dil_i_k.mat[m][n]=k*mat[m][n];}
  44.             }       
  45.        }
  46. return (dil_i_k);
  47.      
  48. }     [cpp]
  49. et pour finir le main.cpp
  50. [cpp]using namespace std;
  51. #include <iostream>
  52. #include <math.h>
  53. #include "Matrix.hh"
  54. main()
  55. {
  56. double mu=2;
  57. matrix<T> M;
  58. int i,j;
  59. for( i=0;i<I*I;i++){for(j=0;j<I*I;j++) { M.mat[i][j]=0 ;}}
  60. for(i=0;i<I*I;i++) { M.mat[i][i]=(1+4*mu);}
  61. for(i=0;i<I*I-1;i++) { M.mat[i][i+1]=-mu;}
  62. for(i=0;i<I*I-1;i++) { M.mat[i+1][i]=-mu;}
  63. for(i=0;i<I*I-I;i++) { M.mat[i][i+I]=-mu;}
  64. for(i=0;i<I*I-I;i++) { M.mat[I+i][i]=-mu;}
  65. matrix<T> Mat(M); 
  66. Mat.dilatation(3,4.);
  67. Mat.affiche();
  68. }


 
 

Reply

Marsh Posté le 28-02-2007 à 23:50:00    

mon probleme est peut etre dans le makefile je n'ai pas réussi à savoir si c'étai la meme chose pour des classe que pour des template de classe pouvez vous me dire!

Reply

Marsh Posté le 01-03-2007 à 00:41:24    

dans Matrix<T>::dilatation,
matrix<T> dil_i_k(); est une déclaration de fonction, à linkage extern, dans la namespace scope directement englobante.
Ce que tu cherches à faire s'écrit probablement comme cela :
matrix<T> dil_i_k;

 

ligne 76: qu'est-ce que T ? remplaceT par float, double, long double, que sais-je ...

 

Après, tu vas avoir des problèmes de link. Soit tu mets tout ton code dans le header et tu utilises l'instanciation implicite, soit tu instancies explicitement :
à la fin du .cpp :
template class Matrix<float>;
template class Matrix<double>;
template class Matrix<long double>;

 

Ou alors, ton compilateur supporte "export" (peu probable), et il faut préfixer ce qui va bien par export.


Message édité par ++fab le 01-03-2007 à 00:42:52
Reply

Marsh Posté le 02-03-2007 à 07:25:43    

Sérieusement, je ne comprends pas pourquoi il se fait chier avec des templates, à l'allure où il avance. Je ne comprends pas pourquoi il met du template<class T> partout.
Mais s'il avait étudié mon code, il aurait mis Matrix<T>& et non Matrix<T> en retour de sa fonction dilatation.
Là, je pense que c'est un problème de compréhension des pointeurs et des templates, alors il faut faire au plus simple: oublier les templates, faire tout en double, faire attention aux types de retour et on n'en parle plus.


Message édité par el muchacho le 02-03-2007 à 08:52:51
Reply

Marsh Posté le    

Reply

Sujets relatifs:

Leave a Replay

Make sure you enter the(*)required information where indicate.HTML code is not allowed