| 1 | /* |
| 2 | * sqrt(a^2 + b^2) |
| 3 | * (but carefully) |
| 4 | */ |
| 5 | |
| 6 | double sqrt(); |
| 7 | double |
| 8 | hypot(a,b) |
| 9 | double a,b; |
| 10 | { |
| 11 | double t; |
| 12 | if(a<0) a = -a; |
| 13 | if(b<0) b = -b; |
| 14 | if(a > b) { |
| 15 | t = a; |
| 16 | a = b; |
| 17 | b = t; |
| 18 | } |
| 19 | if(b==0) return(0.); |
| 20 | a /= b; |
| 21 | /* |
| 22 | * pathological overflow possible |
| 23 | * in the next line. |
| 24 | */ |
| 25 | return(b*sqrt(1. + a*a)); |
| 26 | } |
| 27 | |
| 28 | struct complex |
| 29 | { |
| 30 | double r; |
| 31 | double i; |
| 32 | }; |
| 33 | |
| 34 | double |
| 35 | cabs(arg) |
| 36 | struct complex arg; |
| 37 | { |
| 38 | double hypot(); |
| 39 | |
| 40 | return(hypot(arg.r, arg.i)); |
| 41 | } |