#include "math.h"
#include "fp_private.h"
#if (__WANT_LONG_DOUBLE_FORMAT__ - 0L == 128L)
#include "DD.h"
static const float k1_8To52 = 6755399441055744.0; static const double kPiBy2 = 1.570796326794896619231322; static const double kPiBy2Tail1 = 6.1232339957367660e-17; static const double kPiBy2Tail2 = -1.4973849048591698e-33; static const double kPiBy2Tail3 = 5.5622711043168264e-50; static const double k2ByPi = 0.636619772367581382; static const float k2To52 = 4503599627370496.0; static const float k2To54 = 18014398509481984.0; static const double k2To45m1 = 35184372088831.0;
static const uint32_t CosCoeff[] = {
0xbFE00000, 0x00000000, 0x39685556, 0x00000000,
0x3FA55555, 0x55555555, 0x3C455555, 0x51E8C72E,
0xbF56C16C, 0x16C16C17, 0x3BF02698, 0xA05BBBBD,
0x3EFA01A0, 0x1A0169EB, 0x3B540FCD, 0xE3AE1BE0,
0xbE927E4D, 0x648A174F, 0xbB361705, 0xD6F2C25B
};
static const uint32_t SinCoeff[] = {
0xbFC55555, 0x55555555, 0xbC655555, 0x55555552,
0x3F811111, 0x11111111, 0x3C011111, 0x0DBB7C30,
0xbF2A01A0, 0x1A01A01A, 0xbB510055, 0x358E4100,
0x3EC71DE3, 0xA556A5C2, 0x3B61A976, 0xF8644324,
0xbE5AE642, 0x52F2B067, 0xbAD356A3, 0xF6278F6C
};
static const uint32_t CosCoeff16[] = {
0xBFE00000, 0x00000000, 0x395DA200, 0x00000000,
0x3FA55555, 0x55555555, 0x3C455555, 0x554FC86B,
0xBF56C16C, 0x16C16C17, 0x3BEF49FE, 0x788C8E18,
0x3EFA01A0, 0x1A01A01A, 0xBB6A022F, 0xCA316020,
0xBE927E4F, 0xB7789A12, 0xBB245398, 0xA829E8F5,
0x3E21EED8, 0xEE008FB4, 0x3AC49E38, 0xCA6157DE,
0xBDA938BA, 0xA889516E, 0x3A46CBC1, 0x4AF6E840
};
static const uint32_t SinCoeff16[] = {
0xBFC55555, 0x55555555, 0xBC655555, 0x55555554,
0x3F811111, 0x11111111, 0x3C011111, 0x10F8EB9C,
0xBF2A01A0, 0x1A01A01A, 0xBB69FF8C, 0xFBC69E00,
0x3EC71DE3, 0xA556C733, 0x3B6EFC3D, 0xB602439A,
0xBE5AE645, 0x67F53972, 0xBAF653B7, 0xA917A91D,
0x3DE61246, 0x1048C848, 0x3A7CFCE5, 0x9E3CEA38,
0xBD6AE6FF, 0xCFCF788E, 0x39DA11EB, 0x64C11650
};
static const uint32_t TrigTable[] = {
0x3FB80000, 0x00054E6C, 0x3FB7F701, 0x032A9959, 0xBBAD11F0, 0x5495F5BE, 0x3FEFDC06, 0xBF7E5BB6, 0xBBDCB529, 0x1917F935, 0x3FBC0000, 0x000A5F0A,
0x3FBBF1B7, 0x8572884A, 0x3BBAADFC, 0xF7E52ACA,
0x3FEFCF0C, 0x800E7577, 0xBBE4670B, 0xB43517A0,
0x3FC00000, 0x000F1278,
0x3FBFEAAE, 0xEEA4D6F0, 0xBB9D8153, 0x7343A887,
0x3FEFC015, 0x527CE390, 0x3BE16998, 0xA0C1E22D,
0x3FC20000, 0x00052790,
0x3FC1F0D3, 0xD7B4E938, 0xBBB94CF6, 0x02B38414,
0x3FEFAF22, 0x263C1D96, 0xBBDE4F5F, 0x508C4D7E,
0x3FC40000, 0x000143F6,
0x3FC3EB31, 0x2C5EA6CF, 0x3BBF0AD1, 0xD6EF035D,
0x3FEF9C34, 0x0A7CB78E, 0xBBE87189, 0x89E0691A,
0x3FC60000, 0x0007E6FF,
0x3FC5E44F, 0xCFA8F023, 0xBBAF148F, 0x61A85B8E,
0x3FEF874C, 0x2E1E9676, 0x3BC38F0F, 0xAB6185B5,
0x3FC80000, 0x000E77F8,
0x3FC7DC10, 0x2FC929C2, 0xBBC3C05A, 0x114DA314,
0x3FEF706B, 0xDF9E2180, 0xBBD8BEA1, 0x93F3F8C7,
0x3FCA0000, 0x00102F36,
0x3FC9D252, 0xD0DE9D1A, 0x3B8BF04E, 0xADF4B6FE,
0x3FEF5794, 0x8CFE96A3, 0xBBE58FFE, 0x7BB60718,
0x3FCC0000, 0x0009CE78,
0x3FCBC6F8, 0x4EE5F43E, 0xBBA7E35C, 0xC27E82D6,
0x3FEF3CC7, 0xC3B3493B, 0xBBCC7CA2, 0xE7033121,
0x3FCE0000, 0x00150434,
0x3FCDB9E1, 0x5FCA16EB, 0x3BC107C0, 0xD41B55C5,
0x3FEF2007, 0x30852C41, 0x3BE70655, 0x75B2D6FB,
0x3FD00000, 0x0061B811,
0x3FCFAAEE, 0xD5B07238, 0xBBCB82BA, 0x01ECEAD1,
0x3FEF0154, 0x9F71D818, 0x3BE95043, 0xCC0B6426,
0x3FD10000, 0x0020111B,
0x3FD0CD00, 0xCF125569, 0xBBC2D015, 0xDA2CFEAE,
0x3FEEE0B1, 0xFBBCBB9D, 0xBBE9AA1F, 0xF00E2374,
0x3FD20000, 0x000C1AD3,
0x3FD1C37D, 0x64D25988, 0x3BC4B236, 0xE14BCAA2,
0x3FEEBE21, 0x4F75419B, 0x3BE52928, 0x2D338F16,
0x3FD30000, 0x000E7F89,
0x3FD2B8DD, 0xC44C91CC, 0xBBC5707F, 0x432D5073,
0x3FEE99A4, 0xC3A5AEA3, 0x3BD23919, 0x2EA440E3,
0x3FD40000, 0x0017A801,
0x3FD3AD12, 0x9780568B, 0xBBC59FD1, 0x16DF09C5,
0x3FEE733E, 0xA0159A53, 0xBBB12BF5, 0xDB2DCDCC,
0x3FD50000, 0x0039F7ED,
0x3FD4A00C, 0x9B461D51, 0xBBCDFE0E, 0x1445AAA9,
0x3FEE4AF1, 0x4B20ED63, 0x3BD374F6, 0x9F86063A,
0x3FD60000, 0x0000A1AA,
0x3FD591BC, 0x9FA38DCC, 0xBBCD07DB, 0x418C1954,
0x3FEE20BF, 0x49ACBB83, 0x3B864FA4, 0xEB634A4E,
0x3FD70000, 0x00011BC1,
0x3FD68213, 0x8A39E197, 0xBBBA6968, 0x3C12FB6D,
0x3FEDF4AB, 0x3EBD5578, 0xBBD6766E, 0x746C48E7,
0x3FD80000, 0x0008ED50,
0x3FD77102, 0x557E9094, 0x3BB488F1, 0x4A8DF658,
0x3FEDC6B7, 0xEB97B68C, 0x3BDCA325, 0xE2462D16,
0x3FD90000, 0x00861E6C,
0x3FD85E7A, 0x12FE6D57, 0x3BA65076, 0x5FCB7841,
0x3FED96E8, 0x2F58212A, 0xBBEC65EC, 0x2CEB0E13,
0x3FDA0000, 0x00021040,
0x3FD94A6B, 0xE9F72C06, 0x3BA04A4B, 0xE0FC13D6,
0x3FED653F, 0x073DD7E0, 0x3BEC1A4A, 0xF1A7037E,
0x3FDB0000, 0x005EF712,
0x3FDA34C9, 0x1D1BB055, 0x3BC42083, 0xD43A5519,
0x3FED31BF, 0x8D7A0AAD, 0xBB9D5F11, 0x5DAECA9B,
0x3FDC0000, 0x000C075C,
0x3FDB1D83, 0x053CFB6A, 0xBBDA72E9, 0x9A1C2FBE,
0x3FECFC6C, 0xFA50214C, 0xBBCB5164, 0xBC69A105,
0x3FDD0000, 0x0056E33E,
0x3FDC048B, 0x17FF5F2B, 0x3BD8AF58, 0xF8525AD2,
0x3FECC54A, 0xA29F9263, 0x3BD22EA3, 0xDBB55D0D,
0x3FDE0000, 0x005CCC54,
0x3FDCE9D2, 0xE4276EF1, 0xBBAAEDF9, 0xDD5DC899,
0x3FEC8C5B, 0xF8B9244D, 0xBBC2D9C6, 0xCE4671F4,
0x3FDF0000, 0x00161118,
0x3FDDCD4C, 0x154623DC, 0xBBD97F5B, 0x8DF7FEBF,
0x3FEC51A4, 0x8B85F41A, 0x3BC354E8, 0x86ABC9CE,
0x3FE00000, 0x0013AF20,
0x3FDEAEE8, 0x746D926F, 0x3BD2B7DF, 0xA4EEC7DD,
0x3FEC1528, 0x06520D6D, 0x3BE9EE87, 0x2601172E,
0x3FE08000, 0x0046DF51,
0x3FDF8E99, 0xE7E60D68, 0x3BAF1859, 0x1DD2BAB4,
0x3FEBD6EA, 0x30DFA2E0, 0x3BE95328, 0x8BAB7ACB,
0x3FE10000, 0x0002BED1,
0x3FE03629, 0x39C8F748, 0x3BCFAF90, 0x9F47E8F6,
0x3FEB96EE, 0xEF572000, 0x3BD89C91, 0x79738697,
0x3FE18000, 0x00404634,
0x3FE0A402, 0x1ED4F66C, 0x3BD4E9FD, 0x2527C4F7,
0x3FEB553A, 0x40EAA3C8, 0xBBE9F185, 0x9F0B73D8,
0x3FE20000, 0x0005ABAA,
0x3FE110D0, 0xC4BB683B, 0xBBB78FD4, 0x4D302311,
0x3FEB11D0, 0x415F9E99, 0x3BDAE453, 0xDBCB3704,
0x3FE28000, 0x002D02D5,
0x3FE17C8E, 0x5F549FEE, 0x3B6368B3, 0xE72D341D,
0x3FEACCB5, 0x26DE0531, 0xBB8CE182, 0xC27CD7B8,
0x3FE30000, 0x0002B1A8,
0x3FE1E734, 0x323892EB, 0x3BC9A07C, 0xF3B89184,
0x3FEA85ED, 0x43725E55, 0x3BE2C483, 0xBE536F66,
0x3FE38000, 0x000730F9,
0x3FE250BB, 0x937E7157, 0xBBC3FA7C, 0xADF21ED4,
0x3FEA3D7D, 0x034EA01E, 0xBBD5786E, 0x86ABB722,
0x3FE40000, 0x00076EEF,
0x3FE2B91D, 0xEA8E4953, 0x3BE1CB0F, 0x13421234,
0x3FE9F368, 0xED8CD61E, 0xBBC0A358, 0x32653793,
0x3FE48000, 0x00119ECA,
0x3FE32054, 0xB156DCB6, 0xBBBF1C5F, 0x95EB1A26,
0x3FE9A7B5, 0xA35FDCFF, 0x3BE6C143, 0x88B5E209,
0x3FE50000, 0x000045B6,
0x3FE38659, 0x74565F66, 0x3BD5144E, 0x65AA2BBD,
0x3FE95A67, 0xE00C8774, 0x3BC1A72A, 0xE0D050C3,
0x3FE58000, 0x00A4E7B8,
0x3FE3EB25, 0xD3EDE59C, 0x3BE1F420, 0x5560DD7C,
0x3FE90B84, 0x77E73599, 0x3B9250E3, 0x5B17FE76,
0x3FE60000, 0x000049F4,
0x3FE44EB3, 0x81CF7192, 0xBBD5B364, 0x7D4D895F,
0x3FE8BB10, 0x5A5D9A12, 0xBBD2BB98, 0xF4F2C91C,
0x3FE68000, 0x0008C197,
0x3FE4B0FC, 0x46B16552, 0x3BC863B5, 0xB735A540,
0x3FE86910, 0x8D71FD5A, 0x3BE5C4E3, 0x402BA3EE,
0x3FE70000, 0x00340EE5,
0x3FE511F9, 0xFDA26352, 0x3BCDC369, 0x16E2A8DB,
0x3FE8158A, 0x316F2658, 0x3BD69C38, 0x4E2E4085,
0x3FE78000, 0x005756FD,
0x3FE571A6, 0x96AE2DA7, 0x3BE24CA0, 0xA2E35E02,
0x3FE7C082, 0x7ECF5E07, 0xBBCBD891, 0x2A83004E,
0x3FE80000, 0x004601DC,
0x3FE5CFFC, 0x16F2C847, 0x3BCA6614, 0x963A9B5C,
0x3FE769FE, 0xC62568E3, 0x3BE746CD, 0xBABDB543,
0x3FE88000, 0x002EF948,
0x3FE62CF4, 0x99438A17, 0x3BA3ED99, 0x63704218,
0x3FE71204, 0x6F86E919, 0xBBE0671D, 0x4098C629,
0x3FE90000, 0x0022F620,
0x3FE6888A, 0x4E2C1E13, 0xBBC0D332, 0x374F8A22,
0x3FE6B898, 0xFA865CFA, 0xBBDE0573, 0x1D87548D,
0x3FE98000, 0x00883EB0,
0x3FE6E2B7, 0x7C9FF82A, 0x3BE241BC, 0xA153A6CC,
0x3FE65DC1, 0xFD8A1C59, 0xBBDBAD2A, 0x038DFB75
};
static const uint32_t TanCoeff16[] = {
0x3FD55555, 0x55555555, 0x3C755555, 0x5555554C,
0x3FC11111, 0x11111111, 0x3C411111, 0x11327FD4,
0x3FABA1BA, 0x1BA1BA1C, 0xBC47917B, 0x8730ECE8,
0x3F9664F4, 0x882C10FA, 0xBC08946C, 0x04337BC0,
0x3F8226E3, 0x55E6C239, 0xBC16F3F8, 0x5EB82A9C,
0x3F6D6D3D, 0x0E159BD4, 0xBC0BF02E, 0x91A39520,
0x3F57DA36, 0x44E2BF4A, 0x3BF4CC81, 0xFDF36B54,
0x3F435582, 0xB4E9E8D9, 0xBBEAAB03, 0xE323ED60,
0x3F2F5712, 0xDCD1E82C, 0x3BA67F16, 0xB4DC7710,
0x3F19C9D2, 0x46D3C681, 0x3BA2B1F2, 0xBC709028
};
static const uint32_t CotCoeff16[] = {
0x3FD55555, 0x55555555, 0x3C755555, 0x55555550,
0x3F96C16C, 0x16C16C17, 0xBC2F49F4, 0x9F46AA0E,
0x3F61566A, 0xBC011567, 0xBC050FFB, 0x0703E386,
0x3F2BBD77, 0x9334EF0B, 0xBBC541A3, 0x5FEBE630,
0x3EF66A8F, 0x2BF70EDF, 0xBB960EEE, 0x58853E10,
0x3EC22805, 0xD642B849, 0x3B49A324, 0x92C33F58,
0x3E8D6DB2, 0xD58EAF82, 0x3B103D9B, 0xCE4CCBB8,
0x3E57DA1C, 0x0FF6D02F, 0x3AFBDB54, 0xE9F15B6C,
0x3E2396A6, 0x6CC1B8E2, 0x3AC04FB7, 0x0AEC86D2
};
static const uint32_t TanCoeff[] = {
0x3ff00000, 0x00000000, 0x00000000, 0x00000000,
0x3FD55555, 0x55555555, 0x3C755555, 0x55555549,
0x3FC11111, 0x11111111, 0x3C411111, 0x17771B40,
0x3FABA1BA, 0x1BA1BA1C, 0xBC47A471, 0x53F3999C,
0x3F9664F4, 0x882C119E, 0x3C2A931E, 0x61CFDED8,
0x3F8226E3, 0x5546FD80, 0xBC22700A, 0xDC5E76C2,
0x3F6D6DCB, 0xC6BCC1FD, 0x3C0EB3DD, 0x1EBC7C6C
};
static const uint32_t TanTable[] = {
0x3FB80000, 0x00095B76,
0x3FB81210, 0x420B0DDD, 0xBBADBC24, 0xD03EE61B,
0x40254552, 0xEE62D182, 0xBC23A778, 0xD4D4C82D,
0x3FBC0000, 0x00F6D8DC,
0x3FBC1CB8, 0x85A84FC9, 0xBBA0734F, 0xD6E646D7,
0x40223676, 0x163ABE75, 0x3C19C932, 0xFC451F09,
0x3FC00000, 0x002177FC,
0x3FC01577, 0xAF3710E9, 0xBBBB1D68, 0x25554EB9,
0x401FD549, 0xF004A93B, 0x3BF1B140, 0x750F85AD,
0x3FC20000, 0x0000C2F1,
0x3FC21E9E, 0x0175E475, 0x3BBD91B3, 0xA31924A6,
0x401C41B6, 0xE169DC8D, 0x3C057A49, 0x2ABEE9CF,
0x3FC40000, 0x004D703E,
0x3FC42A13, 0xDFCB15A0, 0xBB94ED3B, 0x3620DA54,
0x4019642D, 0xFDBA3F8F, 0x3BBB47D8, 0x80CDA3C6,
0x3FC60000, 0x000CAC24,
0x3FC6381F, 0x200F2AED, 0x3BBF3967, 0x0381AC58,
0x40170B09, 0x205DE76C, 0xBBFEFF3A, 0x14B71AB8,
0x3FC80000, 0x00473296,
0x3FC84906, 0xF15CE818, 0x3BB03E11, 0x36E0F2D3,
0x4015152E, 0xCDA71D37, 0xBBE77628, 0xDB83D025,
0x3FCA0000, 0x00B2D7A8,
0x3FCA5D14, 0x0081E4D7, 0xBBC0188C, 0x9C6AB2DE,
0x40136BB4, 0xBA046E82, 0xBC0F80BA, 0xA34DE4A8,
0x3FCC0000, 0x0066CDE2,
0x3FCC7490, 0xA23DC3B4, 0xBBCB362F, 0xA69FE4DA,
0x4011FE3C, 0xA58F0345, 0xBBFCB31C, 0x2F2152B4,
0x3FCE0000, 0x00252978,
0x3FCE8FC9, 0x01177F3A, 0xBBA5E8BA, 0x10809B66,
0x4010C0C5, 0xABFB3524, 0xBBFA678B, 0xCA50DEEC,
0x3FD00000, 0x00638D21,
0x3FD05785, 0xA4A65715, 0x3BC29365, 0x7FDA2145,
0x400F549E, 0x3203FCA2, 0x3BE2971C, 0x11C2A8A3,
0x3FD10000, 0x002C35A8,
0x3FD16953, 0xEACF2DA7, 0x3BC40A91, 0xBA0FE60F,
0x400D67EC, 0xF327EFBB, 0x3C05D1EC, 0x4B5FCFCB,
0x3FD20000, 0x001C3AB2,
0x3FD27D78, 0xB42A0CE9, 0x3BCCED21, 0x016F79E8,
0x400BB0C1, 0xF1319BD2, 0xBC051F68, 0x2CBAE56E,
0x3FD30000, 0x000F9C12,
0x3FD3941E, 0xADA8C533, 0xBBBB6B83, 0x08722CFA,
0x400A26A8, 0xA4D1FD5D, 0xBBF11692, 0xFACCF74A,
0x3FD40000, 0x000B6423,
0x3FD4AD71, 0xED5E62C2, 0xBBAD5590, 0x3A1538DA,
0x4008C2DD, 0x5F133FBD, 0xBBDDCA83, 0x1B096361,
0x3FD50000, 0x0017ECB5,
0x3FD5C9A0, 0x105DB3D0, 0x3BCE9A3C, 0xD13B3D17,
0x40077FE6, 0x3A7531D9, 0x3BE036AE, 0x48AA7EC4,
0x3FD60000, 0x00785229,
0x3FD6E8D8, 0x5AEC50EC, 0xBBC83B79, 0x6A9E4244,
0x40065948, 0x26E6F5E4, 0xBBFE36C3, 0x52F1803B,
0x3FD70000, 0x001C8263,
0x3FD80B4B, 0xD8D44659, 0xBBB499F2, 0xDC83C2D8,
0x40054B4F, 0x852854B9, 0x3BFB7522, 0xCB8BEDBE,
0x3FD80000, 0x0092A4EC,
0x3FD9312D, 0x86455645, 0x3B93DC76, 0x9413BABA,
0x400452E6, 0x95D5E9C0, 0xBBDBF037, 0x35429FA9,
0x3FD90000, 0x00353578,
0x3FDA5AB2, 0x70323F27, 0xBBD5E11A, 0x4F5FE60A,
0x40036D75, 0xEB84FD4A, 0xBBE3A384, 0xFD52C107,
0x3FDA0000, 0x002177A0,
0x3FDB8811, 0xE4F7B2DE, 0x3BC5B70A, 0x7D15AD3F,
0x400298CC, 0x1A74BDDF, 0x3BD41900, 0x83E8A76D,
0x3FDB0000, 0x0071F1F6,
0x3FDCB985, 0x9CFB26A0, 0x3BC8BA6E, 0xF094701B,
0x4001D30A, 0xD7B1694B, 0xBBFECEAB, 0x555C2BC9,
0x3FDC0000, 0x0071BD74,
0x3FDDEF49, 0xEB35D726, 0xBBCCC6A5, 0xDDB5912C,
0x40011A98, 0x2075779C, 0xBBF95852, 0x2DBF4C5B,
0x3FDD0000, 0x00A18446,
0x3FDF299D, 0xF3CB9E1D, 0xBBDAB47B, 0x6F4E4DBE,
0x40006E12, 0x72BBF45D, 0x3BE9A288, 0x2C8CD3D9,
0x3FDE0000, 0x0009F5E0,
0x3FE03461, 0xF090AA3D, 0xBBB6411F, 0x2E55E33D,
0x3FFF988E, 0xC8064412, 0xBBE51DE0, 0x7EC2ADBC,
0x3FDF0000, 0x00263EAE,
0x3FE0D680, 0x92D62104, 0x3BDD6620, 0x8142D974,
0x3FFE6858, 0x0C88DF44, 0xBBDD8335, 0xD4838950,
0x3FE00000, 0x002FD5FE,
0x3FE17B4F, 0x5C31640E, 0xBBCF7FD9, 0x4387A68E,
0x3FFD49AD, 0x7DDFB168, 0x3BF71CBF, 0x13808D41,
0x3FE08000, 0x00760E88,
0x3FE222F4, 0xAFFFC590, 0x3BDD3CE8, 0x22236C45,
0x3FFC3AF4, 0x70DCC003, 0xBBFBC120, 0xA2F987AC,
0x3FE10000, 0x000F8806,
0x3FE2CD98, 0xFEB5905C, 0x3BE169C8, 0x3A579A2C,
0x3FFB3AC2, 0x7605E88C, 0xBBC11C38, 0xB96450D3,
0x3FE18000, 0x002B5E5D,
0x3FE37B66, 0xF43D000E, 0x3BD2830E, 0xBAC9646C,
0x3FFA47D6, 0x6F7CCB9D, 0xBBEB41B3, 0x171B37E1,
0x3FE20000, 0x000A3C10,
0x3FE42C8B, 0xA0F7A0E3, 0x3BE38DF3, 0x66A12B85,
0x3FF96112, 0xDACBFCBC, 0xBBEFEF64, 0xEF1EA7D1,
0x3FE28000, 0x00595422,
0x3FE4E136, 0xB0CFA783, 0xBBBA398C, 0xEE061521,
0x3FF88578, 0xFAA0A3DE, 0x3BDA1845, 0x4AE391EC,
0x3FE30000, 0x0026966D,
0x3FE5999A, 0x9E477C73, 0xBBCE4AC5, 0x97E84259,
0x3FF7B424, 0xCEF83049, 0xBBF35B98, 0x9FE75083,
0x3FE38000, 0x0015B3FD,
0x3FE655EC, 0xF397B5BE, 0xBBE51A33, 0x112B4404,
0x3FF6EC49, 0xA1E2DBE5, 0x3BB86EE4, 0x93E45292,
0x3FE40000, 0x00146EE3,
0x3FE71666, 0x89F330B5, 0xBBE35378, 0x5A4867EE,
0x3FF62D2F, 0x215A3E9B, 0x3BEA0C9E, 0xF65CC2CF,
0x3FE48000, 0x000F229A,
0x3FE7DB43, 0xD3A2EEDF, 0xBBE0A454, 0xC31E28C0,
0x3FF5762E, 0xE1176AE6, 0xBBE1C3BD, 0x493469EA,
0x3FE50000, 0x002E680F,
0x3FE8A4C5, 0x2CF1484A, 0xBB9432B2, 0x063946DC,
0x3FF4C6B2, 0x384B4991, 0x3BD8E9D0, 0xF949DFE6,
0x3FE58000, 0x00672096,
0x3FE9732F, 0x3459A187, 0x3BDDADD7, 0x84F3E87C,
0x3FF41E30, 0x6C72DBD4, 0xBBEC67AE, 0xAE1EBB1F,
0x3FE60000, 0x00B6C23F,
0x3FEA46CB, 0x2D189DAF, 0x3BB739F4, 0x75A8271C,
0x3FF37C2D, 0x1B3BE6EB, 0xBBC85203, 0x8CECB4C0,
0x3FE68000, 0x0028A725,
0x3FEB1FE7, 0x6A3CF1FF, 0xBBEA7D29, 0x7EE6B9D1,
0x3FF2E036, 0xDC04F43A, 0x3BED69F5, 0xAF37688F,
0x3FE70000, 0x004CA480,
0x3FEBFED7, 0xCD2DB8E2, 0xBBCD123B, 0x0537EE2D,
0x3FF249E6, 0x09D6492E, 0xBBD38032, 0x80A8D415,
0x3FE78000, 0x000FFDA8,
0x3FECE3F6, 0x42FE6158, 0x3BD569F8, 0x793799F7,
0x3FF1B8DB, 0xBED5846E, 0x3BD9797E, 0x39EAD48F,
0x3FE80000, 0x001E1B0A,
0x3FEDCFA3, 0x61492AAF, 0xBBDA3DEC, 0xB61278D5,
0x3FF12CC0, 0xE58AA7DC, 0x3BB3088E, 0x1168567F,
0x3FE88000, 0x000B208A,
0x3FEEC247, 0x07D4CEEF, 0xBBE3FF76, 0x1F807DFF,
0x3FF0A545, 0x6F3C80B0, 0xBBCF3891, 0xE44E8437,
0x3FE90000, 0x000657F8,
0x3FEFBC51, 0x1E0226B6, 0xBBEC6FF5, 0x5F02EDF6,
0x3FF0221F, 0x9DADF743, 0xBBD42807, 0xFD02BF0B,
0x3FE98000, 0x00655AE2,
0x3FF05F1D, 0x31752EE7, 0xBBD8BE76, 0xB0E002E9,
0x3FEF4616, 0xC9018897, 0xBBDB4E55, 0xDEAF57E4
};
struct TrigTblEntry {
double arg;
double sinHead;
double sinTail;
double cosHead;
double cosTail;
} TrigTblEntry;
struct CosCoeffEntry {
double Head;
double Tail;
} CosCoeffEntry;
struct SinCoeffEntry {
double Head;
double Tail;
} SinCoeffEntry;
struct CosCoeff16Entry {
double Head;
double Tail;
} CosCoeff16Entry;
struct SinCoeff16Entry {
double Head;
double Tail;
} SinCoeff16Entry;
struct TanTblEntry {
double arg;
double TanHead;
double TanTail;
double CotHead;
double CotTail;
} TanTblEntry;
struct TanCoeffEntry {
double Head;
double Tail;
} TanCoeffEntry;
struct CotCoeff16Entry {
double Head;
double Tail;
} CotCoeff16Entry;
struct TanCoeff16Entry {
double Head;
double Tail;
} TanCoeff16Entry;
static void argReduce( double xHead, double xTail, double *dLeft, int *M, double *redHead, double *redTail )
{
double amid, alow, carry, high, q, q1;
double b, blow, b1, c, clow, c1, c2, c3, c4;
double d, dlow, d1, d2, d3, d4, d5, e1, e2, e3;
Double DeN;
amid = 0.0;
alow = 0.0;
while (fabs(xHead) > k2To45m1) { q = xHead*k2ByPi;
if (fabs(q) < k2To54) { if (fabs(q) < k2To52) { if (q > 0)
q = __FADD( q, k2To52 ) - k2To52; else
q = __FSUB( q, k2To52 ) + k2To52;
} q1 = q*0.25; if (q > 0)
q1 = __FADD( q1, k2To52 ) - k2To52; else
q1 = __FSUB( q1, k2To52 ) + k2To52; DeN.f = __FNMSUB( 4.0, q1, q ) + k1_8To52;
*M = *M + DeN.i[1];
}
high = __FNMSUB( q, kPiBy2, xHead ); b = q*kPiBy2Tail1;
blow = __FNMSUB( q, kPiBy2Tail1, b ); c = q*kPiBy2Tail2;
clow = __FNMSUB( q, kPiBy2Tail2, c ); d = q*kPiBy2Tail3;
dlow = __FNMSUB( q, kPiBy2Tail3, d ); b1 = __FADD( xTail, high );
if (fabs(xTail) > fabs(high)) c1 = xTail - b1 + high;
else
c1 = high - b1 + xTail;
xHead = __FSUB( b1, b );
if (fabs(b1) > fabs(b))
c2 = b1 - xHead - b + c1;
else
c2 = b1 - (xHead+b) + c1;
c3 = __FADD( amid, blow);
if (fabs(amid) > fabs(blow))
d1 = amid - c3 + blow;
else
d1 = blow - c3 + amid;
c4 = __FSUB( c2, c );
if (fabs(c2) > fabs(c))
d2 = c2 - c4 - c + d1;
else
d2 = c2 - __FADD( c4, c ) + d1;
xTail = __FADD( c3, c4 );
if (fabs(c3) > fabs(c4))
d3 = c3 - xTail + c4 + d2;
else
d3 = c4 - xTail + c3 + d2;
d4 = __FADD( alow, clow );
if (fabs(alow) > fabs(clow))
e1 = alow - d4 + clow;
else
e1 = clow - d4 + alow;
d5 = __FSUB( d3, d );
if (fabs(d3) > fabs(d))
e2 = d3 - d5 - d + e1;
else
e2 = d3 - (d5+d) + e1;
amid = __FADD( d4, d5 );
if (fabs(d4) > fabs(d5))
e3 = d4 - amid + d5 + e2;
else
e3 = d5 - amid + d4 + e2;
alow = e3 + dlow;
}
high = __FADD( xHead, xTail ); carry = xHead - high + xTail;
xTail = __FADD( carry, __FADD( amid, alow ) );
xHead = high;
*dLeft = __FSUB( carry, xTail ) + __FADD( amid, alow );
*redHead = xHead;
*redTail = xTail;
}
static void SinCosCore(int M, double xHead, double xTail, double dleft, double *yHead, double *yTail )
{
double intquo;
double temp, temps, tempc;
double b, c, ca, d;
double redarg, redarg1, res, reslow, resmid, resbot, restop, resin, resextra;
double ressup, restiny, resinf;
double p1, p2, p3, p4, p5, p6, p6a, p7, p8, p9, p10;
double p11, p12, p13, p14, p15, p16, p17, p18, p19, p20;
double p21, p22, p23, p23a, p23b, p24, p25, p26, p27, p27a, p28, p29, p30, p31;
double arg, arg2, arg2a, argsq, argsq2, argmid, argnew, argres;
double sin, sin2, sint, sinl, sinargtiny;
double sintab, sintab2, arg2b;
double cos, cos2, cosl, cost;
double prods, prodsl, prodc, prodcl;
double blow, carry, cerr;
int i, index;
Double DeN;
struct TrigTblEntry *TrigPtr = (struct TrigTblEntry *)TrigTable - 6;
struct CosCoeffEntry *CosPtr = (struct CosCoeffEntry *) CosCoeff;
struct SinCoeffEntry *SinPtr = (struct SinCoeffEntry *) SinCoeff;
struct CosCoeff16Entry *Cos16Ptr = (struct CosCoeff16Entry *) CosCoeff16;
struct SinCoeff16Entry *Sin16Ptr = (struct SinCoeff16Entry *) SinCoeff16;
DeN.f = __FMADD( xHead, k2ByPi, k1_8To52 );
intquo = DeN.f - k1_8To52; redarg = __FNMSUB( intquo, kPiBy2, xHead );
M = M + DeN.i[1];
b = kPiBy2Tail1*intquo; blow = __FMSUB( kPiBy2Tail1, intquo, b ); redarg1 = __FSUB( xTail, b ); DeN.f = __FMADD( fabs(redarg), 64.0, k1_8To52 ); index = DeN.i[1];
ca = kPiBy2Tail2*intquo;
if (fabs(b) > fabs(xTail))
carry = xTail - (b+redarg1); else
carry = (xTail-redarg1) - b;
if (! (index < 6)) {
if (redarg > 0.0) {
redarg = redarg - TrigPtr[index].arg ; sintab = TrigPtr[index].sinHead ; sintab2 = TrigPtr[index].sinTail ;
}
else { redarg = redarg + TrigPtr[index].arg;
sintab =- TrigPtr[index].sinHead ; sintab2 =- TrigPtr[index].sinTail ;
}
}
arg = __FADD( redarg, redarg1 ); arg2a = (redarg - arg + redarg1);
arg2b = __FSUB( arg2a, blow ); c = __FSUB( carry, ca );
argmid = __FADD( arg2b , c );
argnew = __FADD( arg, argmid );
argres = arg - argnew + argmid;
arg = argnew;
cerr = (carry-c) - ca;
if (fabs(arg2b) > fabs(c))
d = arg2b - argmid + c + dleft;
else
d = c - argmid + arg2b + dleft;
d = __FNMSUB( kPiBy2Tail3, intquo, d + cerr + __FNMSUB( kPiBy2Tail2, intquo, ca ) );
arg2 = __FADD( argres, d ); temp = 2.0*arg*arg2;
argsq = __FMADD( arg, arg, temp ); argsq2 = __FMADD( arg2, arg2, __FMSUB( arg, arg, argsq ) + temp );
if (index < 6) {
if ((M & 1) == 0) { sin = Sin16Ptr[6].Head ; sinl = Sin16Ptr[6].Tail ; sinargtiny = argres - arg2 + d;
for(i = 5; i > -1; i--) {
temps = __FMADD( sin, argsq2, sinl*argsq ); prods = __FMADD( sin, argsq, temps ); prodsl = __FMSUB( sin, argsq, prods ) + temps; sin = __FADD( Sin16Ptr[i].Head, prods );
sin2 = Sin16Ptr[i].Tail + prodsl;
sinl = Sin16Ptr[i].Head - sin+prods + sin2;
}
}
else { cos = Cos16Ptr[6].Head;
cosl = Cos16Ptr[6].Tail ;
for(i = 5; i > -1; i--) {
temps = __FMADD( cos, argsq2, cosl*argsq ); prods = __FMADD( cos, argsq, temps ); prodsl = __FMSUB( cos, argsq, prods) + temps; cos = __FADD( Cos16Ptr[i].Head, prods );
cos2 = Cos16Ptr[i].Tail + prodsl;
cosl = Cos16Ptr[i].Head - cos + prods + cos2;
}
}
}
else {
cos = CosPtr[4].Head; sin = SinPtr[4].Head;
cosl = CosPtr[4].Tail;
sinl = SinPtr[4].Tail;
for (i = 3; i > -1; i--) { tempc = __FMADD( cos, argsq2, cosl*argsq ); prodc = __FMADD( cos, argsq, tempc ); prodcl = __FMSUB( cos, argsq, prodc) + tempc; cos = __FADD( CosPtr[i].Head, prodc );
cos2 = CosPtr[i].Tail + prodcl;
cosl = CosPtr[i].Head - cos + prodc + cos2;
temps = __FMADD( sin, argsq2, sinl*argsq ); prods = __FMADD( sin, argsq, temps ); prodsl = __FMSUB( sin, argsq, prods ) + temps; sin = __FADD( SinPtr[i].Head, prods );
sin2 = SinPtr[i].Tail + prodsl;
sinl = SinPtr[i].Head - sin + prods + sin2;
}
}
tempc = __FMADD( cos, argsq2, cosl*argsq ); cost = __FMADD( cos, argsq, tempc ); cosl = __FMSUB( cos, argsq, cost ) + tempc; cos = cost; temps = __FMADD( sin, argsq2, sinl*argsq ); sint = __FMADD( sin, argsq, temps ); sinl = __FMSUB( sin, argsq, sint ) + temps; temps = __FMADD( sint, arg2, sinl * arg ); sin = __FMADD( sint, arg, temps ); sinl = __FMSUB( sint, arg, sin ) + temps;
if (index < 6) {
if ((M & 1) == 0) {
res = __FADD( arg, sin ); reslow = __FADD( arg2, sinl ); resmid = __FADD( arg - res, sin );
resbot = arg2 - reslow + sinl;
restop = __FADD( res, (reslow+resmid) );
resin = __FADD( res - restop, (reslow+resmid) );
resextra = resmid - (reslow+resmid) + reslow;
ressup = __FADD( restop, (resin+resbot) );
restiny = __FSUB( resin, (resin+resbot) + resbot );
resinf = restop - ressup + (resin+resbot) + (resextra+restiny+sinargtiny);
*yHead = ressup;
*yTail = resinf;
}
else {
res = __FADD( 1.0, cos ); reslow = cosl;
resmid = __FADD( 1.0 - res, cos );
restop = __FADD( res, (reslow+resmid) );
resin = __FADD( res - restop, (reslow+resmid) );
resextra = resmid - (reslow+resmid) + reslow;
ressup = __FADD( restop, resin );
resinf = restop - ressup + resin + resextra;
*yHead = ressup;
*yTail = resinf;
}
} else {
if ((M & 1) == 0) { p1 = cos*sintab; p2 = __FMSUB( cos, sintab, p1 ); p3 = __FMADD( cosl, sintab, cos*sintab2 ); p4 = arg* TrigPtr[index].cosHead ; p5 = __FMSUB(arg, TrigPtr[index].cosHead, p4 ); p6 = arg2* TrigPtr[index].cosHead ; p6a = __FMSUB( arg2, TrigPtr[index].cosHead, p6 ); p7 = sin* TrigPtr[index].cosHead ; p8 = __FMSUB( sin, TrigPtr[index].cosHead, p7 ); p9 = __FMADD( sinl, TrigPtr[index].cosHead, (sin* TrigPtr[index].cosTail ) ); p10 = __FMADD( arg, TrigPtr[index].cosTail, p9 + p8); p11 = __FADD( p2, p5 ); p12 = p5 - p11 + p2 + p6a; p13 = __FADD( p6, p11 ); p15 = __FADD( sintab, p4 ); p16 = sintab - p15 + p4; p17 = sintab2 + p10 + p12 + p3; p18 = __FADD( p15, p1 ); p20 = __FADD( p18, p7 ); if (fabs(p6) > fabs(p11))
p14 = p6 - p13 + p11; else
p14 = p11 - p13 + p6;
p19 = p15 - p18 + p1 + p16; p21 = p18 - p20 + p7; p22 = p13 + p17; p23 = __FADD( p19, p21 ); p25 = __FADD( p20, p23 ); p26 = p20 - p25 + p23; if (fabs(p19) > fabs(p21))
p24 = p19 - p23 + p21; else
p24 = p21 - p23 + p19;
p27 = __FADD( (p14+p24), p22 + p26 ); p27a = p26 - p27 + ((p14+p24)+p22);
p28 = __FADD( p25, p27 );
*yHead = p28;
*yTail = p25 - p28 + p27 + p27a;
}
else { p1 = TrigPtr[index].cosHead *cos; p2 = __FMSUB( TrigPtr[index].cosHead, cos, p1 ); p5 = sintab*arg; p6 = __FMSUB( sintab, arg, p5 ); p7 = sintab*arg2; p8 = sintab*sin; p9 = __FMSUB( sintab, sin, p8 ); p12 = __FMSUB( TrigPtr[index].cosTail, cos, (sintab2*arg) ); p13 = __FNMSUB( sintab2, sin, p12 ); p14 = p13 + p2; p15 = __FMADD( TrigPtr[index].cosHead, cosl, p14 ); p16 = __FNMSUB( sintab, sinl, p15 - p9 ); p17 = __FADD( TrigPtr[index].cosHead, p1 ); p18 = TrigPtr[index].cosHead - p17 + p1; p19 = __FSUB( p17, p5 ); p20 = p17 - p19 - p5; p21 = p18 + p20 ; p22 = __FSUB( p19, p8 ); p23 = p19 - p22 - p8; p24 = TrigPtr[index].cosTail + p16; p23a = __FADD( p23, p21 ); p25 = __FSUB( p23a, p6 ); p26 = p23a - p25 - p6; p27 = __FSUB( p25, p7 ); if (fabs(p23) > fabs(p21))
p23b = p23 - p23a + p21; else
p23b = p21 - p23a + p23;
p28 = p25 - p27 - p7 + p26 + p23b; p29 = __FADD( p27, p24 ); p30 = p27 - p29 + p24 + p28; p31 = __FADD( p22, p29 );
*yHead = p31;
*yTail = p22 - p31 + p29 + p30;
}
}
if((M & 2) != 0) { *yHead = - *yHead;
*yTail = - *yTail;
}
}
long double sinl( long double x )
{
uint32_t hibits;
double xHead, xTail, dleft, fenv, yHead, yTail;
DblDbl ddx;
int M;
FEGETENVD(fenv);
FESETENVD(0.0);
ddx.f = x;
xHead = ddx.d[0];
xTail = ddx.d[1];
hibits = ddx.i[0] & 0x7FFFFFFFu;
if (hibits >= 0x7ff00000u) { if (xHead != xHead) { FESETENVD(fenv);
return x;
}
else {
ddx.d[0] = xHead - xHead; ddx.d[1] = ddx.d[0];
FESETENVD(fenv);
return (ddx.f);
}
}
if ((hibits | ddx.i[1]) == 0) {
FESETENVD(fenv);
return x;
}
if (hibits < 0x29900000u) { FESETENVD(fenv);
return x;
}
M = 0; dleft = 0.0;
if (hibits >= 0x42d00000u) { argReduce(xHead, xTail, &dleft, &M, &xHead, &xTail);
}
SinCosCore( M, xHead, xTail, dleft, &yHead, &yTail );
ddx.d[0] = yHead + yTail; ddx.d[1] = (yHead - ddx.d[0]) + yTail;
FESETENVD(fenv);
return (ddx.f);
}
long double cosl( long double x )
{
uint32_t hibits;
double xHead, xTail, dleft, fenv, yHead, yTail;
DblDbl ddx;
int M;
FEGETENVD(fenv);
FESETENVD(0.0);
ddx.f = x;
xHead = ddx.d[0];
xTail = ddx.d[1];
hibits = ddx.i[0] & 0x7FFFFFFFu;
if (hibits >= 0x7ff00000u) { if (xHead != xHead) { FESETENVD(fenv);
return x;
}
else {
ddx.d[0] = xHead - xHead; ddx.d[1] = ddx.d[0];
FESETENVD(fenv);
return (ddx.f);
}
}
if ((hibits | ddx.i[1]) == 0) {
FESETENVD(fenv);
return (long double) 1.0;
}
if (hibits < 0x1e600000u) { FESETENVD(fenv);
return 1.0;
}
M = 1; dleft = 0.0;
if (hibits >= 0x42d00000u) { argReduce(xHead, xTail, &dleft, &M, &xHead, &xTail);
}
SinCosCore( M, xHead, xTail, dleft, &yHead, &yTail );
ddx.d[0] = __FADD( yHead, yTail ); ddx.d[1] = (yHead - ddx.d[0]) + yTail;
FESETENVD(fenv);
return (ddx.f);
}
long double tanl( long double x )
{
double temp, temps;
double b, c, c1, c2, ca, cerr, d;
double frac, frac2, frac3, t1, t2;
double rec, res, reslow, resmid, resbot, restop, resin, resextra;
double ressup, restiny, resinf;
double rarg, rarg0, rarg2, rarg3, rarg2a, redarg, redarg1, redargold;
double p1, p2, p3, p4, p5, p6, p7, p8, p9, p10;
double p11, p12, p13, p14, p15, p16, p17, p18, p19, p20;
double p21, p22, p23, p24, p25, p26, p27, p28, p29, p30;
double p31, p32, p33, p34, p35, p36, p37, p38, p39, p40;
double p41, p42, p43, p44, p45, p46;
double arg, arg2, arg2a, arg2b, argsq, argsq2, argres, argmid, argnew;
double cot, cot2, cotl, cott, tan, tan2, tant, tanl;
double cotargtiny, tanargtiny;
double prods, prodsl, q;
double blow, carry, dleft;
int i, index, M;
uint32_t hibits;
Double DeN;
double fenv;
double xHead, xTail;
DblDbl ddx;
struct TanTblEntry *TanTablePtr = (struct TanTblEntry *)TanTable - 6;
struct TanCoeffEntry *TanPtr = (struct TanCoeffEntry *) TanCoeff - 1;
struct TanCoeff16Entry *Tan16Ptr = (struct TanCoeff16Entry *) TanCoeff16 - 1;
struct CotCoeff16Entry *Cot16Ptr = (struct CotCoeff16Entry *) CotCoeff16 - 1;
FEGETENVD(fenv);
FESETENVD(0.0);
ddx.f = x;
xHead = ddx.d[0];
xTail = ddx.d[1];
hibits = ddx.i[0] & 0x7FFFFFFFu;
if (hibits >= 0x7ff00000u){ if (xHead != xHead) { FESETENVD(fenv);
return x;
}
else {
ddx.d[0] = xHead - xHead; ddx.d[1] = ddx.d[0];
FESETENVD(fenv);
return (ddx.f);
}
}
if ((hibits | ddx.i[1]) == 0) {
FESETENVD(fenv);
return x;
}
if (hibits < 0x29900000u) { FESETENVD(fenv);
return x;
}
M=0; dleft=0.0;
if (hibits >= 0x42b00000u){ argReduce (xHead, xTail, &dleft, &M, &xHead, &xTail);
}
DeN.f = __FMADD( xHead, k2ByPi, k1_8To52 ); q = DeN.f - k1_8To52; redarg = __FNMSUB( q, kPiBy2, xHead );
M = M + DeN.i[1]; b = kPiBy2Tail1*q; blow = __FMSUB( kPiBy2Tail1, q, b ); redarg1 = __FSUB( xTail, b ); DeN.f = __FMADD( fabs(redarg), 64.0, k1_8To52 ); ca = kPiBy2Tail2*q;
if (fabs(b) > fabs(xTail))
carry = xTail - (b+redarg1); else
carry = (xTail-redarg1) - b;
redargold = redarg;
index = DeN.i[1];
if (!(index < 6)) {
if (redargold > 0.0)
redarg = redarg - TanTablePtr[index].arg; else
redarg = redarg + TanTablePtr[index].arg;
}
arg = __FADD( redarg, redarg1 ); arg2a = (redarg - arg + redarg1);
arg2b = __FSUB( arg2a, blow ); c = __FSUB( carry, ca );
argmid = __FADD( arg2b , c );
argnew = __FADD( arg, argmid );
argres = arg - argnew + argmid;
arg = argnew;
cerr = (carry-c) - ca;
if (fabs(arg2b) > fabs(c))
d = arg2b - argmid + c + dleft;
else
d = c - argmid + arg2b + dleft;
d = __FNMSUB( kPiBy2Tail3, q, d + cerr + __FNMSUB( kPiBy2Tail2, q, ca ) );
arg2 = __FADD( argres, d ); temp = 2.0*arg*arg2;
argsq = __FMADD( arg, arg, temp ); argsq2 = __FMADD( arg2, arg2, __FMSUB( arg, arg, argsq ) + temp );
if (index < 6) { if ((M & 1) == 0) { tan = Tan16Ptr[10].Head; tanl = Tan16Ptr[10].Tail; tanargtiny = argres - arg2 + d;
for(i = 9; i > 0; i--) {
temps = __FMADD( tan, argsq2, tanl*argsq ); prods = __FMADD( tan, argsq, temps ); prodsl = __FMSUB( tan, argsq, prods ) + temps; tan = __FADD( Tan16Ptr[i].Head, prods );
tan2 = Tan16Ptr[i].Tail + prodsl;
tanl = Tan16Ptr[i].Head - tan + prods + tan2;
}
temps = __FMADD( tan, argsq2, tanl*argsq ); tant = __FMADD( tan, argsq, temps ); tanl = __FMSUB( tan, argsq, tant ) + temps; temps = __FMADD( tant, arg2, tanl*arg ); tan = __FMADD( tant, arg, temps ); tanl = __FMSUB( tant, arg, tan ) + temps; res = __FADD( arg, tan ); reslow = __FADD( arg2, tanl ); resmid = arg - res + tan;
resbot = arg2 - reslow + tanl;
restop = __FADD( res, (reslow+resmid) );
resin = res - restop + (reslow+resmid);
resextra = resmid - (reslow+resmid) + reslow;
ressup = __FADD( restop, (resin+resbot) );
restiny = resin - (resin+resbot) + resbot;
resinf = restop - ressup + (resin+resbot) + ( resextra+restiny+tanargtiny);
ddx.d[0] = __FADD( ressup, resinf );
ddx.d[1] = (ressup - ddx.d[0]) + resinf;
FESETENVD(fenv);
return (ddx.f);
}
else { cot = Cot16Ptr[9].Head; cotl = Cot16Ptr[9].Tail; cotargtiny = argres - arg2 + d;
rarg0 = 1.0/arg;
for(i = 8; i > 0; i--){ temps = __FMADD( cot, argsq2, cotl*argsq ); prods = __FMADD( cot, argsq, temps ); prodsl = __FMSUB( cot, argsq, prods ) + temps; cot = __FADD( Cot16Ptr[i].Head, prods );
cot2 = Cot16Ptr[i].Tail + prodsl;
cotl = Cot16Ptr[i].Head - cot + prods + cot2;
}
cott = cot;
rarg = rarg0;
temps = __FMADD( cott, arg2, cotl*arg ); cot = __FMADD( cott, arg, temps ); cotl = __FMSUB( cott, arg, cot ) + temps; p1 = __FNMSUB( arg, rarg, 1.0 ); p2 = arg2*rarg; p3 = __FSUB( p1, p2 ); rarg2a = p3*rarg; rarg2 = rarg2a + __FNMSUB( arg, rarg2a, p3 )*rarg;
if (fabs(p1) > fabs(p2))
p4 = p1 - p3 - p2;
else
p4 = p1 - (p3+p2);
p5 = p4 + __FNMSUB( arg2, rarg, p2 ); p6 = __FNMSUB( cotargtiny, rarg, p5 ); p7 = __FNMSUB( rarg2, arg, p3 ); rarg3 = __FNMSUB( arg2, rarg2, p6+p7 )*rarg;
p8 = __FSUB( cot, rarg );
p9 = cot - (rarg+p8);
p10 = __FSUB( cotl, rarg2 );
p11 = cotl - (rarg2+p10);
p12 = p11 - rarg3;
p13 = __FADD( p9, p10 );
if (fabs(p9) > fabs(p10))
p14 = p9 - p13 + p10;
else
p14 = p10 - p13 + p9;
p15 = p12 + p14;
p16 = __FADD( p8, p13 ); p17 = p8 - p16 + p13; p18 = p17 + p15;
ddx.d[0] = __FADD( p16, p18 );
ddx.d[1] = (p16 - ddx.d[0]) + p18;
FESETENVD(fenv);
return (ddx.f);
}
}
else { tant = TanPtr[7].Head; tant = __FMADD( tant, argsq, TanPtr[6].Head ); tan = __FMADD( tant, argsq, TanPtr[5].Head ); tanl = __FMADD( tant, argsq, TanPtr[5].Head - tan ); tanargtiny = argres - arg2 + d;
for(i = 4; i > 0; i--) {
temps = __FMADD( tan, argsq2, tanl*argsq ); prods = __FMADD( tan, argsq, temps ); prodsl = __FMSUB( tan, argsq, prods ) + temps; tan = __FADD( TanPtr[i].Head, prods );
tan2 = TanPtr[i].Tail + prodsl;
tanl = TanPtr[i].Head -tan + prods + tan2;
}
tant = tan;
temps = __FMADD( tant, arg2, tanl*arg ); tan = __FMADD( tant, arg, temps ); tanl = __FMSUB( tant, arg, tan ) + temps; if (redargold > 0) { t1 = TanTablePtr[index].TanHead; t2 = TanTablePtr[index].TanTail;
c1 = TanTablePtr[index].CotHead; c2 = TanTablePtr[index].CotTail;
}
else {
t1 = -TanTablePtr[index].TanHead;
t2 = -TanTablePtr[index].TanTail;
c1 = -TanTablePtr[index].CotHead;
c2 = -TanTablePtr[index].CotTail;
}
p1 = __FADD( t1, c1 ); p2 = __FADD( t2, c2 ); p3 = c1 - p1 + t1;
p4 = c2 - p2 + t2;
p5 = __FADD( p3, p2 );
p6 = p3 - p5 + p2 + p4; p7 = p1*tan;
p8 = __FMSUB( p1, tan, p7 ); p9 = p5*tan;
p10 = __FMSUB( p5, tan, p9 ); p11 = p1*tanl;
p12 = __FMSUB( p1, tanl, p11 ); p13 = __FMADD( p5, tanl, p10 + p12 ); p14 = __FMADD( p6, tan, p13 ); p15 = __FADD( p8, p9 );
p17 = __FADD( p15, p11 );
if (fabs(p8) > fabs(p9))
p16 = p8 - p15 + p9;
else
p16 = p9 - p15 + p8;
if (fabs(p11) > fabs(p15))
p18 = p11 - p17 + p15 + (p14+p16);
else
p18 = p15 - p17 + p11 + (p14+p16);
p19 = __FADD( p7, p17 );
p20 = p7 - p19 + p17;
p21 = __FADD( p20, p18 );
p22 = p20 - p21 + p18;
if ((M & 1) == 0) {
p23 = __FSUB( c1, tan ); p24 = __FSUB( c2, tanl );
rec = 1.0/p23;
p25 = c1 - p23 - tan;
frac = p19*rec;
if (fabs(c2) > fabs(tanl))
p26 = c2 - p24 - tanl;
else
p26 = c2 - (p24+tanl);
p27 = __FADD( p25, p24 );
p28 = p25 - p27 + p24 + p26; frac = frac + __FNMSUB( p23, frac, p19 )*rec; p29 = __FNMSUB( p23, frac, p19 ); p30 = frac*p27; p31 = __FNMSUB( frac, p27, p30 ) + p22; p32 = __FNMSUB( frac, p28, p31 ); p33 = __FSUB( p29, p30 );
p35 = p33 + p21;
frac2 = p35*rec;
frac2 = frac2 + __FNMSUB( frac2, p23, p35)*rec;
if (fabs(p29) > fabs(p30))
p34 = p29 - p33 - p30 + p32;
else
p34 = p29 - (p33+p30) + p32;
p36 = __FNMSUB( frac2, p23, p35 ) + __FNMSUB( frac2, p27, p34 );
if (fabs(p33) > fabs(p21))
p37 = p33 - p35 + p21;
else
p37 = p21 - p35 + p33;
p38 = p37 + p36;
frac3 = p38*rec;
p39 = __FADD( t1, frac );
p40 = t1 - p39 + frac;
p41 = __FADD( t2, frac2 );
if (fabs(t2) > fabs(frac2))
p42 = t2 - p41 + frac2 + frac3;
else
p42 = frac2 - p41 + t2 + frac3;
p43 = __FADD( p40, p41 );
p44 = p40 - p43 + p41 + p42;
p45 = __FADD( p39, p43 );
p46 = p39 - p45 + p43 + p44;
ddx.d[0] = __FADD( p45, p46 );
ddx.d[1] = (p45 - ddx.d[0]) + p46;
FESETENVD(fenv);
return (ddx.f);
}
else {
p23 = __FADD( t1, tan ); p24 = __FADD( t2, tanl );
rec = 1.0/p23;
p25 = t1 - p23 + tan;
frac = p19*rec;
if (fabs(t2) > fabs(tanl))
p26 = t2 - p24 + tanl;
else
p26 = tanl - p24 + t2;
p27 = __FADD( p25, p24 );
p28 = p25 - p27 + p24 + p26; frac = frac + __FNMSUB( p23, frac, p19 )*rec; p29 = __FNMSUB( p23, frac, p19 ); p30 = frac*p27; p31 = __FNMSUB( frac, p27, p30 ) + p22; p32 = __FNMSUB( frac, p28, p31 ); p33 = __FSUB( p29, p30 );
p35 = __FADD( p33, p21 );
frac2 = p35*rec;
frac2 = frac2, __FNMSUB( frac2, p23, p35 )*rec;
if (fabs(p29) > fabs(p30))
p34 = p29 - p33 - p30 + p32;
else
p34 = p29 - (p33+p30) + p32;
p36 = __FNMSUB( frac2, p23, p35 ) + __FNMSUB( frac2, p27, p34 );
if (fabs(p33) > fabs(p21))
p37 = p33 - p35 + p21;
else
p37 = p21 - p35 + p33;
p38 = p37 + p36;
frac3 = p38*rec;
p39 = __FSUB( frac, c1 );
p40 = frac - (c1+p39);
p41 = __FSUB( frac2, c2 );
if (fabs(c2) > fabs(frac2))
p42 = frac2 - (c2+p41) + frac3;
else
p42 = frac2 - p41 - c2 + frac3;
p43 = __FADD( p40, p41 );
p44 = p40 - p43 + p41 + p42;
p45 = __FADD( p39, p43 );
p46 = p39 - p45 + p43 + p44;
ddx.d[0] = __FADD( p45, p46 );
ddx.d[1] = (p45 - ddx.d[0]) + p46;
FESETENVD(fenv);
return (ddx.f);
}
}
}
#endif