Arcane  4.2.3.0
Documentation développeur
Chargement...
Recherche...
Aucune correspondance
HypreComparer.cc
1// -*- tab-width: 2; indent-tabs-mode: nil; coding: utf-8-with-signature -*-
2//-----------------------------------------------------------------------------
3// Copyright 2000-2026 CEA (www.cea.fr) IFPEN (www.ifpenergiesnouvelles.com)
4// See the top-level COPYRIGHT file for details.
5// SPDX-License-Identifier: Apache-2.0
6//-----------------------------------------------------------------------------
7/*---------------------------------------------------------------------------*/
8/*---------------------------------------------------------------------------*/
9/*
10 * Ce fichier est basé sur le travail sur la bibliothèque AMGCL (version mars 2026)
11 * qui peut être trouvée à https://github.com/ddemidov/amgcl.
12 *
13 * Copyright (c) 2012-2022 Denis Demidov <dennis.demidov@gmail.com>
14 * SPDX-License-Identifier: MIT
15 */
16/*---------------------------------------------------------------------------*/
17/*---------------------------------------------------------------------------*/
18
19/******************************************************************************
20 * Copyright (c) 1998 Lawrence Livermore National Security, LLC and other
21 * HYPRE Project Developers. See the top-level COPYRIGHT file for details.
22 *
23 * SPDX-License-Identifier: (Apache-2.0 OR MIT)
24 ******************************************************************************/
25
26#include "./HypreComparer.h"
27
28#include <stdio.h>
29#include <stdlib.h>
30#include <string.h>
31#include <math.h>
32#include "HYPRE_krylov.h"
33#include "HYPRE.h"
34#include "HYPRE_parcsr_ls.h"
35
36#include <HYPRE_config.h>
37
38#if defined(HYPRE_EXAMPLE_USING_CUDA)
39
40#include <cuda_runtime.h>
41
42#ifndef HYPRE_USING_UNIFIED_MEMORY
43#error *** L'exécution des exemples sur des GPU nécessite une mémoire unifiée. Veuillez reconfigurer et reconstruire avec --enable-unified-memory ***
44#endif
45
46static inline void*
47gpu_malloc(size_t size)
48{
49 void* ptr = NULL;
50 cudaMallocManaged(&ptr, size, cudaMemAttachGlobal);
51 return ptr;
52}
53
54static inline void*
55gpu_calloc(size_t num, size_t size)
56{
57 void *ptr = NULL;
58 cudaMallocManaged(&ptr, num * size, cudaMemAttachGlobal);
59 cudaMemset(ptr, 0, num * size);
60 return ptr;
61}
62
63#define malloc(size) gpu_malloc(size)
64#define calloc(num, size) gpu_calloc(num, size)
65#define free(ptr) ( cudaFree(ptr), ptr = NULL )
66#endif /* #if defined(HYPRE_EXAMPLE_USING_CUDA) */
67
68#include <memory>
69#include <vector>
70#include <iostream>
71
72#include "arccore/base/Convert.h"
73#include "arccore/base/FatalErrorException.h"
74#include "arccore/alina/Profiler.h"
75
76using namespace Arcane;
77
78int hypre_FlexGMRESModifyPCAMGExample(void *precond_data, int iterations,
79 double rel_residual_norm);
80
81#define my_min(a,b) (((a)<(b)) ? (a) : (b))
82
83/*---------------------------------------------------------------------------*/
84/*---------------------------------------------------------------------------*/
85
86HypreComparer::
87~HypreComparer()
88{
89 if (m_need_finalize)
90 MPI_Finalize();
91}
92
93/*---------------------------------------------------------------------------*/
94/*---------------------------------------------------------------------------*/
95
96void HypreComparer::
97solve(int nb_row,
98 std::vector<ptrdiff_t> const& _ptr,
99 std::vector<ptrdiff_t> const& _col,
100 std::vector<double> const& _val,
101 std::vector<double> const& _rhs,
102 std::vector<double>& _x,
103 int argc, char* argv[])
104{
105 auto& prof = Alina::Profiler::globalProfiler();
106 auto t = prof.scoped_tic("Hypre");
107
108 std::cout << "DO_HYPRE nb_row=" << nb_row << "\n";
109 int i;
110 int myid, num_procs;
111
112 int ilower, iupper;
113 int local_size, extra;
114
115 int solver_id;
116 int print_system;
117
118 HYPRE_IJMatrix A;
119 HYPRE_ParCSRMatrix parcsr_A;
120 HYPRE_IJVector b;
121 HYPRE_ParVector par_b;
122 HYPRE_IJVector x;
123 HYPRE_ParVector par_x;
124
125 HYPRE_Solver solver, precond;
126
127 // Initialisation de MPI si demandé
128 if (m_do_mpi_init_and_finalize) {
129 MPI_Init(&argc, &argv);
130 m_need_finalize = true;
131 }
132
133 MPI_Comm_rank(MPI_COMM_WORLD, &myid);
134 MPI_Comm_size(MPI_COMM_WORLD, &num_procs);
135
136 /* Initialisation de HYPRE */
137 HYPRE_Initialize();
138
139 /* Affichage des informations GPU */
140 /* HYPRE_PrintDeviceInfo(); */
141#if defined(HYPRE_USING_GPU)
142 /* Utilise l'implémentation du fournisseur pour SpGEMM */
143 HYPRE_SetSpGemmUseVendor(0);
144#endif
145
146 /* Paramètres par défaut du problème */
147 const int n = nb_row;
148 solver_id = 0;
149 print_system = 0;
150
151 /* Analyse de la ligne de commande */
152 {
153 int arg_index = 0;
154 int print_usage = 0;
155
156 while (arg_index < argc) {
157 if (strcmp(argv[arg_index], "-n") == 0) {
158 arg_index++;
159 //n = atoi(argv[arg_index++]);
160 }
161 else if (strcmp(argv[arg_index], "-solver") == 0) {
162 arg_index++;
163 solver_id = atoi(argv[arg_index++]);
164 }
165 else if (strcmp(argv[arg_index], "-print_system") == 0) {
166 arg_index++;
167 print_system = 1;
168 }
169 else if (strcmp(argv[arg_index], "-help") == 0) {
170 print_usage = 1;
171 break;
172 }
173 else {
174 arg_index++;
175 }
176 }
177
178 if ((print_usage) && (myid == 0)) {
179 printf("\n");
180 printf("Utilisation: %s [<options>]\n", argv[0]);
181 printf("\n");
182 printf(" -n <n> : taille du problème dans chaque direction (par défaut: 33)\n");
183 printf(" -solver <ID> : ID du solveur\n");
184 printf(" 0 - AMG (par défaut) \n");
185 printf(" 1 - AMG-PCG\n");
186 printf(" 8 - ParaSails-PCG\n");
187 printf(" 50 - PCG\n");
188 printf(" 61 - AMG-FlexGMRES\n");
189 printf(" -vis : enregistrer la solution pour la visualisation GLVis\n");
190 printf(" -print_system : imprimer la matrice et le rhs\n");
191 printf("\n");
192 }
193
194 if (print_usage) {
195 return;
196 }
197 }
198
199 // Remplit la valeur nb par ligne.
200 std::vector<int> nb_value_per_row(n);
201 for (int i = 0; i < n; ++i)
202 nb_value_per_row[i] = static_cast<HYPRE_BigInt>(_ptr[i + 1] - _ptr[i]);
203
204 // L'index de colonne est le même que '_col' de la matrice CSR
205 // mais nous faisons une copie si la taille de l'index est différente entre Hypre et Alina.
206 std::vector<HYPRE_BigInt> hypre_column_index(_col.begin(), _col.end());
207
208 // ID de chaque ligne (en séquentiel, c'est le même que l'index)
209 std::vector<HYPRE_Int> hypre_row_index(n);
210 for (int i = 0; i < n; ++i) {
211 hypre_row_index[i] = i;
212 }
213
214 /* Chaque processeur ne connaît que ses propres lignes - la plage est indiquée par ilower
215 et upper. Ici, nous partitionnons les lignes. Nous prenons en compte le fait que
216 N peut ne pas être divisible uniformément par le nombre de processeurs. */
217 local_size = nb_row / num_procs;
218 extra = nb_row - local_size * num_procs;
219
220 ilower = local_size * myid;
221 ilower += my_min(myid, extra);
222
223 iupper = local_size * (myid + 1);
224 iupper += my_min(myid + 1, extra);
225 iupper = iupper - 1;
226 std::cout << "LOWER=" << ilower << " UPPER=" << iupper << "\n";
227 /* Combien de lignes ai-je ? */
228 local_size = iupper - ilower + 1;
229
230 {
231 auto t = prof.scoped_tic("Création de IJMatrix");
232 /* Crée la matrice.
233 Notez qu'il s'agit d'une matrice carrée, nous indiquons donc la taille de la partition
234 de ligne deux fois (car le nombre de lignes = nombre de colonnes) */
235 HYPRE_IJMatrixCreate(MPI_COMM_WORLD, ilower, iupper, ilower, iupper, &A);
236
237 /* Choisis un stockage de format csr parallèle (voir le Manuel de l'utilisateur) */
238 HYPRE_IJMatrixSetObjectType(A, HYPRE_PARCSR);
239
240 /* Initialise avant de définir les coefficients */
241 HYPRE_IJMatrixInitialize(A);
242 }
243
244 // Remplis la matrice.
245 {
246 auto t = prof.scoped_tic("Définition des valeurs de IJMatrix");
247
248 HYPRE_IJMatrixSetValues(A, n,
249 nb_value_per_row.data(),
250 hypre_row_index.data(),
251 hypre_column_index.data(),
252 _val.data());
253 }
254
255 {
256 auto t = prof.scoped_tic("Assemblage de IJMatrix");
257 /* Assemble après avoir défini les coefficients */
258 HYPRE_IJMatrixAssemble(A);
259 }
260
261 /* Remarque : pour le test de petits problèmes, on peut souhaiter lire
262 une matrice au format IJ (pour le format, voir les fichiers de sortie
263 de l'option -print_system).
264 Dans ce cas, on utiliserait la routine suivante:
265 HYPRE_IJMatrixRead( <nom_de_fichier>, MPI_COMM_WORLD,
266 HYPRE_PARCSR, &A );
267 <nom_de_fichier> = IJ.A.out pour lire ce qui a été imprimé par
268 -print_system (les numéros de processeur sont omis).
269 Un appel à HYPRE_IJMatrixRead est une *alternative* à la
270 séquence suivante d'appels HYPRE_IJMatrix:
271 Create, SetObjectType, Initialize, SetValues, et Assemble
272 */
273
274 /* Récupérer l'objet matrice parcsr pour l'utiliser */
275 HYPRE_IJMatrixGetObject(A, (void**)&parcsr_A);
276
277 /* Crée le rhs et la solution */
278 HYPRE_IJVectorCreate(MPI_COMM_WORLD, ilower, iupper, &b);
279 HYPRE_IJVectorSetObjectType(b, HYPRE_PARCSR);
280 HYPRE_IJVectorInitialize(b);
281
282 HYPRE_IJVectorCreate(MPI_COMM_WORLD, ilower, iupper, &x);
283 HYPRE_IJVectorSetObjectType(x, HYPRE_PARCSR);
284 HYPRE_IJVectorInitialize(x);
285
286 /* Définit les valeurs du rhs à h^2 et la solution à zéro */
287 {
288 //double *rhs_values, *x_values;
289 int* rows;
290
291 //rhs_values = (double*)calloc(local_size, sizeof(double));
292 //x_values = (double*)calloc(local_size, sizeof(double));
293 rows = (int*)calloc(local_size, sizeof(int));
294
295 for (i = 0; i < local_size; i++) {
296 //rhs_values[i] = h2;
297 //x_values[i] = 0.0;
298 rows[i] = ilower + i;
299 }
300
301 HYPRE_IJVectorSetValues(b, local_size, rows, _rhs.data());
302 HYPRE_IJVectorSetValues(x, local_size, rows, _x.data());
303
304 //free(x_values);
305 //free(rhs_values);
306 free(rows);
307 }
308
309 HYPRE_IJVectorAssemble(b);
310 /* Tout comme pour la matrice, pour des tests, on peut souhaiter lire un rhs:
311 HYPRE_IJVectorRead( <nom_de_fichier>, MPI_COMM_WORLD,
312 HYPRE_PARCSR, &b );
313 comme alternative à la
314 séquence suivante d'appels HYPRE_IJVectors:
315 Create, SetObjectType, Initialize, SetValues, et Assemble
316 */
317 HYPRE_IJVectorGetObject(b, (void**)&par_b);
318
319 HYPRE_IJVectorAssemble(x);
320 HYPRE_IJVectorGetObject(x, (void**)&par_x);
321
322 /* Imprime le système - les noms de fichiers seront IJ.out.A.XXXXX
323 et IJ.out.b.XXXXX, où XXXXX = ID du processeur */
324 if (print_system) {
325 HYPRE_IJMatrixPrint(A, "IJ.out.A");
326 HYPRE_IJVectorPrint(b, "IJ.out.b");
327 }
328 solver_id = 0;
329 if (auto v = Convert::Type<Int32>::tryParseFromEnvironment("ALINA_HYPRE_SOLVER", true))
330 solver_id = v.value();
331
332 double solver_tolerance = 1.0e-8;
333
334 /* Choisit un solveur et résoudre le système */
335 std::cout << "FIN DE L'ASSEMBLAGE solver_id=" << solver_id << "\n";
336 /* AMG */
337 if (solver_id == 0) {
338 auto t = prof.scoped_tic("HypreSolver AMG");
339 int num_iterations;
340 double final_res_norm;
341
342 /* Crée le solveur */
343 HYPRE_BoomerAMGCreate(&solver);
344
345 /* Définit certains paramètres (Voir le Manuel de Référence pour plus de paramètres) */
346 HYPRE_BoomerAMGSetPrintLevel(solver, 3); /* imprimer les informations de résolution + paramètres */
347 HYPRE_BoomerAMGSetOldDefault(solver); /* Raffinement Falgout avec interpolation classique modifiée */
348 HYPRE_BoomerAMGSetRelaxType(solver, 3); /* Relaxation hybride G-S/Jacobi */
349 HYPRE_BoomerAMGSetRelaxOrder(solver, 1); /* utilise la relaxation C/F */
350 HYPRE_BoomerAMGSetNumSweeps(solver, 1); /* Balayages à chaque niveau */
351 HYPRE_BoomerAMGSetMaxLevels(solver, 20); /* nombre maximum de niveaux */
352 HYPRE_BoomerAMGSetTol(solver, solver_tolerance); /* tolérance de convergence */
353
354 /* Maintenant, configurer et résoudre ! */
355 {
356 auto t = prof.scoped_tic("Configuration AMG");
357 HYPRE_BoomerAMGSetup(solver, parcsr_A, par_b, par_x);
358 }
359 {
360 auto t = prof.scoped_tic("Résolution AMG");
361 HYPRE_BoomerAMGSolve(solver, parcsr_A, par_b, par_x);
362 }
363
364 /* Informations d'exécution - nécessaire si le journal est activé */
365 HYPRE_BoomerAMGGetNumIterations(solver, &num_iterations);
366 HYPRE_BoomerAMGGetFinalRelativeResidualNorm(solver, &final_res_norm);
367 if (myid == 0) {
368 printf("\n");
369 printf("Itérations = %d\n", num_iterations);
370 printf("Norme résiduelle relative finale = %e\n", final_res_norm);
371 printf("\n");
372 }
373
374 /* Détruit le solveur */
375 HYPRE_BoomerAMGDestroy(solver);
376 }
377 /* PCG */
378 else if (solver_id == 50) {
379 auto t = prof.scoped_tic("HypreSolver PCG");
380 int num_iterations;
381 double final_res_norm;
382
383 /* Crée le solveur */
384 HYPRE_ParCSRPCGCreate(MPI_COMM_WORLD, &solver);
385
386 /* Définit certains paramètres (Voir le Manuel de Référence pour plus de paramètres) */
387 HYPRE_PCGSetMaxIter(solver, 1000); /* itérations max */
388 HYPRE_PCGSetTol(solver, solver_tolerance); /* tolérance de convergence */
389 HYPRE_PCGSetTwoNorm(solver, 1); /* utiliser la double norme comme critère d'arrêt */
390 HYPRE_PCGSetPrintLevel(solver, 2); /* affiche les informations d'itération */
391 HYPRE_PCGSetLogging(solver, 1); /* nécessaire pour obtenir les informations d'exécution plus tard */
392
393 /* Maintenant, configurer et résoudre ! */
394 HYPRE_ParCSRPCGSetup(solver, parcsr_A, par_b, par_x);
395 HYPRE_ParCSRPCGSolve(solver, parcsr_A, par_b, par_x);
396
397 /* Informations d'exécution - nécessaire si le journal est activé */
398 HYPRE_PCGGetNumIterations(solver, &num_iterations);
399 HYPRE_PCGGetFinalRelativeResidualNorm(solver, &final_res_norm);
400 if (myid == 0) {
401 printf("\n");
402 printf("Itérations = %d\n", num_iterations);
403 printf("Norme résiduelle relative finale = %e\n", final_res_norm);
404 printf("\n");
405 }
406
407 /* Détruit le solveur */
408 HYPRE_ParCSRPCGDestroy(solver);
409 }
410 /* PCG avec préconditionneur AMG */
411 else if (solver_id == 1) {
412 auto t = prof.scoped_tic("HypreSolver PCG-AMG");
413 int num_iterations;
414 double final_res_norm;
415
416 /* Crée le solveur */
417 HYPRE_ParCSRPCGCreate(MPI_COMM_WORLD, &solver);
418
419 /* Définit certains paramètres (Voir le Manuel de Référence pour plus de paramètres) */
420 HYPRE_PCGSetMaxIter(solver, 1000); /* itérations max */
421 HYPRE_PCGSetTol(solver, solver_tolerance); /* tolérance de convergence */
422 HYPRE_PCGSetTwoNorm(solver, 1); /* utiliser la double norme comme critère d'arrêt */
423 HYPRE_PCGSetPrintLevel(solver, 2); /* affiche les informations de résolution */
424 HYPRE_PCGSetLogging(solver, 1); /* nécessaire pour obtenir les informations d'exécution plus tard */
425
426 /* Maintenant, configurer le préconditionneur AMG et spécifier les paramètres */
427 HYPRE_BoomerAMGCreate(&precond);
428 HYPRE_BoomerAMGSetPrintLevel(precond, 1); /* imprimer les informations de résolution AMG */
429 HYPRE_BoomerAMGSetCoarsenType(precond, 6);
430 HYPRE_BoomerAMGSetOldDefault(precond);
431 HYPRE_BoomerAMGSetRelaxType(precond, 6); /* Hybride G.S./Jacobi symétrique */
432 HYPRE_BoomerAMGSetNumSweeps(precond, 1);
433 HYPRE_BoomerAMGSetTol(precond, 0.0); /* tolérance de convergence zéro */
434 HYPRE_BoomerAMGSetMaxIter(precond, 1); /* faire seulement une itération ! */
435
436 /* Définit le préconditionneur PCG */
437 HYPRE_PCGSetPrecond(solver, (HYPRE_PtrToSolverFcn)HYPRE_BoomerAMGSolve,
438 (HYPRE_PtrToSolverFcn)HYPRE_BoomerAMGSetup, precond);
439
440 /* Maintenant, configurer et résoudre ! */
441 prof.tic("Configuration");
442 HYPRE_ParCSRPCGSetup(solver, parcsr_A, par_b, par_x);
443 prof.toc("Configuration");
444 prof.tic("Résolution");
445 HYPRE_ParCSRPCGSolve(solver, parcsr_A, par_b, par_x);
446 prof.toc("Résolution");
447
448 /* Informations d'exécution - nécessaire si le journal est activé */
449 HYPRE_PCGGetNumIterations(solver, &num_iterations);
450 HYPRE_PCGGetFinalRelativeResidualNorm(solver, &final_res_norm);
451 if (myid == 0) {
452 printf("\n");
453 printf("Itérations = %d\n", num_iterations);
454 printf("Norme résiduelle relative finale = %e\n", final_res_norm);
455 printf("\n");
456 }
457
458 /* Détruit le solveur et le préconditionneur */
459 HYPRE_ParCSRPCGDestroy(solver);
460 HYPRE_BoomerAMGDestroy(precond);
461 }
462 /* PCG avec préconditionneur Parasails */
463 else if (solver_id == 8) {
464 auto t = prof.scoped_tic("HypreSolver PCG - Parasails");
465 int num_iterations;
466 double final_res_norm;
467
468 int sai_max_levels = 1;
469 double sai_threshold = 0.1;
470 double sai_filter = 0.05;
471 int sai_sym = 1;
472
473 /* Crée le solveur */
474 HYPRE_ParCSRPCGCreate(MPI_COMM_WORLD, &solver);
475
476 /* Définit certains paramètres (Voir le Manuel de Référence pour plus de paramètres) */
477 HYPRE_PCGSetMaxIter(solver, 1000); /* itérations max */
478 HYPRE_PCGSetTol(solver, solver_tolerance); /* tolérance de convergence */
479 HYPRE_PCGSetTwoNorm(solver, 1); /* utiliser la double norme comme critère d'arrêt */
480 HYPRE_PCGSetPrintLevel(solver, 2); /* affiche les informations de résolution */
481 HYPRE_PCGSetLogging(solver, 1); /* nécessaire pour obtenir les informations d'exécution plus tard */
482
483 /* Maintenant, configurer le préconditionneur ParaSails et spécifier les paramètres */
484 HYPRE_ParaSailsCreate(MPI_COMM_WORLD, &precond);
485
486 /* Définir certains paramètres (Voir le Manuel de Référence pour plus de paramètres) */
487 HYPRE_ParaSailsSetParams(precond, sai_threshold, sai_max_levels);
488 HYPRE_ParaSailsSetFilter(precond, sai_filter);
489 HYPRE_ParaSailsSetSym(precond, sai_sym);
490 HYPRE_ParaSailsSetLogging(precond, 3);
491
492 /* Définit le préconditionneur PCG */
493 HYPRE_PCGSetPrecond(solver, (HYPRE_PtrToSolverFcn)HYPRE_ParaSailsSolve,
494 (HYPRE_PtrToSolverFcn)HYPRE_ParaSailsSetup, precond);
495
496 /* Maintenant, configurer et résoudre ! */
497 HYPRE_ParCSRPCGSetup(solver, parcsr_A, par_b, par_x);
498 HYPRE_ParCSRPCGSolve(solver, parcsr_A, par_b, par_x);
499
500 /* Informations d'exécution - nécessaire si le journal est activé */
501 HYPRE_PCGGetNumIterations(solver, &num_iterations);
502 HYPRE_PCGGetFinalRelativeResidualNorm(solver, &final_res_norm);
503 if (myid == 0) {
504 printf("\n");
505 printf("Itérations = %d\n", num_iterations);
506 printf("Norme résiduelle relative finale = %e\n", final_res_norm);
507 printf("\n");
508 }
509
510 /* Détruit le solveur et le préconditionneur */
511 HYPRE_ParCSRPCGDestroy(solver);
512 HYPRE_ParaSailsDestroy(precond);
513 }
514 /* Flexible GMRES avec préconditionneur AMG */
515 else if (solver_id == 61) {
516 auto t = prof.scoped_tic("HypreSolver Flexible GMRES - AMG");
517 int num_iterations;
518 double final_res_norm;
519 int restart = 30;
520 int modify = 1;
521
522 /* Crée le solveur */
523 HYPRE_ParCSRFlexGMRESCreate(MPI_COMM_WORLD, &solver);
524
525 /* Définit certains paramètres (Voir le Manuel de Référence pour plus de paramètres) */
526 HYPRE_FlexGMRESSetKDim(solver, restart);
527 HYPRE_FlexGMRESSetMaxIter(solver, 1000); /* itérations max */
528 HYPRE_FlexGMRESSetTol(solver, solver_tolerance); /* tolérance de convergence */
529 HYPRE_FlexGMRESSetPrintLevel(solver, 2); /* affiche les informations de résolution */
530 HYPRE_FlexGMRESSetLogging(solver, 1); /* nécessaire pour obtenir les informations d'exécution plus tard */
531
532 /* Maintenant, configurer le préconditionneur AMG et spécifier les paramètres */
533 HYPRE_BoomerAMGCreate(&precond);
534 HYPRE_BoomerAMGSetPrintLevel(precond, 1); /* imprimer les informations de résolution AMG */
535 HYPRE_BoomerAMGSetCoarsenType(precond, 6);
536 HYPRE_BoomerAMGSetOldDefault(precond);
537 HYPRE_BoomerAMGSetRelaxType(precond, 6); /* Hybride G.S./Jacobi symétrique */
538 HYPRE_BoomerAMGSetNumSweeps(precond, 1);
539 HYPRE_BoomerAMGSetTol(precond, 0.0); /* tolérance de convergence zéro */
540 HYPRE_BoomerAMGSetMaxIter(precond, 1); /* faire seulement une itération ! */
541
542 /* Définit le préconditionneur FlexGMRES */
543 HYPRE_FlexGMRESSetPrecond(solver, (HYPRE_PtrToSolverFcn)HYPRE_BoomerAMGSolve,
544 (HYPRE_PtrToSolverFcn)HYPRE_BoomerAMGSetup, precond);
545
546 if (modify) {
547 /* ceci est un appel optionnel - si vous ne l'appelez pas, hypre_FlexGMRESModifyPCDefault
548 est utilisé - ce qui ne fait rien. Sinon, vous pouvez en définir un propre, similaire à
549 celui utilisé ici */
550 HYPRE_FlexGMRESSetModifyPC(solver, (HYPRE_PtrToModifyPCFcn)hypre_FlexGMRESModifyPCAMGExample);
551 }
552
553 /* Maintenant, configurer et résoudre ! */
554 {
555 auto t = prof.scoped_tic("Configuration FlexGMRES");
556 HYPRE_ParCSRFlexGMRESSetup(solver, parcsr_A, par_b, par_x);
557 }
558 {
559 auto t = prof.scoped_tic("Résolution FlexGMRES");
560 HYPRE_ParCSRFlexGMRESSolve(solver, parcsr_A, par_b, par_x);
561 }
562
563 /* Informations d'exécution - nécessaire si le journal est activé */
564 HYPRE_FlexGMRESGetNumIterations(solver, &num_iterations);
565 HYPRE_FlexGMRESGetFinalRelativeResidualNorm(solver, &final_res_norm);
566 if (myid == 0) {
567 printf("\n");
568 printf("Itérations = %d\n", num_iterations);
569 printf("Norme résiduelle relative finale = %e\n", final_res_norm);
570 printf("\n");
571 }
572
573 /* Détruit le solveur et le préconditionneur */
574 HYPRE_ParCSRFlexGMRESDestroy(solver);
575 HYPRE_BoomerAMGDestroy(precond);
576 }
577 else {
578 if (myid == 0) {
579 ARCCORE_FATAL("ID de solveur '{0}' spécifié invalide.", solver_id);
580 }
581 }
582
583 if (print_system)
584 HYPRE_IJVectorPrint(x, "IJ.out.x");
585
586 /* Nettoyage */
587 HYPRE_IJMatrixDestroy(A);
588 HYPRE_IJVectorDestroy(b);
589 HYPRE_IJVectorDestroy(x);
590
591 /* Finalise HYPRE */
592 HYPRE_Finalize();
593
594 return;
595}
596
597/*---------------------------------------------------------------------------*/
598/*---------------------------------------------------------------------------*/
599/*--------------------------------------------------------------------------
600 hypre_FlexGMRESModifyPCAMGExample -
601
602 Ceci est un exemple (non recommandé)
603 de la façon dont nous pouvons modifier les choses concernant AMG qui
604 affectent la phase de résolution en fonction de la manière dont FlexGMRES progresse...
605 Pour un autre préconditionneur, il pourrait être judicieux de modifier la tolérance..
606 *--------------------------------------------------------------------------*/
607
608int hypre_FlexGMRESModifyPCAMGExample(void* precond_data, [[maybe_unused]] int iterations,
609 double rel_residual_norm)
610{
611
612 if (rel_residual_norm > .1) {
613 HYPRE_BoomerAMGSetNumSweeps((HYPRE_Solver)precond_data, 10);
614 }
615 else {
616 HYPRE_BoomerAMGSetNumSweeps((HYPRE_Solver)precond_data, 1);
617 }
618
619 return 0;
620}
#define ARCCORE_FATAL(...)
Macro envoyant une exception FatalErrorException.
-- tab-width: 2; indent-tabs-mode: nil; coding: utf-8-with-signature --