1 /*
  2     Copyright 2010
  3         Matthias Ehmann,
  4         Michael Gerhaeuser,
  5         Carsten Miller,
  6         Bianca Valentin,
  7         Alfred Wassermann,
  8         Peter Wilfahrt
  9 
 10     This file is part of JSXGraph.
 11 
 12     JSXGraph is free software: you can redistribute it and/or modify
 13     it under the terms of the GNU Lesser General Public License as published by
 14     the Free Software Foundation, either version 3 of the License, or
 15     (at your option) any later version.
 16 
 17     JSXGraph is distributed in the hope that it will be useful,
 18     but WITHOUT ANY WARRANTY; without even the implied warranty of
 19     MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 20     GNU Lesser General Public License for more details.
 21 
 22     You should have received a copy of the GNU Lesser General Public License
 23     along with JSXGraph. If not, see <http://www.gnu.org/licenses/>.
 24 */    
 25 
 26 /**
 27  * @fileoverview In this file the conic sections defined.
 28  */
 29 
 30 /**
 31  * @class This element is used to provide a constructor for an ellipse. An ellipse is given by two points (the foci) and a third point on the the ellipse or 
 32  * the length of the major axis.
 33  * @pseudo
 34  * @description
 35  * @name Ellipse
 36  * @augments JXG.Curve
 37  * @constructor
 38  * @type JXG.Curve
 39  * @throws {Exception} If the element cannot be constructed with the given parent objects an exception is thrown.
 40  * @param {JXG.Point,array_JXG.Point,array_JXG.Point,array} point1,point2,point3 Parent elements can be three elements either of type {@link JXG.Point} or array of
 41  * numbers describing the coordinates of a point. In the latter case the point will be constructed automatically as a fixed invisible point.
 42  * @param {JXG.Point,array_JXG.Point,array_number,function} point1,point2,number Parent elements can be two elements either of type {@link JXG.Point} or array of
 43  * numbers describing the coordinates of a point. The third parameter is a number/function which defines the length of the major axis
 44  * Optional parameters four and five are numbers which define the curve length (e.g. start/end). Default values are -pi and pi.
 45  * @example
 46  * // Create an Ellipse by three points
 47  * var A = board.create('point', [-1,4]);
 48  * var B = board.create('point', [-1,-4]);
 49  * var C = board.create('point', [1,1]);
 50  * var el = board.create('ellipse',[A,B,C]);
 51  * </pre><div id="a4d7fb6f-8708-4e45-87f2-2379ae2bd2c0" style="width: 300px; height: 300px;"></div>
 52  * <script type="text/javascript">
 53  *   var glex1_board = JXG.JSXGraph.initBoard('a4d7fb6f-8708-4e45-87f2-2379ae2bd2c0', {boundingbox:[-6,6,6,-6], keepaspectratio:true, showcopyright: false, shownavigation: false});
 54  *   var A = glex1_board.create('point', [-1,4]);
 55  *   var B = glex1_board.create('point', [-1,-4]);
 56  *   var C = glex1_board.create('point', [1,1]);
 57  *   var el = glex1_board.create('ellipse',[A,B,C]);
 58  * </script><pre>
 59  */
 60 JXG.createEllipse = function(board, parents, attributes) {
 61     var F = [],  // focus 1 and focus 2
 62         C, majorAxis, i,
 63         rotationMatrix,
 64         attr_foci = JXG.copyAttributes(attributes, board.options, 'conic', 'foci'),
 65         attr_curve = JXG.copyAttributes(attributes, board.options, 'conic');
 66 
 67     // The foci and the third point are either points or coordinate arrays.
 68     for (i = 0; i < 2; i++) {
 69         if (parents[i].length > 1) { // focus i given by coordinates
 70             F[i] = board.create('point', parents[i], attr_foci);
 71         } else if (JXG.isPoint(parents[i])) { // focus i given by point
 72             F[i] = JXG.getReference(board,parents[i]);
 73         } else if ((typeof parents[i] == 'function') && (parents[i]().elementClass === JXG.OBJECT_CLASS_POINT)) {  // given by function
 74             F[i] = parents[i]();
 75         } else if (JXG.isString(parents[i])) { // focus i given by point name
 76             F[i] = JXG.getReference(board,parents[i]);
 77         } else
 78             throw new Error("JSXGraph: Can't create Ellipse with parent types '" + 
 79                             (typeof parents[0]) + "' and '" + (typeof parents[1]) + "'." +
 80                             "\nPossible parent types: [point,point,point], [point,point,number|function]");
 81     }
 82     if (JXG.isNumber(parents[2])) { // length of major axis
 83         majorAxis = JXG.createFunction(parents[2],board);
 84     } else if ((typeof parents[2] == 'function') && (JXG.isNumber(parents[2]()))){
 85         majorAxis = parents[2];
 86     } else {
 87         if (JXG.isPoint(parents[2])) {                                               // point on ellipse
 88             C = JXG.getReference(board,parents[2]);
 89         } else if (parents[2].length>1) {                                            // point on ellipse given by coordinates
 90             C = board.create('point', parents[2], attr_foci);
 91         } else if ((typeof parents[2] == 'function') && (parents[2]().elementClass == JXG.OBJECT_CLASS_POINT)) {  // given by function
 92             C = parents[2]();
 93         } else if (JXG.isString(parents[2])) {                                      // focus i given by point name
 94             C = JXG.getReference(board,parents[2]);
 95         } else {
 96             throw new Error("JSXGraph: Can't create Ellipse with parent types '" + 
 97                             (typeof parents[0]) + "' and '" + (typeof parents[1]) + "' and '" + (typeof parents[2]) +"'." +
 98                             "\nPossible parent types: [point,point,point], [point,point,number|function]");
 99         }
100         majorAxis = function(){ return C.Dist(F[0])+C.Dist(F[1]);};
101     }
102 
103     if (typeof parents[4]=='undefined') parents[4] = 1.0001*Math.PI;   // to
104     if (typeof parents[3]=='undefined') parents[3] = -1.0001*Math.PI;  // from
105 
106     var M = board.create('point', [
107                 function(){return (F[0].X()+F[1].X())*0.5;},
108                 function(){return (F[0].Y()+F[1].Y())*0.5;}
109             ], attr_foci);
110 
111     var transformFunc = function() {
112             var ax = F[0].X(),
113                 ay = F[0].Y(),
114                 bx = F[1].X(),
115                 by = F[1].Y(),
116                 beta, co, si;
117 
118             // Rotate by the slope of the line [F[0],F[1]]
119             var sgn = (bx-ax>0)?1:-1;
120             if (Math.abs(bx-ax)>0.0000001) {
121                 beta = Math.atan2(by-ay,bx-ax)+ ((sgn<0)?Math.PI:0);
122             } else {
123                 beta = ((by-ay>0)?0.5:-0.5)*Math.PI;
124             }
125             co = Math.cos(beta);
126             si = Math.sin(beta);
127             var m = [
128                         [1,     0,  0],
129                         [M.X(),co,-si],
130                         [M.Y(),si, co]
131                     ];
132             return m;
133         };
134 
135     var curve = board.create('curve', [function(x) {return 0;}, function(x) {return 0;}, parents[3], parents[4]], attr_curve);
136 
137     var polarForm = function(phi,suspendUpdate) {
138                 var a = majorAxis()*0.5,
139                     aa = a*a,
140                     e = F[1].Dist(F[0])*0.5,
141                     bb = aa-e*e,
142                     b = Math.sqrt(bb), 
143                     transformMat = [[1,0,0],[0,1,0],[0,0,1]],
144                     mx, my;
145                     
146                 if (!suspendUpdate) {
147                     rotationMatrix = transformFunc();
148                     mx = M.X();
149                     my = M.Y();
150                     transformMat[0][0] = rotationMatrix[0][0];
151                     transformMat[0][1] = 0.0;
152                     transformMat[0][2] = 0.0;
153                     transformMat[1][0] = mx*(1-rotationMatrix[1][1])+my*rotationMatrix[1][2];
154                     transformMat[1][1] = rotationMatrix[1][1];
155                     transformMat[1][2] = rotationMatrix[2][1];
156                     transformMat[2][0] = my*(1-rotationMatrix[1][1])-mx*rotationMatrix[1][2];
157                     transformMat[2][1] = rotationMatrix[1][2];
158                     transformMat[2][2] = rotationMatrix[2][2];
159                     curve.quadraticform = 
160                         JXG.Math.matMatMult(JXG.Math.transpose(transformMat),
161                         JXG.Math.matMatMult(
162                             [
163                                 [-1+mx*mx/(a*a)+my*my/bb, -mx/aa , -mx/bb],
164                                 [-mx/aa                 ,   1/aa ,  0    ],
165                                 [-my/bb                 ,      0 ,  1/bb ]
166                             ],
167                         transformMat)); 
168                 }
169                 return JXG.Math.matVecMult(rotationMatrix,[1,a*Math.cos(phi),b*Math.sin(phi)]);
170         };
171 
172     curve.X = function(phi,suspendUpdate) {return polarForm(phi,suspendUpdate)[1];};
173     curve.Y = function(phi,suspendUpdate) {return polarForm(phi,suspendUpdate)[2];};
174     curve.midpoint = M;
175     curve.type = JXG.OBJECT_TYPE_CONIC;
176     return curve;
177 };
178 
179 /**
180  * @class This element is used to provide a constructor for an hyperbola. An hyperbola is given by two points (the foci) and a third point on the the hyperbola or 
181  * the length of the major axis.
182  * @pseudo
183  * @description
184  * @name Hyperbola
185  * @augments JXG.Curve
186  * @constructor
187  * @type JXG.Curve
188  * @throws {Exception} If the element cannot be constructed with the given parent objects an exception is thrown.
189  * @param {JXG.Point,array_JXG.Point,array_JXG.Point,array} point1,point2,point3 Parent elements can be three elements either of type {@link JXG.Point} or array of
190  * numbers describing the coordinates of a point. In the latter case the point will be constructed automatically as a fixed invisible point.
191  * @param {JXG.Point,array_JXG.Point,array_number,function} point1,point2,number Parent elements can be two elements either of type {@link JXG.Point} or array of
192  * numbers describing the coordinates of a point. The third parameter is a number/function which defines the length of the major axis
193  * Optional parameters four and five are numbers which define the curve length (e.g. start/end). Default values are -pi and pi.
194  * @example
195  * // Create an Hyperbola by three points
196  * var A = board.create('point', [-1,4]);
197  * var B = board.create('point', [-1,-4]);
198  * var C = board.create('point', [1,1]);
199  * var el = board.create('hyperbola',[A,B,C]);
200  * </pre><div id="cf99049d-a3fe-407f-b936-27d76550f8c4" style="width: 300px; height: 300px;"></div>
201  * <script type="text/javascript">
202  *   var glex1_board = JXG.JSXGraph.initBoard('cf99049d-a3fe-407f-b936-27d76550f8c4', {boundingbox:[-6,6,6,-6], keepaspectratio:true, showcopyright: false, shownavigation: false});
203  *   var A = glex1_board.create('point', [-1,4]);
204  *   var B = glex1_board.create('point', [-1,-4]);
205  *   var C = glex1_board.create('point', [1,1]);
206  *   var el = glex1_board.create('hyperbola',[A,B,C]);
207  * </script><pre>
208  */
209 JXG.createHyperbola = function(board, parents, attributes) {
210     var F = [],  // focus 1 and focus 2
211         C, 
212         majorAxis,
213         i,
214         rotationMatrix,
215         attr_foci = JXG.copyAttributes(attributes, board.options, 'conic', 'foci'),
216         attr_curve = JXG.copyAttributes(attributes, board.options, 'conic');
217 
218     // The foci and the third point are either points or coordinate arrays.
219     for (i=0;i<2;i++) {
220         if (parents[i].length>1) { // focus i given by coordinates
221             F[i] = board.create('point', parents[i], attr_focu);
222         } else if (JXG.isPoint(parents[i])) { // focus i given by point
223             F[i] = JXG.getReference(board,parents[i]);
224         } else if ((typeof parents[i] == 'function') && (parents[i]().elementClass == JXG.OBJECT_CLASS_POINT)) {  // given by function
225             F[i] = parents[i]();
226         } else if (JXG.isString(parents[i])) { // focus i given by point name
227             F[i] = JXG.getReference(board,parents[i]);
228         } else
229             throw new Error("JSXGraph: Can't create Hyperbola with parent types '" + 
230                             (typeof parents[0]) + "' and '" + (typeof parents[1]) + "'." +
231                             "\nPossible parent types: [point,point,point], [point,point,number|function]");
232     }
233     if (JXG.isNumber(parents[2])) { // length of major axis
234         majorAxis = JXG.createFunction(parents[2],board);
235     } else if ((typeof parents[2] == 'function') && (JXG.isNumber(parents[2]()))){
236         majorAxis = parents[2];
237     } else {
238         if (JXG.isPoint(parents[2])) {                                               // point on ellipse
239             C = JXG.getReference(board,parents[2]);
240         } else if (parents[2].length>1) {                                            // point on ellipse given by coordinates
241             C = board.create('point', parents[2], attr_foci);
242         } else if ((typeof parents[2] == 'function') && (parents[2]().elementClass == JXG.OBJECT_CLASS_POINT)) {  // given by function
243             C = parents[2]();
244         } else if (JXG.isString(parents[2])) {                                      // focus i given by point name
245             C = JXG.getReference(board,parents[2]);
246         } else {
247             throw new Error("JSXGraph: Can't create Hyperbola with parent types '" + 
248                             (typeof parents[0]) + "' and '" + (typeof parents[1]) + "' and '" + (typeof parents[2]) +"'." +
249                             "\nPossible parent types: [point,point,point], [point,point,number|function]");
250         }
251         majorAxis = function(){ return C.Dist(F[0])-C.Dist(F[1]);};
252     }
253 
254     if (typeof parents[4]=='undefined') parents[4] = 1.0001*Math.PI;   // to
255     if (typeof parents[3]=='undefined') parents[3] = -1.0001*Math.PI;  // from
256 
257     var M = board.create('point', [
258                 function(){return (F[0].X()+F[1].X())*0.5;},
259                 function(){return (F[0].Y()+F[1].Y())*0.5;}
260             ], attr_foci);
261 
262     var transformFunc = function() {
263             var ax = F[0].X(),
264                 ay = F[0].Y(),
265                 bx = F[1].X(),
266                 by = F[1].Y(),
267                 beta;
268 
269             // Rotate by the slope of the line [F[0],F[1]]
270             var sgn = (bx-ax>0)?1:-1;
271             if (Math.abs(bx-ax)>0.0000001) {
272                 beta = Math.atan2(by-ay,bx-ax)+ ((sgn<0)?Math.PI:0);
273             } else {
274                 beta = ((by-ay>0)?0.5:-0.5)*Math.PI;
275             }
276             var m = [
277                         [1,    0,             0],
278                         [M.X(),Math.cos(beta),-Math.sin(beta)],
279                         [M.Y(),Math.sin(beta), Math.cos(beta)]
280                     ];
281             return m;
282         };
283 
284     var curve = board.create('curve',[function(x) {return 0;},function(x) {return 0;},parents[3],parents[4]], attr_curve);
285     /*
286     * Hyperbola is defined by (a*sec(t),b*tan(t)) and sec(t) = 1/cos(t)
287     */
288     var polarForm = function(phi,suspendUpdate) {
289                 var a = majorAxis()*0.5,
290                     aa = a*a,
291                     e = F[1].Dist(F[0])*0.5,
292                     b = Math.sqrt(-a*a+e*e), 
293                     bb = b*b,
294                     transformMat = [[1,0,0],[0,1,0],[0,0,1]],
295                     mx, my;
296                     
297                 if (!suspendUpdate) {
298                     rotationMatrix = transformFunc();
299                     mx = M.X();
300                     my = M.Y();
301                     transformMat[0][0] = rotationMatrix[0][0];
302                     transformMat[0][1] = 0.0;
303                     transformMat[0][2] = 0.0;
304                     transformMat[1][0] = mx*(1-rotationMatrix[1][1])+my*rotationMatrix[1][2];
305                     transformMat[1][1] = rotationMatrix[1][1];
306                     transformMat[1][2] = rotationMatrix[2][1];
307                     transformMat[2][0] = my*(1-rotationMatrix[1][1])-mx*rotationMatrix[1][2];
308                     transformMat[2][1] = rotationMatrix[1][2];
309                     transformMat[2][2] = rotationMatrix[2][2];
310                     curve.quadraticform = 
311                         JXG.Math.matMatMult(JXG.Math.transpose(transformMat),
312                         JXG.Math.matMatMult(
313                             [
314                                 [-1+mx*mx/aa+my*my/bb, -mx/aa , my/bb],
315                                 [-mx/aa              ,    1/aa,  0   ],
316                                 [my/bb               ,      0 , -1/bb]
317                             ],
318                         transformMat)); 
319                 }
320                 return JXG.Math.matVecMult(rotationMatrix,[1,a/Math.cos(phi),b*Math.tan(phi)]);
321         };
322 
323     curve.X = function(phi,suspendUpdate) {return polarForm(phi,suspendUpdate)[1];};
324     curve.Y = function(phi,suspendUpdate) {return polarForm(phi,suspendUpdate)[2];};
325     curve.midpoint = M;
326     curve.type = JXG.OBJECT_TYPE_CONIC;
327     return curve;
328 };
329 
330 /**
331  * @class This element is used to provide a constructor for a parabola. A parabola is given by one point (the focus) and a line (the directrix).
332  * @pseudo
333  * @description
334  * @name Parabola
335  * @augments JXG.Curve
336  * @constructor
337  * @type JXG.Curve
338  * @throws {Exception} If the element cannot be constructed with the given parent objects an exception is thrown.
339  * @param {JXG.Point,array_JXG.Line} point,line Parent elements are a point and a line.
340  * Optional parameters three and four are numbers which define the curve length (e.g. start/end). Default values are -pi and pi.
341  * @example
342  * // Create a parabola by a point C and a line l.
343  * var A = board.create('point', [-1,4]);
344  * var B = board.create('point', [-1,-4]);
345  * var l = board.create('line', [A,B]);
346  * var C = board.create('point', [1,1]);
347  * var el = board.create('parabola',[C,l]);
348  * </pre><div id="524d1aae-217d-44d4-ac58-a19c7ab1de36" style="width: 300px; height: 300px;"></div>
349  * <script type="text/javascript">
350  *   var glex1_board = JXG.JSXGraph.initBoard('524d1aae-217d-44d4-ac58-a19c7ab1de36', {boundingbox:[-6,6,6,-6], keepaspectratio:true, showcopyright: false, shownavigation: false});
351  *   var A = glex1_board.create('point', [-1,4]);
352  *   var B = glex1_board.create('point', [-1,-4]);
353  *   var l = glex1_board.create('line', [A,B]);
354  *   var C = glex1_board.create('point', [1,1]);
355  *   var el = glex1_board.create('parabola',[C,l]);
356  * </script><pre>
357  */
358 JXG.createParabola = function(board, parents, attributes) {
359     var F1 = parents[0], // focus
360         l = parents[1],  // directrix
361         rotationMatrix,
362         attr_foci = JXG.copyAttributes(attributes, board.options, 'conic', 'foci'),
363         attr_curve = JXG.copyAttributes(attributes, board.options, 'conic');
364 
365     if (parents[0].length>1) { // focus 1 given by coordinates
366         F1 = board.create('point', parents[0], attr_foci);
367     } else if (JXG.isPoint(parents[0])) { // focus i given by point
368         F1 = JXG.getReference(board,parents[0]);
369     } else if ((typeof parents[0] == 'function') && (parents[0]().elementClass == JXG.OBJECT_CLASS_POINT)) {  // given by function
370         F1 = parents[0]();
371     } else if (JXG.isString(parents[0])) { // focus i given by point name
372         F1 = JXG.getReference(board,parents[0]);
373     } else
374         throw new Error("JSXGraph: Can't create Parabola with parent types '" + 
375                         (typeof parents[0]) + "' and '" + (typeof parents[1]) + "'." +
376                         "\nPossible parent types: [point,line]");
377 
378     if (typeof parents[3]=='undefined') parents[3] = 10.0;   // to
379     if (typeof parents[2]=='undefined') parents[2] = -10.0;  // from
380 
381     var M = board.create('point', [
382                 function() {
383                     var v = [0,l.stdform[1],l.stdform[2]];
384                     v = JXG.Math.crossProduct(v,F1.coords.usrCoords);
385                     return JXG.Math.Geometry.meetLineLine(v,l.stdform,0,board).usrCoords;
386                 }
387             ], attr_foci);
388 
389     var transformFunc = function() {
390             var beta = Math.atan(l.getSlope()),                
391                 x = (M.X()+F1.X())*0.5,
392                 y = (M.Y()+F1.Y())*0.5;
393             beta += (F1.Y()-M.Y()<0 || (F1.Y()==M.Y() && F1.X()>M.X()) ) ? Math.PI : 0;
394 
395             // Rotate by the slope of the line l (Leitlinie = directrix)
396             var m = [
397                         [1,    0,             0],
398                         [x*(1-Math.cos(beta))+y*Math.sin(beta),Math.cos(beta),-Math.sin(beta)],
399                         [y*(1-Math.cos(beta))-x*Math.sin(beta),Math.sin(beta), Math.cos(beta)]
400                     ];
401             return m;
402         };
403 
404     var curve = board.create('curve',[function(x) {return 0;},function(x) {return 0;},parents[2],parents[3]], attr_curve);
405 
406     var polarForm = function(t,suspendUpdate) {
407                 var e = M.Dist(F1)*0.5,
408                     transformMat = [[1,0,0],[0,1,0],[0,0,1]],
409                     a = (M.X()+F1.X())*0.5, 
410                     b = (M.Y()+F1.Y())*0.5;
411                 
412                 if (!suspendUpdate) {
413                     rotationMatrix = transformFunc();
414                     transformMat[0][0] = rotationMatrix[0][0];
415                     transformMat[0][1] = 0.0;
416                     transformMat[0][2] = 0.0;
417                     transformMat[1][0] = a*(1-rotationMatrix[1][1])+b*rotationMatrix[1][2];
418                     transformMat[1][1] = rotationMatrix[1][1];
419                     transformMat[1][2] = rotationMatrix[2][1];
420                     transformMat[2][0] = b*(1-rotationMatrix[1][1])-a*rotationMatrix[1][2];
421                     transformMat[2][1] = rotationMatrix[1][2];
422                     transformMat[2][2] = rotationMatrix[2][2];
423                     curve.quadraticform = 
424                         JXG.Math.matMatMult(JXG.Math.transpose(transformMat),
425                         JXG.Math.matMatMult(
426                             [
427                                 [-b*4*e-a*a, a, 2*e],
428                                 [a,       -1, 0],
429                                 [2*e,      0, 0]
430                             ],
431                         transformMat)); 
432                 }
433                 return JXG.Math.matVecMult(rotationMatrix,[1,t+a,t*t/(e*4)+b]);
434         };
435     curve.X = function(phi,suspendUpdate) {return polarForm(phi,suspendUpdate)[1];};
436     curve.Y = function(phi,suspendUpdate) {return polarForm(phi,suspendUpdate)[2];};
437     curve.type = JXG.OBJECT_TYPE_CONIC;
438     return curve;
439 };
440 
441 /**
442  * 
443  * @class This element is used to provide a constructor for a generic conic section uniquely defined by five points.
444  * @pseudo
445  * @description
446  * @name Conic
447  * @augments JXG.Curve
448  * @constructor
449  * @type JXG.Conic
450  * @throws {Exception} If the element cannot be constructed with the given parent objects an exception is thrown.
451  * @param {JXG.Point,array_JXG.Point,array_JXG.Point,array_JXG.Point,array_JXG.Point,array_} point,point,point,point,point Parent elements are five points.
452  * @param {number_number_number_number_number_number} 6 numbers (a_00,a_11,a_22,a_01,a_12,a_22)
453  * @example
454  * // Create a conic section through the points A, B, C, D, and E.
455  *  var A = board.create('point', [1,5]);
456  *  var B = board.create('point', [1,2]);
457  *  var C = board.create('point', [2,0]);
458  *  var D = board.create('point', [0,0]);
459  *  var E = board.create('point', [-1,5]);
460  *  var conic = board.create('conic',[A,B,C,D,E]);
461  * </pre><div id="2d79bd6a-db9b-423c-9cba-2497f0b06320" style="width: 300px; height: 300px;"></div>
462  * <script type="text/javascript">
463  *   var glex1_board = JXG.JSXGraph.initBoard('2d79bd6a-db9b-423c-9cba-2497f0b06320', {boundingbox:[-6,6,6,-6], keepaspectratio:true, showcopyright: false, shownavigation: false});
464  *   var A = glex1_board.create('point', [1,5]);
465  *   var B = glex1_board.create('point', [1,2]);
466  *   var C = glex1_board.create('point', [2,0]);
467  *   var D = glex1_board.create('point', [0,0]);
468  *   var E = glex1_board.create('point', [-1,5]);
469  *   var conic = glex1_board.create('conic',[A,B,C,D,E]);
470  * </script><pre>
471  */
472 JXG.createConic = function(board, parents, attributes) {
473     var rotationMatrix = [[1,0,0],[0,1,0],[0,0,1]], 
474         eigen, a, b, c, M = [[1,0,0],[0,1,0],[0,0,1]],
475         c1, c2, points = [], i, definingMat, 
476         givenByPoints, 
477         p = [],
478         attr_foci = JXG.copyAttributes(attributes, board.options, 'conic', 'foci'),
479         attr_curve = JXG.copyAttributes(attributes, board.options, 'conic');
480 
481     if (parents.length==5) {
482         givenByPoints = true;
483     } else if (parents.length==6) {
484         givenByPoints = false;
485     } else 
486         throw new Error("JSXGraph: Can't create generic Conic with " + parent.length + " parameters.");  
487 
488     if (givenByPoints) {
489         for (i=0;i<5;i++) {
490             if (parents[i].length>1) { // point i given by coordinates
491                 points[i] = board.create('point', parents[i], attr_foci);
492             } else if (JXG.isPoint(parents[i])) { // point i given by point
493                 points[i] = JXG.getReference(board,parents[i]);
494             } else if ((typeof parents[i] == 'function') && (parents[i]().elementClass == JXG.OBJECT_CLASS_POINT)) {  // given by function
495                 points[i] = parents[i]();
496             } else if (JXG.isString(parents[i])) { // point i given by point name
497                 points[i] = JXG.getReference(board,parents[i]);
498             } else
499                 throw new Error("JSXGraph: Can't create Conic section with parent types '" + (typeof parents[i]) + "'." +
500                                 "\nPossible parent types: [point,point,point,point,point], [a00,a11,a22,a01,a02,a12]");
501         }
502     } else {
503         /* Usual notation (x,y,z):
504          *  [[A0,A3,A4],
505          *   [A3,A1,A5],
506          *   [A4,A5,A2]]. 
507          * Our notation (z,x,y): 
508          *  [[-A2   , A4*2.0, A5*0.5],
509          *   [A4*2.0,    -A0, A3*0.5],
510          *   [A5*0.5, A3*0.5,    -A1]] 
511          * New: (z,x,y): 
512          *  [[A2, A4, A5],
513          *   [A4, A0, A3],
514          *   [A5, A3, A1]] 
515         */
516         definingMat = [[0,0,0],[0,0,0],[0,0,0]];
517         definingMat[0][0] = (JXG.isFunction(parents[2])) ? function(){ return    parents[2]();} : function(){ return    parents[2];};
518         definingMat[0][1] = (JXG.isFunction(parents[4])) ? function(){ return    parents[4]();} : function(){ return    parents[4];};
519         definingMat[0][2] = (JXG.isFunction(parents[5])) ? function(){ return    parents[5]();} : function(){ return    parents[5];};
520         definingMat[1][1] = (JXG.isFunction(parents[0])) ? function(){ return    parents[0]();} : function(){ return    parents[0];};
521         definingMat[1][2] = (JXG.isFunction(parents[3])) ? function(){ return    parents[3]();} : function(){ return    parents[3];};
522         definingMat[2][2] = (JXG.isFunction(parents[1])) ? function(){ return    parents[1]();} : function(){ return    parents[1];};
523     }
524 
525     // sym(A) = A + A^t . Manipulates A in place.
526     var sym = function(A) {
527         var i, j;
528         for (i=0;i<3;i++) {
529             for (j=i;j<3;j++) {
530                 A[i][j] += A[j][i];
531             }
532         }
533         for (i=0;i<3;i++) {
534             for (j=0;j<i;j++) {
535                 A[i][j] = A[j][i];
536             }
537         }
538         return A;
539     };
540 
541     // degconic(v,w) = sym(v*w^t)
542     var degconic = function(v,w) {
543         var i, j, mat = [[0,0,0],[0,0,0],[0,0,0]];
544         for (i=0;i<3;i++) {
545             for (j=0;j<3;j++) {
546                 mat[i][j] = v[i]*w[j];
547             }
548         }
549         return sym(mat);
550     };
551 
552     // (p^t*B*p)*A-(p^t*A*p)*B
553     var fitConic = function(A,B,p)  {
554         var pBp, pAp, Mv, mat = [[0,0,0],[0,0,0],[0,0,0]], i, j;
555         Mv = JXG.Math.matVecMult(B,p);
556         pBp = JXG.Math.innerProduct(p,Mv);
557         Mv = JXG.Math.matVecMult(A,p);
558         pAp = JXG.Math.innerProduct(p,Mv);
559         for (i=0;i<3;i++) {
560             for (j=0;j<3;j++) {
561                 mat[i][j] = pBp*A[i][j]-pAp*B[i][j];
562             }
563         }
564         return mat;
565     };
566  
567     // Here, the defining functions for the curve are just dummy functions.
568     // In polarForm there is a reference to curve.quadraticform.
569     var curve = board.create('curve',[function(x) {return 0;},function(x) {return 0;},0,2*Math.PI], attr_curve);
570 
571     var polarForm = function(phi,suspendUpdate) {
572         var i, j, len, v;
573         if (!suspendUpdate) {
574             if (givenByPoints) {
575                 // Copy the point coordinate vectors
576                 for (i=0;i<5;i++) { 
577                     p[i] = points[i].coords.usrCoords; 
578                 }
579                 // Compute the quadratic form
580                 c1 = degconic(JXG.Math.crossProduct(p[0],p[1]),JXG.Math.crossProduct(p[2],p[3]));
581                 c2 = degconic(JXG.Math.crossProduct(p[0],p[2]),JXG.Math.crossProduct(p[1],p[3]));
582                 M = fitConic(c1,c2,p[4]);
583             } else {
584                 for (i=0;i<3;i++) {
585                     for (j=i;j<3;j++) {
586                         M[i][j] = definingMat[i][j]();
587                         if (j>i) M[j][i] = M[i][j];     
588                     }
589                 }
590             }
591             curve.quadraticform = M;    // Here is the reference back to the curve.
592             
593             // Compute Eigenvalues and Eigenvectors
594             eigen = JXG.Math.Numerics.Jacobi(M);
595             // Scale the Eigenvalues such that the first Eigenvalue is positive
596             if (eigen[0][0][0]<0) {
597                 eigen[0][0][0] *= (-1);
598                 eigen[0][1][1] *= (-1);
599                 eigen[0][2][2] *= (-1);
600             }
601             // Normalize the Eigenvectors
602             for (i=0;i<3;i++) {
603                 len = 0.0;
604                 for (j=0;j<3;j++) {
605                     len += eigen[1][j][i]*eigen[1][j][i];
606                 }
607                 len = Math.sqrt(len);
608                 for (j=0;j<3;j++) {
609                     //eigen[1][j][i] /= len;
610                 }
611             }
612             rotationMatrix = eigen[1];
613             c = Math.sqrt(Math.abs(eigen[0][0][0]));
614             a = Math.sqrt(Math.abs(eigen[0][1][1]));
615             b = Math.sqrt(Math.abs(eigen[0][2][2]));
616 
617         }
618         // The degenerate cases with eigen[0][i][i]==0 are not handled correct yet.
619         if (eigen[0][1][1]<=0.0 && eigen[0][2][2]<=0.0) {
620             v = JXG.Math.matVecMult(rotationMatrix,[1/c,Math.cos(phi)/a,Math.sin(phi)/b]);
621         } else if (eigen[0][1][1]<=0.0 && eigen[0][2][2]>0.0) {
622             v = JXG.Math.matVecMult(rotationMatrix,[Math.cos(phi)/c,1/a,Math.sin(phi)/b]);
623         } else if (eigen[0][2][2]<0.0) {
624             v = JXG.Math.matVecMult(rotationMatrix,[Math.sin(phi)/c,Math.cos(phi)/a,1/b]);
625         } 
626         // Normalize
627         v[1] /= v[0];
628         v[2] /= v[0];
629         v[0] = 1.0;
630         return v;
631     };
632 
633     curve.X = function(phi,suspendUpdate) {return polarForm(phi,suspendUpdate)[1];};
634     curve.Y = function(phi,suspendUpdate) {return polarForm(phi,suspendUpdate)[2];};
635 
636 
637     // Center coordinates see http://en.wikipedia.org/wiki/Matrix_representation_of_conic_sections
638     curve.midpoint = board.create('point',
639         [
640         function(){ 
641             var m = curve.quadraticform;
642             return [
643                 m[1][1]*m[2][2]-m[1][2]*m[1][2],
644                 m[1][2]*m[0][2]-m[2][2]*m[0][1],
645                 m[0][1]*m[1][2]-m[1][1]*m[0][2]
646             ];
647         }
648         ], attr_foci);
649 
650     curve.type = JXG.OBJECT_TYPE_CONIC;
651     return curve;
652 };
653 
654 JXG.JSXGraph.registerElement('ellipse', JXG.createEllipse);
655 JXG.JSXGraph.registerElement('hyperbola', JXG.createHyperbola);
656 JXG.JSXGraph.registerElement('parabola', JXG.createParabola);
657 JXG.JSXGraph.registerElement('conic', JXG.createConic);
658 
659