1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
|
subroutine wwdiv(ar,ai,br,bi,cr,ci,ierr)
c!but
c
c This subroutine wwdiv computes c=a/b where a and b are complex numbers
c
c!Calling sequence
c
c subroutine wwdiv(ar,ai,br,bi,cr,ci,ierr)
c
c ar, ai: double precision, a real and complex parts.
c
c br, bi: double precision, b real and complex parts.
c
c cr, ci: double precision, c real and complex parts.
c
c!author
c
c Serge Steer
c
c!
c Copyright INRIA
double precision ar,ai,br,bi,cr,ci
c c = a/b
double precision s,d,ars,ais,brs,bis
c
ierr=0
if(bi.eq.0.0d0) then
cr=ar/br
ci=ai/br
elseif(br.eq.0.0d0) then
cr=ai/bi
ci=-ar/bi
else
s = abs(br) + abs(bi)
if (s .eq. 0.0d+0) then
ierr=1
cr=ar/s
ci=ai/s
return
endif
ars = ar/s
ais = ai/s
brs = br/s
bis = bi/s
d = brs**2 + bis**2
cr = (ars*brs + ais*bis)/d
ci = (ais*brs - ars*bis)/d
endif
return
end
|
