← Back to branch summary

~ubuntu-branches/ubuntu/hoary/scilab/hoary

« back to all changes in this revision

Viewing changes to routines/calelm/wwdiv.f

  • Committer: Bazaar Package Importer
  • Author(s): Torsten Werner
  • Date: 2005-01-09 22:58:21 UTC
  • mfrom: (1.1.1 upstream)
  • Revision ID: james.westby@ubuntu.com-20050109225821-473xr8vhgugxxx5j
Tags: 3.0-12
http://bugs.debian.org/255869
changed configure.in to build scilab's own malloc.o, closes: #255869

Show diffs side-by-side

added added

removed removed

Lines of Context:
routines/calelm/wwdiv.f
1
 
      subroutine wwdiv(ar,ai,br,bi,cr,ci,ierr)
2
 
c!but
3
 
c
4
 
c     This subroutine wwdiv computes c=a/b where a and b are complex numbers
5
 
c
6
 
c!Calling sequence
7
 
c
8
 
c     subroutine wwdiv(ar,ai,br,bi,cr,ci,ierr)
9
 
c
10
 
c     ar, ai: double precision, a real and complex parts.
11
 
c
12
 
c     br, bi: double precision, b real and complex parts.
13
 
c
14
 
c     cr, ci: double precision, c real and complex parts.
15
 
c
16
 
c!author
17
 
c
18
 
c     Serge Steer
19
 
c
20
 
c!
21
 
c     Copyright INRIA
22
 
      double precision ar,ai,br,bi,cr,ci
23
 
c     c = a/b
24
 
      double precision s,d,ars,ais,brs,bis
25
 
c
26
 
      ierr=0
27
 
      if(bi.eq.0.0d0) then
28
 
         cr=ar/br
29
 
         ci=ai/br
30
 
      elseif(br.eq.0.0d0) then
31
 
         cr=ai/bi
32
 
         ci=-ar/bi
 
1
      subroutine wwdiv(ar, ai, br, bi, cr, ci, ierr)
 
2
*
 
3
*     PURPOSE
 
4
*        complex division algorithm : compute c := a / b
 
5
*        where :
 
6
*       
 
7
*        a = ar + i ai
 
8
*        b = br + i bi
 
9
*        c = cr + i ci
 
10
*
 
11
*        inputs  : ar, ai, br, bi  (double precision)
 
12
*        outputs : cr, ci          (double precision)
 
13
*                  ierr  (integer) ierr = 1 if b = 0 (else 0)
 
14
*
 
15
*     IMPLEMENTATION NOTES
 
16
*        1/ Use scaling with ||b||_oo; the original wwdiv.f used a scaling
 
17
*           with ||b||_1.  It results fewer operations.  From the famous
 
18
*           Golberg paper.  This is known as Smith's method.
 
19
*        2/ Currently set c = NaN + i NaN in case of a division by 0 ;
 
20
*           is that the good choice ?
 
21
*
 
22
*     AUTHOR
 
23
*        Bruno Pincon <Bruno.Pincon@iecn.u-nancy.fr>
 
24
*
 
25
      implicit none
 
26
      
 
27
*     PARAMETERS
 
28
      double precision ar, ai, br, bi, cr, ci
 
29
      integer ierr
 
30
 
 
31
*     LOCAL VARIABLES
 
32
      double precision r, d
 
33
 
 
34
*     TEXT
 
35
      ierr = 0
 
36
 
 
37
*     Treat special cases
 
38
      if (bi .eq. 0d0) then
 
39
         if (br .eq. 0d0) then  
 
40
            ierr = 1    
 
41
*           got NaN + i NaN
 
42
            cr = bi / br
 
43
            ci = cr
 
44
         else
 
45
            cr = ar / br
 
46
            ci = ai / br
 
47
         endif
 
48
      elseif (br .eq. 0d0) then
 
49
         cr = ai / bi
 
50
         ci = (-ar) / bi
33
51
      else
34
 
         s = abs(br) + abs(bi)
35
 
         if (s .eq. 0.0d+0) then
36
 
            ierr=1
37
 
            cr=ar/s
38
 
            ci=ai/s
39
 
            return
 
52
*     Generic division algorithm
 
53
         if (abs(br) .ge. abs(bi)) then
 
54
            r = bi / br
 
55
            d = br + r*bi
 
56
            cr = (ar + ai*r) / d
 
57
            ci = (ai - ar*r) / d
 
58
         else
 
59
            r = br / bi
 
60
            d = bi + r*br
 
61
            cr = (ar*r + ai) / d
 
62
            ci = (ai*r - ar) / d
40
63
         endif
41
 
         ars = ar/s
42
 
         ais = ai/s
43
 
         brs = br/s
44
 
         bis = bi/s
45
 
         d = brs**2 + bis**2
46
 
         cr = (ars*brs + ais*bis)/d
47
 
         ci = (ais*brs - ars*bis)/d
48
64
      endif
49
 
      return
 
65
      
50
66
      end
 
67
 
 
68