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

Last change on this file since 3226 was 3226, checked in by Eric.Larour, 15 years ago

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