[C++] Factory d'une classe qui spécialise une classe mère template

Salut,
je développe un framework C++ pour l’optimisation mathématique et je m’arrache un peu les cheveux pour trouver un design qui tient la route au niveau des matrices creuses.

En gros :

• j’ai deux types de représentations de matrices creuses : appelons les A et B ;
• j’ai deux solveurs linéaires, machin et truc. machin travaille exclusivement avec des matrices de type A et truc exclusivement avec des matrices de type B.
• du coup, j’aimerais templater mon code pour que le choix d’un solveur machin ou truc (à l’exécution) fixe le type des toutes les matrices (A ou B) dans tout le reste du code.

J’ai essayé le code (simplifié) suivant :

template <class MatrixType>
class LinearSolver {
   virtual void factorize(const MatrixType& matrix) = 0;
}

class LinearSolverMachin: public LinearSolver<MatrixTypeA> {
   void factorize(const MatrixTypeA& matrix) override;
}

class LinearSolverTruc: public LinearSolver<MatrixTypeB> {
   void factorize(const MatrixTypeB& matrix) override;
}

Le supertype commun à LinearSolverMachin et LinearSolverTruc étant le type template LinearSolver, j’ai templaté la classe LinearSolverFactory elle-même :

template<class MatrixType>
class LinearSolverFactory {
public:
   std::unique_ptr<LinearSolver<MatrixType> > create(const std::string& solver_name) {
     if (solver_name == "machin") {
        return std::make_unique<LinearSolverMachin>();
     }
     else {
        return std::make_unique<LinearSolverTruc>();
     }
   }

Puis à plus haut niveau, j’appelle LinearSolverFactory::create() avec le bon nom de solveur et le bon paramètre de template.

Malheureusement le code ne compile pas (il me dit que le std::unique_ptr<LinearSolverMachin> ne peut pas être converti en std::unique_ptr<LinearSolver<MatrixTypeA> >). Je ne suis pas assez calé en templates pour voir si le problème peut être corrigé.

Merci d’avance pour votre aide

Charlie

Voilà la solution que j’ai reçue sur StackOverflow : il faut spécialiser la factory.

template <typename MatrixType>
struct Factory;

template <> 
struct Factory<MatrixTypeA> {
   std::unique_ptr<LinearSolver<MatrixTypeA> > create(const std::string& name) {
         // select a LinearSolve<MatrixTypeA> and return it
   }
};

template <> 
struct Factory<MatrixTypeB> {
   std::unique_ptr<LinearSolver<MatrixTypeB> > create(const std::string& name) {
         // select a LinearSolver<MatrixTypeB> and return it
   }
};

Salut,
Tu te fais des noeuds au cerveau, tu mélanges 2 approches pour résoudre le même problème de généricité.
Tu n’as pas besoin du template pour redéfinir une méthode virtuelle, ou tu n’as pas besoin d’un héritage pour spécialiser ton template.

Quant aux unique_ptr qui ne comprennent pas l’affectation d’un unique_ptr d’une classe dérivée … c’est un peu dommage que ça ne marche pas, mais effectivement, à cause du template ce sont deux types distincts et donc ce n’est pas possible. En revanche il est possible de construire un "unique_ptr<Base>" autour d’un pointeur "Derive*"

Pour ton design, il faudrait que tu donnes plus d’infos, pourquoi tu veux utiliser de l’héritage, pourquoi tu as besoin d’un objet solver ? etc ..

Tu pourrais très bien avoir linear_solver::machin(A) et linear_solver::truc(B)

Ça veut dire quoi "fixer le type des matrices" à l’exécution ?

Merci pour vos réponses ! J’essaie effectivement de faire un truc compliqué, je vais récapituler ce qu’il se passe dans mon framework : en gros,

• j’ai une classe InteriorPointMethod qui résout des systèmes linéaires ;
• le système linéaire est formé autour de la Hessienne du problème (matrice des dérivées secondes) au point courant. Il est représenté sous forme de matrice creuse ;
• la Hessienne est évaluée par un membre de InteriorPointMethod de type HessianEvaluation. On peut choisir entre Hessienne exact, convexifiée et quasi-Newton ;
• il existe plusieurs formats pour représenter des matrices creuses. J’en ai (au moins) deux (COO et CSC/R) ;
• le système linéaire est résolu par des codes (solveurs) externes (par exemple MA57 ou PARDISO) pour lesquels j’ai écrit des interfaces MA57Solver et PardisoSolver qui héritent de LinearSolver ;
• ces solveurs linéaires travaillent chacun avec un seul format de matrice creuse : COO pour MA57 et CSC/R pour PARDISO ;
• le choix d’un solveur linéaire (MA57 ou PARDISO) détermine donc le format de matrice creuse qui va être utilisé dans la séquence suivante : évaluation de la Hessienne, formation du système linéaire, résolution du système linéaire par le solveur linéaire ;
• le solveur linéaire est choisi à l’exécution lorsque le fichier de configuration du framework est lu.

J’ai donc templaté :

• ma classe LinearSolver afin de spécifier quel format de matrice le solveur peut résoudre (cf premier code du premier message) ;
• la classe HessianEvaluation qui crée la matrice Hessienne directement dans le bon format creux ;
• la classe InteriorPointMethod qui contient comme membres une instance de LinearSolver et une instance de HessianEvaluation.

LinearSolver<MatrixType>, HessianEvaluation<MatrixType> et InteriorPointMethod<MatrixType> doivent être créés avec le même paramètre de template, ce qui permet de "fixer un unique format de matrice creuse dans tout le code". J’ai donc dû templater mes factories (cf deuxième code du premier message) puisque les 3 classes LinearSolver, HessianEvaluation et InteriorPointMethod sont templatées. Je me retrouve donc avec ce drôle de code.

Ca vous inspire un autre design ?

Ca ne répond toujours pas aux questions

Il faudrait commencer par savoir ce que tu fais runtime et compile time, pourquoi tu as besoin d’avoir des objets pour tes solvers, quel est le code commun etc .. Tu dis que t’as un héritage pour tes solvers, mais qu’ils sont externes. On peut choisir entre 3 types de HessianEvaluation, mais tu mets un paramètre template lié au format de matrice.

Si tu fais une classe template mais que tu spécialises tous les types, autant faire des surcharges ou des classes nommées.

Ecris d’abord l’interface que tu veux pour différents cas, et ensuite tu fais ton design à partir de la. La tu pars du principe que tu veux X qui hérite de Y, et Z qui est templaté, mais on ne sait pas pourquoi.

Je lis le choix de solveur et le choix de Hessienne depuis le fichier de config, donc tout se fait à runtime.

COOMatrix et CSCMatrix héritent de Matrix.

Les solveurs sont des codes externes (en C ou Fortran) et j’implémente une interface par solveur (MA57LinearSolver, etc.). La classe mère LinearSolver ressemble à ça :

class LinearSolver {
   virtual void factorize(const SomeMatrixFormat& matrix) = 0;
   virtual void solve(const SomeMatrixFormat& matrix, const std::vector<double>& rhs) = 0;
}

Evidemment, je n’ai pas envie que tous les formats de matrice apparaissent, puisqu’en pratique, un solveur ne travaille que sur un format. Donc pour éviter ça :

class LinearSolver {
   virtual void factorize(const COOMatrix& matrix) = 0;
   virtual void factorize(const CSCMatrix& matrix) = 0;
   virtual void solve(const COOMatrix& matrix, const std::vector<double>& rhs) = 0;
   virtual void solve(const CSCMatrix& matrix, const std::vector<double>& rhs) = 0;
}

ou des signatures sur la classe mère (ce qui oblige à faire des casts dans les méthodes) :

class LinearSolver {
   virtual void factorize(const Matrix& matrix) = 0;
   virtual void solve(const Matrix& matrix, const std::vector<double>& rhs) = 0;
}

je template la classe :

template <class MatrixFormat>
class LinearSolver {
   virtual void factorize(const MatrixFormat& matrix) = 0;
   virtual void solve(const MatrixFormat& matrix, const std::vector<double>& rhs) = 0;
}

Mes classes MA57Solver et PardisoSolver héritent de LinearSolver et spécialisent le paramètre du template avec le format de matrice qu’ils utilisent :

class MA57Solver: public LinearSolver<COOMatrix>

ce qui me permet de définir les méthodes factorize et solve directement sur le bon format de matrice :

class MA57Solver: public LinearSolver<COOMatrix> {
   // allocation de trucs, appels à Fortran
   void factorize(const COOMatrix& matrix) override;
   void solve(const COOMatrix& matrix, const std::vector<double>& rhs) override;
}

Idem pour HessianEvaluation : j’ai une classe mère abstraite, et autant de classes filles que de stratégies (exact, convexifiée, quasi-Newton). Cette classe stocke une Hessienne (une matrice creuse) qu’elle évalue et stocke directement dans le bon format. Donc j’ai templaté la classe mère et les classes filles :

template <class MatrixFormat>
class HessianEvaluation {
   MatrixFormat hessian;
}

template <class MatrixFormat>
class ExactHessianEvaluation : public HessianEvaluation<MatrixFormat> {
   ...
}
...

Ma classe InteriorPointMethod contient des unique_ptr vers le LinearSolver et la HessianEvaluation, donc je continue de templater :

template <class MatrixFormat>
class InteriorPointMethod {
public:
   void compute_direction(...) {
      // evaluate the Hessian at the current point
      this->hessian_evaluation->evaluate_hessian(...);
      
      // form the linear system around this->hessian_evaluation.get_hessian()
      MatrixFormat linear_system = ...

      // compute symbolic factorization of the linear system
      this->linear_solver->factorize(linear_system);
      ...
   }

private:
   std::unique_ptr<LinearSolver<MatrixFormat> > linear_solver;
   std::unique_ptr<HessianEvaluation<MatrixFormat> > hessian_evaluation;
}

Au niveau des combinaisons possibles : InteriorPointMethod doit pouvoir tourner sur le produit cartésiens des stratégies :
MA57 + Hessienne exacte, MA57 + Hessienne convexifiée, MA57 + Hessienne quasi-Newton (tout dans le format COO)
Pardiso + Hessienne exacte, Pardiso + Hessienne convexifiée, Pardiso + Hessienne quasi-Newton (tout dans le format CSC).

Je suis conscient que c’est bizarre de tout templater alors que tout se fait au runtime Ceci dit, si le paramètre du template a 2 valeurs possibles, ça se gère très bien à la volée dans une factory.

J’espère que le design est plus clair

J’ai de plus en plus l’impression que tu essaies de mettre de la généricité là où il n’y en a pas, puisque tu dis toi-même :

LinearSolver<MatrixType>, HessianEvaluation<MatrixType> et InteriorPointMethod<MatrixType> doivent être créés avec le même paramètre de template, ce qui permet de "fixer un unique format de matrice creuse dans tout le code".

Et je ne fais pas de maths assez pointues pour manipuler ce genre d’objets mathématiques tous les jours, mais je me demande aussi si les classes solver ne peuvent pas être réduit à une fonction, car si leur seule responsabilité est de résoudre un système, je crois que ce n’est qu’une séquence d’opération qui s’applique sur les arguments d’entrée, il n’y a pas besoin d’état interne ou je ne sais quoi, si ?

Je rejoins ads00, tu nous expliques comment tu l’as conçu mais pas pourquoi, du coup difficile de faire des choix de design différent
Et de ce que je comprends SolvA n’est capable de traiter que les matriceA et solvB que les matriceB mais tu veux les rendre generiques parce que tu vas avoir un code en fonction de ta config et tu veux faire

switch(config)
{
case configA:
  MatA mat{...};
  SolvA solv{...};
  solv.factorize(A);
  [...]
  break;
case configB:
  MatB mat{...};
  SolvB solv{...};
  solv.factorize(mat);
  [...] // soit exactement le même code que le case configA mais avec des types différents
  break;
}

Si c’est ce cas, alors je pense qu’il faut garder les spécificités des classes et faire une fonction générique qui les associe, du style

class MatA{
[...]
};

class MatB{
[...]
};

class SolvA{
public:
  void factorize(const MatA & mat); // n'est capable de traiter que des MatA
[...]
};

class SolvB{
public:
  void factorize(const MatB & mat); // n'est capable de traiter que des MatB
[...]
};


class SolvC{
public:
  // pourquoi pas un solveur qui est capable de traiter les deux
  void factorize(const MatA & mat);
  void factorize(const MatB & mat);
[...]
};

template<class Matrice, class Solveur>
void solve_system()
{
  Matrice mat{...};
  Solveur solv{...};
  solv.factorize(mat);
  [...] // le code qu'on avait dans le switch tout à l'heure
}

void traiter_config(Config config)
{
  switch(config)
  {
  case configA:
    solve_system<MatA, SolvA>();
    break;
  case configB:
    solve_system<MatB, SolvB>();
    break;
  case configC:
    solve_system<MatA, SolvC>();
    break;
  case configD:
    solve_system<MatB, SolvC>();
    break;
  case configE:
    solve_system<MatB, SolvA>(); // ======> ERREUR à la compilation
    break;
}

EDIT : je n’avais pas lu le dernier message, je comprends ce qui t’as conduit à mettre du template par dessus ton héritage et je pense effectivement que c’était une bonne idée, mais je ne suis pas tout à fait convaincu que l’héritage est utile.

Cool, ta fonction solve_system() paramétrée par la matrice et le solveur est intéressante. Par contre, j’ai oublié de préciser un truc : on va résoudre une séquence de systèmes linéaires, pas un seul. La méthode InteriorPointMethod::compute_direction() va être appelée plusieurs fois. C’est pour ça que le solveur est un membre de InteriorPointMethod.

Quoi qu’il en soit, ta suggestion est intéressante. Je voulais que la classe abstraite LinearSolver déclare factorize(), mais comme tu l’as montré, chaque sous-classe peut déclarer sa propre factorize(mon_type_de_matrice). Je vais y réfléchir.

Mes classes Solver allouent des tableaux et manipulent pas mal de trucs en interne, donc non, elles ne peuvent pas être réduites à une fonction.

Je fais de la recherche dans le domaine des méthodes d’optimisation et je développe un framework qui implémente des composants génériques qu’on peut combiner à la volée. L’idée est d’obtenir automatiquement les combinaisons :

{solveurs MA57, PARDISO} x {stratégies Hessienne} x {interior point, autres méthodes} x {autres trucs}

Si ça vous intéresse (mais ça rentre dans les détails mathématiques), j’ai donné une présentation dans une conf en juillet pour décrire les possibilités du framework : https://www.researchgate.net/publication/353418272

Clairement le projet est ambitieux et je suis un peu limité parfois Mais je suis presque au bout.

Ca fonctionne, merci
C’est futé, d’exposer le type using matrix_type = MatrixType;.