897616f073
* Update existing 2020.0 version to use tag * Add versions 2018.2 and master * Add patches for GCC/Intel * Use MPI compiler wrappers when +mpi * Constrain CMake build dependency (need >= 3.1) * Add variants for optional components (e.g QFIT library)
94 lines
3.3 KiB
Diff
94 lines
3.3 KiB
Diff
--- a/DALTON/soppa/so_stoppw.F
|
|
+++ b/DALTON/soppa/so_stoppw.F
|
|
@@ -1,19 +1,22 @@
|
|
C
|
|
-C /* Deck so_stoppw */
|
|
- SUBROUTINE SO_STOPPW(STOPP,TRGOS,ISYMTR,IEXCI,EXENG,QVAL)
|
|
+ SUBROUTINE SO_STOPPW(RSTOPP,TRGOS,ISYMTR,IEXCI,EXENG,QVAL)
|
|
C
|
|
C This routine is part of the atomic integral directSOPPA program.
|
|
-C
|
|
+C The charge Z of the incoming ions is set to 1 here.
|
|
C Zhiwen Shi (Clark), Stephan P. A. Sauer, January 2016
|
|
-C
|
|
C PURPOSE: Calculate Stopping Power.
|
|
C
|
|
-#include "implicit.h"
|
|
-#include "cbiexc.h"
|
|
-#include "ccorb.h"
|
|
+ implicit none
|
|
+#include "cbiexc.h" ! LVEL, MXNEXI
|
|
+#include "ccorb.h" ! NSYM
|
|
+#include "pi.h" ! PI
|
|
C
|
|
- DIMENSION STOPP(3,LVEL),TRGOS(3),EXENG(NSYM,MXNEXI)
|
|
- REAL*8 QVAL,VELOC,QMAXV,QMINV
|
|
+ REAL*8, INTENT(INOUT) :: RSTOPP(3,LVEL,2)
|
|
+ REAL*8, INTENT(IN) :: TRGOS(3), EXENG(NSYM,MXNEXI), QVAL
|
|
+ INTEGER, INTENT(IN) :: ISYMTR, IEXCI
|
|
+ REAL*8 :: VELOC, QMAXV, QMINV
|
|
+ REAL*8, PARAMETER :: D4 = 4.0D0
|
|
+ INTEGER :: IVEL
|
|
C
|
|
C------------------
|
|
C Add to trace.
|
|
@@ -28,19 +31,52 @@ C
|
|
DO IVEL=1, LVEL
|
|
C
|
|
VELOC = VMIN+(IVEL-1)*VSTEP
|
|
-C
|
|
QMAXV = VELOC*2
|
|
QMINV = EXENG(ISYMTR,IEXCI)/VELOC
|
|
C
|
|
IF (QMINV .LE. QMAXV) THEN
|
|
C
|
|
- IF ((QVAL .GE. QMINV) .AND. (QVAL .LE. QMAXV)) THEN
|
|
+ IF ((QVAL .GE. QMINV) .AND. (QVAL .LE. QMAXV)) THEN
|
|
+C
|
|
+ RSTOPP(1,IVEL,1) = RSTOPP(1,IVEL,1) + TRGOS(1)*QSTEP*
|
|
+ & PI*D4/(QVAL*VELOC*VELOC)
|
|
+ RSTOPP(2,IVEL,1) = RSTOPP(2,IVEL,1) + TRGOS(2)*QSTEP*
|
|
+ & PI*D4/(QVAL*VELOC*VELOC)
|
|
+ RSTOPP(3,IVEL,1) = RSTOPP(3,IVEL,1) + TRGOS(3)*QSTEP*
|
|
+ & PI*D4/(QVAL*VELOC*VELOC)
|
|
+
|
|
+ ENDIF
|
|
+C
|
|
+C The next line makes sure that the split integration will be applied
|
|
+C after the velocity is larger than half of Q value.
|
|
+C This velocity as half of chosen Q value corresponds to
|
|
+C the highest excitation energy for given basis set.
|
|
+C i.e. integration can be split after this velocity.
|
|
+C
|
|
+ IF (VELOC .GE. QINP/2.0d0) THEN
|
|
+C
|
|
+ IF ((QVAL .GE. QMINV) .AND. (QVAL .LE. QINP)
|
|
+ & .AND. (QVAL .LE. QMAXV)) THEN
|
|
+C
|
|
+ RSTOPP(1,IVEL,2) = RSTOPP(1,IVEL,2) + TRGOS(1)*QSTEP*
|
|
+ & PI*D4/(QVAL*VELOC*VELOC)
|
|
+ RSTOPP(2,IVEL,2) = RSTOPP(2,IVEL,2) + TRGOS(2)*QSTEP*
|
|
+ & PI*D4/(QVAL*VELOC*VELOC)
|
|
+ RSTOPP(3,IVEL,2) = RSTOPP(3,IVEL,2) + TRGOS(3)*QSTEP*
|
|
+ & PI*D4/(QVAL*VELOC*VELOC)
|
|
+C
|
|
+ ENDIF
|
|
+C
|
|
+ ELSEIF ((QVAL .GE. QMINV) .AND. (QVAL .LE. QMAXV)) THEN
|
|
C
|
|
- STOPP(1,IVEL) = STOPP(1,IVEL) + TRGOS(1)*QSTEP
|
|
- STOPP(2,IVEL) = STOPP(2,IVEL) + TRGOS(2)*QSTEP
|
|
- STOPP(3,IVEL) = STOPP(3,IVEL) + TRGOS(3)*QSTEP
|
|
+ RSTOPP(1,IVEL,2) = RSTOPP(1,IVEL,2) + TRGOS(1)*QSTEP*
|
|
+ & PI*D4/(QVAL*VELOC*VELOC)
|
|
+ RSTOPP(2,IVEL,2) = RSTOPP(2,IVEL,2) + TRGOS(2)*QSTEP*
|
|
+ & PI*D4/(QVAL*VELOC*VELOC)
|
|
+ RSTOPP(3,IVEL,2) = RSTOPP(3,IVEL,2) + TRGOS(3)*QSTEP*
|
|
+ & PI*D4/(QVAL*VELOC*VELOC)
|
|
C
|
|
- ENDIF
|
|
+ ENDIF
|
|
C
|
|
ENDIF
|
|
C
|