Commit fe81cb1c authored by Over_score's avatar Over_score

Merge branch '403-implement-agm' into 'master'

Implement AGM

Closes #403

See merge request !494
parents 30ccaa11 94f912a4
No preview for this file type
......@@ -71,7 +71,7 @@ SRC_WP43S = \
stats.c statusBar.c timer.c \
wp43s.c memory.c) \
$(addprefix src/wp43s/mathematics/, \
10pow.c 2pow.c addition.c arccos.c arccosh.c arcsin.c arcsinh.c arctan.c arctanh.c \
10pow.c 2pow.c addition.c agm.c arccos.c arccosh.c arcsin.c arcsinh.c arctan.c arctanh.c \
ceil.c changeSign.c comparisonReals.c conjugate.c cos.c cosh.c cube.c cubeRoot.c \
cxToRe.c idiv.c idivr.c \
division.c exp.c expt.c factorial.c floor.c fractionalPart.c gamma.c gcd.c \
......
;*************************************************************
;*************************************************************
;** **
;** agm **
;** **
;*************************************************************
;*************************************************************
In: FD=0 FI=0 SD=0 RM=0 IM=2compl SS=4 WS=64
Func: fnAgm
;=======================================================================================================================
; Error case: reals and x<0 or y<0 or NaN or wrong data types
;=======================================================================================================================
In: AM=DEG SL=0 FI=0 RX=ShoI:"2#5" RY=LonI:"5"
Out: EC=24 FI=0 SL=0 RX=ShoI:"2#5" RY=LonI:"5"
In: AM=DEG SL=0 FI=0 RX=LonI:"5" RY=ShoI:"2#5"
Out: EC=24 FI=0 SL=0 RX=LonI:"5" RY=ShoI:"2#5"
In: AM=DEG SL=0 FI=0 RX=LonI:"-2" RY=LonI:"5"
Out: EC=1 FI=0 SL=0 RX=LonI:"-2" RY=LonI:"5"
In: AM=DEG SL=0 FI=0 RX=LonI:"5" RY=LonI:"-2"
Out: EC=1 FI=0 SL=0 RX=LonI:"5" RY=LonI:"-2"
In: AM=DEG SL=0 FI=0 RX=Re16:"-2" RY=Re16:"5"
Out: EC=1 FI=0 SL=0 RX=Re16:"-2" RY=Re16:"5"
In: AM=DEG SL=0 FI=0 RX=Re16:"5" RY=Re16:"-2"
Out: EC=1 FI=0 SL=0 RX=Re16:"5" RY=Re16:"-2"
In: AM=DEG SL=0 FI=0 RX=Re34:"-2" RY=Re34:"5"
Out: EC=1 FI=0 SL=0 RX=Re34:"-2" RY=Re34:"5"
In: AM=DEG SL=0 FI=0 RX=Re34:"5" RY=Re34:"-2"
Out: EC=1 FI=0 SL=0 RX=Re34:"5" RY=Re34:"-2"
In: AM=DEG SL=0 FI=0 RX=Re16:"NaN" RY=Re16:"2"
Out: EC=1 FI=0 SL=0 RX=Re16:"NaN" RY=Re16:"2"
In: AM=DEG SL=0 FI=0 RX=Re16:"2" RY=Re16:"NaN"
Out: EC=1 FI=0 SL=0 RX=Re16:"2" RY=Re16:"NaN"
;=======================================================================================================================
; REAL16 Domain
;=======================================================================================================================
In: AM=DEG SL=0 FI=0 NCSD=16 RX=LonI:"5" RY=LonI:"1"
Out: EC=0 FI=0 SL=1 RL=LonI:"5" RX=Re16:"2.60400819053094028869642744872503433092706525563764948774841325446"
In: AM=DEG SL=0 FI=0 RX=Re16:"1" RY=Re16:"5"
Out: EC=0 FI=0 SL=1 RL=Re16:"1" RX=Re16:"2.60400819053094028869642744872503433092706525563764948774841325446"
;=======================================================================================================================
; REAL34 Domain
;=======================================================================================================================
In: AM=DEG SL=0 FI=0 NCSD=34 RX=Re34:"1" RY=Re34:"5"
Out: EC=0 FI=0 SL=1 RL=Re34:"1" RX=Re34:"2.60400819053094028869642744872503433092706525563764948774841325446"
In: AM=DEG SL=0 FI=0 RX=LonI:"1000" RY=Re34:"1000000"
Out: EC=0 FI=0 SL=1 RL=LonI:"1000" RX=Re34:"189388.302409950875557668383899002284150743147627494539244754636078"
In: AM=DEG SL=0 FI=0 RX=Re34:"0.1" RY=Re16:"0.00000001"
Out: EC=0 FI=0 SL=1 RL=Re34:"0.1" RX=Re34:"0.00897372788032617635832975712477846117577814780155175563194533968106"
;=======================================================================================================================
; COMPLEX16 Domain
;=======================================================================================================================
In: AM=DEG SL=0 FI=0 NCSD=16 RX=LonI:"1" RY=Co16:"5 i 1"
Out: EC=0 FI=1 SL=1 RL=LonI:"1" RX=Co16:"2.61069957265504884786303570933994530375554191028876184973777910267 i 0.352968798877474658674656438810498154770461336822115251643519498869"
In: AM=DEG SL=0 FI=0 RX=Co16:"5 i 1" RY=LonI:"1"
Out: EC=0 FI=1 SL=1 RL=Co16:"5 i 1" RX=Co16:"2.61069957265504884786303570933994530375554191028876184973777910267 i 0.352968798877474658674656438810498154770461336822115251643519498869"
In: AM=DEG SL=0 FI=0 RX=Co16:"-9 i 7" RY=Co16:"5 i -3"
Out: EC=0 FI=1 SL=1 RL=Co16:"-9 i 7" RX=Co16:"2.37813658945313198208624760042329802905766383691888005657579693907 i -0.340754010801751188775511862533531607156373857455206703446959965931"
;=======================================================================================================================
; COMPLEX34 Domain
;=======================================================================================================================
In: AM=DEG SL=0 FI=0 NCSD=34 RX=LonI:"1" RY=Co34:"5 i 1"
Out: EC=0 FI=1 SL=1 RL=LonI:"1" RX=Co34:"2.61069957265504884786303570933994530375554191028876184973777910267 i 0.352968798877474658674656438810498154770461336822115251643519498869"
In: AM=DEG SL=0 FI=0 RX=Co34:"5 i 1" RY=LonI:"1"
Out: EC=0 FI=1 SL=1 RL=Co34:"5 i 1" RX=Co34:"2.61069957265504884786303570933994530375554191028876184973777910267 i 0.352968798877474658674656438810498154770461336822115251643519498869"
In: AM=DEG SL=0 FI=0 RX=Co34:"-9 i 7" RY=Co16:"5 i -3"
Out: EC=0 FI=1 SL=1 RL=Co34:"-9 i 7" RX=Co34:"2.37813658945313198208624760042329802905766383691888005657579693907 i -0.340754010801751188775511862533531607156373857455206703446959965931"
In: AM=DEG SL=0 FI=0 RX=Co16:"-9 i 7" RY=Co34:"5 i -3"
Out: EC=0 FI=1 SL=1 RL=Co16:"-9 i 7" RX=Co34:"2.37813658945313198208624760042329802905766383691888005657579693907 i -0.340754010801751188775511862533531607156373857455206703446959965931"
In: AM=DEG SL=0 FI=0 RX=Co34:"-9 i 7" RY=Co34:"5 i -3"
Out: EC=0 FI=1 SL=1 RL=Co34:"-9 i 7" RX=Co34:"2.37813658945313198208624760042329802905766383691888005657579693907 i -0.340754010801751188775511862533531607156373857455206703446959965931"
......@@ -37,6 +37,7 @@ const funcTest_t funcTestNoParam[] = {
{"fn2Pow", fn2Pow },
{"fnAdd", fnAdd },
{"fnAim", fnAim },
{"fnAgm", fnAgm },
{"fnArccos", fnArccos },
{"fnArccosh", fnArccosh },
{"fnArcsin", fnArcsin },
......
......@@ -203,3 +203,4 @@ magnitude.txt
unitVector.txt
slvq.txt
agm.txt
......@@ -139,7 +139,7 @@ const item_t indexOfItems[] = {
/* 15 */ { fnCvtAcreM2, multiply, "ac" STD_RIGHT_ARROW "m" STD_SUP_2, "acre", CAT_FNCT, SLS_ENABLED },
/* 16 */ { fnCvtAcreusM2, multiply, "ac" STD_US STD_RIGHT_ARROW "m" STD_SUP_2, "acre" STD_US, CAT_FNCT, SLS_ENABLED },
/* 17 */ { itemToBeCoded, NOPARAM, "ADV", "ADV", CAT_MENU, SLS_UNCHANGED},
/* 18 */ { itemToBeCoded, NOPARAM, "AGM", "AGM", CAT_FNCT, SLS_UNCHANGED},
/* 18 */ { fnAgm, NOPARAM, "AGM", "AGM", CAT_FNCT, SLS_ENABLED },
/* 19 */ { itemToBeCoded, NOPARAM, "AGRAPH", "AGRAPH", CAT_FNCT, SLS_UNCHANGED},
/* 20 */ { fnDisplayFormatAll, TM_VALUE, "ALL" , "ALL", CAT_FNCT, SLS_UNCHANGED},
/* 21 */ { fnConstant, 3, "a" STD_SUB_M STD_SUB_o STD_SUB_o STD_SUB_n, "a" STD_SUB_M STD_SUB_o STD_SUB_o STD_SUB_n, CAT_CNST, SLS_ENABLED },
......
/* This file is part of 43S.
*
* 43S is free software: you can redistribute it and/or modify
* it under the terms of the GNU General Public License as published by
* the Free Software Foundation, either version 3 of the License, or
* (at your option) any later version.
*
* 43S 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 General Public License for more details.
*
* You should have received a copy of the GNU General Public License
* along with 43S. If not, see <http://www.gnu.org/licenses/>.
*/
/********************************************//**
* \file agm.c
***********************************************/
#include "wp43s.h"
/********************************************//**
* \brief (c, b, a) ==> (x1, x2, r) c ==> regL
* enables stack lift and refreshes the stack
*
* \param[in] unusedParamButMandatory uint16_t
* \return void
***********************************************/
void fnAgm(uint16_t unusedParamButMandatory) {
bool_t result16=true, realInput=true;
complex39_t a, b, c;
switch(getRegisterDataType(REGISTER_X)) {
case dtLongInteger: convertLongIntegerRegisterToReal(REGISTER_X, (real_t *)&a.real, &ctxtReal39);
realZero(&a.imag);
break;
case dtReal16: real16ToReal(REGISTER_REAL16_DATA(REGISTER_X), &a.real);
realZero(&a.imag);
break;
case dtComplex16: real16ToReal(REGISTER_REAL16_DATA(REGISTER_X), &a.real);
real16ToReal(REGISTER_IMAG16_DATA(REGISTER_X), &a.imag);
realInput = false;
break;
case dtReal34: real34ToReal(REGISTER_REAL34_DATA(REGISTER_X), &a.real);
realZero(&a.imag);
result16 = false;
break;
case dtComplex34: real34ToReal(REGISTER_REAL34_DATA(REGISTER_X), &a.real);
real34ToReal(REGISTER_IMAG34_DATA(REGISTER_X), &a.imag);
realInput = false;
result16 = false;
break;
default: displayCalcErrorMessage(ERROR_INVALID_DATA_INPUT_FOR_OP, ERR_REGISTER_LINE, REGISTER_X);
#if (EXTRA_INFO_ON_CALC_ERROR == 1)
sprintf(errorMessage, "cannot calculate AGM with %s in X", getRegisterDataTypeName(REGISTER_X, true, false));
showInfoDialog("In function fnAgm:", errorMessage, NULL, NULL);
#endif
return;
}
switch(getRegisterDataType(REGISTER_Y)) {
case dtLongInteger: convertLongIntegerRegisterToReal(REGISTER_Y, (real_t *)&b.real, &ctxtReal39);
realZero(&b.imag);
break;
case dtReal16: real16ToReal(REGISTER_REAL16_DATA(REGISTER_Y), &b.real);
realZero(&b.imag);
break;
case dtComplex16: real16ToReal(REGISTER_REAL16_DATA(REGISTER_Y), &b.real);
real16ToReal(REGISTER_IMAG16_DATA(REGISTER_Y), &b.imag);
realInput = false;
break;
case dtReal34: real34ToReal(REGISTER_REAL34_DATA(REGISTER_Y), &b.real);
realZero(&b.imag);
result16 = false;
break;
case dtComplex34: real34ToReal(REGISTER_REAL34_DATA(REGISTER_Y), &b.real);
real34ToReal(REGISTER_IMAG34_DATA(REGISTER_Y), &b.imag);
realInput = false;
result16 = false;
break;
default: displayCalcErrorMessage(ERROR_INVALID_DATA_INPUT_FOR_OP, ERR_REGISTER_LINE, REGISTER_Y);
#if (EXTRA_INFO_ON_CALC_ERROR == 1)
sprintf(errorMessage, "cannot calculate AGM with %s in Y", getRegisterDataTypeName(REGISTER_Y, true, false));
showInfoDialog("In function fnAgm:", errorMessage, NULL, NULL);
#endif
return;
}
if( realIsNaN(&a.real) || realIsNaN(&a.imag)
|| realIsNaN(&b.real) || realIsNaN(&b.imag)) {
displayCalcErrorMessage(ERROR_ARG_EXCEEDS_FUNCTION_DOMAIN, ERR_REGISTER_LINE, REGISTER_X);
#if (EXTRA_INFO_ON_CALC_ERROR == 1)
showInfoDialog("In function fnAgm:", "cannot use NaN as X or Y input of AGM", NULL, NULL);
#endif
return;
}
if(realInput && (realIsNegative(&a.real) || realIsNegative(&b.real))) {
displayCalcErrorMessage(ERROR_ARG_EXCEEDS_FUNCTION_DOMAIN, ERR_REGISTER_LINE, REGISTER_X);
#if (EXTRA_INFO_ON_CALC_ERROR == 1)
showInfoDialog("In function fnAgm:", "cannot use negative X and Y as input of AGM", NULL, NULL);
#endif
return;
}
saveStack();
copySourceRegisterToDestRegister(REGISTER_X, REGISTER_L);
if(realInput) {
while(realIdenticalDigits(&a.real, &b.real) <= (result16 ? 16 : 34)) {
realAdd(&a.real, &b.real, &c.real, &ctxtReal39); // c = a + b
realMultiply(&a.real, &b.real, &b.real, &ctxtReal39); // b = a * b
realSquareRoot(&b.real, &b.real, &ctxtReal39); // b = sqrt(a * b)
realMultiply(&c.real, const_1on2, &a.real, &ctxtReal39); // a = (a + b) / 2
}
if(result16) {
reallocateRegister(REGISTER_X, dtReal16, REAL16_SIZE, AM_NONE);
realToReal16(&a.real, REGISTER_REAL16_DATA(REGISTER_X));
}
else {
reallocateRegister(REGISTER_X, dtReal34, REAL34_SIZE, AM_NONE);
realToReal34(&a.real, REGISTER_REAL34_DATA(REGISTER_X));
}
}
else { // Complex input
while(realIdenticalDigits(&a.real, &b.real) <= (result16 ? 16 : 34) || realIdenticalDigits(&a.imag, &b.imag) <= (result16 ? 16 : 34)) {
realAdd(&a.real, &b.real, &c.real, &ctxtReal39); // c = a + b real part
realAdd(&a.imag, &b.imag, &c.imag, &ctxtReal39); // c = a + b imag part
mulCo39Co39(&a, &b, &b); // b = a * b
// b = sqrt(a * b)
real39RectangularToPolar(&b.real, &b.imag, &b.real, &b.imag);
realSquareRoot(&b.real, &b.real, &ctxtReal39);
realMultiply(&b.imag, const_0_5, &b.imag, &ctxtReal39);
real39PolarToRectangular(&b.real, &b.imag, &b.real, &b.imag);
realMultiply(&c.real, const_1on2, &a.real, &ctxtReal39); // a = (a + b) / 2 real part
realMultiply(&c.imag, const_1on2, &a.imag, &ctxtReal39); // a = (a + b) / 2 imag part
}
if(result16) {
reallocateRegister(REGISTER_X, dtComplex16, COMPLEX16_SIZE, AM_NONE);
realToReal16(&a.real, REGISTER_REAL16_DATA(REGISTER_X));
realToReal16(&a.imag, REGISTER_IMAG16_DATA(REGISTER_X));
}
else {
reallocateRegister(REGISTER_X, dtComplex34, COMPLEX34_SIZE, AM_NONE);
realToReal34(&a.real, REGISTER_REAL34_DATA(REGISTER_X));
realToReal34(&a.imag, REGISTER_IMAG34_DATA(REGISTER_X));
}
}
adjustResult(REGISTER_X, true, true, REGISTER_X, -1, -1);
}
/* This file is part of 43S.
*
* 43S is free software: you can redistribute it and/or modify
* it under the terms of the GNU General Public License as published by
* the Free Software Foundation, either version 3 of the License, or
* (at your option) any later version.
*
* 43S 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 General Public License for more details.
*
* You should have received a copy of the GNU General Public License
* along with 43S. If not, see <http://www.gnu.org/licenses/>.
*/
/********************************************//**
* \file agm.h
***********************************************/
void fnAgm (uint16_t unusedButMandatoryParameter);
......@@ -304,3 +304,24 @@ bool_t real39IsAnInteger(const real_t *x) {
return real39CompareEqual(&r, const_0);
}
int16_t realIdenticalDigits(real_t *a, real_t *b) {
int16_t counter, smallest;
if(realGetExponent(a) != realGetExponent(b)) {
return 0;
}
realGetCoefficient(a, tmpStr3000);
realGetCoefficient(b, tmpStr3000 + TMP_STR_LENGTH/2);
smallest = min(a->digits, b->digits);
counter = 0;
while(counter < smallest && tmpStr3000[counter] == tmpStr3000[TMP_STR_LENGTH/2 + counter]) {
counter++;
}
return counter;
}
......@@ -18,31 +18,32 @@
* \file comparisonReals.h
***********************************************/
bool_t real16CompareAbsGreaterThan (const real16_t *number1, const real16_t *number2);
bool_t real16CompareAbsLessThan (const real16_t *number1, const real16_t *number2);
bool_t real16CompareEqual (const real16_t *number1, const real16_t *number2);
bool_t real16CompareGreaterEqual (const real16_t *number1, const real16_t *number2);
bool_t real16CompareGreaterThan (const real16_t *number1, const real16_t *number2);
bool_t real16CompareLessEqual (const real16_t *number1, const real16_t *number2);
bool_t real16CompareLessThan (const real16_t *number1, const real16_t *number2);
bool_t real16CompareAbsGreaterThan (const real16_t *number1, const real16_t *number2);
bool_t real16CompareAbsLessThan (const real16_t *number1, const real16_t *number2);
bool_t real16CompareEqual (const real16_t *number1, const real16_t *number2);
bool_t real16CompareGreaterEqual (const real16_t *number1, const real16_t *number2);
bool_t real16CompareGreaterThan (const real16_t *number1, const real16_t *number2);
bool_t real16CompareLessEqual (const real16_t *number1, const real16_t *number2);
bool_t real16CompareLessThan (const real16_t *number1, const real16_t *number2);
bool_t real34CompareAbsGreaterThan (const real34_t *number1, const real34_t *number2);
bool_t real34CompareAbsGreaterEqual(const real34_t *number1, const real34_t *number2);
bool_t real34CompareAbsLessThan (const real34_t *number1, const real34_t *number2);
bool_t real34CompareEqual (const real34_t *number1, const real34_t *number2);
bool_t real34CompareGreaterEqual (const real34_t *number1, const real34_t *number2);
bool_t real34CompareGreaterThan (const real34_t *number1, const real34_t *number2);
bool_t real34CompareLessEqual (const real34_t *number1, const real34_t *number2);
bool_t real34CompareLessThan (const real34_t *number1, const real34_t *number2);
bool_t real34CompareAbsGreaterThan (const real34_t *number1, const real34_t *number2);
bool_t real34CompareAbsGreaterEqual(const real34_t *number1, const real34_t *number2);
bool_t real34CompareAbsLessThan (const real34_t *number1, const real34_t *number2);
bool_t real34CompareEqual (const real34_t *number1, const real34_t *number2);
bool_t real34CompareGreaterEqual (const real34_t *number1, const real34_t *number2);
bool_t real34CompareGreaterThan (const real34_t *number1, const real34_t *number2);
bool_t real34CompareLessEqual (const real34_t *number1, const real34_t *number2);
bool_t real34CompareLessThan (const real34_t *number1, const real34_t *number2);
bool_t real39CompareAbsGreaterThan (const real_t *number1, const real_t *number2);
bool_t real39CompareAbsLessThan (const real_t *number1, const real_t *number2);
bool_t real39CompareEqual (const real_t *number1, const real_t *number2);
bool_t real39CompareGreaterEqual (const real_t *number1, const real_t *number2);
bool_t real39CompareGreaterThan (const real_t *number1, const real_t *number2);
bool_t real39CompareLessEqual (const real_t *number1, const real_t *number2);
bool_t real39CompareLessThan (const real_t *number1, const real_t *number2);
bool_t real39CompareAbsGreaterThan (const real_t *number1, const real_t *number2);
bool_t real39CompareAbsLessThan (const real_t *number1, const real_t *number2);
bool_t real39CompareEqual (const real_t *number1, const real_t *number2);
bool_t real39CompareGreaterEqual (const real_t *number1, const real_t *number2);
bool_t real39CompareGreaterThan (const real_t *number1, const real_t *number2);
bool_t real39CompareLessEqual (const real_t *number1, const real_t *number2);
bool_t real39CompareLessThan (const real_t *number1, const real_t *number2);
bool_t real16IsAnInteger (const real16_t *x);
bool_t real34IsAnInteger (const real34_t *x);
bool_t real39IsAnInteger (const real_t *x);
bool_t real16IsAnInteger (const real16_t *x);
bool_t real34IsAnInteger (const real34_t *x);
bool_t real39IsAnInteger (const real_t *x);
int16_t realIdenticalDigits (real_t *a, real_t *b);
......@@ -23,6 +23,7 @@
#include "2pow.h"
#include "10pow.h"
#include "addition.h"
#include "agm.h"
#include "arccos.h"
#include "arccosh.h"
#include "arcsin.h"
......
......@@ -19,5 +19,3 @@
***********************************************/
void fnSlvq (uint16_t unusedButMandatoryParameter);
void slvqCo51(void);
void slvqRe51(void);
......@@ -265,6 +265,9 @@
<Unit filename="src/testSuite/addition_time.txt">
<Option target="testSuite" />
</Unit>
<Unit filename="src/testSuite/agm.txt">
<Option target="testSuite" />
</Unit>
<Unit filename="src/testSuite/arccos.txt">
<Option target="testSuite" />
</Unit>
......@@ -982,6 +985,19 @@
<Option target="testSuite" />
<Option target="generateCatalogs" />
</Unit>
<Unit filename="src/wp43s/mathematics/agm.c">
<Option compilerVar="CC" />
<Option target="wp43s" />
<Option target="wp43s_debug" />
<Option target="generateCatalogs" />
<Option target="testSuite" />
</Unit>
<Unit filename="src/wp43s/mathematics/agm.h">
<Option target="wp43s" />
<Option target="wp43s_debug" />
<Option target="generateCatalogs" />
<Option target="testSuite" />
</Unit>
<Unit filename="src/wp43s/mathematics/arccos.c">
<Option compilerVar="CC" />
<Option target="wp43s" />
......
Markdown is supported
0%
or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment