GNU Octave  4.0.0
A high-level interpreted language, primarily intended for numerical computations, mostly compatible with Matlab
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
ellipj.cc
Go to the documentation of this file.
1 /*
2 
3 Copyright (C) 2013-2015 Leopoldo Cerbaro <redbliss@libero.it>
4 
5 This file is part of Octave.
6 
7 Octave is free software; you can redistribute it and/or modify it
8 under the terms of the GNU General Public License as published by the
9 Free Software Foundation; either version 3 of the License, or (at your
10 option) any later version.
11 
12 Octave is distributed in the hope that it will be useful, but WITHOUT
13 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
14 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
15 for more details.
16 
17 You should have received a copy of the GNU General Public License
18 along with Octave; see the file COPYING. If not, see
19 <http://www.gnu.org/licenses/>.
20 
21 */
22 
23 #ifdef HAVE_CONFIG_H
24 #include <config.h>
25 #endif
26 
27 #include "defun.h"
28 #include "error.h"
29 #include "lo-specfun.h"
30 
31 static void
32 gripe_ellipj_arg (const char *arg)
33 {
34  error ("ellipj: expecting scalar or matrix as %s argument", arg);
35 }
36 
37 DEFUN (ellipj, args, nargout,
38  "-*- texinfo -*-\n\
39 @deftypefn {Built-in Function} {[@var{sn}, @var{cn}, @var{dn}, @var{err}] =} ellipj (@var{u}, @var{m})\n\
40 @deftypefnx {Built-in Function} {[@var{sn}, @var{cn}, @var{dn}, @var{err}] =} ellipj (@var{u}, @var{m}, @var{tol})\n\
41 Compute the Jacobi elliptic functions @var{sn}, @var{cn}, and @var{dn}\n\
42 of complex argument @var{u} and real parameter @var{m}.\n\
43 \n\
44 If @var{m} is a scalar, the results are the same size as @var{u}.\n\
45 If @var{u} is a scalar, the results are the same size as @var{m}.\n\
46 If @var{u} is a column vector and @var{m} is a row vector, the\n\
47 results are matrices with @code{length (@var{u})} rows and\n\
48 @code{length (@var{m})} columns. Otherwise, @var{u} and\n\
49 @var{m} must conform in size and the results will be the same size as the\n\
50 inputs.\n\
51 \n\
52 The value of @var{u} may be complex.\n\
53 The value of @var{m} must be 0 @leq{} @var{m} @leq{} 1.\n\
54 \n\
55 The optional input @var{tol} is currently ignored (@sc{matlab} uses this to\n\
56 allow faster, less accurate approximation).\n\
57 \n\
58 If requested, @var{err} contains the following status information\n\
59 and is the same size as the result.\n\
60 \n\
61 @enumerate 0\n\
62 @item\n\
63 Normal return.\n\
64 \n\
65 @item\n\
66 Error---no computation, algorithm termination condition not met,\n\
67 return @code{NaN}.\n\
68 @end enumerate\n\
69 \n\
70 Reference: Milton @nospell{Abramowitz} and Irene A @nospell{Stegun},\n\
71 @cite{Handbook of Mathematical Functions}, Chapter 16 (Sections 16.4, 16.13,\n\
72 and 16.15), Dover, 1965.\n\
73 \n\
74 @seealso{ellipke}\n\
75 @end deftypefn")
76 {
77  octave_value_list retval;
78 
79  int nargin = args.length ();
80 
81  if (nargin < 2 || nargin > 3)
82  {
83  print_usage ();
84  return retval;
85  }
86 
87  octave_value u_arg = args(0);
88  octave_value m_arg = args(1);
89 
90  if (m_arg.is_scalar_type ())
91  {
92  double m = args(1).double_value ();
93 
94  if (error_state)
95  {
96  gripe_ellipj_arg ("second");
97  return retval;
98  }
99 
100  if (u_arg.is_scalar_type ())
101  {
102  if (u_arg.is_real_type ())
103  {
104  // u real, m scalar
105  double u = args(0).double_value ();
106 
107  if (error_state)
108  {
109  gripe_ellipj_arg ("first");
110  return retval;
111  }
112 
113  double sn, cn, dn;
114  double err = 0;
115 
116  ellipj (u, m, sn, cn, dn, err);
117 
118  if (nargout > 3)
119  retval(3) = err;
120  retval(2) = dn;
121  retval(1) = cn;
122  retval(0) = sn;
123  }
124  else
125  {
126  // u complex, m scalar
127  Complex u = u_arg.complex_value ();
128 
129  if (error_state)
130  {
131  gripe_ellipj_arg ("first");
132  return retval;
133  }
134 
135  Complex sn, cn, dn;
136  double err = 0;
137 
138  ellipj (u, m, sn, cn, dn, err);
139 
140  if (nargout > 3)
141  retval(3) = err;
142  retval(2) = dn;
143  retval(1) = cn;
144  retval(0) = sn;
145  }
146  }
147  else
148  {
149  // u is matrix, m is scalar
150  ComplexNDArray u = u_arg.complex_array_value ();
151 
152  if (error_state)
153  {
154  gripe_ellipj_arg ("first");
155  return retval;
156  }
157 
158  dim_vector sz_u = u.dims ();
159 
160  ComplexNDArray sn (sz_u), cn (sz_u), dn (sz_u);
161  NDArray err (sz_u);
162 
163  const Complex *pu = u.data ();
164  Complex *psn = sn.fortran_vec ();
165  Complex *pcn = cn.fortran_vec ();
166  Complex *pdn = dn.fortran_vec ();
167  double *perr = err.fortran_vec ();
168  octave_idx_type nel = u.numel ();
169 
170  for (octave_idx_type i = 0; i < nel; i++)
171  ellipj (pu[i], m, psn[i], pcn[i], pdn[i], perr[i]);
172 
173  if (nargout > 3)
174  retval(3) = err;
175  retval(2) = dn;
176  retval(1) = cn;
177  retval(0) = sn;
178  }
179  }
180  else
181  {
182  NDArray m = args(1).array_value ();
183 
184  if (error_state)
185  {
186  gripe_ellipj_arg ("second");
187  return retval;
188  }
189 
190  dim_vector sz_m = m.dims ();
191 
192  if (u_arg.is_scalar_type ())
193  {
194  // u is scalar, m is array
195  if (u_arg.is_real_type ())
196  {
197  // u is real scalar, m is array
198  double u = u_arg.double_value ();
199 
200  if (error_state)
201  {
202  gripe_ellipj_arg ("first");
203  return retval;
204  }
205 
206  NDArray sn (sz_m), cn (sz_m), dn (sz_m);
207  NDArray err (sz_m);
208 
209  const double *pm = m.data ();
210  double *psn = sn.fortran_vec ();
211  double *pcn = cn.fortran_vec ();
212  double *pdn = dn.fortran_vec ();
213  double *perr = err.fortran_vec ();
214  octave_idx_type nel = m.numel ();
215 
216  for (octave_idx_type i = 0; i < nel; i++)
217  ellipj (u, pm[i], psn[i], pcn[i], pdn[i], perr[i]);
218 
219  if (nargout > 3)
220  retval(3) = err;
221  retval(2) = dn;
222  retval(1) = cn;
223  retval(0) = sn;
224  }
225  else
226  {
227  // u is complex scalar, m is array
228  Complex u = u_arg.complex_value ();
229 
230  if (error_state)
231  {
232  gripe_ellipj_arg ("first");
233  return retval;
234  }
235 
236  ComplexNDArray sn (sz_m), cn (sz_m), dn (sz_m);
237  NDArray err (sz_m);
238 
239  const double *pm = m.data ();
240  Complex *psn = sn.fortran_vec ();
241  Complex *pcn = cn.fortran_vec ();
242  Complex *pdn = dn.fortran_vec ();
243  double *perr = err.fortran_vec ();
244  octave_idx_type nel = m.numel ();
245 
246  for (octave_idx_type i = 0; i < nel; i++)
247  ellipj (u, pm[i], psn[i], pcn[i], pdn[i], perr[i]);
248 
249  if (nargout > 3)
250  retval(3) = err;
251  retval(2) = dn;
252  retval(1) = cn;
253  retval(0) = sn;
254  }
255  }
256  else
257  {
258  // u is array, m is array
259  if (u_arg.is_real_type ())
260  {
261  // u is real array, m is array
262  NDArray u = u_arg.array_value ();
263 
264  if (error_state)
265  {
266  gripe_ellipj_arg ("first");
267  return retval;
268  }
269 
270  dim_vector sz_u = u.dims ();
271 
272  if (sz_u.length () == 2 && sz_m.length () == 2
273  && sz_u(1) == 1 && sz_m(0) == 1)
274  {
275  // u is real column vector, m is row vector
276  octave_idx_type ur = sz_u(0);
277  octave_idx_type mc = sz_m(1);
278  dim_vector sz_out (ur, mc);
279 
280  NDArray sn (sz_out), cn (sz_out), dn (sz_out);
281  NDArray err (sz_out);
282 
283  const double *pu = u.data ();
284  const double *pm = m.data ();
285 
286  for (octave_idx_type j = 0; j < mc; j++)
287  for (octave_idx_type i = 0; i < ur; i++)
288  ellipj (pu[i], pm[j], sn(i,j), cn(i,j), dn(i,j),
289  err(i,j));
290 
291  if (nargout > 3)
292  retval(3) = err;
293  retval(2) = dn;
294  retval(1) = cn;
295  retval(0) = sn;
296  }
297  else if (sz_m == sz_u)
298  {
299  NDArray sn (sz_m), cn (sz_m), dn (sz_m);
300  NDArray err (sz_m);
301 
302  const double *pu = u.data ();
303  const double *pm = m.data ();
304  double *psn = sn.fortran_vec ();
305  double *pcn = cn.fortran_vec ();
306  double *pdn = dn.fortran_vec ();
307  double *perr = err.fortran_vec ();
308  octave_idx_type nel = m.numel ();
309 
310  for (octave_idx_type i = 0; i < nel; i++)
311  ellipj (pu[i], pm[i], psn[i], pcn[i], pdn[i], perr[i]);
312 
313  if (nargout > 3)
314  retval(3) = err;
315  retval(2) = dn;
316  retval(1) = cn;
317  retval(0) = sn;
318  }
319  else
320  error ("ellipj: Invalid size combination for U and M");
321  }
322  else
323  {
324  // u is complex array, m is array
325  ComplexNDArray u = u_arg.complex_array_value ();
326 
327  if (error_state)
328  {
329  gripe_ellipj_arg ("second");
330  return retval;
331  }
332 
333  dim_vector sz_u = u.dims ();
334 
335  if (sz_u.length () == 2 && sz_m.length () == 2
336  && sz_u(1) == 1 && sz_m(0) == 1)
337  {
338  // u is complex column vector, m is row vector
339  octave_idx_type ur = sz_u(0);
340  octave_idx_type mc = sz_m(1);
341  dim_vector sz_out (ur, mc);
342 
343  ComplexNDArray sn (sz_out), cn (sz_out), dn (sz_out);
344  NDArray err (sz_out);
345 
346  const Complex *pu = u.data ();
347  const double *pm = m.data ();
348 
349  for (octave_idx_type j = 0; j < mc; j++)
350  for (octave_idx_type i = 0; i < ur; i++)
351  ellipj (pu[i], pm[j], sn(i,j), cn(i,j), dn(i,j),
352  err(i,j));
353 
354  if (nargout > 3)
355  retval(3) = err;
356  retval(2) = dn;
357  retval(1) = cn;
358  retval(0) = sn;
359  }
360  else if (sz_m == sz_u)
361  {
362  ComplexNDArray sn (sz_m), cn (sz_m), dn (sz_m);
363  NDArray err (sz_m);
364 
365  const Complex *pu = u.data ();
366  const double *pm = m.data ();
367  Complex *psn = sn.fortran_vec ();
368  Complex *pcn = cn.fortran_vec ();
369  Complex *pdn = dn.fortran_vec ();
370  double *perr = err.fortran_vec ();
371  octave_idx_type nel = m.numel ();
372 
373  for (octave_idx_type i = 0; i < nel; i++)
374  ellipj (pu[i], pm[i], psn[i], pcn[i], pdn[i], perr[i]);
375 
376  if (nargout > 3)
377  retval(3) = err;
378  retval(2) = dn;
379  retval(1) = cn;
380  retval(0) = sn;
381  }
382  else
383  error ("ellipj: Invalid size combination for U and M");
384  }
385  }
386  } // m matrix
387 
388  return retval;
389 }
390 
391 /*
392 ## demos taken from inst/ellipj.m
393 
394 %!demo
395 %! N = 150;
396 %! # m = [1-logspace(0,log(eps),N-1), 1]; # m near 1
397 %! # m = [0, logspace(log(eps),0,N-1)]; # m near 0
398 %! m = linspace (0,1,N); # m equally spaced
399 %! u = linspace (-20, 20, N);
400 %! M = ones (length (u), 1) * m;
401 %! U = u' * ones (1, length (m));
402 %! [sn, cn, dn] = ellipj (U,M);
403 %!
404 %! ## Plotting
405 %! data = {sn,cn,dn};
406 %! dname = {"sn","cn","dn"};
407 %! for i=1:3
408 %! subplot (1,3,i);
409 %! data{i}(data{i} > 1) = 1;
410 %! data{i}(data{i} < -1) = -1;
411 %! image (m,u,32*data{i}+32);
412 %! title (dname{i});
413 %! endfor
414 %! colormap (hot (64));
415 
416 %!demo
417 %! N = 200;
418 %! # m = [1-logspace(0,log(eps),N-1), 1]; # m near 1
419 %! # m = [0, logspace(log(eps),0,N-1)]; # m near 0
420 %! m = linspace (0,1,N); # m equally spaced
421 %! u = linspace (0,20,5);
422 %! M = ones (length (u), 1) * m;
423 %! U = u' * ones (1, length (m));
424 %! [sn, cn, dn] = ellipj (U,M);
425 %!
426 %! ## Plotting
427 %! data = {sn,cn,dn};
428 %! dname = {"sn","cn","dn"};
429 %! for i=1:3
430 %! subplot (1,3,i);
431 %! plot (m, data{i});
432 %! title (dname{i});
433 %! grid on;
434 %! endfor
435 */
436 
437 /*
438 ## tests taken from inst/test_sncndn.m
439 
440 %!test
441 %! k = (tan(pi/8.))^2; m = k*k;
442 %! SN = [
443 %! -1. + I * 0. , -0.8392965923 + 0. * I
444 %! -1. + I * 0.2 , -0.8559363407 + 0.108250955 * I
445 %! -1. + I * 0.4 , -0.906529758 + 0.2204040232 * I
446 %! -1. + I * 0.6 , -0.9931306727 + 0.3403783409 * I
447 %! -1. + I * 0.8 , -1.119268095 + 0.4720784944 * I
448 %! -1. + I * 1. , -1.29010951 + 0.6192468708 * I
449 %! -1. + I * 1.2 , -1.512691987 + 0.7850890595 * I
450 %! -1. + I * 1.4 , -1.796200374 + 0.9714821804 * I
451 %! -1. + I * 1.6 , -2.152201882 + 1.177446413 * I
452 %! -1. + I * 1.8 , -2.594547417 + 1.396378892 * I
453 %! -1. + I * 2. , -3.138145339 + 1.611394819 * I
454 %! -0.8 + I * 0. , -0.7158157937 + 0. * I
455 %! -0.8 + I * 0.2 , -0.7301746722 + 0.1394690862 * I
456 %! -0.8 + I * 0.4 , -0.7738940898 + 0.2841710966 * I
457 %! -0.8 + I * 0.6 , -0.8489542135 + 0.4394411376 * I
458 %! -0.8 + I * 0.8 , -0.9588386397 + 0.6107824358 * I
459 %! -0.8 + I * 1. , -1.108848724 + 0.8038415767 * I
460 %! -0.8 + I * 1.2 , -1.306629972 + 1.024193359 * I
461 %! -0.8 + I * 1.4 , -1.563010199 + 1.276740951 * I
462 %! -0.8 + I * 1.6 , -1.893274688 + 1.564345558 * I
463 %! -0.8 + I * 1.8 , -2.318944084 + 1.88491973 * I
464 %! -0.8 + I * 2. , -2.869716809 + 2.225506523 * I
465 %! -0.6 + I * 0. , -0.5638287208 + 0. * I
466 %! -0.6 + I * 0.2 , -0.5752723012 + 0.1654722474 * I
467 %! -0.6 + I * 0.4 , -0.610164314 + 0.3374004736 * I
468 %! -0.6 + I * 0.6 , -0.6702507087 + 0.5224614298 * I
469 %! -0.6 + I * 0.8 , -0.7586657365 + 0.7277663879 * I
470 %! -0.6 + I * 1. , -0.8803349115 + 0.9610513652 * I
471 %! -0.6 + I * 1.2 , -1.042696526 + 1.230800819 * I
472 %! -0.6 + I * 1.4 , -1.256964505 + 1.546195843 * I
473 %! -0.6 + I * 1.6 , -1.540333527 + 1.916612621 * I
474 %! -0.6 + I * 1.8 , -1.919816065 + 2.349972151 * I
475 %! -0.6 + I * 2. , -2.438761841 + 2.848129496 * I
476 %! -0.4 + I * 0. , -0.3891382858 + 0. * I
477 %! -0.4 + I * 0.2 , -0.3971152026 + 0.1850563793 * I
478 %! -0.4 + I * 0.4 , -0.4214662882 + 0.3775700801 * I
479 %! -0.4 + I * 0.6 , -0.4635087491 + 0.5853434119 * I
480 %! -0.4 + I * 0.8 , -0.5256432877 + 0.8168992398 * I
481 %! -0.4 + I * 1. , -0.611733177 + 1.081923504 * I
482 %! -0.4 + I * 1.2 , -0.7278102331 + 1.391822501 * I
483 %! -0.4 + I * 1.4 , -0.8833807998 + 1.760456461 * I
484 %! -0.4 + I * 1.6 , -1.093891878 + 2.205107766 * I
485 %! -0.4 + I * 1.8 , -1.385545188 + 2.747638761 * I
486 %! -0.4 + I * 2. , -1.805081271 + 3.41525351 * I
487 %! -0.2 + I * 0. , -0.1986311721 + 0. * I
488 %! -0.2 + I * 0.2 , -0.2027299916 + 0.1972398665 * I
489 %! -0.2 + I * 0.4 , -0.2152524522 + 0.402598347 * I
490 %! -0.2 + I * 0.6 , -0.2369100139 + 0.6246336356 * I
491 %! -0.2 + I * 0.8 , -0.2690115146 + 0.8728455227 * I
492 %! -0.2 + I * 1. , -0.3136938773 + 1.158323088 * I
493 %! -0.2 + I * 1.2 , -0.3743615191 + 1.494672508 * I
494 %! -0.2 + I * 1.4 , -0.4565255082 + 1.899466033 * I
495 %! -0.2 + I * 1.6 , -0.5694611346 + 2.39667232 * I
496 %! -0.2 + I * 1.8 , -0.7296612675 + 3.020990664 * I
497 %! -0.2 + I * 2. , -0.9685726188 + 3.826022536 * I
498 %! 0. + I * 0. , 0. + 0. * I
499 %! 0. + I * 0.2 , 0. + 0.201376364 * I
500 %! 0. + I * 0.4 , 0. + 0.4111029248 * I
501 %! 0. + I * 0.6 , 0. + 0.6380048435 * I
502 %! 0. + I * 0.8 , 0. + 0.8919321473 * I
503 %! 0. + I * 1. , 0. + 1.184486615 * I
504 %! 0. + I * 1.2 , 0. + 1.530096023 * I
505 %! 0. + I * 1.4 , 0. + 1.947754612 * I
506 %! 0. + I * 1.6 , 0. + 2.464074356 * I
507 %! 0. + I * 1.8 , 0. + 3.119049475 * I
508 %! 0. + I * 2. , 0. + 3.97786237 * I
509 %! 0.2 + I * 0. , 0.1986311721 + 0. * I
510 %! 0.2 + I * 0.2 , 0.2027299916 + 0.1972398665 * I
511 %! 0.2 + I * 0.4 , 0.2152524522 + 0.402598347 * I
512 %! 0.2 + I * 0.6 , 0.2369100139 + 0.6246336356 * I
513 %! 0.2 + I * 0.8 , 0.2690115146 + 0.8728455227 * I
514 %! 0.2 + I * 1. , 0.3136938773 + 1.158323088 * I
515 %! 0.2 + I * 1.2 , 0.3743615191 + 1.494672508 * I
516 %! 0.2 + I * 1.4 , 0.4565255082 + 1.899466033 * I
517 %! 0.2 + I * 1.6 , 0.5694611346 + 2.39667232 * I
518 %! 0.2 + I * 1.8 , 0.7296612675 + 3.020990664 * I
519 %! 0.2 + I * 2. , 0.9685726188 + 3.826022536 * I
520 %! 0.4 + I * 0. , 0.3891382858 + 0. * I
521 %! 0.4 + I * 0.2 , 0.3971152026 + 0.1850563793 * I
522 %! 0.4 + I * 0.4 , 0.4214662882 + 0.3775700801 * I
523 %! 0.4 + I * 0.6 , 0.4635087491 + 0.5853434119 * I
524 %! 0.4 + I * 0.8 , 0.5256432877 + 0.8168992398 * I
525 %! 0.4 + I * 1. , 0.611733177 + 1.081923504 * I
526 %! 0.4 + I * 1.2 , 0.7278102331 + 1.391822501 * I
527 %! 0.4 + I * 1.4 , 0.8833807998 + 1.760456461 * I
528 %! 0.4 + I * 1.6 , 1.093891878 + 2.205107766 * I
529 %! 0.4 + I * 1.8 , 1.385545188 + 2.747638761 * I
530 %! 0.4 + I * 2. , 1.805081271 + 3.41525351 * I
531 %! 0.6 + I * 0. , 0.5638287208 + 0. * I
532 %! 0.6 + I * 0.2 , 0.5752723012 + 0.1654722474 * I
533 %! 0.6 + I * 0.4 , 0.610164314 + 0.3374004736 * I
534 %! 0.6 + I * 0.6 , 0.6702507087 + 0.5224614298 * I
535 %! 0.6 + I * 0.8 , 0.7586657365 + 0.7277663879 * I
536 %! 0.6 + I * 1. , 0.8803349115 + 0.9610513652 * I
537 %! 0.6 + I * 1.2 , 1.042696526 + 1.230800819 * I
538 %! 0.6 + I * 1.4 , 1.256964505 + 1.546195843 * I
539 %! 0.6 + I * 1.6 , 1.540333527 + 1.916612621 * I
540 %! 0.6 + I * 1.8 , 1.919816065 + 2.349972151 * I
541 %! 0.6 + I * 2. , 2.438761841 + 2.848129496 * I
542 %! 0.8 + I * 0. , 0.7158157937 + 0. * I
543 %! 0.8 + I * 0.2 , 0.7301746722 + 0.1394690862 * I
544 %! 0.8 + I * 0.4 , 0.7738940898 + 0.2841710966 * I
545 %! 0.8 + I * 0.6 , 0.8489542135 + 0.4394411376 * I
546 %! 0.8 + I * 0.8 , 0.9588386397 + 0.6107824358 * I
547 %! 0.8 + I * 1. , 1.108848724 + 0.8038415767 * I
548 %! 0.8 + I * 1.2 , 1.306629972 + 1.024193359 * I
549 %! 0.8 + I * 1.4 , 1.563010199 + 1.276740951 * I
550 %! 0.8 + I * 1.6 , 1.893274688 + 1.564345558 * I
551 %! 0.8 + I * 1.8 , 2.318944084 + 1.88491973 * I
552 %! 0.8 + I * 2. , 2.869716809 + 2.225506523 * I
553 %! 1. + I * 0. , 0.8392965923 + 0. * I
554 %! 1. + I * 0.2 , 0.8559363407 + 0.108250955 * I
555 %! 1. + I * 0.4 , 0.906529758 + 0.2204040232 * I
556 %! 1. + I * 0.6 , 0.9931306727 + 0.3403783409 * I
557 %! 1. + I * 0.8 , 1.119268095 + 0.4720784944 * I
558 %! 1. + I * 1. , 1.29010951 + 0.6192468708 * I
559 %! 1. + I * 1.2 , 1.512691987 + 0.7850890595 * I
560 %! 1. + I * 1.4 , 1.796200374 + 0.9714821804 * I
561 %! 1. + I * 1.6 , 2.152201882 + 1.177446413 * I
562 %! 1. + I * 1.8 , 2.594547417 + 1.396378892 * I
563 %! 1. + I * 2. , 3.138145339 + 1.611394819 * I
564 %! ];
565 %! CN = [
566 %! -1. + I * 0. , 0.5436738271 + 0. * I
567 %! -1. + I * 0.2 , 0.5541219664 + 0.1672121517 * I
568 %! -1. + I * 0.4 , 0.5857703552 + 0.3410940893 * I
569 %! -1. + I * 0.6 , 0.6395034233 + 0.5285979063 * I
570 %! -1. + I * 0.8 , 0.716688504 + 0.7372552987 * I
571 %! -1. + I * 1. , 0.8189576795 + 0.9755037374 * I
572 %! -1. + I * 1.2 , 0.9477661951 + 1.253049471 * I
573 %! -1. + I * 1.4 , 1.103540657 + 1.581252712 * I
574 %! -1. + I * 1.6 , 1.284098214 + 1.973449038 * I
575 %! -1. + I * 1.8 , 1.481835651 + 2.4449211 * I
576 %! -1. + I * 2. , 1.679032464 + 3.011729224 * I
577 %! -0.8 + I * 0. , 0.6982891589 + 0. * I
578 %! -0.8 + I * 0.2 , 0.71187169 + 0.1430549855 * I
579 %! -0.8 + I * 0.4 , 0.7530744458 + 0.2920273465 * I
580 %! -0.8 + I * 0.6 , 0.8232501212 + 0.4531616768 * I
581 %! -0.8 + I * 0.8 , 0.9245978896 + 0.6334016187 * I
582 %! -0.8 + I * 1. , 1.060030206 + 0.8408616109 * I
583 %! -0.8 + I * 1.2 , 1.232861756 + 1.085475913 * I
584 %! -0.8 + I * 1.4 , 1.446126965 + 1.379933558 * I
585 %! -0.8 + I * 1.6 , 1.701139468 + 1.741030588 * I
586 %! -0.8 + I * 1.8 , 1.994526268 + 2.191509596 * I
587 %! -0.8 + I * 2. , 2.312257188 + 2.762051518 * I
588 %! -0.6 + I * 0. , 0.8258917445 + 0. * I
589 %! -0.6 + I * 0.2 , 0.842151698 + 0.1130337928 * I
590 %! -0.6 + I * 0.4 , 0.8915487431 + 0.2309124769 * I
591 %! -0.6 + I * 0.6 , 0.975948103 + 0.3588102098 * I
592 %! -0.6 + I * 0.8 , 1.098499209 + 0.5026234141 * I
593 %! -0.6 + I * 1. , 1.263676101 + 0.6695125973 * I
594 %! -0.6 + I * 1.2 , 1.477275851 + 0.8687285705 * I
595 %! -0.6 + I * 1.4 , 1.746262523 + 1.112955966 * I
596 %! -0.6 + I * 1.6 , 2.078179075 + 1.420581466 * I
597 %! -0.6 + I * 1.8 , 2.479425208 + 1.819580713 * I
598 %! -0.6 + I * 2. , 2.950586798 + 2.354077344 * I
599 %! -0.4 + I * 0. , 0.9211793498 + 0. * I
600 %! -0.4 + I * 0.2 , 0.9395019377 + 0.07822091534 * I
601 %! -0.4 + I * 0.4 , 0.9952345231 + 0.1598950363 * I
602 %! -0.4 + I * 0.6 , 1.090715991 + 0.2487465067 * I
603 %! -0.4 + I * 0.8 , 1.229998843 + 0.34910407 * I
604 %! -0.4 + I * 1. , 1.419103868 + 0.4663848201 * I
605 %! -0.4 + I * 1.2 , 1.666426377 + 0.607877235 * I
606 %! -0.4 + I * 1.4 , 1.983347336 + 0.7841054404 * I
607 %! -0.4 + I * 1.6 , 2.385101684 + 1.01134031 * I
608 %! -0.4 + I * 1.8 , 2.89185416 + 1.316448705 * I
609 %! -0.4 + I * 2. , 3.529393374 + 1.74670531 * I
610 %! -0.2 + I * 0. , 0.9800743122 + 0. * I
611 %! -0.2 + I * 0.2 , 0.9997019476 + 0.03999835809 * I
612 %! -0.2 + I * 0.4 , 1.059453907 + 0.08179712295 * I
613 %! -0.2 + I * 0.6 , 1.16200643 + 0.1273503824 * I
614 %! -0.2 + I * 0.8 , 1.312066413 + 0.1789585449 * I
615 %! -0.2 + I * 1. , 1.516804331 + 0.2395555269 * I
616 %! -0.2 + I * 1.2 , 1.786613221 + 0.313189147 * I
617 %! -0.2 + I * 1.4 , 2.136422971 + 0.405890925 * I
618 %! -0.2 + I * 1.6 , 2.588021972 + 0.527357091 * I
619 %! -0.2 + I * 1.8 , 3.174302819 + 0.6944201617 * I
620 %! -0.2 + I * 2. , 3.947361147 + 0.9387994989 * I
621 %! 0. + I * 0. , 1. + 0. * I
622 %! 0. + I * 0.2 , 1.020074723 + 0. * I
623 %! 0. + I * 0.4 , 1.08120563 + 0. * I
624 %! 0. + I * 0.6 , 1.18619146 + 0. * I
625 %! 0. + I * 0.8 , 1.339978715 + 0. * I
626 %! 0. + I * 1. , 1.550164037 + 0. * I
627 %! 0. + I * 1.2 , 1.827893279 + 0. * I
628 %! 0. + I * 1.4 , 2.189462954 + 0. * I
629 %! 0. + I * 1.6 , 2.659259752 + 0. * I
630 %! 0. + I * 1.8 , 3.275434266 + 0. * I
631 %! 0. + I * 2. , 4.101632484 + 0. * I
632 %! 0.2 + I * 0. , 0.9800743122 + 0. * I
633 %! 0.2 + I * 0.2 , 0.9997019476 - 0.03999835809 * I
634 %! 0.2 + I * 0.4 , 1.059453907 - 0.08179712295 * I
635 %! 0.2 + I * 0.6 , 1.16200643 - 0.1273503824 * I
636 %! 0.2 + I * 0.8 , 1.312066413 - 0.1789585449 * I
637 %! 0.2 + I * 1. , 1.516804331 - 0.2395555269 * I
638 %! 0.2 + I * 1.2 , 1.786613221 - 0.313189147 * I
639 %! 0.2 + I * 1.4 , 2.136422971 - 0.405890925 * I
640 %! 0.2 + I * 1.6 , 2.588021972 - 0.527357091 * I
641 %! 0.2 + I * 1.8 , 3.174302819 - 0.6944201617 * I
642 %! 0.2 + I * 2. , 3.947361147 - 0.9387994989 * I
643 %! 0.4 + I * 0. , 0.9211793498 + 0. * I
644 %! 0.4 + I * 0.2 , 0.9395019377 - 0.07822091534 * I
645 %! 0.4 + I * 0.4 , 0.9952345231 - 0.1598950363 * I
646 %! 0.4 + I * 0.6 , 1.090715991 - 0.2487465067 * I
647 %! 0.4 + I * 0.8 , 1.229998843 - 0.34910407 * I
648 %! 0.4 + I * 1. , 1.419103868 - 0.4663848201 * I
649 %! 0.4 + I * 1.2 , 1.666426377 - 0.607877235 * I
650 %! 0.4 + I * 1.4 , 1.983347336 - 0.7841054404 * I
651 %! 0.4 + I * 1.6 , 2.385101684 - 1.01134031 * I
652 %! 0.4 + I * 1.8 , 2.89185416 - 1.316448705 * I
653 %! 0.4 + I * 2. , 3.529393374 - 1.74670531 * I
654 %! 0.6 + I * 0. , 0.8258917445 + 0. * I
655 %! 0.6 + I * 0.2 , 0.842151698 - 0.1130337928 * I
656 %! 0.6 + I * 0.4 , 0.8915487431 - 0.2309124769 * I
657 %! 0.6 + I * 0.6 , 0.975948103 - 0.3588102098 * I
658 %! 0.6 + I * 0.8 , 1.098499209 - 0.5026234141 * I
659 %! 0.6 + I * 1. , 1.263676101 - 0.6695125973 * I
660 %! 0.6 + I * 1.2 , 1.477275851 - 0.8687285705 * I
661 %! 0.6 + I * 1.4 , 1.746262523 - 1.112955966 * I
662 %! 0.6 + I * 1.6 , 2.078179075 - 1.420581466 * I
663 %! 0.6 + I * 1.8 , 2.479425208 - 1.819580713 * I
664 %! 0.6 + I * 2. , 2.950586798 - 2.354077344 * I
665 %! 0.8 + I * 0. , 0.6982891589 + 0. * I
666 %! 0.8 + I * 0.2 , 0.71187169 - 0.1430549855 * I
667 %! 0.8 + I * 0.4 , 0.7530744458 - 0.2920273465 * I
668 %! 0.8 + I * 0.6 , 0.8232501212 - 0.4531616768 * I
669 %! 0.8 + I * 0.8 , 0.9245978896 - 0.6334016187 * I
670 %! 0.8 + I * 1. , 1.060030206 - 0.8408616109 * I
671 %! 0.8 + I * 1.2 , 1.232861756 - 1.085475913 * I
672 %! 0.8 + I * 1.4 , 1.446126965 - 1.379933558 * I
673 %! 0.8 + I * 1.6 , 1.701139468 - 1.741030588 * I
674 %! 0.8 + I * 1.8 , 1.994526268 - 2.191509596 * I
675 %! 0.8 + I * 2. , 2.312257188 - 2.762051518 * I
676 %! 1. + I * 0. , 0.5436738271 + 0. * I
677 %! 1. + I * 0.2 , 0.5541219664 - 0.1672121517 * I
678 %! 1. + I * 0.4 , 0.5857703552 - 0.3410940893 * I
679 %! 1. + I * 0.6 , 0.6395034233 - 0.5285979063 * I
680 %! 1. + I * 0.8 , 0.716688504 - 0.7372552987 * I
681 %! 1. + I * 1. , 0.8189576795 - 0.9755037374 * I
682 %! 1. + I * 1.2 , 0.9477661951 - 1.253049471 * I
683 %! 1. + I * 1.4 , 1.103540657 - 1.581252712 * I
684 %! 1. + I * 1.6 , 1.284098214 - 1.973449038 * I
685 %! 1. + I * 1.8 , 1.481835651 - 2.4449211 * I
686 %! 1. + I * 2. , 1.679032464 - 3.011729224 * I
687 %! ];
688 %! DN = [
689 %! -1. + I * 0. , 0.9895776106 + 0. * I
690 %! -1. + I * 0.2 , 0.9893361555 + 0.002756935338 * I
691 %! -1. + I * 0.4 , 0.9885716856 + 0.005949639805 * I
692 %! -1. + I * 0.6 , 0.9871564855 + 0.01008044183 * I
693 %! -1. + I * 0.8 , 0.9848512162 + 0.01579337596 * I
694 %! -1. + I * 1. , 0.9812582484 + 0.02396648455 * I
695 %! -1. + I * 1.2 , 0.9757399152 + 0.0358288294 * I
696 %! -1. + I * 1.4 , 0.9672786056 + 0.0531049859 * I
697 %! -1. + I * 1.6 , 0.954237868 + 0.0781744383 * I
698 %! -1. + I * 1.8 , 0.933957524 + 0.1141918269 * I
699 %! -1. + I * 2. , 0.9020917489 + 0.1650142936 * I
700 %! -0.8 + I * 0. , 0.992429635 + 0. * I
701 %! -0.8 + I * 0.2 , 0.9924147861 + 0.003020708044 * I
702 %! -0.8 + I * 0.4 , 0.99236555 + 0.00652359532 * I
703 %! -0.8 + I * 0.6 , 0.9922655715 + 0.0110676219 * I
704 %! -0.8 + I * 0.8 , 0.9920785856 + 0.01737733806 * I
705 %! -0.8 + I * 1. , 0.9917291795 + 0.02645738598 * I
706 %! -0.8 + I * 1.2 , 0.9910606387 + 0.03974949378 * I
707 %! -0.8 + I * 1.4 , 0.9897435004 + 0.05935252515 * I
708 %! -0.8 + I * 1.6 , 0.987077644 + 0.08832675281 * I
709 %! -0.8 + I * 1.8 , 0.9815667458 + 0.1310872821 * I
710 %! -0.8 + I * 2. , 0.970020127 + 0.1938136793 * I
711 %! -0.6 + I * 0. , 0.9953099088 + 0. * I
712 %! -0.6 + I * 0.2 , 0.995526009 + 0.002814772354 * I
713 %! -0.6 + I * 0.4 , 0.9962071136 + 0.006083312292 * I
714 %! -0.6 + I * 0.6 , 0.9974557125 + 0.01033463525 * I
715 %! -0.6 + I * 0.8 , 0.9994560563 + 0.01626207722 * I
716 %! -0.6 + I * 1. , 1.00249312 + 0.02484336286 * I
717 %! -0.6 + I * 1.2 , 1.006973922 + 0.0375167093 * I
718 %! -0.6 + I * 1.4 , 1.013436509 + 0.05645315628 * I
719 %! -0.6 + I * 1.6 , 1.022504295 + 0.08499262247 * I
720 %! -0.6 + I * 1.8 , 1.034670023 + 0.1283564595 * I
721 %! -0.6 + I * 2. , 1.049599899 + 0.194806122 * I
722 %! -0.4 + I * 0. , 0.9977686897 + 0. * I
723 %! -0.4 + I * 0.2 , 0.9981836165 + 0.002167241934 * I
724 %! -0.4 + I * 0.4 , 0.9994946045 + 0.004686808612 * I
725 %! -0.4 + I * 0.6 , 1.001910789 + 0.00797144174 * I
726 %! -0.4 + I * 0.8 , 1.005817375 + 0.01256717724 * I
727 %! -0.4 + I * 1. , 1.011836374 + 0.01925509038 * I
728 %! -0.4 + I * 1.2 , 1.020923572 + 0.02920828367 * I
729 %! -0.4 + I * 1.4 , 1.034513743 + 0.04425213602 * I
730 %! -0.4 + I * 1.6 , 1.054725746 + 0.06732276244 * I
731 %! -0.4 + I * 1.8 , 1.08462027 + 0.1033236812 * I
732 %! -0.4 + I * 2. , 1.128407402 + 0.1608240664 * I
733 %! -0.2 + I * 0. , 0.9994191176 + 0. * I
734 %! -0.2 + I * 0.2 , 0.9999683719 + 0.001177128019 * I
735 %! -0.2 + I * 0.4 , 1.001705496 + 0.00254669712 * I
736 %! -0.2 + I * 0.6 , 1.004913944 + 0.004334880912 * I
737 %! -0.2 + I * 0.8 , 1.010120575 + 0.006842775622 * I
738 %! -0.2 + I * 1. , 1.018189543 + 0.01050520136 * I
739 %! -0.2 + I * 1.2 , 1.030482479 + 0.01598431001 * I
740 %! -0.2 + I * 1.4 , 1.049126108 + 0.02433134655 * I
741 %! -0.2 + I * 1.6 , 1.077466003 + 0.0372877718 * I
742 %! -0.2 + I * 1.8 , 1.120863308 + 0.05789156398 * I
743 %! -0.2 + I * 2. , 1.188162088 + 0.09181238708 * I
744 %! 0. + I * 0. , 1. + 0. * I
745 %! 0. + I * 0.2 , 1.000596698 + 0. * I
746 %! 0. + I * 0.4 , 1.002484444 + 0. * I
747 %! 0. + I * 0.6 , 1.005973379 + 0. * I
748 %! 0. + I * 0.8 , 1.011641536 + 0. * I
749 %! 0. + I * 1. , 1.020441432 + 0. * I
750 %! 0. + I * 1.2 , 1.033885057 + 0. * I
751 %! 0. + I * 1.4 , 1.054361188 + 0. * I
752 %! 0. + I * 1.6 , 1.085694733 + 0. * I
753 %! 0. + I * 1.8 , 1.134186672 + 0. * I
754 %! 0. + I * 2. , 1.210701071 + 0. * I
755 %! 0.2 + I * 0. , 0.9994191176 + 0. * I
756 %! 0.2 + I * 0.2 , 0.9999683719 - 0.001177128019 * I
757 %! 0.2 + I * 0.4 , 1.001705496 - 0.00254669712 * I
758 %! 0.2 + I * 0.6 , 1.004913944 - 0.004334880912 * I
759 %! 0.2 + I * 0.8 , 1.010120575 - 0.006842775622 * I
760 %! 0.2 + I * 1. , 1.018189543 - 0.01050520136 * I
761 %! 0.2 + I * 1.2 , 1.030482479 - 0.01598431001 * I
762 %! 0.2 + I * 1.4 , 1.049126108 - 0.02433134655 * I
763 %! 0.2 + I * 1.6 , 1.077466003 - 0.0372877718 * I
764 %! 0.2 + I * 1.8 , 1.120863308 - 0.05789156398 * I
765 %! 0.2 + I * 2. , 1.188162088 - 0.09181238708 * I
766 %! 0.4 + I * 0. , 0.9977686897 + 0. * I
767 %! 0.4 + I * 0.2 , 0.9981836165 - 0.002167241934 * I
768 %! 0.4 + I * 0.4 , 0.9994946045 - 0.004686808612 * I
769 %! 0.4 + I * 0.6 , 1.001910789 - 0.00797144174 * I
770 %! 0.4 + I * 0.8 , 1.005817375 - 0.01256717724 * I
771 %! 0.4 + I * 1. , 1.011836374 - 0.01925509038 * I
772 %! 0.4 + I * 1.2 , 1.020923572 - 0.02920828367 * I
773 %! 0.4 + I * 1.4 , 1.034513743 - 0.04425213602 * I
774 %! 0.4 + I * 1.6 , 1.054725746 - 0.06732276244 * I
775 %! 0.4 + I * 1.8 , 1.08462027 - 0.1033236812 * I
776 %! 0.4 + I * 2. , 1.128407402 - 0.1608240664 * I
777 %! 0.6 + I * 0. , 0.9953099088 + 0. * I
778 %! 0.6 + I * 0.2 , 0.995526009 - 0.002814772354 * I
779 %! 0.6 + I * 0.4 , 0.9962071136 - 0.006083312292 * I
780 %! 0.6 + I * 0.6 , 0.9974557125 - 0.01033463525 * I
781 %! 0.6 + I * 0.8 , 0.9994560563 - 0.01626207722 * I
782 %! 0.6 + I * 1. , 1.00249312 - 0.02484336286 * I
783 %! 0.6 + I * 1.2 , 1.006973922 - 0.0375167093 * I
784 %! 0.6 + I * 1.4 , 1.013436509 - 0.05645315628 * I
785 %! 0.6 + I * 1.6 , 1.022504295 - 0.08499262247 * I
786 %! 0.6 + I * 1.8 , 1.034670023 - 0.1283564595 * I
787 %! 0.6 + I * 2. , 1.049599899 - 0.194806122 * I
788 %! 0.8 + I * 0. , 0.992429635 + 0. * I
789 %! 0.8 + I * 0.2 , 0.9924147861 - 0.003020708044 * I
790 %! 0.8 + I * 0.4 , 0.99236555 - 0.00652359532 * I
791 %! 0.8 + I * 0.6 , 0.9922655715 - 0.0110676219 * I
792 %! 0.8 + I * 0.8 , 0.9920785856 - 0.01737733806 * I
793 %! 0.8 + I * 1. , 0.9917291795 - 0.02645738598 * I
794 %! 0.8 + I * 1.2 , 0.9910606387 - 0.03974949378 * I
795 %! 0.8 + I * 1.4 , 0.9897435004 - 0.05935252515 * I
796 %! 0.8 + I * 1.6 , 0.987077644 - 0.08832675281 * I
797 %! 0.8 + I * 1.8 , 0.9815667458 - 0.1310872821 * I
798 %! 0.8 + I * 2. , 0.970020127 - 0.1938136793 * I
799 %! 1. + I * 0. , 0.9895776106 + 0. * I
800 %! 1. + I * 0.2 , 0.9893361555 - 0.002756935338 * I
801 %! 1. + I * 0.4 , 0.9885716856 - 0.005949639805 * I
802 %! 1. + I * 0.6 , 0.9871564855 - 0.01008044183 * I
803 %! 1. + I * 0.8 , 0.9848512162 - 0.01579337596 * I
804 %! 1. + I * 1. , 0.9812582484 - 0.02396648455 * I
805 %! 1. + I * 1.2 , 0.9757399152 - 0.0358288294 * I
806 %! 1. + I * 1.4 , 0.9672786056 - 0.0531049859 * I
807 %! 1. + I * 1.6 , 0.954237868 - 0.0781744383 * I
808 %! 1. + I * 1.8 , 0.933957524 - 0.1141918269 * I
809 %! 1. + I * 2. , 0.9020917489 - 0.1650142936 * I
810 %! ];
811 %! tol = 1e-9;
812 %! for x = 0:10
813 %! for y = 0:10
814 %! ur = -1 + x * 0.2;
815 %! ui = y * 0.2;
816 %! ii = 1 + y + x*11;
817 %! [sn, cn, dn] = ellipj (ur + I * ui, m);
818 %! assert (sn, SN(ii, 2), tol);
819 %! assert (cn, CN(ii, 2), tol);
820 %! assert (dn, DN(ii, 2), tol);
821 %! endfor
822 %! endfor
823 
824 ## tests taken from test_ellipj.m
825 %!test
826 %! u1 = pi/3; m1 = 0;
827 %! res1 = [sin(pi/3), cos(pi/3), 1];
828 %! [sn,cn,dn] = ellipj (u1,m1);
829 %! assert ([sn,cn,dn], res1, 10*eps);
830 
831 %!test
832 %! u2 = log(2); m2 = 1;
833 %! res2 = [ 3/5, 4/5, 4/5 ];
834 %! [sn,cn,dn] = ellipj (u2,m2);
835 %! assert ([sn,cn,dn], res2, 10*eps);
836 
837 %!test
838 %! u3 = log(2)*1i; m3 = 0;
839 %! res3 = [3i/4,5/4,1];
840 %! [sn,cn,dn] = ellipj (u3,m3);
841 %! assert ([sn,cn,dn], res3, 10*eps);
842 
843 %!test
844 %! u4 = -1; m4 = tan (pi/8)^4;
845 %! res4 = [-0.8392965923,0.5436738271,0.9895776106];
846 %! [sn,cn,dn] = ellipj (u4, m4);
847 %! assert ([sn,cn,dn], res4, 1e-10);
848 
849 %!test
850 %! u5 = -0.2 + 0.4i; m5 = tan(pi/8)^4;
851 %! res5 = [ -0.2152524522 + 0.402598347i, ...
852 %! 1.059453907 + 0.08179712295i, ...
853 %! 1.001705496 + 0.00254669712i ];
854 %! [sn,cn,dn] = ellipj (u5,m5);
855 %! assert ([sn,cn,dn], res5, 1e-9);
856 
857 %!test
858 %! u6 = 0.2 + 0.6i; m6 = tan(pi/8)^4;
859 %! res6 = [ 0.2369100139 + 0.624633635i, ...
860 %! 1.16200643 - 0.1273503824i, ...
861 %! 1.004913944 - 0.004334880912i ];
862 %! [sn,cn,dn] = ellipj (u6,m6);
863 %! assert ([sn,cn,dn], res6, 1e-8);
864 
865 %!test
866 %! u7 = 0.8 + 0.8i; m7 = tan (pi/8)^4;
867 %! res7 = [0.9588386397 + 0.6107824358i, ...
868 %! 0.9245978896 - 0.6334016187i, ...
869 %! 0.9920785856 - 0.01737733806i ];
870 %! [sn,cn,dn] = ellipj (u7,m7);
871 %! assert ([sn,cn,dn], res7, 1e-10);
872 
873 %!test
874 %! u = [0,pi/6,pi/4,pi/2]; m=0;
875 %! res = [0,1/2,1/sqrt(2),1;1,cos(pi/6),1/sqrt(2),0;1,1,1,1];
876 %! [sn,cn,dn] = ellipj (u,m);
877 %! assert ([sn;cn;dn], res, 100*eps);
878 %! [sn,cn,dn] = ellipj (u',0);
879 %! assert ([sn,cn,dn], res', 100*eps);
880 
881 ## FIXME: need to check [real,complex]x[scalar,rowvec,colvec,matrix]x[u,m]
882 
883 ## One test for u column vector x m row vector
884 %!test
885 %! u = [0,pi/6,pi/4,pi/2]'; m = [0 0 0 0];
886 %! res = [0,1/2,1/sqrt(2),1;1,cos(pi/6),1/sqrt(2),0;1,1,1,1]';
887 %! [sn,cn,dn] = ellipj (u,m);
888 %! assert (sn, repmat (res(:,1), [1,4]), 100*eps);
889 %! assert (cn, repmat (res(:,2), [1,4]), 100*eps);
890 %! assert (dn, repmat (res(:,3), [1,4]), 100*eps);
891 
892 %!test
893 %! ## Test Jacobi elliptic functions
894 %! ## against "exact" solution from Mathematica 3.0
895 %! ## David Billinghurst <David.Billinghurst@riotinto.com>
896 %! ## 1 February 2001
897 %! u = [ 0.25; 0.25; 0.20; 0.20; 0.672; 0.5];
898 %! m = [ 0.0; 1.0; 0.19; 0.81; 0.36; 0.9999999999];
899 %! S = [ sin(0.25);
900 %! tanh(0.25);
901 %! 0.19842311013970879516;
902 %! 0.19762082367187648571;
903 %! 0.6095196917919021945;
904 %! 0.4621171572617320908 ];
905 %! C = [ cos(0.25);
906 %! sech(0.25);
907 %! 0.9801164570409401062;
908 %! 0.9802785369736752032;
909 %! 0.7927709286533560550;
910 %! 0.8868188839691764094 ];
911 %! D = [ 1.0;
912 %! sech(0.25);
913 %! 0.9962526643271134302;
914 %! 0.9840560289645665155;
915 %! 0.9307281387786906491;
916 %! 0.8868188839812167635 ];
917 %! [sn,cn,dn] = ellipj (u,m);
918 %! assert (sn, S, 8*eps);
919 %! assert (cn, C, 8*eps);
920 %! assert (dn, D, 8*eps);
921 
922 %!test
923 %! ## Test continuity of dn when cn is near zero (bug #43344)
924 %! m = 0.5;
925 %! u = ellipke (0.5);
926 %! x = [-1e-3, -1e-12, 0, 1e-12, 1e-3];
927 %! [~, ~, dn] = ellipj (u + x, m);
928 %! D = 1/sqrt (2) * ones (size (x));
929 %! assert (dn, D, 1e-6);
930 
931 %!error ellipj ()
932 %!error ellipj (1)
933 %!error ellipj (1,2,3,4)
934 %!warning <expecting 0 <= M <= 1> ellipj (1,2);
935 ## FIXME: errors commented out untill lasterr() truly returns the last error.
936 %!#error <expecting scalar or matrix as second argument> ellipj (1, "1")
937 %!#error <expecting scalar or matrix as first argument> ellipj ("1", 1)
938 %!#error <expecting scalar or matrix as first argument> ellipj ({1}, 1)
939 %!#error <expecting scalar or matrix as first argument> ellipj ({1, 2}, 1)
940 %!#error <expecting scalar or matrix as second argument> ellipj (1, {1, 2})
941 %!#error <expecting scalar or matrix as first argument> ellipj ("1", [1, 2])
942 %!#error <expecting scalar or matrix as first argument> ellipj ({1}, [1, 2])
943 %!#error <expecting scalar or matrix as first argument> ellipj ({1}, [1, 2])
944 %!#error <expecting scalar or matrix as first argument> ellipj ("1,2", [1, 2])
945 %!#error <expecting scalar or matrix as first argument> ellipj ({1, 2}, [1, 2])
946 %!error <Invalid size combination for U and M> ellipj ([1:4], [1:3])
947 %!error <Invalid size combination for U and M> ellipj (complex (1:4,1:4), [1:3])
948 
949 */
950 
ComplexNDArray complex_array_value(bool frc_str_conv=false) const
Definition: ov.h:798
bool is_real_type(void) const
Definition: ov.h:651
OCTINTERP_API void print_usage(void)
Definition: defun.cc:51
octave_idx_type numel(void) const
Number of elements in the array.
Definition: Array.h:275
octave_idx_type length(void) const
Definition: oct-obj.h:89
bool is_scalar_type(void) const
Definition: ov.h:657
#define DEFUN(name, args_name, nargout_name, doc)
Definition: defun.h:44
void error(const char *fmt,...)
Definition: error.cc:476
const dim_vector & dims(void) const
Return a const-reference so that dims ()(i) works efficiently.
Definition: Array.h:337
const T * data(void) const
Definition: Array.h:479
int error_state
Definition: error.cc:101
static void gripe_ellipj_arg(const char *arg)
Definition: ellipj.cc:32
double arg(double x)
Definition: lo-mappers.h:37
NDArray array_value(bool frc_str_conv=false) const
Definition: ov.h:779
Complex complex_value(bool frc_str_conv=false) const
Definition: ov.h:785
void ellipj(double u, double m, double &sn, double &cn, double &dn, double &err)
Definition: lo-specfun.cc:3630
std::complex< double > Complex
Definition: oct-cmplx.h:29
const T * fortran_vec(void) const
Definition: Array.h:481
double double_value(bool frc_str_conv=false) const
Definition: ov.h:759
int length(void) const
Definition: dim-vector.h:281