source: issm/trunk/src/c/Solverx/Solverx.cpp@ 3570

Last change on this file since 3570 was 3570, checked in by Mathieu Morlighem, 15 years ago

removed exprintf of all ISSMERROR (not needed anymore)
Fixed Discontinuous Galerkin (missing thickness in inflow BC)"

File size: 3.0 KB
Line 
1/*!\file Solverx
2 * \brief solver
3 */
4
5#include "./Solverx.h"
6
7#include "../shared/shared.h"
8#include "../include/macros.h"
9
10#ifdef HAVE_CONFIG_H
11 #include "config.h"
12#else
13#error "Cannot compile with HAVE_CONFIG_H symbol! run configure first!"
14#endif
15
16
17void Solverx( Vec* puf, Mat Kff, Vec pf, Vec uf0, char* solver_string){
18
19 /*output: */
20 Vec uf=NULL;
21
22 /*intermediary: */
23 int local_m,local_n;
24
25 /*Solver*/
26 KSP ksp=NULL;
27 PC pc=NULL;
28 int iteration_number;
29 PetscTruth flag;
30 int solver_type;
31
32 /*If initial guess for solution exists, use it to create uf, otherwise,
33 * duplicate right hand side so that solution vector has same structure*/
34 if(uf0){
35 VecDuplicate(uf0,&uf); VecCopy(uf0,uf);
36 }
37 else{
38 MatGetLocalSize(Kff,&local_m,&local_n);uf=NewVecFromLocalSize(local_n);
39 }
40
41 /*Before preparing the solver, add options to the options database*/
42 PetscOptionsInsertMultipleString(solver_string);
43
44 /*Process solver_string to see if we are not using special types of external solvers: */
45 PetscOptionsDetermineSolverType(&solver_type,solver_string);
46
47 #if _PETSC_VERSION_ == 2
48 if (solver_type==MUMPSPACKAGE_LU){
49 /*Convert Kff to MATTAIJMUMPS: */
50 MatConvert(Kff,MATAIJMUMPS,MAT_REUSE_MATRIX,&Kff);
51 }
52 if (solver_type==MUMPSPACKAGE_CHOL){
53 /*Convert Kff to MATTSBAIJMUMPS: */
54 MatConvert(Kff,MATSBAIJMUMPS,MAT_REUSE_MATRIX,&Kff);
55 }
56 if (solver_type==SPOOLESPACKAGE_LU){
57 /*Convert Kff to MATTSBAIJMUMPS: */
58 MatConvert(Kff,MATAIJSPOOLES,MAT_REUSE_MATRIX,&Kff);
59 }
60 if (solver_type==SPOOLESPACKAGE_CHOL){
61 /*Convert Kff to MATTSBAIJMUMPS: */
62 MatConvert(Kff,MATSBAIJSPOOLES,MAT_REUSE_MATRIX,&Kff);
63 }
64 if (solver_type==SUPERLUDISTPACKAGE){
65 /*Convert Kff to MATTSBAIJMUMPS: */
66 MatConvert(Kff,MATSUPERLU_DIST,MAT_REUSE_MATRIX,&Kff);
67 }
68 #endif
69
70 /*Prepare solver*/
71 KSPCreate(MPI_COMM_WORLD,&ksp);
72 KSPSetOperators(ksp,Kff,Kff,DIFFERENT_NONZERO_PATTERN);
73 KSPSetFromOptions(ksp);
74
75 #if _PETSC_VERSION_ == 3
76 /*specific solver?: */
77 KSPGetPC(ksp,&pc);
78 if (solver_type==MUMPSPACKAGE_LU){
79 PCFactorSetMatSolverPackage(pc,MAT_SOLVER_MUMPS);
80 }
81 #endif
82
83
84 /*If initial guess for solution, use it, except if we are using the MUMPS direct solver, where any initial
85 * guess will crash Petsc: */
86 if (uf0){
87 if( (solver_type!=MUMPSPACKAGE_LU) && (solver_type!=MUMPSPACKAGE_CHOL) && (solver_type!=SPOOLESPACKAGE_LU)&& (solver_type!=SPOOLESPACKAGE_CHOL)&& (solver_type!=SUPERLUDISTPACKAGE)){
88 KSPSetInitialGuessNonzero(ksp,PETSC_TRUE);
89 }
90 }
91
92
93 #ifdef _ISSM_DEBUG_
94 KSPView(ksp,PETSC_VIEWER_STDOUT_WORLD);
95 KSPGetInitialGuessNonzero(ksp,&flag);
96 _printf_("User provided initial guess solution: %i\n",flag);
97 #endif
98
99 KSPSolve(ksp,pf,uf);
100
101 /*Check convergence*/
102 KSPGetIterationNumber(ksp,&iteration_number);
103 if (iteration_number<0){
104 ISSMERROR("%s%i"," Solver diverged at iteration number: ",-iteration_number);
105 }
106 else{
107 #ifdef _ISSM_DEBUG_
108 _printf_("%s%i%s\n","Solver converged after ",iteration_number," iterations");
109 #endif
110 }
111
112 /*Free ressources:*/
113 KSPFree(&ksp);
114
115
116 /*Assign output pointers:*/
117 *puf=uf;
118}
Note: See TracBrowser for help on using the repository browser.