I have done another fundamental change to GCS. I'm very proud of it.

The change allows to do vector math in error functions without calculating partials - the code to calculate partials now coincides with the code to calculate the error function.

It was done by converting GCS::Vector2D into GCS::DeriVector2. The new vector class has two additional fields (dx,dy) to store vector's derivative on current parameter.

So, the calculation of ellipse's normal vector and its derivatives used to be like this:

Code: Select all

```
//--------------ellipse
171 Vector2D Ellipse::CalculateNormal(Point &p, double* derivparam)
172 {
173 Vector2D cv (*center.x, *center.y);
174 Vector2D f1v (*focus1X, *focus1Y);
175 Vector2D pv (*p.x, *p.y);
176
177 Vector2D ret(0.0, 0.0);
178
179 Vector2D f2v ( 2*cv.x - f1v.x, 2*cv.y - f1v.y ); //position of focus2
180 if(derivparam){
181
182 Vector2D dc;
183 Vector2D df1;
184 Vector2D dp;
185
186 if (derivparam == center.x) dc.x = 1.0;
187 if (derivparam == center.y) dc.y = 1.0;
188 if (derivparam == focus1X) df1.x = 1.0;
189 if (derivparam == focus1Y) df1.y = 1.0;
190 if (derivparam == p.x) dp.x = 1.0;
191 if (derivparam == p.y) dp.y = 1.0;
192 //todo: exit if all are zero
193
194 Vector2D pf1 = Vector2D(pv.x - f1v.x, pv.y - f1v.y);//same as during error function calculation. I reuse the values during derivative calculation
195 Vector2D pf2 = Vector2D(pv.x - f2v.x, pv.y - f2v.y);
196 Vector2D pf1e = pf1.getNormalized();
197 Vector2D pf2e = pf2.getNormalized();
198
199 Vector2D df2 (2*dc.x - df1.x, 2*dc.y - df1.y );
200
201 Vector2D dpf1 (dp.x - df1.x, dp.y - df1.y);//derivative before normalization
202 Vector2D dpf1e (dpf1.x/pf1.length(), dpf1.y/pf1.length());//first portion of normalization derivative (normalized' = unnormalized'/len + unnormalized*(1/len)')
203 dpf1e.x += -pf1.x/pow(pf1.length(),2)*(dpf1.x*pf1e.x + dpf1.y*pf1e.y);//second part of normalization dreivative
204 dpf1e.y += -pf1.y/pow(pf1.length(),2)*(dpf1.x*pf1e.x + dpf1.y*pf1e.y);
205
206 Vector2D dpf2 (dp.x - df2.x, dp.y - df2.y);//same stuff for pf2
207 Vector2D dpf2e (dpf2.x/pf2.length(), dpf2.y/pf2.length());//first portion of normalization derivative (normalized' = unnormalized'/len + unnormalized*(1/len)')
208 dpf2e.x += -pf2.x/pow(pf2.length(),2)*(dpf2.x*pf2e.x + dpf2.y*pf2e.y);//second part of normalization dreivative
209 dpf2e.y += -pf2.y/pow(pf2.length(),2)*(dpf2.x*pf2e.x + dpf2.y*pf2e.y);
210
211 ret.x = -(dpf1e.x + dpf2e.x);
212 ret.y = -(dpf1e.y + dpf2e.y);//DeepSOIC: derivative calculated manually... error-prone =) Tested, fixed, looks good.
.........(numeric derivative testing code skipped)
232
233 } else {
234 Vector2D pf1 = Vector2D(pv.x - f1v.x, pv.y - f1v.y);
235 Vector2D pf2 = Vector2D(pv.x - f2v.x, pv.y - f2v.y);
236 Vector2D pf1e = pf1.getNormalized();
237 Vector2D pf2e = pf2.getNormalized();
238 ret.x = -(pf1e.x + pf2e.x);
239 ret.y = -(pf1e.y + pf2e.y);
240 };
241
242 return ret;
243 }
```

Now, this is simply:

(both the normal vector and derivative are calculated at the same time)

Code: Select all

```
DeriVector2 Ellipse::CalculateNormal(Point &p, double* derivparam)
182 {
183 //fill some vectors in
184 DeriVector2 cv (center, derivparam);
185 DeriVector2 f1v (focus1, derivparam);
186 DeriVector2 pv (p, derivparam);
187
188 //calculation.
189 //focus2:
190 DeriVector2 f2v = cv.linCombi(2.0, f1v, -1.0); // 2*cv - f1v
191
192 //pf1, pf2 = vectors from p to focus1,focus2
193 DeriVector2 pf1 = f1v.subtr(pv);
194 DeriVector2 pf2 = f2v.subtr(pv);
195 //return sum of normalized pf2, pf2
196 DeriVector2 ret = pf1.getNormalized().sum(pf2.getNormalized());
197
.........(numeric derivative testing code skipped)
218
219 return ret;
220 }
```

This was done to:

+ simplify understanding of the math.

+ eliminate hassle of derivative calculation, thus simplifying adding new constraints.

+ hugely simplify testing of different error function formulations

+ reduce potential for mistakes

The drawbacks are:

- for newcomers, it may be not trivial to understand how the derivatives are actually calculated.

- the derivatives are less simplified, so there is performance drop.

- Derivative calculation code gets executed even when there is no need to do so (e.g. when calculating error() ), but since error calculations are executed O(n) times compared to O(n^2) for derivatives, this isn't much of a penalty.