summaryrefslogtreecommitdiff
path: root/rtl/m68k
diff options
context:
space:
mode:
authorfpc <fpc@3ad0048d-3df7-0310-abae-a5850022a9f2>2005-05-16 18:37:41 +0000
committerfpc <fpc@3ad0048d-3df7-0310-abae-a5850022a9f2>2005-05-16 18:37:41 +0000
commitf206a9c2b1ae1d8727ca27a96d448b61fdb4c766 (patch)
treef28256ff9964c1fc7c0f7fb00891268a117b745d /rtl/m68k
downloadfpc-f206a9c2b1ae1d8727ca27a96d448b61fdb4c766.tar.gz
initial import
git-svn-id: http://svn.freepascal.org/svn/fpc/trunk@1 3ad0048d-3df7-0310-abae-a5850022a9f2
Diffstat (limited to 'rtl/m68k')
-rw-r--r--rtl/m68k/lowmath.inc920
-rw-r--r--rtl/m68k/m68k.inc329
-rw-r--r--rtl/m68k/makefile.cpu7
-rw-r--r--rtl/m68k/math.inc949
-rw-r--r--rtl/m68k/readme7
-rw-r--r--rtl/m68k/set.inc428
-rw-r--r--rtl/m68k/setjump.inc34
-rw-r--r--rtl/m68k/setjumph.inc37
8 files changed, 2711 insertions, 0 deletions
diff --git a/rtl/m68k/lowmath.inc b/rtl/m68k/lowmath.inc
new file mode 100644
index 0000000000..f5a0c3fd92
--- /dev/null
+++ b/rtl/m68k/lowmath.inc
@@ -0,0 +1,920 @@
+{
+ $Id: lowmath.inc,v 1.5 2005/02/14 17:13:30 peter Exp $
+ This file is part of the Free Pascal run time library.
+ Copyright (c) 1999-2000 by Carl-Eric Codere,
+ member of the Free Pascal development team
+
+ See the file COPYING.FPC, included in this distribution,
+ for details about the copyright.
+
+ This program 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.
+
+ **********************************************************************}
+{*************************************************************************}
+{ lowmath.inc }
+{ Ported to FPC-Pascal by Carl Eric Codere }
+{ Terms of use: This source code is freeware. }
+{*************************************************************************}
+{ This inc. implements low-level mathemtical routines for the motorola }
+{ 68000 family of processors. }
+{*************************************************************************}
+{ single floating point routines taken from GCC 2.5.2 Atari compiler }
+{ library source. }
+{ Original credits: }
+{ written by Kai-Uwe Bloem (I5110401@dbstu1.bitnet). }
+{ Based on a 80x86 floating point packet from comp.os.minix, }
+{ written by P.Housel }
+{ Patched by Olaf Flebbe (flebbe@tat.physik.uni-tuebingen.de) }
+{ Revision by michal 05-93 (ntomczak@vm.ucs.ualberta.ca) }
+{*************************************************************************}
+{--------------------------------------------------------------------}
+{ LEFT TO DO: }
+{ o Add support for FPU if present. }
+{ o Verify if single comparison works in all cases. }
+{ o Add support for NANs in SINGLE_CMP }
+{ o Add comp (80-bit) multiplication,addition,substract,division, }
+{ shift. }
+{ o Add stack checking for the routines which use the stack. }
+{ (This will probably have to be done in the code generator). }
+{--------------------------------------------------------------------}
+
+
+
+Procedure Single_Norm;[alias : 'FPC_SINGLE_NORM'];Assembler;
+{--------------------------------------------}
+{ Low-level routine to normalize single e }
+{ IEEE floating point values. Never called }
+{ directly. }
+{ On Exit: }
+{ d0 = result. }
+{ Registers destroyed: d0,d1 }
+{--------------------------------------------}
+Asm
+ tst.l d4 { rounding and u.mant == 0 ? }
+ bne @normlab1
+ tst.b d1
+ beq @retzok
+@normlab1:
+ clr.b d2 { "sticky byte" }
+@normlab3:
+ move.l #$ff000000,d5
+@normlab4:
+ tst.w d0 { divide (shift) }
+ ble @normlab2 { denormalized number }
+ move.l d4,d3
+ and.l d5,d3 { or until no bits above 23 }
+ beq @normlab5
+@normlab2:
+ addq.w #1,d0 { increment exponent }
+ lsr.l #1,d4
+ or.b d1,d2 { set "sticky" }
+ roxr.b #1,d1 { shift into rounding bits }
+ bra @normlab4
+@normlab5:
+ and.b #1,d2
+ or.b d2,d1 { make least sig bit "sticky" }
+ asr.l #1,d5 { #0xff800000 -> d5 }
+@normlab6:
+ move.l d4,d3 { multiply (shift) until }
+ and.l d5,d3 { one in "implied" position }
+ bne @normlab7
+ subq.w #1,d0 { decrement exponent }
+ beq @normlab7 { too small. store as denormalized number }
+ add.b d1,d1 { some doubt about this one * }
+ addx.l d4,d4
+ bra @normlab6
+@normlab7:
+ tst.b d1 { check rounding bits }
+ bge @normlab9 { round down - no action neccessary }
+ neg.b d1
+ bvc @normlab8 { round up }
+ move.w d4,d1 { tie case - round to even }
+ { dont need rounding bits any more }
+ and.w #1,d1 { check if even }
+ beq @normlab9 { mantissa is even - no action necessary }
+ { fall through }
+@normlab8:
+ clr.w d1 { zero rounding bits }
+ add.l #1,d4
+ tst.w d0
+ bne @normlab10 { renormalize if number was denormalized }
+ add.w #1,d0 { correct exponent for denormalized numbers }
+ bra @normlab3
+@normlab10:
+ move.l d4,d3 { check for rounding overflow }
+ asl.l #1,d5 { #0xff000000 -> d5 }
+ and.l d5,d3
+ bne @normlab4 { go back and renormalize }
+@normlab9:
+ tst.l d4 { check if normalization caused an underflow }
+ beq @retz
+ tst.w d0 { check for exponent overflow or underflow }
+ blt @retz
+ cmp.w #255,d0
+ bge @oflow
+
+ lsl.w #8,d0 { re-position exponent - one bit too high }
+ lsl.w #1,d2 { get X bit }
+ roxr.w #1,d0 { shift it into sign position }
+ swap d0 { map to upper word }
+ clr.w d0
+ and.l #$7fffff,d4 { top mantissa bits }
+ or.l d4,d0 { insert exponent and sign }
+ movem.l (sp)+,d2-d5
+ rts
+
+@retz:
+ { handling underflow should be done here... }
+ { by default simply return 0 as retzok... }
+@retzok:
+ moveq.l #0,d0
+ lsl.w #1,d2
+ roxr.l #1,d0 { sign of 0 is the same as of d2 }
+ movem.l (sp)+,d2-d5
+ rts
+
+@oflow:
+ move.l #$7f800000,d0 { +infinity as proposed by IEEE }
+
+ tst.w d2 { transfer sign }
+ bge @ofl_clear { (mjr++) }
+ bset #31,d0 { }
+@ofl_clear:
+ or.b #2,ccr { set overflow flag. }
+ movem.l (sp)+,d2-d5
+ rts
+end;
+
+
+Procedure Single_AddSub; Assembler;
+{--------------------------------------------}
+{ Low-level routine to add/subtract single }
+{ IEEE floating point values. Never called }
+{ directly. }
+{ On Exit: }
+{ d0 = result -- from normalize routine }
+{ Flags : V set if overflow. }
+{ on underflow d0 = 0 }
+{ Registers destroyed: d0,d1 }
+{--------------------------------------------}
+Asm
+{--------------------------------------------}
+{ On Entry: }
+{ d1-d0 = single values to subtract. }
+{--------------------------------------------}
+XDEF SINGLE_SUB
+ eor.l #$80000000,d0 { reverse sign of v }
+{--------------------------------------------}
+{ On Entry: }
+{ d0, d1 = single values to add. }
+{--------------------------------------------}
+XDEF SINGLE_ADD
+ movem.l d2-d5,-(sp) { save registers }
+ move.l d0,d4 { d4 = d0 = v }
+ move.l d1,d5 { d5 = d1 = u }
+
+ move.l #$7fffff,d3
+ move.l d5,d0 { d0 = u.exp }
+ move.l d5,d2 { d2.h = u.sign }
+ swap d0
+ move.w d0,d2 { d2 = u.sign }
+ and.l d3,d5 { remove exponent from u.mantissa }
+
+ move.l d4,d1 { d1 = v.exp }
+ and.l d3,d4 { remove exponent from v.mantissa }
+ swap d1
+ eor.w d1,d2 { d2 = u.sign ^ v.sign (in bit 15)}
+ clr.b d2 { we will use the lowest byte as a flag }
+ moveq.l #15,d3
+ bclr d3,d1 { kill sign bit u.exp }
+ bclr d3,d0 { kill sign bit u.exp }
+ btst d3,d2 { same sign for u and v? }
+ beq @slabel1
+ cmp.l d0,d1 { different signs - maybe x - x ? }
+ seq d2 { set 'cancellation' flag }
+@slabel1:
+ lsr.w #7,d0 { keep here exponents only }
+ lsr.w #7,d1
+{--------------------------------------------------------------------}
+{ Now perform testing of NaN and infinities }
+{--------------------------------------------------------------------}
+ moveq.l #-1,d3
+ cmp.b d3,d0
+ beq @alabel1
+ cmp.b d3,d1
+ bne @nospec
+ bra @alabel2
+{--------------------------------------------------------------------}
+{ u is special. }
+{--------------------------------------------------------------------}
+@alabel1:
+ tst.b d2
+ bne @retnan { cancellation of specials -> NaN }
+ tst.l d5
+ bne @retnan { arith with Nan gives always NaN }
+
+ addq.w #4,a0 { here is an infinity }
+ cmp.b d3,d1
+ bne @alabel3 { skip check for NaN if v not special }
+{--------------------------------------------------------------------}
+{ v is special. }
+{--------------------------------------------------------------------}
+@alabel2:
+ tst.l d4
+ bne @retnan
+@alabel3:
+ move.l (a0),d0
+ bra @return
+{--------------------------------------------------------------------}
+{ Return a quiet nan }
+{--------------------------------------------------------------------}
+@retnan:
+ moveq.l #-1,d0
+ lsr.l #1,d0 { 0x7fffffff -> d0 }
+ bra @return
+{ Ok, no inifinty or NaN involved.. }
+@nospec:
+ tst.b d2
+ beq @alabel4
+ moveq.l #0,d0 { x - x hence we always return +0 }
+@return:
+ movem.l (sp)+,d2-d5
+ rts
+
+@alabel4:
+ moveq.l #23,d3
+ bset d3,d5 { restore implied leading "1" }
+ tst.w d0 { check for zero exponent - no leading "1" }
+ bne @alabel5
+ bclr d3,d5 { remove it }
+ addq.w #1,d0 { "normalize" exponent }
+@alabel5:
+ bset d3,d4 { restore implied leading "1" }
+ tst.w d1 { check for zero exponent - no leading "1" }
+ bne @alabel6
+ bclr d3,d4 { remove it }
+ addq.w #1,d1 { "normalize" exponent }
+@alabel6:
+ moveq.l #0,d3 { (put initial zero rounding bits in d3) }
+ neg.w d1 { d1 = u.exp - v.exp }
+ add.w d0,d1
+ beq @alabel8 { exponents are equal - no shifting neccessary }
+ bgt @alabel7 { not equal but no exchange neccessary }
+ exg d4,d5 { exchange u and v }
+ sub.w d1,d0 { d0 = u.exp - (u.exp - v.exp) = v.exp }
+ neg.w d1
+ tst.w d2 { d2.h = u.sign ^ (u.sign ^ v.sign) = v.sign }
+ bpl @alabel7
+ bchg #31,d2
+@alabel7:
+ cmp.w #26,d1 { is u so much bigger that v is not }
+ bge @alabel9 { significant ? }
+{--------------------------------------------------------------------}
+{ shift mantissa left two digits, to allow cancellation of }
+{ most significant digit, while gaining an additional digit for }
+{ rounding. }
+{--------------------------------------------------------------------}
+ moveq.l #1,d3
+@alabel10:
+ add.l d5,d5
+ subq.w #1,d0 { decrement exponent }
+ subq.w #1,d1 { done shifting altogether ? }
+ dbeq d3,@alabel10 { loop if still can shift u.mant more }
+ moveq.l #0,d3
+
+ cmp.w #16,d1 { see if fast rotate possible }
+ blt @alabel11
+ or.w d4,d3 { set rounding bits }
+ clr.w d4
+ swap d4
+ subq.w #8,d1
+ subq.w #8,d1
+ bra @alabel11
+
+@alabel12:
+ move.b d4,d2
+ and.b #1,d2
+ or.b d2,d3
+ lsr.l #1,d4 { shift v.mant right the rest of the way }
+@alabel11:
+ dbra d1,@alabel12 { loop }
+
+@alabel8:
+ tst.w d2 { are the signs equal ? }
+ bpl @alabel13 { yes, no negate necessary }
+
+
+ tst.w d3 { negate rounding bits and v.mant }
+ beq @alabel14
+ addq.l #1,d4
+@alabel14:
+ neg.l d4
+
+@alabel13:
+ add.l d4,d5 { u.mant = u.mant + v.mant }
+ bcs @alabel9 { needn not negate }
+ tst.w d2 { opposite signs ? }
+ bpl @alabel9 { do not need to negate result }
+
+ neg.l d5
+ not.l d2 { switch sign }
+@alabel9:
+ move.l d5,d4 { move result for normalization }
+ clr.l d1
+ tst.l d3
+ beq @alabel15
+ moveq.l #-1,d1
+@alabel15:
+ swap d2 { put sign into d2 (exponent is in d0) }
+ jmp FPC_SINGLE_NORM { leave registers on stack for norm_sf }
+end;
+
+
+Procedure Single_Mul;Assembler;
+{--------------------------------------------}
+{ Low-level routine to multiply two single }
+{ IEEE floating point values. Never called }
+{ directly. }
+{ Om Entry: }
+{ d0,d1 = values to multiply }
+{ On Exit: }
+{ d0 = result. }
+{ Registers destroyed: d0,d1 }
+{ stack space used (and restored): 8 bytes. }
+{--------------------------------------------}
+Asm
+XDEF SINGLE_MUL
+ movem.l d2-d5,-(sp)
+ move.l d0,d4 { d4 = v }
+ move.l d1,d5 { d5 = u }
+
+ move.l #$7fffff,d3
+ move.l d5,d0 { d0 = u.exp }
+ and.l d3,d5 { remove exponent from u.mantissa }
+ swap d0
+ move.w d0,d2 { d2 = u.sign }
+
+ move.l d4,d1 { d1 = v.exp }
+ and.l d3,d4 { remove exponent from v.mantissa }
+ swap d1
+ eor.w d1,d2 { d2 = u.sign ^ v.sign (in bit 15)}
+
+ moveq.l #15,d3
+ bclr d3,d0 { kill sign bit }
+ bclr d3,d1 { kill sign bit }
+ tst.l d0 { test if one of factors is 0 }
+ beq @mlabel1
+ tst.l d1
+@mlabel1:
+ seq d2 { 'one of factors is 0' flag in the lowest byte }
+ lsr.w #7,d0 { keep here exponents only }
+ lsr.w #7,d1
+
+{--------------------------------------------------------------------}
+{ Now perform testing of NaN and infinities }
+{--------------------------------------------------------------------}
+ moveq.l #-1,d3
+ cmp.b d3,d0
+ beq @mlabel2
+ cmp.b d3,d1
+ bne @mnospec
+ bra @mlabel3
+{--------------------------------------------------------------------}
+{ first operand is special }
+{--------------------------------------------------------------------}
+@mlabel2:
+ tst.l d5 { is it NaN? }
+ bne @mretnan
+@mlabel3:
+ tst.b d2 { 0 times special or special times 0 ? }
+ bne @mretnan { yes -> NaN }
+ cmp.b d3,d1 { is the other special ? }
+ beq @mlabel4 { maybe it is NaN }
+{--------------------------------------------------------------------}
+{ Return infiny with correct sign }
+{--------------------------------------------------------------------}
+@mretinf:
+ move.l #$ff000000,d0 { we will return #0xff800000 or #0x7f800000 }
+ lsl.w #1,d2
+ roxr.l #1,d0 { shift in high bit as given by d2 }
+@mreturn:
+ movem.l (sp)+,d2-d5
+ rts
+
+{--------------------------------------------------------------------}
+{ v is special. }
+{--------------------------------------------------------------------}
+@mlabel4:
+ tst.l d4 { is this NaN? }
+ beq @mretinf { we know that the other is not zero }
+@mretnan:
+ moveq.l #-1,d0
+ lsr.l #1,d0 { 0x7fffffff -> d0 }
+ bra @mreturn
+{--------------------------------------------------------------------}
+{ End of NaN and Inf }
+{--------------------------------------------------------------------}
+@mnospec:
+ tst.b d2 { not needed - but we can waste two instr. }
+ bne @mretzz { return signed 0 if one of factors is 0 }
+ moveq.l #23,d3
+ bset d3,d5 { restore implied leading "1" }
+ subq.w #8,sp { multiplication accumulator }
+ tst.w d0 { check for zero exponent - no leading "1" }
+ bne @mlabel5
+ bclr d3,d5 { remove it }
+ addq.w #1,d0 { "normalize" exponent }
+@mlabel5:
+ tst.l d5
+ beq @mretz { multiplying zero }
+
+ moveq.l #23,d3
+ bset d3,d4 { restore implied leading "1" }
+ tst.w d1 { check for zero exponent - no leading "1" }
+ bne @mlabel6
+ bclr d3,d4 { remove it }
+ addq.w #1,d1 { "normalize" exponent }
+@mlabel6:
+ tst.l d4
+ beq @mretz { multiply by zero }
+
+ add.w d1,d0 { add exponents, }
+ sub.w #BIAS4+16-8,d0 { remove excess bias, acnt for repositioning }
+
+ clr.l (sp) { initialize 64-bit product to zero }
+ clr.l 4(sp)
+{--------------------------------------------------------------------}
+{ see Knuth, Seminumerical Algorithms, section 4.3. algorithm M }
+{--------------------------------------------------------------------}
+ move.w d4,d3
+ mulu.w d5,d3 { mulitply with bigit from multiplier }
+ move.l d3,4(sp) { store into result }
+
+ move.l d4,d3
+ swap d3
+ mulu.w d5,d3
+ add.l d3,2(sp) { add to result }
+
+ swap d5 { [TOP 8 BITS SHOULD BE ZERO !] }
+
+ move.w d4,d3
+ mulu.w d5,d3 { mulitply with bigit from multiplier }
+ add.l d3,2(sp) { store into result (no carry can occur here) }
+
+ move.l d4,d3
+ swap d3
+ mulu.w d5,d3
+ add.l d3,(sp) { add to result }
+ { [TOP 16 BITS SHOULD BE ZERO !] }
+ movem.l 2(sp),d4-d5 { get the 48 valid mantissa bits }
+ clr.w d5 { (pad to 64) }
+
+ move.l #$0000ffff,d3
+@mlabel7:
+ cmp.l d3,d4 { multiply (shift) until }
+ bhi @mlabel8 { 1 in upper 16 result bits }
+ cmp.w #9,d0 { give up for denormalized numbers }
+ ble @mlabel8
+ swap d4 { (we''re getting here only when multiplying }
+ swap d5 { with a denormalized number; there''s an }
+ move.w d5,d4 { eventual loss of 4 bits in the rounding }
+ clr.w d5 { byte -- what a pity 8-) }
+ subq.w #8,d0 { decrement exponent }
+ subq.w #8,d0
+ bra @mlabel7
+@mlabel8:
+ move.l d5,d1 { get rounding bits }
+ rol.l #8,d1
+ move.l d1,d3 { see if sticky bit should be set }
+ and.l #$ffffff00,d3
+ beq @mlabel9
+ or.b #1,d1 { set "sticky bit" if any low-order set }
+@mlabel9:
+ addq.w #8,sp { remove accumulator from stack }
+ jmp FPC_SINGLE_NORM{ (result in d4) }
+
+@mretz:
+ addq.w #8,sp { release accumulator space }
+@mretzz:
+ moveq.l #0,d0 { save zero as result }
+ lsl.w #1,d2 { and set it sign as for d2 }
+ roxr.l #1,d0
+ movem.l (sp)+,d2-d5
+ rts { no normalizing neccessary }
+end;
+
+
+Procedure Single_Div;Assembler;
+{--------------------------------------------}
+{ Low-level routine to dividr two single }
+{ IEEE floating point values. Never called }
+{ directly. }
+{ Om Entry: }
+{ d1/d0 = u/v = operation to perform. }
+{ On Exit: }
+{ d0 = result. }
+{ Registers destroyed: d0,d1 }
+{ stack space used (and restored): 8 bytes. }
+{--------------------------------------------}
+ASM
+XDEF SINGLE_DIV
+ { u = d1 = dividend }
+ { v = d0 = divisor }
+ tst.l d0 { check if divisor is 0 }
+ bne @dno_exception
+
+ move.l #$7f800000,d0
+ btst #31,d1 { transfer sign of dividend }
+ beq @dclear
+ bset #31,d0
+@dclear:
+ rts
+
+@dno_exception:
+ move.l d1,d4 { d4 = u, d5 = v }
+ move.l d0,d5
+ movem.l d2-d5,-(sp) { save registers }
+
+ move.l #$7fffff,d3
+ move.l d4,d0 { d0 = u.exp }
+ and.l d3,d4 { remove exponent from u.mantissa }
+ swap d0
+ move.w d0,d2 { d2 = u.sign }
+
+ move.l d5,d1 { d1 = v.exp }
+ and.l d3,d5 { remove exponent from v.mantissa }
+ swap d1
+ eor.w d1,d2 { d2 = u.sign ^ v.sign (in bit 15) }
+
+ moveq.l #15,d3
+ bclr d3,d0 { kill sign bit }
+ bclr d3,d1 { kill sign bit }
+ lsr.w #7,d0
+ lsr.w #7,d1
+
+ moveq.l #-1,d3
+ cmp.b d3,d0 { comparison with #0xff }
+ beq @dlabel1 { u == NaN ;; u== Inf }
+ cmp.b d3,d1
+ beq @dlabel2 { v == NaN ;; v == Inf }
+ tst.b d0
+ bne @dlabel4 { u not zero nor denorm }
+ tst.l d4
+ beq @dlabel3 { 0/ ? }
+
+@dlabel4:
+ tst.w d1
+ bne @dnospec
+
+ tst.l d5
+ bne @dnospec
+ bra @dretinf { x/0 -> +/- Inf }
+
+@dlabel1:
+ tst.l d4 { u == NaN ? }
+ bne @dretnan { NaN/ x }
+ cmp.b d3,d1
+ beq @dretnan { Inf/Inf or Inf/NaN }
+{ bra dretinf ; Inf/x ; x != Inf && x != NaN }
+{--------------------------------------------------------------------}
+{ Return infinity with correct sign. }
+{--------------------------------------------------------------------}
+@dretinf:
+ move.l #$ff000000,d0
+ lsl.w #1,d2
+ roxr.l #1,d0 { shift in high bit as given by d2 }
+@dreturn:
+ movem.l (sp)+,d2-d5
+ rts
+
+@dlabel2:
+ tst.l d5
+ bne @dretnan { x/NaN }
+{ bra dretzero ; x/Inf -> +/- 0 }
+{--------------------------------------------------------------------}
+{ Return correct signed zero. }
+{--------------------------------------------------------------------}
+@dretzero:
+ moveq.l #0,d0 { zero destination }
+ lsl.w #1,d2 { set X bit accordingly }
+ roxr.l #1,d0
+ bra @dreturn
+
+@dlabel3:
+ tst.w d1
+ bne @dretzero { 0/x ->+/- 0 }
+ tst.l d4
+ bne @dretzero { 0/x }
+{ bra dretnan 0/0 }
+{--------------------------------------------------------------------}
+{ Return NotANumber }
+{--------------------------------------------------------------------}
+@dretnan:
+ move.l d3,d0 { d3 contains 0xffffffff }
+ lsr.l #1,d0
+ bra @dreturn
+{--------------------------------------------------------------------}
+{ End of Special Handling }
+{--------------------------------------------------------------------}
+@dnospec:
+ moveq.l #23,d3
+ bset d3,d4 { restore implied leading "1" }
+ tst.w d0 { check for zero exponent - no leading "1" }
+ bne @dlabel5
+ bclr d3,d4 { remove it }
+ add.w #1,d0 { "normalize" exponent }
+@dlabel5:
+ tst.l d4
+ beq @dretzero { dividing zero }
+
+ bset d3,d5 { restore implied leading "1" }
+ tst.w d1 { check for zero exponent - no leading "1"}
+ bne @dlabel6
+ bclr d3,d5 { remove it }
+ add.w #1,d1 { "normalize" exponent }
+@dlabel6:
+
+ sub.w d1,d0 { subtract exponents, }
+ add.w #BIAS4-8+1,d0 { add bias back in, account for shift }
+ add.w #34,d0 { add loop offset, +2 for extra rounding bits}
+ { for denormalized numbers (2 implied by dbra)}
+ move.w #27,d1 { bit number for "implied" pos (+4 for rounding)}
+ moveq.l #-1,d3 { zero quotient (for speed a one''s complement) }
+ sub.l d5,d4 { initial subtraction, u = u - v }
+@dlabel7:
+ btst d1,d3 { divide until 1 in implied position }
+ beq @dlabel9
+
+ add.l d4,d4
+ bcs @dlabel8 { if carry is set, add, else subtract }
+
+ addx.l d3,d3 { shift quotient and set bit zero }
+ sub.l d5,d4 { subtract, u = u - v }
+ dbra d0,@dlabel7 { give up if result is denormalized }
+ bra @dlabel9
+@dlabel8:
+ addx.l d3,d3 { shift quotient and clear bit zero }
+ add.l d5,d4 { add (restore), u = u + v }
+ dbra d0,@dlabel7 { give up if result is denormalized }
+@dlabel9:
+ subq.w #2,d0 { remove rounding offset for denormalized nums }
+ not.l d3 { invert quotient to get it right }
+
+ clr.l d1 { zero rounding bits }
+ tst.l d4 { check for exact result }
+ beq @dlabel10
+ moveq.l #-1,d1 { prevent tie case }
+@dlabel10:
+ move.l d3,d4 { save quotient mantissa }
+ jmp FPC_SINGLE_NORM{ (registers on stack removed by norm_sf) }
+end;
+
+
+Procedure Single_Cmp; Assembler;
+{--------------------------------------------}
+{ Low-level routine to compare single two }
+{ single point values.. }
+{ Never called directly. }
+{ On Entry: }
+{ d1 and d0 Values to compare }
+{ d0 = first operand }
+{ On Exit: }
+{ Flags according to result }
+{ Registers destroyed: d0,d1 }
+{--------------------------------------------}
+Asm
+XDEF SINGLE_CMP
+ tst.l d0 { check sign bit }
+ bpl @cmplab1
+ neg.l d0 { negate }
+ bchg #31,d0 { toggle sign bit }
+@cmplab1:
+ tst.l d1 { check sign bit }
+ bpl @cmplab2
+ neg.l d1 { negate }
+ bchg #31,d1 { toggle sign bit }
+@cmplab2:
+ cmp.l d0,d1 { compare... }
+ rts
+end;
+
+
+
+Procedure LongMul;Assembler;
+{--------------------------------------------}
+{ Low-level routine to multiply two signed }
+{ 32-bit values. Never called directly. }
+{ On entry: d1,d0 = 32-bit signed values to }
+{ multiply. }
+{ On Exit: }
+{ d0 = result. }
+{ Registers destroyed: d0,d1 }
+{ stack space used and restored: 10 bytes }
+{--------------------------------------------}
+Asm
+XDEF LONGMUL
+ cmp.b #2,Test68000 { Are we on a 68020+ cpu }
+ blt @Lmulcontinue
+ muls.l d1,d0 { yes, then directly mul... }
+ rts { return... result in d0 }
+@Lmulcontinue:
+ move.l d2,a0 { save registers }
+ move.l d3,a1
+
+ move.l d0,-(sp)
+ move.l d1,-(sp)
+
+ movem.w (sp)+,d0-d3 { u = d0-d1, v = d2-d3 }
+
+ move.w d0,-(sp) { sign flag }
+ bpl @LMul1 { is u negative ? }
+ neg.w d1 { yes, force it positive }
+ negx.w d0
+@LMul1:
+ tst.w d2 { is v negative ? }
+ bpl @LMul2
+ neg.w d3 { yes, force it positive ... }
+ negx.w d2
+ not.w (sp) { ... and modify flag word }
+@LMul2:
+ ext.l d0 { u.h <> 0 ? }
+ beq @LMul3
+ mulu.w d3,d0 { r = v.l * u.h }
+@LMul3:
+ tst.w d2 { v.h <> 0 ? }
+ beq @LMul4
+ mulu.w d1,d2 { r += v.h * u.l }
+ add.w d2,d0
+@LMul4:
+ swap d0
+ clr.w d0
+ mulu.w d3,d1 { r += v.l * u.l }
+ add.l d1,d0
+ move.l a1,d3
+ move.l a0,d2
+ tst.w (sp)+ { should the result be negated ? }
+ bpl @LMul5 { no, just return }
+ neg.l d0 { else r = -r }
+@LMul5:
+ rts
+end;
+
+
+
+Procedure Long2Single;Assembler;
+{--------------------------------------------}
+{ Low-level routine to convert a longint }
+{ to a single floating point value. }
+{ On entry: d0 = longint value to convert. }
+{ On Exit: }
+{ d0 = single IEEE value }
+{ Registers destroyed: d0,d1 }
+{ stack space used and restored: 8 bytes }
+{--------------------------------------------}
+Asm
+XDEF LONG2SINGLE
+ movem.l d2-d5,-(sp) { save registers to make norm_sf happy}
+
+ move.l d0,d4 { prepare result mantissa }
+ move.w #BIAS4+32-8,d0 { radix point after 32 bits }
+ move.l d4,d2 { set sign flag }
+ bge @l2slabel1 { nonnegative }
+ neg.l d4 { take absolute value }
+@l2slabel1:
+ swap d2 { follow SINGLE_NORM conventions }
+ clr.w d1 { set rounding = 0 }
+ jmp FPC_SINGLE_NORM
+end;
+
+
+Procedure LongDiv; [alias : 'FPC_LONGDIV'];Assembler;
+{--------------------------------------------}
+{ Low-level routine to do signed long }
+{ division. }
+{ On entry: d0/d1 operation to perform }
+{ On Exit: }
+{ d0 = quotient }
+{ d1 = remainder }
+{ Registers destroyed: d0,d1,d6 }
+{ stack space used and restored: 10 bytes }
+{--------------------------------------------}
+asm
+XDEF LONGDIV
+ cmp.b #2,Test68000 { can we use divs ? }
+ blt @continue
+ tst.l d1
+ beq @zerodiv2
+ move.l d1,d6
+ clr.l d1 { clr }
+ tst.l d0 { check sign of d0 }
+ bpl @posdiv
+ move.l #$ffffffff,d1{ sign extend into d1 }
+@posdiv:
+ divsl.l d6,d1:d0
+ rts
+@continue:
+
+ move.l d2,a0 { save registers }
+ move.l d3,a1
+ move.l d4,-(sp) { divisor = d1 = d4 }
+ move.l d5,-(sp) { divident = d0 = d5 }
+
+ move.l d1,d4 { save divisor }
+ move.l d0,d5 { save dividend }
+
+ clr.w -(sp) { sign flag }
+
+ clr.l d0 { prepare result }
+ move.l d4,d2 { get divisor }
+ beq @zerodiv { divisor = 0 ? }
+ bpl @LDiv1 { divisor < 0 ? }
+ neg.l d2 { negate it }
+ not.w (sp) { remember sign }
+@LDiv1:
+ move.l d5,d1 { get dividend }
+ bpl @LDiv2 { dividend < 0 ? }
+ neg.l d1 { negate it }
+ not.w (sp) { remember sign }
+@LDiv2:
+{;== case 1) divident < divisor}
+ cmp.l d2,d1 { is divident smaller then divisor ? }
+ bcs @LDiv7 { yes, return immediately }
+{;== case 2) divisor has <= 16 significant bits}
+ move.l d4,d6 { put divisor in d6 register }
+ lsr.l #8,d6 { rotate into low word }
+ lsr.l #8,d6
+ tst.l d6
+ bne @LDiv3 { divisor has only 16 bits }
+ move.w d1,d3 { save dividend }
+ clr.w d1 { divide dvd.h by dvs }
+ swap d1
+ beq @LDiv4 { (no division necessary if dividend zero)}
+ divu d2,d1
+@LDiv4:
+ move.w d1,d0 { save quotient.h }
+ swap d0
+ move.w d3,d1 { (d0.h = remainder of prev divu) }
+ divu d2,d1 { divide dvd.l by dvs }
+ move.w d1,d0 { save quotient.l }
+ clr.w d1 { get remainder }
+ swap d1
+ bra @LDiv7 { and return }
+{;== case 3) divisor > 16 bits (corollary is dividend > 16 bits, see case 1)}
+@LDiv3:
+ moveq.l #31,d3 { loop count }
+@LDiv5:
+ add.l d1,d1 { shift divident ... }
+ addx.l d0,d0 { ... into d0 }
+ cmp.l d2,d0 { compare with divisor }
+ bcs @LDiv6
+ sub.l d2,d0 { big enough, subtract }
+ addq.w #1,d1 { and note bit into result }
+@LDiv6:
+ dbra d3,@LDiv5
+ exg d0,d1 { put quotient and remainder in their registers}
+@LDiv7:
+ tst.l d5 { must the remainder be corrected ? }
+ bpl @LDiv8
+ neg.l d1 { yes, apply sign }
+{ the following line would be correct if modulus is defined as in algebra}
+{; add.l sp@(6),d1 ; algebraic correction: modulus can only be >= 0}
+@LDiv8:
+ tst.w (sp)+ { result should be negative ? }
+ bpl @LDiv9
+ neg.l d0 { yes, negate it }
+@LDiv9:
+ move.l a1,d3
+ move.l a0,d2
+ move.l (sp)+,d5
+ move.l (sp)+,d4
+ rts { en exit : remainder = d1, quotient = d0 }
+@zerodiv:
+ move.l a1,d3 { restore stack... }
+ move.l a0,d2
+ move.w (sp)+,d1 { remove sign word }
+ move.l (sp)+,d5
+ move.l (sp)+,d4
+@zerodiv2:
+ move.l #200,d0
+ jsr FPC_HALT_ERROR { RunError(200) }
+ rts { this should never occur... }
+end;
+
+
+Procedure LongMod;[alias : 'FPC_LONGMOD'];Assembler;
+{ see longdiv for info on calling convention }
+asm
+XDEF LONGMOD
+ jsr FPC_LONGDIV
+ move.l d1,d0 { return the remainder in d0 }
+ rts
+end;
+
+{
+ $Log: lowmath.inc,v $
+ Revision 1.5 2005/02/14 17:13:30 peter
+ * truncate log
+
+}
diff --git a/rtl/m68k/m68k.inc b/rtl/m68k/m68k.inc
new file mode 100644
index 0000000000..4e98dcde26
--- /dev/null
+++ b/rtl/m68k/m68k.inc
@@ -0,0 +1,329 @@
+{
+ $Id: m68k.inc,v 1.6 2005/02/14 17:13:30 peter Exp $
+ This file is part of the Free Pascal run time library.
+ Copyright (c) 1999-2000 by Carl-Eric Codere,
+ member of the Free Pascal development team.
+
+ See the file COPYING.FPC, included in this distribution,
+ for details about the copyright.
+
+ This program 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.
+
+ **********************************************************************}
+{****************************************************************************
+
+ m68k.inc : Processor dependent implementation of system unit
+ For Motorola 680x0 Processor.
+
+*****************************************************************************}
+
+{****************************************************************************}
+{ Credit where credit is due: }
+{ -Some of the copy routines taken from the Atari dlib source code: }
+{ Dale Schumacher (alias: Dalnefre') dal@syntel.uucp }
+{ 399 Beacon Ave. St. Paul, MN 55104,USA }
+{ -Some of the routines taken from the freeware ATARI Sozobon C compiler }
+{ 1988 by Sozobon, Limited. Author: Johann Ruegg (freeware) }
+{ Thanks to all these people wherever they maybe today! }
+{****************************************************************************}
+
+
+procedure fpc_cpuinit;
+ begin
+ end;
+
+{$define FPC_SYSTEM_HAS_GET_FRAME}
+function get_frame : pointer; assembler;
+ asm
+ move.l a6,d0
+ end;
+
+
+{$define FPC_SYSTEM_HAS_GET_CALLER_ADDR}
+function get_caller_addr(framebp : pointer) : pointer;
+ begin
+ asm
+ move.l FRAMEBP,a0
+ cmp.l #0,a0
+ beq @Lnul_address
+ move.l 4(a0),a0
+ @Lnul_address:
+ move.l a0,@RESULT
+ end ['a0'];
+ end;
+
+
+{$define FPC_SYSTEM_HAS_GET_CALLER_FRAME}
+function get_caller_frame(framebp : pointer) : pointer;
+ begin
+ asm
+ move.l FRAMEBP,a0
+ cmp.l #0,a0
+ beq @Lnul_frame
+ move.l (a0),a0
+ @Lnul_frame:
+ move.l a0,@RESULT
+ end ['a0'];
+ end;
+
+
+{$define FPC_SYSTEM_HAS_SPTR}
+function Sptr : Longint;
+ begin
+ asm
+ move.l sp,d0
+ add.l #8,d0
+ move.l d0,@RESULT
+ end ['d0'];
+ end;
+
+
+{$define FPC_SYSTEM_HAS_FILLCHAR}
+procedure FillChar(var x;count:longint;value:byte);[public,alias: 'FPC_FILL_OBJECT'];
+ begin
+ asm
+ move.l 8(a6), a0 { destination }
+ move.l 12(a6), d1 { number of bytes to fill }
+ move.b 16(a6),d0 { fill data }
+ cmpi.l #65535, d1 { check, if this is a word move }
+ ble @LMEMSET3 { use fast dbra mode }
+ bra @LMEMSET2
+ @LMEMSET1:
+ move.b d0,(a0)+
+ @LMEMSET2:
+ subq.l #1,d1
+ cmp.l #-1,d1
+ bne @LMEMSET1
+ bra @LMEMSET5 { finished slow mode , exit }
+
+ @LMEMSET4: { fast loop mode section 68010+ }
+ move.b d0,(a0)+
+ @LMEMSET3:
+ dbra d1,@LMEMSET4
+
+ @LMEMSET5:
+ end ['d0','d1','a0'];
+ end;
+
+
+{$ifdef dummy}
+{ procedure strcopy(dstr,sstr : pointer;len : longint);[public,alias: 'STRCOPY'];}
+procedure strcopy; assembler;[public,alias: 'FPC_STRCOPY'];
+{---------------------------------------------------}
+{ Low-level routine to copy a string to another }
+{ string with maximum length. Never call directly! }
+{ On Entry: }
+{ a1.l = string to copy to }
+{ a0.l = source string }
+{ d0.l = maximum length of copy }
+{ registers destroyed: a0,a1,d0,d1 }
+{---------------------------------------------------}
+asm
+{ move.l 12(a6),a0
+ move.l 16(a6),a1
+ move.l 8(a6),d1 }
+ move.l d0,d1
+
+ move.b (a0)+,d0 { Get source length }
+ and.w #$ff,d0
+ cmp.w d1,d0 { This is a signed comparison! }
+ ble @LM4
+ move.b d1,d0 { If longer than maximum size of target, cut
+ source length }
+@LM4:
+ andi.l #$ff,d0 { zero extend d0-byte }
+ move.l d0,d1 { save length to copy }
+ move.b d0,(a1)+ { save new length }
+ { Check if copying length is zero - if so then }
+ { exit without copying anything. }
+ tst.b d1
+ beq @Lend
+ bra @LMSTRCOPY55
+@LMSTRCOPY56: { 68010 Fast loop mode }
+ move.b (a0)+,(a1)+
+@LMSTRCOPY55:
+ dbra d1,@LMSTRCOPY56
+@Lend:
+end;
+
+
+{ Concatenate Strings }
+{ PARAMETERS ARE REVERSED COMPARED TO NORMAL! }
+{ therefore online assembler may not parse the params as normal }
+procedure strconcat(s1,s2 : pointer);[public,alias: 'STRCONCAT'];
+ begin
+ asm
+ move.b #255,d0
+ move.l s1,a0 { a0 = destination }
+ move.l s2,a1 { a1 = source }
+ sub.b (a0),d0 { copyl:= 255 -length(s1) }
+ move.b (a1),d6
+ and.w #$ff,d0 { Sign flags are checked! }
+ and.w #$ff,d6
+ cmp.w d6,d0 { if copyl > length(s2) then }
+ ble @Lcontinue
+ move.b (a1),d0 { copyl:=length(s2) }
+@Lcontinue:
+ move.b (a0),d6
+ and.l #$ff,d6
+ lea 1(a0,d6),a0 { s1[length(s1)+1] }
+ add.l #1,a1 { s2[1] }
+ move.b d0,d6
+ { Check if copying length is zero - if so then }
+ { exit without copying anything. }
+ tst.b d6
+ beq @Lend
+ bra @ALoop
+@Loop:
+ move.b (a1)+,(a0)+ { s1[i] := s2[i]; }
+@ALoop:
+ dbra d6,@Loop
+ move.l s1,a0
+ add.b d0,(a0) { change to new string length }
+@Lend:
+ end ['d0','d1','a0','a1','d6'];
+ end;
+
+{ Compares strings }
+{ DO NOT CALL directly. }
+{ a0 = pointer to first string to compare }
+{ a1 = pointer to second string to compare }
+{ ALL FLAGS are set appropriately. }
+{ ZF = strings are equal }
+{ REGISTERS DESTROYED: a0, a1, d0, d1, d6 }
+procedure strcmp; assembler;[public,alias:'FPC_STRCMP'];
+asm
+ move.b (a0)+,d0 { Get length of first string }
+ move.b (a1)+,d6 { Get length of 2nd string }
+
+ move.b d6,d1 { Save length of string for final compare }
+
+ cmp.b d0,d6 { Get shortest string length }
+ ble @LSTRCONCAT1
+ move.b d0,d6 { Set length to shortest string }
+
+ @LSTRCONCAT1:
+ tst.b d6 { Both strings have a length of zero, exit }
+ beq @LSTRCONCAT2
+
+ andi.l #$ff,d6
+
+
+ subq.l #1,d6 { subtract first attempt }
+ { if value is -1 then don't loop and just compare lengths of }
+ { both strings before exiting. }
+ bmi @LSTRCONCAT2
+ or.l d0,d0 { Make sure to set Zerfo flag to 0 }
+ @LSTRCONCAT5:
+ { Workaroung for GAS v.134 bug }
+ { old: cmp.b (a1)+,(a0)+ }
+ cmpm.b (a1)+,(a0)+
+ @LSTRCONCAT4:
+ dbne d6,@LSTRCONCAT5 { Repeat until not equal }
+ bne @LSTRCONCAT3
+ @LSTRCONCAT2:
+ { If length of both string are equal }
+ { Then set zero flag }
+ cmp.b d1,d0 { Compare length - set flag if equal length strings }
+ @LSTRCONCAT3:
+end;
+{$endif dummy}
+
+
+{$define FPC_SYSTEM_HAS_MOVE}
+procedure move(var source;var dest;count : longint);
+{ base pointer+8 = source }
+{ base pointer+12 = destination }
+{ base pointer+16 = number of bytes to move}
+begin
+ asm
+ clr.l d0
+ move.l 16(a6),d0 { number of bytes }
+ @LMOVE0:
+ move.l 12(a6),a1 { destination }
+ move.l 8(a6),a0 { source }
+
+ cmpi.l #65535, d0 { check, if this is a word move }
+ ble @LMEMSET00 { use fast dbra mode 68010+ }
+
+ cmp.l a0,a1 { check copy direction }
+ bls @LMOVE4
+ add.l d0,a0 { move pointers to end }
+ add.l d0,a1
+ bra @LMOVE2
+ @LMOVE1:
+ move.b -(a0),-(a1) { (s < d) copy loop }
+ @LMOVE2:
+ subq.l #1,d0
+ cmpi.l #-1,d0
+ bne @LMOVE1
+ bra @LMOVE5
+ @LMOVE3:
+ move.b (a0)+,(a1)+ { (s >= d) copy loop }
+ @LMOVE4:
+ subq.l #1,d0
+ cmpi.l #-1,d0
+ bne @LMOVE3
+ bra @LMOVE5
+
+ @LMEMSET00: { use fast loop mode 68010+ }
+ cmp.l a0,a1 { check copy direction }
+ bls @LMOVE04
+ add.l d0,a0 { move pointers to end }
+ add.l d0,a1
+ bra @LMOVE02
+ @LMOVE01:
+ move.b -(a0),-(a1) { (s < d) copy loop }
+ @LMOVE02:
+ dbra d0,@LMOVE01
+ bra @LMOVE5
+ @LMOVE03:
+ move.b (a0)+,(a1)+ { (s >= d) copy loop }
+ @LMOVE04:
+ dbra d0,@LMOVE03
+ { end fast loop mode }
+ @LMOVE5:
+ end ['d0','a0','a1'];
+end;
+
+
+{$define FPC_SYSTEM_HAS_FILLWORD}
+procedure fillword(var x;count : longint;value : word);
+ begin
+ asm
+ move.l 8(a6), a0 { destination }
+ move.l 12(a6), d1 { number of bytes to fill }
+ move.w 16(a6),d0 { fill data }
+ bra @LMEMSET21
+ @LMEMSET11:
+ move.w d0,(a0)+
+ @LMEMSET21:
+ subq.l #1,d1
+ cmp.b #-1,d1
+ bne @LMEMSET11
+ end ['d0','d1','a0'];
+ end;
+
+
+{$define FPC_SYSTEM_HAS_ABS_LONGINT}
+function abs(l : longint) : longint;
+ begin
+ asm
+ move.l 8(a6),d0
+ tst.l d0
+ bpl @LMABS1
+ neg.l d0
+ @LMABS1:
+ move.l d0,@RESULT
+ end ['d0'];
+ end;
+
+
+{
+ $Log: m68k.inc,v $
+ Revision 1.6 2005/02/14 17:13:30 peter
+ * truncate log
+
+}
diff --git a/rtl/m68k/makefile.cpu b/rtl/m68k/makefile.cpu
new file mode 100644
index 0000000000..be1724efeb
--- /dev/null
+++ b/rtl/m68k/makefile.cpu
@@ -0,0 +1,7 @@
+#
+# Here we set processor dependent include file names.
+#
+
+CPUNAMES=m68k lowmath math set
+CPUINCNAMES=$(addsuffix .inc,$(CPUNAMES))
+
diff --git a/rtl/m68k/math.inc b/rtl/m68k/math.inc
new file mode 100644
index 0000000000..c2608289ee
--- /dev/null
+++ b/rtl/m68k/math.inc
@@ -0,0 +1,949 @@
+{
+ $Id: math.inc,v 1.4 2005/02/14 17:13:30 peter Exp $
+ This file is part of the Free Pascal run time library.
+ Copyright (c) 1999-2000 by several people
+ member of the Free Pascal development team.
+
+ See the file COPYING.FPC, included in this distribution,
+ for details about the copyright.
+
+ This program 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.
+
+ **********************************************************************}
+{*************************************************************************}
+{ math.inc }
+{*************************************************************************}
+{ Copyright Abandoned, 1987, Fred Fish }
+{ }
+{ This previously copyrighted work has been placed into the }
+{ public domain by the author (Fred Fish) and may be freely used }
+{ for any purpose, private or commercial. I would appreciate }
+{ it, as a courtesy, if this notice is left in all copies and }
+{ derivative works. Thank you, and enjoy... }
+{ }
+{ The author makes no warranty of any kind with respect to this }
+{ product and explicitly disclaims any implied warranties of }
+{ merchantability or fitness for any particular purpose. }
+{-------------------------------------------------------------------------}
+{ Copyright (c) 1992 Odent Jean Philippe }
+{ }
+{ The source can be modified as long as my name appears and some }
+{ notes explaining the modifications done are included in the file. }
+{-------------------------------------------------------------------------}
+{ Copyright (c) 1997 Carl Eric Codere }
+{ }
+{*************************************************************************}
+{ This is the Motorola 680x0 specific port of the math include. }
+{*************************************************************************}
+{ }
+{ o all reals are mapped to the single type under the motorola version }
+{ }
+{ What is left to do: }
+{ o add support for sqrt with fixed. }
+
+type
+ TabCoef = array[0..6] of Real;
+
+
+const
+ PIO2 = 1.57079632679489661923; { pi/2 }
+ PIO4 = 7.85398163397448309616E-1; { pi/4 }
+ SQRT2 = 1.41421356237309504880; { sqrt(2) }
+ SQRTH = 7.07106781186547524401E-1; { sqrt(2)/2 }
+ LOG2E = 1.4426950408889634073599; { 1/log(2) }
+ SQ2OPI = 7.9788456080286535587989E-1; { sqrt( 2/pi )}
+ LOGE2 = 6.93147180559945309417E-1; { log(2) }
+ LOGSQ2 = 3.46573590279972654709E-1; { log(2)/2 }
+ THPIO4 = 2.35619449019234492885; { 3*pi/4 }
+ TWOOPI = 6.36619772367581343075535E-1; { 2/pi }
+ lossth = 1.073741824e9;
+ MAXLOG = 8.8029691931113054295988E1; { log(2**127) }
+ MINLOG = -8.872283911167299960540E1; { log(2**-128) }
+
+ DP1 = 7.85398125648498535156E-1;
+ DP2 = 3.77489470793079817668E-8;
+ DP3 = 2.69515142907905952645E-15;
+
+const sincof : TabCoef = (
+ 1.58962301576546568060E-10,
+ -2.50507477628578072866E-8,
+ 2.75573136213857245213E-6,
+ -1.98412698295895385996E-4,
+ 8.33333333332211858878E-3,
+ -1.66666666666666307295E-1, 0);
+ coscof : TabCoef = (
+ -1.13585365213876817300E-11,
+ 2.08757008419747316778E-9,
+ -2.75573141792967388112E-7,
+ 2.48015872888517045348E-5,
+ -1.38888888888730564116E-3,
+ 4.16666666666665929218E-2, 0);
+
+
+
+
+ function int(d : real) : real;
+ begin
+ { this will be correct since real = single in the case of }
+ { the motorola version of the compiler... }
+ int:=real(trunc(d));
+ end;
+
+ function trunc(d : real) : longint;
+ var
+ l: longint;
+ Begin
+ asm
+ move.l d,d0 { get number }
+ move.l d2,-(sp) { save register }
+ move.l d0,d1
+ swap d1 { extract exp }
+ move.w d1,d2 { extract sign }
+ bclr #15,d1 { kill sign bit }
+ lsr.w #7,d1
+
+ and.l #$7fffff,d0 { remove exponent from mantissa }
+ bset #23,d0 { restore implied leading "1" }
+
+ cmp.w #BIAS4,d1 { check exponent }
+ blt @zero { strictly factional, no integer part ? }
+ cmp.w #BIAS4+32,d1 { is it too big to fit in a 32-bit integer ? }
+ bgt @toobig
+
+ sub.w #BIAS4+24,d1 { adjust exponent }
+ bgt @trunclab2 { shift up }
+ beq @trunclab7 { no shift (never too big) }
+
+ neg.w d1
+ lsr.l d1,d0 { shift down to align radix point; }
+ { extra bits fall off the end (no rounding) }
+ bra @trunclab7 { never too big }
+ @trunclab2:
+ lsl.l d1,d0 { shift up to align radix point }
+ @trunclab3:
+ cmp.l #$80000000,d0 { -2147483648 is a nasty evil special case }
+ bne @trunclab6
+ tst.w d2 { this had better be -2^31 and not 2^31 }
+ bpl @toobig
+ bra @trunclab8
+ @trunclab6:
+ tst.l d0 { sign bit set ? (i.e. too big) }
+ bmi @toobig
+ @trunclab7:
+ tst.w d2 { is it negative ? }
+ bpl @trunclab8
+ neg.l d0 { negate }
+ bra @trunclab8
+ @zero:
+ clr.l d0 { make the whole thing zero }
+ bra @trunclab8
+ @toobig:
+ moveq #-1,d0 { ugh. Should cause a trap here. }
+ bclr #31,d0 { make it #0x7fffffff }
+ @trunclab8:
+ move.l (sp)+,d2
+ move.l d0,l
+ end;
+ if l = $7fffffff then
+ RunError(207)
+ else
+ trunc := l
+ end;
+
+
+
+
+
+ function abs(d : Real) : Real;
+ begin
+ if( d < 0.0 ) then
+ abs := -d
+ else
+ abs := d ;
+ end;
+
+
+ function frexp(x:Real; var e:Integer ):Real;
+ {* frexp() extracts the exponent from x. It returns an integer *}
+ {* power of two to expnt and the significand between 0.5 and 1 *}
+ {* to y. Thus x = y * 2**expn. *}
+ begin
+ e :=0;
+ if (abs(x)<0.5) then
+ While (abs(x)<0.5) do
+ begin
+ x := x*2;
+ Dec(e);
+ end
+ else
+ While (abs(x)>1) do
+ begin
+ x := x/2;
+ Inc(e);
+ end;
+ frexp := x;
+ end;
+
+
+ function ldexp( x: Real; N: Integer):Real;
+ {* ldexp() multiplies x by 2**n. *}
+ var r : Real;
+ begin
+ R := 1;
+ if N>0 then
+ while N>0 do
+ begin
+ R:=R*2;
+ Dec(N);
+ end
+ else
+ while N<0 do
+ begin
+ R:=R/2;
+ Inc(N);
+ end;
+ ldexp := x * R;
+ end;
+
+
+ function polevl(var x:Real; var Coef:TabCoef; N:Integer):Real;
+ {*****************************************************************}
+ { Evaluate polynomial }
+ {*****************************************************************}
+ { }
+ { SYNOPSIS: }
+ { }
+ { int N; }
+ { double x, y, coef[N+1], polevl[]; }
+ { }
+ { y = polevl( x, coef, N ); }
+ { }
+ { DESCRIPTION: }
+ { }
+ { Evaluates polynomial of degree N: }
+ { }
+ { 2 N }
+ { y = C + C x + C x +...+ C x }
+ { 0 1 2 N }
+ { }
+ { Coefficients are stored in reverse order: }
+ { }
+ { coef[0] = C , ..., coef[N] = C . }
+ { N 0 }
+ { }
+ { The function p1evl() assumes that coef[N] = 1.0 and is }
+ { omitted from the array. Its calling arguments are }
+ { otherwise the same as polevl(). }
+ { }
+ { SPEED: }
+ { }
+ { In the interest of speed, there are no checks for out }
+ { of bounds arithmetic. This routine is used by most of }
+ { the functions in the library. Depending on available }
+ { equipment features, the user may wish to rewrite the }
+ { program in microcode or assembly language. }
+ {*****************************************************************}
+ var ans : Real;
+ i : Integer;
+
+ begin
+ ans := Coef[0];
+ for i:=1 to N do
+ ans := ans * x + Coef[i];
+ polevl:=ans;
+ end;
+
+
+ function p1evl(var x:Real; var Coef:TabCoef; N:Integer):Real;
+ { }
+ { Evaluate polynomial when coefficient of x is 1.0. }
+ { Otherwise same as polevl. }
+ { }
+ var
+ ans : Real;
+ i : Integer;
+ begin
+ ans := x + Coef[0];
+ for i:=1 to N-1 do
+ ans := ans * x + Coef[i];
+ p1evl := ans;
+ end;
+
+
+
+
+
+ function sqr(d : Real) : Real;
+ begin
+ sqr := d*d;
+ end;
+
+
+ function pi : Real;
+ begin
+ pi := 3.1415926535897932385;
+ end;
+
+
+ function sqrt(d:Real):Real;
+ {*****************************************************************}
+ { Square root }
+ {*****************************************************************}
+ { }
+ { SYNOPSIS: }
+ { }
+ { double x, y, sqrt(); }
+ { }
+ { y = sqrt( x ); }
+ { }
+ { DESCRIPTION: }
+ { }
+ { Returns the square root of x. }
+ { }
+ { Range reduction involves isolating the power of two of the }
+ { argument and using a polynomial approximation to obtain }
+ { a rough value for the square root. Then Heron's iteration }
+ { is used three times to converge to an accurate value. }
+ {*****************************************************************}
+ var e : Integer;
+ w,z : Real;
+ begin
+ if( d <= 0.0 ) then
+ begin
+ if( d < 0.0 ) then
+ RunError(207);
+ sqrt := 0.0;
+ end
+ else
+ begin
+ w := d;
+ { separate exponent and significand }
+ z := frexp( d, e );
+
+ { approximate square root of number between 0.5 and 1 }
+ { relative error of approximation = 7.47e-3 }
+ d := 4.173075996388649989089E-1 + 5.9016206709064458299663E-1 * z;
+
+ { adjust for odd powers of 2 }
+ if odd(e) then
+ d := d*SQRT2;
+
+ { re-insert exponent }
+ d := ldexp( d, (e div 2) );
+
+ { Newton iterations: }
+ d := 0.5*(d + w/d);
+ d := 0.5*(d + w/d);
+ d := 0.5*(d + w/d);
+ d := 0.5*(d + w/d);
+ d := 0.5*(d + w/d);
+ d := 0.5*(d + w/d);
+ sqrt := d;
+ end;
+ end;
+
+
+
+
+ function Exp(d:Real):Real;
+ {*****************************************************************}
+ { Exponential Function }
+ {*****************************************************************}
+ { }
+ { SYNOPSIS: }
+ { }
+ { double x, y, exp(); }
+ { }
+ { y = exp( x ); }
+ { }
+ { DESCRIPTION: }
+ { }
+ { Returns e (2.71828...) raised to the x power. }
+ { }
+ { Range reduction is accomplished by separating the argument }
+ { into an integer k and fraction f such that }
+ { }
+ { x k f }
+ { e = 2 e. }
+ { }
+ { A Pade' form of degree 2/3 is used to approximate exp(f)- 1 }
+ { in the basic range [-0.5 ln 2, 0.5 ln 2]. }
+ {*****************************************************************}
+ const P : TabCoef = (
+ 1.26183092834458542160E-4,
+ 3.02996887658430129200E-2,
+ 1.00000000000000000000E0, 0, 0, 0, 0);
+ Q : TabCoef = (
+ 3.00227947279887615146E-6,
+ 2.52453653553222894311E-3,
+ 2.27266044198352679519E-1,
+ 2.00000000000000000005E0, 0 ,0 ,0);
+
+ C1 = 6.9335937500000000000E-1;
+ C2 = 2.1219444005469058277E-4;
+ var n : Integer;
+ px, qx, xx : Real;
+ begin
+ if( d > MAXLOG) then
+ RunError(205)
+ else
+ if( d < MINLOG ) then
+ begin
+ Runerror(205);
+ end
+ else
+ begin
+
+ { Express e**x = e**g 2**n }
+ { = e**g e**( n loge(2) ) }
+ { = e**( g + n loge(2) ) }
+
+ px := d * LOG2E;
+ qx := Trunc( px + 0.5 ); { Trunc() truncates toward -infinity. }
+ n := Trunc(qx);
+ d := d - qx * C1;
+ d := d + qx * C2;
+
+ { rational approximation for exponential }
+ { of the fractional part: }
+ { e**x - 1 = 2x P(x**2)/( Q(x**2) - P(x**2) ) }
+ xx := d * d;
+ px := d * polevl( xx, P, 2 );
+ d := px/( polevl( xx, Q, 3 ) - px );
+ d := ldexp( d, 1 );
+ d := d + 1.0;
+ d := ldexp( d, n );
+ Exp := d;
+ end;
+ end;
+
+
+ function Round(d: Real): longint;
+ var
+ fr: Real;
+ tr: Real;
+ Begin
+ fr := Frac(d);
+ tr := Trunc(d);
+ if fr > 0.5 then
+ Round:=Trunc(d)+1
+ else
+ if fr < 0.5 then
+ Round:=Trunc(d)
+ else { fr = 0.5 }
+ { check sign to decide ... }
+ { as in Turbo Pascal... }
+ if d >= 0.0 then
+ Round := Trunc(d)+1
+ else
+ Round := Trunc(d);
+ end;
+
+
+ function Ln(d:Real):Real;
+ {*****************************************************************}
+ { Natural Logarithm }
+ {*****************************************************************}
+ { }
+ { SYNOPSIS: }
+ { }
+ { double x, y, log(); }
+ { }
+ { y = ln( x ); }
+ { }
+ { DESCRIPTION: }
+ { }
+ { Returns the base e (2.718...) logarithm of x. }
+ { }
+ { The argument is separated into its exponent and fractional }
+ { parts. If the exponent is between -1 and +1, the logarithm }
+ { of the fraction is approximated by }
+ { }
+ { log(1+x) = x - 0.5 x**2 + x**3 P(x)/Q(x). }
+ { }
+ { Otherwise, setting z = 2(x-1)/x+1), }
+ { }
+ { log(x) = z + z**3 P(z)/Q(z). }
+ { }
+ {*****************************************************************}
+ const P : TabCoef = (
+ { Coefficients for log(1+x) = x - x**2/2 + x**3 P(x)/Q(x)
+ 1/sqrt(2) <= x < sqrt(2) }
+
+ 4.58482948458143443514E-5,
+ 4.98531067254050724270E-1,
+ 6.56312093769992875930E0,
+ 2.97877425097986925891E1,
+ 6.06127134467767258030E1,
+ 5.67349287391754285487E1,
+ 1.98892446572874072159E1);
+ Q : TabCoef = (
+ 1.50314182634250003249E1,
+ 8.27410449222435217021E1,
+ 2.20664384982121929218E2,
+ 3.07254189979530058263E2,
+ 2.14955586696422947765E2,
+ 5.96677339718622216300E1, 0);
+
+ { Coefficients for log(x) = z + z**3 P(z)/Q(z),
+ where z = 2(x-1)/(x+1)
+ 1/sqrt(2) <= x < sqrt(2) }
+
+ R : TabCoef = (
+ -7.89580278884799154124E-1,
+ 1.63866645699558079767E1,
+ -6.41409952958715622951E1, 0, 0, 0, 0);
+ S : TabCoef = (
+ -3.56722798256324312549E1,
+ 3.12093766372244180303E2,
+ -7.69691943550460008604E2, 0, 0, 0, 0);
+
+ var e : Integer;
+ z, y : Real;
+
+ Label Ldone;
+ begin
+ if( d <= 0.0 ) then
+ RunError(207);
+ d := frexp( d, e );
+
+ { logarithm using log(x) = z + z**3 P(z)/Q(z),
+ where z = 2(x-1)/x+1) }
+
+ if( (e > 2) or (e < -2) ) then
+ begin
+ if( d < SQRTH ) then
+ begin
+ { 2( 2x-1 )/( 2x+1 ) }
+ Dec(e, 1);
+ z := d - 0.5;
+ y := 0.5 * z + 0.5;
+ end
+ else
+ begin
+ { 2 (x-1)/(x+1) }
+ z := d - 0.5;
+ z := z - 0.5;
+ y := 0.5 * d + 0.5;
+ end;
+ d := z / y;
+ { /* rational form */ }
+ z := d*d;
+ z := d + d * ( z * polevl( z, R, 2 ) / p1evl( z, S, 3 ) );
+ goto ldone;
+ end;
+
+ { logarithm using log(1+x) = x - .5x**2 + x**3 P(x)/Q(x) }
+
+ if( d < SQRTH ) then
+ begin
+ Dec(e, 1);
+ d := ldexp( d, 1 ) - 1.0; { 2x - 1 }
+ end
+ else
+ d := d - 1.0;
+
+ { rational form }
+ z := d*d;
+ y := d * ( z * polevl( d, P, 6 ) / p1evl( d, Q, 6 ) );
+ y := y - ldexp( z, -1 ); { y - 0.5 * z }
+ z := d + y;
+
+ ldone:
+ { recombine with exponent term }
+ if( e <> 0 ) then
+ begin
+ y := e;
+ z := z - y * 2.121944400546905827679e-4;
+ z := z + y * 0.693359375;
+ end;
+
+ Ln:= z;
+ end;
+
+
+
+ function Sin(d:Real):Real;
+ {*****************************************************************}
+ { Circular Sine }
+ {*****************************************************************}
+ { }
+ { SYNOPSIS: }
+ { }
+ { double x, y, sin(); }
+ { }
+ { y = sin( x ); }
+ { }
+ { DESCRIPTION: }
+ { }
+ { Range reduction is into intervals of pi/4. The reduction }
+ { error is nearly eliminated by contriving an extended }
+ { precision modular arithmetic. }
+ { }
+ { Two polynomial approximating functions are employed. }
+ { Between 0 and pi/4 the sine is approximated by }
+ { x + x**3 P(x**2). }
+ { Between pi/4 and pi/2 the cosine is represented as }
+ { 1 - x**2 Q(x**2). }
+ {*****************************************************************}
+ var y, z, zz : Real;
+ j, sign : Integer;
+
+ begin
+ { make argument positive but save the sign }
+ sign := 1;
+ if( d < 0 ) then
+ begin
+ d := -d;
+ sign := -1;
+ end;
+
+ { above this value, approximate towards 0 }
+ if( d > lossth ) then
+ begin
+ sin := 0.0;
+ exit;
+ end;
+
+ y := Trunc( d/PIO4 ); { integer part of x/PIO4 }
+
+ { strip high bits of integer part to prevent integer overflow }
+ z := ldexp( y, -4 );
+ z := Trunc(z); { integer part of y/8 }
+ z := y - ldexp( z, 4 ); { y - 16 * (y/16) }
+
+ j := Trunc(z); { convert to integer for tests on the phase angle }
+ { map zeros to origin }
+ if odd( j ) then
+ begin
+ inc(j);
+ y := y + 1.0;
+ end;
+ j := j and 7; { octant modulo 360 degrees }
+ { reflect in x axis }
+ if( j > 3) then
+ begin
+ sign := -sign;
+ dec(j, 4);
+ end;
+
+ { Extended precision modular arithmetic }
+ z := ((d - y * DP1) - y * DP2) - y * DP3;
+
+ zz := z * z;
+
+ if( (j=1) or (j=2) ) then
+ y := 1.0 - ldexp(zz,-1) + zz * zz * polevl( zz, coscof, 5 )
+ else
+ { y = z + z * (zz * polevl( zz, sincof, 5 )); }
+ y := z + z * z * z * polevl( zz, sincof, 5 );
+
+ if(sign < 0) then
+ y := -y;
+ sin := y;
+ end;
+
+
+
+
+ function Cos(d:Real):Real;
+ {*****************************************************************}
+ { Circular cosine }
+ {*****************************************************************}
+ { }
+ { Circular cosine }
+ { }
+ { SYNOPSIS: }
+ { }
+ { double x, y, cos(); }
+ { }
+ { y = cos( x ); }
+ { }
+ { DESCRIPTION: }
+ { }
+ { Range reduction is into intervals of pi/4. The reduction }
+ { error is nearly eliminated by contriving an extended }
+ { precision modular arithmetic. }
+ { }
+ { Two polynomial approximating functions are employed. }
+ { Between 0 and pi/4 the cosine is approximated by }
+ { 1 - x**2 Q(x**2). }
+ { Between pi/4 and pi/2 the sine is represented as }
+ { x + x**3 P(x**2). }
+ {*****************************************************************}
+ var y, z, zz : Real;
+ j, sign : Integer;
+ i : LongInt;
+ begin
+ { make argument positive }
+ sign := 1;
+ if( d < 0 ) then
+ d := -d;
+
+ { above this value, round towards zero }
+ if( d > lossth ) then
+ begin
+ cos := 0.0;
+ exit;
+ end;
+
+ y := Trunc( d/PIO4 );
+ z := ldexp( y, -4 );
+ z := Trunc(z); { integer part of y/8 }
+ z := y - ldexp( z, 4 ); { y - 16 * (y/16) }
+
+ { integer and fractional part modulo one octant }
+ i := Trunc(z);
+ if odd( i ) then { map zeros to origin }
+ begin
+ inc(i);
+ y := y + 1.0;
+ end;
+ j := i and 07;
+ if( j > 3) then
+ begin
+ dec(j,4);
+ sign := -sign;
+ end;
+ if( j > 1 ) then
+ sign := -sign;
+
+ { Extended precision modular arithmetic }
+ z := ((d - y * DP1) - y * DP2) - y * DP3;
+
+ zz := z * z;
+
+ if( (j=1) or (j=2) ) then
+ { y = z + z * (zz * polevl( zz, sincof, 5 )); }
+ y := z + z * z * z * polevl( zz, sincof, 5 )
+ else
+ y := 1.0 - ldexp(zz,-1) + zz * zz * polevl( zz, coscof, 5 );
+
+ if(sign < 0) then
+ y := -y;
+
+ cos := y ;
+ end;
+
+
+
+ function ArcTan(d:Real):Real;
+ {*****************************************************************}
+ { Inverse circular tangent (arctangent) }
+ {*****************************************************************}
+ { }
+ { SYNOPSIS: }
+ { }
+ { double x, y, atan(); }
+ { }
+ { y = atan( x ); }
+ { }
+ { DESCRIPTION: }
+ { }
+ { Returns radian angle between -pi/2 and +pi/2 whose tangent }
+ { is x. }
+ { }
+ { Range reduction is from four intervals into the interval }
+ { from zero to tan( pi/8 ). The approximant uses a rational }
+ { function of degree 3/4 of the form x + x**3 P(x)/Q(x). }
+ {*****************************************************************}
+ const P : TabCoef = (
+ -8.40980878064499716001E-1,
+ -8.83860837023772394279E0,
+ -2.18476213081316705724E1,
+ -1.48307050340438946993E1, 0, 0, 0);
+ Q : TabCoef = (
+ 1.54974124675307267552E1,
+ 6.27906555762653017263E1,
+ 9.22381329856214406485E1,
+ 4.44921151021319438465E1, 0, 0, 0);
+
+ { tan( 3*pi/8 ) }
+ T3P8 = 2.41421356237309504880;
+ { tan( pi/8 ) }
+ TP8 = 0.41421356237309504880;
+
+ var y,z : Real;
+ Sign : Integer;
+
+ begin
+ { make argument positive and save the sign }
+ sign := 1;
+ if( d < 0.0 ) then
+ begin
+ sign := -1;
+ d := -d;
+ end;
+
+ { range reduction }
+ if( d > T3P8 ) then
+ begin
+ y := PIO2;
+ d := -( 1.0/d );
+ end
+ else if( d > TP8 ) then
+ begin
+ y := PIO4;
+ d := (d-1.0)/(d+1.0);
+ end
+ else
+ y := 0.0;
+
+ { rational form in x**2 }
+
+ z := d * d;
+ y := y + ( polevl( z, P, 3 ) / p1evl( z, Q, 4 ) ) * z * d + d;
+
+ if( sign < 0 ) then
+ y := -y;
+ Arctan := y;
+ end;
+
+ function frac(d : Real) : Real;
+ begin
+ frac := d - Int(d);
+ end;
+
+{$ifdef fixed}
+
+
+ function sqrt(d : fixed) : fixed;
+ begin
+ end;
+
+ function int(d : fixed) : fixed; assembler;
+ {*****************************************************************}
+ { Returns the integral part of d }
+ {*****************************************************************}
+ asm
+ move.l d,d0
+ and.l #$ffff0000,d0 { keep only upper bits .. }
+ end;
+
+
+ function trunc(d : fixed) : longint;
+ {*****************************************************************}
+ { Returns the Truncated integral part of d }
+ {*****************************************************************}
+ begin
+ trunc:=longint(integer(d shr 16)); { keep only upper 16 bits }
+ end;
+
+ function frac(d : fixed) : fixed; assembler;
+ {*****************************************************************}
+ { Returns the Fractional part of d }
+ {*****************************************************************}
+ asm
+ move.l d,d0
+ and.l #$ffff,d0 { keep only decimal parts - lower 16 bits }
+ end;
+
+ function abs(d : fixed) : fixed;
+ {*****************************************************************}
+ { Returns the Absolute value of d }
+ {*****************************************************************}
+ var
+ w: integer;
+ begin
+ w:=integer(d shr 16);
+ if w < 0 then
+ begin
+ w:=-w; { invert sign ... }
+ d:=d and $ffff;
+ d:=d or (fixed(w) shl 16); { add this to fixed number ... }
+ abs:=d;
+ end
+ else
+ abs:=d; { already positive... }
+ end;
+
+
+ function sqr(d : fixed) : fixed;
+ {*****************************************************************}
+ { Returns the Absolute squared value of d }
+ {*****************************************************************}
+ begin
+ {16-bit precision needed, not 32 =)}
+ sqr := d*d;
+{ sqr := (d SHR 8 * d) SHR 8; }
+ end;
+
+
+ function Round(x: fixed): longint;
+ {*****************************************************************}
+ { Returns the Rounded value of d as a longint }
+ {*****************************************************************}
+ var
+ lowf:integer;
+ highf:integer;
+ begin
+ lowf:=x and $ffff; { keep decimal part ... }
+ highf :=integer(x shr 16);
+ if lowf > 5 then
+ highf:=highf+1
+ else
+ if lowf = 5 then
+ begin
+ { here we must check the sign ... }
+ { if greater or equal to zero, then }
+ { greater value will be found by adding }
+ { one... }
+ if highf >= 0 then
+ Highf:=Highf+1;
+ end;
+ Round:= longint(highf);
+ end;
+{$endif fixed}
+
+ function power(bas,expo : real) : real;
+ begin
+ if bas=0.0 then
+ begin
+ if expo<>0.0 then
+ power:=0.0
+ else
+ HandleError(207);
+ end
+ else if expo=0.0 then
+ power:=1
+ else
+ { bas < 0 is not allowed }
+ if bas<0.0 then
+ handleerror(207)
+ else
+ power:=exp(ln(bas)*expo);
+ end;
+
+ function power(bas,expo : longint) : longint;
+ begin
+ if bas=0 then
+ begin
+ if expo<>0 then
+ power:=0
+ else
+ HandleError(207);
+ end
+ else if expo=0 then
+ power:=1
+ else
+ begin
+ if bas<0 then
+ begin
+ if odd(expo) then
+ power:=-round(exp(ln(-bas)*expo))
+ else
+ power:=round(exp(ln(-bas)*expo));
+ end
+ else
+ power:=round(exp(ln(bas)*expo));
+ end;
+ end;
+
+{
+ $Log: math.inc,v $
+ Revision 1.4 2005/02/14 17:13:30 peter
+ * truncate log
+
+}
diff --git a/rtl/m68k/readme b/rtl/m68k/readme
new file mode 100644
index 0000000000..0d769d90fd
--- /dev/null
+++ b/rtl/m68k/readme
@@ -0,0 +1,7 @@
+This directory contains only RTL parts specific
+to the processor M68K family.
+
+
+Missing units are :
+ strings.pp (you can the strings unit in the template directory)
+ getopts.pp
diff --git a/rtl/m68k/set.inc b/rtl/m68k/set.inc
new file mode 100644
index 0000000000..72275a016d
--- /dev/null
+++ b/rtl/m68k/set.inc
@@ -0,0 +1,428 @@
+{
+ $Id: set.inc,v 1.4 2005/02/14 17:13:30 peter Exp $
+ This file is part of the Free Pascal run time library.
+ Copyright (c) 1999-2000 by Carl-Eric Codere,
+ member of the Free Pascal development team.
+
+ See the file COPYING.FPC, included in this distribution,
+ for details about the copyright.
+
+ This program 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.
+
+ **********************************************************************}
+{*************************************************************************}
+{ Converted by Carl Eric Codere }
+{*************************************************************************}
+{ This inc. implements low-level set operations for the motorola }
+{ 68000 familiy of processors. }
+{ Based on original code bt Florian Kl„mpfl for the 80x86. }
+{*************************************************************************}
+
+
+ { add the element b to the set pointed by p }
+ { On entry }
+ { a0 = pointer to set }
+ { d0.b = element to add to the set }
+ { Registers destroyed: d0,a1,d6 }
+ procedure do_set;assembler;
+ asm
+ XDEF SET_SET_BYTE
+ move.l d0,d6
+ { correct long position: }
+ { -> (value div 32)*4 = longint }
+ { (value shr 5)*shl 2 }
+ lsr.l #5,d6
+ lsl.l #2,d6
+ adda.l d6,a0 { correct offset from start address of set }
+
+ move.l d0,d6 { bit is now in here }
+ andi.l #31,d0 { bit number is = value mod 32 }
+
+ { now bit set the value }
+ move.l (a0),d0 { we must put bits into register }
+ bset.l d6,d0 { otherwise btst will be a byte }
+ { put result in carry flag } { operation. }
+ bne @LDOSET1
+ andi.b #$fe,ccr { clear carry flag }
+ bra @LDOSET2
+ @LDOSET1:
+ ori.b #$01,ccr { set carry flag }
+ @LDOSET2:
+ move.l d0,(a0) { restore the value at that location }
+ { of the set. }
+ end;
+
+ { Finds an element in a set }
+ { a0 = address of set }
+ { d0.b = value to compare with }
+ { CARRY SET IF FOUND ON EXIT }
+ { Registers destroyed: d0,a0,d6 }
+ procedure do_in; assembler;
+ { Returns Carry set then = in set , otherwise carry is cleared }
+ { (D0) }
+ asm
+ XDEF SET_IN_BYTE
+ move.l d0,d6
+ { correct long position: }
+ { -> (value div 32)*4 = longint }
+ { (value shr 5)*shl 2 }
+ lsr.l #5,d6
+ lsl.l #2,d6
+ adda.l d6,a0 { correct offset from start address of set }
+
+ move.l d0,d6 { bit is now in here }
+ andi.l #31,d0 { bit number is = value mod 32 }
+
+ move.l (a0),d0 { we must put bits into register }
+ btst.l d6,d0 { otherwise btst will be a byte }
+ { put result in carry flag } { operation. }
+ bne @LDOIN1
+ andi.b #$fe,ccr { clear carry flag }
+ bra @LDOIN2
+ @LDOIN1:
+ ori.b #$01,ccr { set carry flag }
+ @LDOIN2:
+ end;
+
+
+
+ { vereinigt set1 und set2 und speichert das Ergebnis in dest }
+
+ procedure add_sets(set1,set2,dest : pointer);[public,alias: 'SET_ADD_SETS'];
+ { PSEUDO-CODE:
+ type
+ destination = array[1..8] of longint;
+ for i:=1 to 8 do
+ destination(dest^)[i] := destination(set1^)[i] OR destination(set2^)[i];
+ }
+ begin
+ asm
+ { saved used register }
+ move.l a2,-(sp)
+
+ move.l 8(a6),a0
+ move.l 12(a6),a1
+ move.l 16(a6),a2
+
+ move.l #32,d6
+
+ @LMADDSETS1:
+
+ move.b (a0)+,d0
+ or.b (a1)+,d0
+ move.b d0,(a2)+
+ subq.b #1,d6
+ bne @LMADDSETS1
+ { restore register }
+ move.l a2,(sp)+
+ end ['d0','d6','a0','a1'];
+ end;
+
+ { computes the symetric diff from set1 to set2 }
+ { result in dest }
+
+ procedure sym_sub_sets(set1,set2,dest : pointer);[public,alias: 'SET_SYMDIF_SETS'];
+
+ begin
+ asm
+ { saved used register }
+ move.l a2,-(sp)
+
+ move.l 8(a6),a0
+ move.l 12(a6),a1
+ move.l 16(a6),a2
+
+ move.l #32,d6
+
+ @LMADDSETS1:
+
+ move.b (a0)+,d0
+ move.b (a1)+,d1
+ eor.b d1,d0
+ move.b d0,(a2)+
+ subq.b #1,d6
+ bne @LMADDSETS1
+ { restore register }
+ move.l a2,(sp)+
+ end;
+ end;
+
+
+ { bad implementation, but it's very seldom used }
+ procedure do_set(p : pointer;l,h : byte);[public,alias: 'SET_SET_RANGE'];
+
+ begin
+ asm
+ move.b h,d0
+ @LSetRLoop:
+ cmp.b l,d0
+ blt @Lend
+ move.w d0,-(sp)
+ { adjust value to correct endian }
+ lsl.w #8,d0
+ pea p
+ jsr SET_SET_BYTE
+ sub.b #1,d0
+ bra @LSetRLoop
+ @Lend:
+ end;
+ end;
+
+
+ { bildet den Durchschnitt von set1 und set2 }
+ { und speichert das Ergebnis in dest }
+
+ procedure mul_sets(set1,set2,dest : pointer);[public,alias: 'SET_MUL_SETS'];
+ { type
+ larray = array[0..7] of longint;
+ for i:=0 to 7 do
+ larray(dest^)[i] := larray(set1^)[i] AND larray(set2^)[i];
+ }
+ begin
+ asm
+ { saved used register }
+ move.l a2,-(sp)
+ move.l 8(a6),a0
+ move.l 12(a6),a1
+ move.l 16(a6),a2
+
+ move.l #32,d6
+
+ @LMMULSETS1:
+
+ move.b (a0)+,d0
+ and.b (a1)+,d0
+ move.b d0,(a2)+
+ subq.b #1,d6
+ bne @LMMULSETS1
+ { restore register }
+ move.l a2,(sp)+
+ end ['d0','d6','a0','a1'];
+ end;
+
+
+ { bildet die Differenz von set1 und set2 }
+ { und speichert das Ergebnis in dest }
+
+ procedure sub_sets(set1,set2,dest : pointer);[public,alias: 'SET_SUB_SETS'];
+ { type
+ larray = array[0..7] of longint;
+ begin
+ for i:=0 to 7 do
+ larray(dest^)[i] := larray(set1^)[i] AND NOT (larray(set2^)[i]);
+ end;
+ }
+ begin
+ asm
+ { saved used register }
+ move.l a2,-(sp)
+ move.l 8(a6),a0
+ move.l 12(a6),a1
+ move.l 16(a6),a2
+
+ move.l #32,d6
+
+ @LSUBSETS1:
+ move.b (a0)+,d0
+ move.b (a1)+,d1
+ not.b d1
+ and.b d1,d0
+ move.b d0,(a2)+
+ sub.b #1,d6
+ bne @LSUBSETS1
+ { restore register }
+ move.l a2,(sp)+
+ end ['d0','d1','d6','a0','a1'];
+ end;
+
+ { compare both sets }
+ { compares set1 and set2 }
+ { zeroflag is set if they are equal }
+ { on entry : a0 = pointer to first set }
+ { : a1 = pointer to second set }
+ procedure comp_sets; assembler;
+
+ asm
+ XDEF SET_COMP_SETS
+ move.l #32,d6
+ @LMCOMPSETS1:
+ move.b (a0)+,d0
+ move.b (a1),d1
+ cmp.b d1,d0
+ bne @LMCOMPSETEND
+ adda.l #1,a1
+ sub.b #1,d6
+ bne @LMCOMPSETS1
+ { we are here only if the two sets are equal }
+ { we have zero flag set, and that what is expected }
+ cmp.b d0,d0
+ @LMCOMPSETEND:
+ end;
+
+ procedure do_set(p : pointer;b : word);[public,alias: 'SET_SET_WORD'];
+ begin
+ asm
+ move.l 8(a6),a0
+ move.w 12(a6),d6
+ andi.l #$fff8,d6
+ lsl.l #3,d6
+ adda.l d6,a0
+ move.b 12(a6),d6
+ andi.l #7,d6
+
+ move.l (a0),d0 { we must put bits into register }
+ btst.l d6,d0 { otherwise btst will be a byte }
+ { put result in carry flag } { operation. }
+ bne @LBIGDOSET1
+ andi.b #$fe,ccr { clear carry flag }
+ bra @LBIGDOSET2
+ @LBIGDOSET1:
+ ori.b #$01,ccr { set carry flag }
+ @LBIGDOSET2:
+ end ['d0','a0','d6'];
+ end;
+
+ { testet, ob das Element b in der Menge p vorhanden ist }
+ { und setzt das Carryflag entsprechend }
+
+ procedure do_in(p : pointer;b : word);[public,alias: 'SET_IN_WORD'];
+ begin
+ asm
+ move.l 8(a6),a0
+ move.w 12(a6),d6
+ andi.l #$fff8,d6
+ lsl.l #3,d6
+ adda.l d6,a0 { correct offset from start address of set }
+
+ move.b 12(a6),d6
+ andi.l #7,d6
+
+ move.l (a0),d0 { we must put bits into register }
+ btst.l d6,d0 { otherwise btst will be a byte }
+ { put result in carry flag } { operation. }
+ bne @LBIGDOIN1
+ andi.b #$fe,ccr { clear carry flag }
+ bra @LBIGDOIN2
+ @LBIGDOIN1:
+ ori.b #$01,ccr { set carry flag }
+ @LBIGDOIN2:
+ end ['d0','a0','d6'];
+ end;
+
+
+ { vereinigt set1 und set2 und speichert das Ergebnis in dest }
+ { size is the number of bytes in the set }
+
+ procedure add_sets(set1,set2,dest : pointer;size : longint);[public,alias: 'SET_ADD_SETS_SIZE'];
+ begin
+ asm
+ { saved used register }
+ move.l a2,-(sp)
+ move.l 8(a6),a0
+ move.l 12(a6),a1
+ move.l 16(a6),a2
+
+ move.l 20(a6),d6
+
+ @LBIGMADDSETS1:
+
+ move.l (a0)+,d0
+ or.l (a1)+,d0
+ move.l d0,(a2)+
+ subq.l #4,d6
+ bne @LBIGMADDSETS1
+ { restore register }
+ move.l a2,(sp)+
+ end ['d0','d6','a0','a1'];
+ end;
+
+
+ procedure mul_sets(set1,set2,dest : pointer;size : longint);[public,alias: 'SET_MUL_SETS_SIZE'];
+ { bildet den Durchschnitt von set1 und set2 }
+ { und speichert das Ergebnis in dest }
+ { size is the number of bytes in the set }
+ begin
+ asm
+ { saved used register }
+ move.l a2,-(sp)
+ move.l 8(a6),a0
+ move.l 12(a6),a1
+ move.l 16(a6),a2
+
+ move.l 20(a6),d6
+
+ @LBIGMMULSETS1:
+
+ move.l (a0)+,d0
+ and.l (a1)+,d0
+ move.l d0,(a2)+
+ subq.l #4,d6
+ bne @LBIGMMULSETS1
+ { restore register }
+ move.l a2,(sp)+
+ end ['d0','d6','a0','a1'];
+ end;
+
+
+ { bildet die Differenz von set1 und set2 }
+ { und speichert das Ergebnis in dest }
+ { size is the number of bytes in the set }
+
+ procedure sub_sets(set1,set2,dest : pointer;size : longint);[public,alias: 'SET_SUB_SETS_SIZE'];
+ begin
+ asm
+ { saved used register }
+ move.l a2,-(sp)
+ move.l 8(a6),a0
+ move.l 12(a6),a1
+ move.l 16(a6),a2
+
+ move.l 20(a6),d6
+
+ @BIGSUBSETS1:
+
+ move.l (a0)+,d0
+ not.l d0
+ and.l (a1)+,d0
+ move.l d0,(a2)+
+ subq.l #4,d6
+ bne @BIGSUBSETS1
+ { restore register }
+ move.l a2,(sp)+
+ end ['d0','d6','a0','a1'];
+ end;
+
+
+ { vergleicht Mengen und setzt die Flags entsprechend }
+
+ procedure comp_sets(set1,set2 : pointer;size : longint);[public,alias: 'SET_COMP_SETS_SIZE'];
+
+
+ begin
+ asm
+ move.l 8(a6),a0 { set1 - esi}
+ move.l 12(a6),a1 { set2 - edi }
+ move.l 16(a6),d6
+ @MCOMPSETS1:
+ move.l (a0)+,d0
+ move.l (a1),d1
+ cmp.l d1,d0
+ bne @BIGMCOMPSETEND
+ add.l #4,a1
+ subq.l #1,d6
+ bne @MCOMPSETS1
+ { we are here only if the two sets are equal }
+ { we have zero flag set, and that what is expected }
+ cmp.l d0,d0
+ @BIGMCOMPSETEND:
+ end;
+ end;
+
+{
+ $Log: set.inc,v $
+ Revision 1.4 2005/02/14 17:13:30 peter
+ * truncate log
+
+}
diff --git a/rtl/m68k/setjump.inc b/rtl/m68k/setjump.inc
new file mode 100644
index 0000000000..baf7a01e9e
--- /dev/null
+++ b/rtl/m68k/setjump.inc
@@ -0,0 +1,34 @@
+{
+ $Id: setjump.inc,v 1.4 2005/02/14 17:13:30 peter Exp $
+ This file is part of the Free Pascal run time library.
+ Copyright (c) 1999-2000 by xxxx
+ member of the Free Pascal development team
+
+ See the file COPYING.FPC, included in this distribution,
+ for details about the copyright.
+
+ This program 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.
+
+ **********************************************************************}
+
+{**********************************************************************
+ Set_Jmp/Long_jmp
+ **********************************************************************}
+
+Function SetJmp (Var S : Jmp_buf) : longint;assembler;[Public, alias : 'FPC_SETJMP'];
+
+asm
+end;
+
+Procedure longJmp (Var S : Jmp_buf; value : longint); assembler;[Public, alias : 'FPC_LONGJMP'];
+
+asm
+end;
+
+ $Log: setjump.inc,v $
+ Revision 1.4 2005/02/14 17:13:30 peter
+ * truncate log
+
+}
diff --git a/rtl/m68k/setjumph.inc b/rtl/m68k/setjumph.inc
new file mode 100644
index 0000000000..1b6884e927
--- /dev/null
+++ b/rtl/m68k/setjumph.inc
@@ -0,0 +1,37 @@
+{
+ $Id: setjumph.inc,v 1.6 2005/02/14 17:13:30 peter Exp $
+ This file is part of the Free Pascal run time library.
+ Copyright (c) 1999-2000 by xxxx
+ member of the Free Pascal development team
+
+ See the file COPYING.FPC, included in this distribution,
+ for details about the copyright.
+
+ This program 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.
+
+ **********************************************************************}
+
+{**********************************************************************
+ Declarations for SetJmp/LongJmp
+ **********************************************************************}
+
+Type
+ jmp_buf = packed record
+ fp : longint; { frame pointer }
+ sp : longint; { stack pointer }
+ pc : longint; { program counter }
+ aregs : array[0..3] of dword; { address registers (a2,a3,a4,a5) }
+ end;
+ PJmp_buf = ^jmp_buf;
+
+Function Setjmp (Var S : Jmp_buf) : longint;
+Procedure longjmp (Var S : Jmp_buf; value : longint);
+
+{
+ $Log: setjumph.inc,v $
+ Revision 1.6 2005/02/14 17:13:30 peter
+ * truncate log
+
+}
View Lorry Controller Status