src/EllipticFunction.cpp
Go to the documentation of this file.
1 
11 
12 #if defined(_MSC_VER)
13 // Squelch warnings about constant conditional expressions
14 # pragma warning (disable: 4127)
15 #endif
16 
17 namespace GeographicLib {
18 
19  using namespace std;
20 
21  /*
22  * Implementation of methods given in
23  *
24  * B. C. Carlson
25  * Computation of elliptic integrals
26  * Numerical Algorithms 10, 13-26 (1995)
27  */
28 
30  // Carlson, eqs 2.2 - 2.7
31  static const real tolRF =
32  pow(3 * numeric_limits<real>::epsilon() * real(0.01), 1/real(8));
33  real
34  A0 = (x + y + z)/3,
35  An = A0,
36  Q = max(max(abs(A0-x), abs(A0-y)), abs(A0-z)) / tolRF,
37  x0 = x,
38  y0 = y,
39  z0 = z,
40  mul = 1;
41  while (Q >= mul * abs(An)) {
42  // Max 6 trips
43  real lam = sqrt(x0)*sqrt(y0) + sqrt(y0)*sqrt(z0) + sqrt(z0)*sqrt(x0);
44  An = (An + lam)/4;
45  x0 = (x0 + lam)/4;
46  y0 = (y0 + lam)/4;
47  z0 = (z0 + lam)/4;
48  mul *= 4;
49  }
50  real
51  X = (A0 - x) / (mul * An),
52  Y = (A0 - y) / (mul * An),
53  Z = - (X + Y),
54  E2 = X*Y - Z*Z,
55  E3 = X*Y*Z;
56  // http://dlmf.nist.gov/19.36.E1
57  // Polynomial is
58  // (1 - E2/10 + E3/14 + E2^2/24 - 3*E2*E3/44
59  // - 5*E2^3/208 + 3*E3^2/104 + E2^2*E3/16)
60  // convert to Horner form...
61  return (E3 * (6930 * E3 + E2 * (15015 * E2 - 16380) + 17160) +
62  E2 * ((10010 - 5775 * E2) * E2 - 24024) + 240240) /
63  (240240 * sqrt(An));
64  }
65 
67  // Carlson, eqs 2.36 - 2.38
68  static const real tolRG0 =
69  real(2.7) * sqrt((numeric_limits<real>::epsilon() * real(0.01)));
70  real xn = sqrt(x), yn = sqrt(y);
71  if (xn < yn) swap(xn, yn);
72  while (abs(xn-yn) > tolRG0 * xn) {
73  // Max 4 trips
74  real t = (xn + yn) /2;
75  yn = sqrt(xn * yn);
76  xn = t;
77  }
78  return Math::pi() / (xn + yn);
79  }
80 
82  // Defined only for y != 0 and x >= 0.
83  return ( !(x >= y) ? // x < y and catch nans
84  // http://dlmf.nist.gov/19.2.E18
85  atan(sqrt((y - x) / x)) / sqrt(y - x) :
86  ( x == y ? 1 / sqrt(y) :
87  Math::asinh( y > 0 ?
88  // http://dlmf.nist.gov/19.2.E19
89  // atanh(sqrt((x - y) / x))
90  sqrt((x - y) / y) :
91  // http://dlmf.nist.gov/19.2.E20
92  // atanh(sqrt(x / (x - y)))
93  sqrt(-x / y) ) / sqrt(x - y) ) );
94  }
95 
97  if (z == 0)
98  swap(y, z);
99  // Carlson, eq 1.7
100  return (z * RF(x, y, z) - (x-z) * (y-z) * RD(x, y, z) / 3
101  + sqrt(x * y / z)) / 2;
102  }
103 
105  // Carlson, eqs 2.36 - 2.39
106  static const real tolRG0 =
107  real(2.7) * sqrt((numeric_limits<real>::epsilon() * real(0.01)));
108  real
109  x0 = sqrt(max(x, y)),
110  y0 = sqrt(min(x, y)),
111  xn = x0,
112  yn = y0,
113  s = 0,
114  mul = real(0.25);
115  while (abs(xn-yn) > tolRG0 * xn) {
116  // Max 4 trips
117  real t = (xn + yn) /2;
118  yn = sqrt(xn * yn);
119  xn = t;
120  mul *= 2;
121  t = xn - yn;
122  s += mul * t * t;
123  }
124  return (Math::sq( (x0 + y0)/2 ) - s) * Math::pi() / (2 * (xn + yn));
125  }
126 
128  // Carlson, eqs 2.17 - 2.25
129  static const real
130  tolRD = pow(real(0.2) * (numeric_limits<real>::epsilon() * real(0.01)),
131  1/real(8));
132  real
133  A0 = (x + y + z + 2*p)/5,
134  An = A0,
135  delta = (p-x) * (p-y) * (p-z),
136  Q = max(max(abs(A0-x), abs(A0-y)), max(abs(A0-z), abs(A0-p))) / tolRD,
137  x0 = x,
138  y0 = y,
139  z0 = z,
140  p0 = p,
141  mul = 1,
142  mul3 = 1,
143  s = 0;
144  while (Q >= mul * abs(An)) {
145  // Max 7 trips
146  real
147  lam = sqrt(x0)*sqrt(y0) + sqrt(y0)*sqrt(z0) + sqrt(z0)*sqrt(x0),
148  d0 = (sqrt(p0)+sqrt(x0)) * (sqrt(p0)+sqrt(y0)) * (sqrt(p0)+sqrt(z0)),
149  e0 = delta/(mul3 * Math::sq(d0));
150  s += RC(1, 1 + e0)/(mul * d0);
151  An = (An + lam)/4;
152  x0 = (x0 + lam)/4;
153  y0 = (y0 + lam)/4;
154  z0 = (z0 + lam)/4;
155  p0 = (p0 + lam)/4;
156  mul *= 4;
157  mul3 *= 64;
158  }
159  real
160  X = (A0 - x) / (mul * An),
161  Y = (A0 - y) / (mul * An),
162  Z = (A0 - z) / (mul * An),
163  P = -(X + Y + Z) / 2,
164  E2 = X*Y + X*Z + Y*Z - 3*P*P,
165  E3 = X*Y*Z + 2*P * (E2 + 2*P*P),
166  E4 = (2*X*Y*Z + P * (E2 + 3*P*P)) * P,
167  E5 = X*Y*Z*P*P;
168  // http://dlmf.nist.gov/19.36.E2
169  // Polynomial is
170  // (1 - 3*E2/14 + E3/6 + 9*E2^2/88 - 3*E4/22 - 9*E2*E3/52 + 3*E5/26
171  // - E2^3/16 + 3*E3^2/40 + 3*E2*E4/20 + 45*E2^2*E3/272
172  // - 9*(E3*E4+E2*E5)/68)
173  return ((471240 - 540540 * E2) * E5 +
174  (612612 * E2 - 540540 * E3 - 556920) * E4 +
175  E3 * (306306 * E3 + E2 * (675675 * E2 - 706860) + 680680) +
176  E2 * ((417690 - 255255 * E2) * E2 - 875160) + 4084080) /
177  (4084080 * mul * An * sqrt(An)) + 6 * s;
178  }
179 
181  // Carlson, eqs 2.28 - 2.34
182  static const real
183  tolRD = pow(real(0.2) * (numeric_limits<real>::epsilon() * real(0.01)),
184  1/real(8));
185  real
186  A0 = (x + y + 3*z)/5,
187  An = A0,
188  Q = max(max(abs(A0-x), abs(A0-y)), abs(A0-z)) / tolRD,
189  x0 = x,
190  y0 = y,
191  z0 = z,
192  mul = 1,
193  s = 0;
194  while (Q >= mul * abs(An)) {
195  // Max 7 trips
196  real lam = sqrt(x0)*sqrt(y0) + sqrt(y0)*sqrt(z0) + sqrt(z0)*sqrt(x0);
197  s += 1/(mul * sqrt(z0) * (z0 + lam));
198  An = (An + lam)/4;
199  x0 = (x0 + lam)/4;
200  y0 = (y0 + lam)/4;
201  z0 = (z0 + lam)/4;
202  mul *= 4;
203  }
204  real
205  X = (A0 - x) / (mul * An),
206  Y = (A0 - y) / (mul * An),
207  Z = -(X + Y) / 3,
208  E2 = X*Y - 6*Z*Z,
209  E3 = (3*X*Y - 8*Z*Z)*Z,
210  E4 = 3 * (X*Y - Z*Z) * Z*Z,
211  E5 = X*Y*Z*Z*Z;
212  // http://dlmf.nist.gov/19.36.E2
213  // Polynomial is
214  // (1 - 3*E2/14 + E3/6 + 9*E2^2/88 - 3*E4/22 - 9*E2*E3/52 + 3*E5/26
215  // - E2^3/16 + 3*E3^2/40 + 3*E2*E4/20 + 45*E2^2*E3/272
216  // - 9*(E3*E4+E2*E5)/68)
217  return ((471240 - 540540 * E2) * E5 +
218  (612612 * E2 - 540540 * E3 - 556920) * E4 +
219  E3 * (306306 * E3 + E2 * (675675 * E2 - 706860) + 680680) +
220  E2 * ((417690 - 255255 * E2) * E2 - 875160) + 4084080) /
221  (4084080 * mul * An * sqrt(An)) + 3 * s;
222  }
223 
225  real kp2, real alphap2) {
226  // Accept nans here (needed for GeodesicExact)
227  if (k2 > 1)
228  throw GeographicErr("Parameter k2 is not in (-inf, 1]");
229  if (alpha2 > 1)
230  throw GeographicErr("Parameter alpha2 is not in (-inf, 1]");
231  if (kp2 < 0)
232  throw GeographicErr("Parameter kp2 is not in [0, inf)");
233  if (alphap2 < 0)
234  throw GeographicErr("Parameter alphap2 is not in [0, inf)");
235  _k2 = k2;
236  _kp2 = kp2;
237  _alpha2 = alpha2;
238  _alphap2 = alphap2;
239  _eps = _k2/Math::sq(sqrt(_kp2) + 1);
240  // Values of complete elliptic integrals for k = 0,1 and alpha = 0,1
241  // K E D
242  // k = 0: pi/2 pi/2 pi/4
243  // k = 1: inf 1 inf
244  // Pi G H
245  // k = 0, alpha = 0: pi/2 pi/2 pi/4
246  // k = 1, alpha = 0: inf 1 1
247  // k = 0, alpha = 1: inf inf pi/2
248  // k = 1, alpha = 1: inf inf inf
249  //
250  // Pi(0, k) = K(k)
251  // G(0, k) = E(k)
252  // H(0, k) = K(k) - D(k)
253  // Pi(0, k) = K(k)
254  // G(0, k) = E(k)
255  // H(0, k) = K(k) - D(k)
256  // Pi(alpha2, 0) = pi/(2*sqrt(1-alpha2))
257  // G(alpha2, 0) = pi/(2*sqrt(1-alpha2))
258  // H(alpha2, 0) = pi/(2*(1 + sqrt(1-alpha2)))
259  // Pi(alpha2, 1) = inf
260  // H(1, k) = K(k)
261  // G(alpha2, 1) = H(alpha2, 1) = RC(1, alphap2)
262  if (_k2 != 0) {
263  // Complete elliptic integral K(k), Carlson eq. 4.1
264  // http://dlmf.nist.gov/19.25.E1
265  _Kc = _kp2 != 0 ? RF(_kp2, 1) : Math::infinity();
266  // Complete elliptic integral E(k), Carlson eq. 4.2
267  // http://dlmf.nist.gov/19.25.E1
268  _Ec = _kp2 != 0 ? 2 * RG(_kp2, 1) : 1;
269  // D(k) = (K(k) - E(k))/k^2, Carlson eq.4.3
270  // http://dlmf.nist.gov/19.25.E1
271  _Dc = _kp2 != 0 ? RD(0, _kp2, 1) / 3 : Math::infinity();
272  } else {
273  _Kc = _Ec = Math::pi()/2; _Dc = _Kc/2;
274  }
275  if (_alpha2 != 0) {
276  // http://dlmf.nist.gov/19.25.E2
277  real rj = (_kp2 != 0 && _alphap2 != 0) ? RJ(0, _kp2, 1, _alphap2) :
278  Math::infinity(),
279  // Only use rc if _kp2 = 0.
280  rc = _kp2 != 0 ? 0 :
281  (_alphap2 != 0 ? RC(1, _alphap2) : Math::infinity());
282  // Pi(alpha^2, k)
283  _Pic = _kp2 != 0 ? _Kc + _alpha2 * rj / 3 : Math::infinity();
284  // G(alpha^2, k)
285  _Gc = _kp2 != 0 ? _Kc + (_alpha2 - _k2) * rj / 3 : rc;
286  // H(alpha^2, k)
287  _Hc = _kp2 != 0 ? _Kc - (_alphap2 != 0 ? _alphap2 * rj : 0) / 3 : rc;
288  } else {
289  _Pic = _Kc; _Gc = _Ec;
290  // Hc = Kc - Dc but this involves large cancellations if k2 is close to
291  // 1. So write (for alpha2 = 0)
292  // Hc = int(cos(phi)^2/sqrt(1-k2*sin(phi)^2),phi,0,pi/2)
293  // = 1/sqrt(1-k2) * int(sin(phi)^2/sqrt(1-k2/kp2*sin(phi)^2,...)
294  // = 1/kp * D(i*k/kp)
295  // and use D(k) = RD(0, kp2, 1) / 3
296  // so Hc = 1/kp * RD(0, 1/kp2, 1) / 3
297  // = kp2 * RD(0, 1, kp2) / 3
298  // using http://dlmf.nist.gov/19.20.E18
299  // Equivalently
300  // RF(x, 1) - RD(0, x, 1)/3 = x * RD(0, 1, x)/3 for x > 0
301  // For k2 = 1 and alpha2 = 0, we have
302  // Hc = int(cos(phi),...) = 1
303  _Hc = _kp2 != 0 ? _kp2 * RD(0, 1, _kp2) / 3 : 1;
304  }
305  }
306 
307  /*
308  * Implementation of methods given in
309  *
310  * R. Bulirsch
311  * Numerical Calculation of Elliptic Integrals and Elliptic Functions
312  * Numericshe Mathematik 7, 78-90 (1965)
313  */
314 
315  void EllipticFunction::sncndn(real x, real& sn, real& cn, real& dn) const {
316  // Bulirsch's sncndn routine, p 89.
317  static const real tolJAC =
319  if (_kp2 != 0) {
320  real mc = _kp2, d = 0;
321  if (_kp2 < 0) {
322  d = 1 - mc;
323  mc /= -d;
324  d = sqrt(d);
325  x *= d;
326  }
327  real c = 0; // To suppress warning about uninitialized variable
328  real m[num_], n[num_];
329  unsigned l = 0;
330  for (real a = 1; l < num_ || GEOGRAPHICLIB_PANIC; ++l) {
331  // This converges quadratically. Max 5 trips
332  m[l] = a;
333  n[l] = mc = sqrt(mc);
334  c = (a + mc) / 2;
335  if (!(abs(a - mc) > tolJAC * a)) {
336  ++l;
337  break;
338  }
339  mc *= a;
340  a = c;
341  }
342  x *= c;
343  sn = sin(x);
344  cn = cos(x);
345  dn = 1;
346  if (sn != 0) {
347  real a = cn / sn;
348  c *= a;
349  while (l--) {
350  real b = m[l];
351  a *= c;
352  c *= dn;
353  dn = (n[l] + a) / (b + a);
354  a = c / b;
355  }
356  a = 1 / sqrt(c*c + 1);
357  sn = sn < 0 ? -a : a;
358  cn = c * sn;
359  if (_kp2 < 0) {
360  swap(cn, dn);
361  sn /= d;
362  }
363  }
364  } else {
365  sn = tanh(x);
366  dn = cn = 1 / cosh(x);
367  }
368  }
369 
371  // Carlson, eq. 4.5 and
372  // http://dlmf.nist.gov/19.25.E5
373  real cn2 = cn*cn, dn2 = dn*dn,
374  fi = cn2 != 0 ? abs(sn) * RF(cn2, dn2, 1) : K();
375  // Enforce usual trig-like symmetries
376  if (cn < 0)
377  fi = 2 * K() - fi;
378  return Math::copysign(fi, sn);
379  }
380 
382  real
383  cn2 = cn*cn, dn2 = dn*dn, sn2 = sn*sn,
384  ei = cn2 != 0 ?
385  abs(sn) * ( _k2 <= 0 ?
386  // Carlson, eq. 4.6 and
387  // http://dlmf.nist.gov/19.25.E9
388  RF(cn2, dn2, 1) - _k2 * sn2 * RD(cn2, dn2, 1) / 3 :
389  ( _kp2 >= 0 ?
390  // http://dlmf.nist.gov/19.25.E10
391  _kp2 * RF(cn2, dn2, 1) +
392  _k2 * _kp2 * sn2 * RD(cn2, 1, dn2) / 3 +
393  _k2 * abs(cn) / dn :
394  // http://dlmf.nist.gov/19.25.E11
395  - _kp2 * sn2 * RD(dn2, 1, cn2) / 3 +
396  dn / abs(cn) ) ) :
397  E();
398  // Enforce usual trig-like symmetries
399  if (cn < 0)
400  ei = 2 * E() - ei;
401  return Math::copysign(ei, sn);
402  }
403 
405  // Carlson, eq. 4.8 and
406  // http://dlmf.nist.gov/19.25.E13
407  real
408  cn2 = cn*cn, dn2 = dn*dn, sn2 = sn*sn,
409  di = cn2 != 0 ? abs(sn) * sn2 * RD(cn2, dn2, 1) / 3 : D();
410  // Enforce usual trig-like symmetries
411  if (cn < 0)
412  di = 2 * D() - di;
413  return Math::copysign(di, sn);
414  }
415 
417  // Carlson, eq. 4.7 and
418  // http://dlmf.nist.gov/19.25.E14
419  real
420  cn2 = cn*cn, dn2 = dn*dn, sn2 = sn*sn,
421  pii = cn2 != 0 ? abs(sn) * (RF(cn2, dn2, 1) +
422  _alpha2 * sn2 *
423  RJ(cn2, dn2, 1, cn2 + _alphap2 * sn2) / 3) :
424  Pi();
425  // Enforce usual trig-like symmetries
426  if (cn < 0)
427  pii = 2 * Pi() - pii;
428  return Math::copysign(pii, sn);
429  }
430 
432  real
433  cn2 = cn*cn, dn2 = dn*dn, sn2 = sn*sn,
434  gi = cn2 != 0 ? abs(sn) * (RF(cn2, dn2, 1) +
435  (_alpha2 - _k2) * sn2 *
436  RJ(cn2, dn2, 1, cn2 + _alphap2 * sn2) / 3) :
437  G();
438  // Enforce usual trig-like symmetries
439  if (cn < 0)
440  gi = 2 * G() - gi;
441  return Math::copysign(gi, sn);
442  }
443 
445  real
446  cn2 = cn*cn, dn2 = dn*dn, sn2 = sn*sn,
447  // WARNING: large cancellation if k2 = 1, alpha2 = 0, and phi near pi/2
448  hi = cn2 != 0 ? abs(sn) * (RF(cn2, dn2, 1) -
449  _alphap2 * sn2 *
450  RJ(cn2, dn2, 1, cn2 + _alphap2 * sn2) / 3) :
451  H();
452  // Enforce usual trig-like symmetries
453  if (cn < 0)
454  hi = 2 * H() - hi;
455  return Math::copysign(hi, sn);
456  }
457 
459  // Function is periodic with period pi
460  if (cn < 0) { cn = -cn; sn = -sn; }
461  return F(sn, cn, dn) * (Math::pi()/2) / K() - atan2(sn, cn);
462  }
463 
465  // Function is periodic with period pi
466  if (cn < 0) { cn = -cn; sn = -sn; }
467  return E(sn, cn, dn) * (Math::pi()/2) / E() - atan2(sn, cn);
468  }
469 
471  // Function is periodic with period pi
472  if (cn < 0) { cn = -cn; sn = -sn; }
473  return Pi(sn, cn, dn) * (Math::pi()/2) / Pi() - atan2(sn, cn);
474  }
475 
477  // Function is periodic with period pi
478  if (cn < 0) { cn = -cn; sn = -sn; }
479  return D(sn, cn, dn) * (Math::pi()/2) / D() - atan2(sn, cn);
480  }
481 
483  // Function is periodic with period pi
484  if (cn < 0) { cn = -cn; sn = -sn; }
485  return G(sn, cn, dn) * (Math::pi()/2) / G() - atan2(sn, cn);
486  }
487 
489  // Function is periodic with period pi
490  if (cn < 0) { cn = -cn; sn = -sn; }
491  return H(sn, cn, dn) * (Math::pi()/2) / H() - atan2(sn, cn);
492  }
493 
495  real sn = sin(phi), cn = cos(phi), dn = Delta(sn, cn);
496  return abs(phi) < Math::pi() ? F(sn, cn, dn) :
497  (deltaF(sn, cn, dn) + phi) * K() / (Math::pi()/2);
498  }
499 
501  real sn = sin(phi), cn = cos(phi), dn = Delta(sn, cn);
502  return abs(phi) < Math::pi() ? E(sn, cn, dn) :
503  (deltaE(sn, cn, dn) + phi) * E() / (Math::pi()/2);
504  }
505 
507  real n = ceil(ang/360 - real(0.5));
508  ang -= 360 * n;
509  real sn, cn;
510  Math::sincosd(ang, sn, cn);
511  return E(sn, cn, Delta(sn, cn)) + 4 * E() * n;
512  }
513 
515  real sn = sin(phi), cn = cos(phi), dn = Delta(sn, cn);
516  return abs(phi) < Math::pi() ? Pi(sn, cn, dn) :
517  (deltaPi(sn, cn, dn) + phi) * Pi() / (Math::pi()/2);
518  }
519 
521  real sn = sin(phi), cn = cos(phi), dn = Delta(sn, cn);
522  return abs(phi) < Math::pi() ? D(sn, cn, dn) :
523  (deltaD(sn, cn, dn) + phi) * D() / (Math::pi()/2);
524  }
525 
527  real sn = sin(phi), cn = cos(phi), dn = Delta(sn, cn);
528  return abs(phi) < Math::pi() ? G(sn, cn, dn) :
529  (deltaG(sn, cn, dn) + phi) * G() / (Math::pi()/2);
530  }
531 
533  real sn = sin(phi), cn = cos(phi), dn = Delta(sn, cn);
534  return abs(phi) < Math::pi() ? H(sn, cn, dn) :
535  (deltaH(sn, cn, dn) + phi) * H() / (Math::pi()/2);
536  }
537 
539  static const real tolJAC =
541  real n = floor(x / (2 * _Ec) + real(0.5));
542  x -= 2 * _Ec * n; // x now in [-ec, ec)
543  // Linear approximation
544  real phi = Math::pi() * x / (2 * _Ec); // phi in [-pi/2, pi/2)
545  // First order correction
546  phi -= _eps * sin(2 * phi) / 2;
547  // For kp2 close to zero use asin(x/_Ec) or
548  // J. P. Boyd, Applied Math. and Computation 218, 7005-7013 (2012)
549  // https://doi.org/10.1016/j.amc.2011.12.021
550  for (int i = 0; i < num_ || GEOGRAPHICLIB_PANIC; ++i) {
551  real
552  sn = sin(phi),
553  cn = cos(phi),
554  dn = Delta(sn, cn),
555  err = (E(sn, cn, dn) - x)/dn;
556  phi -= err;
557  if (abs(err) < tolJAC)
558  break;
559  }
560  return n * Math::pi() + phi;
561  }
562 
564  // Function is periodic with period pi
565  if (ctau < 0) { ctau = -ctau; stau = -stau; }
566  real tau = atan2(stau, ctau);
567  return Einv( tau * E() / (Math::pi()/2) ) - tau;
568  }
569 
570 } // namespace GeographicLib
Matrix3f m
const char Y
#define max(a, b)
Definition: datatypes.h:20
static T pi()
Definition: Math.hpp:202
void sncndn(real x, real &sn, real &cn, real &dn) const
Jet< T, N > cos(const Jet< T, N > &f)
Definition: jet.h:426
float real
Definition: datatypes.h:10
Scalar * y
Key F(std::uint64_t j)
Scalar * b
Definition: benchVecAdd.cpp:17
void Reset(real k2=0, real alpha2=0)
static T infinity()
Definition: Math.hpp:867
JacobiRotation< float > G
#define min(a, b)
Definition: datatypes.h:19
double mul(const double &a, const double &b)
Math::real deltaPi(real sn, real cn, real dn) const
int n
Scalar Scalar * c
Definition: benchVecAdd.cpp:17
Math::real deltaE(real sn, real cn, real dn) const
EIGEN_DEVICE_FUNC const TanhReturnType tanh() const
Jet< T, N > sin(const Jet< T, N > &f)
Definition: jet.h:439
static Cal3_S2 K(500, 500, 0.1, 640/2, 480/2)
Definition: BFloat16.h:88
static void sincosd(T x, T &sinx, T &cosx)
Definition: Math.hpp:558
Vector3f p0
set noclip points set clip one set noclip two set bar set border lt lw set xdata set ydata set zdata set x2data set y2data set boxwidth set dummy y set format x g set format y g set format x2 g set format y2 g set format z g set angles radians set nogrid set key title set key left top Right noreverse box linetype linewidth samplen spacing width set nolabel set noarrow set nologscale set logscale x set set pointsize set encoding default set nopolar set noparametric set set set set surface set nocontour set clabel set mapping cartesian set nohidden3d set cntrparam order set cntrparam linear set cntrparam levels auto set cntrparam points set size set set xzeroaxis lt lw set x2zeroaxis lt lw set yzeroaxis lt lw set y2zeroaxis lt lw set tics in set ticslevel set tics set mxtics default set mytics default set mx2tics default set my2tics default set xtics border mirror norotate autofreq set ytics border mirror norotate autofreq set ztics border nomirror norotate autofreq set nox2tics set noy2tics set timestamp bottom norotate set rrange [*:*] noreverse nowriteback set trange [*:*] noreverse nowriteback set urange [*:*] noreverse nowriteback set vrange [*:*] noreverse nowriteback set xlabel matrix size set x2label set timefmt d m y n H
static double epsilon
Definition: testRot3.cpp:37
static T asinh(T x)
Definition: Math.hpp:311
EIGEN_DEVICE_FUNC const AtanReturnType atan() const
static real RG(real x, real y, real z)
EIGEN_DEVICE_FUNC const FloorReturnType floor() const
#define Z
Definition: icosphere.cpp:21
static const Line3 l(Rot3(), 1, 1)
int EIGEN_BLAS_FUNC() swap(int *n, RealScalar *px, int *incx, RealScalar *py, int *incy)
Definition: level1_impl.h:130
static T sq(T x)
Definition: Math.hpp:232
static real RF(real x, real y, real z)
static real RC(real x, real y)
Definition: main.h:100
Namespace for GeographicLib.
RealScalar s
Header for GeographicLib::EllipticFunction class.
Math::real deltaH(real sn, real cn, real dn) const
AnnoyingScalar atan2(const AnnoyingScalar &y, const AnnoyingScalar &x)
Math::real deltaG(real sn, real cn, real dn) const
static Symbol x0('x', 0)
static T copysign(T x, T y)
Definition: Math.hpp:751
DiscreteKey E(5, 2)
Exception handling for GeographicLib.
Definition: Constants.hpp:389
The quaternion class used to represent 3D orientations and rotations.
float * p
Math::real deltaD(real sn, real cn, real dn) const
static real RD(real x, real y, real z)
Jet< T, N > sqrt(const Jet< T, N > &f)
Definition: jet.h:418
Jet< T, N > pow(const Jet< T, N > &f, double g)
Definition: jet.h:570
Math::real deltaEinv(real stau, real ctau) const
static real RJ(real x, real y, real z, real p)
EIGEN_DEVICE_FUNC const CoshReturnType cosh() const
#define X
Definition: icosphere.cpp:20
set noclip points set clip one set noclip two set bar set border lt lw set xdata set ydata set zdata set x2data set y2data set boxwidth set dummy x
#define abs(x)
Definition: datatypes.h:17
EIGEN_DEVICE_FUNC const CeilReturnType ceil() const
Point2 t(10, 10)
#define GEOGRAPHICLIB_PANIC
Definition: Math.hpp:87
Math::real deltaF(real sn, real cn, real dn) const


gtsam
Author(s):
autogenerated on Tue Jul 4 2023 02:34:12