*ID RSPMSL */ New stashcode for new pmsl calculation */ *// *//Built by R Stratton (hadas) on Tue Feb 16 09:39:58 GMT 1999 *// *//-------------- *DECLARE ST_DIA21 *//-------------- *I GDR4F405.11 !LL 4.6 16/02/99 New pmsl calculation 16255 R A Stratton *D ASW3F304.8 & PT220,PT221,PT222,PT223,PT224,PT225,PT255 *I GDR4F305.257 PT255=SI(255,16,im_index) *D ST_DIA21.278 & SEC_U_LATITUDE,COS_P_LATITUDE, *D ADP0F401.121 & STASHWORK,INT16,PT_TRACER,STASHWORK(PT255), *D ADP0F401.123 & SF_TRACER,SF(255,16), *// *//Built by R Stratton (hadas) on Tue Feb 16 09:40:01 GMT 1999 *// *//-------------- *DECLARE PHYDIA1A *//-------------- *I GDR4F405.17 !LL 4.6 16/02/99 Call new PMSL diagnostic 16255 R A Stratton *I PHYDIA1A.32 & COS_P_LATITUDE, *D ADP0F401.11 & TRACERWORK,INT16,PT_TRACER,PMSL_NEW, *D ADP0F401.13 & SF_TRACER, QPMSL_NEW, *I GDR4F405.22 *,QPMSL_NEW !IN Logical flag for new pmsl diagnostic *I PHYDIA1A.125 *,PMSL_NEW(P_FIELD) ! new pmsl diagnostic *I PHYDIA1A.160 *,COS_P_LATITUDE(P_FIELD) !IN COS(LAT) AT P POINTS *I GSM1F405.760 & ,ROW_LENGTH,P_ROWS *D GSM1F405.762 & ,LAST_P_FLD_PT,EW_SPACE,NS_SPACE,COS_P_LATITUDE,.false. & ,ICODE,CMESSAGE) ENDIF ! 255 new method of PMSL calculation IF(QPMSL_NEW) THEN CALL PMSL(PMSL_NEW,P(1,Z_REF),PSTAR,P_EXNER_HALF,THETA,Q & ,PHI_STAR,ROW_LENGTH,P_ROWS & ,P_FIELD,P_LEVELS,Q_LEVELS,Z_REF,AKH,BKH,FIRST_FLD_PT & ,LAST_P_FLD_PT,EW_SPACE,NS_SPACE,COS_P_LATITUDE,.true. & ,ICODE,CMESSAGE) *// *//Built by R Stratton (hadas) on Tue Feb 16 09:40:03 GMT 1999 *// *//-------------- *DECLARE ST_MEAN1 *//-------------- *I GDR4F405.7 !LL 4.6 16/02/99 New PMSL diagnostic in call to phy_diag !LL R A Stratton. *D ST_MEAN1.424 & EW_SPACE,NS_SPACE,SEC_U_LATITUDE,COS_P_LATITUDE, *I ADP0F401.141 & STASHWORK , *D ADP0F401.143 & .FALSE.,.FALSE.,SF_TRACER,.FALSE., *// *//Built by R Stratton (hadas) on Tue Feb 16 09:40:05 GMT 1999 *// *//------------ *DECLARE PMSL1A *//------------ *I GDR8F405.12 !LL 4.6 16/02/98 New PMSL calculation R A Stratton *I GSM1F405.558 * ,ROW_LENGTH,P_ROWS *D GSM1F405.560 & ,START,END & ,EW_SPACE,NS_SPACE,COS_P_LATITUDE, L_PMSL_NEW & ,ICODE,CMESSAGE ) *D GSS5F402.18 * ROW_LENGTH !IN *,P_ROWS !IN *,POINTS !IN Number of points per level *I GSS5F402.20 *, ICODE !OUT RETURN CODE : IRET=0 NORMAL EXIT *I PMSL1A.47 *,EW_SPACE !IN LONGITUDE GRID SPACING *,NS_SPACE !IN LATITUDE GRID SPACING *,COS_P_LATITUDE(POINTS) !IN Cosine of latitudes on pressure grid *I PMSL1A.48 CHARACTER*80 CMESSAGE LOGICAL L_PMSL_NEW ! IN ! PARAMETER (L_PMSL_NEW=.TRUE.) now passed in *D PMSL1A.60 REAL TS(points) ! Surface Temperature *I PMSL1A.61 *IF DEF,MPP *CALL PARVARS INTEGER ROW_LENGTH_GLOBAL INTEGER P_ROWS_GLOBAL INTEGER POINTS_GLOBAL *ENDIF *D GSS5F402.38,GSS5F402.39 TS(I)=THETA(I,L)*P_EXNER_FULL*TEMP(I) TS(I)=TS(I)*(1.0+C_VIRTUAL*Q(I,1)) TEMP(I)=(TS(I)+LAPSE*ONE_OVER_G*PHI_STAR(I))/TS(I) *I GSS5F402.54 *IF DEF,MPP *I PMSL1A.105 CAV*IF -DEF,T3E CAV ICODE=1 CAV CMESSAGE='PMSL: MPP version for T3E only using SHMEM' CAV*ENDIF ROW_LENGTH_GLOBAL=GLSIZE(1) P_ROWS_GLOBAL=GLSIZE(2) POINTS_GLOBAL=ROW_LENGTH_GLOBAL*P_ROWS_GLOBAL *ENDIF IF (L_PMSL_NEW) THEN CALL NPMSL(P_MSL,PSTAR *, PHI_STAR,THETA,TS *, COS_P_LATITUDE,EW_SPACE,NS_SPACE *, ROW_LENGTH,P_ROWS *, POINTS,P_LEVELS *IF DEF,MPP *, ROW_LENGTH_GLOBAL,P_ROWS_GLOBAL,POINTS_GLOBAL *ENDIF * ) ENDIF *I PMSL1A.108 *IF -DEF,MPP SUBROUTINE NPMSL(P_MSL,PSTAR *, PHI_STAR,THETA,TS *, COS_P_LATITUDE,EW_SPACE_DEG,NS_SPACE_DEG *, ROW_LENGTH,P_ROWS *, POINTS,P_LEVELS * ) *ELSE SUBROUTINE NPMSL(P_MSL_LOCAL,PSTAR_LOCAL *, PHI_STAR_LOCAL,THETA_LOCAL,TS_LOCAL *, COS_P_LATITUDE_LOCAL,EW_SPACE_DEG,NS_SPACE_DEG *, ROW_LENGTH_LOCAL,P_ROWS_LOCAL *, POINTS_LOCAL,P_LEVELS *, ROW_LENGTH,P_ROWS,POINTS) IMPLICIT NONE C C DEFINE GLOBAL VARIABLES C INTEGER * ROW_LENGTH_LOCAL !IN NUMBER OF POINTS PER ROW *,P_ROWS_LOCAL !IN NUMBER OF ROWS *,POINTS_LOCAL !IN NUMBER OF POINTS PER LEVEL REAL * P_MSL_LOCAL(POINTS_LOCAL) !IN/OUT MEAN SEA LEVEL PRESSURE *,PSTAR_LOCAL(POINTS_LOCAL) !IN SURFACE PRESSURE *,PHI_STAR_LOCAL(POINTS_LOCAL) !IN SURFACE GEOPOTENTIAL *,THETA_LOCAL(POINTS_LOCAL) !IN POTENTIAL TEMPERATURE ON MODEL LEVELS *,TS_LOCAL(POINTS_LOCAL) !IN SURFACE TEMPERATURE *,COS_P_LATITUDE_LOCAL(POINTS_LOCAL) !IN COSINE OF LATITUDES ON PRESSURE GRID *,EXNERT_LOCAL(POINTS_LOCAL) ! NEW EXNER PRESSURE C C DEFINE PE0'S VARIABLES C *ENDIF INTEGER * ROW_LENGTH !IN NUMBER OF POINTS PER ROW *,P_ROWS !IN NUMBER OF ROWS *,POINTS !IN NUMBER OF POINTS PER LEVEL *,P_LEVELS !IN NUMBER OF MODEL LEVELS REAL * P_MSL(POINTS) !IN/OUT MEAN SEA LEVEL PRESSURE *,PSTAR(POINTS) !IN SURFACE PRESSURE *,PHI_STAR(POINTS) !IN SURFACE GEOPOTENTIAL *,THETA(POINTS) !IN POTENTIAL TEMPERATURE ON MODEL LEVELS *,TS(POINTS) !IN SURFACE TEMPERATURE *,COS_P_LATITUDE(POINTS) !IN COSINE OF LATITUDES ON PRESSURE G *,EW_SPACE_DEG !IN EAST-WEST GRIDSPACE IN DEGREES *,NS_SPACE_DEG !IN NORTH-SOUTH GRIDSPACE IN DEGREES C C DEFINE LOCAL VARIABLES C REAL * EXNERS(POINTS) ! EXNER PRESSURE CALCULATED FROM SURFACE PRESSURE *,EXNERM(POINTS) ! EXNER PRESSURE CALCUALTED FROM MSLP *,EXNERT(POINTS,2) ! NEW EXNER PRESSURE *,AVEXNER(POINTS) ! EXNER PRESSURE REQUIRED IN OVERRELAXTION EQN() *,OROG(POINTS) ! LAND SURFACE OROGRAPHY *,TEMPS ! TEMPERATURE AT MEAN SEA LEVEL *,THETAS(POINTS) ! POTENTIAL TEMPERATURE AT MEAN SEA LEVEL *,FUG(POINTS) ! CORIOLIS PARAMETER TIMES GEOSTROPHIC U WIND *,FVG(POINTS) ! CORIOLIS PARAMETER TIMES GEOSTROPHIC V WIND *,F(POINTS) ! FORCING FUNCTION RHS OF EQN() *,SPHLO1(POINTS) ! LONGITUNDINAL SPHERICAL TERM *,SPHLO2(POINTS) ! LONGITUNDINAL SPHERICAL TERM SQUARED *,SPHLA1 ! LATITUNDINAL SPHERICAL TERM *,SPHLA2 ! LATITUNDINAL SPHERICAL TERM SQUARED *,CHGT ! OROGRAPHIC HEIGHT ABOVE WHICH RELAXATION OCCURS *,NITER ! NUMBER OF ITERATIONS *,ORR ! OVER RELAXATION FACTOR *,OORR ! OVER RELAXATION FACTOR AT PREVIOUS ITERATION *,RHO ! SPHERICAL RADIUS OF CONVERGENCE *,N ! TYPICAL NUMBER OF POINTS ABOVE CHGT *,EW_SPACE ! EAST-WEST GRIDSPACE IN RADIANS *,NS_SPACE ! NORTH-SOUTH GRIDSPACE IN RADIANS REAL * TRIGS(ROW_LENGTH) ! HOLDS TRIGONOMETRIC FUNCTIONS NEEDED BY * ! FOURIER, FFT99 AND FFT991 INTEGER * IFAX(10) !OUT HOLDS FACTORS OF ROW_LENGTH *, FILTER_WAVE_NUMBER(P_ROWS) & ! LAST WAVE NUMBER ON EACH ROW WHICH IS NOT FILTERED *, NORTHERN_FILTERED_ROW !IN LAST ROW, MOVING EQUATORWARDS, * ! IN NORTHERN HEMISPHERE * ! ON WHICH FILTERING IS PERFORMED. *, SOUTHERN_FILTERED_ROW !IN LAST ROW, MOVING EQUATORWARDS, *, FILTER_SPACE C C WORKSPACE USAGE C REAL AA(POINTS),BB(POINTS),PRE(POINTS),FAC INTEGER INDEX1,INDEX2 C *IF DEF,MPP *CALL PARVARS *ENDIF C VARIABLES REQUIRED FOR LOOPS C INTEGER I,JJ,II,KK,START,END,MM C C PARAMETERS AND CONSTANTS C REAL * ONE_OVER_G ! *CALL C_G *CALL C_A *CALL C_R_CP *CALL C_LAPSE *CALL C_PI PARAMETER(ONE_OVER_G=1./G,N=15.) REAL FILTER_LAT_LIMIT REAL COS_FILTER_LAT_LIMIT REAL MIN_RETAINED_WAVE *IF DEF,MPP *IF DEF,T3E INTEGER PE_DATASTART,J INTEGER ADD,REM_ADD COMMON/ADDRESSES/ ADD,REM_ADD REAL PE0_ARRAY(POINTS) POINTER(PTR,PE0_ARRAY) *ELSE INTEGER info *ENDIF INTEGER START_LOCAL,END_LOCAL *ENDIF *IF DEF,T3E REAL POWER *ENDIF REAL P_MSL_AVERAGE C C CONSTANTS C START=1 END=POINTS RHO=COS(PI/N) CHGT=500. NITER=20 *IF DEF,GLOBAL FILTER_LAT_LIMIT=60 MIN_RETAINED_WAVE=250000 *ENDIF *IF DEF,MPP *IF DEF,T3E C C SET UP GLOBAL DIMENSIONS C PE_DATASTART=(DATASTART(2)-1)*ROW_LENGTH+ 1 DATASTART(1)-1 C C C COLLECT FROM PSTAR_LOCAL TO PSTAR C ADD=LOC(PSTAR) CALL BARRIER() CALL SHMEM_GET(REM_ADD,ADD,1,0) PTR=REM_ADD CALL BARRIER() DO J=2,P_ROWS_LOCAL-1 CALL SHMEM_PUT( 1 PE0_ARRAY(1+PE_DATASTART+(J-2)*ROW_LENGTH), 1 PSTAR_LOCAL(2+(J-1)*ROW_LENGTH_LOCAL),ROW_LENGTH_LOCAL-2,0) ENDDO CALL BARRIER() ADD=LOC(P_MSL) CALL BARRIER() CALL SHMEM_GET(REM_ADD,ADD,1,0) PTR=REM_ADD CALL BARRIER() DO J=2,P_ROWS_LOCAL-1 CALL SHMEM_PUT( 1 PE0_ARRAY(1+PE_DATASTART+(J-2)*ROW_LENGTH), 1 P_MSL_LOCAL(2+(J-1)*ROW_LENGTH_LOCAL),ROW_LENGTH_LOCAL-2,0) ENDDO CALL BARRIER() ADD=LOC(THETA) CALL BARRIER() CALL SHMEM_GET(REM_ADD,ADD,1,0) PTR=REM_ADD CALL BARRIER() DO J=2,P_ROWS_LOCAL-1 CALL SHMEM_PUT( 1 PE0_ARRAY(1+PE_DATASTART+(J-2)*ROW_LENGTH), 1 THETA_LOCAL(2+(J-1)*ROW_LENGTH_LOCAL),ROW_LENGTH_LOCAL-2,0) ENDDO CALL BARRIER() ADD=LOC(TS) CALL BARRIER() CALL SHMEM_GET(REM_ADD,ADD,1,0) PTR=REM_ADD CALL BARRIER() DO J=2,P_ROWS_LOCAL-1 CALL SHMEM_PUT( 1 PE0_ARRAY(1+PE_DATASTART+(J-2)*ROW_LENGTH), 1 TS_LOCAL(2+(J-1)*ROW_LENGTH_LOCAL),ROW_LENGTH_LOCAL-2,0) ENDDO CALL BARRIER() ADD=LOC(PHI_STAR) CALL BARRIER() CALL SHMEM_GET(REM_ADD,ADD,1,0) PTR=REM_ADD CALL BARRIER() DO J=2,P_ROWS_LOCAL-1 CALL SHMEM_PUT( 1 PE0_ARRAY(1+PE_DATASTART+(J-2)*ROW_LENGTH), 1 PHI_STAR_LOCAL(2+(J-1)*ROW_LENGTH_LOCAL), 1 ROW_LENGTH_LOCAL-2,0) ENDDO CALL BARRIER() ADD=LOC(COS_P_LATITUDE) CALL BARRIER() CALL SHMEM_GET(REM_ADD,ADD,1,0) PTR=REM_ADD CALL BARRIER() DO J=2,P_ROWS_LOCAL-1 CALL SHMEM_PUT( 1 PE0_ARRAY(1+PE_DATASTART+(J-2)*ROW_LENGTH), 1 COS_P_LATITUDE_LOCAL(2+(J-1)*ROW_LENGTH_LOCAL), 1 ROW_LENGTH_LOCAL-2,0) ENDDO CALL BARRIER() *ELSE CALL GATHER_FIELD(PSTAR_LOCAL,PSTAR, & ROW_LENGTH_LOCAL,P_ROWS_LOCAL, & ROW_LENGTH,P_ROWS, & 0,gc_all_proc_group,info) CALL GATHER_FIELD(P_MSL_LOCAL,P_MSL, & ROW_LENGTH_LOCAL,P_ROWS_LOCAL, & ROW_LENGTH,P_ROWS, & 0,gc_all_proc_group,info) CALL GATHER_FIELD(THETA_LOCAL,THETA, & ROW_LENGTH_LOCAL,P_ROWS_LOCAL, & ROW_LENGTH,P_ROWS, & 0,gc_all_proc_group,info) CALL GATHER_FIELD(TS_LOCAL,TS, & ROW_LENGTH_LOCAL,P_ROWS_LOCAL, & ROW_LENGTH,P_ROWS, & 0,gc_all_proc_group,info) CALL GATHER_FIELD(PHI_STAR_LOCAL,PHI_STAR, & ROW_LENGTH_LOCAL,P_ROWS_LOCAL, & ROW_LENGTH,P_ROWS, & 0,gc_all_proc_group,info) *ENDIF START_LOCAL=1 END_LOCAL=POINTS_LOCAL IF(MYPE.EQ.0)THEN *ENDIF EW_SPACE=PI_OVER_180*EW_SPACE_DEG NS_SPACE=PI_OVER_180*NS_SPACE_DEG C C SET UP SOME LOCAL ARRAYS C SPHLA1=1./(NS_SPACE*A) SPHLA2=(NS_SPACE*A)**2 *IF DEF,GLOBAL C C NORTH POLE SET P_MSL AT ROWS 1 AND 2 TO AVERAGE OF ROW 2 C P_MSL_AVERAGE=0. DO II=1,ROW_LENGTH I=II+ROW_LENGTH P_MSL_AVERAGE=P_MSL_AVERAGE+P_MSL(I) ENDDO P_MSL_AVERAGE=P_MSL_AVERAGE/ROW_LENGTH DO JJ=1,2 DO II=1,ROW_LENGTH I=II+ROW_LENGTH*(JJ-1) P_MSL(I)= P_MSL_AVERAGE ENDDO ENDDO C C SOUTH POLE SET P_MSL AT ROWS P_ROWS & P_ROWS-1 C TO AVERAGE OF ROW P_ROWS-1 C P_MSL_AVERAGE=0. DO II=1,ROW_LENGTH I=II+ROW_LENGTH*(P_ROWS-2) P_MSL_AVERAGE=P_MSL_AVERAGE+P_MSL(I) ENDDO P_MSL_AVERAGE=P_MSL_AVERAGE/ROW_LENGTH DO JJ=P_ROWS-1,P_ROWS DO II=1,ROW_LENGTH I=II+ROW_LENGTH*(JJ-1) P_MSL(I)= P_MSL_AVERAGE ENDDO ENDDO *ENDIF C Anette *IF DEF,VECTLIB DO I=START,END EXNERS(I)=(PSTAR(I)/100000.) EXNERM(I)=(P_MSL(I)/100000.) ENDDO CALL POWR_V(END-START+1,EXNERS(START),KAPPA,EXNERS(START)) CALL POWR_V(END-START+1,EXNERM(START),KAPPA,EXNERM(START)) *ELSE DO I=START,END EXNERS(I)=(PSTAR(I)/100000.)**KAPPA EXNERM(I)=(P_MSL(I)/100000.)**KAPPA ENDDO *ENDIF DO I=START,END OROG(I)=ONE_OVER_G*PHI_STAR(I) TEMPS =TS(I)+LAPSE*OROG(I) THETAS(I)=TEMPS /EXNERM(I) SPHLO1(I)=1./(EW_SPACE*A*COS_P_LATITUDE(I)) SPHLO2(I)=(EW_SPACE*A*COS_P_LATITUDE(I))**2. PRE(I)=(SPHLO2(I)*SPHLA2)/((2*SPHLO2(I))+(2*SPHLA2)) EXNERT(I,1)=EXNERM(I) EXNERT(I,2)=EXNERM(I) AVEXNER(I)=EXNERM(I) ENDDO C C CALCULATE GEOSTROPHIC WIND C DO JJ=1,P_ROWS DO II=1,ROW_LENGTH I=II+ROW_LENGTH*(JJ-1) IF (II.EQ.1) THEN INDEX1=I+1 INDEX2=JJ*ROW_LENGTH ELSEIF(II.EQ.ROW_LENGTH) THEN INDEX1=1+((JJ-1)*ROW_LENGTH) INDEX2=I-1 ELSE INDEX1=(I+1) INDEX2=(I-1) ENDIF FVG(I)=0.5*SPHLO1(I)*((PHI_STAR(INDEX1)-PHI_STAR(INDEX2)) + +CP*THETA(I)*(EXNERS(INDEX1)-EXNERS(INDEX2))) IF (JJ.EQ.1)THEN INDEX1=I+ROW_LENGTH INDEX2=I FAC=1.0 ELSEIF (JJ.EQ.P_ROWS)THEN INDEX1=I INDEX2=I-ROW_LENGTH FAC=1.0 ELSE INDEX1=I+ROW_LENGTH INDEX2=I-ROW_LENGTH FAC=0.5 ENDIF FUG(I)=-FAC*SPHLA1*((PHI_STAR(INDEX1)-PHI_STAR(INDEX2)) + +CP*THETA(I)*(EXNERS(INDEX1)-EXNERS(INDEX2))) ENDDO ENDDO C C CALCULATE RHS OF EQUATION F C DO I=START,END AA(I)=FVG(I)/THETAS(I) BB(I)=FUG(I)/THETAS(I) ENDDO DO JJ=2,P_ROWS-1 DO II=1,ROW_LENGTH I=II+ROW_LENGTH*(JJ-1) IF (II.EQ.1)THEN INDEX1=I+1 INDEX2=JJ*ROW_LENGTH ELSEIF (II.EQ.ROW_LENGTH)THEN INDEX1=1+((JJ-1)*ROW_LENGTH) INDEX2=I-1 ELSE INDEX1=I+1 INDEX2=I-1 ENDIF F(I)=(0.5/CP)*SPHLO1(I)*(AA(INDEX1)-AA(INDEX2)) + -(0.5/CP)*SPHLA1*(BB(I+ROW_LENGTH)-BB(I-ROW_LENGTH)) ENDDO ENDDO *IF DEF,GLOBAL C C FILTER F C CALL PMSL_SETFILT(TRIGS,IFAX,ROW_LENGTH,P_ROWS, & POINTS, COS_P_LATITUDE) COS_FILTER_LAT_LIMIT=COS(PI_OVER_180*FILTER_LAT_LIMIT) DO JJ=1,P_ROWS/2 I=JJ*ROW_LENGTH IF(COS_P_LATITUDE(I).GT.COS_FILTER_LAT_LIMIT)THEN NORTHERN_FILTERED_ROW=JJ SOUTHERN_FILTERED_ROW=P_ROWS-JJ GOTO 9876 ENDIF ENDDO 9876 CONTINUE DO JJ=1,P_ROWS I=JJ*ROW_LENGTH FILTER_WAVE_NUMBER(JJ) = ROW_LENGTH*A*EW_SPACE 1 *COS_P_LATITUDE(I)/MIN_RETAINED_WAVE ENDDO FILTER_SPACE = (ROW_LENGTH+2)*(NORTHERN_FILTERED_ROW-1+ * P_ROWS-SOUTHERN_FILTERED_ROW) CALL PMSL_FILTER (F,POINTS,FILTER_SPACE, 2 ROW_LENGTH,P_ROWS, 3 NORTHERN_FILTERED_ROW, 4 SOUTHERN_FILTERED_ROW, 5 FILTER_WAVE_NUMBER, & TRIGS,IFAX) *ENDIF C C RELAXATION TO FIND EXNER TWIDEL C DO KK=1,NITER IF(KK.EQ.1)THEN ORR=1. ELSEIF(KK.EQ.2)THEN ORR=1./(1.-((RHO**2.)/2.)) ELSE ORR=1./(1.-((OORR*RHO**2.)/4.)) ENDIF DO JJ=3,P_ROWS-2 DO II=1,ROW_LENGTH I=II+ROW_LENGTH*(JJ-1) IF(OROG(I).GT.CHGT)THEN IF(II.EQ.1)THEN INDEX1=I+1 INDEX2=JJ*ROW_LENGTH MM=1 ELSEIF(II.EQ.ROW_LENGTH)THEN INDEX1=1+((JJ-1)*ROW_LENGTH) INDEX2=I-1 MM=2 ELSE INDEX1=I+1 INDEX2=I-1 MM=2 ENDIF AVEXNER(I)=PRE(I)*((1./SPHLA2*(EXNERT(INDEX1,1) + +EXNERT(INDEX2,MM))) + +(1./SPHLO2(I)*(EXNERT(I+ROW_LENGTH,1) + +EXNERT(I-ROW_LENGTH,2)))-F(I)) EXNERT(I,2)=EXNERT(I,1)+ORR*(AVEXNER(I)-EXNERT(I,1)) ENDIF ENDDO ENDDO DO I=START,END EXNERT(I,1)=EXNERT(I,2) ENDDO OORR=ORR ENDDO *IF DEF,MPP ENDIF C C DISTRIBUTE FROM THETA_INC_GLOBAL TO THETA_INC C *IF DEF,T3E CALL BARRIER() ADD=LOC(EXNERT) CALL BARRIER() CALL SHMEM_GET(REM_ADD,ADD,1,0) PTR=REM_ADD CALL BARRIER() DO J=2,P_ROWS_LOCAL-1 CALL SHMEM_GET( 1 EXNERT_LOCAL(2+(J-1)*ROW_LENGTH_LOCAL), 1 PE0_ARRAY(1+PE_DATASTART+(J-2)*ROW_LENGTH), 1 ROW_LENGTH_LOCAL-2,0) ENDDO CALL BARRIER() *ELSE CALL SCATTER_FIELD(EXNERT_LOCAL,EXNERT, & ROW_LENGTH_LOCAL,P_ROWS_LOCAL, & ROW_LENGTH,P_ROWS, & 0,gc_all_proc_group,info) *ENDIF CALL SWAPBOUNDS 1 (EXNERT_LOCAL,ROW_LENGTH_LOCAL,P_ROWS_LOCAL,1,1,1) C C CALCULATE NEW PMSL C C Anette *ENDIF ! DO I=START_LOCAL,END_LOCAL ! P_MSL_LOCAL(I)=(EXNERT_LOCAL(I)**(1./KAPPA)) ! ENDDO C Anette *IF DEF,VECTLIB POWER=1./KAPPA CALL POWR_V(END_LOCAL-START_LOCAL+1, & EXNERT_LOCAL(START_LOCAL),POWER,EXNERT_LOCAL(START_LOCAL)) DO I=START_LOCAL,END_LOCAL P_MSL_LOCAL(I)=100000.*EXNERT_LOCAL(I) ENDDO *ELSE DO I=START_LOCAL,END_LOCAL P_MSL_LOCAL(I)=100000.*(EXNERT_LOCAL(I)**(1./KAPPA)) ENDDO *ENDIF *ELSE *IF DEF,VECTLIB POWER=1./KAPPA CALL POWR_V(END-START+1, & EXNERT(I,1),POWER,EXNERT(I,1)) DO I=START,END P_MSL(I)=100000.*EXNERT(I,1) ENDDO *ELSE DO I=START,END P_MSL(I)=100000.*(EXNERT(I,1)**(1./KAPPA)) ENDDO *ENDIF *ENDIF RETURN END *IF DEF,GLOBAL CLL SUBROUTINE PMSL_FILTER --------------------------------------- CLL CLL PURPOSE: FOURIER TRUNCATES ALL WAVES AFTER THE WAVE-NUMBER HELD IN CLL FILTER_WAVE_NUMBER FOR ALL ROWS WHERE FILTERING IS CLL NECESSARY. CALCULATED ON ONE PROCESSOR FOR PMSL CLL CALCULATION. CLL CLL MODEL MODIFICATION HISTORY FROM MODEL VERSION 4.5 CLL VERSION DATE SUBROUTINE PMSL_FILTER 1 (FIELD,FIELD_LENGTH,FILTER_SPACE, 2 ROW_LENGTH,ROWS, 3 NORTHERN_FILTERED_ROW, 4 SOUTHERN_FILTERED_ROW, 5 FILTER_WAVE_NUMBER, & TRIGS,IFAX) INTEGER * FIELD_LENGTH !IN HORIZONTAL DIMENSION OF FIELD TO BE * ! FILTERED. *, ROW_LENGTH !IN NUMBER OF POINTS ON A ROW. *, ROWS *, NORTHERN_FILTERED_ROW !IN LAST ROW, MOVING EQUATORWARDS, * ! IN NORTHERN HEMISPHERE * ! ON WHICH FILTERING IS PERFORMED. *, SOUTHERN_FILTERED_ROW !IN LAST ROW, MOVING EQUATORWARDS, * ! IN SOUTHERN HEMISPHERE * ! ON WHICH FILTERING IS PERFORMED. *, FILTER_SPACE !IN HORIZONTAL DIMENSION OF ARRAY NEEDED TO * ! HOLD DATA TO BE FILTERED TO PASS TO FFT'S. *, IFAX(10) !IN HOLDS FACTORS OF ROW_LENGTH USED IN FFT'S *, FILTER_WAVE_NUMBER(ROWS) & ! LAST WAVE NUMBER ON EACH ROW WHICH IS NOT FILTERED REAL * FIELD(FIELD_LENGTH) !INOUT HOLDS FIELD TO BE FILTERED. REAL * TRIGS(ROW_LENGTH) !IN HOLDS TRIGONOMETRIC TERMS USED IN FFT'S C DEFINE LOCAL VARIABLES INTEGER & I,J ,IJ &, LOT ! NUMBER OF DATA VECTORS PASSED TO FFT'S. &, JUMP ! NUMBER OF STORAGE LOCATIONS BETWEEN THE & ! START OF CONSECUTIVE DATA VECTOR. &, INCREMENT ! NUMBER OF STORAGE LOCATIONS BETWEEN EACH & ! ELEMENT OF THE SAME DATA VECTOR, 1, IF & ! CONSECUTIVE. &, FFTSIGN ! PARAMETER DETERMINING WHETHER SPECTRAL & ! TO GRID-POINT (FFTSIGN = 1) OR GRID-POINT & ! TO SPECTRAL (FFTSIGN = -1) FFT'S ARE & ! REQUIRED. REAL * FILTER_DATA(FILTER_SPACE) ! HOLDS DATA PASSED TO FFT'S C*L EXTERNAL SUBROUTINE CALLS:------------------------------------ EXTERNAL FOURIER C*--------------------------------------------------------------------- DO I=2,NORTHERN_FILTERED_ROW IJ = (I-2)*(ROW_LENGTH+2) IK = (I-1)*ROW_LENGTH DO J=1,ROW_LENGTH FILTER_DATA(IJ+J)=FIELD(IK+J) ENDDO ENDDO CL COLLECT DATA FROM SOUTHERN HEMISPHERE. IL = (NORTHERN_FILTERED_ROW-1)*(ROW_LENGTH+2) DO I=SOUTHERN_FILTERED_ROW,ROWS-1 IJ = IL+(I-SOUTHERN_FILTERED_ROW)*(ROW_LENGTH+2) IK = (I-1)*ROW_LENGTH DO J=1,ROW_LENGTH FILTER_DATA(IJ+J)=FIELD(IK+J) ENDDO ENDDO INCREMENT = 1 JUMP = ROW_LENGTH+2 FFTSIGN=-1 LOT = ((NORTHERN_FILTERED_ROW-1)+(ROWS-SOUTHERN_FILTERED_ROW)) CALL FOURIER(FILTER_DATA,TRIGS,IFAX,INCREMENT,JUMP,ROW_LENGTH, * LOT,FFTSIGN) DO I=2,NORTHERN_FILTERED_ROW IJ = (I-2)*(ROW_LENGTH+2) C 3 IS USED IN LOOP 320 COUNTER AS WAVE-NUMBER 0 HAS TO BE RETAINED. DO J=3+2*FILTER_WAVE_NUMBER(I),ROW_LENGTH+2 FILTER_DATA(IJ+J)=0.0 ENDDO ENDDO CL SOUTHERN HEMISPHERE DATA IL = (NORTHERN_FILTERED_ROW-1)*(ROW_LENGTH+2) DO I=SOUTHERN_FILTERED_ROW,ROWS-1 IJ = IL+(I-SOUTHERN_FILTERED_ROW)*(ROW_LENGTH+2) C 3 IS USED IN LOOP 320 COUNTER AS WAVE-NUMBER 0 HAS TO BE RETAINED. DO J=3+2*FILTER_WAVE_NUMBER(I),ROW_LENGTH+2 FILTER_DATA(IJ+J)=0.0 ENDDO ENDDO INCREMENT = 1 JUMP = ROW_LENGTH+2 FFTSIGN=1 LOT = ((NORTHERN_FILTERED_ROW-1)+(ROWS-SOUTHERN_FILTERED_ROW)) CALL FOURIER(FILTER_DATA,TRIGS,IFAX,INCREMENT,JUMP,ROW_LENGTH, * LOT,FFTSIGN) CL RESTORE DATA IN NORTHERN HEMISPHERE. DO I=2,NORTHERN_FILTERED_ROW IJ = (I-2)*(ROW_LENGTH+2) IK = (I-1)*ROW_LENGTH DO J=1,ROW_LENGTH FIELD(IK+J) = FILTER_DATA(IJ+J) ENDDO ENDDO CL RESTORE DATA IN SOUTHERN HEMISPHERE. IL = (NORTHERN_FILTERED_ROW-1)*(ROW_LENGTH+2) DO I=SOUTHERN_FILTERED_ROW,ROWS-1 IJ = IL+(I-SOUTHERN_FILTERED_ROW)*(ROW_LENGTH+2) IK = (I-1)*ROW_LENGTH DO J=1,ROW_LENGTH FIELD(IK+J) = FILTER_DATA(IJ+J) ENDDO ENDDO RETURN END C CLL SUBROUTINE PMSL_SETFILT ----------------------------------- CLL CLL PURPOSE: CLL SUBROUTINE 'PMSL_SETFILT' - COMPUTES FACTORS OF N & TRIGONOMETRIC CLL FUNCTIONS REQUIRED BY PMSL_FILTER CLL CLL CLL MODEL MODIFICATION HISTORY FROM MODEL VERSION 4.5 SUBROUTINE PMSL_SETFILT(TRIGS,IFAX,N,ROWS, & P_FIELD, COS_P_LATITUDE) IMPLICIT NONE INTEGER * IFAX(10) !OUT HOLDS FACTORS OF N INTEGER * N !IN NUMBER OF POINTS OF DATA ON A SLICE *, ROWS !IN NUMBER OF ROWS IN P_FIELD *, P_FIELD !IN HOLDS NUMBER OF POINTS IN THE FIELD. REAL * COS_P_LATITUDE(P_FIELD) !IN COSINE OF LATITUDE AT P POINTS. REAL * TRIGS(N) !OUT HOLDS TRIGONOMETRIC FUNCTIONS NEEDED BY * ! FOURIER, FFT99 AND FFT991 C* -------------------------------------------------------------- C*L WORKSPACE. 2 ARRAYS ARE REQUIRED ---------------------------- INTEGER * JFAX(10) ! HOLDS FACTORS OF N BEFORE THEY ARE ORDERED * ! AND STORED IN IFAX * ,LFAX(8) ! HOLDS LIST OF VALID FACTORS. LAST VALUE IS * ! MINUS PENULTIMATE ONE. THIS IS USED TO * ! SEND CODE TO THE ERROR MESSAGE AS N CONTAINS * ! AN INVALID FACTOR. C* -------------------------------------------------------------- C*L NO EXTERNAL ROUTINES --------------------------------------- C* --------------------------------------------------------------- C LOCAL VARIABLES REAL * DEL ! HOLDS ANGLE IN RADIANS BETWEEN POINTS ON LATITUDE * ! CIRCLE * ,ANGLE ! HOLDS ANGLE IN RADIANS AT A POINT ON THE LATITUDE * ! CIRCLE INTEGER * K ! LOOP COUNTER * ,I ! LOOP COUNTER * ,NU ! HOLDS CURRENT FACTORISED VALUE OF N * ,IFAC ! HOLDS CURRENT FACTOR BEING TESTED FOR * ,L ! HOLDS CURRENT POSITION IN LFAX ARRAY * ,NFAX ! HOLDS TOTAL NUMBER OF FACTORS OF N CL CL ---------------------------------------------------------------- CL SECTION 1. INITIALISATION. CL ---------------------------------------------------------------- LFAX(1)=6 LFAX(2)=8 LFAX(3)=5 LFAX(4)=4 LFAX(5)=3 LFAX(6)=2 LFAX(7)=1 LFAX(8)=-1 CL CL ---------------------------------------------------------------- CL SECTION 2. SET TRIGONOMETRIC FUNCTIONS CL ---------------------------------------------------------------- DEL=4.0*ASIN(1.0)/N DO 200 K=0,N/2-1 ANGLE=K*DEL TRIGS(2*K+1)=COS(ANGLE) TRIGS(2*K+2)=SIN(ANGLE) 200 CONTINUE CL CL ---------------------------------------------------------------- CL SECTION 3. FIND FACTORS OF N (8,6,5,4,3,2; ONLY ONE 8 ALLOWED) CL STORE FACTORS IN DESCENDING ORDER. CL ---------------------------------------------------------------- C LOOK FOR SIXES FIRST AND STORE FACTORS IN DESCENDING ORDER NU=N IFAC=6 C K HOLDS NUMBER OF FACTORS FOUND. L IS USED TO MOVE THROUGH LIST OF C VALID FACTORS. K=0 L=1 300 CONTINUE IF (MOD(NU,IFAC).EQ.0) THEN C IF IFAC IS A FACTOR OF N. K=K+1 JFAX(K)=IFAC IF (IFAC.EQ.8.AND.K.NE.1) THEN JFAX(1)=8 JFAX(K)=6 ENDIF NU=NU/IFAC C IF N FACTORISES COMPLETELY JUMP PAST ERROR MESSAGE IF(NU.EQ.1) GO TO 310 C IF FACTOR IS NOT AN 8 SEE IF IT APPEARS MORE THAN ONCE. IF(IFAC.NE.8) GO TO 300 ENDIF L=L+1 IFAC=LFAX(L) C IF NOT COMPLETELY FACTORISED GO BACK TO TOP OF PROCEDURE IF(IFAC.GT.1) GO TO 300 C ILLEGAL FACTOR IN N. PRINT ERROR MESSAGE THEN JUMP TO END OF C ROUTINE. WRITE(6,1300) N 1300 FORMAT(' ERROR IN SETFILT. N = ',I4,' CONTAINS AN ILLEGAL FACTOR') GO TO 330 C NOW REVERSE ORDER OF FACTORS SO THAT THEY ARE IN ASCENDING ORDER 310 CONTINUE NFAX=K IFAX(1)=NFAX DO 320 I=1,NFAX IFAX(NFAX+2-I)=JFAX(I) 320 CONTINUE IFAX(10)=N 330 CONTINUE CL END OF ROUTINE SETFILT RETURN END *ENDIF