SchurComplementSolver.h

Go to the documentation of this file.

/* This file is part of MoDELib, the Mechanics Of Defects Evolution Library.
 *
 *
 * MoDELib is distributed without any warranty under the
 * GNU General Public License (GPL) v2 <http://www.gnu.org/licenses/>.
 */
 
#ifndef model_SchurComplementSolver_H_
#define model_SchurComplementSolver_H_
 
#include <Eigen/Dense>
 
namespace model
{

    class SchurComplementSolver
    {
 
        /* constructor ********************************************************/
        void checkDecomposition() const;
 
      public:
        const Eigen::LLT<Eigen::MatrixXd> llt;

        const Eigen::PartialPivLU<Eigen::MatrixXd> pplu;

        Eigen::VectorXd L;

        Eigen::VectorXd X;
 
        /* constructor with Fp ************************************************/
        SchurComplementSolver(const Eigen::MatrixXd &KQQ, const Eigen::MatrixXd &KPQ, const Eigen::VectorXd &Fq,
                              const Eigen::VectorXd &Fp);
 
        /* constructor without Fp *********************************************/
        SchurComplementSolver(const Eigen::MatrixXd &KQQ, const Eigen::MatrixXd &KPQ, const Eigen::VectorXd &Fq);
 
        //        /* constructor without Fp *********************************************/
        //              SchurComplementSolver(const Eigen::MatrixXd& KQQ,
        //        /*                 */ const Eigen::MatrixXd& KPQ) :
        //        /* init list                                     */ llt(KQQ),
        //              /* init list                                     */ pplu(KPQ*llt.solve(KPQ.transpose())){
        //            checkDecomposition();
        //              }
 
        //        /* constructor without Fp *********************************************/
        //              void solve(const Eigen::VectorXd& Fq) {
        //            L=pplu.solve(KPQ*llt.solve(Fq));
        //            X=llt.solve(Fq-KPQ.transpose()*L);
        //              }
        //
        //        /* constructor without Fp *********************************************/
        //        void solve(const Eigen::VectorXd& KPQ, const Eigen::VectorXd& Fq, const Eigen::VectorXd& Fp) {
        //                L=pplu.solve(KPQ*llt.solve(Fq)-Fp);
        //                X=llt.solve(Fq-KPQ.transpose()*L);
        //        }
    };
 
} // namespace model
#endif