[Sketch] placecgs: remove unused includes
- also sort includes - also fix too long lines etc. (done by clang formatter)
This commit is contained in:
@@ -21,11 +21,13 @@
|
||||
***************************************************************************/
|
||||
|
||||
#include <algorithm>
|
||||
#include <cmath>
|
||||
#define DEBUG_DERIVS 0
|
||||
#if DEBUG_DERIVS
|
||||
# include <cassert>
|
||||
#endif
|
||||
#include <cmath>
|
||||
|
||||
#include <boost/graph/graph_concepts.hpp>
|
||||
|
||||
#include "Constraints.h"
|
||||
|
||||
|
||||
@@ -24,8 +24,6 @@
|
||||
#define PLANEGCS_CONSTRAINTS_H
|
||||
|
||||
#include "Geo.h"
|
||||
#include "Util.h"
|
||||
#include <boost/graph/graph_concepts.hpp>
|
||||
|
||||
//#define _GCS_EXTRACT_SOLVER_SUBSYSTEM_ // This enables debugging code intended to extract information to file bug reports against Eigen, not for production code
|
||||
|
||||
@@ -35,6 +33,7 @@
|
||||
#define _PROTECTED_UNLESS_EXTRACT_MODE_ protected
|
||||
#endif
|
||||
|
||||
|
||||
namespace GCS
|
||||
{
|
||||
|
||||
|
||||
File diff suppressed because it is too large
Load Diff
@@ -23,12 +23,11 @@
|
||||
#ifndef PLANEGCS_GCS_H
|
||||
#define PLANEGCS_GCS_H
|
||||
|
||||
#include "SubSystem.h"
|
||||
#include <boost/concept_check.hpp>
|
||||
#include <boost/graph/graph_concepts.hpp>
|
||||
|
||||
#include <Eigen/QR>
|
||||
|
||||
#include "SubSystem.h"
|
||||
|
||||
|
||||
#define EIGEN_VERSION (EIGEN_WORLD_VERSION * 10000 \
|
||||
+ EIGEN_MAJOR_VERSION * 100 \
|
||||
+ EIGEN_MINOR_VERSION)
|
||||
@@ -50,12 +49,14 @@ namespace GCS
|
||||
// Solver
|
||||
///////////////////////////////////////
|
||||
|
||||
enum SolveStatus {
|
||||
Success = 0, // Found a solution zeroing the error function
|
||||
Converged = 1, // Found a solution minimizing the error function
|
||||
Failed = 2, // Failed to find any solution
|
||||
SuccessfulSolutionInvalid = 3, // This is a solution where the solver succeeded, but the resulting geometry is OCE-invalid
|
||||
};
|
||||
enum SolveStatus
|
||||
{
|
||||
Success = 0, // Found a solution zeroing the error function
|
||||
Converged = 1, // Found a solution minimizing the error function
|
||||
Failed = 2, // Failed to find any solution
|
||||
SuccessfulSolutionInvalid = 3,// This is a solution where the solver succeeded, but the
|
||||
// resulting geometry is OCE-invalid
|
||||
};
|
||||
|
||||
enum Algorithm {
|
||||
BFGS = 0,
|
||||
@@ -81,15 +82,19 @@ namespace GCS
|
||||
};
|
||||
|
||||
// Magic numbers for Constraint tags
|
||||
// - Positive Tags identify a higher level constraint form which the solver constraint originates
|
||||
// - Positive Tags identify a higher level constraint form which the solver constraint
|
||||
// originates
|
||||
// - Negative Tags represent temporary constraints, used for example in moving operations, these
|
||||
// have a different handling in component splitting, see GCS::initSolution. Lifetime is defined by
|
||||
// the container object via GCS::clearByTag.
|
||||
// - -1 is typically used as tag for these temporary constraints, its parameters are enforced with
|
||||
// a lower priority than the main system (real sketcher constraints). It gives a nice effect when
|
||||
// dragging the edge of an unconstrained circle, that the center won't move if the edge can be dragged,
|
||||
// and only when/if the edge cannot be dragged, e.g. radius constraint, the center is moved).
|
||||
enum SpecialTag {
|
||||
// have a different handling in component splitting, see GCS::initSolution. Lifetime is defined
|
||||
// by the container object via GCS::clearByTag.
|
||||
// - -1 is typically used as tag for these temporary constraints, its parameters are
|
||||
// enforced with
|
||||
// a lower priority than the main system (real sketcher constraints). It gives a nice
|
||||
// effect when dragging the edge of an unconstrained circle, that the center won't move
|
||||
// if the edge can be dragged, and only when/if the edge cannot be dragged, e.g. radius
|
||||
// constraint, the center is moved).
|
||||
enum SpecialTag
|
||||
{
|
||||
DefaultTemporaryConstraint = -1
|
||||
};
|
||||
|
||||
@@ -119,9 +124,10 @@ namespace GCS
|
||||
void setReference(); // copies the current parameter values to reference
|
||||
void resetToReference(); // reverts all parameter values to the stored reference
|
||||
|
||||
std::vector< VEC_pD > plists; // partitioned plist except equality constraints
|
||||
std::vector< std::vector<Constraint *> > clists; // partitioned clist except equality constraints
|
||||
std::vector< MAP_pD_pD > reductionmaps; // for simplification of equality constraints
|
||||
std::vector<VEC_pD> plists;// partitioned plist except equality constraints
|
||||
std::vector<std::vector<Constraint*>>
|
||||
clists; // partitioned clist except equality constraints
|
||||
std::vector<MAP_pD_pD> reductionmaps;// for simplification of equality constraints
|
||||
|
||||
int dofs;
|
||||
std::set<Constraint *> redundant;
|
||||
@@ -137,59 +143,51 @@ namespace GCS
|
||||
int solve_LM(SubSystem *subsys, bool isRedundantsolving=false);
|
||||
int solve_DL(SubSystem *subsys, bool isRedundantsolving=false);
|
||||
|
||||
void makeReducedJacobian(Eigen::MatrixXd &J, std::map<int,int> &jacobianconstraintmap, GCS::VEC_pD &pdiagnoselist, std::map< int , int> &tagmultiplicity);
|
||||
void makeReducedJacobian(Eigen::MatrixXd& J, std::map<int, int>& jacobianconstraintmap,
|
||||
GCS::VEC_pD& pdiagnoselist, std::map<int, int>& tagmultiplicity);
|
||||
|
||||
void makeDenseQRDecomposition( const Eigen::MatrixXd &J,
|
||||
const std::map<int,int> &jacobianconstraintmap,
|
||||
Eigen::FullPivHouseholderQR<Eigen::MatrixXd>& qrJT,
|
||||
int &rank, Eigen::MatrixXd &R, bool transposeJ = true, bool silent = false);
|
||||
void makeDenseQRDecomposition(const Eigen::MatrixXd& J,
|
||||
const std::map<int, int>& jacobianconstraintmap,
|
||||
Eigen::FullPivHouseholderQR<Eigen::MatrixXd>& qrJT, int& rank,
|
||||
Eigen::MatrixXd& R, bool transposeJ = true,
|
||||
bool silent = false);
|
||||
|
||||
#ifdef EIGEN_SPARSEQR_COMPATIBLE
|
||||
void makeSparseQRDecomposition( const Eigen::MatrixXd &J,
|
||||
const std::map<int,int> &jacobianconstraintmap,
|
||||
Eigen::SparseQR<Eigen::SparseMatrix<double>, Eigen::COLAMDOrdering<int> > &SqrJT,
|
||||
int &rank, Eigen::MatrixXd &R, bool transposeJ = true, bool silent = false);
|
||||
void makeSparseQRDecomposition(
|
||||
const Eigen::MatrixXd& J, const std::map<int, int>& jacobianconstraintmap,
|
||||
Eigen::SparseQR<Eigen::SparseMatrix<double>, Eigen::COLAMDOrdering<int>>& SqrJT,
|
||||
int& rank, Eigen::MatrixXd& R, bool transposeJ = true, bool silent = false);
|
||||
#endif
|
||||
// This function name is long for a reason:
|
||||
// - Only for DenseQR
|
||||
// - Only for Transposed Jacobian QR decomposition
|
||||
void identifyDependentGeometryParametersInTransposedJacobianDenseQRDecomposition(
|
||||
const Eigen::FullPivHouseholderQR<Eigen::MatrixXd>& qrJT,
|
||||
const GCS::VEC_pD &pdiagnoselist,
|
||||
int paramsNum, int rank
|
||||
);
|
||||
const Eigen::FullPivHouseholderQR<Eigen::MatrixXd>& qrJT,
|
||||
const GCS::VEC_pD& pdiagnoselist, int paramsNum, int rank);
|
||||
|
||||
template <typename T>
|
||||
void identifyConflictingRedundantConstraints( Algorithm alg,
|
||||
const T & qrJT,
|
||||
const std::map<int,int> &jacobianconstraintmap,
|
||||
const std::map< int , int> &tagmultiplicity,
|
||||
GCS::VEC_pD &pdiagnoselist,
|
||||
Eigen::MatrixXd &R,
|
||||
int constrNum, int rank,
|
||||
int &nonredundantconstrNum
|
||||
);
|
||||
template<typename T>
|
||||
void identifyConflictingRedundantConstraints(
|
||||
Algorithm alg, const T& qrJT, const std::map<int, int>& jacobianconstraintmap,
|
||||
const std::map<int, int>& tagmultiplicity, GCS::VEC_pD& pdiagnoselist,
|
||||
Eigen::MatrixXd& R, int constrNum, int rank, int& nonredundantconstrNum);
|
||||
|
||||
void eliminateNonZerosOverPivotInUpperTriangularMatrix(Eigen::MatrixXd &R, int rank);
|
||||
void eliminateNonZerosOverPivotInUpperTriangularMatrix(Eigen::MatrixXd& R, int rank);
|
||||
|
||||
#ifdef EIGEN_SPARSEQR_COMPATIBLE
|
||||
void identifyDependentParametersSparseQR( const Eigen::MatrixXd &J,
|
||||
const std::map<int,int> &jacobianconstraintmap,
|
||||
const GCS::VEC_pD &pdiagnoselist,
|
||||
bool silent=true);
|
||||
void identifyDependentParametersSparseQR(const Eigen::MatrixXd& J,
|
||||
const std::map<int, int>& jacobianconstraintmap,
|
||||
const GCS::VEC_pD& pdiagnoselist,
|
||||
bool silent = true);
|
||||
#endif
|
||||
|
||||
void identifyDependentParametersDenseQR( const Eigen::MatrixXd &J,
|
||||
const std::map<int,int> &jacobianconstraintmap,
|
||||
const GCS::VEC_pD &pdiagnoselist,
|
||||
bool silent=true);
|
||||
void identifyDependentParametersDenseQR(const Eigen::MatrixXd& J,
|
||||
const std::map<int, int>& jacobianconstraintmap,
|
||||
const GCS::VEC_pD& pdiagnoselist,
|
||||
bool silent = true);
|
||||
|
||||
template <typename T>
|
||||
void identifyDependentParameters( T & qrJ,
|
||||
Eigen::MatrixXd &Rparams,
|
||||
int rank,
|
||||
const GCS::VEC_pD &pdiagnoselist,
|
||||
bool silent=true);
|
||||
template<typename T>
|
||||
void identifyDependentParameters(T& qrJ, Eigen::MatrixXd& Rparams, int rank,
|
||||
const GCS::VEC_pD& pdiagnoselist, bool silent = true);
|
||||
|
||||
#ifdef _GCS_EXTRACT_SOLVER_SUBSYSTEM_
|
||||
void extractSubsystem(SubSystem *subsys, bool isRedundantsolving);
|
||||
@@ -197,7 +195,8 @@ namespace GCS
|
||||
public:
|
||||
int maxIter;
|
||||
int maxIterRedundant;
|
||||
bool sketchSizeMultiplier; // if true note that the total number of iterations allowed is MaxIterations *xLength
|
||||
bool sketchSizeMultiplier;// if true note that the total number of iterations allowed is
|
||||
// MaxIterations *xLength
|
||||
bool sketchSizeMultiplierRedundant;
|
||||
double convergence;
|
||||
double convergenceRedundant;
|
||||
@@ -230,110 +229,147 @@ namespace GCS
|
||||
void removeConstraint(Constraint *constr);
|
||||
|
||||
// basic constraints
|
||||
int addConstraintEqual(double *param1, double *param2, int tagId=0, bool driving = true, Constraint::Alignment internalalignment = Constraint::Alignment::NoInternalAlignment);
|
||||
int addConstraintProportional(double *param1, double *param2, double ratio, int tagId, bool driving = true);
|
||||
int addConstraintDifference(double *param1, double *param2,
|
||||
double *difference, int tagId=0, bool driving = true);
|
||||
int addConstraintP2PDistance(Point &p1, Point &p2, double *distance, int tagId=0, bool driving = true);
|
||||
int addConstraintP2PAngle(Point &p1, Point &p2, double *angle,
|
||||
double incrAngle, int tagId=0, bool driving = true);
|
||||
int addConstraintP2PAngle(Point &p1, Point &p2, double *angle, int tagId=0, bool driving = true);
|
||||
int addConstraintP2LDistance(Point &p, Line &l, double *distance, int tagId=0, bool driving = true);
|
||||
int addConstraintPointOnLine(Point &p, Line &l, int tagId=0, bool driving = true);
|
||||
int addConstraintPointOnLine(Point &p, Point &lp1, Point &lp2, int tagId=0, bool driving = true);
|
||||
int addConstraintPointOnPerpBisector(Point &p, Line &l, int tagId=0, bool driving = true);
|
||||
int addConstraintPointOnPerpBisector(Point &p, Point &lp1, Point &lp2, int tagId=0, bool driving = true);
|
||||
int addConstraintParallel(Line &l1, Line &l2, int tagId=0, bool driving = true);
|
||||
int addConstraintPerpendicular(Line &l1, Line &l2, int tagId=0, bool driving = true);
|
||||
int addConstraintPerpendicular(Point &l1p1, Point &l1p2,
|
||||
Point &l2p1, Point &l2p2, int tagId=0, bool driving = true);
|
||||
int addConstraintL2LAngle(Line &l1, Line &l2, double *angle, int tagId=0, bool driving = true);
|
||||
int addConstraintL2LAngle(Point &l1p1, Point &l1p2, Point &l2p1, Point &l2p2,
|
||||
double *angle, int tagId=0, bool driving = true);
|
||||
int addConstraintAngleViaPoint(Curve &crv1, Curve &crv2, Point &p,
|
||||
double *angle, int tagId=0, bool driving = true);
|
||||
int addConstraintMidpointOnLine(Line &l1, Line &l2, int tagId=0, bool driving = true);
|
||||
int addConstraintMidpointOnLine(Point &l1p1, Point &l1p2, Point &l2p1, Point &l2p2,
|
||||
int tagId=0, bool driving = true);
|
||||
int addConstraintTangentCircumf(Point &p1, Point &p2, double *rd1, double *rd2,
|
||||
bool internal=false, int tagId=0, bool driving = true);
|
||||
int addConstraintTangentAtBSplineKnot(BSpline &b, Line &l, unsigned int knotindex,
|
||||
int tagId=0, bool driving = true);
|
||||
int addConstraintEqual(
|
||||
double* param1, double* param2, int tagId = 0, bool driving = true,
|
||||
Constraint::Alignment internalalignment = Constraint::Alignment::NoInternalAlignment);
|
||||
int addConstraintProportional(double* param1, double* param2, double ratio, int tagId,
|
||||
bool driving = true);
|
||||
int addConstraintDifference(double* param1, double* param2, double* difference,
|
||||
int tagId = 0, bool driving = true);
|
||||
int addConstraintP2PDistance(Point& p1, Point& p2, double* distance, int tagId = 0,
|
||||
bool driving = true);
|
||||
int addConstraintP2PAngle(Point& p1, Point& p2, double* angle, double incrAngle,
|
||||
int tagId = 0, bool driving = true);
|
||||
int addConstraintP2PAngle(Point& p1, Point& p2, double* angle, int tagId = 0,
|
||||
bool driving = true);
|
||||
int addConstraintP2LDistance(Point& p, Line& l, double* distance, int tagId = 0,
|
||||
bool driving = true);
|
||||
int addConstraintPointOnLine(Point& p, Line& l, int tagId = 0, bool driving = true);
|
||||
int addConstraintPointOnLine(Point& p, Point& lp1, Point& lp2, int tagId = 0,
|
||||
bool driving = true);
|
||||
int addConstraintPointOnPerpBisector(Point& p, Line& l, int tagId = 0, bool driving = true);
|
||||
int addConstraintPointOnPerpBisector(Point& p, Point& lp1, Point& lp2, int tagId = 0,
|
||||
bool driving = true);
|
||||
int addConstraintParallel(Line& l1, Line& l2, int tagId = 0, bool driving = true);
|
||||
int addConstraintPerpendicular(Line& l1, Line& l2, int tagId = 0, bool driving = true);
|
||||
int addConstraintPerpendicular(Point& l1p1, Point& l1p2, Point& l2p1, Point& l2p2,
|
||||
int tagId = 0, bool driving = true);
|
||||
int addConstraintL2LAngle(Line& l1, Line& l2, double* angle, int tagId = 0,
|
||||
bool driving = true);
|
||||
int addConstraintL2LAngle(Point& l1p1, Point& l1p2, Point& l2p1, Point& l2p2, double* angle,
|
||||
int tagId = 0, bool driving = true);
|
||||
int addConstraintAngleViaPoint(Curve& crv1, Curve& crv2, Point& p, double* angle,
|
||||
int tagId = 0, bool driving = true);
|
||||
int addConstraintMidpointOnLine(Line& l1, Line& l2, int tagId = 0, bool driving = true);
|
||||
int addConstraintMidpointOnLine(Point& l1p1, Point& l1p2, Point& l2p1, Point& l2p2,
|
||||
int tagId = 0, bool driving = true);
|
||||
int addConstraintTangentCircumf(Point& p1, Point& p2, double* rd1, double* rd2,
|
||||
bool internal = false, int tagId = 0, bool driving = true);
|
||||
int addConstraintTangentAtBSplineKnot(BSpline& b, Line& l, unsigned int knotindex,
|
||||
int tagId = 0, bool driving = true);
|
||||
|
||||
// derived constraints
|
||||
int addConstraintP2PCoincident(Point &p1, Point &p2, int tagId=0, bool driving = true);
|
||||
int addConstraintHorizontal(Line &l, int tagId=0, bool driving = true);
|
||||
int addConstraintHorizontal(Point &p1, Point &p2, int tagId=0, bool driving = true);
|
||||
int addConstraintVertical(Line &l, int tagId=0, bool driving = true);
|
||||
int addConstraintVertical(Point &p1, Point &p2, int tagId=0, bool driving = true);
|
||||
int addConstraintCoordinateX(Point &p, double *x, int tagId=0, bool driving = true);
|
||||
int addConstraintCoordinateY(Point &p, double *y, int tagId=0, bool driving = true);
|
||||
int addConstraintArcRules(Arc &a, int tagId=0, bool driving = true);
|
||||
int addConstraintPointOnCircle(Point &p, Circle &c, int tagId=0, bool driving = true);
|
||||
int addConstraintPointOnEllipse(Point &p, Ellipse &e, int tagId=0, bool driving = true);
|
||||
int addConstraintPointOnHyperbolicArc(Point &p, ArcOfHyperbola &e, int tagId=0, bool driving = true);
|
||||
int addConstraintPointOnParabolicArc(Point &p, ArcOfParabola &e, int tagId=0, bool driving = true);
|
||||
int addConstraintPointOnBSpline(Point &p, BSpline &b, double* pointparam, int tagId, bool driving = true);
|
||||
int addConstraintArcOfEllipseRules(ArcOfEllipse &a, int tagId=0, bool driving = true);
|
||||
int addConstraintCurveValue(Point &p, Curve &a, double *u, int tagId=0, bool driving = true);
|
||||
int addConstraintArcOfHyperbolaRules(ArcOfHyperbola &a, int tagId=0, bool driving = true);
|
||||
int addConstraintArcOfParabolaRules(ArcOfParabola &a, int tagId=0, bool driving = true);
|
||||
int addConstraintPointOnArc(Point &p, Arc &a, int tagId=0, bool driving = true);
|
||||
int addConstraintPerpendicularLine2Arc(Point &p1, Point &p2, Arc &a,
|
||||
int tagId=0, bool driving = true);
|
||||
int addConstraintPerpendicularArc2Line(Arc &a, Point &p1, Point &p2,
|
||||
int tagId=0, bool driving = true);
|
||||
int addConstraintPerpendicularCircle2Arc(Point ¢er, double *radius, Arc &a,
|
||||
int tagId=0, bool driving = true);
|
||||
int addConstraintPerpendicularArc2Circle(Arc &a, Point ¢er, double *radius,
|
||||
int tagId=0, bool driving = true);
|
||||
int addConstraintPerpendicularArc2Arc(Arc &a1, bool reverse1,
|
||||
Arc &a2, bool reverse2, int tagId=0, bool driving = true);
|
||||
int addConstraintTangent(Line &l, Circle &c, int tagId=0, bool driving = true);
|
||||
int addConstraintTangent(Line &l, Ellipse &e, int tagId=0, bool driving = true);
|
||||
int addConstraintTangent(Line &l, Arc &a, int tagId=0, bool driving = true);
|
||||
int addConstraintTangent(Circle &c1, Circle &c2, int tagId=0, bool driving = true);
|
||||
int addConstraintTangent(Arc &a1, Arc &a2, int tagId=0, bool driving = true);
|
||||
int addConstraintTangent(Circle &c, Arc &a, int tagId=0, bool driving = true);
|
||||
int addConstraintP2PCoincident(Point& p1, Point& p2, int tagId = 0, bool driving = true);
|
||||
int addConstraintHorizontal(Line& l, int tagId = 0, bool driving = true);
|
||||
int addConstraintHorizontal(Point& p1, Point& p2, int tagId = 0, bool driving = true);
|
||||
int addConstraintVertical(Line& l, int tagId = 0, bool driving = true);
|
||||
int addConstraintVertical(Point& p1, Point& p2, int tagId = 0, bool driving = true);
|
||||
int addConstraintCoordinateX(Point& p, double* x, int tagId = 0, bool driving = true);
|
||||
int addConstraintCoordinateY(Point& p, double* y, int tagId = 0, bool driving = true);
|
||||
int addConstraintArcRules(Arc& a, int tagId = 0, bool driving = true);
|
||||
int addConstraintPointOnCircle(Point& p, Circle& c, int tagId = 0, bool driving = true);
|
||||
int addConstraintPointOnEllipse(Point& p, Ellipse& e, int tagId = 0, bool driving = true);
|
||||
int addConstraintPointOnHyperbolicArc(Point& p, ArcOfHyperbola& e, int tagId = 0,
|
||||
bool driving = true);
|
||||
int addConstraintPointOnParabolicArc(Point& p, ArcOfParabola& e, int tagId = 0,
|
||||
bool driving = true);
|
||||
int addConstraintPointOnBSpline(Point& p, BSpline& b, double* pointparam, int tagId,
|
||||
bool driving = true);
|
||||
int addConstraintArcOfEllipseRules(ArcOfEllipse& a, int tagId = 0, bool driving = true);
|
||||
int addConstraintCurveValue(Point& p, Curve& a, double* u, int tagId = 0,
|
||||
bool driving = true);
|
||||
int addConstraintArcOfHyperbolaRules(ArcOfHyperbola& a, int tagId = 0, bool driving = true);
|
||||
int addConstraintArcOfParabolaRules(ArcOfParabola& a, int tagId = 0, bool driving = true);
|
||||
int addConstraintPointOnArc(Point& p, Arc& a, int tagId = 0, bool driving = true);
|
||||
int addConstraintPerpendicularLine2Arc(Point& p1, Point& p2, Arc& a, int tagId = 0,
|
||||
bool driving = true);
|
||||
int addConstraintPerpendicularArc2Line(Arc& a, Point& p1, Point& p2, int tagId = 0,
|
||||
bool driving = true);
|
||||
int addConstraintPerpendicularCircle2Arc(Point& center, double* radius, Arc& a,
|
||||
int tagId = 0, bool driving = true);
|
||||
int addConstraintPerpendicularArc2Circle(Arc& a, Point& center, double* radius,
|
||||
int tagId = 0, bool driving = true);
|
||||
int addConstraintPerpendicularArc2Arc(Arc& a1, bool reverse1, Arc& a2, bool reverse2,
|
||||
int tagId = 0, bool driving = true);
|
||||
int addConstraintTangent(Line& l, Circle& c, int tagId = 0, bool driving = true);
|
||||
int addConstraintTangent(Line& l, Ellipse& e, int tagId = 0, bool driving = true);
|
||||
int addConstraintTangent(Line& l, Arc& a, int tagId = 0, bool driving = true);
|
||||
int addConstraintTangent(Circle& c1, Circle& c2, int tagId = 0, bool driving = true);
|
||||
int addConstraintTangent(Arc& a1, Arc& a2, int tagId = 0, bool driving = true);
|
||||
int addConstraintTangent(Circle& c, Arc& a, int tagId = 0, bool driving = true);
|
||||
|
||||
int addConstraintCircleRadius(Circle &c, double *radius, int tagId=0, bool driving = true);
|
||||
int addConstraintArcRadius(Arc &a, double *radius, int tagId=0, bool driving = true);
|
||||
int addConstraintCircleDiameter(Circle &c, double *diameter, int tagId=0, bool driving = true);
|
||||
int addConstraintArcDiameter(Arc &a, double *diameter, int tagId=0, bool driving = true);
|
||||
int addConstraintEqualLength(Line &l1, Line &l2, int tagId=0, bool driving = true);
|
||||
int addConstraintEqualRadius(Circle &c1, Circle &c2, int tagId=0, bool driving = true);
|
||||
int addConstraintEqualRadii(Ellipse &e1, Ellipse &e2, int tagId=0, bool driving = true);
|
||||
int addConstraintEqualRadii(ArcOfHyperbola &a1, ArcOfHyperbola &a2, int tagId=0, bool driving = true);
|
||||
int addConstraintEqualRadius(Circle &c1, Arc &a2, int tagId=0, bool driving = true);
|
||||
int addConstraintEqualRadius(Arc &a1, Arc &a2, int tagId=0, bool driving = true);
|
||||
int addConstraintEqualFocus(ArcOfParabola &a1, ArcOfParabola &a2, int tagId=0, bool driving = true);
|
||||
int addConstraintP2PSymmetric(Point &p1, Point &p2, Line &l, int tagId=0, bool driving = true);
|
||||
int addConstraintP2PSymmetric(Point &p1, Point &p2, Point &p, int tagId=0, bool driving = true);
|
||||
int addConstraintSnellsLaw(Curve &ray1, Curve &ray2,
|
||||
Curve &boundary, Point p,
|
||||
double* n1, double* n2,
|
||||
bool flipn1, bool flipn2,
|
||||
int tagId, bool driving = true);
|
||||
int addConstraintCircleRadius(Circle& c, double* radius, int tagId = 0,
|
||||
bool driving = true);
|
||||
int addConstraintArcRadius(Arc& a, double* radius, int tagId = 0, bool driving = true);
|
||||
int addConstraintCircleDiameter(Circle& c, double* diameter, int tagId = 0,
|
||||
bool driving = true);
|
||||
int addConstraintArcDiameter(Arc& a, double* diameter, int tagId = 0, bool driving = true);
|
||||
int addConstraintEqualLength(Line& l1, Line& l2, int tagId = 0, bool driving = true);
|
||||
int addConstraintEqualRadius(Circle& c1, Circle& c2, int tagId = 0, bool driving = true);
|
||||
int addConstraintEqualRadii(Ellipse& e1, Ellipse& e2, int tagId = 0, bool driving = true);
|
||||
int addConstraintEqualRadii(ArcOfHyperbola& a1, ArcOfHyperbola& a2, int tagId = 0,
|
||||
bool driving = true);
|
||||
int addConstraintEqualRadius(Circle& c1, Arc& a2, int tagId = 0, bool driving = true);
|
||||
int addConstraintEqualRadius(Arc& a1, Arc& a2, int tagId = 0, bool driving = true);
|
||||
int addConstraintEqualFocus(ArcOfParabola& a1, ArcOfParabola& a2, int tagId = 0,
|
||||
bool driving = true);
|
||||
int addConstraintP2PSymmetric(Point& p1, Point& p2, Line& l, int tagId = 0,
|
||||
bool driving = true);
|
||||
int addConstraintP2PSymmetric(Point& p1, Point& p2, Point& p, int tagId = 0,
|
||||
bool driving = true);
|
||||
int addConstraintSnellsLaw(Curve& ray1, Curve& ray2, Curve& boundary, Point p, double* n1,
|
||||
double* n2, bool flipn1, bool flipn2, int tagId,
|
||||
bool driving = true);
|
||||
|
||||
int addConstraintC2CDistance(Circle &c1, Circle &c2, double *dist, int tagId, bool driving = true);
|
||||
int addConstraintC2CDistance(Circle& c1, Circle& c2, double* dist, int tagId,
|
||||
bool driving = true);
|
||||
|
||||
// internal alignment constraints
|
||||
int addConstraintInternalAlignmentPoint2Ellipse(Ellipse &e, Point &p1, InternalAlignmentType alignmentType, int tagId=0, bool driving = true);
|
||||
int addConstraintInternalAlignmentEllipseMajorDiameter(Ellipse &e, Point &p1, Point &p2, int tagId=0, bool driving = true);
|
||||
int addConstraintInternalAlignmentEllipseMinorDiameter(Ellipse &e, Point &p1, Point &p2, int tagId=0, bool driving = true);
|
||||
int addConstraintInternalAlignmentEllipseFocus1(Ellipse &e, Point &p1, int tagId=0, bool driving = true);
|
||||
int addConstraintInternalAlignmentEllipseFocus2(Ellipse &e, Point &p1, int tagId=0, bool driving = true);
|
||||
int addConstraintInternalAlignmentPoint2Hyperbola(Hyperbola &e, Point &p1, InternalAlignmentType alignmentType, int tagId=0, bool driving = true);
|
||||
int addConstraintInternalAlignmentHyperbolaMajorDiameter(Hyperbola &e, Point &p1, Point &p2, int tagId=0, bool driving = true);
|
||||
int addConstraintInternalAlignmentHyperbolaMinorDiameter(Hyperbola &e, Point &p1, Point &p2, int tagId=0, bool driving = true);
|
||||
int addConstraintInternalAlignmentHyperbolaFocus(Hyperbola &e, Point &p1, int tagId=0, bool driving = true);
|
||||
int addConstraintInternalAlignmentParabolaFocus(Parabola &e, Point &p1, int tagId=0, bool driving = true);
|
||||
int addConstraintInternalAlignmentBSplineControlPoint(BSpline &b, Circle &c, unsigned int poleindex, int tag=0, bool driving = true);
|
||||
int addConstraintInternalAlignmentKnotPoint(BSpline &b, Point &p, unsigned int knotindex, int tagId=0, bool driving=true);
|
||||
int addConstraintInternalAlignmentPoint2Ellipse(Ellipse& e, Point& p1,
|
||||
InternalAlignmentType alignmentType,
|
||||
int tagId = 0, bool driving = true);
|
||||
int addConstraintInternalAlignmentEllipseMajorDiameter(Ellipse& e, Point& p1, Point& p2,
|
||||
int tagId = 0, bool driving = true);
|
||||
int addConstraintInternalAlignmentEllipseMinorDiameter(Ellipse& e, Point& p1, Point& p2,
|
||||
int tagId = 0, bool driving = true);
|
||||
int addConstraintInternalAlignmentEllipseFocus1(Ellipse& e, Point& p1, int tagId = 0,
|
||||
bool driving = true);
|
||||
int addConstraintInternalAlignmentEllipseFocus2(Ellipse& e, Point& p1, int tagId = 0,
|
||||
bool driving = true);
|
||||
int addConstraintInternalAlignmentPoint2Hyperbola(Hyperbola& e, Point& p1,
|
||||
InternalAlignmentType alignmentType,
|
||||
int tagId = 0, bool driving = true);
|
||||
int addConstraintInternalAlignmentHyperbolaMajorDiameter(Hyperbola& e, Point& p1, Point& p2,
|
||||
int tagId = 0,
|
||||
bool driving = true);
|
||||
int addConstraintInternalAlignmentHyperbolaMinorDiameter(Hyperbola& e, Point& p1, Point& p2,
|
||||
int tagId = 0,
|
||||
bool driving = true);
|
||||
int addConstraintInternalAlignmentHyperbolaFocus(Hyperbola& e, Point& p1, int tagId = 0,
|
||||
bool driving = true);
|
||||
int addConstraintInternalAlignmentParabolaFocus(Parabola& e, Point& p1, int tagId = 0,
|
||||
bool driving = true);
|
||||
int addConstraintInternalAlignmentBSplineControlPoint(BSpline& b, Circle& c,
|
||||
unsigned int poleindex, int tag = 0,
|
||||
bool driving = true);
|
||||
int addConstraintInternalAlignmentKnotPoint(BSpline& b, Point& p, unsigned int knotindex,
|
||||
int tagId = 0, bool driving = true);
|
||||
|
||||
double calculateAngleViaPoint(const Curve &crv1, const Curve &crv2, Point &p) const;
|
||||
double calculateAngleViaPoint(const Curve &crv1, const Curve &crv2, Point &p1, Point &p2) const;
|
||||
void calculateNormalAtPoint(const Curve &crv, const Point &p, double &rtnX, double &rtnY) const;
|
||||
double calculateAngleViaPoint(const Curve& crv1, const Curve& crv2, Point& p) const;
|
||||
double calculateAngleViaPoint(const Curve& crv1, const Curve& crv2, Point& p1,
|
||||
Point& p2) const;
|
||||
void calculateNormalAtPoint(const Curve& crv, const Point& p, double& rtnX,
|
||||
double& rtnY) const;
|
||||
|
||||
// Calculates errors of all constraints which have a tag equal to
|
||||
// the one supplied. Individual errors are summed up using RMS.
|
||||
@@ -348,36 +384,69 @@ namespace GCS
|
||||
void declareDrivenParams(VEC_pD ¶ms);
|
||||
void initSolution(Algorithm alg=DogLeg);
|
||||
|
||||
int solve(bool isFine=true, Algorithm alg=DogLeg, bool isRedundantsolving=false);
|
||||
int solve(VEC_pD ¶ms, bool isFine=true, Algorithm alg=DogLeg, bool isRedundantsolving=false);
|
||||
int solve(SubSystem *subsys, bool isFine=true, Algorithm alg=DogLeg, bool isRedundantsolving=false);
|
||||
int solve(SubSystem *subsysA, SubSystem *subsysB, bool isFine=true, bool isRedundantsolving=false);
|
||||
int solve(bool isFine = true, Algorithm alg = DogLeg, bool isRedundantsolving = false);
|
||||
int solve(VEC_pD& params, bool isFine = true, Algorithm alg = DogLeg,
|
||||
bool isRedundantsolving = false);
|
||||
int solve(SubSystem* subsys, bool isFine = true, Algorithm alg = DogLeg,
|
||||
bool isRedundantsolving = false);
|
||||
int solve(SubSystem* subsysA, SubSystem* subsysB, bool isFine = true,
|
||||
bool isRedundantsolving = false);
|
||||
|
||||
void applySolution();
|
||||
void undoSolution();
|
||||
//FIXME: looks like XconvergenceFine is not the solver precision, at least in DogLeg solver.
|
||||
// Note: Yes, every solver has a different way of interpreting precision
|
||||
// but one has to study what is this needed for in order to decide
|
||||
// what to return (this is unchanged from previous versions)
|
||||
double getFinePrecision(){ return convergence;}
|
||||
// FIXME: looks like XconvergenceFine is not the solver precision, at least in DogLeg
|
||||
// solver.
|
||||
// Note: Yes, every solver has a different way of interpreting precision
|
||||
// but one has to study what is this needed for in order to decide
|
||||
// what to return (this is unchanged from previous versions)
|
||||
double getFinePrecision()
|
||||
{
|
||||
return convergence;
|
||||
}
|
||||
|
||||
int diagnose(Algorithm alg=DogLeg);
|
||||
int dofsNumber() const { return hasDiagnosis ? dofs : -1; }
|
||||
void getConflicting(VEC_I &conflictingOut) const
|
||||
{ conflictingOut = hasDiagnosis ? conflictingTags : VEC_I(0); }
|
||||
void getRedundant(VEC_I &redundantOut) const
|
||||
{ redundantOut = hasDiagnosis ? redundantTags : VEC_I(0); }
|
||||
void getPartiallyRedundant (VEC_I &partiallyredundantOut) const
|
||||
{ partiallyredundantOut = hasDiagnosis ? partiallyRedundantTags : VEC_I(0); }
|
||||
void getDependentParams(VEC_pD &pdependentparameterlist) const
|
||||
{ pdependentparameterlist = pDependentParameters;}
|
||||
void getDependentParamsGroups(std::vector<std::vector<double *>> &pdependentparametergroups) const
|
||||
{ pdependentparametergroups = pDependentParametersGroups;}
|
||||
bool isEmptyDiagnoseMatrix() const {return emptyDiagnoseMatrix;}
|
||||
int diagnose(Algorithm alg = DogLeg);
|
||||
int dofsNumber() const
|
||||
{
|
||||
return hasDiagnosis ? dofs : -1;
|
||||
}
|
||||
void getConflicting(VEC_I& conflictingOut) const
|
||||
{
|
||||
conflictingOut = hasDiagnosis ? conflictingTags : VEC_I(0);
|
||||
}
|
||||
void getRedundant(VEC_I& redundantOut) const
|
||||
{
|
||||
redundantOut = hasDiagnosis ? redundantTags : VEC_I(0);
|
||||
}
|
||||
void getPartiallyRedundant(VEC_I& partiallyredundantOut) const
|
||||
{
|
||||
partiallyredundantOut = hasDiagnosis ? partiallyRedundantTags : VEC_I(0);
|
||||
}
|
||||
void getDependentParams(VEC_pD& pdependentparameterlist) const
|
||||
{
|
||||
pdependentparameterlist = pDependentParameters;
|
||||
}
|
||||
void
|
||||
getDependentParamsGroups(std::vector<std::vector<double*>>& pdependentparametergroups) const
|
||||
{
|
||||
pdependentparametergroups = pDependentParametersGroups;
|
||||
}
|
||||
bool isEmptyDiagnoseMatrix() const
|
||||
{
|
||||
return emptyDiagnoseMatrix;
|
||||
}
|
||||
|
||||
bool hasConflicting() const {return !(hasDiagnosis && conflictingTags.empty());}
|
||||
bool hasRedundant() const {return !(hasDiagnosis && redundantTags.empty());}
|
||||
bool hasPartiallyRedundant() const {return !(hasDiagnosis && partiallyRedundantTags.empty());}
|
||||
bool hasConflicting() const
|
||||
{
|
||||
return !(hasDiagnosis && conflictingTags.empty());
|
||||
}
|
||||
bool hasRedundant() const
|
||||
{
|
||||
return !(hasDiagnosis && redundantTags.empty());
|
||||
}
|
||||
bool hasPartiallyRedundant() const
|
||||
{
|
||||
return !(hasDiagnosis && partiallyRedundantTags.empty());
|
||||
}
|
||||
|
||||
void invalidatedDiagnosis();
|
||||
};
|
||||
|
||||
@@ -24,16 +24,19 @@
|
||||
#if DEBUG_DERIVS
|
||||
#endif
|
||||
|
||||
#include <cassert>
|
||||
|
||||
#include "Geo.h"
|
||||
|
||||
#include <cassert>
|
||||
|
||||
namespace GCS{
|
||||
|
||||
DeriVector2::DeriVector2(const Point &p, const double *derivparam)
|
||||
{
|
||||
x=*p.x; y=*p.y;
|
||||
dx=0.0; dy=0.0;
|
||||
x = *p.x;
|
||||
y = *p.y;
|
||||
dx = 0.0;
|
||||
dy = 0.0;
|
||||
if (derivparam == p.x)
|
||||
dx = 1.0;
|
||||
if (derivparam == p.y)
|
||||
@@ -43,31 +46,33 @@ DeriVector2::DeriVector2(const Point &p, const double *derivparam)
|
||||
double DeriVector2::length(double &dlength) const
|
||||
{
|
||||
double l = length();
|
||||
if(l==0){
|
||||
if (l == 0) {
|
||||
dlength = 1.0;
|
||||
return l;
|
||||
} else {
|
||||
dlength = (x*dx + y*dy)/l;
|
||||
}
|
||||
else {
|
||||
dlength = (x * dx + y * dy) / l;
|
||||
return l;
|
||||
}
|
||||
}
|
||||
|
||||
DeriVector2 DeriVector2::getNormalized() const
|
||||
{
|
||||
double l=length();
|
||||
if(l==0.0) {
|
||||
double l = length();
|
||||
if (l == 0.0) {
|
||||
return DeriVector2(0, 0, dx, dy);
|
||||
} else {
|
||||
}
|
||||
else {
|
||||
DeriVector2 rtn;
|
||||
rtn.x = x/l;
|
||||
rtn.y = y/l;
|
||||
//first, simply scale the derivative accordingly.
|
||||
rtn.dx = dx/l;
|
||||
rtn.dy = dy/l;
|
||||
//next, remove the collinear part of dx,dy (make a projection onto a normal)
|
||||
double dsc = rtn.dx*rtn.x + rtn.dy*rtn.y;//scalar product d*v
|
||||
rtn.dx -= dsc*rtn.x;//subtract the projection
|
||||
rtn.dy -= dsc*rtn.y;
|
||||
rtn.x = x / l;
|
||||
rtn.y = y / l;
|
||||
// first, simply scale the derivative accordingly.
|
||||
rtn.dx = dx / l;
|
||||
rtn.dy = dy / l;
|
||||
// next, remove the collinear part of dx,dy (make a projection onto a normal)
|
||||
double dsc = rtn.dx * rtn.x + rtn.dy * rtn.y;// scalar product d*v
|
||||
rtn.dx -= dsc * rtn.x; // subtract the projection
|
||||
rtn.dy -= dsc * rtn.y;
|
||||
return rtn;
|
||||
}
|
||||
}
|
||||
@@ -75,17 +80,15 @@ DeriVector2 DeriVector2::getNormalized() const
|
||||
double DeriVector2::scalarProd(const DeriVector2 &v2, double *dprd) const
|
||||
{
|
||||
if (dprd) {
|
||||
*dprd = dx*v2.x + x*v2.dx + dy*v2.y + y*v2.dy;
|
||||
*dprd = dx * v2.x + x * v2.dx + dy * v2.y + y * v2.dy;
|
||||
};
|
||||
return x*v2.x + y*v2.y;
|
||||
return x * v2.x + y * v2.y;
|
||||
}
|
||||
|
||||
DeriVector2 DeriVector2::divD(double val, double dval) const
|
||||
{
|
||||
return DeriVector2(x/val,y/val,
|
||||
dx/val - x*dval/(val*val),
|
||||
dy/val - y*dval/(val*val)
|
||||
);
|
||||
return DeriVector2(
|
||||
x / val, y / val, dx / val - x * dval / (val * val), dy / val - y * dval / (val * val));
|
||||
}
|
||||
|
||||
DeriVector2 Curve::Value(double /*u*/, double /*du*/, const double* /*derivparam*/) const
|
||||
@@ -149,17 +152,18 @@ DeriVector2 Circle::CalculateNormal(const Point &p, const double* derivparam) co
|
||||
|
||||
DeriVector2 Circle::Value(double u, double du, const double* derivparam) const
|
||||
{
|
||||
//(x,y) = center + cos(u)*(r,0) + sin(u)*(0,r)
|
||||
|
||||
DeriVector2 cv (center, derivparam);
|
||||
DeriVector2 cv(center, derivparam);
|
||||
double r, dr;
|
||||
r = *(this->rad); dr = (derivparam == this->rad) ? 1.0 : 0.0;
|
||||
DeriVector2 ex (r,0.0,dr,0.0);
|
||||
r = *(this->rad);
|
||||
dr = (derivparam == this->rad) ? 1.0 : 0.0;
|
||||
DeriVector2 ex(r, 0.0, dr, 0.0);
|
||||
DeriVector2 ey = ex.rotate90ccw();
|
||||
double si, dsi, co, dco;
|
||||
si = std::sin(u); dsi = du*std::cos(u);
|
||||
co = std::cos(u); dco = du*(-std::sin(u));
|
||||
return cv.sum(ex.multD(co,dco).sum(ey.multD(si,dsi)));
|
||||
si = std::sin(u);
|
||||
dsi = du * std::cos(u);
|
||||
co = std::cos(u);
|
||||
dco = du * (-std::sin(u));
|
||||
return cv.sum(ex.multD(co, dco).sum(ey.multD(si, dsi)));
|
||||
}
|
||||
|
||||
int Circle::PushOwnParams(VEC_pD &pvec)
|
||||
@@ -214,13 +218,18 @@ Arc* Arc::Copy()
|
||||
|
||||
//--------------ellipse
|
||||
|
||||
//this function is exposed to allow reusing pre-filled derivectors in constraints code
|
||||
double Ellipse::getRadMaj(const DeriVector2 ¢er, const DeriVector2 &f1, double b, double db, double &ret_dRadMaj) const
|
||||
// this function is exposed to allow reusing pre-filled derivectors in constraints code
|
||||
double Ellipse::getRadMaj(const DeriVector2& center, const DeriVector2& f1, double b, double db,
|
||||
double& ret_dRadMaj) const
|
||||
{
|
||||
double cf, dcf;
|
||||
cf = f1.subtr(center).length(dcf);
|
||||
DeriVector2 hack (b, cf,
|
||||
db, dcf);//hack = a nonsense vector to calculate major radius with derivatives, useful just because the calculation formula is the same as vector length formula
|
||||
DeriVector2 hack(
|
||||
b,
|
||||
cf,
|
||||
db,
|
||||
dcf);// hack = a nonsense vector to calculate major radius with derivatives, useful just
|
||||
// because the calculation formula is the same as vector length formula
|
||||
return hack.length(ret_dRadMaj);
|
||||
}
|
||||
|
||||
@@ -370,14 +379,15 @@ ArcOfEllipse* ArcOfEllipse::Copy()
|
||||
|
||||
//---------------hyperbola
|
||||
|
||||
//this function is exposed to allow reusing pre-filled derivectors in constraints code
|
||||
double Hyperbola::getRadMaj(const DeriVector2 ¢er, const DeriVector2 &f1, double b, double db, double &ret_dRadMaj) const
|
||||
// this function is exposed to allow reusing pre-filled derivectors in constraints code
|
||||
double Hyperbola::getRadMaj(const DeriVector2& center, const DeriVector2& f1, double b, double db,
|
||||
double& ret_dRadMaj) const
|
||||
{
|
||||
double cf, dcf;
|
||||
cf = f1.subtr(center).length(dcf);
|
||||
double a, da;
|
||||
a = sqrt(cf*cf - b*b);
|
||||
da = (dcf*cf - db*b)/a;
|
||||
a = sqrt(cf * cf - b * b);
|
||||
da = (dcf * cf - db * b) / a;
|
||||
ret_dRadMaj = da;
|
||||
return a;
|
||||
}
|
||||
@@ -397,21 +407,22 @@ double Hyperbola::getRadMaj() const
|
||||
return getRadMaj(nullptr,dradmaj);
|
||||
}
|
||||
|
||||
DeriVector2 Hyperbola::CalculateNormal(const Point &p, const double* derivparam) const
|
||||
DeriVector2 Hyperbola::CalculateNormal(const Point& p, const double* derivparam) const
|
||||
{
|
||||
//fill some vectors in
|
||||
DeriVector2 cv (center, derivparam);
|
||||
DeriVector2 f1v (focus1, derivparam);
|
||||
DeriVector2 pv (p, derivparam);
|
||||
// fill some vectors in
|
||||
DeriVector2 cv(center, derivparam);
|
||||
DeriVector2 f1v(focus1, derivparam);
|
||||
DeriVector2 pv(p, derivparam);
|
||||
|
||||
//calculation.
|
||||
//focus2:
|
||||
DeriVector2 f2v = cv.linCombi(2.0, f1v, -1.0); // 2*cv - f1v
|
||||
// calculation.
|
||||
// focus2:
|
||||
DeriVector2 f2v = cv.linCombi(2.0, f1v, -1.0);// 2*cv - f1v
|
||||
|
||||
//pf1, pf2 = vectors from p to focus1,focus2
|
||||
DeriVector2 pf1 = f1v.subtr(pv).mult(-1.0); // <--- differs from ellipse normal calculation code by inverting this vector
|
||||
// pf1, pf2 = vectors from p to focus1,focus2
|
||||
DeriVector2 pf1 = f1v.subtr(pv).mult(
|
||||
-1.0);// <--- differs from ellipse normal calculation code by inverting this vector
|
||||
DeriVector2 pf2 = f2v.subtr(pv);
|
||||
//return sum of normalized pf2, pf2
|
||||
// return sum of normalized pf2, pf2
|
||||
DeriVector2 ret = pf1.getNormalized().sum(pf2.getNormalized());
|
||||
|
||||
return ret;
|
||||
@@ -506,17 +517,17 @@ ArcOfHyperbola* ArcOfHyperbola::Copy()
|
||||
|
||||
//---------------parabola
|
||||
|
||||
DeriVector2 Parabola::CalculateNormal(const Point &p, const double* derivparam) const
|
||||
DeriVector2 Parabola::CalculateNormal(const Point& p, const double* derivparam) const
|
||||
{
|
||||
//fill some vectors in
|
||||
DeriVector2 cv (vertex, derivparam);
|
||||
DeriVector2 f1v (focus1, derivparam);
|
||||
DeriVector2 pv (p, derivparam);
|
||||
// fill some vectors in
|
||||
DeriVector2 cv(vertex, derivparam);
|
||||
DeriVector2 f1v(focus1, derivparam);
|
||||
DeriVector2 pv(p, derivparam);
|
||||
|
||||
// the normal is the vector from the focus to the intersection of ano thru the point p and direction
|
||||
// of the symmetry axis of the parabola with the directrix.
|
||||
// As both point to directrix and point to focus are of equal magnitude, we can work with unitary vectors
|
||||
// to calculate the normal, substraction of those vectors.
|
||||
// the normal is the vector from the focus to the intersection of ano thru the point p and
|
||||
// direction of the symmetry axis of the parabola with the directrix. As both point to directrix
|
||||
// and point to focus are of equal magnitude, we can work with unitary vectors to calculate the
|
||||
// normal, substraction of those vectors.
|
||||
|
||||
DeriVector2 ret = cv.subtr(f1v).getNormalized().subtr(f1v.subtr(pv).getNormalized());
|
||||
|
||||
@@ -613,9 +624,9 @@ DeriVector2 BSpline::CalculateNormal(const Point &p, const double* derivparam) c
|
||||
// place holder
|
||||
DeriVector2 ret;
|
||||
|
||||
// even if this method is call CalculateNormal, the returned vector is not the normal strictu sensus
|
||||
// but a normal vector, where the vector should point to the left when one walks along the curve from
|
||||
// start to end.
|
||||
// even if this method is call CalculateNormal, the returned vector is not the normal strictu
|
||||
// sensus but a normal vector, where the vector should point to the left when one walks along
|
||||
// the curve from start to end.
|
||||
//
|
||||
// https://forum.freecadweb.org/viewtopic.php?f=10&t=26312#p209486
|
||||
|
||||
@@ -736,10 +747,9 @@ double BSpline::getLinCombFactor(double x, size_t k, size_t i, unsigned int p)
|
||||
|
||||
for (size_t r = 1; r < p + 1; ++r) {
|
||||
for (size_t j = p; j > r - 1; --j) {
|
||||
double alpha =
|
||||
(x - flattenedknots[j + k - p]) /
|
||||
(flattenedknots[j + 1 + k - r] - flattenedknots[j + k - p]);
|
||||
d[j] = (1.0 - alpha) * d[j-1] + alpha * d[j];
|
||||
double alpha = (x - flattenedknots[j + k - p])
|
||||
/ (flattenedknots[j + 1 + k - r] - flattenedknots[j + k - p]);
|
||||
d[j] = (1.0 - alpha) * d[j - 1] + alpha * d[j];
|
||||
}
|
||||
}
|
||||
|
||||
@@ -751,9 +761,8 @@ double BSpline::splineValue(double x, size_t k, unsigned int p, VEC_D& d, const
|
||||
for (size_t r = 1; r < p + 1; ++r) {
|
||||
for (size_t j = p; j > r - 1; --j) {
|
||||
double alpha =
|
||||
(x - flatknots[j + k - p]) /
|
||||
(flatknots[j + 1 + k - r] - flatknots[j + k - p]);
|
||||
d[j] = (1.0 - alpha) * d[j-1] + alpha * d[j];
|
||||
(x - flatknots[j + k - p]) / (flatknots[j + 1 + k - r] - flatknots[j + k - p]);
|
||||
d[j] = (1.0 - alpha) * d[j - 1] + alpha * d[j];
|
||||
}
|
||||
}
|
||||
|
||||
|
||||
@@ -24,17 +24,27 @@
|
||||
#define PLANEGCS_GEO_H
|
||||
|
||||
#include <cmath>
|
||||
|
||||
#include "Util.h"
|
||||
|
||||
|
||||
namespace GCS
|
||||
{
|
||||
class Point
|
||||
{
|
||||
public:
|
||||
Point(){x = nullptr; y = nullptr;}
|
||||
Point(double *px, double *py) {x=px; y=py;}
|
||||
double *x;
|
||||
double *y;
|
||||
Point()
|
||||
{
|
||||
x = nullptr;
|
||||
y = nullptr;
|
||||
}
|
||||
Point(double* px, double* py)
|
||||
{
|
||||
x = px;
|
||||
y = py;
|
||||
}
|
||||
double* x;
|
||||
double* y;
|
||||
};
|
||||
|
||||
using VEC_P = std::vector<Point>;
|
||||
@@ -52,36 +62,76 @@ namespace GCS
|
||||
class DeriVector2
|
||||
{
|
||||
public:
|
||||
DeriVector2(){x=0; y=0; dx=0; dy=0;}
|
||||
DeriVector2(double x, double y) {this->x = x; this->y = y; this->dx = 0; this->dy = 0;}
|
||||
DeriVector2(double x, double y, double dx, double dy) {this->x = x; this->y = y; this->dx = dx; this->dy = dy;}
|
||||
DeriVector2(const Point &p, const double* derivparam);
|
||||
DeriVector2()
|
||||
{
|
||||
x = 0;
|
||||
y = 0;
|
||||
dx = 0;
|
||||
dy = 0;
|
||||
}
|
||||
DeriVector2(double x, double y)
|
||||
{
|
||||
this->x = x;
|
||||
this->y = y;
|
||||
this->dx = 0;
|
||||
this->dy = 0;
|
||||
}
|
||||
DeriVector2(double x, double y, double dx, double dy)
|
||||
{
|
||||
this->x = x;
|
||||
this->y = y;
|
||||
this->dx = dx;
|
||||
this->dy = dy;
|
||||
}
|
||||
DeriVector2(const Point& p, const double* derivparam);
|
||||
double x, dx;
|
||||
double y, dy;
|
||||
|
||||
double length() const {return sqrt(x*x + y*y);}
|
||||
double length(double &dlength) const; //returns length and writes length deriv into the dlength argument.
|
||||
double length() const
|
||||
{
|
||||
return sqrt(x * x + y * y);
|
||||
}
|
||||
double length(double& dlength)
|
||||
const;// returns length and writes length deriv into the dlength argument.
|
||||
|
||||
|
||||
//unlike other vectors in FreeCAD, this normalization creates a new vector instead of modifying existing one.
|
||||
DeriVector2 getNormalized() const; //returns zero vector if the original is zero.
|
||||
double scalarProd(const DeriVector2 &v2, double* dprd=nullptr) const;//calculates scalar product of two vectors and returns the result. The derivative of the result is written into argument dprd.
|
||||
DeriVector2 sum(const DeriVector2 &v2) const {//adds two vectors and returns result
|
||||
return DeriVector2(x + v2.x, y + v2.y,
|
||||
dx + v2.dx, dy + v2.dy);}
|
||||
DeriVector2 subtr(const DeriVector2 &v2) const {//subtracts two vectors and returns result
|
||||
return DeriVector2(x - v2.x, y - v2.y,
|
||||
dx - v2.dx, dy - v2.dy);}
|
||||
DeriVector2 mult(double val) const {
|
||||
return DeriVector2(x*val, y*val, dx*val, dy*val);}//multiplies the vector by a number. Derivatives are scaled.
|
||||
DeriVector2 multD(double val, double dval) const {//multiply vector by a variable with a derivative.
|
||||
return DeriVector2(x*val, y*val, dx*val+x*dval, dy*val+y*dval);}
|
||||
DeriVector2 divD(double val, double dval) const;//divide vector by a variable with a derivative
|
||||
DeriVector2 rotate90ccw() const {return DeriVector2(-y,x,-dy,dx);}
|
||||
DeriVector2 rotate90cw() const {return DeriVector2(y,-x,dy,-dx);}
|
||||
DeriVector2 linCombi(double m1, const DeriVector2 &v2, double m2) const {//linear combination of two vectors
|
||||
return DeriVector2(x*m1 + v2.x*m2, y*m1 + v2.y*m2,
|
||||
dx*m1 + v2.dx*m2, dy*m1 + v2.dy*m2);}
|
||||
// unlike other vectors in FreeCAD, this normalization creates a new vector instead of
|
||||
// modifying existing one.
|
||||
DeriVector2 getNormalized() const;// returns zero vector if the original is zero.
|
||||
double scalarProd(const DeriVector2& v2, double* dprd = nullptr)
|
||||
const;// calculates scalar product of two vectors and returns the result. The derivative
|
||||
// of the result is written into argument dprd.
|
||||
DeriVector2 sum(const DeriVector2& v2) const
|
||||
{// adds two vectors and returns result
|
||||
return DeriVector2(x + v2.x, y + v2.y, dx + v2.dx, dy + v2.dy);
|
||||
}
|
||||
DeriVector2 subtr(const DeriVector2& v2) const
|
||||
{// subtracts two vectors and returns result
|
||||
return DeriVector2(x - v2.x, y - v2.y, dx - v2.dx, dy - v2.dy);
|
||||
}
|
||||
DeriVector2 mult(double val) const
|
||||
{
|
||||
return DeriVector2(x * val, y * val, dx * val, dy * val);
|
||||
}// multiplies the vector by a number. Derivatives are scaled.
|
||||
DeriVector2 multD(double val, double dval) const
|
||||
{// multiply vector by a variable with a derivative.
|
||||
return DeriVector2(x * val, y * val, dx * val + x * dval, dy * val + y * dval);
|
||||
}
|
||||
DeriVector2 divD(double val,
|
||||
double dval) const;// divide vector by a variable with a derivative
|
||||
DeriVector2 rotate90ccw() const
|
||||
{
|
||||
return DeriVector2(-y, x, -dy, dx);
|
||||
}
|
||||
DeriVector2 rotate90cw() const
|
||||
{
|
||||
return DeriVector2(y, -x, dy, -dx);
|
||||
}
|
||||
DeriVector2 linCombi(double m1, const DeriVector2& v2, double m2) const
|
||||
{// linear combination of two vectors
|
||||
return DeriVector2(
|
||||
x * m1 + v2.x * m2, y * m1 + v2.y * m2, dx * m1 + v2.dx * m2, dy * m1 + v2.dy * m2);
|
||||
}
|
||||
|
||||
};
|
||||
|
||||
|
||||
@@ -22,8 +22,10 @@
|
||||
|
||||
#include <iostream>
|
||||
#include <iterator>
|
||||
|
||||
#include "SubSystem.h"
|
||||
|
||||
|
||||
namespace GCS
|
||||
{
|
||||
|
||||
@@ -56,44 +58,45 @@ void SubSystem::initialize(VEC_pD ¶ms, MAP_pD_pD &reductionmap)
|
||||
{
|
||||
SET_pD s1(params.begin(), params.end());
|
||||
SET_pD s2;
|
||||
for (std::vector<Constraint *>::iterator constr=clist.begin();
|
||||
constr != clist.end(); ++constr) {
|
||||
(*constr)->revertParams(); // ensure that the constraint points to the original parameters
|
||||
for (std::vector<Constraint*>::iterator constr = clist.begin(); constr != clist.end();
|
||||
++constr) {
|
||||
(*constr)
|
||||
->revertParams();// ensure that the constraint points to the original parameters
|
||||
VEC_pD constr_params = (*constr)->params();
|
||||
s2.insert(constr_params.begin(), constr_params.end());
|
||||
}
|
||||
std::set_intersection(s1.begin(), s1.end(), s2.begin(), s2.end(),
|
||||
std::back_inserter(tmpplist) );
|
||||
std::set_intersection(
|
||||
s1.begin(), s1.end(), s2.begin(), s2.end(), std::back_inserter(tmpplist));
|
||||
}
|
||||
|
||||
plist.clear();
|
||||
MAP_pD_I rindex;
|
||||
if (!reductionmap.empty()) {
|
||||
int i=0;
|
||||
int i = 0;
|
||||
MAP_pD_I pindex;
|
||||
for (VEC_pD::const_iterator itt=tmpplist.begin();
|
||||
itt != tmpplist.end(); ++itt) {
|
||||
for (VEC_pD::const_iterator itt = tmpplist.begin(); itt != tmpplist.end(); ++itt) {
|
||||
MAP_pD_pD::const_iterator itr = reductionmap.find(*itt);
|
||||
if (itr != reductionmap.end()) {
|
||||
MAP_pD_I::const_iterator itp = pindex.find(itr->second);
|
||||
if (itp == pindex.end()) { // the reduction target is not in plist yet, so add it now
|
||||
if (itp == pindex.end()) {// the reduction target is not in plist yet, so add it now
|
||||
plist.push_back(itr->second);
|
||||
rindex[itr->first] = i;
|
||||
pindex[itr->second] = i;
|
||||
i++;
|
||||
}
|
||||
else // the reduction target is already in plist, just inform rindex
|
||||
else// the reduction target is already in plist, just inform rindex
|
||||
rindex[itr->first] = itp->second;
|
||||
}
|
||||
else if (pindex.find(*itt) == pindex.end()) { // not in plist yet, so add it now
|
||||
else if (pindex.find(*itt) == pindex.end()) {// not in plist yet, so add it now
|
||||
plist.push_back(*itt);
|
||||
pindex[*itt] = i;
|
||||
i++;
|
||||
}
|
||||
}
|
||||
}
|
||||
else
|
||||
else {
|
||||
plist = tmpplist;
|
||||
}
|
||||
|
||||
psize = static_cast<int>(plist.size());
|
||||
pvals.resize(psize);
|
||||
@@ -329,9 +332,8 @@ void SubSystem::applySolution()
|
||||
*(it->first) = *(it->second);
|
||||
}
|
||||
|
||||
void SubSystem::analyse(Eigen::MatrixXd & /*J*/, Eigen::MatrixXd & /*ker*/, Eigen::MatrixXd & /*img*/)
|
||||
{
|
||||
}
|
||||
void SubSystem::analyse(Eigen::MatrixXd& /*J*/, Eigen::MatrixXd& /*ker*/, Eigen::MatrixXd& /*img*/)
|
||||
{}
|
||||
|
||||
void SubSystem::report()
|
||||
{
|
||||
|
||||
@@ -27,8 +27,10 @@
|
||||
#undef max
|
||||
|
||||
#include <Eigen/Core>
|
||||
|
||||
#include "Constraints.h"
|
||||
|
||||
|
||||
namespace GCS
|
||||
{
|
||||
|
||||
|
||||
@@ -23,9 +23,10 @@
|
||||
#ifndef PLANEGCS_UTIL_H
|
||||
#define PLANEGCS_UTIL_H
|
||||
|
||||
#include <vector>
|
||||
#include <map>
|
||||
#include <set>
|
||||
#include <vector>
|
||||
|
||||
|
||||
namespace GCS
|
||||
{
|
||||
|
||||
Reference in New Issue
Block a user