Discussion :

[Asohan]

Bonjour a tous et a toutes,

Excusez moi pour les accent, je n'ai qu'un clavier qwerty, et excusez moi aussi pour mon faible niveau en language C++, c'est pour cette raison que je solicite votre aide et votre immense experience.

Je travaille dans la physique subatomique, et j'ai concu un programme de simulation pour effectuer des interactions a N corps (genre : gestion des effets PhysX de nvidia). Ce programme utilise des calculs assez complexes, mais pour plus de facilite, je stocke, a chaque pas de temps, des tableaux dynamiques de donnees necessaires aux calculs.

Sauf que : mon nombre de particules augmente au fur et a mesure par fission (c'est ce que je veux) : quand cela se produit j'obtiens une erreur :

***glibc detected*** malloc(): memory corruption (fast): 0x.........

Cela veux dire que malgre la quasi-perfection de la STL, la fonction resize que j'utilise pour mes vector< vector<TYPE> > n'utilise pas correctement l'allocation memoire.

Voila une "ebauche" du code que j'utilise :

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
 
Class SYSTEM
{
    public :
         SYSTEM() {};
        ~SYSTEM() {};
 
        ... // quelques reglages de parametres.
 
    private :
        int n_particules;
 
        vector< vector<double> > A_ij;
        vector< vector<double> > B_ij;
 
        ...
 
        void calcul();
        void etape_1();
        void etape_2();
}
 
void SYSTEM::calcul()
{
    n_particules = ...; // calcul.
}
 
void SYSTEM::etape_1()
{
    A_ij.resize(n_particules, vector<double>(n_particules));
}
 
void SYSTEM::etape_2()
{
    B_ij.resize(n_particules, vector<double>(n_particules));
}
 
int main()
{
    SYSTEM system;
 
    system.calcul();
    system.etape_1();
    system.etape_2();
 
    return EXIT_SUCCESS;
}

En gros, si mon programme etait cette ebauche : des que n_particules changerais il effectuerait correctement le resize de A_ij, et il indiquerait le message d'erreur au niveau du resize de B_ij. C'est etrange, non ?

Pourquoi le resize d'un vector contenu dans un objet comme tableau dynamique ferait de la corruption de memoire ?

Si vous avez de l'experience vis-a-vis des memory corruption avec des vector.resize, je suis tres interesse par votre aide. Merci d'avoir lu ce message assez long.

[Joel F]

que vaut ton n_particules ?
Poste un code minimale compilable reproduisant l'erreur.

[Aleph69]

Bonjour,

Cela veux dire que malgre la quasi-perfection de la STL, la fonction resize que j'utilise pour mes vector< vector<TYPE> > n'utilise pas correctement l'allocation memoire.

... ou que tu as fait une erreur dans ton code... les paris sont ouverts!

Il est surtout important que tu envoies la partie du code qui manipule Aij et Bij.
C'est a priori de là que vient le problème... ah oui j'ai parié!

[Asohan]

Eh ben oui !! Vous avez toujours raison c'est pas juste

Bon j'ai perdu combien pour le pari ?

Alors, mon probleme commence a s'eclaircir :
- Il ne s'agit pas de la STL (ca paraissait evident, mais bon, on ne sait jamais),
- il ne s'agit pas de la valeur de n-particule ou d'une autre variable (le calcul est bon et le debugger me donne bien les bonnes valeurs pour mon gros programme)

J'ai teste un code minimal et aussi mon programme apres une grosse cure de Slim Fast, mais du coup je n'ai plus mon erreur

Neanmoins j'ai toujours une erreur glibc memory corruption pour mon code complet. De par votre experience qu'est-ce que ca peut etre autrement ?

Lorsqu'il s'agit d'une segmentation fault ou bien d'une double free or corruption la je vois de quoi il s'agit et je corrige mon erreur, mais pour mon probleme, je ne vois pas du tout. Une idee ?

PS : Je n'ai pas l'autorisation par mon labo de poster le code global librement, mais en revanche, je peux tres bien l'envoyer a quelqu'un par email, car je ne dispose pas d'outils de debuggage avances chez moi.

[Aleph69]

Bonjour,

bon allez on annule les paris : je compatis pour le temps que tu vas passer à chercher l'erreur!

Sinon, en "googlisant" glibc memory corruption, je suis tombé sur un début d'explication :
http://ubuntuforums.org/showthread.php?t=774897

En particulier, ce passage devrait te donner des pistes :

A glibc memory corruption error usually means that your code overwrote some of the memory
between the allocated ranges handed to you by malloc or calloc (or realloc).
It can also be caused by continuing to access memory that you have freed,
or freeing the same address twice, or freeing an address you never allocated.

[Asohan]

Merci pour ta reponse, j'ai deja googler avant de poster mon message et j'ai vu la meme chose. C'est pour ca que j'accusais le resize.

Je n'utilise aucun new ou delete. Je n'utilise que des variables temporaires toutes simples (float, int, unsigned int, ...). Par contre j'utilise des tableaux.
Mais ce ne sont que des matrices carrees N x N qui, apres changement de N-->N', sont "resizees" a N' x N'.

Sinon, en principe, lorsque l'on cree un objet tres simple dans une boucle, celui ci est detruit et desaloue a la fin de la boucle ? Je dis ca parce que j'ai cree l'objet particule qui me sert a trimbaler les caracteristiques et les fonction propres (methodes) de facon plus intuitives. J'utilise parfois des particule temporaires dans des boucles afin d'eviter de faire n'importe quoi sur mon vrai vector<particule>.
Genre :

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
 
for(int i=0 ; i<n_particules ; i++)
{
    PARTICULE test_particule = particule[i];
 
    ... // calculs.
}

[Aleph69]

Sinon, en principe, lorsque l'on cree un objet tres simple dans une boucle, celui ci est detruit et desaloue a la fin de la boucle ?

Ca dépend.
Pour une création statique entre deux accolades,

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
 
{
// ... du code ...
Objet instance;
// ... du code ...
}

instance sera détruit automatiquement à la fin des accolades.

Par contre, pour une création dynamique,

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
 
{
// ... du code ...
Objet* instance=new Objet();
// ... du code ...
}

instance ne sera pas détruit : c'est à toi de le faire en écrivant

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
 
delete instance;

à un endroit du code, à l'intérieur ou à l'extérieur des accolades.

[seeme]

Merci pour ta reponse, j'ai deja googler avant de poster mon message et j'ai vu la meme chose. C'est pour ca que j'accusais le resize.

Je n'utilise aucun new ou delete. Je n'utilise que des variables temporaires toutes simples (float, int, unsigned int, ...). Par contre j'utilise des tableaux.
Mais ce ne sont que des matrices carrees N x N qui, apres changement de N-->N', sont "resizees" a N' x N'.

Sinon, en principe, lorsque l'on cree un objet tres simple dans une boucle, celui ci est detruit et desaloue a la fin de la boucle ? Je dis ca parce que j'ai cree l'objet particule qui me sert a trimbaler les caracteristiques et les fonction propres (methodes) de facon plus intuitives. J'utilise parfois des particule temporaires dans des boucles afin d'eviter de faire n'importe quoi sur mon vrai vector<particule>.
Genre :

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
 
for(int i=0 ; i<n_particules ; i++)
{
    PARTICULE test_particule = particule[i];
 
    ... // calculs.
}

Si tu as utilisé ces mémoires de manières statiques (sans new/alloc), elles seront détruite en sortie de bloc d'instruction.

[Aleph69]

Une question en passant : quelle est la valeur maximale de n_particules?
Est-ce qu'elle entre bien dans le champ de valeurs d'un int?

Pour ta boucle, j'aurais tendance à te conseiller d'utiliser des types d'entiers non signés.
Etant donné que ton indice de boucle fait référence à un indice de table, je te suggère d'utiliser des variables size_t (c'est fait pour cela) :

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
8
9
 
#include <cstddef> // définition de size_t
// ... du code ...
for(size_t i=0 ; i<n_particules ; i++)
{
    PARTICULE test_particule = particule[i];
 
    ... // calculs.
}

[JolyLoic]

Mais ce ne sont que des matrices carrees N x N qui, apres changement de N-->N', sont "resizees" a N' x N'.

Tu es au courant que d'après le code que tu nous as montré, ton "resize" transforme des matrices N*N en "matrices" dont les N premières lignes ont une taille N, et les N'-N dernières une taille N' :

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
8
9
10
11
12
Ex : Si N == 3 et N' == 5
Avant :
***
***
***
 
Après :
***
***
***
*****
*****

[Aleph69]

Ah oui effectivement!
Tu devrais peut-être utiliser reserve au lieu de resize pour ce que tu veux faire?

[Asohan]

Non pas de soucis pour mes valeurs, j'ai des tests de sécurité assez fréquents.

Par exemple, n_particules est forcement positif (unsigned) et compris entre 1 et 5000 (donc on est loin de la limite des entiers).

Sinon, il me semble que resize effectue d'office un reserve lorsque l'on augmente la taille, non ?

Pour la transformation de tableau, de fait :

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
 
tableau.resize(n_particules, vector<double>(n_particules));

Ce qui, en principe me fait passer par exemple de :
***
***
***
a :
****
****
****
****

Je me trompe ? Si c'est le cas, c'est peut être la raison de ma corruption de mémoire.

[JolyLoic]

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
 
tableau.resize(n_particules, vector<double>(n_particules));

Ce qui, en principe me fait passer par exemple de :
***
***
***
a :
****
****
****
****

Je me trompe ?

Oui. La ligne que tu écris signifie :
Change la taille du tableau (externe) pour qu'il contienne n_particules données. Les nouvelles données qu'il créera seront initialisée à partir d'un vector<double>(n_particules). D'où mon :
***
***
***
****

Beaucoup de STL ont un mode de debug où les erreurs de dépassement d'indice sur les tableaux sont détectées, c'est assez pratique.


Mais pour ton cas, une vraie bibliothèque de matrice te permettrait probablement :
- D'éviter les erreurs de ce genre
- D'avoir du code plus simple à écrire
- D'avoir du code plus performant

Regarde par exemple eigen.

[Asohan]

... JolyLoic, je n'ai qu'un seul mot à dire : MERCI !!!!!!!!!

Je crois que ta solution va enfin résoudre mon problème. D'autant plus que la gestion de mes vecteurs et matrices était assez lourde à force. J'avais essayé la STL et Boost, mais c'était du pareil au même niveau facilité et perf'.

Là avec eigen, je crois vraiment que tous mes problème vont être résolus. Merci infiniment. Je vous tiens au courant de la bonne résolution de mon problème.

J'ai encore deux petite question par rapport à cette solution :
- J'ai essayé valgrind pour verifier mes problème de mémoire, connaissez vous une autre solution pour le genre de problème que j'ai rencontré ?
- Pour ce qui est des performances, de part votre expérience, qu'est-ce qui est le mieux (pour l'algèbre linéaire, les matrices, vecteurs, ...) : STL, Boost, eigen, ou un autre ? Sachant que mes principales opération sont : redimensionnement de matrices et vecteurs, inversion de matrices, et bien sur lecture et écriture.

[Aleph69]

Bonjour,

Pour ce qui est des performances, de part votre expérience, qu'est-ce qui est le mieux (pour l'algèbre linéaire, les matrices, vecteurs, ...) : STL, Boost, eigen, ou un autre ? Sachant que mes principales opération sont : redimensionnement de matrices et vecteurs, inversion de matrices, et bien sur lecture et écriture.

Tu peux aussi allez voir du côté de Blitz++.
Cette bibliothèque utilise une technique (expression templates) pour avoir des calculs plus rapides.
La technique est expliquée dans ce tutoriel :
http://loulou.developpez.com/tutoriels/cpp/metaprog/

Lorsqu'on me parle d'inverser des matrices, j'ai toujours un peu peur...
Tu peux m'expliquer pourquoi tu fais une telle inversion (n'hésite pas à être technique)?
Sinon, pour les opérations plus avancées comme celle-ci, j'aurais tendance à te conseiller de quitter le monde du C++ et d'aller voir du côté de LAPACK.

[Goten]

Tu peux aussi allez voir du côté de Blitz++.

Blitz++ était une pionnière du genre, mais n'est plus maintenue depuis un bon moment.

Depuis y'a eu pas mal d'autres choses sorties (dont eigen) basé aussi sur des ET .

[Aleph69]

Blitz++ était une pionnière du genre, mais n'est plus maintenue depuis un bon moment.

Je n'étais pas au courant, merci pour l'information.

Depuis y'a eu pas mal d'autres choses sorties (dont eigen) basé aussi sur des ET .

Je viens de regarder eigen, ça a l'air pas mal.
Je pense que ça va un peu plus loin que les ET étant donnés les perfs.
La compilation de la bibliothèque est dépendante de l'architecture (taille de cache,etc), non?
Sinon, concernant les benchmarks, à mon avis, ils ont foiré la génération de BLAS et LAPACK avec ATLAS.
Les vitesses mesurées sont beaucoup trop lentes par rapport à GOTO.
D'ailleurs, les vitesses d'ATLAS pour les opérations plus lourdes (cholesky,etc) n'apparaissent pas bizarrement.

Sinon, je viens de penser à un truc.
Je me demande si c'est bien utile de faire des resize à tour de bras si on connaît la valeur max de n_particules?
Il me semble plus malin (rapide) de réserver une fois pour toutes tes matrices avec la valeur max et de mettre à jour la taille courante.

[Asohan]

Merci pour vos réponses. Effectivement j'ai fait quelques recherches moi aussi et eigen semble le plus adapté pour mon cas.

Pourquoi rester en C++ ? parce que dans 2 semaines je vais porter mon code sous OpenCL, et qu'il faut être en C++ ou fortran (à choisir je préfère le C++ car les opérations de calcul physique que je fais sont très complexes au niveau hiérarchique et c'est plus lisible en C++).

Je vois mal Lapack rivaliser avec une inversion de matrice en calcul parallèle avec OpenCL mais si vous voulez je peux vous expliquer : il s'agit d'algèbre linéaire. Je dois calculer des paramètres à chaque pas de temps d'un système dynamique. Ces paramètres sont le Nième colonne de la matrice inversée que je cherche. Est-ce qu'on peut faire autrement ? Oui et non. Oui dans le sens ou j'ai résolu analytiquement l'inversion de matrice N x N dans un cas particulier, non dans le sens ou c'est le cas général qui m'intéresse Et pour corser l'affaire : la matrice n'a aucune symétrie ou structure particulière et aucune valeur nulle.

Donc voila, pour faire simple j'ai besoin du C++ et de librairie plutôt "standards" pour simplifier mon code. Par exemple j'ai écrit la routine d'inversion de matrice moi même avec Gauss-Jordan de manière très optimisée, et franchement, j'arrive à rivaliser niveau perfs avec la plupart des références en C++ (en général je suis même un peu meilleur).

Sinon pour finir, pourquoi faire des tableaux dynamiques alors que de toute façon je me limite à environ 5000x5000 ? Bah parce que j'ai plusieurs tableaux et que malheureusement pour moi, tous les GPUs que j'utilise n'ont pas 1Go de GDDR, mais plutôt 512Mo. Et là dans ce cas : je suis sûr de saturer la mémoire ce qui ralentirait tellement le GPU que dans ce cas OpenCL ne servirait plus à rien.

Voila pour l'histoire. Si des discussion sur OpenCL et le C++ vous intéresse, je vous invite à en discuter sur un autre topic (je m'y colle depuis un petit moment maintenant).

Réponse générale de mon topic : j'ai utilisé la solution de JolyLoic pour l'instant avec mes vector :

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
 
int N = /*valeur désirée*/;
vector< vector<TYPE> > tableau;
 
tableau.resize(N);
for(int i=0 ; i<N ; i++)
    tableau[i].resize(N);

Cela marche parfaitement.

Merci à tous pour votre aide ! Sujet résolu

[Aleph69]

Bonsoir,

Pourquoi rester en C++ ? parce que dans 2 semaines je vais porter mon code sous OpenCL, et qu'il faut être en C++ ou fortran (à choisir je préfère le C++ car les opérations de calcul physique que je fais sont très complexes au niveau hiérarchique et c'est plus lisible en C++).

Je parlais bien sûr d'appeler LAPACK (fortran) depuis ton code C++, pas de changer de langage de programmation.

Je dois calculer des paramètres à chaque pas de temps d'un système dynamique. Ces paramètres sont le Nième colonne de la matrice inversée que je cherche. Est-ce qu'on peut faire autrement ? Oui et non.

Est-ce que tu DOIS faire autrement? Oui!

Il faut résoudre un système linéaire par une factorisation LU.
Le vecteur second membre de ce système est le dernier vecteur de la base canonique de R^N (ne contenant que des zéros sauf le dernier coefficient qui vaut 1).

La solution est la Nième colonne de l'inverse de ta matrice.

[Asohan]

Nan mais de toute façon les interactions à N corps nécessite d'être traiter en parallèle si possible, car le temps de calcul est proportionnel en O(N^2).

L'inversion de matrice pure et simple se fait mieux par Gauss-Jordan, mais dans mon cas tu as tout à fait raison : la décomposition LU me donnera juste le résultat qui m'intéresse, et en difficulté O(N^2) également (au lieu de O(N^3) du pivot de Gauss), je suis vraiment débile de ne pas avoir vu ça avant

Merci encore une fois pour ce très très bon conseil. La rapidité de calcul est cruciale dans mon cas, parce qu'avec mon code de cet après midi, il me fallait environ 1200 jours (oui oui vous avez bien lu) pour achever une seule simulation (avec O(N^3)). Du coup avec LU et en plus une implémentation sous OpenCL (GPU) je devrais gagner un facteur d'environ 100(GPU) x 10(LU) = 1000 ! Ce qui me ramène à environ une seule journée de simulation.

Encore un grand merci