464 lines
13 KiB
Fortran
Executable File
464 lines
13 KiB
Fortran
Executable File
SUBROUTINE T_RLP ! File [MAD.SRC]T_RLP.FOR
|
||
c ================
|
||
c
|
||
cdec$ Ident 'V01A'
|
||
c------------------------------------------------------------------
|
||
c Updates:
|
||
c V01A 7-May-1996 DM. Put error output to IOLUN, use IMPLICIT NONE and
|
||
c get the code indented so that it is readable!
|
||
c------------------------------------------------------------------
|
||
c Routines to deal with the reciprocical lattice PB
|
||
c------------------------------------------------------------------
|
||
c Entry points in this file:
|
||
c
|
||
c SETRLP : CALCULATION OF S AND INVS , ORIENTATION MATRIX
|
||
C RL2SPV : TRANSFO FROM RECIP LAT TO SCAT PLANE
|
||
C SP2RLV : TRANSFO FROM SCAT PLANE TO RECIP LAT
|
||
C INVS : INVERT MATRIX S, GENERATED BY SETRLP.
|
||
C ERRESO : DEAL ITH ERROR MESSAGES FOR ALL MODULES
|
||
C
|
||
C SUBROUTINE SETRLP(SAM,IER)
|
||
C SUBROUTINE RL2SPV(QHKL,QT,QM,QS,IER)
|
||
C SUBROUTINE SP2RLV(QHKL,QT,QM,QS,IER)
|
||
C SUBROUTINE INVS(S,SINV,IER)
|
||
C SUBROUTINE ERRESO(MODULE,IER)
|
||
c------------------------------------------------------------------
|
||
implicit none
|
||
end
|
||
c
|
||
SUBROUTINE SETRLP (SAM, IER) ! Part of [MAD.SRC]T_RLP.FOR
|
||
c ============================
|
||
c
|
||
C SETRLP: Computation of matrix S which transforms (QH,QK,QL) to
|
||
C vector (Q1,Q2) in scattering plane (defined by vectors A1,A2)
|
||
C and SINV matrix for the inverse transformation
|
||
C
|
||
C INPUT SAM SAMPLE CHARACTERISTICS
|
||
C SAM(1)=AS LATTICE PARAMETERS
|
||
C SAM(2)=BS ------------------
|
||
C SAM(3)=CS ------------------
|
||
C SAM(4)=AA LATTICE ANGLES
|
||
C SAM(5)=BB --------------
|
||
C SAM(6)=CC --------------
|
||
C SAM(7)=AX VECTOR A IN SCATTERING PLANE
|
||
C SAM(8)=AY ----------------------------
|
||
C SAM(9)=AZ ----------------------------
|
||
C SAM(10)=BX VECTOR B IN SCATTERING PLANE
|
||
C SAM(11)=BY ----------------------------
|
||
C SAM(12)=BZ ----------------------------
|
||
C OUTPUT IER ERROR RETURN TO BE TREATED BY ERRESO
|
||
C IER=1 ERROR ON LATTICE PARAMETERS
|
||
C IER=2 ERROR ON LATTICE ANGLES
|
||
C IER=3 ERROR ON VECTORS A1, A2
|
||
c------------------------------------------------------------------
|
||
IMPLICIT NONE
|
||
c
|
||
DOUBLE PRECISION PI
|
||
PARAMETER (PI=3.14159265358979323846264338327950D0)
|
||
DOUBLE PRECISION RD
|
||
PARAMETER (RD=57.29577951308232087679815481410517D0)
|
||
DOUBLE PRECISION EPS
|
||
PARAMETER (EPS=1.D-1)
|
||
DOUBLE PRECISION EPS6
|
||
PARAMETER (EPS6=1.D-6)
|
||
DOUBLE PRECISION EPS8
|
||
PARAMETER (EPS8=1.D-8)
|
||
c------------------------------------------------------------------
|
||
C Define the dummy arguments
|
||
real*4 sam(64)
|
||
integer*4 ier
|
||
c------------------------------------------------------------------
|
||
C DO NOT EXPORT THE FOLLOWING COMMON !
|
||
C IT IS JUST FOR PERMANENT STORAGE USE
|
||
C
|
||
DOUBLE PRECISION S, SINV
|
||
integer*4 iok
|
||
COMMON /OSOLEM/ S(4,4), SINV(4,4), iok
|
||
C-----------------------------------------------------------------------
|
||
double precision a(3), aspv(3,2), rlb(3,2)
|
||
double precision alfa(3), sina(3), sinb(3), cosa(3), cosb(3)
|
||
double precision b(3), c(3)
|
||
double precision vv(3,3), bb(3,3)
|
||
C
|
||
double precision zp, cc
|
||
integer*4 imod
|
||
integer*4 id, ie
|
||
integer*4 jd, je, jf
|
||
integer*4 kg
|
||
integer*4 lh, lf
|
||
integer*4 md, me
|
||
integer*4 ne
|
||
C-----------------------------------------------------------------------
|
||
C SOME TESTS AND INIT OF CALCUALTION
|
||
C
|
||
ier = 0
|
||
IMOD = 1
|
||
ZP = 2.D0*PI
|
||
IOK = 0
|
||
DO ID = 1,3
|
||
A(ID) = SAM(ID)
|
||
ALFA(ID) = SAM(ID+3)
|
||
ASPV(ID,1) = SAM(6+ID)
|
||
ASPV(ID,2) = SAM(9+ID)
|
||
ENDDO
|
||
C
|
||
DO ID = 1,3
|
||
IER = 1
|
||
IF(ABS(A(ID)).LE.EPS8) GOTO 999
|
||
IER = 0
|
||
ENDDO
|
||
DO ID = 1,3
|
||
A(ID) = A(ID)/ZP
|
||
ALFA(ID) = ALFA(ID)/RD
|
||
COSA(ID) = COS(ALFA(ID))
|
||
SINA(ID) = SIN(ALFA(ID))
|
||
ENDDO
|
||
CC = COSA(1)*COSA(1)+COSA(2)*COSA(2)+COSA(3)*COSA(3)
|
||
CC = (1.D0+2.D0*COSA(1)*COSA(2)*COSA(3)-CC)
|
||
IER = 2
|
||
IF(CC.LE.EPS) GOTO 999
|
||
IER = 0
|
||
CC = SQRT(CC)
|
||
JE = 2
|
||
KG = 3
|
||
DO ID = 1,3
|
||
B(ID) = SINA(ID)/(A(ID)*CC)
|
||
COSB(ID) = (COSA(JE)*COSA(KG)-COSA(ID))/(SINA(JE)*SINA(KG))
|
||
SINB(ID) = SQRT(1.D0-COSB(ID)*COSB(ID))
|
||
RLB(ID,2) = ABS(ATAN(SINB(ID)/COSB(ID)))*RD
|
||
JE = KG
|
||
KG = ID
|
||
ENDDO
|
||
BB(1,1) = B(1)
|
||
BB(2,1) = 0.D0
|
||
BB(3,1) = 0.D0
|
||
BB(1,2) = B(2)*COSB(3)
|
||
BB(2,2) = B(2)*SINB(3)
|
||
BB(3,2) = 0.D0
|
||
BB(1,3) = B(3)*COSB(2)
|
||
BB(2,3) = -B(3)*SINB(2)*COSA(1)
|
||
BB(3,3) = 1.D0/A(3)
|
||
C
|
||
DO ID = 1,3
|
||
RLB(ID,1) = 0.D0
|
||
DO JE = 1,3
|
||
RLB(ID,1) = RLB(ID,1)+BB(JE,ID)**2
|
||
ENDDO
|
||
IER = 1
|
||
IF (ABS(RLB(ID,1)).LE.EPS8) GOTO 999
|
||
IER = 0
|
||
RLB(ID,1) = SQRT(RLB(ID,1))
|
||
ENDDO
|
||
C-----------------------------------------------------------------------
|
||
C GENERATION OF S ORIENTATION MATRIX REC. LATTICE TO SCATTERING PLANE
|
||
C
|
||
DO KG = 1,2
|
||
DO IE = 1,3
|
||
VV(KG,IE) = 0.D0
|
||
DO JF = 1,3
|
||
VV(KG,IE) = VV(KG,IE)+BB(IE,JF)*ASPV(JF,KG)
|
||
ENDDO
|
||
ENDDO
|
||
ENDDO
|
||
DO MD = 3,2,-1
|
||
DO NE = 1,3
|
||
ID = MOD(MD,3)+1
|
||
JE = MOD(MD+1,3)+1
|
||
KG = MOD(NE,3)+1
|
||
LH = MOD(NE+1,3)+1
|
||
VV(MD,NE) = VV(ID,KG)*VV(JE,LH)-VV(ID,LH)*VV(JE,KG)
|
||
ENDDO
|
||
ENDDO
|
||
C
|
||
DO ID = 1,3
|
||
C(ID) = 0.D0
|
||
DO JE = 1,3
|
||
C(ID) = C(ID)+VV(ID,JE)**2
|
||
ENDDO
|
||
IER = 3
|
||
IF (ABS(C(ID)).LE.EPS6) GOTO 999
|
||
IER = 0
|
||
C(ID) = SQRT(C(ID))
|
||
ENDDO
|
||
C
|
||
DO ID = 1,3
|
||
DO JE = 1,3
|
||
VV(JE,ID) = VV(JE,ID)/C(JE)
|
||
ENDDO
|
||
ENDDO
|
||
DO KG = 1,3
|
||
DO ME = 1,3
|
||
S(KG,ME) = 0.D0
|
||
DO LF = 1,3
|
||
S(KG,ME) = S(KG,ME)+VV(KG,LF)*BB(LF,ME)
|
||
ENDDO
|
||
ENDDO
|
||
ENDDO
|
||
S(4,4) = 1.D0
|
||
DO JD = 1,3
|
||
S(4,JD) = 0.D0
|
||
S(JD,4) = 0.D0
|
||
ENDDO
|
||
C-----------------------------------------------------------------------
|
||
C INVERT TRANSFORMATION MATRIX S AND PU RESULT IN SINV
|
||
C
|
||
IER = 3
|
||
CALL INVS(S,SINV,IER)
|
||
IER = 0
|
||
IF (IER.NE.0)GOTO 999
|
||
IOK = 123
|
||
C---------------------------------------------------------------------------
|
||
C SORTIE
|
||
C
|
||
999 CONTINUE
|
||
IF (IER.NE.0) CALL ERRESO(IMOD,IER)
|
||
RETURN
|
||
END
|
||
c
|
||
SUBROUTINE RL2SPV (QHKL, QT, QM, QS, IER) ! Part of [MAD.SRC]T_RLP.FOR
|
||
c =========================================
|
||
c
|
||
c------------------------------------------------------------------
|
||
C INPUT QHKL QHKL -> QT
|
||
C A Q-VECTOR TO BE TRANSFORM FROM RECIP LATTICE TO SCATTERING PLANE
|
||
C CHECK THAT Q-VECTOR IS IN THE PLANE
|
||
C
|
||
C INPUT Q-VECTOR QHKL(3) Q-VECTOR IN RECIPROCICAL LATTICVE
|
||
C
|
||
C OUTPUT Q-VECTOR QT(3) Q-VECTOR IN SCATTERING PLANE
|
||
C OUTPUT QM AND QS QMODULUS AND ITS SQUARE ( TO BE VERIFIED )
|
||
C OUTPUT ERROR IER
|
||
C IER=1 MATRIX S NOT OK
|
||
C IER=2 Q NOT IN SCATTERING PLANE
|
||
C IER=3 Q MODULUS TOO SMALL
|
||
c------------------------------------------------------------------
|
||
IMPLICIT NONE
|
||
c
|
||
DOUBLE PRECISION DEPS4
|
||
PARAMETER (DEPS4=1.D-4)
|
||
DOUBLE PRECISION DEPS8
|
||
PARAMETER (DEPS8=1.D-8)
|
||
c------------------------------------------------------------------
|
||
c Define the dummy arguments
|
||
DOUBLE PRECISION QHKL(3)
|
||
DOUBLE PRECISION QT(3)
|
||
DOUBLE PRECISION QM
|
||
DOUBLE PRECISION QS
|
||
integer*4 IER
|
||
c------------------------------------------------------------------
|
||
C DO NOT EXPORT THE FOLLWING COOMON !
|
||
C IT IS JUST FOR PERMANENT STORAGE USE
|
||
C
|
||
DOUBLE PRECISION S, SINV
|
||
integer*4 iok
|
||
COMMON /OSOLEM/ S(4,4), SINV(4,4), iok
|
||
C---------------------------------------------------------------------------
|
||
integer*4 id, je
|
||
C---------------------------------------------------------------------------
|
||
C INIT AND TEST IF TRANSFO MATRICES ARE OK
|
||
IER = 1
|
||
IF (IOK.NE.123) GOTO 999
|
||
IER = 0
|
||
C-----------------------------------------------------------------------
|
||
C
|
||
DO ID = 1,3
|
||
QT(ID) = 0.D0
|
||
DO JE = 1,3
|
||
QT(ID) = QT(ID) + QHKL(JE)*S(ID,JE)
|
||
ENDDO
|
||
ENDDO
|
||
IER = 2
|
||
IF(ABS(QT(3)).GT.DEPS4) GOTO 999
|
||
IER = 0
|
||
QS = 0.D0
|
||
DO ID = 1,3
|
||
QS = QS+QT(ID)**2
|
||
ENDDO
|
||
IF(QS.LT.DEPS8) THEN
|
||
IER = 3
|
||
ELSE
|
||
QM = SQRT(QS)
|
||
ENDIF
|
||
C---------------------------------------------------------------------------
|
||
C
|
||
999 CONTINUE
|
||
RETURN
|
||
END
|
||
c
|
||
SUBROUTINE SP2RLV (QHKL, QT, QM, QS, IER) ! Part of [MAD.SRC]T_RLP.FOR
|
||
c =========================================
|
||
c
|
||
c------------------------------------------------------------------
|
||
C INPUT QT QHKL <- QT
|
||
C A Q-VECTOR TO BE TRANSFORM FROM SCATTERING PLANE TO RECIP LATTICE
|
||
C CHECK THAT Q, D & G VECTORS ARE IN THE SCATTERING PLANE
|
||
C
|
||
C INPUT Q-VECTOR QT(3) Q-VECTOR IN SCATTERING PLANE
|
||
C
|
||
C OUTPUT Q-VECTOR QHKL(3) Q-VECTOR IN RECIPROCICAL LATTICVE
|
||
C OUTPUT QM AND QS QMODULUS AND ITS SQUARE ( TO BE VERIFIED )
|
||
C OUTPUT ERROR IER
|
||
C IER=1 MATRIX S NOT OK
|
||
C IER=2 Q NOT IN SCATTERING PLANE
|
||
C IER=3 Q MODULUS TOO SMALL
|
||
c------------------------------------------------------------------
|
||
IMPLICIT NONE
|
||
C
|
||
DOUBLE PRECISION EPS4
|
||
PARAMETER (EPS4=1.D-4)
|
||
DOUBLE PRECISION EPS8
|
||
PARAMETER (EPS8=1.D-8)
|
||
c------------------------------------------------------------------
|
||
c Define the dummy arguments
|
||
DOUBLE PRECISION QHKL(3)
|
||
DOUBLE PRECISION QT(3)
|
||
DOUBLE PRECISION QM
|
||
DOUBLE PRECISION QS
|
||
integer*4 IER
|
||
c------------------------------------------------------------------
|
||
C DO NOT EXPORT THE FOLLWING COOMON !
|
||
C IT IS JUST FOR PERMANENT STORAGE USE
|
||
C
|
||
DOUBLE PRECISION S, SINV
|
||
integer*4 iok
|
||
COMMON /OSOLEM/ S(4,4), SINV(4,4), iok
|
||
C---------------------------------------------------------------------------
|
||
integer*4 id, je
|
||
C---------------------------------------------------------------------------
|
||
C INIT AND TEST IF TRANSFO MATRICES ARE OK
|
||
IER = 1
|
||
IF (IOK.NE.123) GOTO 999
|
||
IER = 2
|
||
IF(ABS(QT(3)).GT.EPS4) GOTO 999
|
||
IER = 0
|
||
C-----------------------------------------------------------------------
|
||
QS = 0.D0
|
||
DO ID = 1,3
|
||
QS = QS+QT(ID)**2
|
||
ENDDO
|
||
IF(QS.LT.EPS8) THEN
|
||
IER = 3
|
||
ELSE
|
||
QM = SQRT(QS)
|
||
ENDIF
|
||
C-----------------------------------------------------------------------
|
||
C
|
||
DO ID = 1,3
|
||
QHKL(ID) = 0.D0
|
||
DO JE = 1,3
|
||
QHKL(ID) = QHKL(ID)+SINV(ID,JE)*QT(JE)
|
||
ENDDO
|
||
ENDDO
|
||
C---------------------------------------------------------------------------
|
||
C
|
||
999 CONTINUE
|
||
RETURN
|
||
END
|
||
c
|
||
SUBROUTINE INVS (S, SINV, IER) ! Part of [MAD.SRC]T_RLP.FOR
|
||
c ==============================
|
||
c
|
||
c------------------------------------------------------------------
|
||
C ROUTINE TO INVERT MATRIX S, GENERATED BY SETRLP, WHICH TRANSFORMS
|
||
C (QH,QK,QL) TO (Q1,Q2) IN THE SCATTERING PLANE
|
||
C INPUT MATRIX DOUBLE PRECISION S(4,4)
|
||
C OUTPUT MATRIX DOUBLE PRECISION SINV(4,4)
|
||
C OUTPUT ERROR IER
|
||
C IER=1 DETERMINANT OF MATRIX S TOO SMALL
|
||
c------------------------------------------------------------------
|
||
IMPLICIT NONE
|
||
c
|
||
DOUBLE PRECISION EPS6
|
||
PARAMETER (EPS6=1.D-6)
|
||
c------------------------------------------------------------------
|
||
c Define the dummy arguments
|
||
DOUBLE PRECISION S(4,4)
|
||
DOUBLE PRECISION SINV(4,4)
|
||
integer*4 IER
|
||
c------------------------------------------------------------------
|
||
integer*4 m(3)
|
||
integer*4 n(3)
|
||
c
|
||
integer*4 id, je
|
||
integer*4 mi, mj, ni, nj
|
||
double precision det
|
||
c------------------------------------------------------------------
|
||
DATA M/2,3,1/
|
||
DATA N/3,1,2/
|
||
IER = 0
|
||
DO ID = 1,4
|
||
DO JE = 1,4
|
||
SINV(ID,JE) = 0.D0
|
||
ENDDO
|
||
ENDDO
|
||
DET = 0.D0
|
||
DO ID = 1,3
|
||
DO JE = 1,3
|
||
MI = M(ID)
|
||
MJ = M(JE)
|
||
NI = N(ID)
|
||
NJ = N(JE)
|
||
SINV(JE,ID) = S(MI,MJ)*S(NI,NJ)-S(NI,MJ)*S(MI,NJ)
|
||
ENDDO
|
||
DET = DET+S(ID,1)*SINV(1,ID)
|
||
ENDDO
|
||
IF(ABS(DET).LT.EPS6) THEN
|
||
IER = 1
|
||
ELSE
|
||
DO ID = 1,3
|
||
DO JE = 1,3
|
||
SINV(ID,JE) = SINV(ID,JE)/DET
|
||
ENDDO
|
||
ENDDO
|
||
ENDIF
|
||
SINV(4,4) = 1.D0
|
||
RETURN
|
||
END
|
||
c
|
||
SUBROUTINE ERRESO(MODULE,IER) ! Part of [MAD.SRC]T_RLP.FOR
|
||
c =============================
|
||
c
|
||
c------------------------------------------------------------------
|
||
C SUBROUTINE TO TREAT ERRORS FROM RESOLUTION CALCULATIONS
|
||
C MODULE = 1 -> SETRLP
|
||
C MODULE = 2 -> RL2SPV
|
||
C MODULE = 3 -> EX_CASE
|
||
c------------------------------------------------------------------
|
||
IMPLICIT NONE
|
||
c
|
||
integer*4 MXER
|
||
PARAMETER (MXER = 4)
|
||
integer*4 MXMOD
|
||
PARAMETER (MXMOD = 3)
|
||
c
|
||
C include 'iolsddef.inc'
|
||
c------------------------------------------------------------------
|
||
c Define the dummy arguments
|
||
integer*4 module
|
||
integer*4 ier
|
||
c------------------------------------------------------------------
|
||
integer*4 lmodule, lier
|
||
CHARACTER*64 MESER(MXER,MXMOD)
|
||
DATA MESER/
|
||
+ ' Check lattice spacings (AS,BS,CS)',
|
||
+ ' Check cell angles (AA,BB,CC)',
|
||
+ ' Check scattering plane (AX....BZ)',
|
||
+ ' Check lattice and scattering plane',
|
||
C
|
||
+ ' Check Lattice and Scattering Plane',
|
||
+ ' Q not in scattering plane',
|
||
+ ' Q modulus too small',
|
||
+ ' KI,KF,Q triangle cannot close',
|
||
C
|
||
+ ' Error in KI or KF. Check D-spacings & units',
|
||
+ ' KI or KF cannot be obtained',
|
||
+ ' KI or KF too small',
|
||
+ ' KI,KF,Q triangle will not close'/
|
||
C---------------------------------------------------------------------------
|
||
LIER = MIN(MAX(IER,1),MXER)
|
||
LMODULE = MIN(MAX(MODULE,1),MXMOD)
|
||
C WRITE(iolun,501) MESER(LIER,LMODULE)
|
||
501 FORMAT(A)
|
||
RETURN
|
||
END
|