diff --git a/CHANGES b/CHANGES index 1245230..2178851 100644 --- a/CHANGES +++ b/CHANGES @@ -99,15 +99,6 @@ The following are the changes from calc version 2.14.3.5 to date: Document in help files for builtin functions that take eps arguments, the LIMIT range for such eps values. - Added new versin(x, [,eps]) for versed sine and coversin(x, [,eps]) - for inverse versed sine. - - Added new aversin(x, [,eps]) for inverse versed sine and acoversin(x, [,eps]) - for inverse coversed sine. - - Improved trig function help files to reference use of complex arguments - that while supported were not documented. - Removed old Makefile testing rules for make dbx and make gdb. Improved "make run" to execute calccalc using shared libraries @@ -176,14 +167,15 @@ The following are the changes from calc version 2.14.3.5 to date: number generator in place of the additive 55 generator for a while now. - Improved help files for sin, cos, tan, cot, sec, csc. In case - of tan, cot, sec, csc corrected help file was corrected to - indicate that complex arguments are allowed. This was a help - file oversight from long ago when those trigonometric functions - were expanded to include complex arguments. + Improved help files trigonometric functions. They were corrected + to indicate that complex arguments are allowed: an oversight + from long ago when those trigonometric functions were expanded + to include complex arguments. The EXAMPLE sections were expanded + and made consistent, where applicable, across the trigonometric + help files. Documented libcalc functions in the SEE ALSO sections. Expanded the calc regression test suite test 34dd to test various - real and complex values for sin, cos, tan, cot, sec, csc. + real and complex values for trigonometric functions. Added complex multiple approximation function to commath.c so that users of libcalc may directly round complex number to @@ -338,10 +330,36 @@ The following are the changes from calc version 2.14.3.5 to date: verify that the errnum calc computation error codes and their E_STRING values have not changed. - Improve the ciarify of calc regression suite (regress.cal) to mostly + Improve the clarity of calc regression suite (regress.cal) to mostly use E_STRING errsym instead of numeric errnum values for error() and errno() related tests. + Fixed SEE ALSO typo in help randperm. + + Added the following new trigonometric functions: + + versin(x [,eps]) versed trigonometric sine + coversin(x [,eps]) coversed trigonometric sine + vercos(x [,eps]) versed trigonometric cosine + covercos(x [,eps]) coversed trigonometric cosine + aversin(x [,eps]) inverse versed trigonometric sine + acoversin(x [,eps]) inverse coversed trigonometric sine + avercos(x [,eps]) inverse versed trigonometric cosine + acovercos(x [,eps]) inverse coversed trigonometric cosine + haversin(x [,eps]) half versed trigonometric sine + hacoversin(x [,eps]) half coversed trigonometric sine + havercos(x [,eps]) half versed trigonometric cosine + hacovercos(x [,eps]) half coversed trigonometric cosine + ahaversin(x [,eps]) inverse half versed trigonometric sine + ahacoversin(x [,eps]) inverse half coversed trigonometric sine + ahavercos(x [,eps]) inverse half versed trigonometric cosine + ahacovercos(x [,eps]) inverse half coversed trigonometric cosine + + Fixed calc regression test 42dd to set the display value back to 20. + + Added test 95dd and test9500.trigeq.cal to the calc regression test + suite to perform extensive test of trigonometric functions. + The following are the changes from calc version 2.14.3.4 to 2.14.3.5: diff --git a/Makefile b/Makefile index 29046fa..8c8a637 100644 --- a/Makefile +++ b/Makefile @@ -3103,7 +3103,7 @@ testfuncsort: ./calc${EXT} @${RM} -f func.show func.sort @${CALC_ENV} ./calc${EXT} -d -u show builtin | grep '^[A-Za-z0-9]' > func.show @${CALC_ENV} ./calc${EXT} -d -u show builtin | grep '^[A-Za-z0-9]' | LANG=C LC_ALL=C ${SORT} -d -u > func.sort - @-if ! cmp -s func.show func.sort; then \ + @if ! cmp -s func.show func.sort; then \ echo 1>&2; \ echo "ERROR: builtins[] arrray in func.c is not in dictionary sorted order" 1>&2; \ echo 1>&2; \ diff --git a/cal/Makefile b/cal/Makefile index a9ae95c..50b86af 100644 --- a/cal/Makefile +++ b/cal/Makefile @@ -141,8 +141,9 @@ CALC_FILES= README alg_config.cal beer.cal bernoulli.cal \ test4100.redc.cal test4600.fileop.cal test5100.newdecl.cal \ test5200.globstat.cal test8000.read.cal test8400.quit.cal \ test8500.divmod.cal test8600.maxargs.cal test8700.dotest.cal \ - test8900.special.cal test9300.frem.cal toomcook.cal unitfrac.cal \ - varargs.cal write2file.cal xx_print.cal zeta2.cal + test8900.special.cal test9300.frem.cal test9500.trigeq.cal \ + toomcook.cal unitfrac.cal varargs.cal write2file.cal xx_print.cal \ + zeta2.cal # These calc files are now obsolete and are removed by the install rule. # diff --git a/cal/regress.cal b/cal/regress.cal index 93a124f..8c4523c 100644 --- a/cal/regress.cal +++ b/cal/regress.cal @@ -89,7 +89,7 @@ define vrfy(test, str) ++prob; } if (test != 1) { - print '**** Non-true result (' : test : '): ' : str; + print '**** Non-true result:', str; ++prob; } else { print str; @@ -3877,6 +3877,114 @@ define test_trig() strcat(str(tnum++), ': round(acovercos(2 + 3i, 1e-10), 10) == 0.3076036495+1.8641615442i')); + /* test half versed trigonometric sine */ + vrfy(haversin(0, 1e-10) == 0, + strcat(str(tnum++), ': haversin(0, 1e-10) == 0')); + vrfy(round(haversin(0.2, 1e-10), 10) == 0.0099667111, + strcat(str(tnum++), + ': round(haversin(0.2, 1e-10), 10) == 0.0099667111')); + vrfy(round(haversin(3/7, 1e-10), 10) == 0.0452198242, + strcat(str(tnum++), + ': round(haversin(3/7, 1e-10), 10) == 0.0452198242')); + vrfy(round(haversin(-31, 1e-10), 10) == 0.0426288211, + strcat(str(tnum++), + ': round(haversin(-31, 1e-10), 10) == 0.0426288211')); + vrfy(haversin(pi/3, 1e-10) == 0.25, + strcat(str(tnum++), ': haversin(pi/3, 1e-10) == 0.25')); + vrfy(haversin(pi/2, 1e-10) == 0.5, + strcat(str(tnum++), ': haversin(pi/2, 1e-10) == 0.5')); + vrfy(haversin(pi, 1e-10) == 1, + strcat(str(tnum++), ': haversin(pi, 1e-10) == 1')); + vrfy(haversin(3*pi/2, 1e-10) == 0.5, + strcat(str(tnum++), ': haversin(3*pi/2, 1e-10) == 0.5')); + vrfy(round(haversin(1, 1e-10), 10) == 0.229848847, + strcat(str(tnum++), + ': round(haversin(1, 1e-10), 10) == 0.229848847')); + vrfy(round(haversin(2 + 3i, 1e-10), 10) == 2.5948128455+4.5546139469i, + strcat(str(tnum++), + ': round(haversin(2 + 3i, 1e-10), 10) == 2.5948128455+4.5546139469i')); + + /* test half coversed trigonometric sine */ + vrfy(hacoversin(0, 1e-10) == 0.5, + strcat(str(tnum++), ': hacoversin(0, 1e-10) == 0.5')); + vrfy(round(hacoversin(0.2, 1e-10), 10) == 0.4006653346, + strcat(str(tnum++), + ': round(hacoversin(0.2, 1e-10), 10) == 0.4006653346')); + vrfy(round(hacoversin(3/7, 1e-10), 10) == 0.2922140725, + strcat(str(tnum++), + ': round(hacoversin(3/7, 1e-10), 10) == 0.2922140725')); + vrfy(round(hacoversin(-31, 1e-10), 10) == 0.2979811774, + strcat(str(tnum++), + ': round(hacoversin(-31, 1e-10), 10) == 0.2979811774')); + vrfy(hacoversin(pi/6, 1e-10) == 0.25, + strcat(str(tnum++), ': hacoversin(pi/6, 1e-10) == 0.25')); + vrfy(hacoversin(pi/2, 1e-10) == 0, + strcat(str(tnum++), ': hacoversin(pi/2, 1e-10) == 0')); + vrfy(hacoversin(pi, 1e-10) == 0.5, + strcat(str(tnum++), ': hacoversin(pi, 1e-10) == 0.5')); + vrfy(hacoversin(3*pi/2, 1e-10) == 1, + strcat(str(tnum++), ': hacoversin(3*pi/2, 1e-10) == 1')); + vrfy(round(hacoversin(1, 1e-10), 10) == 0.0792645076, + strcat(str(tnum++), + ': round(hacoversin(1, 1e-10), 10) == 0.0792645076')); + vrfy(round(hacoversin(2 + 3i, 1e-10), 10) == -4.0772495734+2.08445348i, + strcat(str(tnum++), + ': round(hacoversin(2 + 3i, 1e-10), 10) == -4.0772495734+2.08445348i')); + + /* test half versed trigonometric cosine */ + vrfy(havercos(0, 1e-10) == 1, + strcat(str(tnum++), ': havercos(0, 1e-10) == 1')); + vrfy(round(havercos(0.2, 1e-10), 10) == 0.9900332889, + strcat(str(tnum++), + ': round(havercos(0.2, 1e-10), 10) == 0.9900332889')); + vrfy(round(havercos(3/7, 1e-10), 10) == 0.9547801758, + strcat(str(tnum++), + ': round(havercos(3/7, 1e-10), 10) == 0.9547801758')); + vrfy(round(havercos(-31, 1e-10), 10) == 0.9573711789, + strcat(str(tnum++), + ': round(havercos(-31, 1e-10), 10) == 0.9573711789')); + vrfy(havercos(pi/3, 1e-10) == 0.75, + strcat(str(tnum++), ': havercos(pi/3, 1e-10) == 0.75')); + vrfy(havercos(pi/2, 1e-10) == 0.5, + strcat(str(tnum++), ': havercos(pi/2, 1e-10) == 0.5')); + vrfy(havercos(pi, 1e-10) == 0, + strcat(str(tnum++), ': havercos(pi, 1e-10) == 0')); + vrfy(havercos(3*pi/2, 1e-10) == 0.5, + strcat(str(tnum++), ': havercos(3*pi/2, 1e-10) == 0.5')); + vrfy(round(havercos(1, 1e-10), 10) == 0.770151153, + strcat(str(tnum++), + ': round(havercos(1, 1e-10), 10) == 0.770151153')); + vrfy(round(havercos(2 + 3i, 1e-10), 10) == -1.5948128455-4.5546139469i, + strcat(str(tnum++), + ': round(havercos(2 + 3i, 1e-10), 10) == -1.5948128455-4.5546139469i')); + + /* test half coversed trigonometric cosine */ + vrfy(hacovercos(0, 1e-10) == 0.5, + strcat(str(tnum++), ': hacovercos(0, 1e-10) == 0.5')); + vrfy(round(hacovercos(0.2, 1e-10), 10) == 0.5993346654, + strcat(str(tnum++), + ': round(hacovercos(0.2, 1e-10), 10) == 0.5993346654')); + vrfy(round(hacovercos(3/7, 1e-10), 10) == 0.7077859275, + strcat(str(tnum++), + ': round(hacovercos(3/7, 1e-10), 10) == 0.7077859275')); + vrfy(round(hacovercos(-31, 1e-10), 10) == 0.7020188226, + strcat(str(tnum++), + ': round(hacovercos(-31, 1e-10), 10) == 0.7020188226')); + vrfy(hacovercos(pi/6, 1e-10) == 0.75, + strcat(str(tnum++), ': hacovercos(pi/6, 1e-10) == 0.75')); + vrfy(hacovercos(pi/2, 1e-10) == 1, + strcat(str(tnum++), ': hacovercos(pi/2, 1e-10) == 1')); + vrfy(hacovercos(pi, 1e-10) == 0.5, + strcat(str(tnum++), ': hacovercos(pi, 1e-10) == 0.5')); + vrfy(hacovercos(3*pi/2, 1e-10) == 0, + strcat(str(tnum++), ': hacovercos(3*pi/2, 1e-10) == 0')); + vrfy(round(hacovercos(1, 1e-10), 10) == 0.9207354924, + strcat(str(tnum++), + ': round(hacovercos(1, 1e-10), 10) == 0.9207354924')); + vrfy(round(hacovercos(2 + 3i, 1e-10), 10) == 5.0772495734-2.08445348i, + strcat(str(tnum++), + ': round(hacovercos(2 + 3i, 1e-10), 10) == 5.0772495734-2.08445348i')); + print strcat(str(tnum++), ': Ending test_trig'); } print '051: parsed test_trig()'; @@ -4413,8 +4521,10 @@ define test_fileops() print '4280: x = rm("junk4200")'; x = rm("tmp4200"); print '4281: x = rm("tmp4200")'; + x = config("display", 20); + print '4282: x = config("display", 20)'; - print '4282: Ending test_fileops'; + print '4283: Ending test_fileops'; } print '071: parsed test_fileops()'; @@ -9767,9 +9877,29 @@ print '9484: skipping read -once zeta2 - will be read via test8900.special'; print '9485: Ending read of selected calc resource files'; -/* ********************************************************* */ -/* NOTE: ==> Reserved for more read tests 9500-9599 here <== */ -/* ********************************************************* */ +/* + * test 95dd: read test9500.trigeq, test of trigonometric identities + */ +print; +print '9500: Starting trigonometric identities test set'; +epsilon(1e-20),; +print '9501: epsilon(1e-20)'; +read -once test9500.trigeq; +print '9502: read -once test9500.trigeq'; +vrfy(verify_tan(9503) == 0, '9503: verify_tan(9503) == 0'); +vrfy(verify_cot(9504) == 0, '9504: verify_cot(9504) == 0'); +vrfy(verify_sec(9505) == 0, '9505: verify_sec(9505) == 0'); +vrfy(verify_csc(9506) == 0, '9506: verify_csc(9506) == 0'); +vrfy(verify_versin(9507) == 0, '9507: verify_versin(9507) == 0'); +vrfy(verify_coversin(9508) == 0, '9508: verify_coversin(9508) == 0'); +vrfy(verify_vercos(9509) == 0, '9509: verify_vercos(9509) == 0'); +vrfy(verify_covercos(9510) == 0, '9510: verify_covercos(9510) == 0'); +vrfy(verify_haversin(9511) == 0, '9511: verify_haversin(9511) == 0'); +vrfy(verify_hacoversin(9512) == 0, '9512: verify_hacoversin(9512) == 0'); +vrfy(verify_havercos(9513) == 0, '9513: verify_havercos(9513) == 0'); +vrfy(verify_hacovercos(9514) == 0, '9514: verify_hacovercos(9514) == 0'); + +print '9515: Ending trigonometric identities test set'; /* @@ -9790,6 +9920,8 @@ print '9607: u_glob as both static and parameter'; vrfy(config("redecl_warn",1)==0, '9608: config("redecl_warn",1)==0'); vrfy(config("dupvar_warn",1)==0, '9609: config("dupvar_warn",1)==0'); +print '9610: Ending test of dupvar_warn and redecl_warn config parameters'; + /* *********************************************** */ /* NOTE: ==> Room for new tests 9700-9899 here <== */ @@ -10511,6 +10643,12 @@ vrfy_errsym(10553, 10553, "E_ERRSYM_2"); vrfy_errsym(10554, 10554, "E_ERRSYM_3"); vrfy_errsym(10555, 10555, "E_ERRSYM_4"); vrfy_errsym(10556, 10556, "E_ERRSYM_5"); +vrfy_errsym(10557, 10557, "E_HAVERSIN_1"); +vrfy_errsym(10558, 10558, "E_HAVERSIN_2"); +vrfy_errsym(10559, 10559, "E_HAVERSIN_3"); +vrfy_errsym(10560, 10560, "E_AHAVERSIN_1"); +vrfy_errsym(10561, 10561, "E_AHAVERSIN_2"); +vrfy_errsym(10562, 10562, "E_AHAVERSIN_3"); /* ************************************************************** */ /* NOTE: Reserve thru test 10998 for calc computation error codes */ diff --git a/cal/test9500.trigeq.cal b/cal/test9500.trigeq.cal new file mode 100644 index 0000000..fef00d6 --- /dev/null +++ b/cal/test9500.trigeq.cal @@ -0,0 +1,1126 @@ +/* + * test9500.trigeq - test of trigonometric identities for test 95dd + * + * We test over a wide variety of real (NUMBER) and complex (COMPLEX) + * values, that various trigonometric identities hold for all calc + * builtin trigonometric functions. + * + * We assume that the value produced by sin() (trigonometric sine) + * and cos() (trigonometric cosine) are correct. We can reasonably + * assume this because of the calc regression test suite: + * + * (such as test 34dd and in particular cal/test3400.trig.cal and + * such as test 89dd and in particular cal/test8900.special.cal) + * + * will test such those function beforehand. + * + * We use various trigonometric identity to verify the non-inverse + * trigonometric functions other than sin() and cos(). + * + * we verify that it is true within the error tolerance of the eps + * (epsilon) argument. + * + * Copyright (C) 2023 Landon Curt Noll + * + * Primary author: Landon Curt Noll + * + * Calc is open software; you can redistribute it and/or modify it under + * the terms of the version 2.1 of the GNU Lesser General Public License + * as published by the Free Software Foundation. + * + * Calc is distributed in the hope that it will be useful, but WITHOUT + * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY + * or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General + * Public License for more details. + * + * A copy of version 2.1 of the GNU Lesser General Public License is + * distributed with calc under the filename COPYING-LGPL. You should have + * received a copy with calc; if not, write to Free Software Foundation, Inc. + * 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA. + * + * Under source code control: 2023/09/10 14:33:06 + * File existed as early as: 2023 + * + * Share and enjoy! :-) http://www.isthe.com/chongo/tech/comp/calc/ + */ + + +/* + * base_tval - base set of test values + * + * Our 16 base test values must be > 0 and <= 1. + * We use these fractions because, when multiplied + * by pi() will give values that are commonly used. + * + * NOTE: There is no reason to sort these values other than ascetics. :-) + */ +static base_tval = list( + 1/12, /* ~ 0.08333333333333333333 */ + 1/8, /* 0.12500000000000000000 */ + 1/6, /* ~ 0.16666666666666666667 */ + 1/5, /* 0.20000000000000000000 */ + 1/4, /* 0.25000000000000000000 */ + 1/3, /* ~ 0.33333333333333333333 */ + 3/8, /* 0.37500000000000000000 */ + 5/12, /* ~ 0.41666666666666666667 */ + 1/2, /* 0.50000000000000000000 */ + 2/3, /* ~ 0.66666666666666666667 */ + 3/4, /* 0.75000000000000000000 */ + 4/5, /* 0.80000000000000000000 */ + 5/6, /* ~ 0.83333333333333333333 */ + 7/8, /* 0.87500000000000000000 */ + 11/12, /* ~ 0.91666666666666666667 */ + 1 /* 1.00000000000000000000 */ +); + + +/* + * tval - full test list + * + * This long list is formed via the form_tval() function from the base_tval list. + * + * When the form_tval() function forms the full test test values + * we will have 1025 values of various real and complex + * values that also include random values added or subtracted. + * + * See the form_tval() function for details. + */ +tval = list(); + + +/* + * Lists that hold precomputed trigonometric functions on the tval full test list + * + * The precompute_trig() function will fill out these lists. + * They will be the same size as size as the full test list. + */ +sin_tval = list(); /* trigonometric sine of each test test value */ +cos_tval = list(); /* trigonometric cosine of each test test value */ + + +/* + * seed_rand - seed the subtractive 100 shuffle pseudo-random number generator + * + * We want to seed the subtractive 100 shuffle pseudo-random number generator + * with a constant so that the result of this test set will be deterministic. + * + * We seed the subtractive 100 shuffle pseudo-random number generator with + * a 128-bit seed value. The 128-bit seed value was produced when a + * Blum-Blum-Shub pseudo-random number generator was seeded with the + * following 160-bit SHA1 hash of 3 calls to seed(): + * + * sha1(sha1(seed(), seed()<<64, seed()<<128)) == 0x732fda50b1bf6e34b604b7c75b993e5c37a4ad97 + * + * This produced the following 160-bit value that was used to seed + * the Blum-Blum-Shub pseudo-random number: + * + * srandom(0x732fda50b1bf6e34b604b7c75b993e5c37a4ad97) + * + * The seeded Blum-Blum-Shub pseudo-random number generator was then + * used to generate a 128-bit value: + * + * randombit(128) == 0x37301cf47204f7fababc2f32e39ab338 + * + * We seed the subtractive 100 shuffle pseudo-random number generator + * with the above mentioned 128-bit value: + * + * srand(0x37301cf47204f7fababc2f32e39ab338) + */ +define seed_rand() +{ + return srand(0x37301cf47204f7fababc2f32e39ab338); +} + + +/* + * epsilon_bits - number of bits in epsilon() + */ +static epsilon_bits = highbit(1/epsilon()) + 1; + + +/* + * epsilon_norm - norm of a epsilon() + */ +static epsilon_norm = norm(epsilon()); + + +/* + * random_rational_value - compute a random rational value > 0 and < 1 + * + * We compute a pseudo-random value in the range (0,1) with a prevision of epsilon(). + * + * We compute two distinct non-zero pseudo-random integers (i.e., > 0), a and b, + * that are epsilon_bits (see above) long and produced by the seeded + * subtractive 100 shuffle pseudo-random number generator (see seed_rand()). + * + * Return the smaller of a,b divided by the larger of a,b. + */ +define random_rational_value() +{ + local a; /* epsilon_bits pseudo-random integer > 0 */ + local b; /* epsilon_bits pseudo-random integer > 0 */ + + /* + * compute 1st non-zero pseudo-random epsilon_bits integer + */ + do { + a = randbit(epsilon_bits); + } while (a == 0); + + /* + * compute 2nd non-zero pseudo-random epsilon_bits integer also != a + */ + do { + b = randbit(epsilon_bits); + } while (b == 0 || a == b); + + /* + * return the min(a,b) / max(a,b): a value > 0 and < 1 + */ + return min(a,b) / max(a,b); +} + + +/* + * form_tval - form a full list of test values + * + * Given an initial list of 16 values, we form 1025 test values. + * + * From the test set starts with 16 values from the base_tval[] list: + * + * tval = base_tval; + * + * We form a list of 32 values by appending all of the above with: + * + * tval[i] + 1 + * + * We form a list of 64 values by appending all of the above with: + * + * -tval[i] + * + * We form a list of 128 values by appending all of the above with: + * + * For values >= 0: tval[i] + random_rational_value + * For values < 0: tval[i] - random_rational_value + * + * We form a list of 256 values by appending all of the above with: + * + * tval[i] * pi() + * + * We form a list of 512 values by appending all of the above with: + * + * tval[i] * 1i + * + * We for a list of 1024 values by appending new complex values where the + * real part is a randomly selected real value (from the first 256 list values), and the + * imaginary part is a randomly selected complex value (from the last 256 list values). + * + * Finally we append 0 to the list to form the full 1025 value list. + */ +define form_tval() +{ + local tval_len; /* current length of the tval[] array */ + local real_cnt; /* number of real values in tval[] array */ + local imag_cnt; /* number of imaginary values in tval[] array */ + local i; + + /* + * seed the subtractive 100 shuffle pseudo-random number generator + */ + seed_rand(); + + /* + * load base_tval[] array into tval[] array + */ + tval = base_tval; + + /* + * expand all tval with tval+1 + */ + tval_len = size(tval); + for (i=0; i < tval_len; ++i) { + append(tval, tval[i] + 1); + } + + /* + * expand all tval with -tval + */ + tval_len = size(tval); + for (i=0; i < tval_len; ++i) { + append(tval, -tval[i]); + } + + /* + * expand all tval with tval +/- random_rational_value + * + * For values we add random_rational_value if > 0 + * or we subtract crandom_rational_value if < 0. + */ + tval_len = size(tval); + for (i=0; i < tval_len; ++i) { + + /* use the sign to determine if we add or subtract */ + if (tval[i] >= 0) { + append(tval, tval[i] + random_rational_value()); + } else { + append(tval, tval[i] - random_rational_value()); + } + } + + /* + * expand all tval with tval * pi() + */ + tval_len = size(tval); + for (i=0; i < tval_len; ++i) { + append(tval, tval[i] * pi()); + } + + /* + * expand all tval with tval * 1i + */ + tval_len = size(tval); + /* tval[0] thru tval[real_cnt-1] are real values */ + real_cnt = tval_len; + for (i=0; i < tval_len; ++i) { + append(tval, tval[i] * 1i); + } + + /* + * expand tval with randomly selected real value + * plus randomly selected complex value. + */ + tval_len = size(tval); + /* tval[real_cnt] thru tval[real_cnt+imag_cnt-1] are complex values */ + imag_cnt = tval_len - real_cnt; + for (i=0; i < tval_len; ++i) { + append(tval, tval[rand(0,real_cnt)] + tval[real_cnt + rand(0,imag_cnt)]); + } + + /* + * last, append our final 0 value + */ + append(tval, 0); +} + + +/* + * compute - form precomputed trigonometric value lists + */ +define precompute_trig() +{ + local tval_len; /* current length of the tval[] array */ + local i; + + /* + * firewall - be sure tval list is setup + */ + if (size(tval) <= 0) { + form_tval(); + } + tval_len = size(tval); + if (tval_len <= 0) { + quit "FATAL: tval_len is empty"; + } + + /* + * compute trigonometric sine of the tval full test list + */ + for (i=0; i < tval_len; ++i) { + append(sin_tval, sin(tval[i])); + } + if (size(tval) != size(sin_tval)) { + print "****: size(tval):", size(tval), "!= size(sin_tval)", size(sin_tval); + quit "FATAL: sin_tval not same size as tval_len"; + } + + /* + * compute trigonometric cosine of the tval full test list + */ + for (i=0; i < tval_len; ++i) { + append(cos_tval, cos(tval[i])); + } + if (size(tval) != size(cos_tval)) { + print "****: size(tval):", size(tval), "!= size(cos_tval)", size(cos_tval); + quit "FATAL: cos_tval not same size as tval_len"; + } +} + + +/* + * compare - compare the trigonometric identity and the trigonometric function value + * + * given: + * ident_val computed trig value trigonometric identity + * trig_val computed value from the trigonometric function + * name trig function name + * index test index value + * testnum regression test number being performed + * + * returns: + * 0 ==> trigonometric identity and the trigonometric function value match + * 1 ==> trigonometric identity and the trigonometric function value do NOT match + */ +define compare(ident_val, trig_val, name, index, testnum) +{ + local abs_diff; /* absolute difference between ident_val and trig_val */ + + /* + * compute absolute difference + */ + abs_diff = abs(ident_val - trig_val); + + /* + * determine if ident_val within epsilon of trig_val + */ + if (near(norm(abs_diff), epsilon_norm) > 0) { + printf("**** trig test %d-%d failed: %s(tval[%d]): ", + testnum, index, name, index); + printf("%d is not near trig identity value: %d\n", + trig_val, ident_val); + return 1; + } + return 0; +} + + +/* + * not_near_zero - determine if a value is within epsilon() of 0 + * + * given: + * value real or complex value + * + * returns: + * 1 ==> value is farther from 0 than epsilon (not near zero) + * 0 ==> value is within epsilon of 0 (near zero) + */ +define not_near_zero(value) +{ + if (near(norm(value), epsilon_norm) > 0) { + /* value is farther from 0 than epsilon */ + return 1; + } + + /* value is within epsilon of 0 */ + return 0; +} + + +/* + * verify_tan - verify trigonometric tangent + * + * We use the following trigonometric identity: + * + * tan(x) = sin(x) / cos(x) + * + * given: + * testnum regression test number being performed + * + * returns: + * number of tests that failed + */ +define verify_tan(testnum) +{ + local tval_len; /* current length of the tval[] array */ + local ident_val; /* computed trig value trigonometric identity */ + local trig_val; /* computed value from the trigonometric function */ + local error_count; /* number of compare errors detected */ + local i; + + /* + * firewall + */ + if (size(sin_tval) <= 0 || size(cos_tval) <= 0) { + precompute_trig(); + } + + /* + * for each test value, verify the trigonometric identity within epsilon + */ + error_count = 0; + tval_len = size(tval); + for (i=0; i < tval_len; ++i) { + + /* skip test when cos(x) within epsilon of 0 */ + if (not_near_zero(cos_tval[i])) { + + /* compute trigonometric identity */ + ident_val = sin_tval[i] / cos_tval[i]; + + /* compute trigonometric function */ + trig_val = tan(tval[i]); + + /* compare trigonometric identity with trigonometric function value */ + if (compare(ident_val, trig_val, "tan", i, testnum)) { + ++error_count; + } + } + } + + /* + * report test results + */ + if (error_count != 0) { + print '**** test', testnum : ': tan test failure count:', error_count; + } + return error_count; +} + + +/* + * verify_cot - verify trigonometric cotangent + * + * We use the following trigonometric identity: + * + * cot(x) = cos(x) / sin(x) + * + * given: + * testnum regression test number being performed + * + * returns: + * number of tests that failed + */ +define verify_cot(testnum) +{ + local tval_len; /* current length of the tval[] array */ + local ident_val; /* computed trig value trigonometric identity */ + local trig_val; /* computed value from the trigonometric function */ + local error_count; /* number of compare errors detected */ + local i; + + /* + * firewall + */ + if (size(sin_tval) <= 0 || size(cos_tval) <= 0) { + precompute_trig(); + } + + /* + * for each test value, verify the trigonometric identity within epsilon + */ + tval_len = size(tval); + for (i=0; i < tval_len; ++i) { + + /* skip test when sin(x) within epsilon of 0 */ + if (not_near_zero(sin_tval[i])) { + + /* compute trigonometric identity */ + ident_val = cos_tval[i] / sin_tval[i]; + + /* compute trigonometric function */ + trig_val = cot(tval[i]); + + /* compare trigonometric identity with trigonometric function value */ + if (compare(ident_val, trig_val, "cot", i, testnum)) { + ++error_count; + } + } + } + + /* + * report test results + */ + if (error_count != 0) { + print '**** test', testnum : ': cot test failure count:', error_count; + } + return error_count; +} + + +/* + * verify_sec - verify trigonometric secant + * + * We use the following trigonometric identity: + * + * sec(x) = 1 / cos(x) + * + * given: + * testnum regression test number being performed + * + * returns: + * number of tests that failed + */ +define verify_sec(testnum) +{ + local tval_len; /* current length of the tval[] array */ + local ident_val; /* computed trig value trigonometric identity */ + local trig_val; /* computed value from the trigonometric function */ + local error_count; /* number of compare errors detected */ + local i; + + /* + * firewall + */ + if (size(cos_tval) <= 0) { + precompute_trig(); + } + + /* + * for each test value, verify the trigonometric identity within epsilon + */ + tval_len = size(tval); + for (i=0; i < tval_len; ++i) { + + /* skip test when cos(x) within epsilon of 0 */ + if (not_near_zero(cos_tval[i])) { + + /* compute trigonometric identity */ + ident_val = 1 / cos_tval[i]; + + /* compute trigonometric function */ + trig_val = sec(tval[i]); + + /* compare trigonometric identity with trigonometric function value */ + if (compare(ident_val, trig_val, "sec", i, testnum)) { + ++error_count; + } + } + } + + /* + * report test results + */ + if (error_count != 0) { + print '**** test', testnum : ': sec test failure count:', error_count; + } + return error_count; +} + + +/* + * verify_csc - verify trigonometric cosecant + * + * We use the following trigonometric identity: + * + * csc(x) = 1 / sin(x) + * + * given: + * testnum regression test number being performed + * + * returns: + * number of tests that failed + */ +define verify_csc(testnum) +{ + local tval_len; /* current length of the tval[] array */ + local ident_val; /* computed trig value trigonometric identity */ + local trig_val; /* computed value from the trigonometric function */ + local error_count; /* number of compare errors detected */ + local i; + + /* + * firewall + */ + if (size(sin_tval) <= 0) { + precompute_trig(); + } + + /* + * for each test value, verify the trigonometric identity within epsilon + */ + tval_len = size(tval); + for (i=0; i < tval_len; ++i) { + + /* skip test when sin(x) within epsilon of 0 */ + if (not_near_zero(sin_tval[i])) { + + /* compute trigonometric identity */ + ident_val = 1 / sin_tval[i]; + + /* compute trigonometric function */ + trig_val = csc(tval[i]); + + /* compare trigonometric identity with trigonometric function value */ + if (compare(ident_val, trig_val, "csc", i, testnum)) { + ++error_count; + } + } + } + + /* + * report test results + */ + if (error_count != 0) { + print '**** test', testnum : ': csc test failure count:', error_count; + } + return error_count; +} + + +/* + * verify_versin - verify versed trigonometric sine + * + * We use the following trigonometric identity: + * + * versin(x) = 1 - cos(x) + * + * given: + * testnum regression test number being performed + * + * returns: + * number of tests that failed + */ +define verify_versin(testnum) +{ + local tval_len; /* current length of the tval[] array */ + local ident_val; /* computed trig value trigonometric identity */ + local trig_val; /* computed value from the trigonometric function */ + local error_count; /* number of compare errors detected */ + local i; + + /* + * firewall + */ + if (size(cos_tval) <= 0) { + precompute_trig(); + } + + /* + * for each test value, verify the trigonometric identity within epsilon + */ + tval_len = size(tval); + for (i=0; i < tval_len; ++i) { + + /* skip test when cos(x) within epsilon of 0 */ + if (not_near_zero(cos_tval[i])) { + + /* compute trigonometric identity */ + ident_val = 1 - cos_tval[i]; + + /* compute trigonometric function */ + trig_val = versin(tval[i]); + + /* compare trigonometric identity with trigonometric function value */ + if (compare(ident_val, trig_val, "versin", i, testnum)) { + ++error_count; + } + } + } + + /* + * report test results + */ + if (error_count != 0) { + print '**** test', testnum : ': versin test failure count:', error_count; + } + return error_count; +} + + +/* + * verify_coversin - verify coversed trigonometric sine + * + * We use the following trigonometric identity: + * + * coversin(x) = 1 - sin(x) + * + * given: + * testnum regression test number being performed + * + * returns: + * number of tests that failed + */ +define verify_coversin(testnum) +{ + local tval_len; /* current length of the tval[] array */ + local ident_val; /* computed trig value trigonometric identity */ + local trig_val; /* computed value from the trigonometric function */ + local error_count; /* number of compare errors detected */ + local i; + + /* + * firewall + */ + if (size(sin_tval) <= 0) { + precompute_trig(); + } + + /* + * for each test value, verify the trigonometric identity within epsilon + */ + tval_len = size(tval); + for (i=0; i < tval_len; ++i) { + + /* skip test when sin(x) within epsilon of 0 */ + if (not_near_zero(sin_tval[i])) { + + /* compute trigonometric identity */ + ident_val = 1 - sin_tval[i]; + + /* compute trigonometric function */ + trig_val = coversin(tval[i]); + + /* compare trigonometric identity with trigonometric function value */ + if (compare(ident_val, trig_val, "coversin", i, testnum)) { + ++error_count; + } + } + } + + /* + * report test results + */ + if (error_count != 0) { + print '**** test', testnum : ': coversin test failure count:', error_count; + } + return error_count; +} + + +/* + * verify_vercos - verify versed trigonometric cosine + * + * We use the following trigonometric identity: + * + * vercos(x) = 1 + cos(x) + * + * given: + * testnum regression test number being performed + * + * returns: + * number of tests that failed + */ +define verify_vercos(testnum) +{ + local tval_len; /* current length of the tval[] array */ + local ident_val; /* computed trig value trigonometric identity */ + local trig_val; /* computed value from the trigonometric function */ + local error_count; /* number of compare errors detected */ + local i; + + /* + * firewall + */ + if (size(cos_tval) <= 0) { + precompute_trig(); + } + + /* + * for each test value, verify the trigonometric identity within epsilon + */ + tval_len = size(tval); + for (i=0; i < tval_len; ++i) { + + /* skip test when cos(x) within epsilon of 0 */ + if (not_near_zero(cos_tval[i])) { + + /* compute trigonometric identity */ + ident_val = 1 + cos_tval[i]; + + /* compute trigonometric function */ + trig_val = vercos(tval[i]); + + /* compare trigonometric identity with trigonometric function value */ + if (compare(ident_val, trig_val, "vercos", i, testnum)) { + ++error_count; + } + } + } + + /* + * report test results + */ + if (error_count != 0) { + print '**** test', testnum : ': vercos test failure count:', error_count; + } + return error_count; +} + + +/* + * verify_covercos - verify coversed trigonometric cosine + * + * We use the following trigonometric identity: + * + * covercos(x) = 1 + sin(x) + * + * given: + * testnum regression test number being performed + * + * returns: + * number of tests that failed + */ +define verify_covercos(testnum) +{ + local tval_len; /* current length of the tval[] array */ + local ident_val; /* computed trig value trigonometric identity */ + local trig_val; /* computed value from the trigonometric function */ + local error_count; /* number of compare errors detected */ + local i; + + /* + * firewall + */ + if (size(sin_tval) <= 0) { + precompute_trig(); + } + + /* + * for each test value, verify the trigonometric identity within epsilon + */ + tval_len = size(tval); + for (i=0; i < tval_len; ++i) { + + /* skip test when sin(x) within epsilon of 0 */ + if (not_near_zero(sin_tval[i])) { + + /* compute trigonometric identity */ + ident_val = 1 + sin_tval[i]; + + /* compute trigonometric function */ + trig_val = covercos(tval[i]); + + /* compare trigonometric identity with trigonometric function value */ + if (compare(ident_val, trig_val, "covercos", i, testnum)) { + ++error_count; + } + } + } + + /* + * report test results + */ + if (error_count != 0) { + print '**** test', testnum : ': covercos test failure count:', error_count; + } + return error_count; +} + + +/* + * verify_haversin - verify half versed trigonometric sine + * + * We use the following trigonometric identity: + * + * haversin(x) = versin(x) / 2 = (1 - cos(x)) / 2 + * + * given: + * testnum regression test number being performed + * + * returns: + * number of tests that failed + */ +define verify_haversin(testnum) +{ + local tval_len; /* current length of the tval[] array */ + local ident_val; /* computed trig value trigonometric identity */ + local trig_val; /* computed value from the trigonometric function */ + local error_count; /* number of compare errors detected */ + local i; + + /* + * firewall + */ + if (size(cos_tval) <= 0) { + precompute_trig(); + } + + /* + * for each test value, verify the trigonometric identity within epsilon + */ + tval_len = size(tval); + for (i=0; i < tval_len; ++i) { + + /* skip test when cos(x) within epsilon of 0 */ + if (not_near_zero(cos_tval[i])) { + + /* compute trigonometric identity */ + ident_val = (1 - cos_tval[i]) / 2; + + /* compute trigonometric function */ + trig_val = haversin(tval[i]); + + /* compare trigonometric identity with trigonometric function value */ + if (compare(ident_val, trig_val, "haversin", i, testnum)) { + ++error_count; + } + } + } + + /* + * report test results + */ + if (error_count != 0) { + print '**** test', testnum : ': haversin test failure count:', error_count; + } + return error_count; +} + + +/* + * verify_hacoversin - verify half versed trigonometric cosine + * + * We use the following trigonometric identity: + * + * hacoversin(x) = coversin(x) / 2 = (1 - sin(x)) / 2 + * + * given: + * testnum regression test number being performed + * + * returns: + * number of tests that failed + */ +define verify_hacoversin(testnum) +{ + local tval_len; /* current length of the tval[] array */ + local ident_val; /* computed trig value trigonometric identity */ + local trig_val; /* computed value from the trigonometric function */ + local error_count; /* number of compare errors detected */ + local i; + + /* + * firewall + */ + if (size(sin_tval) <= 0) { + precompute_trig(); + } + + /* + * for each test value, verify the trigonometric identity within epsilon + */ + tval_len = size(tval); + for (i=0; i < tval_len; ++i) { + + /* skip test when cos(x) within epsilon of 0 */ + if (not_near_zero(sin_tval[i])) { + + /* compute trigonometric identity */ + ident_val = (1 - sin_tval[i]) / 2; + + /* compute trigonometric function */ + trig_val = hacoversin(tval[i]); + + /* compare trigonometric identity with trigonometric function value */ + if (compare(ident_val, trig_val, "hacoversin", i, testnum)) { + ++error_count; + } + } + } + + /* + * report test results + */ + if (error_count != 0) { + print '**** test', testnum : ': hacoversin test failure count:', error_count; + } + return error_count; +} + + +/* + * verify_havercos - verify half versed trigonometric cosine + * + * We use the following trigonometric identity: + * + * havercos(x) = vercos(x) / 2 = (1 + cos(x)) / 2 + * + * given: + * testnum regression test number being performed + * + * returns: + * number of tests that failed + */ +define verify_havercos(testnum) +{ + local tval_len; /* current length of the tval[] array */ + local ident_val; /* computed trig value trigonometric identity */ + local trig_val; /* computed value from the trigonometric function */ + local error_count; /* number of compare errors detected */ + local i; + + /* + * firewall + */ + if (size(cos_tval) <= 0) { + precompute_trig(); + } + + /* + * for each test value, verify the trigonometric identity within epsilon + */ + tval_len = size(tval); + for (i=0; i < tval_len; ++i) { + + /* skip test when cos(x) within epsilon of 0 */ + if (not_near_zero(cos_tval[i])) { + + /* compute trigonometric identity */ + ident_val = (1 + cos_tval[i]) / 2; + + /* compute trigonometric function */ + trig_val = havercos(tval[i]); + + /* compare trigonometric identity with trigonometric function value */ + if (compare(ident_val, trig_val, "havercos", i, testnum)) { + ++error_count; + } + } + } + + /* + * report test results + */ + if (error_count != 0) { + print '**** test', testnum : ': havercos test failure count:', error_count; + } + return error_count; +} + + +/* + * verify_hacovercos - verify half coversed trigonometric cosine + * + * We use the following trigonometric identity: + * + * hacovercos(x) = covercos(x) / 2 = (1 + sin(x)) / 2 + * + * given: + * testnum regression test number being performed + * + * returns: + * number of tests that failed + */ +define verify_hacovercos(testnum) +{ + local tval_len; /* current length of the tval[] array */ + local ident_val; /* computed trig value trigonometric identity */ + local trig_val; /* computed value from the trigonometric function */ + local error_count; /* number of compare errors detected */ + local i; + + /* + * firewall + */ + if (size(sin_tval) <= 0) { + precompute_trig(); + } + + /* + * for each test value, verify the trigonometric identity within epsilon + */ + tval_len = size(tval); + for (i=0; i < tval_len; ++i) { + + /* skip test when sin(x) within epsilon of 0 */ + if (not_near_zero(sin_tval[i])) { + + /* compute trigonometric identity */ + ident_val = (1 + sin_tval[i]) / 2; + + /* compute trigonometric function */ + trig_val = hacovercos(tval[i]); + + /* compare trigonometric identity with trigonometric function value */ + if (compare(ident_val, trig_val, "hacovercos", i, testnum)) { + ++error_count; + } + } + } + + /* + * report test results + */ + if (error_count != 0) { + print '**** test', testnum : ': hacovercos test failure count:', error_count; + } + return error_count; +} diff --git a/cmath.h b/cmath.h index a1829f6..e069435 100644 --- a/cmath.h +++ b/cmath.h @@ -136,6 +136,14 @@ E_FUNC COMPLEX *c_vercos(COMPLEX *c, NUMBER *epsilon); E_FUNC COMPLEX *c_avercos(COMPLEX *c, NUMBER *epsilon); E_FUNC COMPLEX *c_covercos(COMPLEX *c, NUMBER *epsilon); E_FUNC COMPLEX *c_acovercos(COMPLEX *c, NUMBER *epsilon); +E_FUNC COMPLEX *c_haversin(COMPLEX *c, NUMBER *epsilon); +E_FUNC COMPLEX *c_ahaversin(COMPLEX *c, NUMBER *epsilon); +E_FUNC COMPLEX *c_hacoversin(COMPLEX *c, NUMBER *epsilon); +E_FUNC COMPLEX *c_ahacoversin(COMPLEX *c, NUMBER *epsilon); +E_FUNC COMPLEX *c_havercos(COMPLEX *c, NUMBER *epsilon); +E_FUNC COMPLEX *c_ahavercos(COMPLEX *c, NUMBER *epsilon); +E_FUNC COMPLEX *c_hacovercos(COMPLEX *c, NUMBER *epsilon); +E_FUNC COMPLEX *c_ahacovercos(COMPLEX *c, NUMBER *epsilon); diff --git a/comfunc.c b/comfunc.c index d2c335f..a709530 100644 --- a/comfunc.c +++ b/comfunc.c @@ -1648,7 +1648,7 @@ c_ilog(COMPLEX *c, ZVALUE base) /* - * c_versin - versed sine for COMPLEX values + * c_versin - COMPLEX valued versed trigonometric sine * * This uses the formula: * @@ -1684,13 +1684,13 @@ c_versin(COMPLEX *c, NUMBER *epsilon) */ ctmp = c_cos(c, epsilon); if (ctmp == NULL) { - math_error("Failed to compute complex cos for complex versin"); + math_error("Failed to compute complex cosine for complex versin"); not_reached(); } r = c_sub(&_cone_, ctmp); /* - * return complex 1 - cos(x) + * return trigonometric result */ comfree(ctmp); return r; @@ -1698,7 +1698,7 @@ c_versin(COMPLEX *c, NUMBER *epsilon) /* - * c_aversin - inverse versed sine for COMPLEX values + * c_aversin - COMPLEX valued inverse versed trigonometric sine * * This uses the formula: * @@ -1715,7 +1715,7 @@ COMPLEX * c_aversin(COMPLEX *c, NUMBER *epsilon) { COMPLEX *r; /* inverse trig value result */ - COMPLEX *x; /* argument to inverse trig function */ + COMPLEX *ctmp; /* argument to inverse trig function */ /* * firewall @@ -1732,19 +1732,19 @@ c_aversin(COMPLEX *c, NUMBER *epsilon) /* * calculate complex inverse trig function value */ - x = c_sub(&_cone_, c); - r = c_acos(x, epsilon); - comfree(x); + ctmp = c_sub(&_cone_, c); + r = c_acos(ctmp, epsilon); + comfree(ctmp); /* - * return complex acos(1 - x) + * return inverse trigonometric result */ return r; } /* - * c_coversin - coversed sine for COMPLEX values + * c_coversin - COMPLEX valued coversed trigonometric sine * * This uses the formula: * @@ -1780,21 +1780,21 @@ c_coversin(COMPLEX *c, NUMBER *epsilon) */ ctmp = c_sin(c, epsilon); if (ctmp == NULL) { - math_error("Failed to compute complex sin for complex coversin"); + math_error("Failed to compute complex sine for complex coversin"); not_reached(); } r = c_sub(&_cone_, ctmp); + comfree(ctmp); /* - * return complex 1 - sin(x) + * return trigonometric result */ - comfree(ctmp); return r; } /* - * c_acoversin - inverse versed sine for COMPLEX values + * c_acoversin - COMPLEX valued inverse coversed trigonometric sine * * This uses the formula: * @@ -1811,7 +1811,7 @@ COMPLEX * c_acoversin(COMPLEX *c, NUMBER *epsilon) { COMPLEX *r; /* inverse trig value result */ - COMPLEX *x; /* argument to inverse trig function */ + COMPLEX *ctmp; /* argument to inverse trig function */ /* * firewall @@ -1828,19 +1828,19 @@ c_acoversin(COMPLEX *c, NUMBER *epsilon) /* * calculate complex inverse trig function value */ - x = c_sub(&_cone_, c); - r = c_asin(x, epsilon); - comfree(x); + ctmp = c_sub(&_cone_, c); + r = c_asin(ctmp, epsilon); + comfree(ctmp); /* - * return complex asin(1 - x) + * return inverse trigonometric result */ return r; } /* - * c_vercos - versed sine for COMPLEX values + * c_vercos - COMPLEX valued versed trigonometric cosine * * This uses the formula: * @@ -1876,21 +1876,21 @@ c_vercos(COMPLEX *c, NUMBER *epsilon) */ ctmp = c_cos(c, epsilon); if (ctmp == NULL) { - math_error("Failed to compute complex cos for complex vercos"); + math_error("Failed to compute complex cosine for complex vercos"); not_reached(); } r = c_add(&_cone_, ctmp); + comfree(ctmp); /* - * return complex 1 + cos(x) + * return trigonometric result */ - comfree(ctmp); return r; } /* - * c_avercos - inverse versed sine for COMPLEX values + * c_avercos - COMPLEX valued inverse versed trigonometric cosine * * This uses the formula: * @@ -1907,7 +1907,7 @@ COMPLEX * c_avercos(COMPLEX *c, NUMBER *epsilon) { COMPLEX *r; /* inverse trig value result */ - COMPLEX *x; /* argument to inverse trig function */ + COMPLEX *ctmp; /* argument to inverse trig function */ /* * firewall @@ -1924,19 +1924,19 @@ c_avercos(COMPLEX *c, NUMBER *epsilon) /* * calculate complex inverse trig function value */ - x = c_sub(c, &_cone_); - r = c_acos(x, epsilon); - comfree(x); + ctmp = c_sub(c, &_cone_); + r = c_acos(ctmp, epsilon); + comfree(ctmp); /* - * return complex acos(1 - x) + * return inverse trigonometric result */ return r; } /* - * c_covercos - coversed sine for COMPLEX values + * c_covercos - COMPLEX valued coversed trigonometric cosine * * This uses the formula: * @@ -1972,21 +1972,21 @@ c_covercos(COMPLEX *c, NUMBER *epsilon) */ ctmp = c_sin(c, epsilon); if (ctmp == NULL) { - math_error("Failed to compute complex sin for complex covercos"); + math_error("Failed to compute complex sine for complex covercos"); not_reached(); } r = c_add(&_cone_, ctmp); + comfree(ctmp); /* - * return complex 1 + sin(x) + * return trigonometric result */ - comfree(ctmp); return r; } /* - * c_acovercos - inverse versed sine for COMPLEX values + * c_acovercos - COMPLEX valued inverse coversed trigonometric cosine * * This uses the formula: * @@ -2003,7 +2003,7 @@ COMPLEX * c_acovercos(COMPLEX *c, NUMBER *epsilon) { COMPLEX *r; /* inverse trig value result */ - COMPLEX *x; /* argument to inverse trig function */ + COMPLEX *ctmp; /* argument to inverse trig function */ /* * firewall @@ -2020,12 +2020,408 @@ c_acovercos(COMPLEX *c, NUMBER *epsilon) /* * calculate complex inverse trig function value */ - x = c_sub(c, &_cone_); - r = c_asin(x, epsilon); - comfree(x); + ctmp = c_sub(c, &_cone_); + r = c_asin(ctmp, epsilon); + comfree(ctmp); /* - * return complex asin(x - 1) + * return inverse trigonometric result + */ + return r; +} + + +/* + * c_haversin - COMPLEX valued half versed trigonometric sine + * + * This uses the formula: + * + * haversin(x) = versin(x) / 2 + * + * given: + * q complex value to pass to the trig function + * epsilon error tolerance / precision for trig calculation + * + * returns: + * complex value result of trig function on q with error epsilon + */ +COMPLEX * +c_haversin(COMPLEX *c, NUMBER *epsilon) +{ + COMPLEX *r; /* return COMPLEX value */ + COMPLEX *ctmp; /* complex cos(c) */ + + /* + * firewall + */ + if (c == NULL) { + math_error("%s: c is NULL", __func__); + not_reached(); + } + if (check_epsilon(epsilon) == false) { + math_error("Invalid epsilon arg for %s", __func__); + not_reached(); + } + + /* + * calculate complex trig function value + */ + ctmp = c_versin(c, epsilon); + if (ctmp == NULL) { + math_error("Failed to compute complex versed sine for complex haversin"); + not_reached(); + } + r = c_divq(ctmp, &_qtwo_); + comfree(ctmp); + + /* + * return trigonometric result + */ + return r; +} + + +/* + * c_ahaversin - COMPLEX valued inverse half versed trigonometric sine + * + * This uses the formula: + * + * ahaversin(x) = acos(1 - 2*x) + * + * given: + * q complex value to pass to the trig function + * epsilon error tolerance / precision for trig calculation + * + * returns: + * complex value result of trig function on q with error epsilon + */ +COMPLEX * +c_ahaversin(COMPLEX *c, NUMBER *epsilon) +{ + COMPLEX *r; /* inverse trig value result */ + COMPLEX *ctmp; /* argument to inverse trig function */ + COMPLEX *x2; /* twice x */ + + /* + * firewall + */ + if (c == NULL) { + math_error("%s: c is NULL", __func__); + not_reached(); + } + if (check_epsilon(epsilon) == false) { + math_error("Invalid epsilon arg for %s", __func__); + not_reached(); + } + + /* + * calculate complex inverse trig function value + */ + x2 = c_mulq(c, &_qtwo_); + ctmp = c_sub(&_cone_, x2); + comfree(x2); + r = c_acos(ctmp, epsilon); + comfree(ctmp); + + /* + * return inverse trigonometric result + */ + return r; +} + + +/* + * c_hacoversin - COMPLEX valued half coversed trigonometric sine + * + * This uses the formula: + * + * hacoversin(x) = coversin(x) / 2 + * + * given: + * q complex value to pass to the trig function + * epsilon error tolerance / precision for trig calculation + * + * returns: + * complex value result of trig function on q with error epsilon + */ +COMPLEX * +c_hacoversin(COMPLEX *c, NUMBER *epsilon) +{ + COMPLEX *r; /* return COMPLEX value */ + COMPLEX *ctmp; /* complex sin(c) */ + + /* + * firewall + */ + if (c == NULL) { + math_error("%s: c is NULL", __func__); + not_reached(); + } + if (check_epsilon(epsilon) == false) { + math_error("Invalid epsilon arg for %s", __func__); + not_reached(); + } + + /* + * calculate complex trig function value + */ + ctmp = c_coversin(c, epsilon); + if (ctmp == NULL) { + math_error("Failed to compute complex coversed sine for complex hacoversin"); + not_reached(); + } + r = c_divq(ctmp, &_qtwo_); + comfree(ctmp); + + /* + * return trigonometric result + */ + return r; +} + + +/* + * c_ahacoversin - COMPLEX valued inverse half coversed trigonometric sine + * + * This uses the formula: + * + * ahacoversin(x) = asin(1 - 2*x) + * + * given: + * q complex value to pass to the trig function + * epsilon error tolerance / precision for trig calculation + * + * returns: + * complex value result of trig function on q with error epsilon + */ +COMPLEX * +c_ahacoversin(COMPLEX *c, NUMBER *epsilon) +{ + COMPLEX *r; /* inverse trig value result */ + COMPLEX *ctmp; /* argument to inverse trig function */ + COMPLEX *x2; /* twice x */ + + /* + * firewall + */ + if (c == NULL) { + math_error("%s: c is NULL", __func__); + not_reached(); + } + if (check_epsilon(epsilon) == false) { + math_error("Invalid epsilon arg for %s", __func__); + not_reached(); + } + + /* + * calculate complex inverse trig function value + */ + x2 = c_mulq(c, &_qtwo_); + ctmp = c_sub(&_cone_, x2); + comfree(x2); + r = c_asin(ctmp, epsilon); + comfree(ctmp); + + /* + * return inverse trigonometric result + */ + return r; +} + + +/* + * c_havercos - COMPLEX valued half coversed trigonometric sine + * + * This uses the formula: + * + * havercos(x) = vercos(x) / 2 + * + * given: + * q complex value to pass to the trig function + * epsilon error tolerance / precision for trig calculation + * + * returns: + * complex value result of trig function on q with error epsilon + */ +COMPLEX * +c_havercos(COMPLEX *c, NUMBER *epsilon) +{ + COMPLEX *r; /* return COMPLEX value */ + COMPLEX *ctmp; /* complex sin(c) */ + + /* + * firewall + */ + if (c == NULL) { + math_error("%s: c is NULL", __func__); + not_reached(); + } + if (check_epsilon(epsilon) == false) { + math_error("Invalid epsilon arg for %s", __func__); + not_reached(); + } + + /* + * calculate complex trig function value + */ + ctmp = c_vercos(c, epsilon); + if (ctmp == NULL) { + math_error("Failed to compute complex versed cosine for complex havercos"); + not_reached(); + } + r = c_divq(ctmp, &_qtwo_); + comfree(ctmp); + + /* + * return trigonometric result + */ + return r; +} + + +/* + * c_ahavercos - COMPLEX valued inverse half coversed trigonometric sine + * + * This uses the formula: + * + * ahavercos(x) = acos(2*x - 1) + * + * given: + * q complex value to pass to the trig function + * epsilon error tolerance / precision for trig calculation + * + * returns: + * complex value result of trig function on q with error epsilon + */ +COMPLEX * +c_ahavercos(COMPLEX *c, NUMBER *epsilon) +{ + COMPLEX *r; /* inverse trig value result */ + COMPLEX *ctmp; /* argument to inverse trig function */ + COMPLEX *x2; /* twice x */ + + /* + * firewall + */ + if (c == NULL) { + math_error("%s: c is NULL", __func__); + not_reached(); + } + if (check_epsilon(epsilon) == false) { + math_error("Invalid epsilon arg for %s", __func__); + not_reached(); + } + + /* + * calculate complex inverse trig function value + */ + x2 = c_mulq(c, &_qtwo_); + ctmp = c_sub(&_cone_, x2); + comfree(x2); + r = c_acos(ctmp, epsilon); + comfree(ctmp); + + /* + * return inverse trigonometric result + */ + return r; +} + + +/* + * c_hacovercos - COMPLEX valued half coversed trigonometric cosine + * + * This uses the formula: + * + * hacovercos(x) = covercos(x) / 2 + * + * given: + * q complex value to pass to the trig function + * epsilon error tolerance / precision for trig calculation + * + * returns: + * complex value result of trig function on q with error epsilon + */ +COMPLEX * +c_hacovercos(COMPLEX *c, NUMBER *epsilon) +{ + COMPLEX *r; /* return COMPLEX value */ + COMPLEX *ctmp; /* complex sin(c) */ + + /* + * firewall + */ + if (c == NULL) { + math_error("%s: c is NULL", __func__); + not_reached(); + } + if (check_epsilon(epsilon) == false) { + math_error("Invalid epsilon arg for %s", __func__); + not_reached(); + } + + /* + * calculate complex trig function value + */ + ctmp = c_covercos(c, epsilon); + if (ctmp == NULL) { + math_error("Failed to compute complex coversed cosine for complex hacovercos"); + not_reached(); + } + r = c_divq(ctmp, &_qtwo_); + comfree(ctmp); + + /* + * return trigonometric result + */ + return r; +} + + +/* + * c_ahacovercos - COMPLEX valued inverse half coversed trigonometric cosine + * + * This uses the formula: + * + * ahacovercos(x) = asin(2*x - 1) + * + * given: + * q complex value to pass to the trig function + * epsilon error tolerance / precision for trig calculation + * + * returns: + * complex value result of trig function on q with error epsilon + */ +COMPLEX * +c_ahacovercos(COMPLEX *c, NUMBER *epsilon) +{ + COMPLEX *r; /* inverse trig value result */ + COMPLEX *ctmp; /* argument to inverse trig function */ + COMPLEX *x2; /* twice x */ + + /* + * firewall + */ + if (c == NULL) { + math_error("%s: c is NULL", __func__); + not_reached(); + } + if (check_epsilon(epsilon) == false) { + math_error("Invalid epsilon arg for %s", __func__); + not_reached(); + } + + /* + * calculate complex inverse trig function value + */ + x2 = c_mulq(c, &_qtwo_); + ctmp = c_sub(&_cone_, x2); + comfree(x2); + r = c_asin(ctmp, epsilon); + comfree(ctmp); + + /* + * return inverse trigonometric result */ return r; } diff --git a/errtbl.c b/errtbl.c index b3a7696..9ef2628 100644 --- a/errtbl.c +++ b/errtbl.c @@ -705,6 +705,30 @@ CONST struct errtbl error_table[] = { { 10554, "E_ERRSYM_3", "String argument is not a valid E_STRING for errsym" }, { 10555, "E_ERRSYM_4", "Numeric argument is not an integer for errsym" }, { 10556, "E_ERRSYM_5", "Unable to create a valid E_STRING from the errnum for errsym" }, + { 10557, "E_HAVERSIN_1", "Bad epsilon for haversin" }, + { 10558, "E_HAVERSIN_2", "Bad first argument for haversin" }, + { 10559, "E_HAVERSIN_3", "Too-large im(argument) for haversin" }, + { 10560, "E_AHAVERSIN_1", "Bad epsilon for ahaversin" }, + { 10561, "E_AHAVERSIN_2", "Bad first argument for ahaversin" }, + { 10562, "E_AHAVERSIN_3", "Too-large im(argument) for ahaversin" }, + { 10563, "E_HACOVERSIN_1", "Bad epsilon for hacoversin" }, + { 10564, "E_HACOVERSIN_2", "Bad first argument for hacoversin" }, + { 10565, "E_HACOVERSIN_3", "Too-large im(argument) for hacoversin" }, + { 10566, "E_AHACOVERSIN_1", "Bad epsilon for ahacoversin" }, + { 10567, "E_AHACOVERSIN_2", "Bad first argument for achaoversin" }, + { 10568, "E_AHACOVERSIN_3", "Too-large im(argument) for ahacoversin" }, + { 10569, "E_HAVERCOS_1", "Bad epsilon for havercos" }, + { 10570, "E_HAVERCOS_2", "Bad first argument for havercos" }, + { 10571, "E_HAVERCOS_3", "Too-large im(argument) for havercos" }, + { 10572, "E_AHAVERCOS_1", "Bad epsilon for ahavercos" }, + { 10573, "E_AHAVERCOS_2", "Bad first argument for ahavercos" }, + { 10574, "E_AHAVERCOS_3", "Too-large im(argument) for ahavercos" }, + { 10575, "E_HACOVERCOS_1", "Bad epsilon for hacovercos" }, + { 10576, "E_HACOVERCOS_2", "Bad first argument for hacovercos" }, + { 10577, "E_HACOVERCOS_3", "Too-large im(argument) for hacovercos" }, + { 10578, "E_AHACOVERCOS_1", "Bad epsilon for ahacovercos" }, + { 10579, "E_AHACOVERCOS_2", "Bad first argument for achaovercos" }, + { 10580, "E_AHACOVERCOS_3", "Too-large im(argument) for ahacovercos" }, /* IMPORTANT NOTE: add new entries above here and be sure their errnum numeric value is consecutive! */ /* The next NULL entry must be last */ diff --git a/func.c b/func.c index 971c828..4d7726e 100644 --- a/func.c +++ b/func.c @@ -10821,7 +10821,7 @@ f_version(void) /* - * f_versin - versed sine + * f_versin - versed trigonometric sine */ S_FUNC VALUE f_versin(int count, VALUE **vals) @@ -10878,7 +10878,7 @@ f_versin(int count, VALUE **vals) /* - * f_aversin - inverse versed sine + * f_aversin - inverse versed trigonometric sine */ S_FUNC VALUE f_aversin(int count, VALUE **vals) @@ -10910,7 +10910,7 @@ f_aversin(int count, VALUE **vals) arg1 = *vals[0]; if (arg1.v_type == V_NUM) { - /* try to compute result using real triv function */ + /* try to compute result using real trig function */ result.v_num = qaversin_or_NULL(arg1.v_num, eps); /* @@ -10961,7 +10961,7 @@ f_aversin(int count, VALUE **vals) /* - * f_coversin - coversed sine + * f_coversin - coversed trigonometric sine */ S_FUNC VALUE f_coversin(int count, VALUE **vals) @@ -11018,7 +11018,7 @@ f_coversin(int count, VALUE **vals) /* - * f_acoversin - inverse coversed sine + * f_acoversin - inverse coversed trigonometric sine */ S_FUNC VALUE f_acoversin(int count, VALUE **vals) @@ -11050,7 +11050,7 @@ f_acoversin(int count, VALUE **vals) arg1 = *vals[0]; if (arg1.v_type == V_NUM) { - /* try to compute result using real triv function */ + /* try to compute result using real trig function */ result.v_num = qacoversin_or_NULL(arg1.v_num, eps); /* @@ -11101,7 +11101,7 @@ f_acoversin(int count, VALUE **vals) /* - * f_vercos - versed cosine + * f_vercos - versed trigonometric cosine */ S_FUNC VALUE f_vercos(int count, VALUE **vals) @@ -11158,7 +11158,7 @@ f_vercos(int count, VALUE **vals) /* - * f_avercos - inverse versed cosine + * f_avercos - inverse versed trigonometric cosine */ S_FUNC VALUE f_avercos(int count, VALUE **vals) @@ -11190,7 +11190,7 @@ f_avercos(int count, VALUE **vals) arg1 = *vals[0]; if (arg1.v_type == V_NUM) { - /* try to compute result using real triv function */ + /* try to compute result using real trig function */ result.v_num = qavercos_or_NULL(arg1.v_num, eps); /* @@ -11241,7 +11241,7 @@ f_avercos(int count, VALUE **vals) /* - * f_covercos - coversed cosine + * f_covercos - coversed trigonometric cosine */ S_FUNC VALUE f_covercos(int count, VALUE **vals) @@ -11298,7 +11298,7 @@ f_covercos(int count, VALUE **vals) /* - * f_acovercos - inverse coversed cosine + * f_acovercos - inverse coversed trigonometric cosine */ S_FUNC VALUE f_acovercos(int count, VALUE **vals) @@ -11330,7 +11330,7 @@ f_acovercos(int count, VALUE **vals) arg1 = *vals[0]; if (arg1.v_type == V_NUM) { - /* try to compute result using real triv function */ + /* try to compute result using real trig function */ result.v_num = qacovercos_or_NULL(arg1.v_num, eps); /* @@ -11380,6 +11380,566 @@ f_acovercos(int count, VALUE **vals) } +/* + * f_haversin - half versed trigonometric sine + */ +S_FUNC VALUE +f_haversin(int count, VALUE **vals) +{ + VALUE result; + COMPLEX *c; + NUMBER *eps; + + /* initialize VALUE */ + result.v_subtype = V_NOSUBTYPE; + + /* + * set error tolerance for builtin function + * + * Use eps VALUE arg if given and value is in a valid range. + */ + eps = conf->epsilon; + if (count == 2) { + if (verify_eps(vals[1]) == false) { + return error_value(E_HAVERSIN_1); + } + eps = vals[1]->v_num; + } + + /* + * compute trig function to a given error tolerance + */ + switch (vals[0]->v_type) { + case V_NUM: + result.v_num = qhaversin(vals[0]->v_num, eps); + result.v_type = V_NUM; + break; + case V_COM: + c = c_haversin(vals[0]->v_com, eps); + if (c == NULL) { + return error_value(E_HAVERSIN_3); + } + result.v_com = c; + result.v_type = V_COM; + + /* + * case: complex trig function returned real, convert result to NUMBER + */ + if (cisreal(c)) { + result.v_num = c_to_q(c, true); + result.v_type = V_NUM; + } + break; + default: + return error_value(E_HAVERSIN_2); + } + return result; +} + + +/* + * f_ahaversin - inverse half versed trigonometric sine + */ +S_FUNC VALUE +f_ahaversin(int count, VALUE **vals) +{ + VALUE arg1; /* 1st arg if it is a COMPLEX value */ + VALUE result; /* value to return */ + COMPLEX *c; /* COMPLEX trig result */ + NUMBER *eps; /* epsilon error tolerance */ + + /* initialize VALUE */ + result.v_subtype = V_NOSUBTYPE; + + /* + * set error tolerance for builtin function + * + * Use eps VALUE arg if given and value is in a valid range. + */ + eps = conf->epsilon; + if (count == 2) { + if (verify_eps(vals[1]) == false) { + return error_value(E_AHAVERSIN_1); + } + eps = vals[1]->v_num; + } + + /* + * compute inverse trig function to a given error tolerance + */ + arg1 = *vals[0]; + if (arg1.v_type == V_NUM) { + + /* try to compute result using real trig function */ + result.v_num = qahaversin_or_NULL(arg1.v_num, eps); + + /* + * case: trig function returned a NUMBER + */ + if (result.v_num != NULL) { + result.v_type = V_NUM; + + /* + * case: trig function returned NULL - need to try COMPLEX trig function + */ + } else { + /* convert NUMBER argument from NUMBER to COMPLEX */ + arg1.v_com = qqtoc(arg1.v_num, &_qzero_); + arg1.v_type = V_COM; + } + } + if (arg1.v_type == V_COM) { + + /* + * case: argument was COMPLEX or + * trig function returned NULL and argument was converted to COMPLEX + */ + c = c_ahaversin(arg1.v_com, eps); + if (c == NULL) { + return error_value(E_AHAVERSIN_3); + } + result.v_com = c; + result.v_type = V_COM; + + /* + * case: complex trig function returned real, convert result to NUMBER + */ + if (cisreal(c)) { + result.v_num = c_to_q(c, true); + result.v_type = V_NUM; + } + } + if (arg1.v_type != V_NUM && arg1.v_type != V_COM) { + + /* + * case: argument type is not valid for this function + */ + return error_value(E_AHAVERSIN_2); + } + return result; +} + + +/* + * f_hacoversin - half coversed trigonometric sine + */ +S_FUNC VALUE +f_hacoversin(int count, VALUE **vals) +{ + VALUE result; + COMPLEX *c; + NUMBER *eps; + + /* initialize VALUE */ + result.v_subtype = V_NOSUBTYPE; + + /* + * set error tolerance for builtin function + * + * Use eps VALUE arg if given and value is in a valid range. + */ + eps = conf->epsilon; + if (count == 2) { + if (verify_eps(vals[1]) == false) { + return error_value(E_HACOVERSIN_1); + } + eps = vals[1]->v_num; + } + + /* + * compute trig function to a given error tolerance + */ + switch (vals[0]->v_type) { + case V_NUM: + result.v_num = qhacoversin(vals[0]->v_num, eps); + result.v_type = V_NUM; + break; + case V_COM: + c = c_hacoversin(vals[0]->v_com, eps); + if (c == NULL) { + return error_value(E_HACOVERSIN_3); + } + result.v_com = c; + result.v_type = V_COM; + + /* + * case: complex trig function returned real, convert result to NUMBER + */ + if (cisreal(c)) { + result.v_num = c_to_q(c, true); + result.v_type = V_NUM; + } + break; + default: + return error_value(E_HACOVERSIN_2); + } + return result; +} + + +/* + * f_ahacoversin - inverse half coversed trigonometric sine + */ +S_FUNC VALUE +f_ahacoversin(int count, VALUE **vals) +{ + VALUE arg1; /* 1st arg if it is a COMPLEX value */ + VALUE result; /* value to return */ + COMPLEX *c; /* COMPLEX trig result */ + NUMBER *eps; /* epsilon error tolerance */ + + /* initialize VALUE */ + result.v_subtype = V_NOSUBTYPE; + + /* + * set error tolerance for builtin function + * + * Use eps VALUE arg if given and value is in a valid range. + */ + eps = conf->epsilon; + if (count == 2) { + if (verify_eps(vals[1]) == false) { + return error_value(E_AHACOVERSIN_1); + } + eps = vals[1]->v_num; + } + + /* + * compute inverse trig function to a given error tolerance + */ + arg1 = *vals[0]; + if (arg1.v_type == V_NUM) { + + /* try to compute result using real trig function */ + result.v_num = qahacoversin_or_NULL(arg1.v_num, eps); + + /* + * case: trig function returned a NUMBER + */ + if (result.v_num != NULL) { + result.v_type = V_NUM; + + /* + * case: trig function returned NULL - need to try COMPLEX trig function + */ + } else { + /* convert NUMBER argument from NUMBER to COMPLEX */ + arg1.v_com = qqtoc(arg1.v_num, &_qzero_); + arg1.v_type = V_COM; + } + } + if (arg1.v_type == V_COM) { + + /* + * case: argument was COMPLEX or + * trig function returned NULL and argument was converted to COMPLEX + */ + c = c_ahacoversin(arg1.v_com, eps); + if (c == NULL) { + return error_value(E_AHACOVERSIN_3); + } + result.v_com = c; + result.v_type = V_COM; + + /* + * case: complex trig function returned real, convert result to NUMBER + */ + if (cisreal(c)) { + result.v_num = c_to_q(c, true); + result.v_type = V_NUM; + } + } + if (arg1.v_type != V_NUM && arg1.v_type != V_COM) { + + /* + * case: argument type is not valid for this function + */ + return error_value(E_AHACOVERSIN_2); + } + return result; +} + + +/* + * f_havercos - half versed trigonometric cosine + */ +S_FUNC VALUE +f_havercos(int count, VALUE **vals) +{ + VALUE result; + COMPLEX *c; + NUMBER *eps; + + /* initialize VALUE */ + result.v_subtype = V_NOSUBTYPE; + + /* + * set error tolerance for builtin function + * + * Use eps VALUE arg if given and value is in a valid range. + */ + eps = conf->epsilon; + if (count == 2) { + if (verify_eps(vals[1]) == false) { + return error_value(E_HAVERCOS_1); + } + eps = vals[1]->v_num; + } + + /* + * compute trig function to a given error tolerance + */ + switch (vals[0]->v_type) { + case V_NUM: + result.v_num = qhavercos(vals[0]->v_num, eps); + result.v_type = V_NUM; + break; + case V_COM: + c = c_havercos(vals[0]->v_com, eps); + if (c == NULL) { + return error_value(E_HAVERCOS_3); + } + result.v_com = c; + result.v_type = V_COM; + + /* + * case: complex trig function returned real, convert result to NUMBER + */ + if (cisreal(c)) { + result.v_num = c_to_q(c, true); + result.v_type = V_NUM; + } + break; + default: + return error_value(E_HAVERCOS_2); + } + return result; +} + + +/* + * f_ahavercos - inverse half versed trigonometric cosine + */ +S_FUNC VALUE +f_ahavercos(int count, VALUE **vals) +{ + VALUE arg1; /* 1st arg if it is a COMPLEX value */ + VALUE result; /* value to return */ + COMPLEX *c; /* COMPLEX trig result */ + NUMBER *eps; /* epsilon error tolerance */ + + /* initialize VALUE */ + result.v_subtype = V_NOSUBTYPE; + + /* + * set error tolerance for builtin function + * + * Use eps VALUE arg if given and value is in a valid range. + */ + eps = conf->epsilon; + if (count == 2) { + if (verify_eps(vals[1]) == false) { + return error_value(E_AHAVERCOS_1); + } + eps = vals[1]->v_num; + } + + /* + * compute inverse trig function to a given error tolerance + */ + arg1 = *vals[0]; + if (arg1.v_type == V_NUM) { + + /* try to compute result using real trig function */ + result.v_num = qahavercos_or_NULL(arg1.v_num, eps); + + /* + * case: trig function returned a NUMBER + */ + if (result.v_num != NULL) { + result.v_type = V_NUM; + + /* + * case: trig function returned NULL - need to try COMPLEX trig function + */ + } else { + /* convert NUMBER argument from NUMBER to COMPLEX */ + arg1.v_com = qqtoc(arg1.v_num, &_qzero_); + arg1.v_type = V_COM; + } + } + if (arg1.v_type == V_COM) { + + /* + * case: argument was COMPLEX or + * trig function returned NULL and argument was converted to COMPLEX + */ + c = c_ahavercos(arg1.v_com, eps); + if (c == NULL) { + return error_value(E_AHAVERCOS_3); + } + result.v_com = c; + result.v_type = V_COM; + + /* + * case: complex trig function returned real, convert result to NUMBER + */ + if (cisreal(c)) { + result.v_num = c_to_q(c, true); + result.v_type = V_NUM; + } + } + if (arg1.v_type != V_NUM && arg1.v_type != V_COM) { + + /* + * case: argument type is not valid for this function + */ + return error_value(E_AHAVERCOS_2); + } + return result; +} + + +/* + * f_hacovercos - half coversed trigonometric cosine + */ +S_FUNC VALUE +f_hacovercos(int count, VALUE **vals) +{ + VALUE result; + COMPLEX *c; + NUMBER *eps; + + /* initialize VALUE */ + result.v_subtype = V_NOSUBTYPE; + + /* + * set error tolerance for builtin function + * + * Use eps VALUE arg if given and value is in a valid range. + */ + eps = conf->epsilon; + if (count == 2) { + if (verify_eps(vals[1]) == false) { + return error_value(E_HACOVERCOS_1); + } + eps = vals[1]->v_num; + } + + /* + * compute trig function to a given error tolerance + */ + switch (vals[0]->v_type) { + case V_NUM: + result.v_num = qhacovercos(vals[0]->v_num, eps); + result.v_type = V_NUM; + break; + case V_COM: + c = c_hacovercos(vals[0]->v_com, eps); + if (c == NULL) { + return error_value(E_HACOVERCOS_3); + } + result.v_com = c; + result.v_type = V_COM; + + /* + * case: complex trig function returned real, convert result to NUMBER + */ + if (cisreal(c)) { + result.v_num = c_to_q(c, true); + result.v_type = V_NUM; + } + break; + default: + return error_value(E_HACOVERCOS_2); + } + return result; +} + + +/* + * f_ahacovercos - inverse half coversed trigonometric cosine + */ +S_FUNC VALUE +f_ahacovercos(int count, VALUE **vals) +{ + VALUE arg1; /* 1st arg if it is a COMPLEX value */ + VALUE result; /* value to return */ + COMPLEX *c; /* COMPLEX trig result */ + NUMBER *eps; /* epsilon error tolerance */ + + /* initialize VALUE */ + result.v_subtype = V_NOSUBTYPE; + + /* + * set error tolerance for builtin function + * + * Use eps VALUE arg if given and value is in a valid range. + */ + eps = conf->epsilon; + if (count == 2) { + if (verify_eps(vals[1]) == false) { + return error_value(E_AHACOVERCOS_1); + } + eps = vals[1]->v_num; + } + + /* + * compute inverse trig function to a given error tolerance + */ + arg1 = *vals[0]; + if (arg1.v_type == V_NUM) { + + /* try to compute result using real trig function */ + result.v_num = qahacovercos_or_NULL(arg1.v_num, eps); + + /* + * case: trig function returned a NUMBER + */ + if (result.v_num != NULL) { + result.v_type = V_NUM; + + /* + * case: trig function returned NULL - need to try COMPLEX trig function + */ + } else { + /* convert NUMBER argument from NUMBER to COMPLEX */ + arg1.v_com = qqtoc(arg1.v_num, &_qzero_); + arg1.v_type = V_COM; + } + } + if (arg1.v_type == V_COM) { + + /* + * case: argument was COMPLEX or + * trig function returned NULL and argument was converted to COMPLEX + */ + c = c_ahacovercos(arg1.v_com, eps); + if (c == NULL) { + return error_value(E_AHACOVERCOS_3); + } + result.v_com = c; + result.v_type = V_COM; + + /* + * case: complex trig function returned real, convert result to NUMBER + */ + if (cisreal(c)) { + result.v_num = c_to_q(c, true); + result.v_type = V_NUM; + } + } + if (arg1.v_type != V_NUM && arg1.v_type != V_COM) { + + /* + * case: argument type is not valid for this function + */ + return error_value(E_AHACOVERCOS_2); + } + return result; +} + + #endif /* !FUNCLIST */ @@ -11442,6 +12002,14 @@ STATIC CONST struct builtin builtins[] = { "inverse csch of a within accuracy b"}, {"agd", 1, 2, 0, OP_NOP, {.null = NULL}, {.valfunc_cnt = f_agd}, "inverse Gudermannian function"}, + {"ahacovercos", 1, 2, 0, OP_NOP, {.null = NULL}, {.valfunc_cnt = f_ahacovercos}, + "inverse half coversed cosine of a within accuracy b"}, + {"ahacoversin", 1, 2, 0, OP_NOP, {.null = NULL}, {.valfunc_cnt = f_ahacoversin}, + "inverse half coversed sine of a within accuracy b"}, + {"ahavercos", 1, 2, 0, OP_NOP, {.null = NULL}, {.valfunc_cnt = f_ahavercos}, + "inverse half versed cosine of a within accuracy b"}, + {"ahaversin", 1, 2, 0, OP_NOP, {.null = NULL}, {.valfunc_cnt = f_ahaversin}, + "inverse half versed sine of a within accuracy b"}, {"append", 1, IN, FA, OP_NOP, {.null = NULL}, {.valfunc_cnt = f_listappend}, "append values to end of list"}, {"appr", 1, 3, 0, OP_NOP, {.null = NULL}, {.valfunc_cnt = f_appr}, @@ -11693,9 +12261,17 @@ STATIC CONST struct builtin builtins[] = { "convert a to b hours, c min, rounding type d\n"}, {"h2hms", 4, 5, FA, OP_NOP, {.null = NULL}, {.valfunc_cnt = f_h2hms}, "convert a to b hours, c min, d sec, rounding type e\n"}, + {"hacovercos", 1, 2, 0, OP_NOP, {.null = NULL}, {.valfunc_cnt = f_hacovercos}, + "half coversed cosine of value a within accuracy b"}, + {"hacoversin", 1, 2, 0, OP_NOP, {.null = NULL}, {.valfunc_cnt = f_hacoversin}, + "half coversed sine of value a within accuracy b"}, {"hash", 1, IN, 0, OP_NOP, {.null = NULL}, {.valfunc_cnt = f_hash}, "return non-negative hash value for one or\n" "\t\t\tmore values"}, + {"havercos", 1, 2, 0, OP_NOP, {.null = NULL}, {.valfunc_cnt = f_havercos}, + "half versed cosine of value a within accuracy b"}, + {"haversin", 1, 2, 0, OP_NOP, {.null = NULL}, {.valfunc_cnt = f_haversin}, + "half versed sine of value a within accuracy b"}, {"head", 2, 2, 0, OP_NOP, {.null = NULL}, {.valfunc_2 = f_head}, "return list of specified number at head of a list"}, {"highbit", 1, 1, 0, OP_HIGHBIT, {.null = NULL}, {.null = NULL}, diff --git a/help/randperm b/help/randperm index bc50ee6..4d07ad0 100644 --- a/help/randperm +++ b/help/randperm @@ -44,9 +44,9 @@ LINK LIBRARY none SEE ALSO - comp, fact, rand, perm + comb, fact, rand, perm -## Copyright (C) 1999-2006 Landon Curt Noll +## Copyright (C) 1999-2006,2023 Landon Curt Noll ## ## Calc is open software; you can redistribute it and/or modify it under ## the terms of the version 2.1 of the GNU Lesser General Public License diff --git a/qmath.h b/qmath.h index 9daa909..7c8e369 100644 --- a/qmath.h +++ b/qmath.h @@ -238,7 +238,18 @@ E_FUNC NUMBER *qavercos(NUMBER *q, NUMBER *epsilon); E_FUNC NUMBER *qcovercos(NUMBER *q, NUMBER *epsilon); E_FUNC NUMBER *qacovercos_or_NULL(NUMBER *q, NUMBER *epsilon); E_FUNC NUMBER *qacovercos(NUMBER *q, NUMBER *epsilon); - +E_FUNC NUMBER *qhaversin(NUMBER *q, NUMBER *epsilon); +E_FUNC NUMBER *qahaversin_or_NULL(NUMBER *q, NUMBER *epsilon); +E_FUNC NUMBER *qahaversin(NUMBER *q, NUMBER *epsilon); +E_FUNC NUMBER *qhacoversin(NUMBER *q, NUMBER *epsilon); +E_FUNC NUMBER *qahacoversin_or_NULL(NUMBER *q, NUMBER *epsilon); +E_FUNC NUMBER *qahacoversin(NUMBER *q, NUMBER *epsilon); +E_FUNC NUMBER *qhavercos(NUMBER *q, NUMBER *epsilon); +E_FUNC NUMBER *qahavercos_or_NULL(NUMBER *q, NUMBER *epsilon); +E_FUNC NUMBER *qahavercos(NUMBER *q, NUMBER *epsilon); +E_FUNC NUMBER *qhacovercos(NUMBER *q, NUMBER *epsilon); +E_FUNC NUMBER *qahacovercos_or_NULL(NUMBER *q, NUMBER *epsilon); +E_FUNC NUMBER *qahacovercos(NUMBER *q, NUMBER *epsilon); /* * pseudo-seed generator diff --git a/qtrans.c b/qtrans.c index 46b2caa..fb73bad 100644 --- a/qtrans.c +++ b/qtrans.c @@ -1979,7 +1979,7 @@ qacoth(NUMBER *q, NUMBER *epsilon) /* - * qversin - versed sine for NUMBER values + * qversin - versed trigonometric sine * * This uses the formula: * @@ -1995,8 +1995,8 @@ qacoth(NUMBER *q, NUMBER *epsilon) NUMBER * qversin(NUMBER *q, NUMBER *epsilon) { - NUMBER *res; - NUMBER *cos; + NUMBER *res; /* inverse trig value result */ + NUMBER *qtmp; /* argument to inverse trig function */ /* * firewall @@ -2013,19 +2013,19 @@ qversin(NUMBER *q, NUMBER *epsilon) /* * calculate trig function value */ - cos = qcos(q, epsilon); - res = qsub(&_qone_, cos); - qfree(cos); + qtmp = qcos(q, epsilon); + res = qsub(&_qone_, qtmp); + qfree(qtmp); /* - * return 1 - cos(x) + * return trigonometric result */ return res; } /* - * qaversin_or_NULL - inverse versed sine for NUMBER values + * qaversin_or_NULL - return NULL or non-complex inverse versed trigonometric sine * * This uses the formula: * @@ -2047,7 +2047,7 @@ NUMBER * qaversin_or_NULL(NUMBER *q, NUMBER *epsilon) { NUMBER *res; /* inverse trig value result */ - NUMBER *x; /* argument to inverse trig function */ + NUMBER *qtmp; /* argument to inverse trig function */ /* * firewall @@ -2064,22 +2064,22 @@ qaversin_or_NULL(NUMBER *q, NUMBER *epsilon) /* * calculate inverse trig function value */ - x = qsub(&_qone_, q); - res = qacos(x, epsilon); - qfree(x); + qtmp = qsub(&_qone_, q); + res = qacos(qtmp, epsilon); + qfree(qtmp); if (res == NULL) { return NULL; } /* - * return acos(1 - x) + * return inverse trigonometric result */ return res; } /* - * qaversin - inverse versed sine for NUMBER values + * qaversin - non-complex inverse versed trigonometric sine * * This uses the formula: * @@ -2114,19 +2114,19 @@ qaversin(NUMBER *q, NUMBER *epsilon) */ res = qaversin_or_NULL(q, epsilon); if (res == NULL) { - math_error("cannot compute inverse cos for aversin"); + math_error("cannot compute inverse cosine for aversin"); not_reached(); } /* - * return acos(1 - x) + * return inverse trigonometric result */ return res; } /* - * qcoversin - coversed sine for NUMBER values + * qcoversin - coversed trigonometric sine * * This uses the formula: * @@ -2142,7 +2142,8 @@ qaversin(NUMBER *q, NUMBER *epsilon) NUMBER * qcoversin(NUMBER *q, NUMBER *epsilon) { - NUMBER *sin, *res; + NUMBER *res; /* inverse trig value result */ + NUMBER *qtmp; /* argument to inverse trig function */ /* * firewall @@ -2159,19 +2160,19 @@ qcoversin(NUMBER *q, NUMBER *epsilon) /* * calculate trig function value */ - sin = qsin(q, epsilon); - res = qsub(&_qone_, sin); - qfree(sin); + qtmp = qsin(q, epsilon); + res = qsub(&_qone_, qtmp); + qfree(qtmp); /* - * return 1 - sin(x) + * return trigonometric result */ return res; } /* - * qacoversin_or_NULL - inverse coversed sine for NUMBER values + * qacoversin_or_NULL - return NULL or non-complex inverse coversed trigonometric sine * * This uses the formula: * @@ -2196,7 +2197,7 @@ NUMBER * qacoversin_or_NULL(NUMBER *q, NUMBER *epsilon) { NUMBER *res; /* inverse trig value result */ - NUMBER *x; /* argument to inverse trig function */ + NUMBER *qtmp; /* argument to inverse trig function */ /* * firewall @@ -2213,22 +2214,22 @@ qacoversin_or_NULL(NUMBER *q, NUMBER *epsilon) /* * calculate inverse trig function value */ - x = qsub(&_qone_, q); - res = qasin(x, epsilon); - qfree(x); + qtmp = qsub(&_qone_, q); + res = qasin(qtmp, epsilon); + qfree(qtmp); if (res == NULL) { return NULL; } /* - * return asin(1 - x) + * return inverse trigonometric result */ return res; } /* - * qacoversin - inverse coversed sine for NUMBER values + * qacoversin - non-complex inverse coversed trigonometric sine * * This uses the formula: * @@ -2263,19 +2264,19 @@ qacoversin(NUMBER *q, NUMBER *epsilon) */ res = qacoversin_or_NULL(q, epsilon); if (res == NULL) { - math_error("cannot compute inverse sin for acoversin"); + math_error("cannot compute inverse sine for acoversin"); not_reached(); } /* - * return asin(1 - x) + * return inverse trigonometric result */ return res; } /* - * qvercos - versed sine for NUMBER values + * qvercos - versed trigonometric cosine * * This uses the formula: * @@ -2291,7 +2292,8 @@ qacoversin(NUMBER *q, NUMBER *epsilon) NUMBER * qvercos(NUMBER *q, NUMBER *epsilon) { - NUMBER *cos, *res; + NUMBER *res; /* inverse trig value result */ + NUMBER *qtmp; /* argument to inverse trig function */ /* * firewall @@ -2308,19 +2310,19 @@ qvercos(NUMBER *q, NUMBER *epsilon) /* * calculate trig function value */ - cos = qcos(q, epsilon); - res = qqadd(&_qone_, cos); - qfree(cos); + qtmp = qcos(q, epsilon); + res = qqadd(&_qone_, qtmp); + qfree(qtmp); /* - * return 1 + cos(x) + * return trigonometric result */ return res; } /* - * qavercos_or_NULL - inverse versed sine for NUMBER values + * qavercos_or_NULL - return NULL or non-complex inverse versed trigonometric cosine * * This uses the formula: * @@ -2342,7 +2344,7 @@ NUMBER * qavercos_or_NULL(NUMBER *q, NUMBER *epsilon) { NUMBER *res; /* inverse trig value result */ - NUMBER *x; /* argument to inverse trig function */ + NUMBER *qtmp; /* argument to inverse trig function */ /* * firewall @@ -2359,22 +2361,22 @@ qavercos_or_NULL(NUMBER *q, NUMBER *epsilon) /* * calculate inverse trig function value */ - x = qsub(q, &_qone_); - res = qacos(x, epsilon); - qfree(x); + qtmp = qsub(q, &_qone_); + res = qacos(qtmp, epsilon); + qfree(qtmp); if (res == NULL) { return NULL; } /* - * return acos(x - 1) + * return inverse trigonometric result */ return res; } /* - * qavercos - inverse versed sine for NUMBER values + * qavercos - non-complex inverse versed trigonometric cosine * * This uses the formula: * @@ -2409,19 +2411,19 @@ qavercos(NUMBER *q, NUMBER *epsilon) */ res = qavercos_or_NULL(q, epsilon); if (res == NULL) { - math_error("cannot compute inverse cos for avercos"); + math_error("cannot compute inverse cosine for avercos"); not_reached(); } /* - * return acos(x - 1) + * return inverse trigonometric result */ return res; } /* - * qcovercos - coversed sine for NUMBER values + * qcovercos - coversed trigonometric cosine * * This uses the formula: * @@ -2437,7 +2439,8 @@ qavercos(NUMBER *q, NUMBER *epsilon) NUMBER * qcovercos(NUMBER *q, NUMBER *epsilon) { - NUMBER *sin, *res; + NUMBER *res; /* inverse trig value result */ + NUMBER *qtmp; /* argument to inverse trig function */ /* * firewall @@ -2454,19 +2457,19 @@ qcovercos(NUMBER *q, NUMBER *epsilon) /* * calculate trig function value */ - sin = qsin(q, epsilon); - res = qqadd(&_qone_, sin); - qfree(sin); + qtmp = qsin(q, epsilon); + res = qqadd(&_qone_, qtmp); + qfree(qtmp); /* - * return 1 + sin(x) + * return trigonometric result */ return res; } /* - * qacovercos_or_NULL - inverse coversed sine for NUMBER values + * qacovercos_or_NULL - return NULL or non-complex inverse coversed trigonometric cosine * * This uses the formula: * @@ -2491,7 +2494,7 @@ NUMBER * qacovercos_or_NULL(NUMBER *q, NUMBER *epsilon) { NUMBER *res; /* inverse trig value result */ - NUMBER *x; /* argument to inverse trig function */ + NUMBER *qtmp; /* argument to inverse trig function */ /* * firewall @@ -2508,22 +2511,22 @@ qacovercos_or_NULL(NUMBER *q, NUMBER *epsilon) /* * calculate inverse trig function value */ - x = qsub(&_qone_, q); - res = qasin(x, epsilon); - qfree(x); + qtmp = qsub(&_qone_, q); + res = qasin(qtmp, epsilon); + qfree(qtmp); if (res == NULL) { return NULL; } /* - * return asin(x - 1) + * return inverse trigonometric result */ return res; } /* - * qacovercos - inverse coversed sine for NUMBER values + * qacovercos - non-complex inverse coversed trigonometric cosine * * This uses the formula: * @@ -2558,12 +2561,612 @@ qacovercos(NUMBER *q, NUMBER *epsilon) */ res = qacovercos_or_NULL(q, epsilon); if (res == NULL) { - math_error("cannot compute inverse sin for acovercos"); + math_error("cannot compute inverse sine for acovercos"); not_reached(); } /* - * return asin(x - 1) + * return inverse trigonometric result + */ + return res; +} + + +/* + * qhaversin - half versed trigonometric sine + * + * This uses the formula: + * + * haversin(x) = versin(x) / 2 + * + * given: + * q real value to pass to the trig function + * epsilon error tolerance / precision for trig calculation + * + * returns: + * real value result of trig function on q with error epsilon + */ +NUMBER * +qhaversin(NUMBER *q, NUMBER *epsilon) +{ + NUMBER *res; /* inverse trig value result */ + NUMBER *qtmp; /* argument to inverse trig function */ + + /* + * firewall + */ + if (q == NULL) { + math_error("q is NULL for %s", __func__); + not_reached(); + } + if (check_epsilon(epsilon) == false) { + math_error("Invalid epsilon arg for %s", __func__); + not_reached(); + } + + /* + * calculate trig function value + */ + qtmp = qversin(q, epsilon); + res = qdivi(qtmp, 2); + qfree(qtmp); + + /* + * return trigonometric result + */ + return res; +} + + +/* + * qahaversin_or_NULL - return NULL or non-complex inverse half versed trigonometric sine + * + * This uses the formula: + * + * ahaversin(x) = acos(1 - 2*x) + * + * given: + * q real value to pass to the trig function + * epsilon error tolerance / precision for trig calculation + * + * returns: + * != NULL ==> real value result of trig function on q with error epsilon, + * NULL ==> trig function value cannot be expressed as a NUMBER + * + * NOTE: If this function returns NULL, consider calling the equivalent + * COMPLEX function from comfunc.c. See the help file for the + * related builtin for details. + */ +NUMBER * +qahaversin_or_NULL(NUMBER *q, NUMBER *epsilon) +{ + NUMBER *res; /* inverse trig value result */ + NUMBER *qtmp; /* argument to inverse trig function */ + NUMBER *x2; /* twice x */ + + /* + * firewall + */ + if (q == NULL) { + math_error("q is NULL for %s", __func__); + not_reached(); + } + if (check_epsilon(epsilon) == false) { + math_error("Invalid epsilon arg for %s", __func__); + not_reached(); + } + + /* + * calculate inverse trig function value + */ + x2 = qmuli(q, 2); + qtmp = qsub(&_qone_, x2); + qfree(x2); + res = qacos(qtmp, epsilon); + qfree(qtmp); + if (res == NULL) { + return NULL; + } + + /* + * return inverse trigonometric result + */ + return res; +} + + +/* + * qahaversin - non-complex inverse half versed trigonometric sine + * + * This uses the formula: + * + * ahaversin(x) = acos(1 - 2*x) + * + * given: + * q real value to pass to the trig function + * epsilon error tolerance / precision for trig calculation + * + * returns: + * real value result of trig function on q with error epsilon + */ +NUMBER * +qahaversin(NUMBER *q, NUMBER *epsilon) +{ + NUMBER *res; /* inverse trig value result */ + + /* + * firewall + */ + if (q == NULL) { + math_error("q is NULL for %s", __func__); + not_reached(); + } + if (check_epsilon(epsilon) == false) { + math_error("Invalid epsilon arg for %s", __func__); + not_reached(); + } + + /* + * calculate inverse trig function value + */ + res = qahaversin_or_NULL(q, epsilon); + if (res == NULL) { + math_error("cannot compute inverse cosine for ahaversin"); + not_reached(); + } + + /* + * return inverse trigonometric result + */ + return res; +} + + +/* + * qhacoversin - coversed trigonometric sine + * + * This uses the formula: + * + * hacoversin((x) = coversin(x) / 2 + * + * given: + * q real value to pass to the trig function + * epsilon error tolerance / precision for trig calculation + * + * returns: + * real value result of trig function on q with error epsilon + */ +NUMBER * +qhacoversin(NUMBER *q, NUMBER *epsilon) +{ + NUMBER *res; /* inverse trig value result */ + NUMBER *qtmp; /* argument to inverse trig function */ + + /* + * firewall + */ + if (q == NULL) { + math_error("q is NULL for %s", __func__); + not_reached(); + } + if (check_epsilon(epsilon) == false) { + math_error("Invalid epsilon arg for %s", __func__); + not_reached(); + } + + /* + * calculate trig function value + */ + qtmp = qcoversin(q, epsilon); + res = qdivi(qtmp, 2); + qfree(qtmp); + + /* + * return trigonometric result + */ + return res; +} + + +/* + * qahacoversin_or_NULL - return NULL or non-complex inverse coversed trigonometric sine + * + * This uses the formula: + * + * ahacoversin(x) = asin(1 - 2*x) + * + * given: + * q real value to pass to the trig function + * epsilon error tolerance / precision for trig calculation + * + * returns: + * real value result of trig function on q with error epsilon + * + * returns: + * != NULL ==> real value result of trig function on q with error epsilon, + * NULL ==> trig function value cannot be expressed as a NUMBER + * + * NOTE: If this function returns NULL, consider calling the equivalent + * COMPLEX function from comfunc.c. See the help file for the + * related builtin for details. + */ +NUMBER * +qahacoversin_or_NULL(NUMBER *q, NUMBER *epsilon) +{ + NUMBER *res; /* inverse trig value result */ + NUMBER *qtmp; /* argument to inverse trig function */ + + /* + * firewall + */ + if (q == NULL) { + math_error("q is NULL for %s", __func__); + not_reached(); + } + if (check_epsilon(epsilon) == false) { + math_error("Invalid epsilon arg for %s", __func__); + not_reached(); + } + + /* + * calculate inverse trig function value + */ + qtmp = qsub(&_qone_, q); + res = qasin(qtmp, epsilon); + qfree(qtmp); + if (res == NULL) { + return NULL; + } + + /* + * return inverse trigonometric result + */ + return res; +} + + +/* + * qahacoversin - non-complex inverse coversed trigonometric sine + * + * This uses the formula: + * + * ahacoversin(x) = asin(1 - 2*x) + * + * given: + * q real value to pass to the trig function + * epsilon error tolerance / precision for trig calculation + * + * returns: + * real value result of trig function on q with error epsilon + */ +NUMBER * +qahacoversin(NUMBER *q, NUMBER *epsilon) +{ + NUMBER *res; /* inverse trig value result */ + + /* + * firewall + */ + if (q == NULL) { + math_error("q is NULL for %s", __func__); + not_reached(); + } + if (check_epsilon(epsilon) == false) { + math_error("Invalid epsilon arg for %s", __func__); + not_reached(); + } + + /* + * calculate inverse trig function value + */ + res = qahacoversin_or_NULL(q, epsilon); + if (res == NULL) { + math_error("cannot compute inverse sine for ahacoversin"); + not_reached(); + } + + /* + * return inverse trigonometric result + */ + return res; +} + + +/* + * qhavercos - half versed trigonometric cosine + * + * This uses the formula: + * + * havercos(x) = vercos(x) / 2 + * + * given: + * q real value to pass to the trig function + * epsilon error tolerance / precision for trig calculation + * + * returns: + * real value result of trig function on q with error epsilon + */ +NUMBER * +qhavercos(NUMBER *q, NUMBER *epsilon) +{ + NUMBER *res; /* inverse trig value result */ + NUMBER *qtmp; /* argument to inverse trig function */ + + /* + * firewall + */ + if (q == NULL) { + math_error("q is NULL for %s", __func__); + not_reached(); + } + if (check_epsilon(epsilon) == false) { + math_error("Invalid epsilon arg for %s", __func__); + not_reached(); + } + + /* + * calculate trig function value + */ + qtmp = qvercos(q, epsilon); + res = qdivi(qtmp, 2); + qfree(qtmp); + + /* + * return trigonometric result + */ + return res; +} + + +/* + * qahavercos_or_NULL - return NULL or non-complex inverse half versed trigonometric cosine + * + * This uses the formula: + * + * ahavercos(x) = acos(1 - 2*x) + * + * given: + * q real value to pass to the trig function + * epsilon error tolerance / precision for trig calculation + * + * returns: + * != NULL ==> real value result of trig function on q with error epsilon, + * NULL ==> trig function value cannot be expressed as a NUMBER + * + * NOTE: If this function returns NULL, consider calling the equivalent + * COMPLEX function from comfunc.c. See the help file for the + * related builtin for details. + */ +NUMBER * +qahavercos_or_NULL(NUMBER *q, NUMBER *epsilon) +{ + NUMBER *res; /* inverse trig value result */ + NUMBER *qtmp; /* argument to inverse trig function */ + NUMBER *x2; /* twice x */ + + /* + * firewall + */ + if (q == NULL) { + math_error("q is NULL for %s", __func__); + not_reached(); + } + if (check_epsilon(epsilon) == false) { + math_error("Invalid epsilon arg for %s", __func__); + not_reached(); + } + + /* + * calculate inverse trig function value + */ + x2 = qmuli(q, 2); + qtmp = qsub(&_qone_, x2); + qfree(x2); + res = qacos(qtmp, epsilon); + qfree(qtmp); + if (res == NULL) { + return NULL; + } + + /* + * return inverse trigonometric result + */ + return res; +} + + +/* + * qahavercos - non-complex inverse half versed trigonometric cosine + * + * This uses the formula: + * + * ahavercos(x) = acos(1 - 2*x) + * + * given: + * q real value to pass to the trig function + * epsilon error tolerance / precision for trig calculation + * + * returns: + * real value result of trig function on q with error epsilon + */ +NUMBER * +qahavercos(NUMBER *q, NUMBER *epsilon) +{ + NUMBER *res; /* inverse trig value result */ + + /* + * firewall + */ + if (q == NULL) { + math_error("q is NULL for %s", __func__); + not_reached(); + } + if (check_epsilon(epsilon) == false) { + math_error("Invalid epsilon arg for %s", __func__); + not_reached(); + } + + /* + * calculate inverse trig function value + */ + res = qahaversin_or_NULL(q, epsilon); + if (res == NULL) { + math_error("cannot compute inverse cocosine for ahavercos"); + not_reached(); + } + + /* + * return inverse trigonometric result + */ + return res; +} + + +/* + * qhacovercos - half coversed trigonometric cosine + * + * This uses the formula: + * + * hacovercos((x) = covercos(x) / 2 + * + * given: + * q real value to pass to the trig function + * epsilon error tolerance / precision for trig calculation + * + * returns: + * real value result of trig function on q with error epsilon + */ +NUMBER * +qhacovercos(NUMBER *q, NUMBER *epsilon) +{ + NUMBER *res; /* inverse trig value result */ + NUMBER *qtmp; /* argument to inverse trig function */ + + /* + * firewall + */ + if (q == NULL) { + math_error("q is NULL for %s", __func__); + not_reached(); + } + if (check_epsilon(epsilon) == false) { + math_error("Invalid epsilon arg for %s", __func__); + not_reached(); + } + + /* + * calculate trig function value + */ + qtmp = qcovercos(q, epsilon); + res = qdivi(qtmp, 2); + qfree(qtmp); + + /* + * return trigonometric result + */ + return res; +} + + +/* + * qahacovercos_or_NULL - return NULL or non-complex inverse half coversed trigonometric cosine + * + * This uses the formula: + * + * ahacovercos(x) = acos(2*x - 1) + * + * given: + * q real value to pass to the trig function + * epsilon error tolerance / precision for trig calculation + * + * returns: + * real value result of trig function on q with error epsilon + * + * returns: + * != NULL ==> real value result of trig function on q with error epsilon, + * NULL ==> trig function value cannot be expressed as a NUMBER + * + * NOTE: If this function returns NULL, consider calling the equivalent + * COMPLEX function from comfunc.c. See the help file for the + * related builtin for details. + */ +NUMBER * +qahacovercos_or_NULL(NUMBER *q, NUMBER *epsilon) +{ + NUMBER *res; /* inverse trig value result */ + NUMBER *qtmp; /* argument to inverse trig function */ + + /* + * firewall + */ + if (q == NULL) { + math_error("q is NULL for %s", __func__); + not_reached(); + } + if (check_epsilon(epsilon) == false) { + math_error("Invalid epsilon arg for %s", __func__); + not_reached(); + } + + /* + * calculate inverse trig function value + */ + qtmp = qsub(&_qone_, q); + res = qacos(qtmp, epsilon); + qfree(qtmp); + if (res == NULL) { + return NULL; + } + + /* + * return inverse trigonometric result + */ + return res; +} + + +/* + * qahacovercos - non-complex inverse half coversed trigonometric cosine + * + * This uses the formula: + * + * ahacovercos(x) = acos(2*x - 1) + * + * given: + * q real value to pass to the trig function + * epsilon error tolerance / precision for trig calculation + * + * returns: + * real value result of trig function on q with error epsilon + */ +NUMBER * +qahacovercos(NUMBER *q, NUMBER *epsilon) +{ + NUMBER *res; /* inverse trig value result */ + + /* + * firewall + */ + if (q == NULL) { + math_error("q is NULL for %s", __func__); + not_reached(); + } + if (check_epsilon(epsilon) == false) { + math_error("Invalid epsilon arg for %s", __func__); + not_reached(); + } + + /* + * calculate inverse trig function value + */ + res = qahacoversin_or_NULL(q, epsilon); + if (res == NULL) { + math_error("cannot compute inverse cosine for ahacovercos"); + not_reached(); + } + + /* + * return inverse trigonometric result */ return res; }