1 /* ******************************************************
 2    *                                                    *
 3    *                                                    *
 4    * Copyright (c) 1972 by Massachusetts Institute of   *
 5    * Technology and Honeywell Information Systems, Inc. *
 6    *                                                    *
 7    *                                                    *
 8    ****************************************************** */
 9 
10 dcsqrt_: proc (number) returns (complex float bin (63));
11 
12 dcl (number, a, b) complex float bin (63);
13 
14           a = number;
15           real (b) = abs (real (a));
16 
17           if imag (a) = 0.0e0
18           then do;
19                real (b) = sqrt (real (b));
20                imag (b) = 0.0e0;
21           end;
22           else do;
23                real (b) = sqrt ((abs (a)+real (b))*0.5e0);
24 
25                if real (a)<0.0e0
26                     & imag (a)<0.0e0
27                     then real (b) = -real (b);
28 
29                imag (b) = imag (a) * 0.5e0/real (b);
30           end;
31 
32           if real (a) >= 0.0e0 then return (b);
33 
34           real (a) = imag (b);
35           imag (a) = real (b);
36 
37           return (a);
38 
39      end dcsqrt_;