Merge branch '8-improve-epsilon-tests-and-error-detection-in-module-salt' into 'main'

#8: Improve epsilon tests and error detection in module SALT:

See merge request dragon/5.1!14

6 files changed, 54 insertions, 7 deletions

diff --git a/Dragon/src/EVOSIG.f b/Dragon/src/EVOSIG.f
index 3e5c172..a9c282e 100644
--- a/Dragon/src/EVOSIG.f
+++ b/Dragon/src/EVOSIG.f
@@ -198,7 +198,7 @@
         SIG(IS,NREAC,IBM)=SIG(IS,NREAC,IBM)+1.0E-3*FACT*REAL(GAR)
         DEALLOCATE(ZKERMA)
       ELSE
-        IF(IGLOB.EQ.-1) THEN
+        IF((IGLOB.EQ.-1).AND.(AWR(IS).GT.210.0)) THEN
           CALL XABORT('EVOSIG: EDP0 OPTION NEEDS H-FACTOR INFORMATION.')
         ENDIF
       ENDIF
diff --git a/Dragon/src/MCCGT.f b/Dragon/src/MCCGT.f
index cb890d9..458a1f4 100644
--- a/Dragon/src/MCCGT.f
+++ b/Dragon/src/MCCGT.f
@@ -158,19 +158,21 @@
       ENDDO
       READ(IFTRAK) NDIM,ISPEC,N2REG,N2SOU,NALBG,NCOR,NANGL,MXSUB,MXSEG
       IF((NDIM.NE.2).AND.(NDIM.NE.3)) 
-     & CALL XABORT('2D OR 3D EXCELT TRACKING EXPECTED')
+     1 CALL XABORT('2D OR 3D EXCELT TRACKING EXPECTED')
 *----
 *  RECOVER TRACKING STATE-VECTOR AND USER INPUT INFORMATION
 *----
       IGP(:NSTATE)=0 
       CALL LCMGET(IPTRK,'STATE-VECTOR',IGP)
-      CALL LCMGET(IPTRK,'ALBEDO',ALBEDO)
+      IF(IGP(7).EQ.5) CALL XABORT('MCCGT: THE SHORT CHARACTERISTICS M'
+     1 //'ETHOD IS NOT IMPLEMENTED. USE NOIC OPTION IN SALT:.')
       NREG=IGP(1)
       NSOU=IGP(5)
       NANIS=IGP(6)
       TRTY=IGP(9)
       CYCLIC=(TRTY.EQ.1)
       ISYMM=IGP(12)
+      CALL LCMGET(IPTRK,'ALBEDO',ALBEDO)
 *
       IMPX=1
       LCACT=IGP(13)
diff --git a/Dragon/src/MUSACG.f90 b/Dragon/src/MUSACG.f90
index 4f68a52..9ceaeb8 100644
--- a/Dragon/src/MUSACG.f90
+++ b/Dragon/src/MUSACG.f90
@@ -57,6 +57,7 @@ SUBROUTINE MUSACG(ITRACK,IFTRK,IPRINT,IMACRO,NBSLIN,RCUTOF,GG,LGINF,NBNODE_MACRO
   INTEGER ELEM, OK, TYPE
   REAL(PDB) :: X1,X2,Y1,Y2,DET1,DET2
   REAL(PDB) :: DGMESHX(2),DGMESHY(2)
+  LOGICAL :: LTEST
   INTEGER, DIMENSION(NSTATE) :: I_STATE,IEDIMG
   CHARACTER(LEN=72) :: TEXT72
   !----
@@ -313,6 +314,11 @@ SUBROUTINE MUSACG(ITRACK,IFTRK,IPRINT,IMACRO,NBSLIN,RCUTOF,GG,LGINF,NBNODE_MACRO
     !----
     !  Unfold macro geometry (many times, if required)
     !----
+    LTEST=.TRUE.
+    DO IB=1,GG_MAC%NBBCDA
+      IF(.NOT.LFOLD(IB)) LTEST=.FALSE.
+    ENDDO
+    IF(LTEST) CALL XABORT('MUSACG: YOU CANNOT UNFOLD ALL PERIMETERS OF A MACROCELL.')
     DO IB=1,GG_MAC%NBBCDA
       IF(LFOLD(IB)) THEN
         ALLOCATE(IFOLD(2*GG_MAC%NB_ELEM))
diff --git a/Dragon/src/SALTLS.f90 b/Dragon/src/SALTLS.f90
index 32a4500..2ce2c4f 100644
--- a/Dragon/src/SALTLS.f90
+++ b/Dragon/src/SALTLS.f90
@@ -77,6 +77,7 @@ SUBROUTINE SALTLS(IFTEMP,IPRINT,IGTRK,NFREG,NBANGL,NQUAD,RENO,NBDR,IFMT,DENUSR,
   REAL(PDB), DIMENSION(2) :: THETA0
   REAL(PDB), ALLOCATABLE, DIMENSION(:)  :: VOLN,SURFN,CURRN
   REAL(PDB), ALLOCATABLE, DIMENSION(:,:) :: FACNRM ! aux for normalisation
+  CHARACTER(LEN=131) :: HSMG
   REAL, PARAMETER :: EPS3  = 1E-3
   INTEGER, PARAMETER :: FOUT =6
   !
@@ -174,7 +175,10 @@ SUBROUTINE SALTLS(IFTEMP,IPRINT,IGTRK,NFREG,NBANGL,NQUAD,RENO,NBDR,IFMT,DENUSR,
                  II0=CNT0+1
                  IF(IPRINT > 3) CALL SAL231(RTRAC2(II0:),ITRAC2(II0:),DELX,EX0,EY0,ANGLE)
                  DO II=1,ITRAC2(II0)
-                   IF(RTRAC2(II0+II+NNN-1) <= 0.0) CALL XABORT('SALTLS: INVALID SEGMENT LENGTH')
+                   IF(RTRAC2(II0+II+NNN-1) <= -EPSILON_PDB) THEN
+                     WRITE(HSMG,'(39HSALTLS: INVALID SEGMENT LENGTH (RTRAC2=,1P,E12.4,1H))') RTRAC2(II0+II+NNN-1)
+                     CALL XABORT(HSMG)
+                   ENDIF
                  ENDDO
                  ! compute volumes
                  CALL SAL232(ITRAC2(II0:),RTRAC2(II0:),FACNRM,GG,SURFN,CURRN)
diff --git a/Dragon/src/SAL_GEOMETRY_MOD.f90 b/Dragon/src/SAL_GEOMETRY_MOD.f90
index 1056016..4bce6e8 100644
--- a/Dragon/src/SAL_GEOMETRY_MOD.f90
+++ b/Dragon/src/SAL_GEOMETRY_MOD.f90
@@ -16,7 +16,7 @@ MODULE SAL_GEOMETRY_MOD
 
   USE SAL_GEOMETRY_TYPES
   USE PRECISION_AND_KINDS, ONLY : PDB, PI,TWOPI,HALFPI
-  USE SAL_NUMERIC_MOD,    ONLY : SAL141
+  USE SAL_NUMERIC_MOD, ONLY : SAL141
   USE SALGET_FUNS_MOD
 
 CONTAINS
@@ -341,7 +341,7 @@ CONTAINS
     !
     !---------------------------------------------------------------------
     !
-    USE SAL_GEOMETRY_TYPES, ONLY : ALLSUR,PREC
+    USE SAL_GEOMETRY_TYPES, ONLY : TYPGEO,ALLSUR,PREC
     !****
     IMPLICIT NONE
     ! in variable
@@ -351,6 +351,7 @@ CONTAINS
     INTEGER, PARAMETER :: N_DATAIN=25
     INTEGER, DIMENSION (N_DATAIN) :: DATAIN
     INTEGER   :: ELEM,I,TYPE,NBER
+    CHARACTER(LEN=131) :: HSMG
     INTEGER, PARAMETER, DIMENSION(0:4) :: READ_BC_LEN=(/1,1,3,3,3/)
     INTEGER, PARAMETER :: FOUT =6
     REAL(PDB) :: ANGLE,BCDATA_TDT(3)
@@ -409,9 +410,43 @@ CONTAINS
             GG%BCDATAREAD(I)%BCDATA(6)=1._PDB
             LBCDIAG=LBCDIAG.OR.((GG%BCDATAREAD(I)%BCDATA(1)==0._PDB).AND.(GG%BCDATAREAD(I)%BCDATA(2)==0._PDB) &
                            .AND.(GG%BCDATAREAD(I)%BCDATA(5)==PI/4._PDB))
+            IF((TYPGEO.EQ.3).OR.(TYPGEO.EQ.6).OR.(TYPGEO.EQ.6)) THEN
+              IF((ANGLE/=0.0).AND.(ANGLE/=PI/2._PDB)) THEN
+                WRITE(HSMG,100) ANGLE,I,TYPGEO
+                CALL XABORT(HSMG)
+              ENDIF
+            ELSE IF(TYPGEO.EQ.7) THEN
+              IF((ANGLE/=0.0).AND.(ANGLE/=PI/2._PDB).AND.(ANGLE/=PI/4._PDB)) THEN
+                WRITE(HSMG,100) ANGLE,I,TYPGEO
+                CALL XABORT(HSMG)
+              ENDIF
+            ELSE IF((TYPGEO.EQ.8).OR.(TYPGEO.EQ.10)) THEN
+              IF((ANGLE/=0.0).AND.(ANGLE/=PI/3._PDB).AND.(ANGLE/=2._PDB*PI/3._PDB)) THEN
+                WRITE(HSMG,100) ANGLE,I,TYPGEO
+                CALL XABORT(HSMG)
+              ENDIF
+            ELSE IF(TYPGEO.EQ.9) THEN
+              IF((ANGLE/=0.0).AND.(ANGLE/=PI/6._PDB).AND.(ANGLE/=2._PDB*PI/3._PDB).AND. &
+              & (ANGLE/=PI).AND.(ANGLE/=7._PDB*PI/6._PDB).AND.(ANGLE/=5._PDB*PI/3._PDB)) THEN
+                WRITE(HSMG,100) ANGLE,I,TYPGEO
+                CALL XABORT(HSMG)
+              ENDIF
+            ELSE IF(TYPGEO.EQ.11) THEN
+              IF((ANGLE/=0.0).AND.(ANGLE/=PI/3._PDB)) THEN
+                WRITE(HSMG,100) ANGLE,I,TYPGEO
+                CALL XABORT(HSMG)
+              ENDIF
+            ELSE IF(TYPGEO.EQ.12) THEN
+              IF((ANGLE/=0.0).AND.(ANGLE/=PI/6._PDB).AND.(ANGLE/=2._PDB*PI/3._PDB)) THEN
+                WRITE(HSMG,100) ANGLE,I,TYPGEO
+                CALL XABORT(HSMG)
+              ENDIF
+            ENDIF
           ENDIF
        ENDDO
     ENDIF
+    100 FORMAT(34HSALINP: FORBIDDEN PERIMETER ANGLE=,1P,E13.4,18H RADIANS FOR SIDE=,I3, &
+    & 12H AND TYPGEO=,I3,1H.)
   END SUBROUTINE SALINP
   !
   SUBROUTINE SAL126(RPAR,IPAR)
diff --git a/Dragon/src/SAL_TRAJECTORY_MOD.f90 b/Dragon/src/SAL_TRAJECTORY_MOD.f90
index e108fcb..ac1714d 100644
--- a/Dragon/src/SAL_TRAJECTORY_MOD.f90
+++ b/Dragon/src/SAL_TRAJECTORY_MOD.f90
@@ -110,7 +110,7 @@ CONTAINS
       ELSE
         WRITE(*,*) 'PPERIM_MAC2(N_AXIS+1),PPERIM_MAC2(N_AXIS) :',PPERIM_MAC2(N_AXIS+1),PPERIM_MAC2(N_AXIS)
         WRITE(*,*) 'DIST_AXIS(PPERIM_MAC2(N_AXIS+1)-1) :',DIST_AXIS(PPERIM_MAC2(N_AXIS+1)-1)
-        WRITE(*,*) 'DELX :',DELX
+        WRITE(*,*) 'DELX :',DELX,' RADIA=',RADIA
         CALL XABORT('SALTRA: Cant find entry point')
       ENDIF
     ENDIF