Actual source code: itcreate.c
1: /*
2: The basic KSP routines, Create, View etc. are here.
3: */
4: #include src/ksp/ksp/kspimpl.h
5: #include petscsys.h
7: /* Logging support */
8: PetscCookie KSP_COOKIE = 0;
9: PetscEvent KSP_GMRESOrthogonalization = 0, KSP_SetUp = 0, KSP_Solve = 0;
12: PetscTruth KSPRegisterAllCalled = PETSC_FALSE;
16: /*@C
17: KSPView - Prints the KSP data structure.
19: Collective on KSP
21: Input Parameters:
22: + ksp - the Krylov space context
23: - viewer - visualization context
25: Options Database Keys:
26: . -ksp_view - print the ksp data structure at the end of a KSPSolve call
28: Note:
29: The available visualization contexts include
30: + PETSC_VIEWER_STDOUT_SELF - standard output (default)
31: - PETSC_VIEWER_STDOUT_WORLD - synchronized standard
32: output where only the first processor opens
33: the file. All other processors send their
34: data to the first processor to print.
36: The user can open an alternative visualization context with
37: PetscViewerASCIIOpen() - output to a specified file.
39: Level: beginner
41: .keywords: KSP, view
43: .seealso: PCView(), PetscViewerASCIIOpen()
44: @*/
45: PetscErrorCode KSPView(KSP ksp,PetscViewer viewer)
46: {
47: char *type;
49: PetscTruth iascii;
53: if (!viewer) viewer = PETSC_VIEWER_STDOUT_(ksp->comm);
57: PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&iascii);
58: if (iascii) {
59: KSPGetType(ksp,&type);
60: if (ksp->prefix) {
61: PetscViewerASCIIPrintf(viewer,"KSP Object:(%s)\n",ksp->prefix);
62: } else {
63: PetscViewerASCIIPrintf(viewer,"KSP Object:\n");
64: }
65: if (type) {
66: PetscViewerASCIIPrintf(viewer," type: %s\n",type);
67: } else {
68: PetscViewerASCIIPrintf(viewer," type: not yet set\n");
69: }
70: if (ksp->ops->view) {
71: PetscViewerASCIIPushTab(viewer);
72: (*ksp->ops->view)(ksp,viewer);
73: PetscViewerASCIIPopTab(viewer);
74: }
75: if (ksp->guess_zero) {PetscViewerASCIIPrintf(viewer," maximum iterations=%D, initial guess is zero\n",ksp->max_it);}
76: else {PetscViewerASCIIPrintf(viewer," maximum iterations=%D\n", ksp->max_it);}
77: if (ksp->guess_knoll) {PetscViewerASCIIPrintf(viewer," using preconditioner applied to right hand side for initial guess\n");}
78: PetscViewerASCIIPrintf(viewer," tolerances: relative=%g, absolute=%g, divergence=%g\n",ksp->rtol,ksp->abstol,ksp->divtol);
79: if (ksp->pc_side == PC_RIGHT) {PetscViewerASCIIPrintf(viewer," right preconditioning\n");}
80: else if (ksp->pc_side == PC_SYMMETRIC) {PetscViewerASCIIPrintf(viewer," symmetric preconditioning\n");}
81: else {PetscViewerASCIIPrintf(viewer," left preconditioning\n");}
82: } else {
83: if (ksp->ops->view) {
84: (*ksp->ops->view)(ksp,viewer);
85: }
86: }
87: PCView(ksp->pc,viewer);
88: return(0);
89: }
91: /*
92: Contains the list of registered KSP routines
93: */
94: PetscFList KSPList = 0;
98: /*@C
99: KSPSetNormType - Sets the norm that is used for convergence testing.
101: Collective on KSP
103: Input Parameter:
104: + ksp - Krylov solver context
105: - normtype - one of
106: $ KSP_NO_NORM - skips computing the norm, this should only be used if you are using
107: $ the Krylov method as a smoother with a fixed small number of iterations.
108: $ You must also call KSPSetConvergenceTest(ksp,KSPSkipConverged,PETSC_NULL);
109: $ supported only by CG, Richardson, Bi-CG-stab, CR, and CGS methods.
110: $ KSP_PRECONDITIONED_NORM - the default for left preconditioned solves, uses the l2 norm
111: $ of the preconditioned residual
112: $ KSP_UNPRECONDITIONED_NORM - uses the l2 norm of the true b - Ax residual, supported only by
113: $ CG, CHEBYCHEV, and RICHARDSON
114: $ KSP_NATURAL_NORM - supported by cg, cr, and cgs
117: Options Database Key:
118: . -ksp_norm_type <none,preconditioned,unpreconditioned,natural>
120: Notes:
121: Currently only works with the CG, Richardson, Bi-CG-stab, CR, and CGS methods.
123: Level: advanced
125: .keywords: KSP, create, context, norms
127: .seealso: KSPSetUp(), KSPSolve(), KSPDestroy(), KSPSkipConverged()
128: @*/
129: PetscErrorCode KSPSetNormType(KSP ksp,KSPNormType normtype)
130: {
133: ksp->normtype = normtype;
134: if (normtype == KSP_NO_NORM) {
135: PetscLogInfo(ksp,"KSPSetNormType:Warning seting KSPNormType to skip computing the norm\n\
136: make sure you set the KSP convergence test to KSPSkipConvergence\n");
137: }
138: return(0);
139: }
143: static PetscErrorCode KSPPublish_Petsc(PetscObject obj)
144: {
146: return(0);
147: }
151: /*@
152: KSPSetOperators - Sets the matrix associated with the linear system
153: and a (possibly) different one associated with the preconditioner.
155: Collective on KSP and Mat
157: Input Parameters:
158: + ksp - the KSP context
159: . Amat - the matrix associated with the linear system
160: . Pmat - the matrix to be used in constructing the preconditioner, usually the
161: same as Amat.
162: - flag - flag indicating information about the preconditioner matrix structure
163: during successive linear solves. This flag is ignored the first time a
164: linear system is solved, and thus is irrelevant when solving just one linear
165: system.
167: Notes:
168: The flag can be used to eliminate unnecessary work in the preconditioner
169: during the repeated solution of linear systems of the same size. The
170: available options are
171: $ SAME_PRECONDITIONER -
172: $ Pmat is identical during successive linear solves.
173: $ This option is intended for folks who are using
174: $ different Amat and Pmat matrices and want to reuse the
175: $ same preconditioner matrix. For example, this option
176: $ saves work by not recomputing incomplete factorization
177: $ for ILU/ICC preconditioners.
178: $ SAME_NONZERO_PATTERN -
179: $ Pmat has the same nonzero structure during
180: $ successive linear solves.
181: $ DIFFERENT_NONZERO_PATTERN -
182: $ Pmat does not have the same nonzero structure.
184: Caution:
185: If you specify SAME_NONZERO_PATTERN, PETSc believes your assertion
186: and does not check the structure of the matrix. If you erroneously
187: claim that the structure is the same when it actually is not, the new
188: preconditioner will not function correctly. Thus, use this optimization
189: feature carefully!
191: If in doubt about whether your preconditioner matrix has changed
192: structure or not, use the flag DIFFERENT_NONZERO_PATTERN.
194: Level: beginner
196: .keywords: KSP, set, operators, matrix, preconditioner, linear system
198: .seealso: KSPSolve(), KSPGetPC(), PCGetOperators(), PCSetOperators(), KSPGetOperators()
199: @*/
200: PetscErrorCode KSPSetOperators(KSP ksp,Mat Amat,Mat Pmat,MatStructure flag)
201: {
208: PCSetOperators(ksp->pc,Amat,Pmat,flag);
209: if (ksp->setupcalled > 1) ksp->setupcalled = 1; /* so that next solve call will call setup */
210: return(0);
211: }
215: /*@
216: KSPGetOperators - Gets the matrix associated with the linear system
217: and a (possibly) different one associated with the preconditioner.
219: Collective on KSP and Mat
221: Input Parameter:
222: . ksp - the KSP context
224: Output Parameters:
225: + Amat - the matrix associated with the linear system
226: . Pmat - the matrix to be used in constructing the preconditioner, usually the
227: same as Amat.
228: - flag - flag indicating information about the preconditioner matrix structure
229: during successive linear solves. This flag is ignored the first time a
230: linear system is solved, and thus is irrelevant when solving just one linear
231: system.
233: Level: intermediate
235: .keywords: KSP, set, get, operators, matrix, preconditioner, linear system
237: .seealso: KSPSolve(), KSPGetPC(), PCGetOperators(), PCSetOperators(), KSPSetOperators()
238: @*/
239: PetscErrorCode KSPGetOperators(KSP ksp,Mat *Amat,Mat *Pmat,MatStructure *flag)
240: {
245: PCGetOperators(ksp->pc,Amat,Pmat,flag);
246: return(0);
247: }
251: /*@C
252: KSPCreate - Creates the default KSP context.
254: Collective on MPI_Comm
256: Input Parameter:
257: . comm - MPI communicator
259: Output Parameter:
260: . ksp - location to put the KSP context
262: Notes:
263: The default KSP type is GMRES with a restart of 30, using modified Gram-Schmidt
264: orthogonalization.
266: Level: beginner
268: .keywords: KSP, create, context
270: .seealso: KSPSetUp(), KSPSolve(), KSPDestroy(), KSP
271: @*/
272: PetscErrorCode KSPCreate(MPI_Comm comm,KSP *inksp)
273: {
274: KSP ksp;
279: *inksp = 0;
280: #ifndef PETSC_USE_DYNAMIC_LIBRARIES
281: KSPInitializePackage(PETSC_NULL);
282: #endif
284: PetscHeaderCreate(ksp,_p_KSP,struct _KSPOps,KSP_COOKIE,-1,"KSP",comm,KSPDestroy,KSPView);
285: PetscLogObjectCreate(ksp);
286: *inksp = ksp;
287: ksp->bops->publish = KSPPublish_Petsc;
289: ksp->type = -1;
290: ksp->max_it = 10000;
291: ksp->pc_side = PC_LEFT;
292: ksp->rtol = 1.e-5;
293: ksp->abstol = 1.e-50;
294: ksp->divtol = 1.e4;
296: ksp->normtype = KSP_PRECONDITIONED_NORM;
297: ksp->rnorm = 0.0;
298: ksp->its = 0;
299: ksp->guess_zero = PETSC_TRUE;
300: ksp->calc_sings = PETSC_FALSE;
301: ksp->res_hist = PETSC_NULL;
302: ksp->res_hist_len = 0;
303: ksp->res_hist_max = 0;
304: ksp->res_hist_reset = PETSC_TRUE;
305: ksp->numbermonitors = 0;
306: ksp->converged = KSPDefaultConverged;
307: ksp->ops->buildsolution = KSPDefaultBuildSolution;
308: ksp->ops->buildresidual = KSPDefaultBuildResidual;
310: ksp->ops->setfromoptions = 0;
312: ksp->vec_sol = 0;
313: ksp->vec_rhs = 0;
314: ksp->pc = 0;
316: ksp->ops->solve = 0;
317: ksp->ops->setup = 0;
318: ksp->ops->destroy = 0;
320: ksp->data = 0;
321: ksp->nwork = 0;
322: ksp->work = 0;
324: ksp->cnvP = 0;
326: ksp->reason = KSP_CONVERGED_ITERATING;
328: ksp->setupcalled = 0;
329: PetscPublishAll(ksp);
330: PCCreate(comm,&ksp->pc);
331: return(0);
332: }
333:
336: /*@C
337: KSPSetType - Builds KSP for a particular solver.
339: Collective on KSP
341: Input Parameters:
342: + ksp - the Krylov space context
343: - type - a known method
345: Options Database Key:
346: . -ksp_type <method> - Sets the method; use -help for a list
347: of available methods (for instance, cg or gmres)
349: Notes:
350: See "petsc/include/petscksp.h" for available methods (for instance,
351: KSPCG or KSPGMRES).
353: Normally, it is best to use the KSPSetFromOptions() command and
354: then set the KSP type from the options database rather than by using
355: this routine. Using the options database provides the user with
356: maximum flexibility in evaluating the many different Krylov methods.
357: The KSPSetType() routine is provided for those situations where it
358: is necessary to set the iterative solver independently of the command
359: line or options database. This might be the case, for example, when
360: the choice of iterative solver changes during the execution of the
361: program, and the user's application is taking responsibility for
362: choosing the appropriate method. In other words, this routine is
363: not for beginners.
365: Level: intermediate
367: .keywords: KSP, set, method
369: .seealso: PCSetType(), KSPType
371: @*/
372: PetscErrorCode KSPSetType(KSP ksp,const KSPType type)
373: {
374: PetscErrorCode ierr,(*r)(KSP);
375: PetscTruth match;
381: PetscTypeCompare((PetscObject)ksp,type,&match);
382: if (match) return(0);
384: if (ksp->data) {
385: /* destroy the old private KSP context */
386: (*ksp->ops->destroy)(ksp);
387: ksp->data = 0;
388: }
389: /* Get the function pointers for the iterative method requested */
390: if (!KSPRegisterAllCalled) {KSPRegisterAll(PETSC_NULL);}
391: PetscFListFind(ksp->comm,KSPList,type,(void (**)(void)) &r);
392: if (!r) SETERRQ1(PETSC_ERR_ARG_UNKNOWN_TYPE,"Unknown KSP type given: %s",type);
393: ksp->setupcalled = 0;
394: (*r)(ksp);
395: PetscObjectChangeTypeName((PetscObject)ksp,type);
396: return(0);
397: }
401: /*@C
402: KSPRegisterDestroy - Frees the list of KSP methods that were
403: registered by KSPRegisterDynamic().
405: Not Collective
407: Level: advanced
409: .keywords: KSP, register, destroy
411: .seealso: KSPRegisterDynamic(), KSPRegisterAll()
412: @*/
413: PetscErrorCode KSPRegisterDestroy(void)
414: {
418: if (KSPList) {
419: PetscFListDestroy(&KSPList);
420: KSPList = 0;
421: }
422: KSPRegisterAllCalled = PETSC_FALSE;
423: return(0);
424: }
428: /*@C
429: KSPGetType - Gets the KSP type as a string from the KSP object.
431: Not Collective
433: Input Parameter:
434: . ksp - Krylov context
436: Output Parameter:
437: . name - name of KSP method
439: Level: intermediate
441: .keywords: KSP, get, method, name
443: .seealso: KSPSetType()
444: @*/
445: PetscErrorCode KSPGetType(KSP ksp,KSPType *type)
446: {
450: *type = ksp->type_name;
451: return(0);
452: }
456: /*@
457: KSPSetFromOptions - Sets KSP options from the options database.
458: This routine must be called before KSPSetUp() if the user is to be
459: allowed to set the Krylov type.
461: Collective on KSP
463: Input Parameters:
464: . ksp - the Krylov space context
466: Options Database Keys:
467: + -ksp_max_it - maximum number of linear iterations
468: . -ksp_rtol rtol - relative tolerance used in default determination of convergence, i.e.
469: if residual norm decreases by this factor than convergence is declared
470: . -ksp_atol abstol - absolute tolerance used in default convergence test, i.e. if residual
471: norm is less than this then convergence is declared
472: . -ksp_divtol tol - if residual norm increases by this factor than divergence is declared
473: . -ksp_norm_type - none - skip norms used in convergence tests (useful only when not using
474: $ convergence test (say you always want to run with 5 iterations) to
475: $ save on communication overhead
476: $ preconditioned - default for left preconditioning
477: $ unpreconditioned - see KSPSetNormType()
478: $ natural - see KSPSetNormType()
479: . -ksp_constant_null_space - assume the operator (matrix) has the constant vector in its null space
480: . -ksp_test_null_space - tests the null space set with KSPSetNullSpace() to see if it truly is a null space
481: . -ksp_knoll - compute initial guess by applying the preconditioner to the right hand side
482: . -ksp_cancelmonitors - cancel all previous convergene monitor routines set
483: . -ksp_monitor - print residual norm at each iteration
484: . -ksp_xmonitor - plot residual norm at each iteration
485: . -ksp_vecmonitor - plot solution at each iteration
486: - -ksp_singmonitor - monitor extremem singular values at each iteration
488: Notes:
489: To see all options, run your program with the -help option
490: or consult the users manual.
492: Level: beginner
494: .keywords: KSP, set, from, options, database
496: .seealso:
497: @*/
498: PetscErrorCode KSPSetFromOptions(KSP ksp)
499: {
501: PetscInt indx;
502: char type[256];
503: const char *stype[] = {"none","preconditioned","unpreconditioned","natural"};
504: PetscTruth flg;
508: PCSetFromOptions(ksp->pc);
510: if (!KSPRegisterAllCalled) {KSPRegisterAll(PETSC_NULL);}
511: PetscOptionsBegin(ksp->comm,ksp->prefix,"Krylov Method (KSP) Options","KSP");
512: PetscOptionsList("-ksp_type","Krylov method","KSPSetType",KSPList,(char*)(ksp->type_name?ksp->type_name:KSPGMRES),type,256,&flg);
513: if (flg) {
514: KSPSetType(ksp,type);
515: }
516: /*
517: Set the type if it was never set.
518: */
519: if (!ksp->type_name) {
520: KSPSetType(ksp,KSPGMRES);
521: }
523: PetscOptionsInt("-ksp_max_it","Maximum number of iterations","KSPSetTolerances",ksp->max_it,&ksp->max_it,PETSC_NULL);
524: PetscOptionsReal("-ksp_rtol","Relative decrease in residual norm","KSPSetTolerances",ksp->rtol,&ksp->rtol,PETSC_NULL);
525: PetscOptionsReal("-ksp_atol","Absolute value of residual norm","KSPSetTolerances",ksp->abstol,&ksp->abstol,PETSC_NULL);
526: PetscOptionsReal("-ksp_divtol","Residual norm increase cause divergence","KSPSetTolerances",ksp->divtol,&ksp->divtol,PETSC_NULL);
527: PetscOptionsLogical("-ksp_knoll","Use preconditioner applied to b for initial guess","KSPSetInitialGuessKnoll",ksp->guess_knoll,
528: &ksp->guess_knoll,PETSC_NULL);
530: PetscOptionsEList("-ksp_norm_type","KSP Norm type","KSPSetNormType",stype,4,"preconditioned",&indx,&flg);
531: if (flg) {
532: switch (indx) {
533: case 0:
534: KSPSetNormType(ksp,KSP_NO_NORM);
535: KSPSetConvergenceTest(ksp,KSPSkipConverged,0);
536: break;
537: case 1:
538: KSPSetNormType(ksp,KSP_PRECONDITIONED_NORM);
539: break;
540: case 2:
541: KSPSetNormType(ksp,KSP_UNPRECONDITIONED_NORM);
542: break;
543: case 3:
544: KSPSetNormType(ksp,KSP_NATURAL_NORM);
545: break;
546: }
547: }
549: PetscOptionsName("-ksp_diagonal_scale","Diagonal scale matrix before building preconditioner","KSPSetDiagonalScale",&flg);
550: if (flg) {
551: KSPSetDiagonalScale(ksp,PETSC_TRUE);
552: }
553: PetscOptionsName("-ksp_diagonal_scale_fix","Fix diagonaled scaled matrix after solve","KSPSetDiagonalScaleFix",&flg);
554: if (flg) {
555: KSPSetDiagonalScaleFix(ksp,PETSC_TRUE);
556: }
559: PetscOptionsName("-ksp_constant_null_space","Add constant null space to Krylov solver","KSPSetNullSpace",&flg);
560: if (flg) {
561: MatNullSpace nsp;
563: MatNullSpaceCreate(ksp->comm,PETSC_TRUE,0,0,&nsp);
564: KSPSetNullSpace(ksp,nsp);
565: MatNullSpaceDestroy(nsp);
566: }
568: /* option is actually checked in KSPSetUp() */
569: if (ksp->nullsp) {
570: PetscOptionsName("-ksp_test_null_space","Is provided null space correct","None",&flg);
571: }
573: /*
574: Prints reason for convergence or divergence of each linear solve
575: */
576: PetscOptionsName("-ksp_converged_reason","Print reason for converged or diverged","KSPSolve",&flg);
577: if (flg) {
578: ksp->printreason = PETSC_TRUE;
579: }
581: PetscOptionsName("-ksp_cancelmonitors","Remove any hardwired monitor routines","KSPClearMonitor",&flg);
582: /* -----------------------------------------------------------------------*/
583: /*
584: Cancels all monitors hardwired into code before call to KSPSetFromOptions()
585: */
586: if (flg) {
587: KSPClearMonitor(ksp);
588: }
589: /*
590: Prints preconditioned residual norm at each iteration
591: */
592: PetscOptionsName("-ksp_monitor","Monitor preconditioned residual norm","KSPSetMonitor",&flg);
593: if (flg) {
594: KSPSetMonitor(ksp,KSPDefaultMonitor,PETSC_NULL,PETSC_NULL);
595: }
596: /*
597: Plots the vector solution
598: */
599: PetscOptionsName("-ksp_vecmonitor","Monitor solution graphically","KSPSetMonitor",&flg);
600: if (flg) {
601: KSPSetMonitor(ksp,KSPVecViewMonitor,PETSC_NULL,PETSC_NULL);
602: }
603: /*
604: Prints preconditioned and true residual norm at each iteration
605: */
606: PetscOptionsName("-ksp_truemonitor","Monitor true (unpreconditioned) residual norm","KSPSetMonitor",&flg);
607: if (flg) {
608: KSPSetMonitor(ksp,KSPTrueMonitor,PETSC_NULL,PETSC_NULL);
609: }
610: /*
611: Prints extreme eigenvalue estimates at each iteration
612: */
613: PetscOptionsName("-ksp_singmonitor","Monitor singular values","KSPSetMonitor",&flg);
614: if (flg) {
615: KSPSetComputeSingularValues(ksp,PETSC_TRUE);
616: KSPSetMonitor(ksp,KSPSingularValueMonitor,PETSC_NULL,PETSC_NULL);
617: }
618: /*
619: Prints preconditioned residual norm with fewer digits
620: */
621: PetscOptionsName("-ksp_smonitor","Monitor preconditioned residual norm with fewer digitis","KSPSetMonitor",&flg);
622: if (flg) {
623: KSPSetMonitor(ksp,KSPDefaultSMonitor,PETSC_NULL,PETSC_NULL);
624: }
625: /*
626: Graphically plots preconditioned residual norm
627: */
628: PetscOptionsName("-ksp_xmonitor","Monitor graphically preconditioned residual norm","KSPSetMonitor",&flg);
629: if (flg) {
630: KSPSetMonitor(ksp,KSPLGMonitor,PETSC_NULL,PETSC_NULL);
631: }
632: /*
633: Graphically plots preconditioned and true residual norm
634: */
635: PetscOptionsName("-ksp_xtruemonitor","Monitor graphically true residual norm","KSPSetMonitor",&flg);
636: if (flg){
637: KSPSetMonitor(ksp,KSPLGTrueMonitor,PETSC_NULL,PETSC_NULL);
638: }
640: /* -----------------------------------------------------------------------*/
642: PetscOptionsLogicalGroupBegin("-ksp_left_pc","Use left preconditioning","KSPSetPreconditionerSide",&flg);
643: if (flg) { KSPSetPreconditionerSide(ksp,PC_LEFT); }
644: PetscOptionsLogicalGroup("-ksp_right_pc","Use right preconditioning","KSPSetPreconditionerSide",&flg);
645: if (flg) { KSPSetPreconditionerSide(ksp,PC_RIGHT);}
646: PetscOptionsLogicalGroupEnd("-ksp_symmetric_pc","Use symmetric (factorized) preconditioning","KSPSetPreconditionerSide",&flg);
647: if (flg) { KSPSetPreconditionerSide(ksp,PC_SYMMETRIC);}
649: PetscOptionsName("-ksp_compute_singularvalues","Compute singular values of preconditioned operator","KSPSetComputeSingularValues",&flg);
650: if (flg) { KSPSetComputeSingularValues(ksp,PETSC_TRUE); }
651: PetscOptionsName("-ksp_compute_eigenvalues","Compute eigenvalues of preconditioned operator","KSPSetComputeSingularValues",&flg);
652: if (flg) { KSPSetComputeSingularValues(ksp,PETSC_TRUE); }
653: PetscOptionsName("-ksp_plot_eigenvalues","Scatter plot extreme eigenvalues","KSPSetComputeSingularValues",&flg);
654: if (flg) { KSPSetComputeSingularValues(ksp,PETSC_TRUE); }
656: if (ksp->ops->setfromoptions) {
657: (*ksp->ops->setfromoptions)(ksp);
658: }
659: PetscOptionsName("-ksp_view","View linear solver parameters","KSPView",&flg);
660: PetscOptionsEnd();
661: return(0);
662: }
666: /*@C
667: KSPRegister - See KSPRegisterDynamic()
669: Level: advanced
670: @*/
671: PetscErrorCode KSPRegister(const char sname[],const char path[],const char name[],PetscErrorCode (*function)(KSP))
672: {
674: char fullname[PETSC_MAX_PATH_LEN];
677: PetscFListConcat(path,name,fullname);
678: PetscFListAdd(&KSPList,sname,fullname,(void (*)(void))function);
679: return(0);
680: }
684: /*@C
685: KSPSetNullSpace - Sets the null space of the operator
687: Collective on KSP
689: Input Parameters:
690: + ksp - the Krylov space object
691: - nullsp - the null space of the operator
693: Level: advanced
695: .seealso: KSPSetOperators(), MatNullSpaceCreate(), KSPGetNullSpace()
696: @*/
697: PetscErrorCode KSPSetNullSpace(KSP ksp,MatNullSpace nullsp)
698: {
702: if (ksp->nullsp) {
703: MatNullSpaceDestroy(ksp->nullsp);
704: }
705: ksp->nullsp = nullsp;
706: PetscObjectReference((PetscObject)ksp->nullsp);
707: return(0);
708: }
712: /*@C
713: KSPGetNullSpace - Gets the null space of the operator
715: Collective on KSP
717: Input Parameters:
718: + ksp - the Krylov space object
719: - nullsp - the null space of the operator
721: Level: advanced
723: .seealso: KSPSetOperators(), MatNullSpaceCreate(), KSPSetNullSpace()
724: @*/
725: PetscErrorCode KSPGetNullSpace(KSP ksp,MatNullSpace *nullsp)
726: {
728: *nullsp = ksp->nullsp;
729: return(0);
730: }