Mercurial > hg > octave-nkf

annotate libcruft/amos/xzsqrt.f @ 9584:0fcbfddaa87f

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
allow abbreviated graphics property names to match, with optional warning
author John W. Eaton <jwe@octave.org>
date Fri, 28 Aug 2009 05:30:29 -0400
parents 8b0cb8f79fdc
children
Ignore whitespace changes - Everywhere: Within whitespace: At end of lines:
rev   line source
3217
8b0cb8f79fdc [project @ 1998-11-11 23:32:20 by jwe]
jwe
parents:
diff changeset
1 SUBROUTINE XZSQRT(AR, AI, BR, BI)
8b0cb8f79fdc [project @ 1998-11-11 23:32:20 by jwe]
jwe
parents:
diff changeset
2 C***BEGIN PROLOGUE XZSQRT
8b0cb8f79fdc [project @ 1998-11-11 23:32:20 by jwe]
jwe
parents:
diff changeset
3 C***REFER TO ZBESH,ZBESI,ZBESJ,ZBESK,ZBESY,ZAIRY,ZBIRY
8b0cb8f79fdc [project @ 1998-11-11 23:32:20 by jwe]
jwe
parents:
diff changeset
4 C
8b0cb8f79fdc [project @ 1998-11-11 23:32:20 by jwe]
jwe
parents:
diff changeset
5 C DOUBLE PRECISION COMPLEX SQUARE ROOT, B=CSQRT(A)
8b0cb8f79fdc [project @ 1998-11-11 23:32:20 by jwe]
jwe
parents:
diff changeset
6 C
8b0cb8f79fdc [project @ 1998-11-11 23:32:20 by jwe]
jwe
parents:
diff changeset
7 C***ROUTINES CALLED XZABS
8b0cb8f79fdc [project @ 1998-11-11 23:32:20 by jwe]
jwe
parents:
diff changeset
8 C***END PROLOGUE XZSQRT
8b0cb8f79fdc [project @ 1998-11-11 23:32:20 by jwe]
jwe
parents:
diff changeset
9 DOUBLE PRECISION AR, AI, BR, BI, ZM, DTHETA, DPI, DRT
8b0cb8f79fdc [project @ 1998-11-11 23:32:20 by jwe]
jwe
parents:
diff changeset
10 DOUBLE PRECISION XZABS
8b0cb8f79fdc [project @ 1998-11-11 23:32:20 by jwe]
jwe
parents:
diff changeset
11 DATA DRT , DPI / 7.071067811865475244008443621D-1,
8b0cb8f79fdc [project @ 1998-11-11 23:32:20 by jwe]
jwe
parents:
diff changeset
12 1 3.141592653589793238462643383D+0/
8b0cb8f79fdc [project @ 1998-11-11 23:32:20 by jwe]
jwe
parents:
diff changeset
13 ZM = XZABS(AR,AI)
8b0cb8f79fdc [project @ 1998-11-11 23:32:20 by jwe]
jwe
parents:
diff changeset
14 ZM = DSQRT(ZM)
8b0cb8f79fdc [project @ 1998-11-11 23:32:20 by jwe]
jwe
parents:
diff changeset
15 IF (AR.EQ.0.0D+0) GO TO 10
8b0cb8f79fdc [project @ 1998-11-11 23:32:20 by jwe]
jwe
parents:
diff changeset
16 IF (AI.EQ.0.0D+0) GO TO 20
8b0cb8f79fdc [project @ 1998-11-11 23:32:20 by jwe]
jwe
parents:
diff changeset
17 DTHETA = DATAN(AI/AR)
8b0cb8f79fdc [project @ 1998-11-11 23:32:20 by jwe]
jwe
parents:
diff changeset
18 IF (DTHETA.LE.0.0D+0) GO TO 40
8b0cb8f79fdc [project @ 1998-11-11 23:32:20 by jwe]
jwe
parents:
diff changeset
19 IF (AR.LT.0.0D+0) DTHETA = DTHETA - DPI
8b0cb8f79fdc [project @ 1998-11-11 23:32:20 by jwe]
jwe
parents:
diff changeset
20 GO TO 50
8b0cb8f79fdc [project @ 1998-11-11 23:32:20 by jwe]
jwe
parents:
diff changeset
21 10 IF (AI.GT.0.0D+0) GO TO 60
8b0cb8f79fdc [project @ 1998-11-11 23:32:20 by jwe]
jwe
parents:
diff changeset
22 IF (AI.LT.0.0D+0) GO TO 70
8b0cb8f79fdc [project @ 1998-11-11 23:32:20 by jwe]
jwe
parents:
diff changeset
23 BR = 0.0D+0
8b0cb8f79fdc [project @ 1998-11-11 23:32:20 by jwe]
jwe
parents:
diff changeset
24 BI = 0.0D+0
8b0cb8f79fdc [project @ 1998-11-11 23:32:20 by jwe]
jwe
parents:
diff changeset
25 RETURN
8b0cb8f79fdc [project @ 1998-11-11 23:32:20 by jwe]
jwe
parents:
diff changeset
26 20 IF (AR.GT.0.0D+0) GO TO 30
8b0cb8f79fdc [project @ 1998-11-11 23:32:20 by jwe]
jwe
parents:
diff changeset
27 BR = 0.0D+0
8b0cb8f79fdc [project @ 1998-11-11 23:32:20 by jwe]
jwe
parents:
diff changeset
28 BI = DSQRT(DABS(AR))
8b0cb8f79fdc [project @ 1998-11-11 23:32:20 by jwe]
jwe
parents:
diff changeset
29 RETURN
8b0cb8f79fdc [project @ 1998-11-11 23:32:20 by jwe]
jwe
parents:
diff changeset
30 30 BR = DSQRT(AR)
8b0cb8f79fdc [project @ 1998-11-11 23:32:20 by jwe]
jwe
parents:
diff changeset
31 BI = 0.0D+0
8b0cb8f79fdc [project @ 1998-11-11 23:32:20 by jwe]
jwe
parents:
diff changeset
32 RETURN
8b0cb8f79fdc [project @ 1998-11-11 23:32:20 by jwe]
jwe
parents:
diff changeset
33 40 IF (AR.LT.0.0D+0) DTHETA = DTHETA + DPI
8b0cb8f79fdc [project @ 1998-11-11 23:32:20 by jwe]
jwe
parents:
diff changeset
34 50 DTHETA = DTHETA*0.5D+0
8b0cb8f79fdc [project @ 1998-11-11 23:32:20 by jwe]
jwe
parents:
diff changeset
35 BR = ZM*DCOS(DTHETA)
8b0cb8f79fdc [project @ 1998-11-11 23:32:20 by jwe]
jwe
parents:
diff changeset
36 BI = ZM*DSIN(DTHETA)
8b0cb8f79fdc [project @ 1998-11-11 23:32:20 by jwe]
jwe
parents:
diff changeset
37 RETURN
8b0cb8f79fdc [project @ 1998-11-11 23:32:20 by jwe]
jwe
parents:
diff changeset
38 60 BR = ZM*DRT
8b0cb8f79fdc [project @ 1998-11-11 23:32:20 by jwe]
jwe
parents:
diff changeset
39 BI = ZM*DRT
8b0cb8f79fdc [project @ 1998-11-11 23:32:20 by jwe]
jwe
parents:
diff changeset
40 RETURN
8b0cb8f79fdc [project @ 1998-11-11 23:32:20 by jwe]
jwe
parents:
diff changeset
41 70 BR = ZM*DRT
8b0cb8f79fdc [project @ 1998-11-11 23:32:20 by jwe]
jwe
parents:
diff changeset
42 BI = -ZM*DRT
8b0cb8f79fdc [project @ 1998-11-11 23:32:20 by jwe]
jwe
parents:
diff changeset
43 RETURN
8b0cb8f79fdc [project @ 1998-11-11 23:32:20 by jwe]
jwe
parents:
diff changeset
44 END