/*  snfunc.h  by K.Tsuru  */
/*********************
SN library
Mathematical functions' sorting list.
**********************/
#ifndef SN_FUNCTIONS_H
#define SN_FUNCTIONS_H
                                                    // function ID
/*****************************************************************
Comparing two values from the position of one, return how many digits
agree each other.In the case of SDouble the exponent is not considerd
and the head zeros are not counted.
If an error ocuured return a negative value as follows.
not DRADIX           : -1
different type       : -2
not in standard form : -3
*******************************************************************/
long AgreeUpTo(const SNumber& a, const SNumber& b); // 1000

// SLong class
SLong LFact(ulong n);   // factorial n! in SLong arithmetric for reference 2001
SLong FactPF(ulong N);  // factorial N! by prime factorization method
                        // 2009 since ver 2.82
SLong TournamentProduct(const SNBlock <SLong>& a, const uint N);// 2010 since ver 2.182
SLong BothSideProduct(const SNBlock <SLong>& a, const uint N, bool sort = true);    // 2011 since ver 2.182

#ifndef STRUCT_PRIMFACTOR
#define STRUCT_PRIMFACTOR
struct primeFactor {
    int prime;
    int power;
};
#endif  // STRUCT_PRIMFACTOR
SLong PrimeFactorProduct(const SNBlock <primeFactor>& pf, int n);// 2012 since ver 2.182

/************* assistant functions for trigonometric ones **************/
/****************************************************************
 arcsin by series
Set "fe" a value of "fixedExp".For the calculation of sum
  y = a + AsinSeries(x);
in the case of |a|>>|x|, set fe = a.RdxExp() to avoid calculating useless precision.
This is a measure to decrease terms, because convergence of arcsin's series is
very slow.
If fe = INT_MAX > DRADIX_EXP_MAX, set fe = x.RdxExp().
*****************************************************************/
SDouble  AsinSeries(const SDouble& x, int fe = INT_MAX);    //3204
SDouble  AsinBSNW(const SDouble& x);
SDouble  CosRN(const SDouble& x); // x/(R^N) |x|< pi/4. The assistant function of CosS_RN(x) 3205
SDouble  SinRN(const SDouble& x); //                  ditto                                  3206
SDouble  TanRN(const SDouble& x); // 3106 Old version's name is "SDouble Tan(const SDouble& x)".
SDouble  CosS_RN(const SDouble& x); // Cos(x)
SDouble  SinS_RN(const SDouble& x); // Sin(x)

SDouble  CosDiv(const SDouble& x);      //x=short+long&small |x| < pi/4 3207(ref)
SDouble  SinDiv(const SDouble& x);      //                   3208(ref)

SDouble  CosSeries(const SDouble& x, int first = 1);    //for first = 1 return cos(x), first=0 cos(x)-1 3209
SDouble  SinSeries(const SDouble& x);                       //3210
/*--------------------------------------------------------------------------
Get a calculation method for sin(x), cos(x) or tan(x) functions.
 y = x - n*(pi/2) ( n is an integer, |y| <= pi/4 ).
Return the sign which attatch to calculation result.
If sign = 0, y has an exact value as y = -1.0, 0.0 or 1.0
Give *func a value COS_CALC, SIN_CALC or TAN_CALC, an appropriate function
is set in it.
SDecimal (*pfBSeries)(const SDouble& x), int maxExp, const char* func);
----------------------------------------------------------------------------*/
enum { COS_CALC = 1, SIN_CALC = 2, TAN_CALC = 3, COT_CALC = 4};
int GetTriCalcMethod(const SDouble& x, SDouble& y, int* func, int* quadrant = NULL);//3201
SDecimal BsinSeries(const SDouble& x);

// exponential and hyperbolic functions
enum { COSH_CALC = 1, SINH_CALC = -1};

// if first != 0 calculate exp(x), else exp(x) - 1 for small x
SDouble ExpSeries(const SDouble& x, int first = 1); //  3306
SDouble CoshSeries(const SDouble& x, int first);    //  3307
SDouble CoshSeries(const SDouble& x);//CoshSeries(x, 1) 3307
SDouble SinhSeries(const SDouble& x);               //  3308

// artan(x) =0.5*log{(1+x)/(1-x)}   |x| << 1.0
SDouble AtanhSeries(const SDouble& x);              // 3312
SDouble ExpDS(const SDouble& x); // divide and series method 3313
SDouble ExpL(const SDouble& x);  // exp(x) using Log(x) 3314
/******************************************************
ExpBSR(n, d) returns exp(n/d) with SLong n and d using binary splitting method. 3315 ver. 2.18
******************************************************/
SDouble ExpBSR(const SLong& n, const SLong& d); //  3316 ver. 2.18
/******************************************************
"ExpBSR(n, d, a, c)" returns the value exp(n/d) as a rational number a/c.
******************************************************/
void ExpBSR(const SLong& n, const SLong& d, SLong& a, SLong& c);//  3317 ver. 2.18
/***************************************************
It evaluates exp(x) using "ExpBSR(...)" above.
****************************************************/
SDouble ExpBS(const SDouble& x);                    //  3318 ver. 2.18
SNManager::SNErrorFlag ExpArgCheck(const SDouble& x, int funcID);// 3319 ver. 2.18
/*
It evaluates cos(x) and sin(x) using the binary splitting method. ver. 2.18
*/
struct TrigFuncValues {
    int quadrant;
    SDouble angle, cos, sin, tan;
    TrigFuncValues() {}
    TrigFuncValues(const SDouble& a, const SDouble& c, const SDouble& s, const SDouble& t)
        : angle(a), cos(c), sin(s), tan(t) {}
    TrigFuncValues(const TrigFuncValues& a);
    TrigFuncValues& operator=(const TrigFuncValues& a);
};
int CosSinBS(const SDouble& x, SDouble& cosVal, SDouble& sinVal); // 3320
TrigFuncValues CosSinBS(const SDouble& x, bool getTan = false); // 3321
SDouble CosBS(const SDouble& x);    //                      3322 ver. 2.18
SDouble SinBS(const SDouble& x);    //                      3323 ver. 2.18
SDouble TanBS(const SDouble& x);    //                      3324 ver. 2.18
/**************************************************************
It evaluates hyperbolic functions using binary splitting method.
***************************************************************/
void Hyperbolic(const SDouble& x, SDouble& ch, SDouble& sh);
SDouble CoshBS(const SDouble& x);   //                      3325 ver. 2.18
SDouble SinhBS(const SDouble& x);   //                      3326 ver. 2.18
SDouble TanhBS(const SDouble& x);   //                      3327 ver. 2.18
SDouble TanhBSLx(const SDouble& x); //for large x rename "TanBS(x)" since version 2.21

SDouble Log1_X(const SDouble& x); // log(1 - x) |x|<<1.0    3403
SDouble LogE(const SDouble& x);   // log(x) using Exp(x) 3404
SDouble LogAGM(const SDouble& x); // log(x) using AGM.   3405
SDouble LogNW(const SDouble& x); //log(x) using Newton's method 3406 ver. 2.18
/*--------------------------------------------------------------------------
Get a calculation method for log(x) functions.
x = X * 10^exp 
log(x) = log(X)+exp*log(10), add = exp*log(10)
It returns "true" if x = 1 or 10, set X = 0 or log(10).
----------------------------------------------------------------------------*/
bool GetLogxCalcMethod(const SDouble& x, SDouble& X, SDouble& add); //3408

// constants
/*************************************************
Entry constants
userV     : a value given by user
snV       : result, snV != NULL
pfCalcFunc: default function of calculation
maxSize   : figure of *snV
Used for Pi(), E(), Log10() functions
**************************************************/
void EntryConst(const SDouble* userV, SDouble* snV,
SDouble (*pfCalcFunc)(),  uint* maxSize); // 3500
SDouble SNE();          // The base of natural logarithm using BSE() below. 3501
SDouble SNLog10();      // log(10.0)        3502
SDouble SNPi();         // pi : 3503 Chudnovskys' formura is used since ver 2.30
//inline SDouble SNPi() { return ChudnovskysPi(); } is not used. Because Pi() needs the function pointer of SNPi().
SDouble GaussPi();      // pi(Gauss's formula)      3504 rename since ver 2.18
SDouble MachinPi();     // pi(Machin's formula)     3505 rename since ver 2.18
//RPi() is fast but out of memory occure in the FFT for the calculation 1/sqrt(2.0)
//and reciprocal.
SDouble RamanujanPi();  //pi by Ramanujan's formula 3508 rename since ver 2.30
SDouble RamanujanPi2(); //pi by Ramanujan's formula 3512 since ver 2.30
SDouble StomerPi();     //pi by Stomer's formula    3509 rename since ver 2.18
SDouble KlingenstiernaPi(); //pi by Klingenstierna  3510 rename since ver 2.18
SDouble MPi2();         // pi/2                     3511

SDouble MPi4();         // pi/4                     3550
SDouble RecLog10();     // 1/log(10.0)              3551
SDouble AGMPi();        //pi by AGM method          3552
//pi by Borweins' method improved by by Takahashi and Kanada.
SDouble BorweinsPi();       //                      3553 rename since ver 2.18
//Evaluation of e by binary splitting algorithm
SDouble BSE();          // default method of e      3554 ver. 2.17
SDouble ChudnovskysPi(); // pi by Chudnovskys' formura  3555 ver. 2.17 rename since ver 2.30

/**************************
arctan functions in DRADIX
***************************/
// arctan2 function(reference)
SDouble Datan2L(ulong num, ulong den);                  //3601
// arctan(x) =  arcsin(x/sqrt(1+x*x))
SDouble AtanA(const SDouble& x);                        //3602
// arctan(x) = arctan(y) + arctan( (x-y)/(1+x*y) )
SDouble AtanT(const SDouble& x);                        //3603
/******************************************************************
 arctan(x) by series
 arctan x = x - x^3/3 + x^5/5 -....
For the meaning of "fe", see "AsinSeries()".
Usually set fe = x.RdxExp().
******************************************************************/
SDouble AtanSeries(const SDouble& x, int fe);           //3604

/***********************************************************************
Set a constant on memory by reading from a text file.
1. fname = file name
2. Set pfunc = Log10, E or Pi. Do not check that the given value is correct.
3. Set "pfCalcFunc" a function pointer for the calculation.If sufficient length
(accuracy)in the file, NULL is ok.If pfCalcFunc == NULL and length is not enough,
a syntax error occures in "EntryConst()".
4. Set "ef" your effective figures in DRADIX. If ef = 0, set current effective
figures given by "EffFig()".
An example:
    SetConstByFile("pi.snc", Pi, NULL, 10000u);
10000*4 = 40000 digits in decimal are read from "pi.snc" file and set in "Pi()".

return 1 : enough length in file
return 0 : not enough, then needs calculation calling (*pfCalcFunc)()
*************************************************************************/
bool SetConstByFile(const char* fname,
SDouble (*pfunc)(const SDouble* c, SDouble (*pf)()),
SDouble (*pfCalcFunc)(), uint ef);                      // 3700

// SLong and SInteger class
SInteger BFact(ulong n);    // n! in BRADIX                 4102
SInteger BdFact(ulong n);// n!(partially multiplication method)BRADIX   and FFT 4106
SInteger Ipow(const SInteger& x, ulong n);  // n-th power of x(x^n) 4106
SLong combPF(ulong n, ulong k);// combination number(binomial) nCk using "PrimeFactor" 4107 since ver 2.182
SLong permL(const ulong n, const ulong r); // permutation number nPr 4108 since ver 2.182

/********** Using SDecimal arithmetric *****************/
// arctan(num/den); for num < den
SDecimal Batan2(ulong num, ulong den);  //              5101
SDecimal BE();          // base of natural logarithms 5102
SDecimal BGaussPi();        // pi by Gauss's fromula    5103 rename since ver 2.18
SDecimal BMachinPi();   // pi by Machin's formula   5104 rename since ver 2.18
// Batan2 convert to DRADIX
SDouble  Atan2L(ulong num, ulong den);  //              5105
SDecimal BStomerPi();         // pi by Stomer's formula 5106 rename since ver 2.18
SDecimal BKlingenstiernaPi();   // pi by Klingenstierna 5107

/**********************************************************
trigonometric functions in BRAIDX
Not so fast. Nearly same as Cos() and Sin().
Cos(), Sin() are faster for short decimal as x = 0.1;
BcosSeries and BsinSeries are used as assistant functions of SDouble ones.
************************************************************/
SDecimal Bcos(const SDouble& x);                        //5201
SDecimal BcosSeries(const SDouble& x);                  //5202
SDecimal Bsin(const SDouble& x);                        //5203
SDecimal BsinSeries(const SDouble& x);                  //5204
// reference function
SDecimal BexpSeries(const SDouble& x);                  //5211

/************************
inline functions group
**************************/
// abusolute value
inline SLong Labs(const SLong& x){
    if(x.Sign() < 0 ) return -x;
    return x;
}
inline SDouble Dabs(const SDouble& x){
    if(x.SNSign() < 0) return -x;
    return x;
}
inline SInteger Iabs(const SInteger& x){
    if(x.Sign() < 0 ) return -x;
    return x;
}
inline SDecimal Xabs(const SDecimal& x){
    if(x.Sign() < 0) return -x;
    return x;
}
/// over loaded function absolute value "abs(x)"
inline SDouble abs(const SDouble& x){ return Dabs(x); } // over load
inline SLong abs(const SLong& x){ return Labs(x); } // over load
inline SInteger abs(const SInteger& x) { return Iabs(x); }
inline SDecimal abs(const SDecimal& x) { return Xabs(x); }

#endif // SN_FUNCTIONS_H