:orphan:
# KSPFETIDP
The FETI-DP method [1] 
## Options Database Keys

- ***-ksp_fetidp_fullyredundant <false>   -*** use a fully redundant set of Lagrange multipliers
- ***-ksp_fetidp_saddlepoint <false>      -*** activates support for saddle point problems, see [2]
- ***-ksp_fetidp_saddlepoint_flip <false> -*** usually, an incompressible Stokes problem is written as
| A B^T | | v | = | f |
| B 0   | | p | = | g |
with B representing -\int_\Omega \nabla \cdot u q.
If -ksp_fetidp_saddlepoint_flip is true, the code assumes that the user provides it as
| A B^T | | v | = | f |
|-B 0   | | p | = |-g |
- ***-ksp_fetidp_pressure_field <-*** 1>      - activates support for saddle point problems, and identifies the pressure field id.
If this information is not provided, the pressure field is detected by using MatFindZeroDiagonals().
- ***-ksp_fetidp_pressure_all <false>     -*** if false, uses the interface pressures, as described in [2]. If true, uses the entire pressure field.





## Notes
The matrix for the `KSP` must be of type `MATIS`.

The FETI-DP linear system (automatically generated constructing an internal `PCBDDC` object) is solved using an internal `KSP` object.

Options for the inner `KSP` and for the customization of the `PCBDDC` object can be specified at command line by using the prefixes `-fetidp_` and `-fetidp_bddc_`. E.g.,
```none
      -fetidp_ksp_type gmres -fetidp_bddc_pc_bddc_symmetric false
```

will use `KSPGMRES` for the solution of the linear system on the Lagrange multipliers, generated using a non-symmetric `PCBDDC`.

For saddle point problems with continuous pressures, the preconditioned operator for the pressure solver can be specified with `KSPFETIDPSetPressureOperator()`.
Alternatively, the pressure operator is extracted from the precondioned matrix (if it is different from the linear solver matrix).
If none of the above, an identity matrix will be created; the user then needs to scale it through a Richardson solver.
Options for the pressure solver can be prefixed with `-fetidp_fielsplit_p_`, E.g.
```none
      -fetidp_fielsplit_p_ksp_type preonly -fetidp_fielsplit_p_pc_type lu -fetidp_fielsplit_p_pc_factor_mat_solver_type mumps
```

In order to use the deluxe version of FETI-DP, you must customize the inner `PCBDDC` operator with -fetidp_bddc_pc_bddc_use_deluxe_scaling -fetidp_bddc_pc_bddc_deluxe_singlemat and use
non-redundant multipliers, i.e. `-ksp_fetidp_fullyredundant false`. Options for the scaling solver are prefixed by `-fetidp_bddelta_`, E.g.
```none
      -fetidp_bddelta_pc_factor_mat_solver_type mumps -fetidp_bddelta_pc_type lu
```


Some of the basic options such as the maximum number of iterations and tolerances are automatically passed from this `KSP` to the inner `KSP` that actually performs the iterations.

The converged reason and number of iterations computed are passed from the inner `KSP` to this `KSP` at the end of the solution.


## Developer Note
Even though this method does not directly use any norms, the user is allowed to set the `KSPNormType` to any value.
This is so users do not have to change `KSPNormType` options when they switch from other `KSP` methods to this one.


## References

- ***[1] -*** C. Farhat, M. Lesoinne, P. LeTallec, K. Pierson, and D. Rixen, FETI-DP: a dual-primal unified FETI method. I. A faster alternative to the two-level FETI method, Internat. J. Numer. Methods Engrg., 50 (2001), pp. 1523--1544
- ***[2] -*** X. Tu, J. Li, A FETI-DP type domain decomposition algorithm for three-dimensional incompressible Stokes equations, SIAM J. Numer. Anal., 53 (2015), pp. 720-742



## See Also
 [](ch_ksp), `MATIS`, `PCBDDC`, `KSPFETIDPSetInnerBDDC()`, `KSPFETIDPGetInnerBDDC()`, `KSPFETIDPGetInnerKSP()`

## Level
Advanced

## Location
<A HREF="PETSC_DOC_OUT_ROOT_PLACEHOLDER/src/ksp/ksp/impls/fetidp/fetidp.c.html#KSPFETIDP">src/ksp/ksp/impls/fetidp/fetidp.c</A>

## Examples
<A HREF="PETSC_DOC_OUT_ROOT_PLACEHOLDER/src/ksp/ksp/tutorials/ex59.c.html">src/ksp/ksp/tutorials/ex59.c</A><BR>


---
[Edit on GitLab](https://gitlab.com/petsc/petsc/-/edit/release/src/ksp/ksp/impls/fetidp/fetidp.c)


[Index of all KSP routines](index.md)  
[Table of Contents for all manual pages](/manualpages/index.md)  
[Index of all manual pages](/manualpages/singleindex.md)  
