OpticConstraints: Adding Snell's law

Fix AngleViaPoint to support new derivative calculation technique.

OpticConstraints: Adding Snell's law. Fix AngleViaPoint to support new derivative calculation technique.

Snell's law constraint added to GCS, but not yet exposed and cannot be
tested.
Since the way CalculateNormal() returns derivatives had changed,
AngleViaPoint constraint needed modifications. Nothing serious.

OpticConstraints: SnellsLaw progress

Addable through python. Fix math. Some quick-and-dirty visual stuff to
get rid of hangs and to see the constraint in action.

OpticConstraints: SnellsLaw: flipping logic fix

OpticConstraints: SnellsLaw progress

Added toolbar button. Allowed editing a datum by doubleclick. New error
message approach during constraint creation.

OpticConstraints: SnellsLaw

OpticConstraints: SnellsLaw: list label improvement

OpticConstraints: SnellsLaw: fix after rebase

OpticConstraints: SnellsLaw: expose helper constraints

Snell's law internally is made of three constraints: point-on-object,
coincident and the Snell's sin/sin. They were all buried under one UI
constraint. Exposing them allows to construct reflection and
birefringence on the point (attempting to do so used to result in
redundant constraints and was often not functional at all).
This commit breaks compatibility with older files.

OpticConstraints: SnellsLaw: small refactor of math

Placing the duplicated code of error and gradient calculation into a
private method.

OpticConstraints: SnellsLaw: fix datum edit unit

OpticConstraints: SnellsLaw: fix datum edit bug

After previous fix, the dimensionless value was not accepted (the
constraint's value did not change, the changes were ignored).

src/Mod/Sketcher/App/Constraint.h

@@ -46,7 +46,8 @@ enum ConstraintType {
     Equal = 12,
     PointOnObject = 13,
     Symmetric = 14,
-    InternalAlignment = 15
+    InternalAlignment = 15,
+    SnellsLaw = 16
 };
 
 enum InternalAlignmentType {

src/Mod/Sketcher/App/ConstraintPyImp.cpp

@@ -378,17 +378,34 @@ int ConstraintPy::PyInit(PyObject* args, PyObject* /*kwd*/)
     }
     PyErr_Clear();
 
-    if (PyArg_ParseTuple(args, "siiiiii", &ConstraintType, &FirstIndex, &FirstPos, &SecondIndex, &SecondPos, &ThirdIndex, &ThirdPos)) {
-        // ConstraintType, GeoIndex1, PosIndex1, GeoIndex2, PosIndex2, GeoIndex3, PosIndex3
-        if (strcmp("Symmetric",ConstraintType) == 0 ) {
-            this->getConstraintPtr()->Type = Symmetric;
-            this->getConstraintPtr()->First     = FirstIndex;
-            this->getConstraintPtr()->FirstPos  = (Sketcher::PointPos) FirstPos;
-            this->getConstraintPtr()->Second    = SecondIndex;
-            this->getConstraintPtr()->SecondPos = (Sketcher::PointPos) SecondPos;
-            this->getConstraintPtr()->Third     = ThirdIndex;
-            this->getConstraintPtr()->ThirdPos  = (Sketcher::PointPos) ThirdPos;
-            return 0;
+    if (PyArg_ParseTuple(args, "siiiiiO", &ConstraintType, &FirstIndex, &FirstPos, &SecondIndex, &SecondPos, &ThirdIndex, &index_or_value)) {
+        if (PyInt_Check(index_or_value)) {
+            ThirdPos = PyInt_AsLong(index_or_value);
+            // ConstraintType, GeoIndex1, PosIndex1, GeoIndex2, PosIndex2, GeoIndex3, PosIndex3
+            if (strcmp("Symmetric",ConstraintType) == 0 ) {
+                this->getConstraintPtr()->Type = Symmetric;
+                this->getConstraintPtr()->First     = FirstIndex;
+                this->getConstraintPtr()->FirstPos  = (Sketcher::PointPos) FirstPos;
+                this->getConstraintPtr()->Second    = SecondIndex;
+                this->getConstraintPtr()->SecondPos = (Sketcher::PointPos) SecondPos;
+                this->getConstraintPtr()->Third     = ThirdIndex;
+                this->getConstraintPtr()->ThirdPos  = (Sketcher::PointPos) ThirdPos;
+                return 0;
+            }
+        }
+        if (PyNumber_Check(index_or_value)) { // can be float or int
+            Value = PyFloat_AsDouble(index_or_value);
+            if (strcmp("SnellsLaw",ConstraintType) == 0 ) {
+                this->getConstraintPtr()->Type = SnellsLaw;
+                this->getConstraintPtr()->First     = FirstIndex;
+                this->getConstraintPtr()->FirstPos  = (Sketcher::PointPos) FirstPos;
+                this->getConstraintPtr()->Second    = SecondIndex;
+                this->getConstraintPtr()->SecondPos = (Sketcher::PointPos) SecondPos;
+                this->getConstraintPtr()->Third     = ThirdIndex;
+                this->getConstraintPtr()->ThirdPos  = none;
+                this->getConstraintPtr()->Value = Value;
+                return 0;
+            }
         }
     }
 
@@ -436,6 +453,7 @@ std::string ConstraintPy::representation(void) const
             else
                 result << "'AngleViaPoint'>";
         break;
+        case SnellsLaw          : result << "'SnellsLaw'>"; break;
         case InternalAlignment  : 
             switch(this->getConstraintPtr()->AlignmentType) {
                 case Undef                  : result << "'InternalAlignment:Undef'>";break;

src/Mod/Sketcher/App/Sketch.cpp

@@ -766,6 +766,13 @@ int Sketch::addConstraint(const Constraint *constraint)
                 break;
         }
         break;
+    case SnellsLaw:
+        assert(constraint->ThirdPos==none);
+        rtn = addSnellsLawConstraint(constraint->First, constraint->FirstPos,
+                                     constraint->Second, constraint->SecondPos,
+                                     constraint->Third,
+                                     constraint->Value);
+        break;
     case None:
         break;
     }
@@ -1642,6 +1649,69 @@ int Sketch::addSymmetricConstraint(int geoId1, PointPos pos1, int geoId2, PointP
     return -1;
 }
 
+int Sketch::addSnellsLawConstraint(int geoIdRay1, PointPos posRay1,
+                                   int geoIdRay2, PointPos posRay2,
+                                   int geoIdBnd,
+                                   double n2divn1)
+{
+
+    geoIdRay1 = checkGeoId(geoIdRay1);
+    geoIdRay2 = checkGeoId(geoIdRay2);
+    geoIdBnd = checkGeoId(geoIdBnd);
+
+    if (Geoms[geoIdRay1].type == Point ||
+        Geoms[geoIdRay2].type == Point){
+        assert(0);//point is not a curve. Not applicable!
+        return -1;
+    }
+
+    GCS::Curve* ray1 =getGCSCurveByGeoId(geoIdRay1);
+    GCS::Curve* ray2 =getGCSCurveByGeoId(geoIdRay2);
+    GCS::Curve* boundary =getGCSCurveByGeoId(geoIdBnd);
+    if (!ray1 || !ray2 || !boundary) {
+        assert(0);
+        Base::Console().Error("addSnellsLawConstraint: getGCSCurveByGeoId returned NULL!\n");
+        return -1;
+    }
+
+    int pointId1 = getPointId(geoIdRay1, posRay1);
+    int pointId2 = getPointId(geoIdRay2, posRay2);
+    if ( pointId1 < 0 || pointId1 >= int(Points.size()) ||
+         pointId2 < 0 || pointId2 >= int(Points.size()) ){
+        assert(0);
+        Base::Console().Error("addSnellsLawConstraint: point index out of range.\n");
+        return -1;
+    }
+    GCS::Point &p1 = Points[pointId1];
+    GCS::Point &p2 = Points[pointId2];
+
+    // add the parameters (refractive indexes)
+    FixParameters.push_back(new double(0.0));
+    double *n1 = FixParameters[FixParameters.size()-1];
+    FixParameters.push_back(new double(0.0));
+    double *n2 = FixParameters[FixParameters.size()-1];
+
+    if ( abs(n2divn1) >= 1.0 ){
+        *n2 = n2divn1;
+        *n1 = 1.0;
+    } else {
+        *n2 = 1.0;
+        *n1 = 1/n2divn1;
+    }
+
+    int tag = -1;
+    //tag = Sketch::addPointOnObjectConstraint(geoIdRay1, posRay1, geoIdBnd);//increases ConstraintsCounter
+    tag = ++ConstraintsCounter;
+    //GCSsys.addConstraintP2PCoincident(p1, p2, tag);
+    GCSsys.addConstraintSnellsLaw(*ray1, *ray2,
+                                  *boundary, p1,
+                                  n1, n2,
+                                  posRay1==start, posRay2 == end,
+                                  tag);
+    return ConstraintsCounter;
+}
+
+
 int Sketch::addInternalAlignmentEllipseMajorDiameter(int geoId1, int geoId2)
 {
     std::swap(geoId1, geoId2);

src/Mod/Sketcher/App/Sketch.h

@@ -187,6 +187,11 @@ class SketcherExport Sketch :public Base::Persistence
     int addSymmetricConstraint(int geoId1, PointPos pos1, int geoId2, PointPos pos2, int geoId3);
     /// add a symmetric constraint between three points, the last point is in the middle of the first two
     int addSymmetricConstraint(int geoId1, PointPos pos1, int geoId2, PointPos pos2, int geoId3, PointPos pos3);
+    /// add a snell's law constraint
+    int addSnellsLawConstraint(int geoIdRay1, PointPos posRay1,
+                               int geoIdRay2, PointPos posRay2,
+                               int geoIdBnd,
+                               double n2divn1);
     //@}
 
     /// Internal Alignment constraints

src/Mod/Sketcher/App/SketchObject.cpp

@@ -196,7 +196,8 @@ int SketchObject::setDatum(int ConstrId, double Datum)
         type != Radius &&
         type != Angle &&
         type != Tangent && //for tangent, value==0 is autodecide, value==Pi/2 is external and value==-Pi/2 is internal
-        type != Perpendicular)
+        type != Perpendicular &&
+        type != SnellsLaw)
         return -1;
 
     if ((type == Distance || type == Radius) && Datum <= 0)