solver_petsc.cpp 14.6 KB
Newer Older
Kirill Terekhov's avatar
Kirill Terekhov committed
1
#define _CRT_SECURE_NO_WARNINGS
Kirill Terekhov's avatar
Kirill Terekhov committed
2 3 4 5 6 7 8 9 10 11 12 13 14
#include "inmost_solver.h"
#if defined(USE_SOLVER)
#if defined(USE_SOLVER_PETSC)
#include <petsc.h>
#include "solver_petsc.h"
#define PETSC_SUCCESS 0


bool SolverIsInitializedPetsc()
{
	PetscBool isInitialized;
	PetscErrorCode ierr = PetscInitialized(&isInitialized);
	if( ierr != PETSC_SUCCESS ) throw INMOST::ErrorInSolver;
Kirill Terekhov's avatar
Kirill Terekhov committed
15
	return (isInitialized == PETSC_TRUE);
Kirill Terekhov's avatar
Kirill Terekhov committed
16 17 18 19 20 21 22
}

void SolverInitializePetsc(int * argc,char *** argv, const char * file_options)
{
	if( !SolverIsInitializedPetsc() )
	{
		PetscErrorCode ierr = PetscInitialize(argc,argv,file_options,"solver");
Kirill Terekhov's avatar
Kirill Terekhov committed
23 24 25
		//prevent petsc from catching signals, petsc error handling is misleading for 
		//unexperienced user and results in the opinion that errors emerge from petsc
		PetscPopSignalHandler(); 
Kirill Terekhov's avatar
Kirill Terekhov committed
26 27 28 29 30 31 32 33 34
		if( ierr != PETSC_SUCCESS ) throw INMOST::ErrorInSolver;
	}
}

bool SolverIsFinalizedPetsc()
{
	PetscBool isFinalized;
	PetscErrorCode ierr = PetscFinalized(&isFinalized);
	if( ierr != PETSC_SUCCESS ) throw INMOST::ErrorInSolver;
Kirill Terekhov's avatar
Kirill Terekhov committed
35
	return (isFinalized == PETSC_TRUE);
Kirill Terekhov's avatar
Kirill Terekhov committed
36 37
}

38 39 40 41 42
void SolverFinalizePetsc() {
    if (!SolverIsFinalizedPetsc()) {
        PetscErrorCode ierr = PetscFinalize();
        if (ierr != PETSC_SUCCESS) throw INMOST::ErrorInSolver;
    }
Kirill Terekhov's avatar
Kirill Terekhov committed
43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58
}

void MatrixDestroyDataPetsc(void ** data)
{
	Mat * m = static_cast<Mat *>(*data);
	PetscErrorCode ierr = MatDestroy(m);
	if( ierr != PETSC_SUCCESS ) throw INMOST::ErrorInSolver;
	delete m;
	*data = NULL;
}


void MatrixInitDataPetsc(void ** data, INMOST_MPI_Comm comm, const char * name)
{
	PetscErrorCode ierr;
	if( data == NULL ) throw INMOST::DataCorruptedInSolver;
59
	if( *data != NULL )
Kirill Terekhov's avatar
Kirill Terekhov committed
60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 167 168 169 170 171 172 173 174 175 176 177 178 179 180 181 182 183 184 185 186 187 188 189 190 191 192 193 194 195 196 197 198 199 200 201 202 203 204 205 206 207 208 209 210 211 212 213 214 215 216 217
		MatrixDestroyDataPetsc(data);
	*data = static_cast<void *>(new Mat);
#if !defined(USE_MPI)
	ierr = MatCreate(PETSC_COMM_WORLD,static_cast<Mat *>(*data));
#else
	ierr = MatCreate(comm,static_cast<Mat *>(*data));
#endif
	if( ierr != PETSC_SUCCESS ) throw INMOST::ErrorInSolver;
	ierr = MatSetOptionsPrefix(*static_cast<Mat *>(*data),name);
	if( ierr != PETSC_SUCCESS ) throw INMOST::ErrorInSolver;
	ierr = MatSetFromOptions(*static_cast<Mat *>(*data));
	if( ierr != PETSC_SUCCESS ) throw INMOST::ErrorInSolver;
}

void MatrixCopyDataPetsc(void ** data, void * other_data)
{
	PetscErrorCode ierr;
	if( data == NULL || other_data == NULL ) throw INMOST::DataCorruptedInSolver;
	if( *data != NULL )
		MatrixDestroyDataPetsc(data);
	*data = static_cast<void *>(new Mat);
	ierr = MatConvert(*static_cast<Mat *>(other_data),MATSAME,MAT_INITIAL_MATRIX,static_cast<Mat *>(*data));
	if( ierr != PETSC_SUCCESS ) throw INMOST::ErrorInSolver;
}

void MatrixAssignDataPetsc(void * data, void * other_data)
{
	PetscErrorCode ierr;
	if( data == NULL || other_data == NULL ) throw INMOST::DataCorruptedInSolver;
	ierr = MatCopy(*static_cast<Mat *>(other_data),*static_cast<Mat *>(data),DIFFERENT_NONZERO_PATTERN);
	if( ierr != PETSC_SUCCESS ) throw INMOST::ErrorInSolver;
}


void MatrixPreallocatePetsc(void * data, int local_size, int global_size, int * diag_nonzeroes, int * off_diag_nonzeroes)
{
	PetscErrorCode ierr;
	Mat * m = static_cast<Mat *>(data);
	ierr = MatSetSizes(*m,local_size,local_size,global_size,global_size);
	if( ierr != PETSC_SUCCESS ) throw INMOST::ErrorInSolver;
	ierr = MatXAIJSetPreallocation(*m,1,diag_nonzeroes,off_diag_nonzeroes,PETSC_NULL,PETSC_NULL);
	if( ierr != PETSC_SUCCESS ) throw INMOST::ErrorInSolver;
}

void MatrixFillPetsc(void * data, int row, int cols, int * col_positions, double * col_values)
{
	PetscErrorCode ierr;
	Mat * m = static_cast<Mat *>(data);
	ierr = MatSetValues(*m,1,&row,cols,col_positions,col_values,INSERT_VALUES);
	if( ierr != PETSC_SUCCESS ) throw INMOST::ErrorInSolver;
}

void MatrixFinalizePetsc(void * data)
{
	PetscErrorCode ierr;
	Mat * m = static_cast<Mat *>(data);
	ierr = MatAssemblyBegin(*m,MAT_FINAL_ASSEMBLY);
	if( ierr != PETSC_SUCCESS ) throw INMOST::ErrorInSolver;
	ierr = MatAssemblyEnd(*m,MAT_FINAL_ASSEMBLY);
	if( ierr != PETSC_SUCCESS ) throw INMOST::ErrorInSolver;
}


void VectorDestroyDataPetsc(void ** data)
{
	Vec * m = static_cast<Vec *>(*data);
	PetscErrorCode ierr = VecDestroy(m);
	if( ierr != PETSC_SUCCESS ) throw INMOST::ErrorInSolver;
	delete m;
	*data = NULL;
}


void VectorInitDataPetsc(void ** data, INMOST_MPI_Comm comm, const char * name)
{
	PetscErrorCode ierr;
	if( data == NULL ) throw INMOST::DataCorruptedInSolver;
	if( *data != NULL ) 
		VectorDestroyDataPetsc(data);
	*data = static_cast<void *>(new Vec);
#if !defined(USE_MPI)
	ierr = VecCreate(PETSC_COMM_WORLD,static_cast<Vec *>(*data));
#else
	ierr = VecCreate(comm,static_cast<Vec *>(*data));
#endif
	if( ierr != PETSC_SUCCESS ) throw INMOST::ErrorInSolver;
	ierr = VecSetOptionsPrefix(*static_cast<Vec *>(*data),name);
	if( ierr != PETSC_SUCCESS ) throw INMOST::ErrorInSolver;
	ierr = VecSetFromOptions(*static_cast<Vec *>(*data));
	if( ierr != PETSC_SUCCESS ) throw INMOST::ErrorInSolver;
}

void VectorCopyDataPetsc(void ** data, void * other_data)
{
	PetscErrorCode ierr;
	if( data == NULL || other_data == NULL ) throw INMOST::DataCorruptedInSolver;
	if( *data != NULL )
		VectorDestroyDataPetsc(data);
	*data = static_cast<void *>(new Vec);
	ierr = VecDuplicate(*static_cast<Vec *>(other_data),static_cast<Vec *>(*data));
	if( ierr != PETSC_SUCCESS ) throw INMOST::ErrorInSolver;
	ierr = VecCopy(*static_cast<Vec *>(other_data),*static_cast<Vec *>(*data));
	if( ierr != PETSC_SUCCESS ) throw INMOST::ErrorInSolver;
}

void VectorAssignDataPetsc(void * data, void * other_data)
{
	PetscErrorCode ierr;
	if( data == NULL || other_data == NULL ) throw INMOST::DataCorruptedInSolver;
	ierr = VecCopy(*static_cast<Vec *>(other_data),*static_cast<Vec *>(data));
	if( ierr != PETSC_SUCCESS ) throw INMOST::ErrorInSolver;
}



void VectorPreallocatePetsc(void * data, int local_size, int global_size)
{
	Vec * m = static_cast<Vec *>(data);
	PetscErrorCode ierr = VecSetSizes(*m,local_size,global_size);
	if( ierr != PETSC_SUCCESS ) throw INMOST::ErrorInSolver;
}
void VectorFillPetsc(void * data, int size, int * positions, double * values)
{
	Vec * m = static_cast<Vec *>(data);
	PetscErrorCode ierr = VecSetValues(*m,size,positions,values,INSERT_VALUES);
	if( ierr != PETSC_SUCCESS ) throw INMOST::ErrorInSolver;
}
void VectorLoadPetsc(void * data, int size, int * positions, double * values)
{
	Vec * m = static_cast<Vec *>(data);
	PetscErrorCode ierr = VecGetValues(*m,size,positions,values);
	if( ierr != PETSC_SUCCESS ) throw INMOST::ErrorInSolver;
}

void VectorFinalizePetsc(void * data)
{
	PetscErrorCode ierr;
	Vec * m = static_cast<Vec *>(data);
	ierr = VecAssemblyBegin(*m);
	if( ierr != PETSC_SUCCESS ) throw INMOST::ErrorInSolver;
	ierr = VecAssemblyEnd(*m);
	if( ierr != PETSC_SUCCESS ) throw INMOST::ErrorInSolver;
}

void SolverDestroyDataPetsc(void ** data)
{
	KSP * m = static_cast<KSP *>(*data);
	PetscErrorCode ierr = KSPDestroy(m);
	if( ierr != PETSC_SUCCESS ) throw INMOST::ErrorInSolver;
	delete m;
	*data = NULL;
}


void SolverInitDataPetsc(void ** data, INMOST_MPI_Comm comm, const char * name)
{
	PetscErrorCode ierr;
	if( data == NULL ) throw INMOST::DataCorruptedInSolver;
218
	if( *data != NULL )
Kirill Terekhov's avatar
Kirill Terekhov committed
219 220 221 222 223 224 225 226 227 228 229 230 231 232 233 234 235 236 237 238 239 240 241 242 243 244 245 246 247 248 249 250 251 252 253 254 255 256 257 258 259 260 261 262 263 264 265 266 267 268 269 270 271 272 273 274 275
		SolverDestroyDataPetsc(data);
	*data = static_cast<void *>(new KSP);
#if !defined(USE_MPI)
	ierr = KSPCreate(PETSC_COMM_WORLD,static_cast<KSP *>(*data));
#else
	ierr = KSPCreate(comm,static_cast<KSP *>(*data));
#endif
	if( ierr != PETSC_SUCCESS ) throw INMOST::ErrorInSolver;
	ierr = KSPSetOptionsPrefix(*static_cast<KSP *>(*data),name);
	if( ierr != PETSC_SUCCESS ) throw INMOST::ErrorInSolver;
	ierr = KSPSetFromOptions(*static_cast<KSP *>(*data));
	if( ierr != PETSC_SUCCESS ) throw INMOST::ErrorInSolver;
}

void SolverCopyDataPetsc(void ** data, void * other_data, INMOST_MPI_Comm comm)
{
	PetscErrorCode ierr;
	if( data == NULL || other_data == NULL ) throw INMOST::DataCorruptedInSolver;
	if( *data != NULL )
		SolverDestroyDataPetsc(data);
	*data = static_cast<void *>(new KSP);
#if !defined(USE_MPI)
	ierr = KSPCreate(PETSC_COMM_WORLD,static_cast<KSP *>(*data));
#else
	ierr = KSPCreate(comm,static_cast<KSP *>(*data));
#endif
	if( ierr != PETSC_SUCCESS ) throw INMOST::ErrorInSolver;
	KSP * mine = static_cast<KSP *>(*data);
	KSP * other = static_cast<KSP *>(other_data);
	char * prefix;
	ierr = KSPGetOptionsPrefix(*other,const_cast<const char **>(&prefix));
	if( ierr != PETSC_SUCCESS ) throw INMOST::ErrorInSolver;
	ierr = KSPSetOptionsPrefix(*mine,prefix);
	if( ierr != PETSC_SUCCESS ) throw INMOST::ErrorInSolver;
	ierr = KSPSetFromOptions(*mine);
	if( ierr != PETSC_SUCCESS ) throw INMOST::ErrorInSolver;
}
void SolverAssignDataPetsc(void * data, void * other_data)
{
	PetscErrorCode ierr;
	if( data == NULL || other_data == NULL ) throw INMOST::DataCorruptedInSolver;
	KSP * mine = static_cast<KSP *>(data);
	KSP * other = static_cast<KSP *>(other_data);
	char * prefix;
	ierr = KSPGetOptionsPrefix(*other,const_cast<const char **>(&prefix));
	if( ierr != PETSC_SUCCESS ) throw INMOST::ErrorInSolver;
	ierr = KSPSetOptionsPrefix(*mine,prefix);
	if( ierr != PETSC_SUCCESS ) throw INMOST::ErrorInSolver;
	ierr = KSPSetFromOptions(*mine);
	if( ierr != PETSC_SUCCESS ) throw INMOST::ErrorInSolver;
}


void SolverSetMatrixPetsc(void * data, void * matrix_data, bool same_pattern, bool reuse_preconditioner)
{
	KSP * m = static_cast<KSP *>(data);
	Mat * mat = static_cast<Mat *>(matrix_data);
Kirill Terekhov's avatar
Kirill Terekhov committed
276 277 278 279 280 281
	PetscErrorCode ierr;
	if( reuse_preconditioner ) 
		ierr = KSPSetReusePreconditioner(*m,PETSC_TRUE);
	else 
		ierr = KSPSetReusePreconditioner(*m,PETSC_FALSE);
	ierr = KSPSetOperators(*m,*mat,*mat);//,reuse_preconditioner? SAME_PRECONDITIONER : (same_pattern? SAME_NONZERO_PATTERN : DIFFERENT_NONZERO_PATTERN));
Kirill Terekhov's avatar
Kirill Terekhov committed
282 283 284 285 286 287 288 289 290 291 292 293 294 295 296 297 298 299 300 301 302 303
	if( ierr != PETSC_SUCCESS ) throw INMOST::ErrorInSolver;
}

bool SolverSolvePetsc(void * data, void * rhs_data, void * sol_data)
{
	PetscErrorCode ierr;
	PetscBool haveksp;
	Vec * rhs = static_cast<Vec *>(rhs_data), * sol = static_cast<Vec *>(sol_data);
	KSP * ksp = static_cast<KSP *>(data);
	bool guess = true;
	char typeksp[2048];
	char * prefix;
	ierr = KSPGetOptionsPrefix(*ksp,const_cast<const char **>(&prefix));
	if( ierr != PETSC_SUCCESS ) throw INMOST::ErrorInSolver;
	ierr = PetscOptionsGetString(prefix,"-ksp_type",typeksp,2048,&haveksp);
	if( ierr != PETSC_SUCCESS ) throw INMOST::ErrorInSolver;
	if( haveksp && !strcmp(typeksp,"preonly") ) guess = false;
	if( guess )
	{
		ierr = KSPSetInitialGuessNonzero(*ksp,PETSC_TRUE);
		if( ierr != PETSC_SUCCESS ) throw INMOST::ErrorInSolver;
	}
Kirill Terekhov's avatar
Kirill Terekhov committed
304 305 306
  //MatNullSpace nullsp;
  //ierr = MatNullSpaceCreate(PETSC_COMM_WORLD, PETSC_TRUE, 0, PETSC_NULL,&nullsp);
  //ierr = KSPSetNullSpace(*ksp,nullsp);
Kirill Terekhov's avatar
Kirill Terekhov committed
307 308 309 310 311 312 313 314 315 316 317 318 319 320 321 322 323 324 325 326 327 328 329 330 331
	ierr = KSPSolve(*ksp,*rhs,*sol);
	if( ierr != PETSC_SUCCESS ) throw INMOST::ErrorInSolver;
	KSPConvergedReason reason;
	ierr = KSPGetConvergedReason(*ksp,&reason);
	if( ierr != PETSC_SUCCESS ) throw INMOST::ErrorInSolver;
	if( reason >= 0 )
		return true;
	return false;
}
int SolverIterationNumberPetsc(void * data)
{
	PetscInt its;
	KSP * ksp = static_cast<KSP *>(data);
	PetscErrorCode ierr = KSPGetIterationNumber(*ksp,&its);
	if( ierr != PETSC_SUCCESS ) throw INMOST::ErrorInSolver;
	return its;
}
double SolverResidualNormPetsc(void * data)
{
	PetscReal norm;
	KSP * ksp = static_cast<KSP *>(data);
	PetscErrorCode ierr = KSPGetResidualNorm(*ksp,&norm);
	if( ierr != PETSC_SUCCESS ) throw INMOST::ErrorInSolver;
	return norm;
}
Kirill Terekhov's avatar
Kirill Terekhov committed
332

333 334 335 336 337 338 339 340 341 342 343 344 345 346 347 348 349 350 351 352 353 354 355 356 357 358 359 360 361 362 363 364 365 366 367 368 369 370 371 372 373 374 375 376 377 378 379 380 381 382 383

void SolverSetTolerancesPetsc(void * data, double rtol, double atol, double divtol, int maxits)
{
	KSP * ksp = static_cast<KSP *>(data);
	PetscErrorCode ierr = KSPSetTolerances(*ksp,rtol,atol,divtol,maxits);
	if( ierr != PETSC_SUCCESS ) throw INMOST::ErrorInSolver;
}


void SolverSetOverlapPetsc(void * data, int levels)
{
	KSP * ksp = static_cast<KSP *>(data);
	PetscErrorCode ierr;
	PC pc; PCType pc_type;
	ierr = KSPGetPC(*ksp,&pc);
	if( ierr != PETSC_SUCCESS ) throw INMOST::ErrorInSolver;
	ierr = PCGetType(pc,&pc_type);
	if( ierr != PETSC_SUCCESS ) throw INMOST::ErrorInSolver;
	if( !strcmp(pc_type,"asm") )
	{
		ierr = PCASMSetOverlap(pc,levels);
	}
	else if( !strcmp(pc_type,"gasm") )
	{
		ierr = PCGASMSetOverlap(pc,levels);
	}
	if( ierr != PETSC_SUCCESS ) throw INMOST::ErrorInSolver;
}

void SolverSetDropTolerancePetsc(void * data, double dtol)
{
	KSP * ksp = static_cast<KSP *>(data);
	PetscErrorCode ierr;
	PC pc;
	ierr = KSPGetPC(*ksp,&pc);
	if( ierr != PETSC_SUCCESS ) throw INMOST::ErrorInSolver;
	ierr = PCFactorSetDropTolerance(pc,dtol,0.01,10000); //2 other parameters are set to extreme
	if( ierr != PETSC_SUCCESS ) throw INMOST::ErrorInSolver;
}

void SolverSetFillLevelPetsc(void * data, double lfill)
{
	KSP * ksp = static_cast<KSP *>(data);
	PetscErrorCode ierr;
	PC pc;
	ierr = KSPGetPC(*ksp,&pc);
	if( ierr != PETSC_SUCCESS ) throw INMOST::ErrorInSolver;
	ierr = PCFactorSetLevels(pc,static_cast<PetscInt>(lfill));
	if( ierr != PETSC_SUCCESS ) throw INMOST::ErrorInSolver;
}

Kirill Terekhov's avatar
Kirill Terekhov committed
384 385
const char * SolverConvergedReasonPetsc(void * data)
{
Kirill Terekhov's avatar
Kirill Terekhov committed
386
	static char reason_str[4096];
Kirill Terekhov's avatar
Kirill Terekhov committed
387 388 389 390 391 392 393 394
	KSPConvergedReason reason;
	KSP * ksp = static_cast<KSP *>(data);
	PetscErrorCode ierr = KSPGetConvergedReason(*ksp,&reason);
	if( ierr != PETSC_SUCCESS ) throw INMOST::ErrorInSolver;
	switch(reason)
	{
	case KSP_CONVERGED_RTOL:
	case KSP_CONVERGED_RTOL_NORMAL:
Kirill Terekhov's avatar
Kirill Terekhov committed
395
		strcpy(reason_str,"norm decreased by a factor of relative tolerance");
Kirill Terekhov's avatar
Kirill Terekhov committed
396 397 398
		break;
	case KSP_CONVERGED_ATOL:
	case KSP_CONVERGED_ATOL_NORMAL:
Kirill Terekhov's avatar
Kirill Terekhov committed
399
		strcpy(reason_str,"norm less then absolute tolerance");
Kirill Terekhov's avatar
Kirill Terekhov committed
400 401
		break;
	case  KSP_CONVERGED_ITS:
Kirill Terekhov's avatar
Kirill Terekhov committed
402
		strcpy(reason_str,"converged by direct solver");
Kirill Terekhov's avatar
Kirill Terekhov committed
403 404 405
		break;
	case KSP_CONVERGED_CG_NEG_CURVE:
	case KSP_CONVERGED_CG_CONSTRAINED:
Kirill Terekhov's avatar
Kirill Terekhov committed
406
		strcpy(reason_str,"converged due to some condition in conjugate gradient method");
Kirill Terekhov's avatar
Kirill Terekhov committed
407 408
		break;
	case KSP_CONVERGED_STEP_LENGTH:
Kirill Terekhov's avatar
Kirill Terekhov committed
409
		strcpy(reason_str,"converged due to step length");
Kirill Terekhov's avatar
Kirill Terekhov committed
410 411
		break;
	case KSP_CONVERGED_HAPPY_BREAKDOWN :
Kirill Terekhov's avatar
Kirill Terekhov committed
412
		strcpy(reason_str,"converged due to happy breakdown");
Kirill Terekhov's avatar
Kirill Terekhov committed
413 414
		break;
	case  KSP_CONVERGED_ITERATING:
Kirill Terekhov's avatar
Kirill Terekhov committed
415
		strcpy(reason_str,"converged");
Kirill Terekhov's avatar
Kirill Terekhov committed
416 417
		break;
	case KSP_DIVERGED_NULL:
Kirill Terekhov's avatar
Kirill Terekhov committed
418
		strcpy(reason_str,"diverged due to null");
Kirill Terekhov's avatar
Kirill Terekhov committed
419 420
		break;
	case KSP_DIVERGED_ITS:
Kirill Terekhov's avatar
Kirill Terekhov committed
421
		strcpy(reason_str,"diverged due to maximum iteration number");
Kirill Terekhov's avatar
Kirill Terekhov committed
422 423
		break;
	case KSP_DIVERGED_DTOL:
Kirill Terekhov's avatar
Kirill Terekhov committed
424
		strcpy(reason_str,"diverged as divergence tolerance was reached");
Kirill Terekhov's avatar
Kirill Terekhov committed
425 426
		break;
	case KSP_DIVERGED_BREAKDOWN:
Kirill Terekhov's avatar
Kirill Terekhov committed
427
		strcpy(reason_str,"diverged due to breakdown in method");
Kirill Terekhov's avatar
Kirill Terekhov committed
428 429
		break;
	case KSP_DIVERGED_BREAKDOWN_BICG:
Kirill Terekhov's avatar
Kirill Terekhov committed
430
		strcpy(reason_str,"diverged due to breakdown in Biconjugate Gradients method");
Kirill Terekhov's avatar
Kirill Terekhov committed
431 432
		break;
	case KSP_DIVERGED_NONSYMMETRIC:
Kirill Terekhov's avatar
Kirill Terekhov committed
433
		strcpy(reason_str,"diverged since matrix is nonsymmetric");
Kirill Terekhov's avatar
Kirill Terekhov committed
434 435
		break;
	case KSP_DIVERGED_INDEFINITE_PC:
Kirill Terekhov's avatar
Kirill Terekhov committed
436
		strcpy(reason_str,"diverged since preconditioner is not defined");
Kirill Terekhov's avatar
Kirill Terekhov committed
437 438
		break;
	case KSP_DIVERGED_NANORINF:
Kirill Terekhov's avatar
Kirill Terekhov committed
439
		strcpy(reason_str,"diverged since not a number or infinite was encountered");
Kirill Terekhov's avatar
Kirill Terekhov committed
440 441
		break;
	case KSP_DIVERGED_INDEFINITE_MAT:
Kirill Terekhov's avatar
Kirill Terekhov committed
442
		strcpy(reason_str,"diverged due to matrix is not definite");
Kirill Terekhov's avatar
Kirill Terekhov committed
443 444
		break;
	default:
Kirill Terekhov's avatar
Kirill Terekhov committed
445
		strcpy(reason_str,"reason code was not defined in manual");
Kirill Terekhov's avatar
Kirill Terekhov committed
446 447
		break;
	}
Kirill Terekhov's avatar
Kirill Terekhov committed
448
	return reason_str;
Kirill Terekhov's avatar
Kirill Terekhov committed
449
}
Kirill Terekhov's avatar
Kirill Terekhov committed
450 451
#endif
#endif