math.inc 9.8 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343
  1. {
  2. $Id$
  3. This file is part of the Free Pascal run time library.
  4. Copyright (c) 2000 by Jonas Maebe and other members of the
  5. Free Pascal development team
  6. Implementation of mathamatical Routines (only for real)
  7. See the file COPYING.FPC, included in this distribution,
  8. for details about the copyright.
  9. This program is distributed in the hope that it will be useful,
  10. but WITHOUT ANY WARRANTY; without even the implied warranty of
  11. MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.
  12. **********************************************************************}
  13. {****************************************************************************
  14. Int to real helpers
  15. ****************************************************************************}
  16. const
  17. longint_to_real_helper: int64 = $4330000080000000;
  18. cardinal_to_real_helper: int64 = $430000000000000;
  19. int_to_real_factor: double = double(high(cardinal))+1.0;
  20. function fpc_int64_to_double(i: int64): double; compilerproc;
  21. {assembler;}
  22. { input: high(i) in r3, low(i) in r4 }
  23. { output: double(i) in f0 }
  24. var
  25. temp: packed record
  26. case byte of
  27. 0: (l1,l2: cardinal);
  28. 1: (d: double);
  29. end;
  30. begin{asm}
  31. (* lis r0,0x4330
  32. stw r0,temp
  33. xoris r3,r3,0x8000
  34. stw r3,4+temp
  35. {$ifndef macos}
  36. lis r3,longint_to_real_helper@ha
  37. lfd f1,longint_to_real_helper@l(r3)
  38. {$else}
  39. lfd f1,longint_to_real_helper(r2)
  40. {$endif}
  41. lfd f0,temp
  42. stw r4,4+temp
  43. fsub f0,f0,f1
  44. {$ifndef macos}
  45. lis r4,cardinal_to_real_helper@ha
  46. lfd f1,cardinal_to_real_helper@l(r4)
  47. lis r3,int_to_real_factor@ha
  48. lfd f3,temp
  49. lfd f2,int_to_real_factor@l(r3)
  50. {$else}
  51. lfd f1,cardinal_to_real_helper(r2)
  52. lfd f3,temp
  53. lfd f2,int_to_real_factor(r2)
  54. {$endif}
  55. fsub f3,f3,f1
  56. fmadd f1,f0,f2,f3*)
  57. end{ ['R0','R3','R4','F0','F1','F2','F3']};
  58. function fpc_qword_to_double(q: qword): double; compilerproc;
  59. {assembler;}
  60. { input: high(q) in r3, low(q) in r4 }
  61. { output: double(q) in f0 }
  62. var
  63. temp: packed record
  64. case byte of
  65. 0: (l1,l2: cardinal);
  66. 1: (d: double);
  67. end;
  68. begin{asm}
  69. (* lis r0,0x4330
  70. stw r0,temp
  71. stw r3,4+temp
  72. lfd f0,temp
  73. {$ifndef macos}
  74. lis r3,cardinal_to_real_helper@ha
  75. lfd f1,cardinal_to_real_helper@l(r3)
  76. {$else}
  77. lfd f1,cardinal_to_real_helper(r2)
  78. {$endif}
  79. stw r4,4+temp
  80. fsub f0,f0,f1
  81. lfd f3,temp
  82. {$ifndef macos}
  83. lis r3,int_to_real_factor@ha
  84. lfd f2,int_to_real_factor@l(r3)
  85. {$else}
  86. lfd f2,int_to_real_factor(r2)
  87. {$endif}
  88. fsub f3,f3,f1
  89. fmadd f1,f0,f2,f3*)
  90. end{ ['R0','R3','F0','F1','F2','F3']};
  91. {****************************************************************************
  92. EXTENDED data type routines
  93. ****************************************************************************}
  94. {$define FPC_SYSTEM_HAS_PI}
  95. function pi : double;[internproc:in_pi];
  96. {$define FPC_SYSTEM_HAS_ABS}
  97. function abs(d : extended) : extended;[internproc:in_abs_extended];
  98. {$define FPC_SYSTEM_HAS_SQR}
  99. function sqr(d : extended) : extended;[internproc:in_sqr_extended];
  100. {$define FPC_SYSTEM_HAS_SQRT}
  101. function sqrt(d : extended) : extended;[internproc:in_sqrt_extended];
  102. {
  103. function arctan(d : extended) : extended;[internconst:in_arctan_extended];
  104. begin
  105. runerror(207);
  106. end;
  107. function ln(d : extended) : extended;[internconst:in_ln_extended];
  108. begin
  109. runerror(207);
  110. end;
  111. function sin(d : extended) : extended;[internconst: in_sin_extended];
  112. begin
  113. runerror(207);
  114. end;
  115. function cos(d : extended) : extended;[internconst:in_cos_extended];
  116. begin
  117. runerror(207);
  118. end;
  119. function exp(d : extended) : extended;[internconst:in_const_exp];
  120. begin
  121. runerror(207);
  122. end;
  123. function frac(d : extended) : extended;[internconst:in_const_frac];
  124. begin
  125. runerror(207);
  126. end;
  127. }
  128. {$define FPC_SYSTEM_HAS_INT}
  129. {$warning FIX ME}
  130. function int(d : extended) : extended;[internconst:in_const_int];
  131. begin
  132. runerror(207);
  133. end;
  134. {$define FPC_SYSTEM_HAS_TRUNC}
  135. {$warning FIX ME}
  136. function trunc(d : extended) : int64;{assembler;}[internconst:in_const_trunc];
  137. { input: d in fr1 }
  138. { output: result in r3 }
  139. {assembler;}
  140. var
  141. temp: packed record
  142. case byte of
  143. 0: (l1,l2: longint);
  144. 1: (d: double);
  145. end;
  146. begin{asm}
  147. { fctiwz f1,f1
  148. stfd f1,temp
  149. lwz r3,temp
  150. lwz r4,4+temp}
  151. end{ ['R3','F1']};
  152. {$define FPC_SYSTEM_HAS_ROUND}
  153. {$ifdef hascompilerproc}
  154. function round(d : extended) : int64;[internconst:in_const_round, external name 'FPC_ROUND'];
  155. function fpc_round(d : extended) : int64;[public, alias:'FPC_ROUND'];{$ifdef hascompilerproc}compilerproc;{$endif hascompilerproc}
  156. {$else}
  157. function round(d : extended) : int64;[internconst:in_const_round];
  158. {$endif hascompilerproc}
  159. { input: d in fr1 }
  160. { output: result in r3 }
  161. {assembler;}
  162. var
  163. temp : packed record
  164. case byte of
  165. 0: (l1,l2: longint);
  166. 1: (d: double);
  167. end;
  168. begin{asm}
  169. { fctiw f1,f1
  170. stfd f1,temp
  171. lwz r3,temp
  172. lwz r4,4+temp}
  173. end{ ['R3','F1']};
  174. {$define FPC_SYSTEM_HAS_POWER}
  175. function power(bas,expo : extended) : extended;
  176. begin
  177. if bas=0 then
  178. begin
  179. if expo<>0 then
  180. power:=0.0
  181. else
  182. HandleError(207);
  183. end
  184. else if expo=0 then
  185. power:=1
  186. else
  187. { bas < 0 is not allowed }
  188. if bas<0 then
  189. handleerror(207)
  190. else
  191. power:=exp(ln(bas)*expo);
  192. end;
  193. {****************************************************************************
  194. Longint data type routines
  195. ****************************************************************************}
  196. {$define FPC_SYSTEM_HAS_POWER_INT64}
  197. function power(bas,expo : Int64) : Int64;
  198. begin
  199. if bas=0 then
  200. begin
  201. if expo<>0 then
  202. power:=0
  203. else
  204. HandleError(207);
  205. end
  206. else if expo=0 then
  207. power:=1
  208. else
  209. begin
  210. if bas<0 then
  211. begin
  212. if odd(expo) then
  213. power:=-round(exp(ln(-bas)*expo))
  214. else
  215. power:=round(exp(ln(-bas)*expo));
  216. end
  217. else
  218. power:=round(exp(ln(bas)*expo));
  219. end;
  220. end;
  221. {****************************************************************************
  222. Helper routines to support old TP styled reals
  223. ****************************************************************************}
  224. { warning: the following converts a little-endian TP-style real }
  225. { to a big-endian double. So don't byte-swap the TP real! }
  226. {$define FPC_SYSTEM_HAS_REAL2DOUBLE}
  227. function real2double(r : real48) : double;
  228. var
  229. res : array[0..7] of byte;
  230. exponent : word;
  231. begin
  232. { copy mantissa }
  233. res[6]:=0;
  234. res[5]:=r[1] shl 5;
  235. res[4]:=(r[1] shr 3) or (r[2] shl 5);
  236. res[3]:=(r[2] shr 3) or (r[3] shl 5);
  237. res[2]:=(r[3] shr 3) or (r[4] shl 5);
  238. res[1]:=(r[4] shr 3) or (r[5] and $7f) shl 5;
  239. res[0]:=(r[5] and $7f) shr 3;
  240. { copy exponent }
  241. { correct exponent: }
  242. exponent:=(word(r[0])+(1023-129));
  243. res[1]:=res[1] or ((exponent and $f) shl 4);
  244. res[0]:=exponent shr 4;
  245. { set sign }
  246. res[0]:=res[0] or (r[5] and $80);
  247. real2double:=double(res);
  248. end;
  249. {
  250. $Log$
  251. Revision 1.4 2003-04-23 21:28:21 peter
  252. * fpc_round added, needed for int64 currency
  253. Revision 1.3 2003/01/22 20:45:15 mazen
  254. * making math code in RTL compiling.
  255. *NB : This does NOT mean necessary that it will generate correct code!
  256. Revision 1.2 2003/01/20 22:21:36 mazen
  257. * many stuff related to RTL fixed
  258. Revision 1.1 2002/12/24 21:30:20 mazen
  259. - some writeln(s) removed in compiler
  260. + many files added to RTL
  261. * some errors fixed in RTL
  262. Revision 1.14 2002/11/28 11:04:16 olle
  263. * macos: refs to globals in begin{asm} adapted to macos
  264. Revision 1.13 2002/10/21 18:08:28 jonas
  265. * round has int64 instead of longint result
  266. Revision 1.12 2002/09/08 13:00:21 jonas
  267. * made pi an internproc instead of internconst
  268. Revision 1.11 2002/09/07 16:01:26 peter
  269. * old logs removed and tabs fixed
  270. Revision 1.10 2002/08/18 22:11:10 florian
  271. * fixed remaining assembler errors
  272. Revision 1.9 2002/08/18 21:37:48 florian
  273. * several errors in inline assembler fixed
  274. Revision 1.8 2002/08/10 17:14:36 jonas
  275. * various fixes, mostly changing the names of the modifies registers to
  276. upper case since that seems to be required by the compiler
  277. Revision 1.7 2002/07/31 16:58:12 jonas
  278. * fixed conversion from int64/qword to double errors
  279. Revision 1.6 2002/07/29 21:28:17 florian
  280. * several fixes to get further with linux/ppc system unit compilation
  281. Revision 1.5 2002/07/28 21:39:29 florian
  282. * made abs a compiler proc if it is generic
  283. Revision 1.4 2002/07/28 20:43:49 florian
  284. * several fixes for linux/powerpc
  285. * several fixes to MT
  286. }