1 #ifndef __CS_SLES_H__ 2 #define __CS_SLES_H__ 3 4 /*============================================================================ 5 * Sparse Linear Equation Solver driver 6 *============================================================================*/ 7 8 /* 9 This file is part of Code_Saturne, a general-purpose CFD tool. 10 11 Copyright (C) 1998-2021 EDF S.A. 12 13 This program is free software; you can redistribute it and/or modify it under 14 the terms of the GNU General Public License as published by the Free Software 15 Foundation; either version 2 of the License, or (at your option) any later 16 version. 17 18 This program is distributed in the hope that it will be useful, but WITHOUT 19 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS 20 FOR A PARTICULAR PURPOSE. See the GNU General Public License for more 21 details. 22 23 You should have received a copy of the GNU General Public License along with 24 this program; if not, write to the Free Software Foundation, Inc., 51 Franklin 25 Street, Fifth Floor, Boston, MA 02110-1301, USA. 26 */ 27 28 /*----------------------------------------------------------------------------*/ 29 30 /*---------------------------------------------------------------------------- 31 * Local headers 32 *----------------------------------------------------------------------------*/ 33 34 #include "cs_base.h" 35 #include "cs_log.h" 36 #include "cs_halo_perio.h" 37 #include "cs_matrix.h" 38 #include "cs_matrix.h" 39 40 /*----------------------------------------------------------------------------*/ 41 42 BEGIN_C_DECLS 43 44 /*============================================================================ 45 * Macro definitions 46 *============================================================================*/ 47 48 /*============================================================================ 49 * Type definitions 50 *============================================================================*/ 51 52 /*---------------------------------------------------------------------------- 53 * Convergence status 54 *----------------------------------------------------------------------------*/ 55 56 typedef enum { 57 58 CS_SLES_DIVERGED = -3, 59 CS_SLES_BREAKDOWN = -2, 60 CS_SLES_MAX_ITERATION = -1, 61 CS_SLES_ITERATING = 0, 62 CS_SLES_CONVERGED = 1 63 64 } cs_sles_convergence_state_t; 65 66 /* General linear solver context (opaque) */ 67 68 typedef struct _cs_sles_t cs_sles_t; 69 70 /*---------------------------------------------------------------------------- 71 * Function pointer for pre-resolution setup of a linear system solvers's 72 * context. 73 * 74 * This setup may include building a multigrid hierarchy, or a preconditioner. 75 * 76 * Use of this type of function is optional: the context is expected to 77 * maintain state, so that if a cs_sles_solve_t function is called before a 78 * cs_sles_setup_t function, the latter will be called automatically. 79 * 80 * parameters: 81 * context <-> pointer to solver context 82 * name <-- pointer to name of linear system 83 * a <-- matrix 84 * verbosity <-- associated verbosity 85 *----------------------------------------------------------------------------*/ 86 87 typedef void 88 (cs_sles_setup_t) (void *context, 89 const char *name, 90 const cs_matrix_t *a, 91 int verbosity); 92 93 /*---------------------------------------------------------------------------- 94 * Function pointer for resolution of a linear system. 95 * 96 * If the associated cs_sles_setup_t function has not been called before 97 * this function, it will be called automatically. 98 * 99 * The solution context setup by this call (or that of the matching setup 100 * function) will be maintained until the matching cs_sles_free_t function 101 * is called. 102 * 103 * The matrix is not expected to change between successive calls, although 104 * the right hand side may. If the matrix changes, the associated 105 * cs_sles_setup_t or cs_sles_free_t function must be called between 106 * solves. 107 * 108 * The system is considered to have converged when 109 * residue/r_norm <= precision, residue being the L2 norm of a.vx-rhs. 110 * 111 * parameters: 112 * context <-> pointer to solver context 113 * name <-- pointer to name of linear system 114 * a <-- matrix 115 * verbosity <-- associated verbosity 116 * precision <-- solver precision 117 * r_norm <-- residue normalization 118 * n_iter --> number of "equivalent" iterations 119 * residue --> residue 120 * rhs <-- right hand side 121 * vx <-- system solution 122 * aux_size <-- number of elements in aux_vectors 123 * aux_vectors <-- optional working area (internal allocation if NULL) 124 * 125 * returns: 126 * convergence status 127 *----------------------------------------------------------------------------*/ 128 129 typedef cs_sles_convergence_state_t 130 (cs_sles_solve_t) (void *context, 131 const char *name, 132 const cs_matrix_t *a, 133 int verbosity, 134 double precision, 135 double r_norm, 136 int *n_iter, 137 double *residue, 138 const cs_real_t *rhs, 139 cs_real_t *vx, 140 size_t aux_size, 141 void *aux_vectors); 142 143 /*---------------------------------------------------------------------------- 144 * Function pointer for freeing of a linear system's context data. 145 * 146 * Note that this function should free resolution-related data, such as 147 * multigrid hierarchy, preconditioning, and any other temporary arrays or 148 * objects required for resolution, but should not free the whole context, 149 * as info used for logging (especially performance data) should be 150 * maintained. 151 * 152 * parameters: 153 * context <-> pointer to solver context 154 *----------------------------------------------------------------------------*/ 155 156 typedef void 157 (cs_sles_free_t) (void *context); 158 159 /*---------------------------------------------------------------------------- 160 * Function pointer for logging of linear solver setup, 161 * history and performance data. 162 * 163 * This function will be called for each solver when cs_sles_finalize() 164 * is called. 165 * 166 * parameters: 167 * context <-- pointer to solver context 168 * log_type <-- log type 169 *----------------------------------------------------------------------------*/ 170 171 typedef void 172 (cs_sles_log_t) (const void *context, 173 cs_log_t log_type); 174 175 /*---------------------------------------------------------------------------- 176 * Function pointer for creation of a solver context based on the copy 177 * of another. 178 * 179 * The new context copies the settings of the copied context, but not 180 * its setup data and logged info, such as performance data. 181 * 182 * This type of function is optional, but enables associating different 183 * solvers to related systems (to differentiate logging) while using 184 * the same settings by default. 185 * 186 * parameters: 187 * context <-- source context 188 * 189 * returns: 190 * pointer to newly created context 191 *----------------------------------------------------------------------------*/ 192 193 typedef void * 194 (cs_sles_copy_t) (const void *context); 195 196 /*---------------------------------------------------------------------------- 197 * Function pointer for destruction of a linear system solver context. 198 * 199 * This function should free all context data, and will be called for each 200 * system when cs_sles_finalize() is called. 201 * 202 * parameters: 203 * context <-> pointer to solver context 204 *----------------------------------------------------------------------------*/ 205 206 typedef void 207 (cs_sles_destroy_t) (void **context); 208 209 /*---------------------------------------------------------------------------- 210 * Function pointer for handling of non-convegence when solving 211 * a linear system. 212 * 213 * Such a function is optional, and may be used for a variety of purposes, 214 * such as logging, postprocessing, re-trying with different parameters, 215 * aborting the run, or any combination thereof. 216 * 217 * An error handler may be associated with a given solver using 218 * cs_sles_set_error_handler(), in which case it will be called whenever 219 * convergence fails. 220 * 221 * parameters: 222 * sles <-> pointer to solver object 223 * state <-- convergence status 224 * a <-- matrix 225 * rhs <-- Right hand side 226 * vx <-- System solution 227 * 228 * returns: 229 * true if solve should be re-executed, false otherwise 230 *----------------------------------------------------------------------------*/ 231 232 typedef bool 233 (cs_sles_error_handler_t) (cs_sles_t *sles, 234 cs_sles_convergence_state_t state, 235 const cs_matrix_t *a, 236 const cs_real_t *rhs, 237 cs_real_t *vx); 238 239 /*---------------------------------------------------------------------------- 240 * Function pointer for the default definition of a sparse 241 * linear equation solver 242 * 243 * The function may be associated using cs_sles_set_default_define(), so 244 * that it may provide a definition that will be used when 245 * cs_sles_setup() or cs_sles_solve() is used for a system for which 246 * no matching call to cs_sles_define() has been done. 247 * 248 * The function should call cs_sles_define() with arguments f_id 249 * and name, and appropriately chosen function pointers. 250 * 251 * A pointer to the matrix of the system to be solved is also provided, 252 * so that the corresponding information may be used to better choose 253 * defaults. 254 * 255 * parameters: 256 * f_id <-- associated field id, or < 0 257 * name <-- associated name if f_id < 0, or NULL 258 * a <-- Matrix 259 *----------------------------------------------------------------------------*/ 260 261 typedef void 262 (cs_sles_define_t) (int f_id, 263 const char *name, 264 const cs_matrix_t *a); 265 266 /*---------------------------------------------------------------------------- 267 * Function pointer for the default definition of a sparse 268 * linear equation solver's verbosity 269 * 270 * The function may be associated using cs_sles_set_default_verbosity(), so 271 * that it may provide a definition that will be used when 272 * cs_sles_default_verbosity() is called. 273 * 274 * parameters: 275 * f_id <-- associated field id, or < 0 276 * name <-- associated name if f_id < 0, or NULL 277 * 278 * returns: 279 * default verbosity value 280 *----------------------------------------------------------------------------*/ 281 282 typedef int 283 (cs_sles_verbosity_t) (int f_id, 284 const char *name); 285 286 /*============================================================================ 287 * Global variables 288 *============================================================================*/ 289 290 /*============================================================================= 291 * Public function prototypes for Fortran API 292 *============================================================================*/ 293 294 /*============================================================================= 295 * Public function prototypes 296 *============================================================================*/ 297 298 /*---------------------------------------------------------------------------- 299 * \brief Initialize sparse linear equation solver API. 300 *----------------------------------------------------------------------------*/ 301 302 void 303 cs_sles_initialize(void); 304 305 /*---------------------------------------------------------------------------- 306 * \brief Finalize sparse linear equation solver API. 307 *----------------------------------------------------------------------------*/ 308 309 void 310 cs_sles_finalize(void); 311 312 /*----------------------------------------------------------------------------*/ 313 /*! 314 * \brief Log sparse linear equation solver info 315 * 316 * \param[in] log_type log type (CS_LOG_SETUP or CS_LOG_PERFORMANCE) 317 */ 318 /*----------------------------------------------------------------------------*/ 319 320 void 321 cs_sles_log(cs_log_t log_type); 322 323 /*----------------------------------------------------------------------------*/ 324 /*! 325 * \brief Return pointer to linear system object, based on matching field id or 326 * system name. 327 * 328 * If this system did not previously exist, NULL is returned. 329 * 330 * \param[in] f_id associated field id, or < 0 331 * \param[in] name associated name if f_id < 0, or NULL 332 * 333 * \return pointer to associated linear system object, or NULL 334 */ 335 /*----------------------------------------------------------------------------*/ 336 337 cs_sles_t * 338 cs_sles_find(int f_id, 339 const char *name); 340 341 /*----------------------------------------------------------------------------*/ 342 /*! 343 * \brief Return pointer to linear system object, based on matching field id or 344 * system name. 345 * 346 * If this system did not previously exist, it is created and added to 347 * to the list of "known" systems. In this case, it will be usable 348 * only if cs_sles_define() is called for the same field id and name 349 * (in which case calling the present function is redundant), or if 350 * cs_sles_set_sefault_define() has been previously used to define 351 * the default solver policy. 352 * 353 * \param[in] f_id associated field id, or < 0 354 * \param[in] name associated name if f_id < 0, or NULL 355 * 356 * \return pointer to associated linear system object, or NULL 357 */ 358 /*----------------------------------------------------------------------------*/ 359 360 cs_sles_t * 361 cs_sles_find_or_add(int f_id, 362 const char *name); 363 364 /*----------------------------------------------------------------------------*/ 365 /*! 366 * \brief Temporarily replace field id with name for matching calls 367 * to \ref cs_sles_setup, \ref cs_sles_solve, \ref cs_sles_free, and other 368 * operations involving access through a field id. 369 * 370 * \deprecated This function is provided to allow some peculiar 371 * calling sequences, in which \ref cs_equation_iterative_solve_scalar 372 * is called with a nonzero \c ivar value, but specific solver options must 373 * still be set. 374 * In the future, a more consistent mechanism (using a zero \c ivar 375 * value or designing a cleaner method to handle those exceptional cases) 376 * is preferred. As such, only a stack depth of 1 is allowed. 377 * 378 * \param[in] f_id associated field id, or < 0 379 * \param[in] name associated name if f_id < 0, or NULL 380 */ 381 /*----------------------------------------------------------------------------*/ 382 383 void 384 cs_sles_push(int f_id, 385 const char *name); 386 387 /*----------------------------------------------------------------------------*/ 388 /*! 389 * \brief Restore behavior temporarily modified by \ref cs_sles_push. 390 * 391 * \deprecated This function matches \ref cs_sles_push, which is deprecated. 392 * 393 * \param[in] f_id associated field id, or < 0 394 */ 395 /*----------------------------------------------------------------------------*/ 396 397 void 398 cs_sles_pop(int f_id); 399 400 /*----------------------------------------------------------------------------*/ 401 /*! 402 * \brief Define sparse linear equation solver for a given field or 403 * equation name. 404 * 405 * If this system did not previously exist, it is added to the list of 406 * "known" systems. 407 * 408 * The context pointer is used to point to a structure adapted to the function 409 * pointers given here, and combined with those functions, allows using 410 * both built-in, external, or user-defined solvers. 411 * 412 * It is recommended the context type name provided here directly relate 413 * to the associated structure type (for example, "cs_sles_it_t" or 414 * "cs_multigrid_t"). 415 * 416 * \param[in] f_id associated field id, or < 0 417 * \param[in] name associated name if f_id < 0, or NULL 418 * \param[in, out] context pointer to solver context management 419 * structure (cs_sles subsystem becomes owner) 420 * \param[in] type_name context structure or object type name 421 * \param[in] setup_func pointer to system setup function 422 * \param[in] solve_func pointer to system solution function (also 423 * calls setup_func if not done yet or free_func 424 * called since last solve) 425 * \param[in] free_func pointer function freeing system setup 426 * \param[in] log_func pointer to system info logging function 427 (optional, but recommended) 428 * \param[in] copy_func pointer to settings copy function (optional) 429 * \param[in] destroy_func pointer to function destroying solver context 430 * (called with \ref cs_sles_finalize or with a 431 * new call to this function for the same system) 432 * 433 * \return pointer to associated linear system object 434 */ 435 /*----------------------------------------------------------------------------*/ 436 437 cs_sles_t * 438 cs_sles_define(int f_id, 439 const char *name, 440 void *context, 441 const char *type_name, 442 cs_sles_setup_t *setup_func, 443 cs_sles_solve_t *solve_func, 444 cs_sles_free_t *free_func, 445 cs_sles_log_t *log_func, 446 cs_sles_copy_t *copy_func, 447 cs_sles_destroy_t *destroy_func); 448 449 /*----------------------------------------------------------------------------*/ 450 /*! 451 * \brief Set the verbosity for a given linear equation solver. 452 * 453 * This verbosity will be used by cs_sles_setup and cs_sles_solve. 454 * 455 * By default, the verbosity is set to 0, or the value returned by the 456 * function set with cs_sles_set_default_define(). 457 * 458 * \param[in, out] sles pointer to solver object 459 * \param[in] verbosity verbosity level 460 */ 461 /*----------------------------------------------------------------------------*/ 462 463 void 464 cs_sles_set_verbosity(cs_sles_t *sles, 465 int verbosity); 466 467 /*----------------------------------------------------------------------------*/ 468 /*! 469 * \brief Get the verbosity for a given linear equation solver. 470 * 471 * This verbosity will be used by cs_sles_setup and cs_sles_solve. 472 * 473 * \param[in, out] sles pointer to solver object 474 * 475 * \return verbosity level 476 */ 477 /*----------------------------------------------------------------------------*/ 478 479 int 480 cs_sles_get_verbosity(cs_sles_t *sles); 481 482 /*----------------------------------------------------------------------------*/ 483 /*! 484 * \brief Activate postprocessing output for a given linear equation solver. 485 * 486 * This allows the output of the residual at the end of each solution 487 * series, using a single postprocessing writer. 488 * By default, no output is activated. 489 * 490 * \param[in, out] sles pointer to solver object 491 * \param[in] writer_id id of the writer 492 */ 493 /*----------------------------------------------------------------------------*/ 494 495 void 496 cs_sles_set_post_output(cs_sles_t *sles, 497 int writer_id); 498 499 /*----------------------------------------------------------------------------*/ 500 /*! 501 * \brief Return the id of the associated writer if postprocessing output 502 * is active for a given linear equation solver. 503 * 504 * \param[in] sles pointer to solver object 505 * 506 * \return id od associated writer, or 0 507 */ 508 /*----------------------------------------------------------------------------*/ 509 510 int 511 cs_sles_get_post_output(cs_sles_t *sles); 512 513 /*----------------------------------------------------------------------------*/ 514 /*! 515 * \brief Return type name of solver context. 516 * 517 * The returned string is intended to help determine which type is associated 518 * with the void * pointer returned by \ref cs_sles_get_context for a given 519 * solver definition, so as to be able to call additional specific functions 520 * beyond the generic functions assigned to a cs_sles_t object. 521 * 522 * If no type_name string was associated to the solver upon its definition by 523 * \ref cs_sles_define, or it has not been defined yet, the string returned 524 * is "<undefined>". It is recommended the type name provided 525 * \ref cs_sles_define directly relate to the associated structure type 526 * (for example, "cs_sles_it_t" or "cs_multigrid_t"). 527 * 528 * \param[in] sles pointer to solver object 529 * 530 * \return pointer to linear system solver specific type name 531 */ 532 /*----------------------------------------------------------------------------*/ 533 534 const char * 535 cs_sles_get_type(cs_sles_t *sles); 536 537 /*----------------------------------------------------------------------------*/ 538 /*! 539 * \brief Return pointer to solver context structure pointer. 540 * 541 * The context structure depends on the type of solver used, which may in 542 * turn be determined by the string returned by cs_sles_get_type(). 543 * If may be used by appropriate functions specific to that type. 544 * 545 * \param[in] sles pointer to solver object 546 * 547 * \return pointer to solver-specific linear system info and context 548 */ 549 /*----------------------------------------------------------------------------*/ 550 551 void * 552 cs_sles_get_context(cs_sles_t *sles); 553 554 /*----------------------------------------------------------------------------*/ 555 /*! 556 * \brief Return field id associated with a given sparse linear equation solver. 557 * 558 * \param[in] sles pointer to solver object 559 * 560 * \return associated field id (or -1 if defined by name) 561 */ 562 /*----------------------------------------------------------------------------*/ 563 564 int 565 cs_sles_get_f_id(const cs_sles_t *sles); 566 567 /*----------------------------------------------------------------------------*/ 568 /*! 569 * \brief Return name associated with a given sparse linear equation solver. 570 * 571 * This is simply a utility function which will return its name argument 572 * if f_id < 0, and the associated field's name or label otherwise. 573 * 574 * \param[in] sles pointer to solver object 575 * 576 * \return pointer to associated linear system object name 577 */ 578 /*----------------------------------------------------------------------------*/ 579 580 const char * 581 cs_sles_get_name(const cs_sles_t *sles); 582 583 /*----------------------------------------------------------------------------*/ 584 /*! 585 * \brief Setup sparse linear equation solver. 586 * 587 * Use of this function is optional: if a \ref cs_sles_solve is called 588 * for the same system before this function is called, the latter will be 589 * called automatically. 590 * 591 * If no options were previously provided for the matching system, 592 * default options will be used. 593 * 594 * \param[in, out] sles pointer to solver object 595 * \param[in] a matrix 596 */ 597 /*----------------------------------------------------------------------------*/ 598 599 void 600 cs_sles_setup(cs_sles_t *sles, 601 const cs_matrix_t *a); 602 603 /*----------------------------------------------------------------------------*/ 604 /*! 605 * \brief General sparse linear system resolution. 606 * 607 * If no options were previously provided for the matching system, 608 * default options will be used. 609 * 610 * Note that if \ref cs_sles_setup was previously called for this 611 * system, and \ref cs_sles_free has not been called since, the matrix 612 * provided should be the same. The optional separation between the 613 * two stages is intended to allow amortizing the cost of setup 614 * over multiple solutions. 615 * 616 * The system is considered to have converged when 617 * residue/r_norm <= precision, residue being the L2 norm of a.vx-rhs. 618 * 619 * \param[in, out] sles pointer to solver object 620 * \param[in] a matrix 621 * \param[in] precision solver precision 622 * \param[in] r_norm residue normalization 623 * \param[out] n_iter number of "equivalent" iterations 624 * \param[out] residue residue 625 * \param[in] rhs right hand side 626 * \param[in, out] vx system solution 627 * \param[in] aux_size size of aux_vectors (in bytes) 628 * \param aux_vectors optional working area 629 * (internal allocation if NULL) 630 * 631 * \return convergence state 632 */ 633 /*----------------------------------------------------------------------------*/ 634 635 cs_sles_convergence_state_t 636 cs_sles_solve(cs_sles_t *sles, 637 const cs_matrix_t *a, 638 double precision, 639 double r_norm, 640 int *n_iter, 641 double *residue, 642 const cs_real_t *rhs, 643 cs_real_t *vx, 644 size_t aux_size, 645 void *aux_vectors); 646 647 /*----------------------------------------------------------------------------*/ 648 /*! 649 * \brief Free sparse linear equation solver setup. 650 * 651 * This function frees resolution-related data, such as multigrid hierarchy, 652 * preconditioning, and any other temporary arrays or objects required for 653 * resolution, but maintains context information such as that used for 654 * logging (especially performance data). 655 * 656 * \param[in, out] sles pointer to solver object 657 */ 658 /*----------------------------------------------------------------------------*/ 659 660 void 661 cs_sles_free(cs_sles_t *sles); 662 663 /*----------------------------------------------------------------------------*/ 664 /*! 665 * \brief Copy the definition of a sparse linear equation solver to another. 666 * 667 * The intended use of this function is to allow associating different 668 * solvers to related systems, so as to differentiate logging, while using 669 * the same settings by default. 670 * 671 * If the source solver does not provide a \ref cs_sles_copy_t function, 672 * No modification is done to the solver. If the copy function is available, 673 * the context is copied, as are the matching function pointers. 674 * 675 * If previous settings have been defined and used, they are saved as 676 * per \ref cs_sles_define. 677 * 678 * \param[in, out] dest pointer to destination solver object 679 * \param[in] src pointer to source solver object 680 * 681 * \return 0 in case of success, 1 in case of failure 682 */ 683 /*----------------------------------------------------------------------------*/ 684 685 int 686 cs_sles_copy(cs_sles_t *dest, 687 const cs_sles_t *src); 688 689 /*----------------------------------------------------------------------------*/ 690 /*! 691 * \brief Associate a convergence error handler to a given sparse linear 692 * equation solver. 693 * 694 * The error will be called whenever convergence fails. To dissassociate 695 * the error handler, this function may be called with \p handler = NULL. 696 * 697 * The association will only be successful if the matching solver 698 * has already been defined. 699 * 700 * \param[in, out] sles pointer to solver object 701 * \param[in] error_handler_func pointer to convergence error 702 * handler function 703 */ 704 /*----------------------------------------------------------------------------*/ 705 706 void 707 cs_sles_set_error_handler(cs_sles_t *sles, 708 cs_sles_error_handler_t *error_handler_func); 709 710 /*----------------------------------------------------------------------------*/ 711 /*! 712 * \brief Return pointer to default sparse linear solver definition function. 713 * 714 * The associated function will be used to provide a definition when 715 * \ref cs_sles_setup or \ref cs_sles_solve is used for a system for which no 716 * matching call to \ref cs_sles_define has been done. 717 * 718 * \return define_func pointer to default definition function 719 */ 720 /*----------------------------------------------------------------------------*/ 721 722 cs_sles_define_t * 723 cs_sles_get_default_define(void); 724 725 /*----------------------------------------------------------------------------*/ 726 /*! 727 * \brief Set default sparse linear solver definition function. 728 * 729 * The provided function will be used to provide a definition when 730 * \ref cs_sles_setup or \ref cs_sles_solve is used for a system for which no 731 * matching call to \ref cs_sles_define has been done. 732 * 733 * \param[in] define_func pointer to default definition function 734 */ 735 /*----------------------------------------------------------------------------*/ 736 737 void 738 cs_sles_set_default_define(cs_sles_define_t *define_func); 739 740 /*----------------------------------------------------------------------------*/ 741 /*! 742 * \brief Set default verbosity definition function. 743 * 744 * The provided function will be used to define the verbosity when 745 * \ref cs_sles_find_or_add is called. 746 * 747 * \param[in] verbosity_func pointer to default verbosity function 748 */ 749 /*----------------------------------------------------------------------------*/ 750 751 void 752 cs_sles_set_default_verbosity(cs_sles_verbosity_t *verbosity_func); 753 754 /*----------------------------------------------------------------------------*/ 755 /*! 756 * \brief Output default post-processing data for failed system convergence. 757 * 758 * \param[in] name variable name 759 * \param[in] mesh_id id of error output mesh, or 0 if none 760 * \param[in] a linear equation matrix 761 * \param[in] rhs right hand side 762 * \param[in, out] vx current system solution 763 */ 764 /*----------------------------------------------------------------------------*/ 765 766 void 767 cs_sles_post_error_output_def(const char *name, 768 int mesh_id, 769 const cs_matrix_t *a, 770 const cs_real_t *rhs, 771 cs_real_t *vx); 772 773 /*----------------------------------------------------------------------------*/ 774 /*! 775 * \brief Output post-processing variable related to system convergence. 776 * 777 * \param[in] name variable name 778 * \param[in] mesh_id id of error output mesh, or 0 if none 779 * \param[in] location_id mesh location id (cells or vertices) 780 * \param[in] writer_id id of specified associated writer, or 781 * \ref CS_POST_WRITER_ALL_ASSOCIATED for all 782 * \param[in] diag_block_size block size for diagonal 783 * \param[in, out] var variable values 784 */ 785 /*----------------------------------------------------------------------------*/ 786 787 void 788 cs_sles_post_output_var(const char *name, 789 int mesh_id, 790 int location_id, 791 int writer_id, 792 int diag_block_size, 793 cs_real_t var[]); 794 795 /*----------------------------------------------------------------------------*/ 796 /*! 797 * \brief Return base name associated to a field id, name couple. 798 * 799 * This is simply a utility function which will return its name argument 800 * if f_id < 0, and the associated field's name or label otherwise. 801 * 802 * \param[in] f_id associated field id, or < 0 803 * \param[in] name associated name if f_id < 0, or NULL 804 * 805 * \return pointer to associated linear system object, or NULL 806 */ 807 /*----------------------------------------------------------------------------*/ 808 809 const char * 810 cs_sles_base_name(int f_id, 811 const char *name); 812 813 /*----------------------------------------------------------------------------*/ 814 /*! 815 * \brief Return name associated to a field id, name couple. 816 * 817 * \param[in] f_id associated field id, or < 0 818 * \param[in] name associated name if f_id < 0, or NULL 819 * 820 * \return pointer to name associated to the field id, name couple 821 */ 822 /*----------------------------------------------------------------------------*/ 823 824 const char * 825 cs_sles_name(int f_id, 826 const char *name); 827 828 /*----------------------------------------------------------------------------*/ 829 830 END_C_DECLS 831 832 #endif /* __CS_SLES_H__ */ 833