diff --git a/grid_gen/ascii_netcdf_packager/Makefile b/grid_gen/ascii_netcdf_packager/Makefile index 154cc06d8..85dd726cb 100644 --- a/grid_gen/ascii_netcdf_packager/Makefile +++ b/grid_gen/ascii_netcdf_packager/Makefile @@ -1,11 +1,31 @@ +#If used, assumes ${NETCDF} is defined to the root of the netcdf library + CXX = g++ -CPPFLAGS = -I${NETCDF}/include CXXFLAGS = -O3 -LIBS = -L${NETCDF}/lib -lnetcdf -lnetcdf_c++ + EXE = AsciiNetCDFPackager.x + +ifneq ($(NETCDF), ) +ifneq ($(shell which ${NETCDF}/bin/nc-config 2> /dev/null), ) + LIBS = $(shell ${NETCDF}/bin/nc-config --libs) -lnetcdf_c++ + INCS = $(shell ${NETCDF}/bin/nc-config --cflags) +else + LIBS = -L${NETCDF}/lib + LIBS += -lnetcdf_c++ -lnetcdf -lstdc++ + INCS = -I${NETCDF}/include +endif +else ifneq ($(shell which nc-config 2> /dev/null), ) + LIBS = $(shell nc-config --libs) -lnetcdf_c++ + INCS = $(shell nc-config --cflags) +else + LIBS = -L${NETCDF}/lib + LIBS += -lnetcdf_c++ -lnetcdf -lstdc++ + INCS = -I${NETCDF}/include +endif + all: - $(CXX) $(CPPFLAGS) ascii_to_netcdf_packager.cpp $(CXXFLAGS) $(LIBS) -o $(EXE) + $(CXX) $(INCS) ascii_to_netcdf_packager.cpp $(CXXFLAGS) $(LIBS) -o $(EXE) clean: rm -f $(EXE) diff --git a/grid_gen/basin/src/Makefile b/grid_gen/basin/src/Makefile index aade64cf0..e7ac4ef85 100644 --- a/grid_gen/basin/src/Makefile +++ b/grid_gen/basin/src/Makefile @@ -38,15 +38,26 @@ LDFLAGS = #-g -traceback -check all CPP = cpp -P -traditional CPPFLAGS = CPPINCLUDES = -INCLUDES = -I$(NETCDF)/include - -LIBS = -L$(NETCDF)/lib -NCLIB = -lnetcdf -NCLIBF = -lnetcdff -ifneq ($(wildcard $(NETCDF)/lib/libnetcdff.*), ) # CHECK FOR NETCDF4 - LIBS += $(NCLIBF) -endif # CHECK FOR NETCDF4 -LIBS += $(NCLIB) + +ifneq ($(NETCDF), ) +ifneq ($(shell which ${NETCDF}/bin/nc-config 2> /dev/null), ) + LIBS = $(shell ${NETCDF}/bin/nc-config --libs) -lnetcdff + INCS = $(shell ${NETCDF}/bin/nc-config --cflags) +else + LIBS= -L${NETCDF}/lib + LIBS += -lnetcdf -lnetcdff + INCS = -I${NETCDF}/include +endif +else ifneq ($(shell which nc-config 2> /dev/null), ) + LIBS = $(shell nc-config --libs) -lnetcdff + INCS = $(shell nc-config --cflags) +else + LIBS= -L${NETCDF}/lib + LIBS += -lnetcdf_c++ -lnetcdf -lstdc++ -lnetcdff + INCS = -I${NETCDF}/include +endif + +INCLUDES = $(INCS) RM = rm -f diff --git a/grid_gen/global_scvt/Makefile b/grid_gen/global_scvt/Makefile index aaef05ca9..dc8e1fe55 100644 --- a/grid_gen/global_scvt/Makefile +++ b/grid_gen/global_scvt/Makefile @@ -1,23 +1,23 @@ #FC = ifort -#FFLAGS = -FR -m64 -O3 -fast -ipo -openmp +#FFLAGS = -FR -m64 -O3 -fast -ipo -openmp $(shell nc-config --fflags) #F77FLAGS = -FI -m64 -O3 -fast -ipo -openmp #CPPFLAGS = -DRKIND=8 #PROMOTION = -r8 -#LDFLAGS = -m64 -O3 -fast -ipo -openmp +#LDFLAGS = -m64 -O3 -fast -ipo -openmp $(shell nc-config --flibs) FC = gfortran -FFLAGS = -ffree-form -O3 -fopenmp -ffree-line-length-none +FFLAGS = -ffree-form -O3 -fopenmp -ffree-line-length-none $(shell nc-config --fflags) F77FLAGS = -ffixed-form -O3 -fopenmp -fsecond-underscore -CPPFLAGS = -DRKIND=8 +CPPFLAGS = -DRKIND=8 -D_RKIND=_8 PROMOTION = -fdefault-real-8 -LDFLAGS = -O3 -fopenmp +LDFLAGS = -O3 -fopenmp $(shell nc-config --flibs) #FC = pgf90 -#FFLAGS = -Mfree -O3 -mp -byteswapio +#FFLAGS = -Mfree -O3 -mp -byteswapio $(shell nc-config --fflags) #F77FLAGS = -O3 -byteswapio -#CPPFLAGS = -DRKIND=8 +#CPPFLAGS = -DRKIND=8 -D_RKIND=_8 #PROMOTION = -r8 -#LDFLAGS = -O3 -mp -byteswapio +#LDFLAGS = -O3 -mp -byteswapio $(shell nc-config --flibs) all: grid_gen grid_ref diff --git a/grid_gen/global_scvt/src/Makefile b/grid_gen/global_scvt/src/Makefile index 9f44d4d9e..c8229ebad 100644 --- a/grid_gen/global_scvt/src/Makefile +++ b/grid_gen/global_scvt/src/Makefile @@ -3,7 +3,7 @@ OBJS = STRIPACK.o module_grid_params.o module_grid_constants.o module_data_types.o module_sphere_utilities.o module_voronoi_utils.o module_grid_gen_utils.o module_scvt.o module_write_netcdf.o module_grid_meta.o grid_gen.o all: $(OBJS) - $(FC) $(PROMOTION) $(LDFLAGS) -o grid_gen $(OBJS) -L$(NETCDF)/lib -lnetcdff -lnetcdf + $(FC) $(PROMOTION) $(LDFLAGS) -o grid_gen $(OBJS) $(shell nc-config --fflags) grid_gen.o: module_grid_params.o module_grid_constants.o module_data_types.o module_grid_gen_utils.o module_voronoi_utils.o STRIPACK.o module_scvt.o module_grid_meta.o @@ -28,12 +28,14 @@ module_voronoi_utils.o: module_grid_constants.o STRIPACK.o .F.o: - cpp -C -P -traditional $(CPPFLAGS) $< > $*.f90 - $(FC) $(FFLAGS) $(PROMOTION) -c $*.f90 -I$(NETCDF)/include + gcc-4.8 -C -P -traditional-cpp $(CPPFLAGS) $< > $*.f90 + $(FC) $(FFLAGS) $(PROMOTION) -c $*.f90 rm -f $*.f90 .f.o: - $(FC) $(F77FLAGS) $(PROMOTION) -c $< + gcc-4.8 -C -P -traditional-cpp $(CPPFLAGS) $< > $*.for + $(FC) $(F77FLAGS) $(PROMOTION) -c $*.for + rm -f $.*.for clean: rm -f *.o *.mod grid_gen diff --git a/grid_gen/global_scvt/src/STRIPACK.f b/grid_gen/global_scvt/src/STRIPACK.f index 968c818a5..d50595a73 100644 --- a/grid_gen/global_scvt/src/STRIPACK.f +++ b/grid_gen/global_scvt/src/STRIPACK.f @@ -5,7 +5,7 @@ MODULE STRIPACK SUBROUTINE ADDNOD (NST,K,X,Y,Z, LIST,LPTR,LEND, . LNEW, IER) INTEGER NST, K, LIST(*), LPTR(*), LEND(K), LNEW, IER - REAL X(K), Y(K), Z(K) + REAL(KIND=RKIND) X(K), Y(K), Z(K) C C*********************************************************** C @@ -82,7 +82,7 @@ SUBROUTINE ADDNOD (NST,K,X,Y,Z, LIST,LPTR,LEND, C INTEGER I1, I2, I3, IO1, IO2, IN1, IST, KK, KM1, L, . LP, LPF, LPO1, LPO1S - REAL B1, B2, B3, P(3) + REAL(KIND=RKIND) B1, B2, B3, P(3) C C Local parameters: C @@ -206,8 +206,8 @@ SUBROUTINE ADDNOD (NST,K,X,Y,Z, LIST,LPTR,LEND, 5 IER = L RETURN END SUBROUTINE - REAL FUNCTION AREAS (V1,V2,V3) - REAL V1(3), V2(3), V3(3) + REAL(KIND=RKIND) FUNCTION AREAS (V1,V2,V3) + REAL(KIND=RKIND) V1(3), V2(3), V3(3) C C*********************************************************** C @@ -554,7 +554,7 @@ SUBROUTINE BNODES (N,LIST,LPTR,LEND, NODES,NB,NA,NT) END SUBROUTINE SUBROUTINE CIRCUM (V1,V2,V3, C,IER) INTEGER IER - REAL V1(3), V2(3), V3(3), C(3) + REAL(KIND=RKIND) V1(3), V2(3), V3(3), C(3) C C*********************************************************** C @@ -601,7 +601,7 @@ SUBROUTINE CIRCUM (V1,V2,V3, C,IER) C*********************************************************** C INTEGER I - REAL CNORM, CU(3), E1(3), E2(3) + REAL(KIND=RKIND) CNORM, CU(3), E1(3), E2(3) C C Local parameters: C @@ -734,7 +734,7 @@ SUBROUTINE CRLIST (N,NCOL,X,Y,Z,LIST,LEND, LPTR,LNEW, . LTRI, LISTC,NB,XC,YC,ZC,RC,IER) INTEGER N, NCOL, LIST(*), LEND(N), LPTR(*), LNEW, . LTRI(6,NCOL), LISTC(*), NB, IER - REAL X(N), Y(N), Z(N), XC(*), YC(*), ZC(*), RC(*) + REAL(KIND=RKIND) X(N), Y(N), Z(N), XC(*), YC(*), ZC(*), RC(*) C C*********************************************************** C @@ -882,7 +882,7 @@ SUBROUTINE CRLIST (N,NCOL,X,Y,Z,LIST,LEND, LPTR,LNEW, . KT12, KT21, KT22, LP, LPL, LPN, N0, N1, N2, . N3, N4, NM2, NN, NT LOGICAL SWP - REAL C(3), T, V1(3), V2(3), V3(3) + REAL(KIND=RKIND) C(3), T, V1(3), V2(3), V3(3) C C Local parameters: C @@ -1560,7 +1560,7 @@ SUBROUTINE DELNOD (K, N,X,Y,Z,LIST,LPTR,LEND,LNEW,LWK, . IWK, IER) INTEGER K, N, LIST(*), LPTR(*), LEND(*), LNEW, LWK, . IWK(2,*), IER - REAL X(*), Y(*), Z(*) + REAL(KIND=RKIND) X(*), Y(*), Z(*) C C*********************************************************** C @@ -1668,7 +1668,8 @@ SUBROUTINE DELNOD (K, N,X,Y,Z,LIST,LPTR,LEND,LNEW,LWK, . LPL2, LPN, LWKL, N1, N2, NFRST, NIT, NL, NN, . NNB, NR LOGICAL BDRY - REAL X1, X2, XL, XR, Y1, Y2, YL, YR, Z1, Z2, ZL, ZR + REAL(KIND=RKIND) X1, X2, XL, XR, Y1, Y2, YL, YR, + . Z1, Z2, ZL, ZR C C Local parameters: C @@ -2001,7 +2002,7 @@ SUBROUTINE EDGE (IN1,IN2,X,Y,Z, LWK,IWK,LIST,LPTR, . LEND, IER) INTEGER IN1, IN2, LWK, IWK(2,*), LIST(*), LPTR(*), . LEND(*), IER - REAL X(*), Y(*), Z(*) + REAL(KIND=RKIND) X(*), Y(*), Z(*) C C*********************************************************** C @@ -2103,8 +2104,8 @@ SUBROUTINE EDGE (IN1,IN2,X,Y,Z, LWK,IWK,LIST,LPTR, INTEGER I, IERR, IWC, IWCP1, IWEND, IWF, IWL, LFT, LP, . LP21, LPL, N0, N1, N1FRST, N1LST, N2, NEXT, . NIT, NL, NR - REAL DP12, DP1L, DP1R, DP2L, DP2R, X0, X1, X2, Y0, - . Y1, Y2, Z0, Z1, Z2 + REAL(KIND=RKIND) DP12, DP1L, DP1R, DP2L, DP2R, + . X0, X1, X2, Y0, Y1, Y2, Z0, Z1, Z2 C C Local parameters: C @@ -2530,7 +2531,7 @@ SUBROUTINE EDGE (IN1,IN2,X,Y,Z, LWK,IWK,LIST,LPTR, SUBROUTINE GETNP (X,Y,Z,LIST,LPTR,LEND,L, NPTS, DF, . IER) INTEGER LIST(*), LPTR(*), LEND(*), L, NPTS(L), IER - REAL X(*), Y(*), Z(*), DF + REAL(KIND=RKIND) X(*), Y(*), Z(*), DF C C*********************************************************** C @@ -2594,7 +2595,7 @@ SUBROUTINE GETNP (X,Y,Z,LIST,LPTR,LEND,L, NPTS, DF, C*********************************************************** C INTEGER I, LM1, LP, LPL, N1, NB, NI, NP - REAL DNB, DNP, X1, Y1, Z1 + REAL(KIND=RKIND) DNB, DNP, X1, Y1, Z1 C C Local parameters: C @@ -2722,7 +2723,7 @@ SUBROUTINE INSERT (K,LP, LIST,LPTR,LNEW ) END SUBROUTINE LOGICAL FUNCTION INSIDE (P,LV,XV,YV,ZV,NV,LISTV, IER) INTEGER LV, NV, LISTV(NV), IER - REAL P(3), XV(LV), YV(LV), ZV(LV) + REAL(KIND=RKIND) P(3), XV(LV), YV(LV), ZV(LV) C C*********************************************************** C @@ -2818,8 +2819,8 @@ LOGICAL FUNCTION INSIDE (P,LV,XV,YV,ZV,NV,LISTV, IER) C INTEGER I1, I2, IERR, IMX, K, K0, N, NI LOGICAL EVEN, LFT1, LFT2, PINR, QINR - REAL B(3), BP, BQ, CN(3), D, EPS, PN(3), Q(3), - . QN(3), QNRM, V1(3), V2(3), VN(3), VNRM + REAL(KIND=RKIND) B(3), BP, BQ, CN(3), D, EPS, PN(3), + . Q(3), QN(3), QNRM, V1(3), V2(3), VN(3), VNRM C C Local parameters: C @@ -3115,7 +3116,7 @@ SUBROUTINE INTADD (KK,I1,I2,I3, LIST,LPTR,LEND,LNEW ) END SUBROUTINE SUBROUTINE INTRSC (P1,P2,CN, P,IER) INTEGER IER - REAL P1(3), P2(3), CN(3), P(3) + REAL(KIND=RKIND) P1(3), P2(3), CN(3), P(3) C C*********************************************************** C @@ -3172,7 +3173,7 @@ SUBROUTINE INTRSC (P1,P2,CN, P,IER) C*********************************************************** C INTEGER I - REAL D1, D2, PP(3), PPN, T + REAL(KIND=RKIND) D1, D2, PP(3), PPN, T C C Local parameters: C @@ -3283,7 +3284,7 @@ INTEGER FUNCTION JRAND (N, IX,IY,IZ ) RETURN END FUNCTION LOGICAL FUNCTION LEFT (X1,Y1,Z1,X2,Y2,Z2,X0,Y0,Z0) - REAL X1, Y1, Z1, X2, Y2, Z2, X0, Y0, Z0 + REAL(KIND=RKIND) X1, Y1, Z1, X2, Y2, Z2, X0, Y0, Z0 C C*********************************************************** C @@ -3439,7 +3440,7 @@ INTEGER FUNCTION NBCNT (LPL,LPTR) INTEGER FUNCTION NEARND (P,IST,N,X,Y,Z,LIST,LPTR, . LEND, AL) INTEGER IST, N, LIST(*), LPTR(*), LEND(N) - REAL P(3), X(N), Y(N), Z(N), AL + REAL(KIND=RKIND) P(3), X(N), Y(N), Z(N), AL C C*********************************************************** C @@ -3506,8 +3507,8 @@ INTEGER FUNCTION NEARND (P,IST,N,X,Y,Z,LIST,LPTR, PARAMETER (LMAX=25) INTEGER I1, I2, I3, L, LISTP(LMAX), LP, LP1, LP2, . LPL, LPTRP(LMAX), N1, N2, N3, NN, NR, NST - REAL B1, B2, B3, DS1, DSR, DX1, DX2, DX3, DY1, - . DY2, DY3, DZ1, DZ2, DZ3 + REAL(KIND=RKIND) B1, B2, B3, DS1, DSR, DX1, DX2, + . DX3, DY1, DY2, DY3, DZ1, DZ2, DZ3 C C Local parameters: C @@ -3679,7 +3680,7 @@ SUBROUTINE OPTIM (X,Y,Z,NA, LIST,LPTR,LEND,NIT, . IWK, IER) INTEGER NA, LIST(*), LPTR(*), LEND(*), NIT, IWK(2,NA), . IER - REAL X(*), Y(*), Z(*) + REAL(KIND=RKIND) X(*), Y(*), Z(*) C C*********************************************************** C @@ -3871,7 +3872,7 @@ SUBROUTINE OPTIM (X,Y,Z,NA, LIST,LPTR,LEND,NIT, RETURN END SUBROUTINE SUBROUTINE SCOORD (PX,PY,PZ, PLAT,PLON,PNRM) - REAL PX, PY, PZ, PLAT, PLON, PNRM + REAL(KIND=RKIND) PX, PY, PZ, PLAT, PLON, PNRM C C*********************************************************** C @@ -3924,8 +3925,8 @@ SUBROUTINE SCOORD (PX,PY,PZ, PLAT,PLON,PNRM) ENDIF RETURN END SUBROUTINE - REAL FUNCTION STORE (X) - REAL X + REAL(KIND=RKIND) FUNCTION STORE (X) + REAL(KIND=RKIND) X C C*********************************************************** C @@ -3959,7 +3960,7 @@ REAL FUNCTION STORE (X) C C*********************************************************** C - REAL Y + REAL(KIND=RKIND) Y COMMON/STCOM/Y Y = X STORE = Y @@ -4077,7 +4078,7 @@ SUBROUTINE SWAP (IN1,IN2,IO1,IO2, LIST,LPTR, END SUBROUTINE LOGICAL FUNCTION SWPTST (N1,N2,N3,N4,X,Y,Z) INTEGER N1, N2, N3, N4 - REAL X(*), Y(*), Z(*) + REAL(KIND=RKIND) X(*), Y(*), Z(*) C C*********************************************************** C @@ -4124,8 +4125,8 @@ LOGICAL FUNCTION SWPTST (N1,N2,N3,N4,X,Y,Z) C C*********************************************************** C - REAL DX1, DX2, DX3, DY1, DY2, DY3, DZ1, DZ2, DZ3, - . X4, Y4, Z4 + REAL(KIND=RKIND) DX1, DX2, DX3, DY1, DY2, DY3, DZ1, + . DZ2, DZ3, X4, Y4, Z4 C C Local parameters: C @@ -4158,7 +4159,7 @@ LOGICAL FUNCTION SWPTST (N1,N2,N3,N4,X,Y,Z) END FUNCTION SUBROUTINE TRANS (N,RLAT,RLON, X,Y,Z) INTEGER N - REAL RLAT(N), RLON(N), X(N), Y(N), Z(N) + REAL(KIND=RKIND) RLAT(N), RLON(N), X(N), Y(N), Z(N) C C*********************************************************** C @@ -4203,7 +4204,7 @@ SUBROUTINE TRANS (N,RLAT,RLON, X,Y,Z) C*********************************************************** C INTEGER I, NN - REAL COSPHI, PHI, THETA + REAL(KIND=RKIND) COSPHI, PHI, THETA C C Local parameters: C @@ -4227,7 +4228,7 @@ SUBROUTINE TRANS (N,RLAT,RLON, X,Y,Z) SUBROUTINE TRFIND (NST,P,N,X,Y,Z,LIST,LPTR,LEND, B1, . B2,B3,I1,I2,I3) INTEGER NST, N, LIST(*), LPTR(*), LEND(N), I1, I2, I3 - REAL P(3), X(N), Y(N), Z(N), B1, B2, B3 + REAL(KIND=RKIND) P(3), X(N), Y(N), Z(N), B1, B2, B3 C C*********************************************************** C @@ -4296,9 +4297,9 @@ SUBROUTINE TRFIND (NST,P,N,X,Y,Z,LIST,LPTR,LEND, B1, C INTEGER IX, IY, IZ, LP, N0, N1, N1S, N2, N2S, N3, N4, . NEXT, NF, NL - REAL DET, EPS, PTN1, PTN2, Q(3), S12, TOL, XP, YP, - . ZP - REAL X0, X1, X2, Y0, Y1, Y2, Z0, Z1, Z2 + REAL(KIND=RKIND) DET, EPS, PTN1, PTN2, Q(3), S12, TOL, + . XP, YP, ZP + REAL(KIND=RKIND) X0, X1, X2, Y0, Y1, Y2, Z0, Z1, Z2 C SAVE IX, IY, IZ DATA IX/1/, IY/2/, IZ/3/ @@ -4915,7 +4916,7 @@ SUBROUTINE TRLIST (N,LIST,LPTR,LEND,NROW, NT,LTRI,IER) END SUBROUTINE SUBROUTINE TRLPRT (N,X,Y,Z,IFLAG,NROW,NT,LTRI,LOUT) INTEGER N, IFLAG, NROW, NT, LTRI(NROW,NT), LOUT - REAL X(N), Y(N), Z(N) + REAL(KIND=RKIND) X(N), Y(N), Z(N) C C*********************************************************** C @@ -5110,7 +5111,7 @@ SUBROUTINE TRMESH (N,X,Y,Z, LIST,LPTR,LEND,LNEW,NEAR, . NEXT,DIST,IER) INTEGER N, LIST(*), LPTR(*), LEND(N), LNEW, NEAR(N), . NEXT(N), IER - REAL X(N), Y(N), Z(N), DIST(N) + REAL(KIND=RKIND) X(N), Y(N), Z(N), DIST(N) C C*********************************************************** C @@ -5311,7 +5312,7 @@ SUBROUTINE TRMESH (N,X,Y,Z, LIST,LPTR,LEND,LNEW,NEAR, C*********************************************************** C INTEGER I, I0, J, K, LP, LPL, NEXTI, NN - REAL D, D1, D2, D3 + REAL(KIND=RKIND) D, D1, D2, D3 C C Local parameters: C @@ -5531,7 +5532,7 @@ SUBROUTINE TRPLOT (LUN,PLTSIZ,ELAT,ELON,A,N,X,Y,Z, CHARACTER*(*) TITLE INTEGER LUN, N, LIST(*), LPTR(*), LEND(N), IER LOGICAL NUMBR - REAL PLTSIZ, ELAT, ELON, A, X(N), Y(N), Z(N) + REAL(KIND=RKIND) PLTSIZ, ELAT, ELON, A, X(N), Y(N), Z(N) C C*********************************************************** C @@ -5638,9 +5639,9 @@ SUBROUTINE TRPLOT (LUN,PLTSIZ,ELAT,ELON,A,N,X,Y,Z, C INTEGER IPX1, IPX2, IPY1, IPY2, IR, LP, LPL, N0, N1 LOGICAL ANNOT - REAL CF, CT, EX, EY, EZ, FSIZN, FSIZT, R11, R12, - . R21, R22, R23, SF, T, TX, TY, WR, WRS, X0, X1, - . Y0, Y1, Z0, Z1 + REAL(KIND=RKIND) CF, CT, EX, EY, EZ, FSIZN, FSIZT, + . R11, R12, R21, R22, R23, SF, T, TX, TY, WR, + . WRS, X0, X1, Y0, Y1, Z0, Z1 C DATA ANNOT/.TRUE./, FSIZN/10.0/, FSIZT/16.0/ C @@ -5958,7 +5959,7 @@ SUBROUTINE TRPLOT (LUN,PLTSIZ,ELAT,ELON,A,N,X,Y,Z, END SUBROUTINE SUBROUTINE TRPRNT (N,X,Y,Z,IFLAG,LIST,LPTR,LEND,LOUT) INTEGER N, IFLAG, LIST(*), LPTR(*), LEND(N), LOUT - REAL X(N), Y(N), Z(N) + REAL(KIND=RKIND) X(N), Y(N), Z(N) C C*********************************************************** C @@ -6221,8 +6222,8 @@ SUBROUTINE VRPLOT (LUN,PLTSIZ,ELAT,ELON,A,N,X,Y,Z, CHARACTER*(*) TITLE INTEGER LUN, N, NT, LISTC(*), LPTR(*), LEND(N), IER LOGICAL NUMBR - REAL PLTSIZ, ELAT, ELON, A, X(N), Y(N), Z(N), - . XC(NT), YC(NT), ZC(NT) + REAL(KIND=RKIND) PLTSIZ, ELAT, ELON, A, X(N), + . Y(N), Z(N), XC(NT), YC(NT), ZC(NT) C C*********************************************************** C @@ -6366,9 +6367,9 @@ SUBROUTINE VRPLOT (LUN,PLTSIZ,ELAT,ELON,A,N,X,Y,Z, INTEGER IPX1, IPX2, IPY1, IPY2, IR, KV1, KV2, LP, LPL, . N0 LOGICAL ANNOT, IN1, IN2 - REAL CF, CT, EX, EY, EZ, FSIZN, FSIZT, R11, R12, - . R21, R22, R23, SF, T, TX, TY, WR, WRS, X0, X1, - . X2, Y0, Y1, Y2, Z1, Z2 + REAL(KIND=RKIND) CF, CT, EX, EY, EZ, FSIZN, FSIZT, + . R11, R12, R21, R22, R23, SF, T, TX, TY, WR, + . WRS, X0, X1, X2, Y0, Y1, Y2, Z1, Z2 C DATA ANNOT/.TRUE./, FSIZN/10.0/, FSIZT/16.0/ C diff --git a/grid_gen/global_scvt/src/grid_gen.F b/grid_gen/global_scvt/src/grid_gen.F index 34d64bb50..48d394708 100644 --- a/grid_gen/global_scvt/src/grid_gen.F +++ b/grid_gen/global_scvt/src/grid_gen.F @@ -11,20 +11,20 @@ program grid_gen implicit none - real :: dlat, dlon - real :: dl - real :: d1, d2, d3, d4 + real(kind=RKIND) :: dlat, dlon + real(kind=RKIND) :: dl + real(kind=RKIND) :: d1, d2, d3, d4 integer :: p1, p2, p3 integer :: if character (len=80) :: frame_name - real :: pi - real :: area_per_sample, nhexs, sum_nhexs, hex_area + real(kind=RKIND) :: pi + real(kind=RKIND) :: area_per_sample, nhexs, sum_nhexs, hex_area type (geo_point) :: p integer :: i, j, k, nb, ier - real, allocatable, dimension(:) :: rlat, rlon, vclat, vclon, x, y, z, xc, yc, zc + real(kind=RKIND), allocatable, dimension(:) :: rlat, rlon, vclat, vclon, x, y, z, xc, yc, zc integer, allocatable, dimension(:) :: list, lptr, listc, lend integer, allocatable, dimension(:,:) :: ltri @@ -34,24 +34,24 @@ program grid_gen type (adjacency_list) :: alist, clist - real :: lat1, lon1, lat2, lon2, lat3, lon3, latc, lonc + real(kind=RKIND) :: lat1, lon1, lat2, lon2, lat3, lon3, latc, lonc call read_namelist() - pi = 4.0*atan(1.0) + pi = 4.0_RKIND*atan(1.0_RKIND) - area_per_sample = 4.0 * pi * 6370000**2.0 / 6000000.0 - sum_nhexs = 0.0 + area_per_sample = 4.0_RKIND * pi * 6370000**2.0_RKIND / 6000000.0_RKIND + sum_nhexs = 0.0_RKIND write(0,'(a,f10.1)') 'Computing an estimate for the required number of cells to reach dx=', min_dx do if = 1,5 - nhexs = 0.0 + nhexs = 0.0_RKIND do i=1,6000000 call random_point(p) d1 = density_for_point(p) dl = min_dx / (d1 ** 0.25) - hex_area = sqrt(3.0) / 2.0 * dl**2.0 + hex_area = sqrt(3.0_RKIND) / 2.0_RKIND * dl**2.0_RKIND nhexs = nhexs + area_per_sample / hex_area end do ! write(0,'(a,i2,a,i)') 'Estimate ',if,' for required # hexs:', nint(nhexs) @@ -59,7 +59,7 @@ program grid_gen write(0,'(a,i3,a)',advance='no') ' ...',if*20,'%' end do write(0,*) ' ' - write(0,*) 'Estimated # hexs:', nint(sum_nhexs/5.0) + write(0,*) 'Estimated # hexs:', nint(sum_nhexs/5.0_RKIND) write(0,*) ' ' @@ -117,7 +117,7 @@ program grid_gen call TRANS (nvc, vclat, vclon, xc, yc, zc) write(frame_name,'(a)') 'scvt_initial.ps' open(32,file=trim(frame_name),form='formatted',status='unknown') - call vrplot(32, 8.0, 0.0, 0.0, 90.0 ,N, X,Y,Z, 2*n-4,LISTC,LPTR,LEND,XC,YC,ZC,'(spherical centroidal voronoi tessellation)',.false.,IER) + call vrplot(32, 8.0_RKIND, 0.0_RKIND, 0.0_RKIND, 90.0_RKIND ,N, X,Y,Z, 2*n-4,LISTC,LPTR,LEND,XC,YC,ZC,'(spherical centroidal voronoi tessellation)',.false.,IER) close(32) call scvt_solve(n, lend, rlat, rlon, nvc, list, lptr, if) @@ -127,7 +127,7 @@ program grid_gen call TRANS (nvc, vclat, vclon, xc, yc, zc) write(frame_name,'(a)') 'scvt_final.ps' open(32,file=trim(frame_name),form='formatted',status='unknown') - call vrplot(32, 8.0, 0.0, 0.0, 90.0 ,N, X,Y,Z, 2*n-4,LISTC,LPTR,LEND,XC,YC,ZC,'(spherical centroidal voronoi tessellation)',.false.,IER) + call vrplot(32, 8.0_RKIND, 0.0_RKIND, 0.0_RKIND, 90.0_RKIND ,N, X,Y,Z, 2*n-4,LISTC,LPTR,LEND,XC,YC,ZC,'(spherical centroidal voronoi tessellation)',.false.,IER) close(32) diff --git a/grid_gen/global_scvt/src/module_data_types.F b/grid_gen/global_scvt/src/module_data_types.F index 011bd8fbf..e9c6816b7 100644 --- a/grid_gen/global_scvt/src/module_data_types.F +++ b/grid_gen/global_scvt/src/module_data_types.F @@ -3,7 +3,7 @@ module data_types integer, parameter :: LESS = -1, EQUAL = 0, GREATER = 1 type geo_point - real :: lat, lon + real(kind=RKIND) :: lat, lon end type geo_point type send_list_ptr @@ -29,7 +29,7 @@ module data_types type binary_tree integer :: node1, node2 integer :: vertex1, vertex2 - real :: lat1, lon1, lat2, lon2 + real(kind=RKIND) :: lat1, lon1, lat2, lon2 type (binary_tree), pointer :: left, right, parent end type binary_tree diff --git a/grid_gen/global_scvt/src/module_grid_gen_utils.F b/grid_gen/global_scvt/src/module_grid_gen_utils.F index 80c0fbd52..49263773b 100644 --- a/grid_gen/global_scvt/src/module_grid_gen_utils.F +++ b/grid_gen/global_scvt/src/module_grid_gen_utils.F @@ -238,7 +238,7 @@ subroutine compute_h_area(corners, centers, areas, n) if(itp1 > 6) itp1 = 1 hex_area = hex_area + triangle_area( centers(i,j), & corners( it,i,j), & - corners(itp1,i,j), 1.) + corners(itp1,i,j), 1.0_RKIND) enddo areas(i,j) = hex_area enddo @@ -269,7 +269,7 @@ subroutine compute_edge_lengths(corners, edge_lengths, n) itp1 = it+1 if(itp1 > 6) itp1 = 1 edge_lengths(it,i,j) = sphere_distance( corners( it,i,j), & - corners(itp1,i,j), 1.) + corners(itp1,i,j), 1.0_RKIND) end do end do end do @@ -296,11 +296,11 @@ subroutine compute_dx( centers, dx, n ) do j=1,n do i=1,2*n-1 dx(1,i,j) = sphere_distance( centers(i ,j ), & - centers(i-1,j ), 1. ) + centers(i-1,j ), 1.0_RKIND ) dx(2,i,j) = sphere_distance( centers(i ,j ), & - centers(i-1,j-1), 1. ) + centers(i-1,j-1), 1.0_RKIND ) dx(3,i,j) = sphere_distance( centers(i ,j ), & - centers(i ,j-1), 1. ) + centers(i ,j-1), 1.0_RKIND ) enddo enddo diff --git a/grid_gen/global_scvt/src/module_grid_meta.F b/grid_gen/global_scvt/src/module_grid_meta.F index 4d513830a..db91c6fda 100644 --- a/grid_gen/global_scvt/src/module_grid_meta.F +++ b/grid_gen/global_scvt/src/module_grid_meta.F @@ -18,8 +18,8 @@ subroutine write_grid(rlat, rlon, n, vclat, vclon, nvc, alist, clist) implicit none integer, intent(in) :: n, nvc - real, dimension(n), intent(inout) :: rlat, rlon - real, dimension(nvc), intent(inout) :: vclat, vclon + real(kind=RKIND), dimension(n), intent(inout) :: rlat, rlon + real(kind=RKIND), dimension(nvc), intent(inout) :: vclat, vclon type (adjacency_list) :: alist, clist integer, parameter :: maxEdges = 10 @@ -38,16 +38,16 @@ subroutine write_grid(rlat, rlon, n, vclat, vclon, nvc, alist, clist) integer, allocatable, dimension(:,:) :: cellsOnEdge, edgesOnCell, verticesOnCell, cellsOnCell, & verticesOnEdge, edgesOnEdge, edgesOnVertex, cellsOnVertex integer, allocatable, dimension(:) :: isObtuse - real, allocatable, dimension(:,:) :: weightsOnEdge, kiteAreasOnVertex + real(kind=RKIND), allocatable, dimension(:,:) :: weightsOnEdge, kiteAreasOnVertex integer :: temp logical :: tdrtest = .true. - real :: sum_r, area, r, s, de, rtmp - real, allocatable, dimension(:) :: latCell, lonCell, latEdge, lonEdge, angleEdge, latVertex, lonVertex, & + real(kind=RKIND) :: sum_r, area, r, s, de, rtmp + real(kind=RKIND), allocatable, dimension(:) :: latCell, lonCell, latEdge, lonEdge, angleEdge, latVertex, lonVertex, & lat1Edge, lon1Edge, lat2Edge, lon2Edge, dvEdge, dv1Edge, dv2Edge, dcEdge, & areaCell, areaTriangle, fEdge, fVertex, h_s, u_sbr - real, allocatable, dimension(:,:,:) :: u, v, h, vh, circulation, vorticity, ke - real, allocatable, dimension(:,:,:,:) :: tracers - real, allocatable, dimension(:) :: xCell, yCell, zCell, xEdge, yEdge, zEdge, xVertex, yVertex, zVertex, meshDensity + real(kind=RKIND), allocatable, dimension(:,:,:) :: u, v, h, vh, circulation, vorticity, ke + real(kind=RKIND), allocatable, dimension(:,:,:,:) :: tracers + real(kind=RKIND), allocatable, dimension(:) :: xCell, yCell, zCell, xEdge, yEdge, zEdge, xVertex, yVertex, zVertex, meshDensity type (geo_point) :: vertex1GP, vertex2GP, cell1GP, cell2GP, cell3GP, edgeGP, edgeGP_prev, edgeGP_next, pCell type (geo_point) :: center type (geo_point), allocatable, dimension(:) :: points @@ -207,7 +207,7 @@ subroutine write_grid(rlat, rlon, n, vclat, vclon, nvc, alist, clist) cell1GP % lon = lonCell(cellsOnEdge(1,i)) cell2GP % lat = latCell(cellsOnEdge(2,i)) cell2GP % lon = lonCell(cellsOnEdge(2,i)) - dcEdge(i) = sphere_distance(cell1GP, cell2GP, 1.0) + dcEdge(i) = sphere_distance(cell1GP, cell2GP, 1.0_RKIND) end do @@ -291,10 +291,10 @@ subroutine write_grid(rlat, rlon, n, vclat, vclon, nvc, alist, clist) cell2GP % lon = lonCell(cellsOnEdge(2,j)) call compute_edge_latlon(cell1GP, cell2GP, vertex1GP, vertex2GP, edgeGP) - dvEdge(j) = sphere_distance(vertex1GP, vertex2GP, 1.0) + dvEdge(j) = sphere_distance(vertex1GP, vertex2GP, 1.0_RKIND) rtmp = (vertex2GP % lat - vertex1GP % lat) / dvEdge(j) - if (rtmp > 1.0) rtmp = 1.0 - if (rtmp < -1.0) rtmp = -1.0 + if (rtmp > 1.0_RKIND) rtmp = 1.0_RKIND + if (rtmp < -1.0_RKIND) rtmp = -1.0_RKIND rtmp = acos(rtmp) angleEdge(j) = meridian_angle(edgeGP, vertex2GP) angleEdge(j) = sign(rtmp, angleEdge(j)) @@ -315,7 +315,7 @@ subroutine write_grid(rlat, rlon, n, vclat, vclon, nvc, alist, clist) temp = verticesOnEdge(1,j) verticesOnEdge(1,j) = verticesOnEdge(2,j) verticesOnEdge(2,j) = temp - u_sbr(j) = -1.0*u_sbr(j) + u_sbr(j) = -1.0_RKIND*u_sbr(j) angleEdge(j) = angleEdge(j) + pii if (angleEdge(j) > pii) angleEdge(j) = angleEdge(j) - 2.0*pii if (angleEdge(j) < -pii) angleEdge(j) = angleEdge(j) + 2.0*pii @@ -457,13 +457,13 @@ subroutine write_grid(rlat, rlon, n, vclat, vclon, nvc, alist, clist) vertex1GP % lon = lonVertex(verticesOnCell(j,i)) vertex2GP % lat = latVertex(verticesOnCell(j+1,i)) vertex2GP % lon = lonVertex(verticesOnCell(j+1,i)) - areaCell(i) = areaCell(i) + triangle_area(cell1GP, vertex1GP, vertex2GP, 1.0) + areaCell(i) = areaCell(i) + triangle_area(cell1GP, vertex1GP, vertex2GP, 1.0_RKIND) end do vertex1GP % lat = latVertex(verticesOnCell(j,i)) vertex1GP % lon = lonVertex(verticesOnCell(j,i)) vertex2GP % lat = latVertex(verticesOnCell(1,i)) vertex2GP % lon = lonVertex(verticesOnCell(1,i)) - areaCell(i) = areaCell(i) + triangle_area(cell1GP, vertex1GP, vertex2GP, 1.0) + areaCell(i) = areaCell(i) + triangle_area(cell1GP, vertex1GP, vertex2GP, 1.0_RKIND) end do @@ -477,7 +477,7 @@ subroutine write_grid(rlat, rlon, n, vclat, vclon, nvc, alist, clist) cell2GP % lon = lonCell(cellsOnVertex(2,i)) cell3GP % lat = latCell(cellsOnVertex(3,i)) cell3GP % lon = lonCell(cellsOnVertex(3,i)) - areaTriangle(i) = triangle_area(cell1GP, cell2GP, cell3GP, 1.0) + areaTriangle(i) = triangle_area(cell1GP, cell2GP, cell3GP, 1.0_RKIND) end do ! @@ -497,7 +497,7 @@ subroutine write_grid(rlat, rlon, n, vclat, vclon, nvc, alist, clist) write(6,*) ' number of obtuse triangles ', nObtuse - kiteAreasOnVertex(:,:) = -1.0 + kiteAreasOnVertex(:,:) = -1.0_RKIND ! ! Compute weights used in tangential velocity reconstruction @@ -544,8 +544,8 @@ subroutine write_grid(rlat, rlon, n, vclat, vclon, nvc, alist, clist) edgeGP_next % lon = lonEdge(next_edge) if(isObtuse(vtx).eq.0) then - area = abs(triangle_area(cell1GP, vertex1GP, edgeGP_prev, 1.0)) - area = area + abs(triangle_area(cell1GP, vertex1GP, edgeGP_next, 1.0)) + area = abs(triangle_area(cell1GP, vertex1GP, edgeGP_prev, 1.0_RKIND)) + area = area + abs(triangle_area(cell1GP, vertex1GP, edgeGP_next, 1.0_RKIND)) else if(cellsOnVertex(isObtuse(vtx),vtx).eq.cellsOnEdge(1,j)) then iFlag = 0 @@ -563,10 +563,10 @@ subroutine write_grid(rlat, rlon, n, vclat, vclon, nvc, alist, clist) enddo edgeGP % lat = latEdge(iEdge) edgeGP % lon = lonEdge(iEdge) - area = abs(triangle_area(cell1GP, edgeGP, edgeGP_prev, 1.0)) - area = area + abs(triangle_area(cell1GP, edgeGP, edgeGP_next, 1.0)) + area = abs(triangle_area(cell1GP, edgeGP, edgeGP_prev, 1.0_RKIND)) + area = area + abs(triangle_area(cell1GP, edgeGP, edgeGP_next, 1.0_RKIND)) else - area = abs(triangle_area(cell1GP, edgeGP_prev, edgeGP_next, 1.0)) + area = abs(triangle_area(cell1GP, edgeGP_prev, edgeGP_next, 1.0_RKIND)) endif endif @@ -582,9 +582,9 @@ subroutine write_grid(rlat, rlon, n, vclat, vclon, nvc, alist, clist) r = area / areaCell(cellsOnEdge(1,j)) sum_r = sum_r + r if (cellsOnEdge(1,j) == cellsOnEdge(1,edgesOnEdge(i,j))) then - s = 1.0 + s = 1.0_RKIND else - s = -1.0 + s = -1.0_RKIND end if weightsOnEdge(i,j) = s*(0.5-sum_r)*dvEdge(edgesOnEdge(i,j))/de endif @@ -621,8 +621,8 @@ subroutine write_grid(rlat, rlon, n, vclat, vclon, nvc, alist, clist) edgeGP_next % lon = lonEdge(next_edge) if(isObtuse(vtx).eq.0) then - area = abs(triangle_area(cell1GP, vertex1GP, edgeGP_prev, 1.0)) - area = area + abs(triangle_area(cell1GP, vertex1GP, edgeGP_next, 1.0)) + area = abs(triangle_area(cell1GP, vertex1GP, edgeGP_prev, 1.0_RKIND)) + area = area + abs(triangle_area(cell1GP, vertex1GP, edgeGP_next, 1.0_RKIND)) else if(cellsOnVertex(isObtuse(vtx),vtx).eq.cellsOnEdge(2,j)) then iFlag = 0 @@ -640,10 +640,10 @@ subroutine write_grid(rlat, rlon, n, vclat, vclon, nvc, alist, clist) enddo edgeGP % lat = latEdge(iEdge) edgeGP % lon = lonEdge(iEdge) - area = abs(triangle_area(cell1GP, edgeGP, edgeGP_prev, 1.0)) - area = area + abs(triangle_area(cell1GP, edgeGP, edgeGP_next, 1.0)) + area = abs(triangle_area(cell1GP, edgeGP, edgeGP_prev, 1.0_RKIND)) + area = area + abs(triangle_area(cell1GP, edgeGP, edgeGP_next, 1.0_RKIND)) else - area = abs(triangle_area(cell1GP, edgeGP_prev, edgeGP_next, 1.0)) + area = abs(triangle_area(cell1GP, edgeGP_prev, edgeGP_next, 1.0_RKIND)) endif endif @@ -659,9 +659,9 @@ subroutine write_grid(rlat, rlon, n, vclat, vclon, nvc, alist, clist) r = area / areaCell(cellsOnEdge(2,j)) sum_r = sum_r + r if (cellsOnEdge(2,j) == cellsOnEdge(1,edgesOnEdge(i,j))) then - s = -1.0 + s = -1.0_RKIND else - s = 1.0 + s = 1.0_RKIND end if weightsOnEdge(i,j) = s*(0.5-sum_r)*dvEdge(edgesOnEdge(i,j))/de endif @@ -712,9 +712,9 @@ subroutine write_grid(rlat, rlon, n, vclat, vclon, nvc, alist, clist) r = area / areaCell(cellsOnEdge(1,j)) sum_r = sum_r + r if (cellsOnEdge(1,j) == cellsOnEdge(1,edgesOnEdge(i,j))) then - s = 1.0 + s = 1.0_RKIND else - s = -1.0 + s = -1.0_RKIND end if weightsOnEdge(i,j) = s*(0.5-sum_r)*dvEdge(edgesOnEdge(i,j))/de prev_edge = next_edge @@ -744,9 +744,9 @@ subroutine write_grid(rlat, rlon, n, vclat, vclon, nvc, alist, clist) r = area / areaCell(cellsOnEdge(2,j)) sum_r = sum_r + r if (cellsOnEdge(2,j) == cellsOnEdge(1,edgesOnEdge(i,j))) then - s = -1.0 + s = -1.0_RKIND else - s = 1.0 + s = 1.0_RKIND end if weightsOnEdge(i,j) = s*(0.5-sum_r)*dvEdge(edgesOnEdge(i,j))/de prev_edge = next_edge @@ -767,8 +767,8 @@ subroutine write_grid(rlat, rlon, n, vclat, vclon, nvc, alist, clist) vertex1GP % lon = lonVertex(verticesOnEdge(1,i)) vertex2GP % lat = latVertex(verticesOnEdge(2,i)) vertex2GP % lon = lonVertex(verticesOnEdge(2,i)) - dv1Edge(i) = sphere_distance(edgeGP, vertex1GP, 1.0) - dv2Edge(i) = sphere_distance(edgeGP, vertex2GP, 1.0) + dv1Edge(i) = sphere_distance(edgeGP, vertex1GP, 1.0_RKIND) + dv2Edge(i) = sphere_distance(edgeGP, vertex2GP, 1.0_RKIND) end do @@ -782,17 +782,17 @@ subroutine write_grid(rlat, rlon, n, vclat, vclon, nvc, alist, clist) do i=1,nCells cell1GP % lat = latCell(i) cell1GP % lon = lonCell(i) - call convert_lx(xCell(i), yCell(i), zCell(i), 1.0, cell1GP) + call convert_lx(xCell(i), yCell(i), zCell(i), 1.0_RKIND, cell1GP) end do do i=1,nVertices vertex1GP % lat = latVertex(i) vertex1GP % lon = lonVertex(i) - call convert_lx(xVertex(i), yVertex(i), zVertex(i), 1.0, vertex1GP) + call convert_lx(xVertex(i), yVertex(i), zVertex(i), 1.0_RKIND, vertex1GP) end do do i=1,nEdges edgeGP % lat = latEdge(i) edgeGP % lon = lonEdge(i) - call convert_lx(xEdge(i), yEdge(i), zEdge(i), 1.0, edgeGP) + call convert_lx(xEdge(i), yEdge(i), zEdge(i), 1.0_RKIND, edgeGP) end do @@ -1034,10 +1034,10 @@ subroutine insert_edge_to_tree(cellID, vertex1ID, vertex2ID, lat1, lon1, lat2, l integer, intent(in) :: cellID integer, intent(in) :: vertex1ID, vertex2ID - real, intent(in) :: lat1, lon1, lat2, lon2 + real(kind=RKIND), intent(in) :: lat1, lon1, lat2, lon2 type (binary_tree), pointer :: root - real :: tLat1, tLat2, tLon1, tLon2 + real(kind=RKIND) :: tLat1, tLat2, tLon1, tLon2 integer :: tID1, tID2 logical :: found type (binary_tree), pointer :: pre_cursor, cursor @@ -1122,7 +1122,7 @@ integer function point_compare(lat1, lon1, lat2, lon2) implicit none - real, intent(in) :: lat1, lon1, lat2, lon2 + real(kind=RKIND), intent(in) :: lat1, lon1, lat2, lon2 point_compare = -1 @@ -1147,7 +1147,7 @@ integer function edge_compare(aLat1, aLon1, aLat2, aLon2, bLat1, bLon1, bLat2, b implicit none - real, intent(in) :: aLat1, aLon1, aLat2, aLon2, bLat1, bLon1, bLat2, bLon2 + real(kind=RKIND), intent(in) :: aLat1, aLon1, aLat2, aLon2, bLat1, bLon1, bLat2, bLon2 edge_compare = -1 @@ -1197,18 +1197,18 @@ logical function is_flipped_vertex_order(latCell1, lonCell1, & implicit none - real, intent(in) :: latCell1, lonCell1, & + real(kind=RKIND), intent(in) :: latCell1, lonCell1, & latCell2, lonCell2, & latVertex1, lonVertex1, & latVertex2, lonVertex2 - real :: xCell1, yCell1, zCell1 - real :: xCell2, yCell2, zCell2 - real :: xVertex1, yVertex1, zVertex1 - real :: xVertex2, yVertex2, zVertex2 - real :: xV1, yV1, zV1 - real :: xV2, yV2, zV2 - real :: ci, cj, ck + real(kind=RKIND) :: xCell1, yCell1, zCell1 + real(kind=RKIND) :: xCell2, yCell2, zCell2 + real(kind=RKIND) :: xVertex1, yVertex1, zVertex1 + real(kind=RKIND) :: xVertex2, yVertex2, zVertex2 + real(kind=RKIND) :: xV1, yV1, zV1 + real(kind=RKIND) :: xV2, yV2, zV2 + real(kind=RKIND) :: ci, cj, ck type (geo_point) :: cell1, cell2, vertex1, vertex2 cell1 % lat = latCell1 @@ -1220,10 +1220,10 @@ logical function is_flipped_vertex_order(latCell1, lonCell1, & vertex2 % lat = latVertex2 vertex2 % lon = lonVertex2 - call convert_lx(xCell1, yCell1, zCell1, 1.0, cell1) - call convert_lx(xCell2, yCell2, zCell2, 1.0, cell2) - call convert_lx(xVertex1, yVertex1, zVertex1, 1.0, vertex1) - call convert_lx(xVertex2, yVertex2, zVertex2, 1.0, vertex2) + call convert_lx(xCell1, yCell1, zCell1, 1.0_RKIND, cell1) + call convert_lx(xCell2, yCell2, zCell2, 1.0_RKIND, cell2) + call convert_lx(xVertex1, yVertex1, zVertex1, 1.0_RKIND, vertex1) + call convert_lx(xVertex2, yVertex2, zVertex2, 1.0_RKIND, vertex2) xV1 = xCell2 - xCell1 yV1 = yCell2 - yCell1 @@ -1259,14 +1259,14 @@ logical function is_flipped_vertex_order2(edge, cell2, vertex2) type (geo_point), intent(in) :: edge, cell2, vertex2 - real :: xEdge, yEdge, zEdge - real :: xCell2, yCell2, zCell2 - real :: xVertex2, yVertex2, zVertex2 - real :: angle + real(kind=RKIND) :: xEdge, yEdge, zEdge + real(kind=RKIND) :: xCell2, yCell2, zCell2 + real(kind=RKIND) :: xVertex2, yVertex2, zVertex2 + real(kind=RKIND) :: angle - call convert_lx(xEdge, yEdge, zEdge, 1.0, edge) - call convert_lx(xCell2, yCell2, zCell2, 1.0, cell2) - call convert_lx(xVertex2, yVertex2, zVertex2, 1.0, vertex2) + call convert_lx(xEdge, yEdge, zEdge, 1.0_RKIND, edge) + call convert_lx(xCell2, yCell2, zCell2, 1.0_RKIND, cell2) + call convert_lx(xVertex2, yVertex2, zVertex2, 1.0_RKIND, vertex2) angle = plane_angle(xEdge, yEdge, zEdge, & xCell2, yCell2, zCell2, & @@ -1305,20 +1305,20 @@ subroutine order_points_ccw(center, npts, points, permutation) integer :: i, j integer :: itemp - real :: rtemp - real :: nx, ny, nz - real :: px, py, pz - real :: p0x, p0y, p0z - real, dimension(npts) :: angle + real(kind=RKIND) :: rtemp + real(kind=RKIND) :: nx, ny, nz + real(kind=RKIND) :: px, py, pz + real(kind=RKIND) :: p0x, p0y, p0z + real(kind=RKIND), dimension(npts) :: angle type (geo_point) :: ptemp - call convert_lx(nx, ny, nz, 1.0, center) - call convert_lx(p0x, p0y, p0z, 1.0, points(1)) + call convert_lx(nx, ny, nz, 1.0_RKIND, center) + call convert_lx(p0x, p0y, p0z, 1.0_RKIND, points(1)) angle(1) = 0.0 do i=2,npts - call convert_lx(px, py, pz, 1.0, points(i)) + call convert_lx(px, py, pz, 1.0_RKIND, points(i)) angle(i) = plane_angle(nx, ny, nz, p0x, p0y, p0z, px, py, pz, nx, ny, nz) if (angle(i) < 0.0) angle(i) = angle(i) + 2.0*pii if (angle(i) > 2.0*pii) angle(i) = angle(i) - 2.0*pii @@ -1358,15 +1358,15 @@ subroutine write_OpenDX( nCells, & integer, intent(in) :: nCells integer, intent(in) :: nVertices - real (kind=RKIND), dimension(:), intent(in) :: xCell - real (kind=RKIND), dimension(:), intent(in) :: yCell - real (kind=RKIND), dimension(:), intent(in) :: zCell - real (kind=RKIND), dimension(:), intent(in) :: xVertex - real (kind=RKIND), dimension(:), intent(in) :: yVertex - real (kind=RKIND), dimension(:), intent(in) :: zVertex + real(kind=RKIND), dimension(:), intent(in) :: xCell + real(kind=RKIND), dimension(:), intent(in) :: yCell + real(kind=RKIND), dimension(:), intent(in) :: zCell + real(kind=RKIND), dimension(:), intent(in) :: xVertex + real(kind=RKIND), dimension(:), intent(in) :: yVertex + real(kind=RKIND), dimension(:), intent(in) :: zVertex integer, dimension(:), intent(in) :: nEdgesOnCell integer, dimension(:,:), intent(in) :: verticesOnCell - real (kind=RKIND), dimension(:), intent(in) :: areaCell + real(kind=RKIND), dimension(:), intent(in) :: areaCell character(len=80) :: a, b, c, d, e, f integer :: i, j, k, nVerticesTotal, iEdge, iLoop diff --git a/grid_gen/global_scvt/src/module_scvt.F b/grid_gen/global_scvt/src/module_scvt.F index e414d56eb..167db4005 100644 --- a/grid_gen/global_scvt/src/module_scvt.F +++ b/grid_gen/global_scvt/src/module_scvt.F @@ -21,21 +21,21 @@ subroutine scvt_solve(n, lend, rlat, rlon, nvc, list, lptr, fn) integer, intent(in) :: n, nvc, fn integer, dimension(n), intent(inout) :: lend integer, dimension(nvc), intent(inout) :: list, lptr - real, dimension(n), intent(inout) :: rlat, rlon + real(kind=RKIND), dimension(n), intent(inout) :: rlat, rlon integer :: maxitr integer :: i, j, k, iter integer :: ntmax, nrow, nptri integer, allocatable, dimension(:) :: listc - real :: area, density, tot_mass - real :: x, y, z, new_ctr_x, new_ctr_y, new_ctr_z - real, allocatable, dimension(:) :: vclat, vclon - real, allocatable, dimension(:) :: rlat_2, rlon_2 + real(kind=RKIND) :: area, density, tot_mass + real(kind=RKIND) :: x, y, z, new_ctr_x, new_ctr_y, new_ctr_z + real(kind=RKIND), allocatable, dimension(:) :: vclat, vclon + real(kind=RKIND), allocatable, dimension(:) :: rlat_2, rlon_2 type (geo_point) :: p1, p2, p3, pc type (geo_point) :: p_n1, p_n2 type (geo_point), dimension(3,64) :: ptri - real :: avg_movement, maxmovement, movement + real(kind=RKIND) :: avg_movement, maxmovement, movement logical converged maxitr = n_scvt_iterations @@ -95,14 +95,14 @@ subroutine scvt_solve(n, lend, rlat, rlon, nvc, list, lptr, fn) call divide_triangle(p1, p2, p3, nptri, ptri) do j=1,nptri - area = triangle_area(ptri(1,j), ptri(2,j), ptri(3,j), 1.0) + area = triangle_area(ptri(1,j), ptri(2,j), ptri(3,j), 1.0_RKIND) call center_of_mass(ptri(1,j), ptri(2,j), ptri(3,j), pc) if (p1%lon - pc%lon > pii) pc%lon = pc%lon + 2.0*pii if (p1%lon - pc%lon < -pii) pc%lon = pc%lon - 2.0*pii density = density_for_point(pc) tot_mass = tot_mass + area * density - call convert_lx(x, y, z, 1.0, pc) + call convert_lx(x, y, z, 1.0_RKIND, pc) new_ctr_x = new_ctr_x + x*area*density new_ctr_y = new_ctr_y + y*area*density new_ctr_z = new_ctr_z + z*area*density @@ -120,14 +120,14 @@ subroutine scvt_solve(n, lend, rlat, rlon, nvc, list, lptr, fn) call divide_triangle(p1, p2, p3, nptri, ptri) do j=1,nptri - area = triangle_area(ptri(1,j), ptri(2,j), ptri(3,j), 1.0) + area = triangle_area(ptri(1,j), ptri(2,j), ptri(3,j), 1.0_RKIND) call center_of_mass(ptri(1,j), ptri(2,j), ptri(3,j), pc) if (p1%lon - pc%lon > pii) pc%lon = pc%lon + 2.0*pii if (p1%lon - pc%lon < -pii) pc%lon = pc%lon - 2.0*pii density = density_for_point(pc) tot_mass = tot_mass + area * density - call convert_lx(x, y, z, 1.0, pc) + call convert_lx(x, y, z, 1.0_RKIND, pc) new_ctr_x = new_ctr_x + x*area*density new_ctr_y = new_ctr_y + y*area*density new_ctr_z = new_ctr_z + z*area*density @@ -156,8 +156,8 @@ subroutine scvt_solve(n, lend, rlat, rlon, nvc, list, lptr, fn) p_n2%lat = rlat_2(i) p_n2%lon = rlon_2(i) - call convert_lx(x,y,z,1.0,p_n1) - call convert_lx(new_ctr_x, new_ctr_y, new_ctr_z,1.0,p_n2) + call convert_lx(x,y,z,1.0_RKIND,p_n1) + call convert_lx(new_ctr_x, new_ctr_y, new_ctr_z,1.0_RKIND,p_n2) !x y z computation movement = sqrt((x - new_ctr_x)**2 + (y - new_ctr_y)**2 + (z - new_ctr_z)**2) @@ -195,27 +195,27 @@ end subroutine scvt_solve subroutine random_point(p) type (geo_point), intent(inout) :: p - real :: x, y, z, m - real :: pi + real(kind=RKIND) :: x, y, z, m + real(kind=RKIND) :: pi - pi = 4.0*atan(1.0) + pi = 4.0*atan(1.0_RKIND) x = 0.0 y = 0.0 z = 0.0 m = 2.0 - do while (m > 1.0 .or. (x == 0.0 .and. y == 0.0 .and. z == 0.0)) + do while (m > 1.0_RKIND .or. (x == 0.0 .and. y == 0.0 .and. z == 0.0)) call random_number(x) call random_number(y) call random_number(z) - x = x * 2.0 - 1.0 - y = y * 2.0 - 1.0 - z = z * 2.0 - 1.0 + x = x * 2.0 - 1.0_RKIND + y = y * 2.0 - 1.0_RKIND + z = z * 2.0 - 1.0_RKIND m = x**2 + y**2 + z**2 end do - m = 1.0 / sqrt(m) + m = 1.0_RKIND / sqrt(m) x = x * m y = y * m z = z * m @@ -229,27 +229,27 @@ end subroutine random_point ! FUNCTION DENSITY_FOR_POINT ! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! - real function density_for_point(p) + real(kind=RKIND) function density_for_point(p) implicit none type (geo_point), intent(in) :: p character (len=256) :: fname - real :: rx, ry, rz, prx, pry, prz + real(kind=RKIND) :: rx, ry, rz, prx, pry, prz type (geo_point) :: p_local - real :: hgt - real :: r, norm, t_cent - real :: r1 - real :: pi - real :: width, trans_center, min_val + real(kind=RKIND) :: hgt + real(kind=RKIND) :: r, norm, t_cent + real(kind=RKIND) :: r1 + real(kind=RKIND) :: pi + real(kind=RKIND) :: width, trans_center, min_val - pi = 4.0*atan(1.0) + pi = 4.0*atan(1.0_RKIND) - !density_for_point = 1.0 + (1.19*cos(p%lat-3.141592654/4.0))**16.0 + !density_for_point = 1.0_RKIND + (1.19*cos(p%lat-3.141592654/4.0))**16.0 ! Uniform Density Function - density_for_point = 1.0 + density_for_point = 1.0_RKIND !Target Density Function based on hyperbolic tangent ! p_local%lat = latitude (radians) center of high-resolution region @@ -257,19 +257,19 @@ real function density_for_point(p) ! width = width of transition zone ! trans_center = width (radians) of high resolution zone ! minval = minimum density value. to have grid spacing vary by a factor of 8 - ! set minval = (1.0 / 8.0)**4 + ! set minval = (1.0_RKIND / 8.0_RKIND)**4 ! p_local%lat = pii/4.0 ! p_local%lon = 1.25*pii - ! call convert_lx(rx, ry, rz, 1.0, p) - ! call convert_lx(prx, pry, prz, 1.0, p_local) + ! call convert_lx(rx, ry, rz, 1.0_RKIND, p) + ! call convert_lx(prx, pry, prz, 1.0_RKIND, p_local) ! r = acos(rx*prx + ry*pry + rz*prz) ! width = 0.15 ! trans_center = pi/6.0 - ! min_val = (1.0/8.0)**4 - ! norm = 1.0/(1.0-min_val) - ! density_for_point = ((tanh((trans_center-r)*(1.0/width))+1.0)/2)/norm + min_val + ! min_val = (1.0_RKIND/8.0_RKIND)**4 + ! norm = 1.0_RKIND/(1.0_RKIND-min_val) + ! density_for_point = ((tanh((trans_center-r)*(1.0_RKIND/width))+1.0_RKIND)/2)/norm + min_val end function density_for_point diff --git a/grid_gen/global_scvt/src/module_sphere_utilities.F b/grid_gen/global_scvt/src/module_sphere_utilities.F index f0026bec7..ad19f9c83 100644 --- a/grid_gen/global_scvt/src/module_sphere_utilities.F +++ b/grid_gen/global_scvt/src/module_sphere_utilities.F @@ -8,18 +8,18 @@ module sphere_utilities ! Given the (latitude, longitude) coordinates of the corners of a triangle, ! plus the radius of the sphere, compute the area of the triangle. !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! -real function triangle_area(p1, p2, p3, radius) +real(kind=RKIND) function triangle_area(p1, p2, p3, radius) use data_types implicit none type (geo_point), intent(in) :: p1, p2, p3 - real, intent(in) :: radius + real(kind=RKIND), intent(in) :: radius - real :: a, b, c, s, e, pii, tanqe + real(kind=RKIND) :: a, b, c, s, e, pii, tanqe - pii = 2.*asin(1.0) + pii = 2.*asin(1.0_RKIND) a = sphere_distance(p1,p2,radius) b = sphere_distance(p2,p3,radius) @@ -51,13 +51,13 @@ integer function obtuse(p1, p2, p3) type (geo_point), intent(in) :: p1, p2, p3 - real :: x1(3), x2(3), x3(3), dot, r(3), s(3), rmag, smag + real(kind=RKIND) :: x1(3), x2(3), x3(3), dot, r(3), s(3), rmag, smag obtuse = 0 - call convert_lx(x1(1), x1(2), x1(3), 1.0, p1) - call convert_lx(x2(1), x2(2), x2(3), 1.0, p2) - call convert_lx(x3(1), x3(2), x3(3), 1.0, p3) + call convert_lx(x1(1), x1(2), x1(3), 1.0_RKIND, p1) + call convert_lx(x2(1), x2(2), x2(3), 1.0_RKIND, p2) + call convert_lx(x3(1), x3(2), x3(3), 1.0_RKIND, p3) ! test angle formed by x3,x1,x2 r(:) = x3(:) - x1(:) @@ -99,16 +99,16 @@ end function obtuse ! plus the radius of the sphere, compute the great circle distance between ! the points. !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! -real function sphere_distance(p1, p2, radius) +real(kind=RKIND) function sphere_distance(p1, p2, radius) use data_types implicit none type (geo_point), intent(in) :: p1, p2 - real, intent(in) :: radius + real(kind=RKIND), intent(in) :: radius - real :: arg1 + real(kind=RKIND) :: arg1 arg1 = sqrt( sin(0.5*(p2%lat-p1%lat))**2 + & cos(p1%lat)*cos(p2%lat)*sin(0.5*(p2%lon-p1%lon))**2 ) @@ -124,16 +124,16 @@ end function sphere_distance ! plus the radius of the sphere, compute the secant distance between ! the points. !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! -real function plane_distance(p1, p2, radius) +real(kind=RKIND) function plane_distance(p1, p2, radius) use data_types implicit none type (geo_point), intent(in) :: p1, p2 - real, intent(in) :: radius + real(kind=RKIND), intent(in) :: radius - real :: x1, x2, y1, y2, z1, z2 + real(kind=RKIND) :: x1, x2, y1, y2, z1, z2 z1 = sin(p1%lat) z2 = sin(p2%lat) @@ -154,7 +154,7 @@ end function plane_distance ! compute the angle between the points as measured from the origin of the ! sphere. !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! -real function arc_angle(p1, p2) +real(kind=RKIND) function arc_angle(p1, p2) use data_types @@ -162,7 +162,7 @@ real function arc_angle(p1, p2) type (geo_point), intent(in) :: p1, p2 - real :: arg1 + real(kind=RKIND) :: arg1 arg1 = sqrt( sin(0.5*(p2%lat-p1%lat))**2 + & cos(p1%lat)*cos(p2%lat)*sin(0.5*(p2%lon-p1%lon))**2 ) @@ -187,15 +187,15 @@ subroutine great_circle_points(p1, p2, pl, n) type (geo_point), intent(in) :: p1, p2 type (geo_point), dimension(n), intent(inout) :: pl - real :: x1, x2, y1, y2, z1, z2 - real :: dx, dl - real :: x, y, z + real(kind=RKIND) :: x1, x2, y1, y2, z1, z2 + real(kind=RKIND) :: dx, dl + real(kind=RKIND) :: x, y, z integer :: i - real :: dtheta, dinc, dt + real(kind=RKIND) :: dtheta, dinc, dt - real :: pii, rtod + real(kind=RKIND) :: pii, rtod - pii = 2.*asin(1.0) + pii = 2.*asin(1.0_RKIND) rtod = 180./pii ! write(6,*) ' in gcp ',rtod*lat1,rtod*lon1,rtod*lat2,rtod*lon2 @@ -213,8 +213,8 @@ subroutine great_circle_points(p1, p2, pl, n) dtheta = arc_angle(p1, p2) dinc = dtheta/float(n-1) - call convert_lx(x1,y1,z1,1.,p1) - call convert_lx(x2,y2,z2,1.,p2) + call convert_lx(x1,y1,z1,1.0_RKIND,p1) + call convert_lx(x2,y2,z2,1.0_RKIND,p2) ! set the end points @@ -230,14 +230,14 @@ subroutine great_circle_points(p1, p2, pl, n) dt = float(i-1)*dinc if (dt <= 0.5*dtheta) then - dx = 1.-tan(0.5*dtheta-dt)/tan(0.5*dtheta) + dx = 1._RKIND-tan(0.5*dtheta-dt)/tan(0.5*dtheta) ! write(6,*) ' case 1 ',dx x = x1+0.5*dx*(x2-x1) y = y1+0.5*dx*(y2-y1) z = z1+0.5*dx*(z2-z1) else dt = dtheta-dt - dx = 1.-tan(0.5*dtheta-dt)/tan(0.5*dtheta) + dx = 1._RKIND-tan(0.5*dtheta-dt)/tan(0.5*dtheta) ! write(6,*) ' case 2 ',dx x = x2+0.5*dx*(x1-x2) y = y2+0.5*dx*(y1-y2) @@ -265,7 +265,7 @@ end subroutine great_circle_points ! type (geo_point), intent(in) :: p1, p2, p3 ! type (geo_point), dimension(6), intent(inout) :: pnew ! -! real :: t_area, area_total, radius +! real(kind=RKIND) :: t_area, area_total, radius ! type (geo_point), dimension(3) :: pts ! type (geo_point) :: c ! @@ -349,9 +349,9 @@ end subroutine great_circle_points ! type (geo_point), intent(in) :: p0, p1, p2 ! type (geo_point), intent(out) :: vc ! -! real :: x0, y0, z0, x1, y1, z1, x2, y2, z2, xc, yc, zc, cabs -! real :: a1, a2, a3, b1, b2, b3 -! real :: dot0 +! real(kind=RKIND) :: x0, y0, z0, x1, y1, z1, x2, y2, z2, xc, yc, zc, cabs +! real(kind=RKIND) :: a1, a2, a3, b1, b2, b3 +! real(kind=RKIND) :: dot0 ! ! z0 = sin(p0%lat) ! z1 = sin(p1%lat) @@ -441,9 +441,9 @@ subroutine convert_lx(x, y, z, radius, latlon) implicit none - real, intent(in) :: radius + real(kind=RKIND), intent(in) :: radius type (geo_point), intent(in) :: latlon - real, intent(out) :: x, y, z + real(kind=RKIND), intent(out) :: x, y, z z = radius * sin(latlon%lat) x = radius * cos(latlon%lon) * cos(latlon%lat) @@ -464,15 +464,15 @@ subroutine convert_xl(x, y, z, latlon) implicit none - real, intent(in) :: x, y, z + real(kind=RKIND), intent(in) :: x, y, z type (geo_point), intent(out) :: latlon - real :: dl, clat, pii, rtod - real :: eps + real(kind=RKIND) :: dl, clat, pii, rtod + real(kind=RKIND) :: eps parameter (eps=1.e-10) - pii = 2.*asin(1.0) - rtod=180./pii + pii = 2.*asin(1.0_RKIND) + rtod=180.0_RKIND/pii dl = sqrt(x*x + y*y + z*z) latlon%lat = asin(z/dl) @@ -533,10 +533,10 @@ subroutine gc_intersect(p0, p1, p2, p3, pc) type (geo_point), intent(in) :: p0, p1, p2, p3 type (geo_point), intent(out) :: pc - real :: x0, y0, z0, x1, y1, z1, x2, y2, z2, x3, y3, z3 - real :: n1, n2, n3, m1, m2, m3 - real :: xc, yc, zc, dot - real, parameter :: radius=1.0 + real(kind=RKIND) :: x0, y0, z0, x1, y1, z1, x2, y2, z2, x3, y3, z3 + real(kind=RKIND) :: n1, n2, n3, m1, m2, m3 + real(kind=RKIND) :: xc, yc, zc, dot + real(kind=RKIND), parameter :: radius=1.0_RKIND call convert_lx(x0,y0,z0,radius,p0) call convert_lx(x1,y1,z1,radius,p1) @@ -573,15 +573,15 @@ end subroutine gc_intersect ! ! Normalize an angle, given in radians, to lie in the interval [0,2*PI]. !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! -real function pos_ang(angle) +real(kind=RKIND) function pos_ang(angle) implicit none - real, intent(in) :: angle + real(kind=RKIND), intent(in) :: angle - real :: pii + real(kind=RKIND) :: pii - pii = 2.*asin(1.0) + pii = 2.*asin(1.0_RKIND) pos_ang = angle if(angle > 2.*pii) then @@ -603,7 +603,7 @@ end function pos_ang ! Convention: zero points north, 90 points west, -90 point east, ! points south 180, -180 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! -real function meridian_angle(p1, p2) +real(kind=RKIND) function meridian_angle(p1, p2) use data_types @@ -613,21 +613,21 @@ real function meridian_angle(p1, p2) type (geo_point) :: np - real :: pii, da, db, dc + real(kind=RKIND) :: pii, da, db, dc type (geo_point) :: p3 - real :: cosa - real :: eps + real(kind=RKIND) :: cosa + real(kind=RKIND) :: eps parameter (eps = 1.e-04) -real :: ax, ay, az -real :: bx, by, bz -real :: cx, cy, cz +real(kind=RKIND) :: ax, ay, az +real(kind=RKIND) :: bx, by, bz +real(kind=RKIND) :: cx, cy, cz np = p1 np%lat = np%lat + 0.05 -call convert_lx(ax, ay, az, 1.0, p1) -call convert_lx(bx, by, bz, 1.0, np) -call convert_lx(cx, cy, cz, 1.0, p2) +call convert_lx(ax, ay, az, 1.0_RKIND, p1) +call convert_lx(bx, by, bz, 1.0_RKIND, np) +call convert_lx(cx, cy, cz, 1.0_RKIND, p2) meridian_angle = plane_angle(ax, ay, az, bx, by, bz, cx, cy, cz, ax, ay, az) return @@ -638,7 +638,7 @@ real function meridian_angle(p1, p2) else - pii = 2.*asin(1.0) + pii = 2.*asin(1.0_RKIND) dc = arc_angle(p1,p2) p3%lon = p1%lon @@ -683,13 +683,13 @@ subroutine center_of_mass(p1, p2, p3, pc) type (geo_point), intent(in) :: p1, p2, p3 type (geo_point), intent(out) :: pc - real :: x1, x2, x3, xc - real :: y1, y2, y3, yc - real :: z1, z2, z3, zc + real(kind=RKIND) :: x1, x2, x3, xc + real(kind=RKIND) :: y1, y2, y3, yc + real(kind=RKIND) :: z1, z2, z3, zc - call convert_lx(x1,y1,z1,1.,p1) - call convert_lx(x2,y2,z2,1.,p2) - call convert_lx(x3,y3,z3,1.,p3) + call convert_lx(x1,y1,z1,1.0_RKIND,p1) + call convert_lx(x2,y2,z2,1.0_RKIND,p2) + call convert_lx(x3,y3,z3,1.0_RKIND,p3) xc = (x1+x2+x3)/3. yc = (y1+y2+y3)/3. @@ -784,14 +784,14 @@ subroutine point_to_plane(a, b, c, d, Qx, Qy, Qz, xp, yp, zp) implicit none - real, intent(in) :: a, b, c, d ! The coefficients in the equation of the plane - real, intent(in) :: Qx, Qy, Qz ! The coordinates of the point Q to be projected to the plane - real, intent(out) :: xp, yp, zp ! The coordinates of the point projected in the plane + real(kind=RKIND), intent(in) :: a, b, c, d ! The coefficients in the equation of the plane + real(kind=RKIND), intent(in) :: Qx, Qy, Qz ! The coordinates of the point Q to be projected to the plane + real(kind=RKIND), intent(out) :: xp, yp, zp ! The coordinates of the point projected in the plane - real :: Px, Py, Pz ! A point P in the plane ax + by + cz + d = 0 - real :: PQx, PQy, PQz ! Components of the vector from P to Q - real :: PQn ! The dot product of PQ and the vector normal to the plane - real :: m2 ! The magnitude and squared magnitude of the vector n normal to the plane + real(kind=RKIND) :: Px, Py, Pz ! A point P in the plane ax + by + cz + d = 0 + real(kind=RKIND) :: PQx, PQy, PQz ! Components of the vector from P to Q + real(kind=RKIND) :: PQn ! The dot product of PQ and the vector normal to the plane + real(kind=RKIND) :: m2 ! The magnitude and squared magnitude of the vector n normal to the plane m2 = (a**2.0 + b**2.0 + c**2.0) @@ -830,13 +830,13 @@ subroutine point_to_sphere(a, b, c, d, r, Qx, Qy, Qz, xp, yp, zp) implicit none - real, intent(in) :: a, b, c, d ! The coefficients in the equation of the plane - real, intent(in) :: r ! The radius of the sphere - real, intent(in) :: Qx, Qy, Qz ! The coordinates of the point Q to be projected to the sphere - real, intent(out) :: xp, yp, zp ! The coordinates of the point projected to the sphere + real(kind=RKIND), intent(in) :: a, b, c, d ! The coefficients in the equation of the plane + real(kind=RKIND), intent(in) :: r ! The radius of the sphere + real(kind=RKIND), intent(in) :: Qx, Qy, Qz ! The coordinates of the point Q to be projected to the sphere + real(kind=RKIND), intent(out) :: xp, yp, zp ! The coordinates of the point projected to the sphere - real :: aa, bb, cc ! Coefficients of quadratic equation - real :: disc, t1, t2 + real(kind=RKIND) :: aa, bb, cc ! Coefficients of quadratic equation + real(kind=RKIND) :: disc, t1, t2 ! Solve for the interesection of the line (Qx - at, Qy - bt, Qz - ct) and the ! sphere x^2 + y^2 + z^2 - r^2 = 0 @@ -884,11 +884,11 @@ subroutine rotate_about_vector(x, y, z, theta, a, b, c, u, v, w, xp, yp, zp) implicit none - real, intent(in) :: x, y, z, theta, a, b, c, u, v, w - real, intent(out) :: xp, yp, zp + real(kind=RKIND), intent(in) :: x, y, z, theta, a, b, c, u, v, w + real(kind=RKIND), intent(out) :: xp, yp, zp - real :: vw2, uw2, uv2 - real :: m + real(kind=RKIND) :: vw2, uw2, uv2 + real(kind=RKIND) :: m vw2 = v**2.0 + w**2.0 uw2 = u**2.0 + w**2.0 @@ -908,22 +908,22 @@ end subroutine rotate_about_vector ! Computes the angle between vectors AB and AC, given points A, B, and C, and ! a vector (u,v,w) normal to the plane. !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! -real function plane_angle(ax, ay, az, bx, by, bz, cx, cy, cz, u, v, w) +real(kind=RKIND) function plane_angle(ax, ay, az, bx, by, bz, cx, cy, cz, u, v, w) implicit none - real, intent(in) :: ax, ay, az, bx, by, bz, cx, cy, cz, u, v, w + real(kind=RKIND), intent(in) :: ax, ay, az, bx, by, bz, cx, cy, cz, u, v, w - real :: ABx, ABy, ABz ! The components of the vector AB - real :: mAB ! The magnitude of AB - real :: ACx, ACy, ACz ! The components of the vector AC - real :: mAC ! The magnitude of AC + real(kind=RKIND) :: ABx, ABy, ABz ! The components of the vector AB + real(kind=RKIND) :: mAB ! The magnitude of AB + real(kind=RKIND) :: ACx, ACy, ACz ! The components of the vector AC + real(kind=RKIND) :: mAC ! The magnitude of AC - real :: Dx ! The i-components of the cross product AB x AC - real :: Dy ! The j-components of the cross product AB x AC - real :: Dz ! The k-components of the cross product AB x AC + real(kind=RKIND) :: Dx ! The i-components of the cross product AB x AC + real(kind=RKIND) :: Dy ! The j-components of the cross product AB x AC + real(kind=RKIND) :: Dz ! The k-components of the cross product AB x AC - real :: cos_angle + real(kind=RKIND) :: cos_angle ABx = bx - ax ABy = by - ay @@ -942,10 +942,10 @@ real function plane_angle(ax, ay, az, bx, by, bz, cx, cy, cz, u, v, w) cos_angle = (ABx*ACx + ABy*ACy + ABz*ACz) / (mAB * mAC) - if (cos_angle < -1.0) then - cos_angle = -1.0 - else if (cos_angle > 1.0) then - cos_angle = 1.0 + if (cos_angle < -1.0_RKIND) then + cos_angle = -1.0_RKIND + else if (cos_angle > 1.0_RKIND) then + cos_angle = 1.0_RKIND end if if ((Dx*u + Dy*v + Dz*w) >= 0.0) then diff --git a/grid_gen/global_scvt/src/module_voronoi_utils.F b/grid_gen/global_scvt/src/module_voronoi_utils.F index 3079cd320..1da42e1e0 100644 --- a/grid_gen/global_scvt/src/module_voronoi_utils.F +++ b/grid_gen/global_scvt/src/module_voronoi_utils.F @@ -17,13 +17,13 @@ subroutine compute_dt(rlat, rlon, n, ltri, nrow, ntmx, nt) integer, intent(in) :: n, nrow, ntmx integer, intent(inout) :: nt integer, dimension(nrow, ntmx), intent(in) :: ltri - real, dimension(n), intent(in) :: rlat, rlon + real(kind=RKIND), dimension(n), intent(in) :: rlat, rlon integer :: ierr, lnew, nscr integer, dimension(n) :: near, next integer, dimension(n) :: lend integer, dimension(6*n+12) :: list, lptr - real, dimension(n) :: x, y, z, dist + real(kind=RKIND), dimension(n) :: x, y, z, dist nscr = 6*n+12 @@ -57,15 +57,15 @@ subroutine compute_vc(rlat, rlon, n, nrow, ntmx, list, lptr, lend, listc, vclat, integer, intent(in) :: n, nrow, ntmx, nvc integer, dimension(nvc), intent(inout) :: list, lptr, listc - real, dimension(nvc), intent(inout) :: vclat, vclon + real(kind=RKIND), dimension(nvc), intent(inout) :: vclat, vclon integer, dimension(n), intent(inout) :: lend - real, dimension(n), intent(in) :: rlat, rlon + real(kind=RKIND), dimension(n), intent(in) :: rlat, rlon integer :: ierr, lnew, nb integer, dimension(n) :: near, next integer, dimension(nrow, ntmx) :: ltri - real, dimension(n) :: x, y, z, dist - real, dimension(nvc) :: xc, yc, zc, rc + real(kind=RKIND), dimension(n) :: x, y, z, dist + real(kind=RKIND), dimension(nvc) :: xc, yc, zc, rc if (nvc < 6*n-12) then write(0,*) 'Error: Argument nvc to COMPUTE_VC must be at least 6*n+12' @@ -98,8 +98,8 @@ subroutine trans_inv(x, y, z, lat, lon, n) implicit none integer, intent(in) :: n - real, dimension(n), intent(in) :: x, y, z - real, dimension(n), intent(out) :: lat, lon + real(kind=RKIND), dimension(n), intent(in) :: x, y, z + real(kind=RKIND), dimension(n), intent(out) :: lat, lon integer :: i diff --git a/grid_gen/grid_rotate/Makefile b/grid_gen/grid_rotate/Makefile index 67520dae6..75b4850a9 100644 --- a/grid_gen/grid_rotate/Makefile +++ b/grid_gen/grid_rotate/Makefile @@ -1,12 +1,28 @@ FC = gfortran FFLAGS = -m64 -O2 -ffree-line-length-none -ffree-form -Wall -FCINCLUDES = $(shell nc-config --fflags) -FCLIBS = $(shell nc-config --flibs) + +ifneq ($(NETCDF), ) +ifneq ($(shell which ${NETCDF}/bin/nc-config 2> /dev/null), ) + LIBS = $(shell ${NETCDF}/bin/nc-config --libs) -lnetcdf_c++ -lnetcdff + INCS = $(shell ${NETCDF}/bin/nc-config --cflags) +else + LIBS= -L${NETCDF}/lib + LIBS += -lnetcdf_c++ -lnetcdf -lstdc++ -lnetcdff + INCS = -I${NETCDF}/include +endif +else ifneq ($(shell which nc-config 2> /dev/null), ) + LIBS = $(shell nc-config --libs) -lnetcdf_c++ + INCS = $(shell nc-config --cflags) +else + LIBS= -L${NETCDF}/lib + LIBS += -lnetcdf_c++ -lnetcdf -lstdc++ -lnetcdff + INCS = -I${NETCDF}/include +endif all: grid_rotate grid_rotate: grid_rotate.f90 - $(FC) grid_rotate.f90 -o grid_rotate $(FFLAGS) $(FCINCLUDES) $(FCLIBS) + $(FC) grid_rotate.f90 -o grid_rotate $(FFLAGS) $(INCS) $(LIBS) clean: rm grid_rotate diff --git a/grid_gen/mesh_conversion_tools/Makefile b/grid_gen/mesh_conversion_tools/Makefile index 186601a0a..1f86b268f 100644 --- a/grid_gen/mesh_conversion_tools/Makefile +++ b/grid_gen/mesh_conversion_tools/Makefile @@ -2,30 +2,45 @@ # gfortran CC=g++ -CFLAGS= -O3 -std=c++0x -fopenmp -DFLAGS= -g -std=c++0x -D_DEBUG -fopenmp +CFLAGS= -O3 -std=c++0x -fopenmp -lstdc++ +DFLAGS= -g -std=c++0x -D_DEBUG -fopenmp -lstdc++ # intel # CC=icpc -# CFLAGS= -O3 -std=c++0x -openmp -# DFLAGS= -g -std=c++0x -D_DEBUG -openmp +# CFLAGS= -O3 -std=c++0x -openmp -lstdc++ +# DFLAGS= -g -std=c++0x -D_DEBUG -openmp -lstdc++ CONV_EXECUTABLE= MpasMeshConverter.x CULL_EXECUTABLE= MpasCellCuller.x MASK_EXECUTABLE= MpasMaskCreator.x -NCINCDIR=${NETCDF}/include -NCLIBDIR=${NETCDF}/lib +ifneq ($(NETCDF), ) +ifneq ($(shell which ${NETCDF}/bin/nc-config 2> /dev/null), ) + LIBS = $(shell ${NETCDF}/bin/nc-config --libs) -lnetcdf_c++ -lnetcdff + INCS = $(shell ${NETCDF}/bin/nc-config --cflags) +else + LIBS= -L${NETCDF}/lib + LIBS += -lnetcdf_c++ -lnetcdf -lstdc++ -lnetcdff + INCS = -I${NETCDF}/include +endif +else ifneq ($(shell which nc-config 2> /dev/null), ) + LIBS = $(shell nc-config --libs) -lnetcdf_c++ + INCS = $(shell nc-config --cflags) +else + LIBS= -L${NETCDF}/lib + LIBS += -lnetcdf_c++ -lnetcdf -lstdc++ -lnetcdff + INCS = -I${NETCDF}/include +endif all: - ${CC} mpas_mesh_converter.cpp netcdf_utils.cpp ${CFLAGS} -o ${CONV_EXECUTABLE} -I${NCINCDIR} -L${NCLIBDIR} -lnetcdf_c++ -lnetcdf -lstdc++ - ${CC} mpas_cell_culler.cpp netcdf_utils.cpp ${CFLAGS} -o ${CULL_EXECUTABLE} -I${NCINCDIR} -L${NCLIBDIR} -lnetcdf_c++ -lnetcdf -lstdc++ - ${CC} mpas_mask_creator.cpp netcdf_utils.cpp jsoncpp.cpp ${CFLAGS} -o ${MASK_EXECUTABLE} -I. -I${NCINCDIR} -L${NCLIBDIR} -lnetcdf_c++ -lnetcdf -lstdc++ + ${CC} mpas_mask_creator.cpp netcdf_utils.cpp jsoncpp.cpp ${CFLAGS} -o ${MASK_EXECUTABLE} -I. ${INCS} ${LIBS} + ${CC} mpas_mesh_converter.cpp netcdf_utils.cpp ${CFLAGS} -o ${CONV_EXECUTABLE} ${INCS} ${LIBS} + ${CC} mpas_cell_culler.cpp netcdf_utils.cpp ${CFLAGS} -o ${CULL_EXECUTABLE} ${INCS} ${LIBS} debug: - ${CC} mpas_mesh_converter.cpp netcdf_utils.cpp ${DFLAGS} -o ${CONV_EXECUTABLE} -I${NCINCDIR} -L${NCLIBDIR} -lnetcdf_c++ -lnetcdf -lstdc++ - ${CC} mpas_cell_culler.cpp netcdf_utils.cpp ${DFLAGS} -o ${CULL_EXECUTABLE} -I${NCINCDIR} -L${NCLIBDIR} -lnetcdf_c++ -lnetcdf -lstdc++ - ${CC} mpas_mask_creator.cpp netcdf_utils.cpp jsoncpp.cpp ${DFLAGS} -o ${MASK_EXECUTABLE} -I. -I${NCINCDIR} -L${NCLIBDIR} -lnetcdf_c++ -lnetcdf -lstdc++ + ${CC} mpas_mesh_converter.cpp netcdf_utils.cpp ${DFLAGS} -o ${CONV_EXECUTABLE} ${INCS} ${LIBS} + ${CC} mpas_cell_culler.cpp netcdf_utils.cpp ${DFLAGS} -o ${CULL_EXECUTABLE} ${INCS} ${LIBS} + ${CC} mpas_mask_creator.cpp netcdf_utils.cpp jsoncpp.cpp ${DFLAGS} -o ${MASK_EXECUTABLE} -I. ${INCS} ${LIBS} clean: rm -f grid.nc diff --git a/grid_gen/mesh_conversion_tools/json/json.h b/grid_gen/mesh_conversion_tools/json/json.h index 42e7e7f4a..01225cac6 100644 --- a/grid_gen/mesh_conversion_tools/json/json.h +++ b/grid_gen/mesh_conversion_tools/json/json.h @@ -425,8 +425,8 @@ namespace Json { class JSON_API Exception : public std::exception { public: Exception(std::string const& msg); - ~Exception() throw() override; - char const* what() const throw() override; + ~Exception() throw() ; + char const* what() const throw() ; protected: std::string msg_; }; @@ -1581,9 +1581,9 @@ class JSON_API CharReaderBuilder : public CharReader::Factory { Json::Value settings_; CharReaderBuilder(); - ~CharReaderBuilder() override; + ~CharReaderBuilder() ; - CharReader* newCharReader() const override; + CharReader* newCharReader() const ; /** \return true if 'settings' are legal and consistent; * otherwise, indicate bad settings via 'invalid'. @@ -1778,12 +1778,12 @@ class JSON_API StreamWriterBuilder : public StreamWriter::Factory { Json::Value settings_; StreamWriterBuilder(); - ~StreamWriterBuilder() override; + ~StreamWriterBuilder() ; /** * \throw std::exception if something goes wrong (e.g. invalid settings) */ - StreamWriter* newStreamWriter() const override; + StreamWriter* newStreamWriter() const ; /** \return true if 'settings' are legal and consistent; * otherwise, indicate bad settings via 'invalid'. @@ -1824,7 +1824,7 @@ class JSON_API FastWriter : public Writer { public: FastWriter(); - ~FastWriter() override {} + ~FastWriter() {} void enableYAMLCompatibility(); @@ -1838,7 +1838,7 @@ class JSON_API FastWriter : public Writer { void omitEndingLineFeed(); public: // overridden from Writer - std::string write(const Value& root) override; + std::string write(const Value& root) ; private: void writeValue(const Value& value); @@ -1876,14 +1876,14 @@ class JSON_API FastWriter : public Writer { class JSON_API StyledWriter : public Writer { public: StyledWriter(); - ~StyledWriter() override {} + ~StyledWriter() {} public: // overridden from Writer /** \brief Serialize a Value in JSON format. * \param root Value to serialize. * \return String containing the JSON document that represents the root value. */ - std::string write(const Value& root) override; + std::string write(const Value& root) ; private: void writeValue(const Value& value); diff --git a/grid_gen/mesh_conversion_tools/jsoncpp.cpp b/grid_gen/mesh_conversion_tools/jsoncpp.cpp index ddb5137f5..985589c82 100644 --- a/grid_gen/mesh_conversion_tools/jsoncpp.cpp +++ b/grid_gen/mesh_conversion_tools/jsoncpp.cpp @@ -2088,7 +2088,7 @@ class OurCharReader : public CharReader { {} bool parse( char const* beginDoc, char const* endDoc, - Value* root, std::string* errs) override { + Value* root, std::string* errs) { bool ok = reader_.parse(beginDoc, endDoc, *root, collectComments_); if (errs) { *errs = reader_.getFormattedErrorMessages(); @@ -2652,7 +2652,7 @@ Value::CZString::CZString(const CZString& other) #if JSON_HAS_RVALUE_REFERENCES Value::CZString::CZString(CZString&& other) : cstr_(other.cstr_), index_(other.index_) { - other.cstr_ = nullptr; + other.cstr_ = NULL; } #endif @@ -4806,7 +4806,7 @@ struct BuiltStyledStreamWriter : public StreamWriter std::string const& endingLineFeedSymbol, bool useSpecialFloats, unsigned int precision); - int write(Value const& root, std::ostream* sout) override; + int write(Value const& root, std::ostream* sout) ; private: void writeValue(Value const& value); void writeArrayValue(Value const& value); diff --git a/grid_gen/mesh_conversion_tools/mpas_cell_culler.cpp b/grid_gen/mesh_conversion_tools/mpas_cell_culler.cpp index 291973ac6..b35568ee6 100644 --- a/grid_gen/mesh_conversion_tools/mpas_cell_culler.cpp +++ b/grid_gen/mesh_conversion_tools/mpas_cell_culler.cpp @@ -139,11 +139,11 @@ int main ( int argc, char *argv[] ) { for ( int i = 3; i < argc; i+=2 ) { if (strcmp(argv[i], "-m") == 0 ) { - mask_ops.push_back(merge); + mask_ops.push_back(static_cast(merge)); } else if ( strcmp(argv[i], "-i") == 0 ){ - mask_ops.push_back(invert); + mask_ops.push_back(static_cast(invert)); } else if ( strcmp(argv[i], "-p") == 0 ){ - mask_ops.push_back(preserve); + mask_ops.push_back(static_cast(preserve)); } else { cout << " ERROR: Invalid option passed on the command line " << argv[i] << ". Exiting..." << endl; exit(1); diff --git a/grid_gen/mesh_conversion_tools/mpas_mask_creator.cpp b/grid_gen/mesh_conversion_tools/mpas_mask_creator.cpp index 288416752..39cff2a1d 100644 --- a/grid_gen/mesh_conversion_tools/mpas_mask_creator.cpp +++ b/grid_gen/mesh_conversion_tools/mpas_mask_creator.cpp @@ -978,7 +978,7 @@ int resetFeatureInfo(){/*{{{*/ return 0; }/*}}}*/ int getFeatureInfo(const string featureFilename){/*{{{*/ - ifstream json_file(featureFilename); + ifstream json_file(featureFilename.c_str()); Json::Value root; string groupName, tempGroupName; vector pointIndices; @@ -1304,7 +1304,7 @@ int getFeatureInfo(const string featureFilename){/*{{{*/ return 0; }/*}}}*/ int getSeedInfo(const string seedFilename){/*{{{*/ - ifstream json_file(seedFilename); + ifstream json_file(seedFilename.c_str()); Json::Value root; json_file >> root; diff --git a/grid_gen/periodic_hex/Makefile b/grid_gen/periodic_hex/Makefile index cd290ca6d..0e4ee4842 100644 --- a/grid_gen/periodic_hex/Makefile +++ b/grid_gen/periodic_hex/Makefile @@ -37,16 +37,26 @@ LDFLAGS = -O3 -m64 CPP = cpp -P -traditional CPPFLAGS = CPPINCLUDES = -INCLUDES = -I$(NETCDF)/include - -# Specify NetCDF libraries, checking if netcdff is required (it will be present in v4 of netCDF) -LIBS = -L$(NETCDF)/lib -NCLIB = -lnetcdf -NCLIBF = -lnetcdff -ifneq ($(wildcard $(NETCDF)/lib/libnetcdff.*), ) # CHECK FOR NETCDF4 - LIBS += $(NCLIBF) -endif # CHECK FOR NETCDF4 -LIBS += $(NCLIB) + +ifneq ($(NETCDF), ) +ifneq ($(shell which ${NETCDF}/bin/nc-config 2> /dev/null), ) + LIBS = $(shell ${NETCDF}/bin/nc-config --libs) -lnetcdff + INCS = $(shell ${NETCDF}/bin/nc-config --cflags) +else + LIBS= -L${NETCDF}/lib + LIBS += -lnetcdf -lnetcdff + INCS = -I${NETCDF}/include +endif +else ifneq ($(shell which nc-config 2> /dev/null), ) + LIBS = $(shell nc-config --libs) -lnetcdff + INCS = $(shell nc-config --cflags) +else + LIBS= -L${NETCDF}/lib + LIBS += -lnetcdf_c++ -lnetcdf -lstdc++ -lnetcdff + INCS = -I${NETCDF}/include +endif + +INCLUDES = $(INCS) diff --git a/grid_gen/periodic_hex_minimal/Makefile b/grid_gen/periodic_hex_minimal/Makefile index 3799439a2..475b1a026 100644 --- a/grid_gen/periodic_hex_minimal/Makefile +++ b/grid_gen/periodic_hex_minimal/Makefile @@ -38,16 +38,26 @@ LDFLAGS = -O3 -m64 CPP = cpp -P -traditional CPPFLAGS = CPPINCLUDES = -INCLUDES = -I$(NETCDF)/include - -# Specify NetCDF libraries, checking if netcdff is required (it will be present in v4 of netCDF) -LIBS = -L$(NETCDF)/lib -NCLIB = -lnetcdf -NCLIBF = -lnetcdff -ifneq ($(wildcard $(NETCDF)/lib/libnetcdff.*), ) # CHECK FOR NETCDF4 - LIBS += $(NCLIBF) -endif # CHECK FOR NETCDF4 -LIBS += $(NCLIB) + +ifneq ($(NETCDF), ) +ifneq ($(shell which ${NETCDF}/bin/nc-config 2> /dev/null), ) + LIBS = $(shell ${NETCDF}/bin/nc-config --libs) -lnetcdff + INCS = $(shell ${NETCDF}/bin/nc-config --cflags) +else + LIBS= -L${NETCDF}/lib + LIBS += -lnetcdf -lnetcdff + INCS = -I${NETCDF}/include +endif +else ifneq ($(shell which nc-config 2> /dev/null), ) + LIBS = $(shell nc-config --libs) -lnetcdff + INCS = $(shell nc-config --cflags) +else + LIBS= -L${NETCDF}/lib + LIBS += -lnetcdf_c++ -lnetcdf -lstdc++ -lnetcdff + INCS = -I${NETCDF}/include +endif + +INCLUDES = $(INCS) diff --git a/grid_gen/periodic_quad/Makefile b/grid_gen/periodic_quad/Makefile index bcae7b09f..c60f922fb 100644 --- a/grid_gen/periodic_quad/Makefile +++ b/grid_gen/periodic_quad/Makefile @@ -37,16 +37,26 @@ LDFLAGS = -O3 -m64 CPP = cpp -P -traditional CPPFLAGS = CPPINCLUDES = -INCLUDES = -I$(NETCDF)/include - -# Specify NetCDF libraries, checking if netcdff is required (it will be present in v4 of netCDF) -LIBS = -L$(NETCDF)/lib -NCLIB = -lnetcdf -NCLIBF = -lnetcdff -ifneq ($(wildcard $(NETCDF)/lib/libnetcdff.*), ) # CHECK FOR NETCDF4 - LIBS += $(NCLIBF) -endif # CHECK FOR NETCDF4 -LIBS += $(NCLIB) + +ifneq ($(NETCDF), ) +ifneq ($(shell which ${NETCDF}/bin/nc-config 2> /dev/null), ) + LIBS = $(shell ${NETCDF}/bin/nc-config --libs) -lnetcdff + INCS = $(shell ${NETCDF}/bin/nc-config --cflags) +else + LIBS= -L${NETCDF}/lib + LIBS += -lnetcdf -lnetcdff + INCS = -I${NETCDF}/include +endif +else ifneq ($(shell which nc-config 2> /dev/null), ) + LIBS = $(shell nc-config --libs) -lnetcdff + INCS = $(shell nc-config --cflags) +else + LIBS= -L${NETCDF}/lib + LIBS += -lnetcdf_c++ -lnetcdf -lstdc++ -lnetcdff + INCS = -I${NETCDF}/include +endif + +INCLUDES = $(INCS) RM = rm -f diff --git a/grid_gen/points-mpas/Makefile b/grid_gen/points-mpas/Makefile index a38e7d139..c4e25ef33 100644 --- a/grid_gen/points-mpas/Makefile +++ b/grid_gen/points-mpas/Makefile @@ -4,14 +4,29 @@ CC=g++ CFLAGS= -O3 EXECUTABLE= PointsToMpas.x -NCINCDIR=${NETCDF}/include -NCLIBDIR=${NETCDF}/lib +ifneq ($(NETCDF), ) +ifneq ($(shell which ${NETCDF}/bin/nc-config 2> /dev/null), ) + LIBS = $(shell ${NETCDF}/bin/nc-config --libs) -lnetcdf_c++ + INCS = $(shell ${NETCDF}/bin/nc-config --cflags) +else + LIBS = -L${NETCDF}/lib + LIBS += -lnetcdf_c++ -lnetcdf -lstdc++ + INCS = -I${NETCDF}/include +endif +else ifneq ($(shell which nc-config 2> /dev/null), ) + LIBS = $(shell nc-config --libs) -lnetcdf_c++ + INCS = $(shell nc-config --cflags) +else + LIBS = -L${NETCDF}/lib + LIBS += -lnetcdf_c++ -lnetcdf -lstdc++ + INCS = -I${NETCDF}/include +endif all: - ${CC} points-mpas.cpp ${CFLAGS} -o ${EXECUTABLE} -I${NCINCDIR} -L${NCLIBDIR} -lnetcdf_c++ -lnetcdf + ${CC} points-mpas.cpp ${CFLAGS} -o ${EXECUTABLE} ${INCS} ${LIBS} -lnetcdf_c++ -lnetcdf debug: - ${CC} points-mpas.cpp -g -o ${EXECUTABLE} -I${NCINCDIR} -L${NCLIBDIR} -lnetcdf_c++ -lnetcdf + ${CC} points-mpas.cpp -g -o ${EXECUTABLE} ${INCS} ${LIBS} -lnetcdf_c++ -lnetcdf clean: rm -f grid.nc