libigl_Martin/include/igl/shapeup.h

// This file is part of libigl, a simple c++ geometry processing library.
//
// Copyright (C) 2017 Amir Vaxman <avaxman@gmail.com>
//
// This Source Code Form is subject to the terms of the Mozilla Public License
// v. 2.0. If a copy of the MPL was not distributed with this file, You can
// obtain one at http://mozilla.org/MPL/2.0/.
#ifndef IGL_SHAPEUP_H
#define IGL_SHAPEUP_H

#include "shapeup_local_projections.h"
#include <igl/igl_inline.h>
#include <igl/setdiff.h>
#include <igl/cat.h>
#include <Eigen/Core>
#include <vector>


//This file implements the following algorithm:

//Boaziz et al.
//Shape-Up: Shaping Discrete Geometry with Projections
//Computer Graphics Forum (Proc. SGP) 31(5), 2012

namespace igl
{
    
    struct ShapeupData{

        //input data
        Eigen::MatrixXd P;
        Eigen::VectorXi SC;
        Eigen::MatrixXi S;
        Eigen::VectorXi b;
        //  maxIterations     referring to number of local-global pairs.
        int maxIterations;
     //  pTolerance        algorithm stops when max(|P_k-P_{k-1}|)<pTolerance.
        double pTolerance;
        double shapeCoeff, closeCoeff, smoothCoeff;
        std::function<bool(const Eigen::MatrixXd&, const Eigen::VectorXi&, const Eigen::MatrixXi&,  Eigen::MatrixXd&)> local_projection;
        
        //Internally-used matrices
        Eigen::SparseMatrix<double> DShape, DClose, DSmooth, Q, A, At;
        
         min_quad_with_fixed_data<double> solver_data;
        
        ShapeupData():
        maxIterations(50),
        pTolerance(10e-6),
        shapeCoeff(1.0),
        closeCoeff(100.0),
        smoothCoeff(0.0),
        local_projection(igl::shapeup_identity_projection),
        {}
    };
    
    
    
    //This function precomputation the necessary matrices for the ShapeUp process, and prefactorizes them.
    
    //input:
    //  P   #P by 3             point positions
    //  SC  #Set by 1           cardinalities of sets in S
    //  S   #Sets by max(SC)    independent sets where the local projection applies. Values beyond column SC(i) in row S(i,:) are "don't care"
    //  E   #E by 2             the "edges" of a mesh; used for the smoothness energy.
    //  b   #b by 1             boundary (fixed) vertices from P.
    //  w   #Set by 1           weight for each set (used in the global step)
    //  local_projection              function pointer taking (P,SC,S,projP),
    // where the first three parameters are as defined, and "projP" is the output, as a #S by 3*max(SC) function in format xyzxyzxyz, and where it returns the projected points corresponding to each set in S in the same order.
    
    // Output:
    //  sudata struct ShapeupData     the data necessary to solve the system in shapeup_solve

    template <
    typename DerivedP,
    typename DerivedSC,
    typename DerivedS,
    typename Derivedb,
    typename Derivedw>
    IGL_INLINE void shapeup_precomputation(
                                           const Eigen::PlainObjectBase<DerivedP>& P,
                                           const Eigen::PlainObjectBase<DerivedSC>& SC,
                                           const Eigen::PlainObjectBase<DerivedS>& S,
                                           const Eigen::PlainObjectBase<DerivedS>& E,
                                           const Eigen::PlainObjectBase<Derivedb>& b,
                                           const Eigen::PlainObjectBase<Derivedw>& w,
                                           const std::function<bool(const Eigen::PlainObjectBase<DerivedP>&, const Eigen::PlainObjectBase<DerivedSX>&, const Eigen::PlainObjectBase<DerivedS>&,  Eigen::PlainObjectBase<Derivedb&)>& local_projection,
                                           struct ShapeupData& sudata);
    
    //This function solve the shapeup project optimization. shapeup_precompute must be called before with the same sudata, or results are unpredictable
    
    //Input:
    //bc                #b by 3 fixed point values corresonding to "b" in sudata
    //NOTE: the input values in P0 don't need to correspond to prescribed values in bc; the iterations will project them automatically (by design).
    //P0                #P by 3 initial solution (point positions)

    //Output:
    //P                 the solution to the problem, corresponding to P0.
    template <
    typename Derivedbc,
    typename DerivedP>
    IGL_INLINE void shapeup_solve(const Eigen::PlainObjectBase<Derivedbc>& bc,
                                   const Eigen::PlainObjectBase<DerivedP>& P0,
                                    const struct ShapeupData& sudata,
                                    Eigen::PlainObjectBase<DerivedP>& P);
}

#ifndef IGL_STATIC_LIBRARY
#include "shapeup.cpp"
#endif

#endif