mbdyn_siconos.h File Reference

This graph shows which files directly or indirectly include this file:

Go to the source code of this file.

Classes
struct	solver_parameters

Enumerations
enum	LCPsolver {   QP, CPG, PGS, RPGS,   PSOR, NSQP, LATIN, LATIN_W,   LEXICO_LEMKE, NEWTON_MIN, NEWTON_FB }

Functions
void	mbdyn_siconos_LCP_call (int size, double M[], double blcp[], double zlem[], double wlem[], solver_parameters &solparam)

Enumeration Type Documentation

enum LCPsolver

Enumerator
QP
CPG
PGS
RPGS
PSOR
NSQP
LATIN
LATIN_W
LEXICO_LEMKE
NEWTON_MIN
NEWTON_FB

Definition at line 40 of file mbdyn_siconos.h.

               {        
        // quadratic programming formulation
        QP,

        // CPG (Conjugated Projected Gradient) solver for LCP based on quadratic minimization.
        CPG,

        // PGS is a basic Projected Gauss-Seidel solver for LCP.
        PGS,

        // Regularized Projected Gauss-Seidel, is a solver for LCP, 
        // able to handle matrices with null diagonal terms
        RPGS,

        // Projected Succesive over relaxation solver for LCP. See cottle, Pang Stone Chap 5
        PSOR,

        // quadratic programm formulation for solving an non symmetric LCP
        NSQP,

        // (LArge Time INcrements) is a basic latin solver for LCP.
        LATIN,

        // (LArge Time INcrements) is a basic latin solver with relaxation for LCP
        LATIN_W,

        // direct solver for LCP based on pivoting method principle for degenerate problem.
        // Choice of pivot variable is performed via lexicographic ordering 
        LEXICO_LEMKE,

        // nonsmooth Newton method based on the min formulation (or max formulation) of the LCP
        NEWTON_MIN,

        // uses a nonsmooth newton method based on the Fischer-Bursmeister convex function
        NEWTON_FB

#if 0
        // Gauss-Seidel solver based on a Sparse-Block storage for the matrix M of the LCP.
        // Can't be used here because Matrix M of the LCP must be formulated as SparseBlockStructuredMatrix. 
        NSGS_SBM
#endif
};

Function Documentation

void mbdyn_siconos_LCP_call	(	int	size,
		double	M[],
		double	blcp[],
		double	zlem[],
		double	wlem[],
		solver_parameters &	solparam
	)

Definition at line 43 of file mbdyn_siconos.cc.

References CPG, solver_parameters::info, LATIN, LATIN_W, LEXICO_LEMKE, NEWTON_FB, NEWTON_MIN, NSQP, PGS, solver_parameters::processed_iterations, PSOR, QP, solver_parameters::resulting_error, RPGS, solver_parameters::solver, solver_parameters::solveritermax, and solver_parameters::solvertol.

Referenced by ModuleNonsmoothNode::mbs_get_force(), and ModuleNonsmoothNode::mbs_get_force_frictional().

{
        solparam.info = 0;

        LCPsolver lcpsol = solparam.solver;
        double tolerance = solparam.solvertol;
        int maxiternum = solparam.solveritermax;

        // LCP problem description:
        LinearComplementarityProblem OSNSProblem;
        OSNSProblem.size = sizep;
        OSNSProblem.q = bLCP;

        NumericsMatrix MM;
        MM.storageType = 0;
        MM.matrix0 = W_NN;              // W_delassus
        MM.size0 = sizep;
        MM.size1 = sizep;
        OSNSProblem.M = &MM;

        SolverOptions numerics_solver_options = { 0 };

        solparam.resulting_error = 0.;

        // Solving LCP problem
        // FIXME: initialize once, keep alive?
        switch (lcpsol) {
        case LEXICO_LEMKE:
                linearComplementarity_lexicolemke_setDefaultSolverOptions(&numerics_solver_options);
                numerics_solver_options.iparam[0] = maxiternum;
                lcp_lexicolemke(&OSNSProblem, Pkp1, wlem, &solparam.info, &numerics_solver_options);
                solparam.processed_iterations = numerics_solver_options.iparam[1];
                break;

        case RPGS:
                linearComplementarity_rpgs_setDefaultSolverOptions(&numerics_solver_options);
                numerics_solver_options.iparam[0] = maxiternum;
                numerics_solver_options.dparam[0] = tolerance;
                lcp_rpgs(&OSNSProblem, Pkp1, wlem, &solparam.info, &numerics_solver_options);
                solparam.processed_iterations = numerics_solver_options.iparam[1];
                solparam.resulting_error = numerics_solver_options.dparam[1];
                break;

        case QP:
                linearComplementarity_qp_setDefaultSolverOptions(&numerics_solver_options);
                numerics_solver_options.dparam[0] = tolerance;
                lcp_qp(&OSNSProblem, Pkp1, wlem, &solparam.info, &numerics_solver_options);
                break;

        case CPG:
                linearComplementarity_cpg_setDefaultSolverOptions(&numerics_solver_options);
                numerics_solver_options.iparam[0] = maxiternum;
                numerics_solver_options.dparam[0] = tolerance;
                lcp_cpg(&OSNSProblem, Pkp1, wlem, &solparam.info, &numerics_solver_options);
                solparam.processed_iterations = numerics_solver_options.iparam[1];
                solparam.resulting_error = numerics_solver_options.dparam[1];
                break;

        case PGS:
                linearComplementarity_pgs_setDefaultSolverOptions(&numerics_solver_options);
                numerics_solver_options.iparam[0] = maxiternum;
                numerics_solver_options.dparam[0] = tolerance;
                lcp_pgs(&OSNSProblem, Pkp1, wlem, &solparam.info, &numerics_solver_options);
                solparam.processed_iterations = numerics_solver_options.iparam[1];
                solparam.resulting_error = numerics_solver_options.dparam[1];
                break;

        case PSOR:
                linearComplementarity_psor_setDefaultSolverOptions(&numerics_solver_options);
                numerics_solver_options.iparam[0] = maxiternum;
                numerics_solver_options.dparam[0] = tolerance;
                lcp_psor(&OSNSProblem, Pkp1, wlem, &solparam.info, &numerics_solver_options);
                solparam.processed_iterations = numerics_solver_options.iparam[1];
                solparam.resulting_error = numerics_solver_options.dparam[1];
                break;

        case NSQP:
                linearComplementarity_nsqp_setDefaultSolverOptions(&numerics_solver_options);
                numerics_solver_options.dparam[0] = tolerance;
                lcp_nsqp(&OSNSProblem, Pkp1, wlem, &solparam.info, &numerics_solver_options);
                break;

        case LATIN:
                linearComplementarity_latin_setDefaultSolverOptions(&numerics_solver_options);
                numerics_solver_options.iparam[0] = maxiternum;
                numerics_solver_options.dparam[0] = tolerance;
                lcp_latin(&OSNSProblem, Pkp1, wlem, &solparam.info, &numerics_solver_options);
                break;

        case LATIN_W:
                linearComplementarity_latin_w_setDefaultSolverOptions(&numerics_solver_options);
                numerics_solver_options.iparam[0] = maxiternum;
                numerics_solver_options.dparam[0] = tolerance;
                lcp_latin_w(&OSNSProblem, Pkp1, wlem, &solparam.info, &numerics_solver_options);
                break;

        case NEWTON_MIN:
                linearComplementarity_newton_min_setDefaultSolverOptions(&numerics_solver_options);
                numerics_solver_options.iparam[0] = maxiternum;
                numerics_solver_options.dparam[0] = tolerance;
                lcp_newton_min(&OSNSProblem, Pkp1, wlem, &solparam.info, &numerics_solver_options);
                solparam.processed_iterations = numerics_solver_options.iparam[1];
                solparam.resulting_error = numerics_solver_options.dparam[1];
                break;

        case NEWTON_FB:
                linearComplementarity_newton_FB_setDefaultSolverOptions(&numerics_solver_options);
                numerics_solver_options.iparam[0] = maxiternum;
                numerics_solver_options.dparam[0] = tolerance;
                lcp_newton_FB(&OSNSProblem, Pkp1, wlem, &solparam.info, &numerics_solver_options);
                solparam.processed_iterations = numerics_solver_options.iparam[1];
                solparam.resulting_error = numerics_solver_options.dparam[1];
                break;
        }

#if 0
        // analyzing the LCP solver final status
        switch (solparam.info) {
        case 0: // convergence
                break;

        case 1: // iter=itermax
                std::cout << std::endl
                        << "loadable element nonsmooth node: max iterations reached in LCP solver "
                        << std::endl;
                std::cout << "processed_iterations "<< solparam.processed_iterations
                        << " resulting_error " << solparam.resulting_error <<std::endl;
                break;

        default: // other problem
                std::cout << std::endl
                        <<"loadable element nonsmooth node: problem in solution of LCP"
                        << std::endl;
                std::cout << "processed_iterations "<< solparam.processed_iterations
                        << " resulting_error " << solparam.resulting_error <<std::endl;
                break;
        }
#endif

        deleteSolverOptions(&numerics_solver_options);
}