← Back to branch summary

~ubuntu-branches/debian/lenny/fpc/lenny

« back to all changes in this revision

Viewing changes to rtl/m68k/lowmath.inc

  • Committer: Bazaar Package Importer
  • Author(s): Mazen Neifer, Torsten Werner, Mazen Neifer
  • Date: 2008-05-17 17:12:11 UTC
  • mfrom: (3.1.9 intrepid)
  • Revision ID: james.westby@ubuntu.com-20080517171211-9qi33xhd9evfa0kg
Tags: 2.2.0-dfsg1-9
http://bugs.debian.org/477743
[ Torsten Werner ]
* Add Mazen Neifer to Uploaders field.

[ Mazen Neifer ]
* Moved FPC sources into a version dependent directory from /usr/share/fpcsrc
  to /usr/share/fpcsrc/${FPCVERSION}. This allow installing more than on FPC
  release.
* Fixed far call issue in compiler preventing building huge binearies.
  (closes: #477743)
* Updated building dependencies, recomennded and suggested packages.
* Moved fppkg to fp-utils as it is just a helper tool and is not required by
  compiler.

Show diffs side-by-side

added added

removed removed

Lines of Context:
rtl/m68k/lowmath.inc
1
 
{
2
 
    $Id: lowmath.inc,v 1.4 2002/09/07 16:01:20 peter Exp $
3
 
    This file is part of the Free Pascal run time library.
4
 
    Copyright (c) 1999-2000 by Carl-Eric Codere,
5
 
    member of the Free Pascal development team
6
 
 
7
 
    See the file COPYING.FPC, included in this distribution,
8
 
    for details about the copyright.
9
 
 
10
 
    This program is distributed in the hope that it will be useful,
11
 
    but WITHOUT ANY WARRANTY; without even the implied warranty of
12
 
    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.
13
 
 
14
 
 **********************************************************************}
15
 
{*************************************************************************}
16
 
{  lowmath.inc                                                            }
17
 
{  Ported to FPC-Pascal by Carl Eric Codere                               }
18
 
{  Terms of use: This source code is freeware.                            }
19
 
{*************************************************************************}
20
 
{  This inc. implements low-level mathemtical routines  for the motorola  }
21
 
{  68000 family of processors.                                            }
22
 
{*************************************************************************}
23
 
{ single floating point routines taken from GCC 2.5.2 Atari compiler      }
24
 
{ library source.                                                         }
25
 
{  Original credits:                                                      }
26
 
{    written by Kai-Uwe Bloem (I5110401@dbstu1.bitnet).                   }
27
 
{      Based on a 80x86 floating point packet from comp.os.minix,         }
28
 
{        written by P.Housel                                              }
29
 
{    Patched by Olaf Flebbe (flebbe@tat.physik.uni-tuebingen.de)          }
30
 
{    Revision by michal 05-93 (ntomczak@vm.ucs.ualberta.ca)               }
31
 
{*************************************************************************}
32
 
{--------------------------------------------------------------------}
33
 
{ LEFT TO DO:                                                        }
34
 
{ o Add support for FPU if present.                                  }
35
 
{ o Verify if single comparison works in all cases.                  }
36
 
{ o Add support for  NANs in SINGLE_CMP                              }
37
 
{ o Add comp (80-bit) multiplication,addition,substract,division,    }
38
 
{    shift.                                                          }
39
 
{ o Add stack checking for the routines which use the stack.         }
40
 
{    (This will probably have to be done in the code generator).     }
41
 
{--------------------------------------------------------------------}
42
 
 
43
 
 
44
 
 
45
 
Procedure Single_Norm;[alias : 'FPC_SINGLE_NORM'];Assembler;
46
 
{--------------------------------------------}
47
 
{ Low-level routine to normalize single  e   }
48
 
{  IEEE floating point values. Never called  }
49
 
{  directly.                                 }
50
 
{ On Exit:                                   }
51
 
{      d0 = result.                          }
52
 
{  Registers destroyed: d0,d1                }
53
 
{--------------------------------------------}
54
 
Asm
55
 
   tst.l    d4             { rounding and u.mant == 0 ?    }
56
 
   bne      @normlab1
57
 
   tst.b    d1
58
 
   beq      @retzok
59
 
@normlab1:
60
 
   clr.b    d2             { "sticky byte"                 }
61
 
@normlab3:
62
 
   move.l   #$ff000000,d5
63
 
@normlab4:
64
 
   tst.w    d0             { divide (shift)                }
65
 
   ble      @normlab2      {  denormalized number          }
66
 
   move.l   d4,d3
67
 
   and.l    d5,d3          {  or until no bits above 23    }
68
 
   beq      @normlab5
69
 
@normlab2:
70
 
   addq.w   #1,d0          { increment exponent            }
71
 
   lsr.l    #1,d4
72
 
   or.b     d1,d2          { set "sticky"                  }
73
 
   roxr.b   #1,d1          { shift into rounding bits      }
74
 
   bra      @normlab4
75
 
@normlab5:
76
 
   and.b    #1,d2
77
 
   or.b     d2,d1          { make least sig bit "sticky"   }
78
 
   asr.l    #1,d5          { #0xff800000 -> d5             }
79
 
@normlab6:
80
 
   move.l   d4,d3          { multiply (shift) until        }
81
 
   and.l    d5,d3          { one in "implied" position     }
82
 
   bne      @normlab7
83
 
   subq.w   #1,d0          { decrement exponent            }
84
 
   beq      @normlab7      {  too small. store as denormalized number }
85
 
   add.b    d1,d1          { some doubt about this one *              }
86
 
   addx.l   d4,d4
87
 
   bra      @normlab6
88
 
@normlab7:
89
 
   tst.b    d1             { check rounding bits          }
90
 
   bge      @normlab9      { round down - no action neccessary }
91
 
   neg.b    d1
92
 
   bvc      @normlab8      { round up                     }
93
 
   move.w   d4,d1          { tie case - round to even     }
94
 
                           {   dont need rounding bits any more }
95
 
   and.w    #1,d1          { check if even                }
96
 
   beq      @normlab9      {   mantissa is even - no action necessary  }
97
 
                           { fall through                 }
98
 
@normlab8:
99
 
   clr.w    d1             { zero rounding bits            }
100
 
   add.l    #1,d4
101
 
   tst.w    d0
102
 
   bne      @normlab10     { renormalize if number was denormalized   }
103
 
   add.w    #1,d0          { correct exponent for denormalized numbers }
104
 
   bra      @normlab3
105
 
@normlab10:
106
 
   move.l   d4,d3          { check for rounding overflow              }
107
 
   asl.l    #1,d5          { #0xff000000 -> d5                        }
108
 
   and.l    d5,d3
109
 
   bne      @normlab4      { go back and renormalize                  }
110
 
@normlab9:
111
 
   tst.l    d4             { check if normalization caused an underflow }
112
 
   beq      @retz
113
 
   tst.w    d0             { check for exponent overflow or underflow  }
114
 
   blt      @retz
115
 
   cmp.w    #255,d0
116
 
   bge      @oflow
117
 
 
118
 
   lsl.w    #8,d0          { re-position exponent - one bit too high }
119
 
   lsl.w    #1,d2          { get X bit                               }
120
 
   roxr.w   #1,d0          { shift it into sign position             }
121
 
   swap     d0             { map to upper word                       }
122
 
   clr.w    d0
123
 
   and.l    #$7fffff,d4    { top mantissa bits                       }
124
 
   or.l     d4,d0          { insert exponent and sign                }
125
 
   movem.l  (sp)+,d2-d5
126
 
   rts
127
 
 
128
 
@retz:
129
 
   { handling underflow should be done here... }
130
 
   { by default simply return 0 as retzok...   }
131
 
@retzok:
132
 
   moveq.l   #0,d0
133
 
   lsl.w     #1,d2
134
 
   roxr.l    #1,d0          { sign of 0 is the same as of d2          }
135
 
   movem.l   (sp)+,d2-d5
136
 
   rts
137
 
 
138
 
@oflow:
139
 
   move.l    #$7f800000,d0  { +infinity as proposed by IEEE         }
140
 
 
141
 
   tst.w     d2             { transfer sign                         }
142
 
   bge       @ofl_clear     { (mjr++)                               }
143
 
   bset      #31,d0         {                                       }
144
 
@ofl_clear:
145
 
   or.b      #2,ccr         { set overflow flag.                    }
146
 
   movem.l   (sp)+,d2-d5
147
 
   rts
148
 
end;
149
 
 
150
 
 
151
 
Procedure Single_AddSub; Assembler;
152
 
{--------------------------------------------}
153
 
{ Low-level routine to add/subtract single   }
154
 
{  IEEE floating point values. Never called  }
155
 
{  directly.                                 }
156
 
{ On Exit:                                   }
157
 
{      d0 = result -- from normalize routine }
158
 
{      Flags : V set if overflow.            }
159
 
{              on underflow d0 = 0           }
160
 
{  Registers destroyed: d0,d1                }
161
 
{--------------------------------------------}
162
 
Asm
163
 
{--------------------------------------------}
164
 
{ On Entry:                                  }
165
 
{      d1-d0 = single values to subtract.    }
166
 
{--------------------------------------------}
167
 
XDEF SINGLE_SUB
168
 
   eor.l      #$80000000,d0         {   reverse sign of v }
169
 
{--------------------------------------------}
170
 
{ On Entry:                                  }
171
 
{      d0, d1 = single values to add.        }
172
 
{--------------------------------------------}
173
 
XDEF SINGLE_ADD
174
 
   movem.l d2-d5,-(sp)              { save registers      }
175
 
   move.l   d0,d4                   { d4 = d0 = v         }
176
 
   move.l   d1,d5                   { d5 = d1 = u         }
177
 
 
178
 
   move.l  #$7fffff,d3
179
 
   move.l  d5,d0                    { d0 = u.exp          }
180
 
   move.l  d5,d2                    { d2.h = u.sign       }
181
 
   swap       d0
182
 
   move.w  d0,d2                    { d2 = u.sign         }
183
 
   and.l   d3,d5                    { remove exponent from u.mantissa }
184
 
 
185
 
   move.l  d4,d1                    { d1 = v.exp          }
186
 
   and.l   d3,d4                    { remove exponent from v.mantissa }
187
 
   swap    d1
188
 
   eor.w   d1,d2                    { d2 = u.sign ^ v.sign (in bit 15)}
189
 
   clr.b   d2                       { we will use the lowest byte as a flag }
190
 
   moveq.l #15,d3
191
 
   bclr    d3,d1                    { kill sign bit u.exp             }
192
 
   bclr    d3,d0                    { kill sign bit u.exp             }
193
 
   btst    d3,d2                    { same sign for u and v?          }
194
 
   beq     @slabel1
195
 
   cmp.l   d0,d1                    { different signs - maybe x - x ? }
196
 
   seq     d2                       { set 'cancellation' flag         }
197
 
@slabel1:
198
 
   lsr.w   #7,d0                    { keep here exponents only        }
199
 
   lsr.w   #7,d1
200
 
{--------------------------------------------------------------------}
201
 
{          Now perform testing of NaN and infinities                 }
202
 
{--------------------------------------------------------------------}
203
 
   moveq.l #-1,d3
204
 
   cmp.b   d3,d0
205
 
   beq     @alabel1
206
 
   cmp.b   d3,d1
207
 
   bne     @nospec
208
 
   bra     @alabel2
209
 
{--------------------------------------------------------------------}
210
 
{                        u is special.                               }
211
 
{--------------------------------------------------------------------}
212
 
@alabel1:
213
 
   tst.b      d2
214
 
   bne        @retnan      { cancellation of specials -> NaN }
215
 
   tst.l      d5
216
 
   bne        @retnan      { arith with Nan gives always NaN }
217
 
 
218
 
   addq.w     #4,a0        { here is an infinity             }
219
 
   cmp.b      d3,d1
220
 
   bne        @alabel3     { skip check for NaN if v not special }
221
 
{--------------------------------------------------------------------}
222
 
{                        v is special.                               }
223
 
{--------------------------------------------------------------------}
224
 
@alabel2:
225
 
   tst.l           d4
226
 
   bne        @retnan
227
 
@alabel3:
228
 
   move.l     (a0),d0
229
 
   bra        @return
230
 
{--------------------------------------------------------------------}
231
 
{                     Return a quiet nan                             }
232
 
{--------------------------------------------------------------------}
233
 
@retnan:
234
 
   moveq.l    #-1,d0
235
 
   lsr.l      #1,d0                { 0x7fffffff -> d0        }
236
 
   bra        @return
237
 
{ Ok, no inifinty or NaN involved.. }
238
 
@nospec:
239
 
   tst.b           d2
240
 
   beq        @alabel4
241
 
   moveq.l    #0,d0                { x - x hence we always return +0 }
242
 
@return:
243
 
   movem.l    (sp)+,d2-d5
244
 
   rts
245
 
 
246
 
@alabel4:
247
 
   moveq.l    #23,d3
248
 
   bset       d3,d5           { restore implied leading "1"   }
249
 
   tst.w      d0              { check for zero exponent - no leading "1" }
250
 
   bne        @alabel5
251
 
   bclr       d3,d5           { remove it                     }
252
 
   addq.w     #1,d0           { "normalize" exponent          }
253
 
@alabel5:
254
 
   bset       d3,d4           { restore implied leading "1"   }
255
 
   tst.w      d1              { check for zero exponent - no leading "1"  }
256
 
   bne        @alabel6
257
 
   bclr       d3,d4           { remove it                     }
258
 
   addq.w  #1,d1              { "normalize" exponent          }
259
 
@alabel6:
260
 
   moveq.l    #0,d3           { (put initial zero rounding bits in d3) }
261
 
   neg.w      d1              { d1 = u.exp - v.exp            }
262
 
   add.w      d0,d1
263
 
   beq        @alabel8      { exponents are equal - no shifting neccessary }
264
 
   bgt        @alabel7      { not equal but no exchange neccessary         }
265
 
   exg        d4,d5                { exchange u and v }
266
 
   sub.w      d1,d0                { d0 = u.exp - (u.exp - v.exp) = v.exp }
267
 
   neg.w      d1
268
 
   tst.w      d2              { d2.h = u.sign ^ (u.sign ^ v.sign) = v.sign }
269
 
   bpl        @alabel7
270
 
   bchg       #31,d2
271
 
@alabel7:
272
 
   cmp.w      #26,d1       { is u so much bigger that v is not }
273
 
   bge        @alabel9      { significant ?                     }
274
 
{--------------------------------------------------------------------}
275
 
{ shift mantissa left two digits, to allow cancellation of           }
276
 
{ most significant digit, while gaining an additional digit for      }
277
 
{ rounding.                                                          }
278
 
{--------------------------------------------------------------------}
279
 
   moveq.l #1,d3
280
 
@alabel10:
281
 
   add.l           d5,d5
282
 
   subq.w  #1,d0              { decrement exponent            }
283
 
   subq.w  #1,d1              { done shifting altogether ?    }
284
 
   dbeq    d3,@alabel10        { loop if still can shift u.mant more }
285
 
   moveq.l    #0,d3
286
 
 
287
 
   cmp.w      #16,d1          { see if fast rotate possible    }
288
 
   blt        @alabel11
289
 
   or.w       d4,d3           { set rounding bits              }
290
 
   clr.w      d4
291
 
   swap       d4
292
 
   subq.w  #8,d1
293
 
   subq.w  #8,d1
294
 
   bra      @alabel11
295
 
 
296
 
@alabel12:
297
 
   move.b    d4,d2
298
 
   and.b      #1,d2
299
 
   or.b       d2,d3
300
 
   lsr.l      #1,d4               { shift v.mant right the rest of the way }
301
 
@alabel11:
302
 
   dbra    d1,@alabel12           { loop                                   }
303
 
 
304
 
@alabel8:
305
 
   tst.w      d2                  { are the signs equal ?                  }
306
 
   bpl        @alabel13           { yes, no negate necessary               }
307
 
 
308
 
 
309
 
   tst.w      d3                  { negate rounding bits and v.mant        }
310
 
   beq        @alabel14
311
 
   addq.l  #1,d4
312
 
@alabel14:
313
 
   neg.l           d4
314
 
 
315
 
@alabel13:
316
 
   add.l      d4,d5                { u.mant = u.mant + v.mant              }
317
 
   bcs        @alabel9             { needn not negate                      }
318
 
   tst.w      d2                   { opposite signs ?                      }
319
 
   bpl        @alabel9             { do not need to negate result          }
320
 
 
321
 
   neg.l      d5
322
 
   not.l      d2                   { switch sign                           }
323
 
@alabel9:
324
 
   move.l  d5,d4                   { move result for normalization         }
325
 
   clr.l      d1
326
 
   tst.l      d3
327
 
   beq        @alabel15
328
 
   moveq.l #-1,d1
329
 
@alabel15:
330
 
   swap    d2                      { put sign into d2 (exponent is in d0)  }
331
 
   jmp     FPC_SINGLE_NORM             { leave registers on stack for norm_sf  }
332
 
end;
333
 
 
334
 
 
335
 
Procedure Single_Mul;Assembler;
336
 
{--------------------------------------------}
337
 
{ Low-level routine to multiply two single   }
338
 
{  IEEE floating point values. Never called  }
339
 
{  directly.                                 }
340
 
{ Om Entry:                                  }
341
 
{      d0,d1 = values to multiply            }
342
 
{ On Exit:                                   }
343
 
{      d0 = result.                          }
344
 
{  Registers destroyed: d0,d1                }
345
 
{ stack space used (and restored): 8 bytes.  }
346
 
{--------------------------------------------}
347
 
Asm
348
 
XDEF SINGLE_MUL
349
 
   movem.l  d2-d5,-(sp)
350
 
   move.l   d0,d4       { d4 = v                   }
351
 
   move.l   d1,d5       { d5 = u                   }
352
 
 
353
 
   move.l   #$7fffff,d3
354
 
   move.l   d5,d0       { d0 = u.exp               }
355
 
   and.l    d3,d5       { remove exponent from u.mantissa }
356
 
   swap     d0
357
 
   move.w   d0,d2       { d2 = u.sign              }
358
 
 
359
 
   move.l   d4,d1       { d1 = v.exp               }
360
 
   and.l    d3,d4       { remove exponent from v.mantissa }
361
 
   swap     d1
362
 
   eor.w    d1,d2       { d2 = u.sign ^ v.sign (in bit 15)}
363
 
 
364
 
   moveq.l  #15,d3
365
 
   bclr     d3,d0       { kill sign bit                   }
366
 
   bclr     d3,d1       { kill sign bit                   }
367
 
   tst.l    d0          { test if one of factors is 0     }
368
 
   beq      @mlabel1
369
 
   tst.l    d1
370
 
@mlabel1:
371
 
   seq      d2         { 'one of factors is 0' flag in the lowest byte }
372
 
   lsr.w    #7,d0      { keep here exponents only         }
373
 
   lsr.w    #7,d1
374
 
 
375
 
{--------------------------------------------------------------------}
376
 
{          Now perform testing of NaN and infinities                 }
377
 
{--------------------------------------------------------------------}
378
 
   moveq.l  #-1,d3
379
 
   cmp.b    d3,d0
380
 
   beq      @mlabel2
381
 
   cmp.b    d3,d1
382
 
   bne      @mnospec
383
 
   bra      @mlabel3
384
 
{--------------------------------------------------------------------}
385
 
{                   first operand is special                         }
386
 
{--------------------------------------------------------------------}
387
 
@mlabel2:
388
 
   tst.l    d5         { is it NaN?                       }
389
 
   bne      @mretnan
390
 
@mlabel3:
391
 
   tst.b    d2         { 0 times special or special times 0 ? }
392
 
   bne      @mretnan   { yes -> NaN                       }
393
 
   cmp.b    d3,d1      { is the other special ?           }
394
 
   beq      @mlabel4   { maybe it is NaN                  }
395
 
{--------------------------------------------------------------------}
396
 
{                   Return infiny with correct sign                  }
397
 
{--------------------------------------------------------------------}
398
 
@mretinf:
399
 
   move.l   #$ff000000,d0  { we will return #0xff800000 or #0x7f800000 }
400
 
   lsl.w    #1,d2
401
 
   roxr.l   #1,d0      { shift in high bit as given by d2 }
402
 
@mreturn:
403
 
   movem.l  (sp)+,d2-d5
404
 
   rts
405
 
 
406
 
{--------------------------------------------------------------------}
407
 
{                        v is special.                               }
408
 
{--------------------------------------------------------------------}
409
 
@mlabel4:
410
 
   tst.l    d4          { is this NaN?                    }
411
 
   beq      @mretinf    { we know that the other is not zero }
412
 
@mretnan:
413
 
   moveq.l  #-1,d0
414
 
   lsr.l    #1,d0       { 0x7fffffff -> d0                }
415
 
   bra      @mreturn
416
 
{--------------------------------------------------------------------}
417
 
{                    End of NaN and Inf                              }
418
 
{--------------------------------------------------------------------}
419
 
@mnospec:
420
 
   tst.b    d2          { not needed - but we can waste two instr. }
421
 
   bne      @mretzz     { return signed 0 if one of factors is 0   }
422
 
   moveq.l  #23,d3
423
 
   bset     d3,d5       { restore implied leading "1"              }
424
 
   subq.w   #8,sp       { multiplication accumulator               }
425
 
   tst.w    d0          { check for zero exponent - no leading "1" }
426
 
   bne      @mlabel5
427
 
   bclr     d3,d5       { remove it                                }
428
 
   addq.w   #1,d0       { "normalize" exponent                     }
429
 
@mlabel5:
430
 
   tst.l   d5
431
 
   beq     @mretz       { multiplying zero                         }
432
 
 
433
 
   moveq.l #23,d3
434
 
   bset    d3,d4        { restore implied leading "1"              }
435
 
   tst.w   d1           { check for zero exponent - no leading "1" }
436
 
   bne     @mlabel6
437
 
   bclr    d3,d4        { remove it                                }
438
 
   addq.w  #1,d1        { "normalize" exponent                     }
439
 
@mlabel6:
440
 
   tst.l   d4
441
 
   beq     @mretz       { multiply by zero                         }
442
 
 
443
 
   add.w   d1,d0        { add exponents,                           }
444
 
   sub.w   #BIAS4+16-8,d0 { remove excess bias, acnt for repositioning }
445
 
 
446
 
   clr.l   (sp)         { initialize 64-bit product to zero        }
447
 
   clr.l   4(sp)
448
 
{--------------------------------------------------------------------}
449
 
{ see Knuth, Seminumerical Algorithms, section 4.3. algorithm M      }
450
 
{--------------------------------------------------------------------}
451
 
   move.w  d4,d3
452
 
   mulu.w  d5,d3        { mulitply with bigit from multiplier      }
453
 
   move.l  d3,4(sp)     { store into result                        }
454
 
 
455
 
   move.l  d4,d3
456
 
   swap    d3
457
 
   mulu.w  d5,d3
458
 
   add.l   d3,2(sp)     { add to result                            }
459
 
 
460
 
   swap    d5           { [TOP 8 BITS SHOULD BE ZERO !]            }
461
 
 
462
 
   move.w  d4,d3
463
 
   mulu.w  d5,d3        { mulitply with bigit from multiplier      }
464
 
   add.l   d3,2(sp)     { store into result (no carry can occur here) }
465
 
 
466
 
   move.l  d4,d3
467
 
   swap    d3
468
 
   mulu.w  d5,d3
469
 
   add.l   d3,(sp)      { add to result                            }
470
 
            { [TOP 16 BITS SHOULD BE ZERO !] }
471
 
   movem.l 2(sp),d4-d5  { get the 48 valid mantissa bits           }
472
 
   clr.w   d5           { (pad to 64)                              }
473
 
 
474
 
   move.l  #$0000ffff,d3
475
 
@mlabel7:
476
 
   cmp.l   d3,d4        { multiply (shift) until                  }
477
 
   bhi     @mlabel8     {  1 in upper 16 result bits              }
478
 
   cmp.w   #9,d0        { give up for denormalized numbers        }
479
 
   ble     @mlabel8
480
 
   swap    d4         { (we''re getting here only when multiplying }
481
 
   swap    d5         {  with a denormalized number; there''s an   }
482
 
   move.w  d5,d4      {  eventual loss of 4 bits in the rounding   }
483
 
   clr.w   d5         {  byte -- what a pity 8-)                   }
484
 
   subq.w  #8,d0      { decrement exponent                         }
485
 
   subq.w  #8,d0
486
 
   bra     @mlabel7
487
 
@mlabel8:
488
 
   move.l  d5,d1      { get rounding bits                          }
489
 
   rol.l   #8,d1
490
 
   move.l  d1,d3      { see if sticky bit should be set            }
491
 
   and.l   #$ffffff00,d3
492
 
   beq     @mlabel9
493
 
   or.b    #1,d1      { set "sticky bit" if any low-order set      }
494
 
@mlabel9:
495
 
   addq.w  #8,sp      { remove accumulator from stack              }
496
 
   jmp     FPC_SINGLE_NORM{ (result in d4)                             }
497
 
 
498
 
@mretz:
499
 
   addq.w  #8,sp      { release accumulator space                  }
500
 
@mretzz:
501
 
   moveq.l #0,d0      { save zero as result                        }
502
 
   lsl.w   #1,d2      { and set it sign as for d2                  }
503
 
   roxr.l  #1,d0
504
 
   movem.l (sp)+,d2-d5
505
 
   rts                { no normalizing neccessary                  }
506
 
end;
507
 
 
508
 
 
509
 
Procedure Single_Div;Assembler;
510
 
{--------------------------------------------}
511
 
{ Low-level routine to dividr two single     }
512
 
{  IEEE floating point values. Never called  }
513
 
{  directly.                                 }
514
 
{ Om Entry:                                  }
515
 
{      d1/d0 = u/v = operation to perform.   }
516
 
{ On Exit:                                   }
517
 
{      d0 = result.                          }
518
 
{  Registers destroyed: d0,d1                }
519
 
{ stack space used (and restored): 8 bytes.  }
520
 
{--------------------------------------------}
521
 
ASM
522
 
XDEF SINGLE_DIV
523
 
   { u = d1 = dividend }
524
 
   { v = d0 = divisor  }
525
 
   tst.l  d0              { check if divisor is 0                  }
526
 
   bne    @dno_exception
527
 
 
528
 
   move.l #$7f800000,d0
529
 
   btst   #31,d1          { transfer sign of dividend              }
530
 
   beq    @dclear
531
 
   bset   #31,d0
532
 
@dclear:
533
 
   rts
534
 
 
535
 
@dno_exception:
536
 
   move.l  d1,d4          { d4 = u, d5 = v                         }
537
 
   move.l  d0,d5
538
 
   movem.l d2-d5,-(sp)    { save registers                         }
539
 
 
540
 
   move.l  #$7fffff,d3
541
 
   move.l  d4,d0          { d0 = u.exp                             }
542
 
   and.l   d3,d4          { remove exponent from u.mantissa        }
543
 
   swap    d0
544
 
   move.w  d0,d2          { d2 = u.sign                            }
545
 
 
546
 
   move.l  d5,d1          { d1 = v.exp                             }
547
 
   and.l   d3,d5          { remove exponent from v.mantissa        }
548
 
   swap    d1
549
 
   eor.w   d1,d2          { d2 = u.sign ^ v.sign (in bit 15)       }
550
 
 
551
 
   moveq.l #15,d3
552
 
   bclr    d3,d0          { kill sign bit                          }
553
 
   bclr    d3,d1          { kill sign bit                          }
554
 
   lsr.w   #7,d0
555
 
   lsr.w   #7,d1
556
 
 
557
 
   moveq.l #-1,d3
558
 
   cmp.b   d3,d0          { comparison with #0xff                  }
559
 
   beq     @dlabel1       { u == NaN ;; u== Inf                    }
560
 
   cmp.b   d3,d1
561
 
   beq     @dlabel2       { v == NaN ;; v == Inf                   }
562
 
   tst.b   d0
563
 
   bne     @dlabel4       { u not zero nor denorm                  }
564
 
   tst.l   d4
565
 
   beq     @dlabel3       { 0/ ?                                   }
566
 
 
567
 
@dlabel4:
568
 
   tst.w   d1
569
 
   bne     @dnospec
570
 
 
571
 
   tst.l   d5
572
 
   bne     @dnospec
573
 
   bra     @dretinf       { x/0 -> +/- Inf                         }
574
 
 
575
 
@dlabel1:
576
 
   tst.l   d4             { u == NaN ?                             }
577
 
   bne     @dretnan       { NaN/ x                                 }
578
 
   cmp.b   d3,d1
579
 
   beq     @dretnan       { Inf/Inf or Inf/NaN                     }
580
 
{  bra     dretinf      ; Inf/x ; x != Inf && x != NaN             }
581
 
{--------------------------------------------------------------------}
582
 
{                  Return infinity with correct sign.                }
583
 
{--------------------------------------------------------------------}
584
 
@dretinf:
585
 
   move.l  #$ff000000,d0
586
 
   lsl.w   #1,d2
587
 
   roxr.l  #1,d0          { shift in high bit as given by d2       }
588
 
@dreturn:
589
 
   movem.l (sp)+,d2-d5
590
 
   rts
591
 
 
592
 
@dlabel2:
593
 
   tst.l   d5
594
 
   bne     @dretnan       { x/NaN                                  }
595
 
{   bra    dretzero   ; x/Inf -> +/- 0                             }
596
 
{--------------------------------------------------------------------}
597
 
{                  Return correct signed zero.                       }
598
 
{--------------------------------------------------------------------}
599
 
@dretzero:
600
 
   moveq.l #0,d0          { zero destination                       }
601
 
   lsl.w   #1,d2          { set X bit accordingly                  }
602
 
   roxr.l  #1,d0
603
 
   bra     @dreturn
604
 
 
605
 
@dlabel3:
606
 
   tst.w   d1
607
 
   bne     @dretzero      { 0/x ->+/- 0                            }
608
 
   tst.l   d4
609
 
   bne     @dretzero      { 0/x                                    }
610
 
{   bra    dretnan         0/0                                     }
611
 
{--------------------------------------------------------------------}
612
 
{                       Return NotANumber                            }
613
 
{--------------------------------------------------------------------}
614
 
@dretnan:
615
 
   move.l  d3,d0          { d3 contains 0xffffffff                 }
616
 
   lsr.l   #1,d0
617
 
   bra     @dreturn
618
 
{--------------------------------------------------------------------}
619
 
{                      End of Special Handling                       }
620
 
{--------------------------------------------------------------------}
621
 
@dnospec:
622
 
   moveq.l #23,d3
623
 
   bset    d3,d4          { restore implied leading "1"            }
624
 
   tst.w   d0             { check for zero exponent - no leading "1" }
625
 
   bne     @dlabel5
626
 
   bclr    d3,d4          { remove it                              }
627
 
   add.w   #1,d0          { "normalize" exponent                   }
628
 
@dlabel5:
629
 
   tst.l   d4
630
 
   beq     @dretzero      { dividing zero                          }
631
 
 
632
 
   bset    d3,d5          { restore implied leading "1"            }
633
 
   tst.w   d1             { check for zero exponent - no leading "1"}
634
 
   bne     @dlabel6
635
 
   bclr    d3,d5          { remove it                              }
636
 
   add.w   #1,d1          { "normalize" exponent                   }
637
 
@dlabel6:
638
 
 
639
 
   sub.w   d1,d0          { subtract exponents,                    }
640
 
   add.w   #BIAS4-8+1,d0  { add bias back in, account for shift    }
641
 
   add.w   #34,d0     { add loop offset, +2 for extra rounding bits}
642
 
                      { for denormalized numbers (2 implied by dbra)}
643
 
   move.w  #27,d1     { bit number for "implied" pos (+4 for rounding)}
644
 
   moveq.l #-1,d3     { zero quotient (for speed a one''s complement) }
645
 
   sub.l   d5,d4      { initial subtraction, u = u - v                }
646
 
@dlabel7:
647
 
   btst    d1,d3      { divide until 1 in implied position            }
648
 
   beq     @dlabel9
649
 
 
650
 
   add.l   d4,d4
651
 
   bcs     @dlabel8   { if carry is set, add, else subtract           }
652
 
 
653
 
   addx.l  d3,d3      { shift quotient and set bit zero               }
654
 
   sub.l   d5,d4      { subtract, u = u - v                           }
655
 
   dbra    d0,@dlabel7      { give up if result is denormalized        }
656
 
   bra     @dlabel9
657
 
@dlabel8:
658
 
   addx.l  d3,d3      { shift quotient and clear bit zero             }
659
 
   add.l   d5,d4      { add (restore), u = u + v                      }
660
 
   dbra    d0,@dlabel7      { give up if result is denormalized        }
661
 
@dlabel9:
662
 
   subq.w  #2,d0      { remove rounding offset for denormalized nums  }
663
 
   not.l   d3         { invert quotient to get it right               }
664
 
 
665
 
   clr.l   d1         { zero rounding bits                            }
666
 
   tst.l   d4         { check for exact result                        }
667
 
   beq     @dlabel10
668
 
   moveq.l #-1,d1     { prevent tie case                              }
669
 
@dlabel10:
670
 
   move.l  d3,d4      { save quotient mantissa                        }
671
 
   jmp     FPC_SINGLE_NORM{ (registers on stack removed by norm_sf)       }
672
 
end;
673
 
 
674
 
 
675
 
Procedure Single_Cmp; Assembler;
676
 
{--------------------------------------------}
677
 
{ Low-level routine to compare single two    }
678
 
{  single point values..                     }
679
 
{  Never called directly.                    }
680
 
{ On Entry:                                  }
681
 
{      d1 and d0 Values to compare           }
682
 
{      d0 = first operand                    }
683
 
{ On Exit:                                   }
684
 
{      Flags according to result             }
685
 
{  Registers destroyed: d0,d1                }
686
 
{--------------------------------------------}
687
 
Asm
688
 
XDEF  SINGLE_CMP
689
 
   tst.l     d0               { check sign bit             }
690
 
   bpl       @cmplab1
691
 
   neg.l     d0               { negate                     }
692
 
   bchg      #31,d0           { toggle sign bit            }
693
 
@cmplab1:
694
 
   tst.l     d1               { check sign bit             }
695
 
   bpl       @cmplab2
696
 
   neg.l     d1               { negate                     }
697
 
   bchg      #31,d1           { toggle sign bit            }
698
 
@cmplab2:
699
 
   cmp.l     d0,d1            { compare...                 }
700
 
   rts
701
 
end;
702
 
 
703
 
 
704
 
 
705
 
Procedure LongMul;Assembler;
706
 
{--------------------------------------------}
707
 
{ Low-level routine to multiply two signed   }
708
 
{  32-bit values. Never called directly.     }
709
 
{ On entry: d1,d0 = 32-bit signed values to  }
710
 
{           multiply.                        }
711
 
{ On Exit:                                   }
712
 
{      d0 = result.                          }
713
 
{  Registers destroyed: d0,d1                }
714
 
{  stack space used and restored: 10 bytes   }
715
 
{--------------------------------------------}
716
 
Asm
717
 
XDEF LONGMUL
718
 
             cmp.b   #2,Test68000    { Are we on a 68020+ cpu                  }
719
 
             blt     @Lmulcontinue
720
 
             muls.l  d1,d0           { yes, then directly mul...               }
721
 
             rts                     { return... result in d0                  }
722
 
@Lmulcontinue:
723
 
             move.l    d2,a0         { save registers                          }
724
 
             move.l    d3,a1
725
 
 
726
 
             move.l    d0,-(sp)
727
 
             move.l    d1,-(sp)
728
 
 
729
 
             movem.w   (sp)+,d0-d3   { u = d0-d1, v = d2-d3                    }
730
 
 
731
 
             move.w    d0,-(sp)      { sign flag                               }
732
 
             bpl       @LMul1        { is u negative ?                         }
733
 
             neg.w     d1            { yes, force it positive                  }
734
 
             negx.w    d0
735
 
@LMul1:
736
 
             tst.w     d2            { is v negative ?                         }
737
 
             bpl       @LMul2
738
 
             neg.w     d3            { yes, force it positive ...              }
739
 
             negx.w    d2
740
 
             not.w     (sp)          {  ... and modify flag word               }
741
 
@LMul2:
742
 
             ext.l     d0            { u.h <> 0 ?                              }
743
 
             beq       @LMul3
744
 
             mulu.w    d3,d0         { r  = v.l * u.h                          }
745
 
@LMul3:
746
 
             tst.w     d2            { v.h <> 0 ?                              }
747
 
             beq       @LMul4
748
 
             mulu.w    d1,d2         { r += v.h * u.l                          }
749
 
             add.w     d2,d0
750
 
@LMul4:
751
 
             swap      d0
752
 
             clr.w     d0
753
 
             mulu.w    d3,d1        { r += v.l * u.l                           }
754
 
             add.l     d1,d0
755
 
             move.l    a1,d3
756
 
             move.l    a0,d2
757
 
             tst.w     (sp)+        { should the result be negated ?           }
758
 
             bpl       @LMul5       { no, just return                          }
759
 
             neg.l     d0           { else r = -r                              }
760
 
@LMul5:
761
 
             rts
762
 
end;
763
 
 
764
 
 
765
 
 
766
 
Procedure Long2Single;Assembler;
767
 
{--------------------------------------------}
768
 
{ Low-level routine to convert a longint     }
769
 
{  to a single floating point value.         }
770
 
{ On entry: d0 = longint value to convert.   }
771
 
{ On Exit:                                   }
772
 
{      d0 = single IEEE value                }
773
 
{  Registers destroyed: d0,d1                }
774
 
{  stack space used and restored:  8 bytes   }
775
 
{--------------------------------------------}
776
 
Asm
777
 
XDEF LONG2SINGLE
778
 
   movem.l  d2-d5,-(sp)  { save registers to make norm_sf happy}
779
 
 
780
 
   move.l   d0,d4        { prepare result mantissa          }
781
 
   move.w   #BIAS4+32-8,d0   { radix point after 32 bits    }
782
 
   move.l   d4,d2        { set sign flag                    }
783
 
   bge      @l2slabel1   { nonnegative                      }
784
 
   neg.l    d4           { take absolute value              }
785
 
@l2slabel1:
786
 
   swap     d2           { follow SINGLE_NORM conventions   }
787
 
   clr.w    d1           { set rounding = 0                 }
788
 
   jmp      FPC_SINGLE_NORM
789
 
end;
790
 
 
791
 
 
792
 
Procedure LongDiv; [alias : 'FPC_LONGDIV'];Assembler;
793
 
{--------------------------------------------}
794
 
{ Low-level routine to do signed long        }
795
 
{ division.                                  }
796
 
{ On entry: d0/d1 operation to perform       }
797
 
{ On Exit:                                   }
798
 
{      d0 = quotient                         }
799
 
{      d1 = remainder                        }
800
 
{  Registers destroyed: d0,d1,d6             }
801
 
{  stack space used and restored: 10 bytes   }
802
 
{--------------------------------------------}
803
 
asm
804
 
XDEF LONGDIV
805
 
   cmp.b       #2,Test68000  { can we use divs ?                  }
806
 
   blt         @continue
807
 
   tst.l       d1
808
 
   beq         @zerodiv2
809
 
   move.l      d1,d6
810
 
   clr.l       d1           { clr                                 }
811
 
   tst.l       d0           { check sign of d0                    }
812
 
   bpl         @posdiv
813
 
   move.l      #$ffffffff,d1{ sign extend into  d1                }
814
 
@posdiv:
815
 
   divsl.l     d6,d1:d0
816
 
   rts
817
 
@continue:
818
 
 
819
 
   move.l      d2,a0      { save registers                        }
820
 
   move.l      d3,a1
821
 
   move.l      d4,-(sp)   { divisor = d1 = d4                     }
822
 
   move.l      d5,-(sp)   { divident = d0 = d5                    }
823
 
 
824
 
   move.l      d1,d4      { save divisor                          }
825
 
   move.l      d0,d5      { save dividend                         }
826
 
 
827
 
   clr.w       -(sp)      { sign flag                             }
828
 
 
829
 
   clr.l       d0         { prepare result                        }
830
 
   move.l      d4,d2      { get divisor                           }
831
 
   beq         @zerodiv   { divisor = 0 ?                         }
832
 
   bpl         @LDiv1     { divisor < 0 ?                         }
833
 
   neg.l       d2         { negate it                             }
834
 
   not.w       (sp)       { remember sign                         }
835
 
@LDiv1:
836
 
   move.l      d5,d1      { get dividend                          }
837
 
   bpl         @LDiv2     { dividend < 0 ?                        }
838
 
   neg.l       d1         { negate it                             }
839
 
   not.w       (sp)       { remember sign                         }
840
 
@LDiv2:
841
 
{;== case 1) divident < divisor}
842
 
   cmp.l       d2,d1      { is divident smaller then divisor ?    }
843
 
   bcs         @LDiv7     { yes, return immediately               }
844
 
{;== case 2) divisor has <= 16 significant bits}
845
 
   move.l      d4,d6      { put divisor in d6 register            }
846
 
   lsr.l       #8,d6      { rotate into low word                  }
847
 
   lsr.l       #8,d6
848
 
   tst.l       d6
849
 
   bne         @LDiv3     { divisor has only 16 bits              }
850
 
   move.w      d1,d3      { save dividend                         }
851
 
   clr.w       d1         { divide dvd.h by dvs                   }
852
 
   swap        d1
853
 
   beq         @LDiv4     { (no division necessary if dividend zero)}
854
 
   divu        d2,d1
855
 
@LDiv4:
856
 
   move.w      d1,d0      { save quotient.h                       }
857
 
   swap        d0
858
 
   move.w      d3,d1      { (d0.h = remainder of prev divu)       }
859
 
   divu        d2,d1      { divide dvd.l by dvs                   }
860
 
   move.w      d1,d0      { save quotient.l                       }
861
 
   clr.w       d1         { get remainder                         }
862
 
   swap        d1
863
 
   bra         @LDiv7     { and return                            }
864
 
{;== case 3) divisor > 16 bits (corollary is dividend > 16 bits, see case 1)}
865
 
@LDiv3:
866
 
   moveq.l     #31,d3     { loop count                            }
867
 
@LDiv5:
868
 
   add.l       d1,d1      { shift divident ...                    }
869
 
   addx.l      d0,d0      {  ... into d0                          }
870
 
   cmp.l       d2,d0      { compare with divisor                  }
871
 
   bcs         @LDiv6
872
 
   sub.l       d2,d0      { big enough, subtract                  }
873
 
   addq.w      #1,d1      { and note bit into result              }
874
 
@LDiv6:
875
 
   dbra        d3,@LDiv5
876
 
   exg         d0,d1      { put quotient and remainder in their registers}
877
 
@LDiv7:
878
 
   tst.l       d5         { must the remainder be corrected ?     }
879
 
   bpl         @LDiv8
880
 
   neg.l       d1         { yes, apply sign                       }
881
 
{ the following line would be correct if modulus is defined as in algebra}
882
 
{;   add.l   sp@(6),d1   ; algebraic correction: modulus can only be >= 0}
883
 
@LDiv8:
884
 
   tst.w       (sp)+      { result should be negative ?          }
885
 
   bpl         @LDiv9
886
 
   neg.l       d0         { yes, negate it                      }
887
 
@LDiv9:
888
 
   move.l      a1,d3
889
 
   move.l      a0,d2
890
 
   move.l      (sp)+,d5
891
 
   move.l      (sp)+,d4
892
 
   rts                    { en exit : remainder = d1, quotient = d0 }
893
 
@zerodiv:
894
 
   move.l      a1,d3      { restore stack...                        }
895
 
   move.l      a0,d2
896
 
   move.w      (sp)+,d1   { remove sign word                        }
897
 
   move.l      (sp)+,d5
898
 
   move.l      (sp)+,d4
899
 
@zerodiv2:
900
 
   move.l      #200,d0
901
 
   jsr         FPC_HALT_ERROR { RunError(200)                          }
902
 
   rts                    { this should never occur...             }
903
 
end;
904
 
 
905
 
 
906
 
Procedure LongMod;[alias : 'FPC_LONGMOD'];Assembler;
907
 
{ see longdiv for info on calling convention }
908
 
asm
909
 
XDEF LONGMOD
910
 
   jsr     FPC_LONGDIV
911
 
   move.l  d1,d0      { return the remainder in d0 }
912
 
   rts
913
 
end;
914
 
 
915
 
{
916
 
  $Log: lowmath.inc,v $
917
 
  Revision 1.4  2002/09/07 16:01:20  peter
918
 
    * old logs removed and tabs fixed
919
 
 
920
 
}