File indexing completed on 2018-03-02 18:41:43 UTC
view on githubraw file Latest commit add29e06 on 2018-01-31 20:35:05 UTC
108a00eab9 Ryan*0001 #include "LAYERS_OPTIONS.h"
0002
0003
0004
0005 SUBROUTINE LAYERS_CHECK( myThid )
0006
0007
0008
0009 IMPLICIT NONE
0010 #include "SIZE.h"
0011 #include "EEPARAMS.h"
0012 #include "PARAMS.h"
0013 #include "EOS.h"
0014 #include "LAYERS_SIZE.h"
0015 #include "LAYERS.h"
0016
0017
0018 INTEGER myThid
0019
0020
c4343d1bc6 Jean*0021
108a00eab9 Ryan*0022 CHARACTER*(MAX_LEN_MBUF) msgBuf
2c2df907b2 Jean*0023 CHARACTER*(40) tmpName
0024 CHARACTER*(1) sfx
0025 INTEGER iLa, k, errCount
44cfc9294d Jean*0026 _RL tmpVar
108a00eab9 Ryan*0027
0028 #ifdef ALLOW_LAYERS
0029 _BEGIN_MASTER(myThid)
0030
0031 WRITE(msgBuf,'(A)') 'LAYERS_CHECK: #define LAYERS'
0032 CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
2c2df907b2 Jean*0033 & SQUEEZE_RIGHT, myThid )
108a00eab9 Ryan*0034
44cfc9294d Jean*0035
108a00eab9 Ryan*0036 CALL WRITE_0D_I( NZZ, INDEX_NONE, 'NZZ =',
0037 & ' /* number of levels in the fine vertical grid */')
f61c9be03a Ryan*0038 CALL WRITE_1D_RL( dZZf, NZZ, INDEX_K, 'dZZf =',
0039 & ' /* fine vertical grid spacing for isopycnal interp */')
0040
406891c1a3 Gael*0041 DO iLa=1,layers_maxNum
44cfc9294d Jean*0042 IF ( layers_num(iLa).NE.0 ) THEN
2c2df907b2 Jean*0043 sfx = '#'
0044 IF ( iLa.LE.9 ) WRITE(sfx,'(I1)') iLa
0045
0046 WRITE(tmpName,'(3A)') 'layers_num(', sfx, ') ='
0047 CALL WRITE_0D_I( layers_num(iLa), INDEX_NONE, tmpName(1:15),
0048 & ' /* averaging field: 1= theta, 2= salt, 3= prho */' )
0049 WRITE(tmpName,'(3A)') 'layers_name(', sfx, ') ='
0050 CALL WRITE_0D_C( layers_name(iLa),-1,INDEX_NONE, tmpName(1:16),
0051 & ' /* averaging field: TH = theta, SLT= salt, RHO= prho */' )
0052 WRITE(tmpName,'(3A)') 'layers_bolus(', sfx, ') ='
0053 IF ( useGMRedi )
0054 & CALL WRITE_0D_L ( layers_bolus(iLa), INDEX_NONE, tmpName(1:17),
0055 & ' /* include potential GM bolus velocity */')
0056 WRITE(tmpName,'(3A)') 'layers_krho(', sfx, ') ='
44cfc9294d Jean*0057 IF ( layers_num(iLa).EQ.3 )
2c2df907b2 Jean*0058 & CALL WRITE_0D_I( layers_krho(iLa), INDEX_NONE, tmpName(1:16),
0059 & ' /* model level to reference potential density to */' )
0060 WRITE(tmpName,'(3A)') 'layers_bounds(*,', sfx, ') ='
44cfc9294d Jean*0061 CALL WRITE_1D_RL( layers_bounds(1,iLa), Nlayers+1, INDEX_K,
2c2df907b2 Jean*0062 & tmpName(1:20), ' /* boundaries of tracer-averaging bins */')
17ce8d85dd Davi*0063
44cfc9294d Jean*0064 ENDIF
406891c1a3 Gael*0065 ENDDO
17ce8d85dd Davi*0066
44cfc9294d Jean*0067
2c2df907b2 Jean*0068 errCount = 0
44cfc9294d Jean*0069 DO iLa=1,layers_maxNum
108a00eab9 Ryan*0070
44cfc9294d Jean*0071 IF ( layers_num(iLa).NE.0 ) THEN
2c2df907b2 Jean*0072
0073
44cfc9294d Jean*0074 DO k=1,Nlayers
0075 IF ( layers_bounds(k,iLa).GE.layers_bounds(k+1,iLa) ) THEN
0076 WRITE(msgBuf,'(A,I2,A,I4)') 'LAYERS_CHECK(iLa=', iLa,
0077 & '): layers_bounds k -> k+1 not increasing at k=', k
0078 CALL PRINT_ERROR( msgBuf, myThid )
2c2df907b2 Jean*0079 errCount = errCount + 1
44cfc9294d Jean*0080 ENDIF
0081 ENDDO
0082 ENDIF
0083
0084 IF ( layers_num(iLa).EQ.3 ) THEN
0085
2c2df907b2 Jean*0086
44cfc9294d Jean*0087 tmpVar = layers_bounds(Nlayers+1,iLa) - layers_bounds(1,iLa)
0088 IF ( tmpVar.LE.50. .AND. layers_bounds(1,iLa).GE.950. ) THEN
0089 WRITE(msgBuf,'(A,I2,A)') 'LAYERS_CHECK(iLa=', iLa,
0090 & '): layers_bounds seems to be expressed as "rho"'
0091 CALL PRINT_ERROR( msgBuf, myThid )
0092 WRITE(msgBuf,'(A,I2,A)') 'LAYERS_CHECK(iLa=', iLa,
0093 & '): while it should be expressed as "rho - 1000"'
0094 CALL PRINT_ERROR( msgBuf, myThid )
2c2df907b2 Jean*0095 errCount = errCount + 1
0096 ENDIF
0097
0098 IF ( layers_krho(iLa).LT.1 .OR. layers_krho(iLa).GT.Nr ) THEN
0099 WRITE(msgBuf,'(2A,I3,A,I9)') 'LAYERS_CHECK: ',
0100 & 'Invalid layer_krho(iLa=', iLa,') =', layers_krho(iLa)
0101 CALL PRINT_ERROR( msgBuf, myThid )
0102 errCount = errCount + 1
44cfc9294d Jean*0103 ENDIF
0104 ENDIF
0105
0106 ENDDO
108a00eab9 Ryan*0107
2c2df907b2 Jean*0108 IF ( errCount.GE.1 ) THEN
0109 WRITE(msgBuf,'(A,I3,A)')
0110 & 'LAYERS_CHECK: detected', errCount,' fatal error(s)'
0111 CALL PRINT_ERROR( msgBuf, myThid )
0112 CALL ALL_PROC_DIE( 0 )
0113 STOP 'ABNORMAL END: S/R LAYERS_CHECK'
0114 ELSE
0115 WRITE(msgBuf,'(A)') 'LAYERS_CHECK: done'
0116 CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
0117 & SQUEEZE_RIGHT, myThid )
0118 ENDIF
0119
108a00eab9 Ryan*0120 _END_MASTER(myThid)
c4343d1bc6 Jean*0121 #endif /* ALLOW_LAYERS */
108a00eab9 Ryan*0122
0123 RETURN
0124 END