Directory of image this file is from
This file as a plain text file

	SUBROUTINE COSOLV(G,H,A,B,BET,SD,RABE,IGA,C)
	DIMENSION G0(21),U0(21),V0(21),ISN(4,21)
	REAL M
	PI=3.1415926
	R=A/B
	DELT=ATAN((H/G)*(1./R))
	M=SQRT((G*A/(B**2))**2 + (H/B)**2)
	IF(ABS(H).GT..1E-5.AND.ABS(G).GT..1E-5) GO TO 20
	IH=1
	II=1
	IF(ABS(H).GT..1E-5) GO TO 22
C------ H = 0
	V0(IH)=0.0
	U0(IH)=G-A
	GO TO 12
C------ G = 0
22	U0(IH)=0.0
	V0(IH)=H-B
	GO TO 12
C------OUTER LOOP ZERO TO TWO PI
20	IH=0
	DO 1 IJ=1,2
	IF(IJ.EQ.2) DELT=-DELT
	BL=0.0
	BU=2.*PI
	DO 1 I=1,41
	GAM0=BL+(I-1)*(BU-BL)/40.
	GAM=GAM0
	Q1=((R*R-1.)/2.)*SIN(2.*GAM)
	Q2=M*SIN(GAM+DELT)
	DPR=Q1+Q2
	IF(ABS(DPR).LT..1E-6) GO TO 3
	SPR=ABS(DPR)/DPR
	IF(I.LT.2) GO TO 11
	IF(SPR.EQ.SPP) GO TO 11
C-------CROSSING DETECTED, GO INTO NEWTON-RAPHSON
	GAM=.5*(GAM0+GAMP)
6	FN=.5*(1.-R*R)*SIN(2.*GAM)-M*SIN(GAM+DELT)
	FNP=(1.-R*R)*COS(2.*GAM)-M*COS(GAM+DELT)
	GNP1=GAM-FN/FNP
	IF(ABS(GNP1-GAM).LT..1E-6) GO TO 3
	GAM=GNP1
	GO TO 6
C------FOUND ONE, WRITE IT
3	IH=IH+1
	G0(IH)=GAM
11	SPP=SPR
1	GAMP=GAM0
C------DETERMINE U0,V0, AND SIGNS
	ISN(1,1)=ABS(G)/G
	ISN(2,1)=ABS(H)/H
	DO 10 I=1,IH
	CS=COS(G0(I))
	SS=SIN(G0(I))
	U0(I)=G+A*COS(G0(I))
	V0(I)=H+B*SIN(G0(I))
	ISN(3,I)=ABS(U0(I))/U0(I)
10	ISN(4,I)=ABS(V0(I))/V0(I)
C------FIND THE CORRECT ROOT
	DO 16 I=1,IH
	II=I
	IF(ISN(1,1).NE.ISN(3,I).OR.ISN(2,1).NE.ISN(4,I)) GO TO 16
18	IF(ABS(U0(I)).LE.ABS(G).AND.ABS(V0(I)).LE.ABS(H)) GO TO 12
16	CONTINUE
	IF(IGA.EQ.0) GO TO 17
C	WRITE(5,15)
15	FORMAT('  NO SOLUTION FOUND')
C	WRITE(5,704) G,H,A,B,(G0(IK),U0(IK),V0(IK),IK=1,IH)
704	FORMAT(' G, H, A, B',/,4(3X,E12.5),/,
	28(3(4X,E16.5),/))
	IGA=2
	RETURN
C------CHECK TO SEE IF REF.ASSMNT. CORRECT
12	UV0=SQRT((U0(II))**2 + (V0(II))**2)
	T0= -U0(II)*(SIN(BET)/COS(BET))*(1.-(SD/(2.*UV0)))
	T=RABE+SQRT((SD/2.)**2-(G-U0(II))**2-(H-V0(II))**2)
	IF(T.LT.T0) GO TO 17
	IGA=0
	C=UV0-SD/2.
	IF(T.GT.0.) RETURN
	C=SQRT(T**2+(U0(II)**2+V0(II)**2)*(1.-(SD/(2.*UV0)))**2)
	RETURN
17	IGA=IGA+1
C	IF(IGA.EQ.2) WRITE(5,21) T0,T,U0(II),V0(II),UV0,BET
21	FORMAT(' NEITHER REF. SOL. WORKED, T0 = ',F10.4,4X,'T = ',
	2F10.4,/,' U0 = ',F10.4,'  V0 = ',F10.4,'  UV0 = ',F10.4
	3,'    BETA = ',F10.4)
	RETURN
	END