Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/GNUmakefile =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/GNUmakefile (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/GNUmakefile (revision 1795) @@ -0,0 +1,14 @@ +# This is the makefile for "Max_diff_nc" with its libraries. For GNU +# make. + +SUBDIRS = Max_diff_nc NetCDF95 NR_util Jumble + +.PHONY: all clean $(SUBDIRS) + +all clean: $(SUBDIRS) + +$(SUBDIRS): + $(MAKE) -C $@ ${MAKECMDGOALS} + +Max_diff_nc: NetCDF95 Jumble +Jumble: NR_util Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/Jumble/GNUmakefile =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/Jumble/GNUmakefile (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/Jumble/GNUmakefile (revision 1795) @@ -0,0 +1,38 @@ +# This is a makefile for GNU make. +# This makefile builds "libjumble.a". + +# 1. Source files + +VPATH = Numerical + +sources := $(sort avg_mag.f90 count_lines.f90 opt_merge.f90 point.f90 compare.f90 csvread.f90 new_unit.f90 read_column.f90 jumble.f90 averge.f90 get_divisors.f90 dtridgl.f90 quadrat.f90 spherical.f90 prep_file.f90 prt_cmp.f90) + +# 2. Objects and library + +objects := $(sources:.f90=.o) +lib = libjumble.a + +# 3. Compiler-dependent part + +override FFLAGS += -I../NR_util + +# 4. Rules + +# Extend known suffixes: +%.o: %.f90 + $(COMPILE.f) $(OUTPUT_OPTION) $< + +.PHONY: all clean depend + +all: ${lib} + +${lib}: ${lib}(${objects}) + +depend depend.mk: + makedepf90 -Wmissing -Wconfused -I${VPATH} -nosrc -u nr_util ${sources} >depend.mk + +clean: + rm -f ${lib} ${objects} + +# Dependencies between object files and include files: +include depend.mk Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/Jumble/Numerical/averge.f90 =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/Jumble/Numerical/averge.f90 (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/Jumble/Numerical/averge.f90 (revision 1795) @@ -0,0 +1,18 @@ +module AVERGE_m + + implicit none + +contains + + pure real FUNCTION AVERGE(X,Y) + + real, intent(in):: x, y + + !------------- + + AVERGE = .5 * SQRT(X) * sqrt(Y) + 0.25 * (X + Y) + ! (the square roots should be separated to avoid underflow) + + END FUNCTION AVERGE + +end module AVERGE_m Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/Jumble/Numerical/dtridgl.f90 =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/Jumble/Numerical/dtridgl.f90 (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/Jumble/Numerical/dtridgl.f90 (revision 1795) @@ -0,0 +1,42 @@ +module DTRIDGL_m + + implicit none + +contains + + SUBROUTINE DTRIDGL(L,AF,BF,CF,DF,XK) + + ! Double precision version of tridgl. This subroutine solves a + ! system of tridiagional matrix equations. The form of the + ! equations are: + ! A(I)*X(I-1) + B(I)*X(I) + C(I)*X(I+1) = D(I) + ! where i=1,l less than 103. Reviewed -CP + + integer, intent(in):: l + double precision, intent(in):: AF(L),BF(L),CF(L),DF(L) + double precision, intent(out):: XK(L) + + ! Variables local to the procedure: + + integer, PARAMETER:: NMAX=201 + double precision AS(NMAX),DS(NMAX), xkb, x + integer i + + !---------------------------- + + AS(L) = AF(L)/BF(L) + DS(L) = DF(L)/BF(L) + DO I=2,L + X=1./(BF(L+1-I) - CF(L+1-I)*AS(L+2-I)) + AS(L+1-I)=AF(L+1-I)*X + DS(L+1-I)=(DF(L+1-I)-CF(L+1-I)*DS(L+2-I))*X + end DO + XK(1)=DS(1) + DO I=2,L + XKB=XK(I-1) + XK(I)=DS(I)-AS(I)*XKB + end DO + + END SUBROUTINE DTRIDGL + +end module DTRIDGL_m Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/Jumble/Numerical/get_divisors.f90 =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/Jumble/Numerical/get_divisors.f90 (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/Jumble/Numerical/get_divisors.f90 (revision 1795) @@ -0,0 +1,61 @@ +module get_divisors_m + + implicit none + +contains + + subroutine get_divisors(n, div) + + ! Find all the divisors of a given integer. + + integer, intent(in):: n + integer, pointer:: div(:) + + ! Variables local to the procedure: + + integer i, i_max + integer n_div ! number of divisors of "n" + integer work(n) ! temporary array to hold divisors + + !------------------------ + + if (n <= 0) then + print *, "get_divisors: n = ", n + stop 1 + end if + + if (n == 1) then + work(1) = 1 + n_div = 1 + else + ! n >= 2, there are at least two divisors: 1 and "n" + work(1) = 1 + work(2) = n + n_div = 2 + + i_max = int(sqrt(real(n))) + + do i = 2, i_max - 1 + if (mod(n, i) == 0) then + work(n_div+1) = i + work(n_div+2) = n / i + n_div = n_div + 2 + end if + end do + + if (i_max >= 2 .and. mod(n, i_max) == 0) then + work(n_div+1) = i_max + n_div = n_div + 1 + if (i_max**2 /= n) then + work(n_div+1) = n / i_max + n_div = n_div + 1 + end if + end if + end if + + allocate(div(n_div)) + div = work(:n_div) + + end subroutine get_divisors + +end module get_divisors_m Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/Jumble/Numerical/quadrat.f90 =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/Jumble/Numerical/quadrat.f90 (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/Jumble/Numerical/quadrat.f90 (revision 1795) @@ -0,0 +1,37 @@ +module quadrat_m + + implicit none + +contains + + subroutine quadrat(a, b, c, delta, root) + + ! This subroutine computes the real roots of a quadratic equation + ! with real coefficients. If there is a double root, it appears in + ! double in the output array. If there are two distinct roots, + ! they are output in increasing order. If there is no real root, + ! the output array is undefined. + + real, intent(in):: a, b, c + real, intent(out):: delta + real, intent(out):: root(:) ! (2) + + ! Variables local to the procedure + real q + + !--------------- + + delta = b**2 - 4. * a * c + if (delta == 0.) then + root = - b / 2 / a + else if (delta > 0.) then + q = - (b + sign(sqrt(delta), b)) / 2 + root(1) = q / a + root(2) = c / q + ! Sort the roots: + if (root(1) > root(2)) root = root(2:1:-1) + end if + + end subroutine quadrat + +end module quadrat_m Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/Jumble/Numerical/spherical.f90 =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/Jumble/Numerical/spherical.f90 (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/Jumble/Numerical/spherical.f90 (revision 1795) @@ -0,0 +1,120 @@ +module spherical + + implicit none + +contains + + + subroutine rectsph(rect,col,azm,r,r2) + + ! Converts rectangular coordinates to spherical coordinates. + ! (Ox) is taken as the polar axis. + ! Azimuth is the angle in the (Oyz) plane, measured from vector y. + ! If neither "r" nor "r2" is present then they are assumed to be + ! equal to 1. + + real, intent(in):: rect(3) + ! (Rectangular coordinates x,y,z. Should be different from (0,0 + ! ,0)) + + real, intent(out):: col ! Colatitude (in radians) + real, intent(out):: azm ! Azimuth (in radians). -pi < azm <= pi + real, intent(out), optional:: r ! Radius + real, intent(out), optional:: r2 ! Square of radius + + ! Local variables: + logical p, p2 + real d ! Radius + real d2 ! Square of radius + + !---------------------------------------------------------------- + + p = present(r) + p2 = present(r2) + + if (.not. p .and. .not. p2) then + ! Distance to the origin is assumed to be 1 + col = acos(rect(1)) + else + ! Compute the distance to the origin: + d2 = dot_product(rect,rect) + d = sqrt(d2) + col = acos(rect(1)/d) + if (p) r = d + if (p2) r2 = d2 + end if + + if (rect(2) == 0. .and. rect(3) == 0.) then + azm = 0. + ! (arbitrary value, azimuth is not well defined since vector + ! position is parallel + ! to polar axis) + else + azm = atan2(rect(3),rect(2)) + end if + + end subroutine rectsph + + !******************************************************************* + + real function sphbase(col,azm) + + ! This function returns the matrix of the spherical base: (radial + ! vector, colatitude vector, azimuthal vector) in the cartesian + ! vector base: (x, y, z). Vector x is assumed to be the polar + ! vector. Colatitude "col" and azimuth "azm" are the angular + ! spherical coordinates of the radial vector (azimuth is measured + ! from vector "y"). Note that if col = 0 (radial vector parallel + ! to x) then the choice of either the colatitude vector or the + ! azimuthal vector is arbitrary. The choice made depends on the + ! input value of azm. If azm is 0 then sphbase will return vector + ! y as the colatitude vector ("sphbase" is then the identity + ! matrix). Uses "sphrect", which converts spherical coordinates to + ! rectangular coordinates. + + dimension sphbase(3,3) ! Transformation matrix (no dimension) + + real, intent(in):: col ! Colatitude (in radians) + real, intent(in):: azm ! Azimuth (in radians) + + ! Local variable: + real pi + + !----------------------------------------------------------------- + + pi = acos(-1.) + + ! Radial vector: + sphbase(:,1) = sphrect(col, azm) + + ! Colatitude vector: + sphbase(:,2) = sphrect(col + pi/2, azm) + + ! Azimuthal vector: + sphbase(:,3) = (/0., - sin(azm), cos(azm)/) + + end function sphbase + + !***************************************************************** + + real function sphrect(col,azm,r) + + ! Converts spherical coordinates to rectangular coordinates. (Ox) + ! is taken as the polar axis. Azimuth "azm" is the angle in the + ! (Oyz) plane, measured from vector y. + ! If r is not present then we assume r = 1. + + dimension sphrect(3) ! Rectangular coordinates = (x,y,z) + + real, intent(in):: col ! Colatitude (in radians) + real, intent(in):: azm ! Azimuth (in radians) + real, intent(in), optional:: r ! Radius + + !--------------------------------------------------------------- + + sphrect = (/cos(col), sin(col) * cos(azm), sin(col) * sin(azm)/) + if (present(r)) sphrect = sphrect * r + + end function sphrect + +end module spherical Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/Jumble/avg_mag.f90 =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/Jumble/avg_mag.f90 (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/Jumble/avg_mag.f90 (revision 1795) @@ -0,0 +1,127 @@ +module avg_mag_m + + ! The generic procedure computes the average magnitude, that is + ! log10 of absolute value of an array. + ! The difference between the specific procedures is the kind and + ! rank of the array. + ! We do not care here about precision so all specific procedures + ! compute and return a default real kind value. + + implicit none + + interface avg_mag + module procedure avg_mag1, avg_mag1_dble, avg_mag2, avg_mag2_dble, & + avg_mag3, avg_mag3_dble, avg_mag4, avg_mag4_dble + end interface + + private + public avg_mag + +contains + + pure real function avg_mag1(a) + + real, intent(in):: a(:) + + ! Variables local to the procedure: + logical not_zero(size(a)) ! not zero in "a" + real magnit(size(a)) ! magnitudes of elements of "a" + + !------------------------------------- + + not_zero = a /= 0. + + if (any(not_zero)) then + where (not_zero) magnit = log10(abs(a)) + avg_mag1 = sum(magnit, mask=not_zero) / count(not_zero) + else + avg_mag1 = - huge(0.) ! minus infinity + end if + + end function avg_mag1 + + !******************************************************************* + + pure real function avg_mag1_dble(a) + + double precision, intent(in):: a(:) + + !------------------------------------- + + avg_mag1_dble = avg_mag1(real(a)) + + end function avg_mag1_dble + + !******************************************************************* + + pure real function avg_mag2(a) + + real, intent(in):: a(:, :) + + !------------------------------------- + + avg_mag2 = avg_mag1(pack(a, .true.)) + + end function avg_mag2 + + !******************************************************************* + + pure real function avg_mag2_dble(a) + + double precision, intent(in):: a(:, :) + + !------------------------------------- + + avg_mag2_dble = avg_mag1(pack(real(a), .true.)) + + end function avg_mag2_dble + + !******************************************************************* + + pure real function avg_mag3(a) + + real, intent(in):: a(:, :, :) + + !------------------------------------- + + avg_mag3 = avg_mag1(pack(a, .true.)) + + end function avg_mag3 + + !******************************************************************* + + pure real function avg_mag3_dble(a) + + double precision, intent(in):: a(:, :, :) + + !------------------------------------- + + avg_mag3_dble = avg_mag1(pack(real(a), .true.)) + + end function avg_mag3_dble + + !******************************************************************* + + pure real function avg_mag4(a) + + real, intent(in):: a(:, :, :, :) + + !------------------------------------- + + avg_mag4 = avg_mag1(pack(a, .true.)) + + end function avg_mag4 + + !******************************************************************* + + pure real function avg_mag4_dble(a) + + double precision, intent(in):: a(:, :, :, :) + + !------------------------------------- + + avg_mag4_dble = avg_mag1(pack(real(a), .true.)) + + end function avg_mag4_dble + +end module avg_mag_m Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/Jumble/compare.f90 =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/Jumble/compare.f90 (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/Jumble/compare.f90 (revision 1795) @@ -0,0 +1,775 @@ +module compare_m + + use prt_cmp_m, only: prt_cmp + use avg_mag_m, only: avg_mag + + implicit none + + interface compare + ! Makes a numerical comparison between two arrays of rank 1 to 4, + ! of types real or double precision. + module procedure compare1, compare1_dble, compare2, compare2_dble, & + compare3, compare3_dble, compare4, compare4_dble + end interface + + character(len=*), parameter:: dashes & + = "---------------------------------------------" + + private + public compare + +contains + + subroutine compare1(data_old, data_new, tag, comp_mag, report_id, quiet) + + ! Rank 1, real + + use nr_util, only: assert_eq + + real, intent(in):: data_old(:), data_new(:) + character(len=*), intent(in):: tag + logical, intent(in):: comp_mag + logical, intent(in):: report_id ! report identical variables + logical, intent(in):: quiet + + ! Variables local to the subprogram: + + logical zero(size(data_old)) + + real rel_diff(size(data_old)) + ! (absolute value of relative difference) + + real abs_diff(size(data_old)) + ! (absolute value of absolute difference) + + integer data_size + integer location(1) + character(len=len(dashes)+len(tag)+20) tag_fmt + + !------------------------------------------------------ + + data_size = assert_eq(size(data_old), size(data_new), & + "compare1: sizes differ -- " // tag) + if (quiet) then + tag_fmt = '(a, ":")' + else + tag_fmt = '(/, "' // dashes // '", /, a, ":")' + end if + + if (data_size == 0) then + print tag_fmt, tag + print *, "Zero-sized array" + else + if (all(data_old == data_new)) then + if (report_id) then + write(unit=*, fmt=tag_fmt, advance="no") tag + print *, "arrays are identical." + end if + else + if (quiet) then + write(unit=*, fmt=tag_fmt, advance="no") tag + print *, "arrays are different." + else + zero = data_old == 0. .or. data_new == 0. + where (.not. zero) rel_diff = abs(data_new / data_old - 1.) + abs_diff = abs(data_new - data_old) + + print tag_fmt, tag + + if (comp_mag) then + print * + print '("Average value of log10(abs(data_old)): ", 1pg8.1)', & + avg_mag(data_old) + end if + + if (any(.not. zero)) then + print * + print *, 'Relative difference for non-zero values:' + location = maxloc(rel_diff, mask=.not. zero) + call prt_cmp(data_old(location(1)), data_new(location(1)), & + rel_diff(location(1)), location) + end if + + if (any(zero)) then + print * + print *, 'Absolute difference when there is a zero:' + location = maxloc(abs_diff, mask=zero) + call prt_cmp(data_old(location(1)), data_new(location(1)), & + abs_diff(location(1)), location) + end if + + print * + print *, 'Absolute difference:' + location = maxloc(abs_diff) + call prt_cmp(data_old(location(1)), data_new(location(1)), & + abs_diff(location(1)), location) + end if + end if + end if + + end subroutine compare1 + + !*********************************************************** + + subroutine compare1_dble(data_old, data_new, tag, comp_mag, report_id, quiet) + + ! Rank 1, double precision + + use nr_util, only: assert_eq + + double precision, intent(in):: data_old(:), data_new(:) + character(len=*), intent(in):: tag + logical, intent(in):: comp_mag + logical, intent(in):: report_id ! report identical variables + logical, intent(in):: quiet + + ! Variables local to the subprogram: + + logical zero(size(data_old)) + + double precision rel_diff(size(data_old)) + ! (absolute value of relative difference) + + double precision abs_diff(size(data_old)) + ! (absolute value of absolute difference) + + integer data_size + integer location(1) + character(len=len(dashes)+len(tag)+20) tag_fmt + + !------------------------------------------------------ + + data_size = assert_eq(size(data_old), size(data_new), & + "compare1_dble: sizes differ -- " // tag) + if (quiet) then + tag_fmt = '(a, ":")' + else + tag_fmt = '(/, "' // dashes // '", /, a, ":")' + end if + + if (data_size == 0) then + print tag_fmt, tag + print *, "Zero-sized array" + else + if (all(data_old == data_new)) then + if (report_id) then + write(unit=*, fmt=tag_fmt, advance="no") tag + print *, "arrays are identical." + end if + else + if (quiet) then + write(unit=*, fmt=tag_fmt, advance="no") tag + print *, "arrays are different." + else + zero = data_old == 0d0 .or. data_new == 0d0 + where (.not. zero) rel_diff = abs(data_new / data_old - 1d0) + abs_diff = abs(data_new - data_old) + + print tag_fmt, tag + + if (comp_mag) then + print * + print '("Average value of log10(abs(data_old)): ", 1pg8.1)', & + avg_mag(data_old) + end if + + if (any(.not. zero)) then + print * + print *, 'Relative difference for non-zero values:' + location = maxloc(rel_diff, mask=.not. zero) + call prt_cmp(data_old(location(1)), data_new(location(1)), & + rel_diff(location(1)), location) + end if + + if (any(zero)) then + print * + print *, 'Absolute difference when there is a zero:' + location = maxloc(abs_diff, mask=zero) + call prt_cmp(data_old(location(1)), data_new(location(1)), & + abs_diff(location(1)), location) + end if + + print * + print *, 'Absolute difference:' + location = maxloc(abs_diff) + call prt_cmp(data_old(location(1)), data_new(location(1)), & + abs_diff(location(1)), location) + end if + end if + end if + + end subroutine compare1_dble + + !*********************************************************** + + subroutine compare2(data_old, data_new, tag, comp_mag, report_id, quiet) + + ! Rank 2, real + + use nr_util, only: assert + use point_m, only: point + + real, intent(in):: data_old(:,:), data_new(:,:) + character(len=*), intent(in):: tag + logical, intent(in):: comp_mag + logical, intent(in):: report_id ! report identical variables + logical, intent(in):: quiet + + ! Variables local to the subprogram: + + logical zero(size(data_old,1), size(data_old,2)) + + real rel_diff(size(data_old,1), size(data_old,2)) + ! (absolute value of relative difference) + + real abs_diff(size(data_old,1), size(data_old,2)) + ! (absolute value of absolute difference) + + integer location(2) + character(len=len(dashes)+len(tag)+20) tag_fmt + + !------------------------------------------------------ + + call assert(shape(data_old) == shape(data_new), & + "compare2: shapes differ -- " // tag) + if (quiet) then + tag_fmt = '(a, ":")' + else + tag_fmt = '(/, "' // dashes // '", /, a, ":")' + end if + + if (size(data_old) == 0) then + print tag_fmt, tag + print *, "Zero-sized array" + else + if (all(data_old == data_new)) then + if (report_id) then + write(unit=*, fmt=tag_fmt, advance="no") tag + print *, "arrays are identical." + end if + else + if (quiet) then + write(unit=*, fmt=tag_fmt, advance="no") tag + print *, "arrays are different." + else + zero = data_old == 0. .or. data_new == 0. + where (.not. zero) rel_diff = abs(data_new / data_old - 1.) + abs_diff = abs(data_new - data_old) + + print tag_fmt, tag + + if (comp_mag) then + print * + print '("Average value of log10(abs(data_old)): ", 1pg8.1)', & + avg_mag(data_old) + end if + + if (any(.not. zero)) then + print * + print *, 'Relative difference for non-zero values:' + location = maxloc(rel_diff, mask=.not. zero) + call prt_cmp(point(data_old, location), & + point(data_new, location), point(rel_diff, location), & + location) + end if + + if (any(zero)) then + print * + print *, 'Absolute difference when there is a zero:' + location = maxloc(abs_diff, mask=zero) + call prt_cmp(point(data_old, location), & + point(data_new, location), point(abs_diff, location), & + location) + end if + + print * + print *, 'Absolute difference:' + location = maxloc(abs_diff) + call prt_cmp(point(data_old, location), & + point(data_new, location), point(abs_diff, location), & + location) + end if + end if + end if + + end subroutine compare2 + + !*********************************************************** + + subroutine compare2_dble(data_old, data_new, tag, comp_mag, report_id, quiet) + + ! Rank 2, double precision + + use nr_util, only: assert + use point_m, only: point + + double precision, intent(in):: data_old(:,:), data_new(:,:) + character(len=*), intent(in):: tag + logical, intent(in):: comp_mag + logical, intent(in):: report_id ! report identical variables + logical, intent(in):: quiet + + ! Variables local to the subprogram: + + logical zero(size(data_old,1), size(data_old,2)) + + double precision rel_diff(size(data_old,1), size(data_old,2)) + ! (absolute value of relative difference) + + double precision abs_diff(size(data_old,1), size(data_old,2)) + ! (absolute value of absolute difference) + + integer location(2) + character(len=len(dashes)+len(tag)+20) tag_fmt + + !------------------------------------------------------ + + call assert(shape(data_old) == shape(data_new), & + "compare2_dble: shapes differ -- " // tag) + if (quiet) then + tag_fmt = '(a, ":")' + else + tag_fmt = '(/, "' // dashes // '", /, a, ":")' + end if + + if (size(data_old) == 0) then + print tag_fmt, tag + print *, "Zero-sized array" + else + if (all(data_old == data_new)) then + if (report_id) then + write(unit=*, fmt=tag_fmt, advance="no") tag + print *, "arrays are identical." + end if + else + if (quiet) then + write(unit=*, fmt=tag_fmt, advance="no") tag + print *, "arrays are different." + else + zero = data_old == 0d0 .or. data_new == 0d0 + where (.not. zero) rel_diff = abs(data_new / data_old - 1d0) + abs_diff = abs(data_new - data_old) + + print tag_fmt, tag + + if (comp_mag) then + print * + print '("Average value of log10(abs(data_old)): ", 1pg8.1)', & + avg_mag(data_old) + end if + + if (any(.not. zero)) then + print * + print *, 'Relative difference for non-zero values:' + location = maxloc(rel_diff, mask=.not. zero) + call prt_cmp(point(data_old, location), & + point(data_new, location), point(rel_diff, location), & + location) + end if + + if (any(zero)) then + print * + print *, 'Absolute difference when there is a zero:' + location = maxloc(abs_diff, mask=zero) + call prt_cmp(point(data_old, location), & + point(data_new, location), point(abs_diff, location), & + location) + end if + + print * + print *, 'Absolute difference:' + location = maxloc(abs_diff) + call prt_cmp(point(data_old, location), & + point(data_new, location), point(abs_diff, location), & + location) + end if + end if + end if + + end subroutine compare2_dble + + !*********************************************************** + + subroutine compare3(data_old, data_new, tag, comp_mag, report_id, quiet) + + ! Rank 3, real + + use nr_util, only: assert + use point_m, only: point + + real, intent(in):: data_old(:, :, :), data_new(:, : ,:) + character(len=*), intent(in):: tag + logical, intent(in):: comp_mag + logical, intent(in):: report_id ! report identical variables + logical, intent(in):: quiet + + ! Variables local to the subprogram: + + logical zero(size(data_old,1), size(data_old,2), size(data_old,3)) + + real rel_diff(size(data_old,1), size(data_old,2), size(data_old,3)) + ! (absolute value of relative difference) + + real abs_diff(size(data_old,1), size(data_old,2), size(data_old,3)) + ! (absolute value of absolute difference) + + integer location(3) + character(len=len(dashes)+len(tag)+20) tag_fmt + + !------------------------------------------------------ + + call assert(shape(data_old) == shape(data_new), & + "compare3: shapes differ -- " // tag) + if (quiet) then + tag_fmt = '(a, ":")' + else + tag_fmt = '(/, "' // dashes // '", /, a, ":")' + end if + + if (size(data_old) == 0) then + print tag_fmt, tag + print *, "Zero-sized array" + else + if (all(data_old == data_new)) then + if (report_id) then + write(unit=*, fmt=tag_fmt, advance="no") tag + print *, "arrays are identical." + end if + else + if (quiet) then + write(unit=*, fmt=tag_fmt, advance="no") tag + print *, "arrays are different." + else + zero = data_old == 0. .or. data_new == 0. + where (.not. zero) rel_diff = abs(data_new / data_old - 1.) + abs_diff = abs(data_new - data_old) + + print tag_fmt, tag + + if (comp_mag) then + print * + print '("Average value of log10(abs(data_old)): ", 1pg8.1)', & + avg_mag(data_old) + end if + + if (any(.not. zero)) then + print * + print *, 'Relative difference for non-zero values:' + location = maxloc(rel_diff, mask=.not. zero) + call prt_cmp(point(data_old, location), & + point(data_new, location), point(rel_diff, location), & + location) + end if + + if (any(zero)) then + print * + print *, 'Absolute difference when there is a zero:' + location = maxloc(abs_diff, mask=zero) + call prt_cmp(point(data_old, location), & + point(data_new, location), point(abs_diff, location), & + location) + end if + + print * + print *, 'Absolute difference:' + location = maxloc(abs_diff) + call prt_cmp(point(data_old, location), & + point(data_new, location), point(abs_diff, location), & + location) + end if + end if + end if + + end subroutine compare3 + + !*********************************************************** + + subroutine compare3_dble(data_old, data_new, tag, comp_mag, report_id, quiet) + + ! Rank 3, double precision + + use nr_util, only: assert + use point_m, only: point + + double precision, intent(in):: data_old(:, :, :), data_new(:, : ,:) + character(len=*), intent(in):: tag + logical, intent(in):: comp_mag + logical, intent(in):: report_id ! report identical variables + logical, intent(in):: quiet + + ! Variables local to the subprogram: + + logical zero(size(data_old,1), size(data_old,2), size(data_old,3)) + + double precision rel_diff(size(data_old,1), size(data_old,2), & + size(data_old,3)) + ! (absolute value of relative difference) + + double precision abs_diff(size(data_old,1), size(data_old,2), & + size(data_old,3)) + ! (absolute value of absolute difference) + + integer location(3) + character(len=len(dashes)+len(tag)+20) tag_fmt + + !------------------------------------------------------ + + call assert(shape(data_old) == shape(data_new), & + "compare3_dble: shapes differ -- " // tag) + if (quiet) then + tag_fmt = '(a, ":")' + else + tag_fmt = '(/, "' // dashes // '", /, a, ":")' + end if + + if (size(data_old) == 0) then + print tag_fmt, tag + print *, "Zero-sized array" + else + if (all(data_old == data_new)) then + if (report_id) then + write(unit=*, fmt=tag_fmt, advance="no") tag + print *, "arrays are identical." + end if + else + if (quiet) then + write(unit=*, fmt=tag_fmt, advance="no") tag + print *, "arrays are different." + else + zero = data_old == 0. .or. data_new == 0. + where (.not. zero) rel_diff = abs(data_new / data_old - 1.) + abs_diff = abs(data_new - data_old) + + print tag_fmt, tag + + if (comp_mag) then + print * + print '("Average value of log10(abs(data_old)): ", 1pg8.1)', & + avg_mag(data_old) + end if + + if (any(.not. zero)) then + print * + print *, 'Relative difference for non-zero values:' + location = maxloc(rel_diff, mask=.not. zero) + call prt_cmp(point(data_old, location), & + point(data_new, location), point(rel_diff, location), & + location) + end if + + if (any(zero)) then + print * + print *, 'Absolute difference when there is a zero:' + location = maxloc(abs_diff, mask=zero) + call prt_cmp(point(data_old, location), & + point(data_new, location), point(abs_diff, location), & + location) + end if + + print * + print *, 'Absolute difference:' + location = maxloc(abs_diff) + call prt_cmp(point(data_old, location), & + point(data_new, location), point(abs_diff, location), & + location) + end if + end if + end if + + end subroutine compare3_dble + + !*********************************************************** + + subroutine compare4(data_old, data_new, tag, comp_mag, report_id, quiet) + + ! Rank 4, real + + use nr_util, only: assert + use point_m, only: point + + real, intent(in):: data_old(:, :, :, :), data_new(:, : ,:, :) + character(len=*), intent(in):: tag + logical, intent(in):: comp_mag + logical, intent(in):: report_id ! report identical variables + logical, intent(in):: quiet + + ! Variables local to the subprogram: + + logical zero(size(data_old,1), size(data_old,2), size(data_old,3), & + size(data_old,4)) + + real rel_diff(size(data_old,1), size(data_old,2), size(data_old,3), & + size(data_old,4)) + ! (absolute value of relative difference) + + real abs_diff(size(data_old,1), size(data_old,2), size(data_old,3), & + size(data_old,4)) + ! (absolute value of absolute difference) + + integer location(4) + character(len=len(dashes)+len(tag)+20) tag_fmt + + !------------------------------------------------------ + + call assert(shape(data_old) == shape(data_new), & + "compare4: shapes differ -- " // tag) + if (quiet) then + tag_fmt = '(a, ":")' + else + tag_fmt = '(/, "' // dashes // '", /, a, ":")' + end if + + if (size(data_old) == 0) then + print tag_fmt, tag + print *, "Zero-sized array" + else + if (all(data_old == data_new)) then + if (report_id) then + write(unit=*, fmt=tag_fmt, advance="no") tag + print *, "arrays are identical." + end if + else + if (quiet) then + write(unit=*, fmt=tag_fmt, advance="no") tag + print *, "arrays are different." + else + zero = data_old == 0. .or. data_new == 0. + where (.not. zero) rel_diff = abs(data_new / data_old - 1.) + abs_diff = abs(data_new - data_old) + + print tag_fmt, tag + + if (comp_mag) then + print * + print '("Average value of log10(abs(data_old)): ", 1pg8.1)', & + avg_mag(data_old) + end if + + if (any(.not. zero)) then + print * + print *, 'Relative difference for non-zero values:' + location = maxloc(rel_diff, mask=.not. zero) + call prt_cmp(point(data_old, location), & + point(data_new, location), point(rel_diff, location), & + location) + end if + + if (any(zero)) then + print * + print *, 'Absolute difference when there is a zero:' + location = maxloc(abs_diff, mask=zero) + call prt_cmp(point(data_old, location), & + point(data_new, location), point(abs_diff, location), & + location) + end if + + print * + print *, 'Absolute difference:' + location = maxloc(abs_diff) + call prt_cmp(point(data_old, location), & + point(data_new, location), point(abs_diff, location), & + location) + end if + end if + end if + + end subroutine compare4 + + !*********************************************************** + + subroutine compare4_dble(data_old, data_new, tag, comp_mag, report_id, quiet) + + ! Rank 4, double precision + + use nr_util, only: assert + use point_m, only: point + + double precision, intent(in):: data_old(:, :, :, :), data_new(:, : ,:, :) + character(len=*), intent(in):: tag + logical, intent(in):: comp_mag + logical, intent(in):: report_id ! report identical variables + logical, intent(in):: quiet + + ! Variables local to the subprogram: + + logical zero(size(data_old,1), size(data_old,2), size(data_old,3), & + size(data_old,4)) + + double precision rel_diff(size(data_old,1), size(data_old,2), & + size(data_old,3), size(data_old,4)) + ! (absolute value of relative difference) + + double precision abs_diff(size(data_old,1), size(data_old,2), & + size(data_old,3), size(data_old,4)) + ! (absolute value of absolute difference) + + integer location(4) + character(len=len(dashes)+len(tag)+20) tag_fmt + + !------------------------------------------------------ + + call assert(shape(data_old) == shape(data_new), & + "compare4_dble: shapes differ -- " // tag) + if (quiet) then + tag_fmt = '(a, ":")' + else + tag_fmt = '(/, "' // dashes // '", /, a, ":")' + end if + + if (size(data_old) == 0) then + print tag_fmt, tag + print *, "Zero-sized array" + else + if (all(data_old == data_new)) then + if (report_id) then + write(unit=*, fmt=tag_fmt, advance="no") tag + print *, "arrays are identical." + end if + else + if (quiet) then + write(unit=*, fmt=tag_fmt, advance="no") tag + print *, "arrays are different." + else + zero = data_old == 0. .or. data_new == 0. + where (.not. zero) rel_diff = abs(data_new / data_old - 1.) + abs_diff = abs(data_new - data_old) + + print tag_fmt, tag + + if (comp_mag) then + print * + print '("Average value of log10(abs(data_old)): ", 1pg8.1)', & + avg_mag(data_old) + end if + + if (any(.not. zero)) then + print * + print *, 'Relative difference for non-zero values:' + location = maxloc(rel_diff, mask=.not. zero) + call prt_cmp(point(data_old, location), & + point(data_new, location), point(rel_diff, location), & + location) + end if + + if (any(zero)) then + print * + print *, 'Absolute difference when there is a zero:' + location = maxloc(abs_diff, mask=zero) + call prt_cmp(point(data_old, location), & + point(data_new, location), point(abs_diff, location), & + location) + end if + + print * + print *, 'Absolute difference:' + location = maxloc(abs_diff) + call prt_cmp(point(data_old, location), & + point(data_new, location), point(abs_diff, location), & + location) + end if + end if + end if + + end subroutine compare4_dble + +end module compare_m Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/Jumble/count_lines.f90 =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/Jumble/count_lines.f90 (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/Jumble/count_lines.f90 (revision 1795) @@ -0,0 +1,33 @@ +module count_lines_m + + implicit none + +contains + + subroutine count_lines(unit, n) + + ! This subroutine counts the number of lines in an external file, + ! from the current position, not necessarily the first record of the + ! file. + ! On return, the position is at the end of the file. + ! The file should be connected for sequential access. + ! The records of the file shoud be formatted. + + integer, intent(in):: unit ! external file unit + integer, intent(out):: n ! number of lines + + ! Variable local to the procedure: + integer iostat + + !------------------------------- + + n = 0 + do + read(unit, fmt=*, iostat=iostat) + if (iostat /= 0) exit + n = n + 1 + end do + + end subroutine count_lines + +end module count_lines_m Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/Jumble/csvread.f90 =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/Jumble/csvread.f90 (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/Jumble/csvread.f90 (revision 1795) @@ -0,0 +1,98 @@ +module csvread_m + + use new_unit_m, only: new_unit + use prep_file_m, only: prep_file + + implicit none + + private + public csvread + + interface csvread + ! Reads numeric values from a file. + ! Values must be separated by comma and/or blanks. + ! Values are read into a default real kind or double precision array. + ! The last column and/or last row parameters may be 0. + ! This is interpreted as "last in the file". + ! The only difference between the interfaces of the specific + ! procedures is the type of "a". + module procedure csvread_sp, csvread_dp + end interface + +contains + + subroutine csvread_sp(file, a, first_r, first_c, last_r, last_c) + + character(len=*), intent(in):: file + real, pointer:: a(:,:) + integer, intent(in), optional:: first_r ! (first row to read) + integer, intent(in), optional:: first_c ! (first column to read) + integer, intent(in), optional:: last_r ! (last row to read) + integer, intent(in), optional:: last_c ! (last column to read) + + ! Variables local to the subprogram: + integer i, j, unit + integer f_r_loc ! (first row to read, local variable) + integer f_c_loc ! (first column to read, local variable) + integer l_r_loc ! (last row to read, local variable) + integer l_c_loc ! (last column to read, local variable) + real trash + + !------------------------------------------------------ + + print *, 'Reading data from file "' // file // '"' + call new_unit(unit) + open(unit, file=file, status='old', action='read', position='rewind') + + call prep_file(unit, first_r, first_c, last_r, last_c, f_r_loc, & + f_c_loc, l_r_loc, l_c_loc) + + allocate(a(l_r_loc - f_r_loc + 1, l_c_loc - f_c_loc + 1)) + + do i = 1, l_r_loc - f_r_loc + 1 + read(unit, fmt=*) (trash, j = 1, f_c_loc - 1), a(i, :) + end do + + close(unit) + + end subroutine csvread_sp + + !*********************************************************** + + subroutine csvread_dp(file, a, first_r, first_c, last_r, last_c) + + character(len=*), intent(in):: file + double precision, pointer:: a(:,:) + integer, intent(in), optional:: first_r ! (first row to read) + integer, intent(in), optional:: first_c ! (first column to read) + integer, intent(in), optional:: last_r ! (last row to read) + integer, intent(in), optional:: last_c ! (last column to read) + + ! Variables local to the subprogram: + integer i, j, unit + integer f_r_loc ! (first row to read, local variable) + integer f_c_loc ! (first column to read, local variable) + integer l_r_loc ! (last row to read, local variable) + integer l_c_loc ! (last column to read, local variable) + double precision trash + + !------------------------------------------------------ + + print *, 'Reading data from file "' // file // '"' + call new_unit(unit) + open(unit, file=file, status='old', action='read', position='rewind') + + call prep_file(unit, first_r, first_c, last_r, last_c, f_r_loc, & + f_c_loc, l_r_loc, l_c_loc) + + allocate(a(l_r_loc - f_r_loc + 1, l_c_loc - f_c_loc + 1)) + + do i = 1, l_r_loc - f_r_loc + 1 + read(unit, fmt=*) (trash, j = 1, f_c_loc - 1), a(i, :) + end do + + close(unit) + + end subroutine csvread_dp + +end module csvread_m Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/Jumble/depend.mk =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/Jumble/depend.mk (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/Jumble/depend.mk (revision 1795) @@ -0,0 +1,5 @@ +compare.o : point.o avg_mag.o prt_cmp.o +csvread.o : prep_file.o new_unit.o +jumble.o : spherical.o quadrat.o dtridgl.o get_divisors.o averge.o read_column.o new_unit.o csvread.o compare.o point.o opt_merge.o count_lines.o avg_mag.o +prep_file.o : opt_merge.o +read_column.o : opt_merge.o new_unit.o Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/Jumble/jumble.f90 =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/Jumble/jumble.f90 (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/Jumble/jumble.f90 (revision 1795) @@ -0,0 +1,17 @@ +module jumble + + use avg_mag_m + use count_lines_m + use opt_merge_m + use point_m + use compare_m + use csvread_m + use new_unit_m + use read_column_m + use averge_m + use get_divisors_m + use dtridgl_m + use quadrat_m + use spherical + +end module jumble Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/Jumble/new_unit.f90 =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/Jumble/new_unit.f90 (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/Jumble/new_unit.f90 (revision 1795) @@ -0,0 +1,25 @@ +module new_unit_m + + implicit none + +contains + + subroutine new_unit(unit) + + integer, intent(out):: unit + + ! Variables local to the procedure: + logical opened, exist + + !------------------------------------------------------ + + unit = 0 + do + inquire(unit=unit, opened=opened, exist=exist) + if (exist .and. .not. opened) exit + unit = unit + 1 + end do + + end subroutine new_unit + +end module new_unit_m Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/Jumble/opt_merge.f90 =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/Jumble/opt_merge.f90 (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/Jumble/opt_merge.f90 (revision 1795) @@ -0,0 +1,26 @@ +module opt_merge_m + + implicit none + +contains + + integer function opt_merge(param, default) + + ! Analogous to the intrinsic procedure "merge" : merges an + ! optional parameter and a default value depending on the + ! presence of the optional parameter. + + integer, intent(in), optional:: param + integer, intent(in):: default + + !-------------- + + if (present(param)) then + opt_merge = param + else + opt_merge = default + end if + + end function opt_merge + +end module opt_merge_m Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/Jumble/point.f90 =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/Jumble/point.f90 (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/Jumble/point.f90 (revision 1795) @@ -0,0 +1,117 @@ +module point_m + + implicit none + + interface point + ! The difference between the procedures is the rank and type of + ! the first argument. + module procedure point_2, point_2_dble, point_3, point_3_dble, point_4, & + point_4_dble + end interface + + private + public point + +contains + + real function point_2(array, location) + + use nr_util, only: assert + + real, intent(in):: array(:, :) + integer, intent(in):: location(:) + + !--------------------------------- + + call assert(size(location) == 2, "point_2") + + point_2 = array(location(1), location(2)) + + end function point_2 + + !*************************************************** + + double precision function point_2_dble(array, location) + + use nr_util, only: assert + + double precision, intent(in):: array(:, :) + integer, intent(in):: location(:) + + !--------------------------------- + + call assert(size(location) == 2, "point_2_dble") + + point_2_dble = array(location(1), location(2)) + + end function point_2_dble + + !*************************************************** + + real function point_3(array, location) + + use nr_util, only: assert + + real, intent(in):: array(:, :, :) + integer, intent(in):: location(:) + + !--------------------------------- + + call assert(size(location) == 3, "point_3") + + point_3 = array(location(1), location(2), location(3)) + + end function point_3 + + !*************************************************** + + double precision function point_3_dble(array, location) + + use nr_util, only: assert + + double precision, intent(in):: array(:, :, :) + integer, intent(in):: location(:) + + !--------------------------------- + + call assert(size(location) == 3, "point_3_dble") + + point_3_dble = array(location(1), location(2), location(3)) + + end function point_3_dble + + !*************************************************** + + real function point_4(array, location) + + use nr_util, only: assert + + real, intent(in):: array(:, :, :, :) + integer, intent(in):: location(:) + + !--------------------------------- + + call assert(size(location) == 4, "point_4") + + point_4 = array(location(1), location(2), location(3), location(4)) + + end function point_4 + + !*************************************************** + + double precision function point_4_dble(array, location) + + use nr_util, only: assert + + double precision, intent(in):: array(:, :, :, :) + integer, intent(in):: location(:) + + !--------------------------------- + + call assert(size(location) == 4, "point_4_dble") + + point_4_dble = array(location(1), location(2), location(3), location(4)) + + end function point_4_dble + +end module point_m Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/Jumble/prep_file.f90 =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/Jumble/prep_file.f90 (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/Jumble/prep_file.f90 (revision 1795) @@ -0,0 +1,79 @@ +module prep_file_m + + implicit none + +contains + + subroutine prep_file(unit, first_r, first_c, last_r, last_c, f_r_not_opt, & + f_c_not_opt, l_r_not_opt, l_c_not_opt) + + ! This subroutine is used by "csvread". It fills non-optional + ! arguments: first and last row, first and last column which will + ! actually be read, taking information from the file itself if + ! necessary. It also positions the input file on the first row to + ! read. + + use opt_merge_m, only: opt_merge + + integer, intent(in):: unit ! logical unit for input file + integer, intent(in), optional:: first_r ! (first row to read) + integer, intent(in), optional:: first_c ! (first column to read) + integer, intent(in), optional:: last_r ! (last row to read) + integer, intent(in), optional:: last_c ! (last column to read) + integer, intent(out):: f_r_not_opt ! (first row to read, not optional) + integer, intent(out):: f_c_not_opt ! (first column to read, not optional) + integer, intent(out):: l_r_not_opt ! (last row to read, not optional) + integer, intent(out):: l_c_not_opt ! (last column to read, not optional) + + ! Variables local to the subprogram: + integer iostat, i + character c + logical prev_value ! previous character was part of a value + logical curr_value ! current character is part of a value + + !------------------------------------------------------ + + f_r_not_opt = opt_merge(first_r, 1) + f_c_not_opt = opt_merge(first_c, 1) + l_r_not_opt = opt_merge(last_r, 0) + l_c_not_opt = opt_merge(last_c, 0) + + if (l_r_not_opt == 0) then + ! Count the number of lines in the file: + l_r_not_opt = 0 + do + read(unit, fmt=*, iostat=iostat) + if (iostat /= 0) exit + l_r_not_opt = l_r_not_opt + 1 + end do + if (l_r_not_opt == 0) stop 'Empty file.' + + rewind(unit) + end if + + ! Go to first row to read: + do i = 1, f_r_not_opt - 1 + read(unit, fmt=*) + end do + + if (l_c_not_opt == 0) then + ! Count the number of values per line: + l_c_not_opt = 0 + curr_value = .false. + do + read(unit, fmt='(a)', advance='no', iostat=iostat) c + if (iostat /= 0) exit + prev_value = curr_value + curr_value = c /= " " .and. c /= "," + if (curr_value .and. .not. prev_value) l_c_not_opt = l_c_not_opt + 1 + end do + + backspace(unit) + end if + + print *, 'Reading column(s) ', f_c_not_opt, ':', l_c_not_opt, & + ', row(s) ', f_r_not_opt, ':', l_r_not_opt + + end subroutine prep_file + +end module prep_file_m Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/Jumble/prt_cmp.f90 =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/Jumble/prt_cmp.f90 (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/Jumble/prt_cmp.f90 (revision 1795) @@ -0,0 +1,83 @@ +module prt_cmp_m + + implicit none + + interface prt_cmp + ! Used by "compare". + module procedure prt_cmp_real, prt_cmp_double + end interface prt_cmp + + private + public prt_cmp + +contains + + subroutine prt_cmp_real(data_old, data_new, maxv, location) + + real, intent(in):: data_old, data_new, maxv + integer, intent(in):: location(:) + + ! Variable local to the procedure: + character(len=30) format_string + + !------------------------------------------------------ + + if (size(location) == 1) then + format_string = "(a, i0, a)" + else + ! size(location) >= 2 + write(unit=format_string, fmt="(a, i0, a)") "(a, ", & + size(location) - 1, "(i0, 1x), i0, a)" + end if + + print * + print '("Maximum difference: ", 1pg8.1)', maxv +!! print '("Maximum difference: ", 1pg9.2)', maxv + print *, 'Occurring at:' + + write(unit=*, fmt=format_string, advance="no") 'data_old(', location, & + ') = ' + print *, data_old + + write(unit=*, fmt=format_string, advance="no") 'data_new(', location, & + ') = ' + print *, data_new + + end subroutine prt_cmp_real + + !*********************************************************** + + subroutine prt_cmp_double(data_old, data_new, maxv, location) + + double precision, intent(in):: data_old, data_new, maxv + integer, intent(in):: location(:) + + ! Variable local to the procedure: + character(len=30) format_string + + !------------------------------------------------------ + + if (size(location) == 1) then + format_string = "(a, i0, a)" + else + ! size(location) >= 2 + write(unit=format_string, fmt="(a, i0, a)") "(a, ", & + size(location) - 1, "(i0, 1x), i0, a)" + end if + + print * + print '("Maximum difference: ", 1pg8.1)', maxv +!! print '("Maximum difference: ", 1pg9.2)', maxv + print *, 'Occurring at:' + + write(unit=*, fmt=format_string, advance="no") 'data_old(', location, & + ') = ' + print *, data_old + + write(unit=*, fmt=format_string, advance="no") 'data_new(', location, & + ') = ' + print *, data_new + + end subroutine prt_cmp_double + +end module prt_cmp_m Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/Jumble/read_column.f90 =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/Jumble/read_column.f90 (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/Jumble/read_column.f90 (revision 1795) @@ -0,0 +1,119 @@ +module read_column_m + + implicit none + + private + public read_column + + interface read_column + ! This generic procedure reads a column of values in an external + ! file. + ! The records of the file shoud be formatted. + ! If the argument "last" is 0 or is absent then the procedure + ! reads to the last line in the file. + ! The difference between the specific procedures is the type of + ! argument "a". + module procedure read_column_real, read_column_char + end interface + +contains + + subroutine read_column_real(file, a, first, last) + + use new_unit_m, only: new_unit + + character(len=*), intent(in):: file + real, pointer:: a(:) + integer, intent(in), optional:: first ! (first line to read) + integer, intent(in), optional:: last ! (last line to read) + + ! Variables local to the subprogram: + integer unit ! external file unit + integer first_not_opt ! first line to read, local variable + integer last_not_opt ! last line to read, local variable + + !------------------------------------------------------ + + call new_unit(unit) + open(unit, file=file, status='old', action='read', position='rewind') + call prep_file(unit, first, last, first_not_opt, last_not_opt) + allocate(a(last_not_opt - first_not_opt + 1)) + read(unit, fmt=*) a + close(unit) + + end subroutine read_column_real + + !*********************************************************** + + subroutine read_column_char(file, a, first, last) + + use new_unit_m, only: new_unit + + character(len=*), intent(in):: file + character(len=*), pointer:: a(:) + integer, intent(in), optional:: first ! (first line to read) + integer, intent(in), optional:: last ! (last line to read) + + ! Variables local to the subprogram: + integer unit + integer first_not_opt ! first line to read, local variable + integer last_not_opt ! last line to read, local variable + + !------------------------------------------------------ + + call new_unit(unit) + open(unit, file=file, status='old', action='read', position='rewind') + call prep_file(unit, first, last, first_not_opt, last_not_opt) + allocate(a(last_not_opt - first_not_opt + 1)) + read(unit, fmt=*) a + close(unit) + + end subroutine read_column_char + + !*********************************************************** + + subroutine prep_file(unit, first, last, first_not_opt, last_not_opt) + + ! This subroutine is used by the various versions of "read_column". + ! It fills non-optional arguments: first and last line which will + ! actually be read, taking information from the file itself if necessary. + ! It also positions the input file on the first line to read. + + use opt_merge_m, only: opt_merge + + integer, intent(in):: unit ! logical unit for input file + integer, intent(in), optional:: first ! (first line to read) + integer, intent(in), optional:: last ! (last line to read) + integer, intent(out):: first_not_opt ! (first line to read, not optional) + integer, intent(out):: last_not_opt ! (last line to read, not optional) + + ! Variables local to the subprogram: + integer iostat, i + + !------------------------------------------------------ + + first_not_opt = opt_merge(first, 1) + last_not_opt = opt_merge(last, 0) + + if (last_not_opt == 0) then + ! Count the number of lines in the file: + i = 0 + do + read(unit, fmt=*, iostat=iostat) + if (iostat /= 0) exit + i = i + 1 + end do + last_not_opt = i + if (last_not_opt == 0) stop 'Empty file.' + + rewind(unit) + end if + + ! Go to first line to read: + do i = 1, first_not_opt - 1 + read(unit, fmt=*) + end do + + end subroutine prep_file + +end module read_column_m Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/Max_diff_nc/GNUmakefile =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/Max_diff_nc/GNUmakefile (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/Max_diff_nc/GNUmakefile (revision 1795) @@ -0,0 +1,21 @@ +# This is a makefile for GNU make. +# This makefile builds "max_diff_nc". + +override FFLAGS += $(addprefix -I, ${NETCDF_INC_DIR} ../NetCDF95 ../Jumble) + +override LDLIBS := $(addprefix -L, ../NetCDF95 ../NR_util ../Jumble) -ljumble -lnr_util -lnetcdf95 ${LDLIBS} + +name = max_diff_nc + +%: %.f90 + $(LINK.f) $^ $(LOADLIBES) $(LDLIBS) -o $@ + +.PHONY: all clean +all: ${name} + +clean: + -rm ${name} + +# Some compilers leave an object file if linking fails, so we must be +# able to process this object: +LINK.o = $(FC) $(LDFLAGS) $(TARGET_ARCH) Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/Max_diff_nc/max_diff_nc.f90 =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/Max_diff_nc/max_diff_nc.f90 (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/Max_diff_nc/max_diff_nc.f90 (revision 1795) @@ -0,0 +1,136 @@ +program max_diff_nc + + ! This is a program in Fortran 95. + ! Author: Lionel GUEZ + ! See description in wrapper script. + + ! Maximum memory used will normally be: + + ! -- without computation of average order of magnitude: about 5 + ! times the memory occupied by the largest variable; + + ! -- with computation of average order of magnitude: about 7 times + ! the memory occupied by the largest variable. + + ! This program is meant to be used with a wrapper script so input + ! statements do not have prompts. + + use netcdf95, only: nf95_close, nf95_gw_var, nf95_inq_varid, nf95_inquire, & + nf95_inquire_variable, nf95_open + use netcdf, only: nf90_noerr, nf90_nowrite, nf90_max_name, NF90_FLOAT, & + NF90_double + use jumble, only: compare + + implicit none + + integer, parameter:: wp = kind(0d0) ! working precision + + integer ncid1, ncid2, ncerr, xtype1, ndims + integer nvariables ! number of variables in the first file + integer nvar_comp ! number of variables which will be compared + integer, allocatable:: varid1(:), varid2(:) + real(wp), pointer:: v1_1d(:), v2_1d(:) + real(wp), pointer:: v1_2d(:, :), v2_2d(:, :) + real(wp), pointer:: v1_3d(:, :, :), v2_3d(:, :, :) + real(wp), pointer:: v1_4d(:, :, :, :), v2_4d(:, :, :, :) + character(len=nf90_max_name) name1 + logical same_varid ! compare variables with same varid + logical report_id ! report identical variables + logical comp_mag ! compute avergage order of magnitude + logical quiet + character(len=30+nf90_max_name), allocatable:: tag(:) + integer i + + !---------------------- + + read *, same_varid + read *, report_id + read *, quiet + read *, comp_mag + read "(a)", name1 + + call nf95_open("in1.nc", nf90_nowrite, ncid1) + call nf95_open("in2.nc", nf90_nowrite, ncid2) + + ! Define "nvar_comp", "varid1(:nvar_comp)", "varid2(:nvar_comp)" and + ! "tag(:nvar_comp)": + if (name1 == "") then + ! We want to compare all the variables + call nf95_inquire(ncid1, nvariables=nvariables) + print *, "Found ", nvariables, " variable(s) in the first file." + allocate(varid1(nvariables), varid2(nvariables), tag(nvariables)) + if (same_varid) then + nvar_comp = nvariables + varid1 = (/(i, i = 1, nvariables)/) + varid2 = varid1 + do i = 1, nvariables + call nf95_inquire_variable(ncid1, varid1(i), name1) + tag(i) = 'Variable "' // trim(name1) // '" (name in the first file)' + end do + else + nvar_comp = 0 + do i = 1, nvariables + call nf95_inquire_variable(ncid1, i, name1) + call nf95_inq_varid(ncid2, trim(name1), varid2(nvar_comp + 1), & + ncerr) + if (ncerr == nf90_noerr) then + varid1(nvar_comp + 1) = i + tag(nvar_comp + 1) = 'Variable "' // trim(name1) // '"' + nvar_comp = nvar_comp + 1 + else + print *, 'Could not find "' // trim(name1) & + // '" in the second file. Comparison will be skipped.' + end if + end do + end if + else + nvar_comp = 1 + allocate(varid1(1), varid2(1), tag(1)) + call nf95_inq_varid(ncid1, trim(name1), varid1(1)) + call nf95_inq_varid(ncid2, trim(name1), varid2(1)) + tag(1) = 'Variable "' // trim(name1) // '"' + end if + + do i = 1, nvar_comp + call nf95_inquire_variable(ncid1, varid1(i), xtype=xtype1, ndims=ndims) + if (xtype1 == nf90_float .or. xtype1 == nf90_double) then + select case (ndims) + case (1) + call nf95_gw_var(ncid1, varid1(i), v1_1d) + call nf95_gw_var(ncid2, varid2(i), v2_1d) + call compare(v1_1d, v2_1d, trim(tag(i)), comp_mag, report_id, quiet) + deallocate(v1_1d, v2_1d) + case (2) + call nf95_gw_var(ncid1, varid1(i), v1_2d) + call nf95_gw_var(ncid2, varid2(i), v2_2d) + call compare(v1_2d, v2_2d, trim(tag(i)), comp_mag, report_id, quiet) + deallocate(v1_2d, v2_2d) + case (3) + call nf95_gw_var(ncid1, varid1(i), v1_3d) + call nf95_gw_var(ncid2, varid2(i), v2_3d) + call compare(v1_3d, v2_3d, trim(tag(i)), comp_mag, report_id, quiet) + deallocate(v1_3d, v2_3d) + case (4) + call nf95_gw_var(ncid1, varid1(i), v1_4d) + call nf95_gw_var(ncid2, varid2(i), v2_4d) + call compare(v1_4d, v2_4d, trim(tag(i)), comp_mag, report_id, quiet) + deallocate(v1_4d, v2_4d) + case default + print * + print *, "******************" + print *, trim(tag(i)) // ":" + print *, "Rank not supported." + print *, "ndims = ", ndims + end select + else + print * + print *, "******************" + print *, trim(tag(i)) // ":" + print *, 'Not of type "nf90_float or "nf90_double".' + end if + end do + + call nf95_close(ncid1) + call nf95_close(ncid2) + +end program max_diff_nc Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/Max_diff_nc/max_diff_nc.sh =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/Max_diff_nc/max_diff_nc.sh (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/Max_diff_nc/max_diff_nc.sh (revision 1795) @@ -0,0 +1,82 @@ +# This is a script in Bash. +# Author: Lionel GUEZ + +# This script is a wrapper for a Fortran program. It compares NetCDF +# variables of type "NF90_FLOAT" or "NF90_DOUBLE" and of rank 1 to 4 +# in two files. NetCDF variables are read into Fortran variables of +# type "double precision". The program computes the maximum of the +# absolute value of the difference and the maximum of the absolute +# value of the relative difference. The program either compares +# variables with the same varid or variables with the same name, or a +# single variable selected on the command line. Compared variables are +# assumed to have the same type. The program checks that compared +# variables have the same shape. + +USAGE="usage: `basename $0` [-s] [-m] [-i] [-q] [-v ] file file + -s: report identical variables + -m: compute average order of magnitude + -i: compare variables with same varid, regardless of variable name + -q: only report names of variables which differ, without maximum difference + -v : only compare variable " + +# Default values: +report_id=False +comp_mag=False +same_varid=False +quiet=False + +while getopts smiqv: name + do + case $name in + s) report_id=True;; + m) comp_mag=True;; + i) same_varid=True;; + q) quiet=True;; + v) name1=$OPTARG;; + \?) echo "$USAGE" + exit 1;; + esac +done + +shift $((OPTIND - 1)) + +if (($# != 2)) + then + echo "$USAGE" >&2 + exit 1 +fi + +##set -x +trap 'exit 1' ERR + +for argument in $* + do + if [[ ! -f $argument ]] + then + echo "$argument not found" + exit 1 + fi +done + +max_diff_nc=... +if [[ ! -x $max_diff_nc ]] + then + echo "Fortran executable not found" >&2 + exit 1 +fi + +ln -s $1 in1.nc +ln -s $2 in2.nc + +trap 'rm in1.nc in2.nc; exit 1' ERR + +# Run the Fortran program: +$max_diff_nc <depend.mk + +clean: + rm -f ${lib} ${objects} + +# Dependencies between object files and include files: +include depend.mk Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/array_copy.f90 =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/array_copy.f90 (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/array_copy.f90 (revision 1795) @@ -0,0 +1,40 @@ +MODULE array_copy_m + + IMPLICIT NONE + + INTERFACE array_copy + MODULE PROCEDURE array_copy_r, array_copy_d, array_copy_i + END INTERFACE + + private array_copy_r, array_copy_d, array_copy_i + +CONTAINS + + SUBROUTINE array_copy_r(src,dest,n_copied,n_not_copied) + REAL, DIMENSION(:), INTENT(IN) :: src + REAL, DIMENSION(:), INTENT(OUT) :: dest + INTEGER, INTENT(OUT) :: n_copied, n_not_copied + n_copied=min(size(src),size(dest)) + n_not_copied=size(src)-n_copied + dest(1:n_copied)=src(1:n_copied) + END SUBROUTINE array_copy_r + !BL + SUBROUTINE array_copy_d(src,dest,n_copied,n_not_copied) + double precision, DIMENSION(:), INTENT(IN) :: src + double precision, DIMENSION(:), INTENT(OUT) :: dest + INTEGER, INTENT(OUT) :: n_copied, n_not_copied + n_copied=min(size(src),size(dest)) + n_not_copied=size(src)-n_copied + dest(1:n_copied)=src(1:n_copied) + END SUBROUTINE array_copy_d + !BL + SUBROUTINE array_copy_i(src,dest,n_copied,n_not_copied) + INTEGER, DIMENSION(:), INTENT(IN) :: src + INTEGER, DIMENSION(:), INTENT(OUT) :: dest + INTEGER, INTENT(OUT) :: n_copied, n_not_copied + n_copied=min(size(src),size(dest)) + n_not_copied=size(src)-n_copied + dest(1:n_copied)=src(1:n_copied) + END SUBROUTINE array_copy_i + +END MODULE array_copy_m Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/arth.f90 =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/arth.f90 (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/arth.f90 (revision 1795) @@ -0,0 +1,111 @@ +MODULE arth_m + + IMPLICIT NONE + + INTEGER, PARAMETER, private:: NPAR_ARTH=16, NPAR2_ARTH=8 + + INTERFACE arth + ! Returns an arithmetic progression, given a first term "first", an + ! increment and a number of terms "n". + + MODULE PROCEDURE arth_r, arth_d, arth_i + ! The difference between the procedures is the kind and type of + ! arguments first and increment and of function result. + END INTERFACE + + private arth_r, arth_d, arth_i + +CONTAINS + + pure FUNCTION arth_r(first,increment,n) + + + REAL, INTENT(IN) :: first,increment + INTEGER, INTENT(IN) :: n + REAL, DIMENSION(n) :: arth_r + + ! Variables local to the procedure: + + INTEGER :: k,k2 + REAL :: temp + + !--------------------------------------- + + if (n > 0) arth_r(1)=first + if (n <= NPAR_ARTH) then + do k=2,n + arth_r(k)=arth_r(k-1)+increment + end do + else + do k=2,NPAR2_ARTH + arth_r(k)=arth_r(k-1)+increment + end do + temp=increment*NPAR2_ARTH + k=NPAR2_ARTH + do + if (k >= n) exit + k2=k+k + arth_r(k+1:min(k2,n)) = temp + arth_r(1:min(k,n-k)) + temp=temp+temp + k=k2 + end do + end if + END FUNCTION arth_r + + !************************************* + + pure FUNCTION arth_d(first,increment,n) + DOUBLE PRECISION, INTENT(IN) :: first,increment + INTEGER, INTENT(IN) :: n + DOUBLE PRECISION, DIMENSION(n) :: arth_d + INTEGER :: k,k2 + DOUBLE PRECISION :: temp + if (n > 0) arth_d(1)=first + if (n <= NPAR_ARTH) then + do k=2,n + arth_d(k)=arth_d(k-1)+increment + end do + else + do k=2,NPAR2_ARTH + arth_d(k)=arth_d(k-1)+increment + end do + temp=increment*NPAR2_ARTH + k=NPAR2_ARTH + do + if (k >= n) exit + k2=k+k + arth_d(k+1:min(k2,n))=temp+arth_d(1:min(k,n-k)) + temp=temp+temp + k=k2 + end do + end if + END FUNCTION arth_d + + !************************************* + + pure FUNCTION arth_i(first,increment,n) + INTEGER, INTENT(IN) :: first,increment,n + INTEGER, DIMENSION(n) :: arth_i + INTEGER :: k,k2,temp + if (n > 0) arth_i(1)=first + if (n <= NPAR_ARTH) then + do k=2,n + arth_i(k)=arth_i(k-1)+increment + end do + else + do k=2,NPAR2_ARTH + arth_i(k)=arth_i(k-1)+increment + end do + temp=increment*NPAR2_ARTH + k=NPAR2_ARTH + do + if (k >= n) exit + k2=k+k + arth_i(k+1:min(k2,n))=temp+arth_i(1:min(k,n-k)) + temp=temp+temp + k=k2 + end do + end if + END FUNCTION arth_i + +END MODULE arth_m Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/assert.f90 =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/assert.f90 (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/assert.f90 (revision 1795) @@ -0,0 +1,71 @@ +MODULE assert_m + + implicit none + + INTERFACE assert + MODULE PROCEDURE assert1,assert2,assert3,assert4,assert_v + END INTERFACE + + private assert1,assert2,assert3,assert4,assert_v + +CONTAINS + + SUBROUTINE assert1(n1,string) + CHARACTER(LEN=*), INTENT(IN) :: string + LOGICAL, INTENT(IN) :: n1 + if (.not. n1) then + print *, 'An assertion failed with this tag: ' // string + print *, 'program terminated by assert1' + stop 1 + end if + END SUBROUTINE assert1 + + !**************************** + + SUBROUTINE assert2(n1,n2,string) + CHARACTER(LEN=*), INTENT(IN) :: string + LOGICAL, INTENT(IN) :: n1,n2 + if (.not. (n1 .and. n2)) then + print *, 'An assertion failed with this tag: ' // string + print *, 'program terminated by assert2' + stop 1 + end if + END SUBROUTINE assert2 + + !**************************** + + SUBROUTINE assert3(n1,n2,n3,string) + CHARACTER(LEN=*), INTENT(IN) :: string + LOGICAL, INTENT(IN) :: n1,n2,n3 + if (.not. (n1 .and. n2 .and. n3)) then + print *, 'An assertion failed with this tag: ' // string + print *, 'program terminated by assert3' + stop 1 + end if + END SUBROUTINE assert3 + + !**************************** + + SUBROUTINE assert4(n1,n2,n3,n4,string) + CHARACTER(LEN=*), INTENT(IN) :: string + LOGICAL, INTENT(IN) :: n1,n2,n3,n4 + if (.not. (n1 .and. n2 .and. n3 .and. n4)) then + print *, 'An assertion failed with this tag: ' // string + print *, 'program terminated by assert4' + stop 1 + end if + END SUBROUTINE assert4 + + !**************************** + + SUBROUTINE assert_v(n,string) + CHARACTER(LEN=*), INTENT(IN) :: string + LOGICAL, DIMENSION(:), INTENT(IN) :: n + if (.not. all(n)) then + print *, 'An assertion failed with this tag: ' // string + print *, 'program terminated by assert_v' + stop 1 + end if + END SUBROUTINE assert_v + +END MODULE assert_m Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/assert_eq.f90 =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/assert_eq.f90 (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/assert_eq.f90 (revision 1795) @@ -0,0 +1,69 @@ +MODULE assert_eq_m + + implicit none + + INTERFACE assert_eq + MODULE PROCEDURE assert_eq2,assert_eq3,assert_eq4,assert_eqn + END INTERFACE + + private assert_eq2,assert_eq3,assert_eq4,assert_eqn + +CONTAINS + + FUNCTION assert_eq2(n1,n2,string) + CHARACTER(LEN=*), INTENT(IN) :: string + INTEGER, INTENT(IN) :: n1,n2 + INTEGER assert_eq2 + if (n1 == n2) then + assert_eq2=n1 + else + write (*,*) 'nrerror: an assert_eq failed with this tag: ', & + string + print *, 'program terminated by assert_eq2' + stop 1 + end if + END FUNCTION assert_eq2 + !BL + FUNCTION assert_eq3(n1,n2,n3,string) + CHARACTER(LEN=*), INTENT(IN) :: string + INTEGER, INTENT(IN) :: n1,n2,n3 + INTEGER assert_eq3 + if (n1 == n2 .and. n2 == n3) then + assert_eq3=n1 + else + write (*,*) 'nrerror: an assert_eq failed with this tag: ', & + string + print *, 'program terminated by assert_eq3' + stop 1 + end if + END FUNCTION assert_eq3 + !BL + FUNCTION assert_eq4(n1,n2,n3,n4,string) + CHARACTER(LEN=*), INTENT(IN) :: string + INTEGER, INTENT(IN) :: n1,n2,n3,n4 + INTEGER assert_eq4 + if (n1 == n2 .and. n2 == n3 .and. n3 == n4) then + assert_eq4=n1 + else + write (*,*) 'nrerror: an assert_eq failed with this tag: ', & + string + print *, 'program terminated by assert_eq4' + stop 1 + end if + END FUNCTION assert_eq4 + !BL + FUNCTION assert_eqn(nn,string) + CHARACTER(LEN=*), INTENT(IN) :: string + INTEGER, DIMENSION(:), INTENT(IN) :: nn + INTEGER assert_eqn + if (all(nn(2:) == nn(1))) then + assert_eqn=nn(1) + else + write (*,*) 'nrerror: an assert_eq failed with this tag: ', & + string + print *, 'program terminated by assert_eqn' + stop 1 + end if + END FUNCTION assert_eqn + +END MODULE assert_eq_m Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/cumprod.f90 =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/cumprod.f90 (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/cumprod.f90 (revision 1795) @@ -0,0 +1,40 @@ +module cumprod_m + + implicit none + +contains + + RECURSIVE FUNCTION cumprod(arr,seed) RESULT(ans) + + ! Cumulative product on an array, with optional multiplicative seed. + + USE nrtype, only: wp + + REAL(WP), DIMENSION(:), INTENT(IN) :: arr + REAL(WP), OPTIONAL, INTENT(IN) :: seed + REAL(WP), DIMENSION(size(arr)) :: ans + + ! Local: + INTEGER n,j + REAL(WP) :: sd + INTEGER, PARAMETER :: NPAR_CUMPROD=8 + + !-------------------------------------------------------- + + n=size(arr) + if (n == 0) RETURN + sd=1.0_wp + if (present(seed)) sd=seed + ans(1)=arr(1)*sd + if (n < NPAR_CUMPROD) then + do j=2,n + ans(j)=ans(j-1)*arr(j) + end do + else + ans(2:n:2)=cumprod(arr(2:n:2)*arr(1:n-1:2),sd) + ans(3:n:2)=ans(2:n-1:2)*arr(3:n:2) + end if + + END FUNCTION cumprod + +end module cumprod_m Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/cumsum.f90 =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/cumsum.f90 (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/cumsum.f90 (revision 1795) @@ -0,0 +1,56 @@ +MODULE cumsum_m + + IMPLICIT NONE + + INTEGER, PARAMETER, private :: NPAR_CUMSUM=16 + + INTERFACE cumsum + MODULE PROCEDURE cumsum_r,cumsum_i + END INTERFACE + + private cumsum_r,cumsum_i + +CONTAINS + + RECURSIVE FUNCTION cumsum_r(arr,seed) RESULT(ans) + REAL, DIMENSION(:), INTENT(IN) :: arr + REAL, OPTIONAL, INTENT(IN) :: seed + REAL, DIMENSION(size(arr)) :: ans + INTEGER :: n,j + REAL :: sd + n=size(arr) + if (n == 0) RETURN + sd=0.0 + if (present(seed)) sd=seed + ans(1)=arr(1)+sd + if (n < NPAR_CUMSUM) then + do j=2,n + ans(j)=ans(j-1)+arr(j) + end do + else + ans(2:n:2)=cumsum_r(arr(2:n:2)+arr(1:n-1:2),sd) + ans(3:n:2)=ans(2:n-1:2)+arr(3:n:2) + end if + END FUNCTION cumsum_r + !BL + RECURSIVE FUNCTION cumsum_i(arr,seed) RESULT(ans) + INTEGER, DIMENSION(:), INTENT(IN) :: arr + INTEGER, OPTIONAL, INTENT(IN) :: seed + INTEGER, DIMENSION(size(arr)) :: ans + INTEGER :: n,j,sd + n=size(arr) + if (n == 0) RETURN + sd=0 + if (present(seed)) sd=seed + ans(1)=arr(1)+sd + if (n < NPAR_CUMSUM) then + do j=2,n + ans(j)=ans(j-1)+arr(j) + end do + else + ans(2:n:2)=cumsum_i(arr(2:n:2)+arr(1:n-1:2),sd) + ans(3:n:2)=ans(2:n-1:2)+arr(3:n:2) + end if + END FUNCTION cumsum_i + +END MODULE cumsum_m Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/depend.mk =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/depend.mk (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/depend.mk (revision 1795) @@ -0,0 +1,20 @@ +array_copy.o : nrtype.o +arth.o : nrtype.o +cumprod.o : nrtype.o +cumsum.o : nrtype.o +diagadd.o : assert_eq.o +diagmult.o : assert_eq.o +geop.o : nrtype.o +get_diag.o : assert_eq.o +lower_triangle.o : outerdiff.o arth.o +nr_util.o : zroots_unity.o vabs.o upper_triangle.o unit_matrix.o swap.o scatter_max.o scatter_add.o reallocate.o put_diag.o poly_term.o poly.o outersum.o outerprod.o outerdiv.o outerdiff.o outerand.o nrtype.o nrerror.o lower_triangle.o iminloc.o imaxloc.o ifirstloc.o get_diag.o geop.o diagmult.o diagadd.o cumsum.o cumprod.o assert.o assert_eq.o arth.o array_copy.o +outerdiv.o : nrtype.o +outersum.o : nrtype.o +put_diag.o : assert_eq.o +scatter_add.o : assert_eq.o +scatter_max.o : assert_eq.o +swap.o : nrtype.o +unit_matrix.o : nrtype.o +upper_triangle.o : outerdiff.o arth.o +vabs.o : nrtype.o +zroots_unity.o : nrtype.o Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/diagadd.f90 =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/diagadd.f90 (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/diagadd.f90 (revision 1795) @@ -0,0 +1,34 @@ +MODULE diagadd_m + + IMPLICIT NONE + + INTERFACE diagadd + MODULE PROCEDURE diagadd_rv,diagadd_r + END INTERFACE + + private diagadd_rv,diagadd_r + +CONTAINS + + SUBROUTINE diagadd_rv(mat,diag) + use assert_eq_m, only: assert_eq + REAL, DIMENSION(:,:), INTENT(INOUT) :: mat + REAL, DIMENSION(:), INTENT(IN) :: diag + INTEGER :: j,n + n = assert_eq(size(diag),min(size(mat,1),size(mat,2)),'diagadd_rv') + do j=1,n + mat(j,j)=mat(j,j)+diag(j) + end do + END SUBROUTINE diagadd_rv + !BL + SUBROUTINE diagadd_r(mat,diag) + REAL, DIMENSION(:,:), INTENT(INOUT) :: mat + REAL, INTENT(IN) :: diag + INTEGER :: j,n + n = min(size(mat,1),size(mat,2)) + do j=1,n + mat(j,j)=mat(j,j)+diag + end do + END SUBROUTINE diagadd_r + +END MODULE diagadd_m Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/diagmult.f90 =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/diagmult.f90 (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/diagmult.f90 (revision 1795) @@ -0,0 +1,34 @@ +MODULE diagmult_m + + IMPLICIT NONE + + INTERFACE diagmult + MODULE PROCEDURE diagmult_rv,diagmult_r + END INTERFACE + + private diagmult_rv,diagmult_r + +CONTAINS + + SUBROUTINE diagmult_rv(mat,diag) + use assert_eq_m, only: assert_eq + REAL, DIMENSION(:,:), INTENT(INOUT) :: mat + REAL, DIMENSION(:), INTENT(IN) :: diag + INTEGER j,n + n = assert_eq(size(diag),min(size(mat,1),size(mat,2)),'diagmult_rv') + do j=1,n + mat(j,j)=mat(j,j)*diag(j) + end do + END SUBROUTINE diagmult_rv + !BL + SUBROUTINE diagmult_r(mat,diag) + REAL, DIMENSION(:,:), INTENT(INOUT) :: mat + REAL, INTENT(IN) :: diag + INTEGER j,n + n = min(size(mat,1),size(mat,2)) + do j=1,n + mat(j,j)=mat(j,j)*diag + end do + END SUBROUTINE diagmult_r + +END MODULE diagmult_m Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/geop.f90 =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/geop.f90 (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/geop.f90 (revision 1795) @@ -0,0 +1,148 @@ +MODULE geop_m + + IMPLICIT NONE + + INTERFACE geop + MODULE PROCEDURE geop_r, geop_d, geop_i, geop_c, geop_dv + END INTERFACE + + INTEGER, PARAMETER, private :: NPAR_GEOP=4, NPAR2_GEOP=2 + private geop_r, geop_d, geop_i, geop_c, geop_dv + +CONTAINS + + FUNCTION geop_r(first,factor,n) + REAL, INTENT(IN) :: first,factor + INTEGER, INTENT(IN) :: n + REAL, DIMENSION(n) :: geop_r + INTEGER :: k,k2 + REAL :: temp + if (n > 0) geop_r(1)=first + if (n <= NPAR_GEOP) then + do k=2,n + geop_r(k)=geop_r(k-1)*factor + end do + else + do k=2,NPAR2_GEOP + geop_r(k)=geop_r(k-1)*factor + end do + temp=factor**NPAR2_GEOP + k=NPAR2_GEOP + do + if (k >= n) exit + k2=k+k + geop_r(k+1:min(k2,n))=temp*geop_r(1:min(k,n-k)) + temp=temp*temp + k=k2 + end do + end if + END FUNCTION geop_r + !BL + FUNCTION geop_d(first,factor,n) + DOUBLE PRECISION, INTENT(IN) :: first,factor + INTEGER, INTENT(IN) :: n + DOUBLE PRECISION, DIMENSION(n) :: geop_d + INTEGER :: k,k2 + DOUBLE PRECISION :: temp + if (n > 0) geop_d(1)=first + if (n <= NPAR_GEOP) then + do k=2,n + geop_d(k)=geop_d(k-1)*factor + end do + else + do k=2,NPAR2_GEOP + geop_d(k)=geop_d(k-1)*factor + end do + temp=factor**NPAR2_GEOP + k=NPAR2_GEOP + do + if (k >= n) exit + k2=k+k + geop_d(k+1:min(k2,n))=temp*geop_d(1:min(k,n-k)) + temp=temp*temp + k=k2 + end do + end if + END FUNCTION geop_d + !BL + FUNCTION geop_i(first,factor,n) + INTEGER, INTENT(IN) :: first,factor,n + INTEGER, DIMENSION(n) :: geop_i + INTEGER :: k,k2,temp + if (n > 0) geop_i(1)=first + if (n <= NPAR_GEOP) then + do k=2,n + geop_i(k)=geop_i(k-1)*factor + end do + else + do k=2,NPAR2_GEOP + geop_i(k)=geop_i(k-1)*factor + end do + temp=factor**NPAR2_GEOP + k=NPAR2_GEOP + do + if (k >= n) exit + k2=k+k + geop_i(k+1:min(k2,n))=temp*geop_i(1:min(k,n-k)) + temp=temp*temp + k=k2 + end do + end if + END FUNCTION geop_i + !BL + FUNCTION geop_c(first,factor,n) + COMPLEX, INTENT(IN) :: first,factor + INTEGER, INTENT(IN) :: n + COMPLEX, DIMENSION(n) :: geop_c + INTEGER :: k,k2 + COMPLEX :: temp + if (n > 0) geop_c(1)=first + if (n <= NPAR_GEOP) then + do k=2,n + geop_c(k)=geop_c(k-1)*factor + end do + else + do k=2,NPAR2_GEOP + geop_c(k)=geop_c(k-1)*factor + end do + temp=factor**NPAR2_GEOP + k=NPAR2_GEOP + do + if (k >= n) exit + k2=k+k + geop_c(k+1:min(k2,n))=temp*geop_c(1:min(k,n-k)) + temp=temp*temp + k=k2 + end do + end if + END FUNCTION geop_c + !BL + FUNCTION geop_dv(first,factor,n) + DOUBLE PRECISION, DIMENSION(:), INTENT(IN) :: first,factor + INTEGER, INTENT(IN) :: n + DOUBLE PRECISION, DIMENSION(size(first),n) :: geop_dv + INTEGER :: k,k2 + DOUBLE PRECISION, DIMENSION(size(first)) :: temp + if (n > 0) geop_dv(:,1)=first(:) + if (n <= NPAR_GEOP) then + do k=2,n + geop_dv(:,k)=geop_dv(:,k-1)*factor(:) + end do + else + do k=2,NPAR2_GEOP + geop_dv(:,k)=geop_dv(:,k-1)*factor(:) + end do + temp=factor**NPAR2_GEOP + k=NPAR2_GEOP + do + if (k >= n) exit + k2=k+k + geop_dv(:,k+1:min(k2,n))=geop_dv(:,1:min(k,n-k))*& + spread(temp,2,size(geop_dv(:,1:min(k,n-k)),2)) + temp=temp*temp + k=k2 + end do + end if + END FUNCTION geop_dv + +END MODULE geop_m Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/get_diag.f90 =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/get_diag.f90 (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/get_diag.f90 (revision 1795) @@ -0,0 +1,35 @@ +MODULE get_diag_m + + IMPLICIT NONE + + INTERFACE get_diag + MODULE PROCEDURE get_diag_rv, get_diag_dv + END INTERFACE + + private get_diag_rv, get_diag_dv + +CONTAINS + + FUNCTION get_diag_rv(mat) + use assert_eq_m, only: assert_eq + REAL, DIMENSION(:,:), INTENT(IN) :: mat + REAL, DIMENSION(size(mat,1)) :: get_diag_rv + INTEGER j, n + n=assert_eq(size(mat,1),size(mat,2),'get_diag_rv') + do j=1, n + get_diag_rv(j)=mat(j,j) + end do + END FUNCTION get_diag_rv + !BL + FUNCTION get_diag_dv(mat) + use assert_eq_m, only: assert_eq + DOUBLE PRECISION, DIMENSION(:,:), INTENT(IN) :: mat + DOUBLE PRECISION, DIMENSION(size(mat,1)) :: get_diag_dv + INTEGER j, n + n=assert_eq(size(mat,1),size(mat,2),'get_diag_dv') + do j=1, n + get_diag_dv(j)=mat(j,j) + end do + END FUNCTION get_diag_dv + +END MODULE get_diag_m Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/ifirstloc.f90 =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/ifirstloc.f90 (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/ifirstloc.f90 (revision 1795) @@ -0,0 +1,35 @@ +module ifirstloc_m + + implicit none + +contains + + INTEGER FUNCTION ifirstloc(mask) + + ! Location of first true value in a logical array, returned as an + ! integer. Returns size(mask)+1 if mask has zero element or all + ! elements of mask are false. So the result is always >= 1. + + ! See notes on programming choices. + + LOGICAL, INTENT(IN):: mask(:) + + ! Local: + integer n + + !------------------------------------------------------- + + n = size(mask) + ifirstloc = 1 + + if (n >= 1) then + do while (ifirstloc <= n - 1 .and. .not. mask(ifirstloc)) + ifirstloc = ifirstloc + 1 + end do + ! {1 <= ifirstloc <= n} + if (.not. mask(ifirstloc)) ifirstloc = n + 1 + end if + + END FUNCTION ifirstloc + +end module ifirstloc_m Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/imaxloc.f90 =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/imaxloc.f90 (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/imaxloc.f90 (revision 1795) @@ -0,0 +1,25 @@ +MODULE imaxloc_m + + ! Useless in Fortran 95. + + implicit none + + INTERFACE imaxloc + MODULE PROCEDURE imaxloc_r,imaxloc_i + END INTERFACE + + private imaxloc_r,imaxloc_i + +CONTAINS + + INTEGER FUNCTION imaxloc_r(arr) + REAL, INTENT(IN) :: arr(:) + imaxloc_r=maxloc(arr, dim=1) + END FUNCTION imaxloc_r + + INTEGER FUNCTION imaxloc_i(iarr) + INTEGER, INTENT(IN) :: iarr(:) + imaxloc_i=maxloc(iarr, dim=1) + END FUNCTION imaxloc_i + +END MODULE imaxloc_m Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/iminloc.f90 =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/iminloc.f90 (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/iminloc.f90 (revision 1795) @@ -0,0 +1,27 @@ +module iminloc_m + + ! Useless in Fortran 95. + + implicit none + + interface iminloc + module procedure iminloc_r, iminloc_d + end interface + + private iminloc_r, iminloc_d + +contains + + pure INTEGER FUNCTION iminloc_r(arr) + REAL, INTENT(IN) :: arr(:) + iminloc_r=minloc(arr, dim=1) + END FUNCTION iminloc_r + + !****************************************************** + + pure INTEGER FUNCTION iminloc_d(arr) + DOUBLE PRECISION, INTENT(IN) :: arr(:) + iminloc_d=minloc(arr, dim=1) + END FUNCTION iminloc_d + +end module iminloc_m Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/lower_triangle.f90 =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/lower_triangle.f90 (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/lower_triangle.f90 (revision 1795) @@ -0,0 +1,19 @@ +module lower_triangle_m + + implicit none + +contains + + FUNCTION lower_triangle(j,k,extra) + use arth_m, only: arth + use outerdiff_m, only: outerdiff + INTEGER, INTENT(IN) :: j,k + INTEGER, OPTIONAL, INTENT(IN) :: extra + LOGICAL, DIMENSION(j,k) :: lower_triangle + INTEGER :: n + n=0 + if (present(extra)) n=extra + lower_triangle=(outerdiff(arth(1,1,j),arth(1,1,k)) > -n) + END FUNCTION lower_triangle + +end module lower_triangle_m Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/nr_util.f90 =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/nr_util.f90 (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/nr_util.f90 (revision 1795) @@ -0,0 +1,38 @@ +module nr_util + + ! Modules used in alphabetical order. + + use array_copy_m + use arth_m + use assert_eq_m + use assert_m + use cumprod_m + use cumsum_m + use diagadd_m + use diagmult_m + use geop_m + use get_diag_m + use ifirstloc_m + use imaxloc_m + use iminloc_m + use lower_triangle_m + use nrerror_m + use nrtype + use outerand_m + use outerdiff_m + use outerdiv_m + use outerprod_m + use outersum_m + use poly_m + use poly_term_m + use put_diag_m + use reallocate_m + use scatter_add_m + use scatter_max_m + use swap_m + use unit_matrix_m + use upper_triangle_m + use vabs_m + use zroots_unity_m + +end module nr_util Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/nrerror.f90 =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/nrerror.f90 (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/nrerror.f90 (revision 1795) @@ -0,0 +1,17 @@ +module nrerror_m + + implicit none + +contains + + SUBROUTINE nrerror(string) + + CHARACTER(LEN=*), INTENT(IN) :: string + + print *, 'nrerror: ', string + print *, 'program terminated by nrerror' + stop 1 + + END SUBROUTINE nrerror + +end module nrerror_m Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/nrtype.f90 =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/nrtype.f90 (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/nrtype.f90 (revision 1795) @@ -0,0 +1,39 @@ +MODULE nrtype + + implicit none + + integer, parameter:: wp = kind(0.) ! working precision for real type + + ! Frequently used mathematical constants (with precision to spare): + + REAL, PARAMETER :: PI=3.141592653589793238462643383279502884197 + REAL, PARAMETER :: PIO2=1.57079632679489661923132169163975144209858 + REAL, PARAMETER :: TWOPI=6.283185307179586476925286766559005768394 + REAL, PARAMETER :: SQRT2=1.41421356237309504880168872420969807856967 + REAL, PARAMETER :: EULER=0.5772156649015328606065120900824024310422 + + DOUBLE PRECISION, PARAMETER:: & + PI_D = 3.141592653589793238462643383279502884197d0 + DOUBLE PRECISION, PARAMETER:: & + PIO2_D=1.57079632679489661923132169163975144209858d0 + DOUBLE PRECISION, PARAMETER:: & + TWOPI_D=6.283185307179586476925286766559005768394d0 + + ! Derived data types for sparse matrices, single and double + ! precision (see use in Chapter B2): + + TYPE sprs2_sp + INTEGER :: n,len + REAL, DIMENSION(:), POINTER :: val + INTEGER, DIMENSION(:), POINTER :: irow + INTEGER, DIMENSION(:), POINTER :: jcol + END TYPE sprs2_sp + + TYPE sprs2_dp + INTEGER :: n,len + DOUBLE PRECISION, DIMENSION(:), POINTER :: val + INTEGER, DIMENSION(:), POINTER :: irow + INTEGER, DIMENSION(:), POINTER :: jcol + END TYPE sprs2_dp + +END MODULE nrtype Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/outerand.f90 =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/outerand.f90 (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/outerand.f90 (revision 1795) @@ -0,0 +1,14 @@ +module outerand_m + + implicit none + +contains + + FUNCTION outerand(a,b) + LOGICAL, DIMENSION(:), INTENT(IN) :: a,b + LOGICAL, DIMENSION(size(a),size(b)) :: outerand + outerand = spread(a,dim=2,ncopies=size(b)) .and. & + spread(b,dim=1,ncopies=size(a)) + END FUNCTION outerand + +end module outerand_m Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/outerdiff.f90 =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/outerdiff.f90 (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/outerdiff.f90 (revision 1795) @@ -0,0 +1,34 @@ +MODULE outerdiff_m + + IMPLICIT NONE + + INTERFACE outerdiff + MODULE PROCEDURE outerdiff_r,outerdiff_d,outerdiff_i + END INTERFACE + + private outerdiff_r,outerdiff_d,outerdiff_i + +CONTAINS + + FUNCTION outerdiff_r(a,b) + REAL, DIMENSION(:), INTENT(IN) :: a,b + REAL, DIMENSION(size(a),size(b)) :: outerdiff_r + outerdiff_r = spread(a,dim=2,ncopies=size(b)) - & + spread(b,dim=1,ncopies=size(a)) + END FUNCTION outerdiff_r + !BL + FUNCTION outerdiff_d(a,b) + DOUBLE PRECISION, DIMENSION(:), INTENT(IN) :: a,b + DOUBLE PRECISION, DIMENSION(size(a),size(b)) :: outerdiff_d + outerdiff_d = spread(a,dim=2,ncopies=size(b)) - & + spread(b,dim=1,ncopies=size(a)) + END FUNCTION outerdiff_d + !BL + FUNCTION outerdiff_i(a,b) + INTEGER, DIMENSION(:), INTENT(IN) :: a,b + INTEGER, DIMENSION(size(a),size(b)) :: outerdiff_i + outerdiff_i = spread(a,dim=2,ncopies=size(b)) - & + spread(b,dim=1,ncopies=size(a)) + END FUNCTION outerdiff_i + +END MODULE outerdiff_m Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/outerdiv.f90 =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/outerdiv.f90 (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/outerdiv.f90 (revision 1795) @@ -0,0 +1,16 @@ +module outerdiv_m + + implicit none + +contains + + FUNCTION outerdiv(a,b) + ! Returns a matrix that is the outer quotient of two vectors. + USE nrtype, only: wp + REAL(WP), DIMENSION(:), INTENT(IN) :: a,b + REAL(WP), DIMENSION(size(a),size(b)) :: outerdiv + outerdiv = spread(a,dim=2,ncopies=size(b)) / & + spread(b,dim=1,ncopies=size(a)) + END FUNCTION outerdiv + +end module outerdiv_m Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/outerprod.f90 =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/outerprod.f90 (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/outerprod.f90 (revision 1795) @@ -0,0 +1,27 @@ +MODULE outerprod_m + + IMPLICIT NONE + + INTERFACE outerprod + MODULE PROCEDURE outerprod_r,outerprod_d + END INTERFACE + + private outerprod_r,outerprod_d + +CONTAINS + + FUNCTION outerprod_r(a,b) + REAL, DIMENSION(:), INTENT(IN) :: a,b + REAL, DIMENSION(size(a),size(b)) :: outerprod_r + outerprod_r = spread(a,dim=2,ncopies=size(b)) * & + spread(b,dim=1,ncopies=size(a)) + END FUNCTION outerprod_r + !BL + FUNCTION outerprod_d(a,b) + DOUBLE PRECISION, DIMENSION(:), INTENT(IN) :: a,b + DOUBLE PRECISION, DIMENSION(size(a),size(b)) :: outerprod_d + outerprod_d = spread(a,dim=2,ncopies=size(b)) * & + spread(b,dim=1,ncopies=size(a)) + END FUNCTION outerprod_d + +END MODULE outerprod_m Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/outersum.f90 =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/outersum.f90 (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/outersum.f90 (revision 1795) @@ -0,0 +1,15 @@ +module outersum_m + + implicit none + +contains + + FUNCTION outersum(a,b) + USE nrtype, only: wp + REAL(WP), DIMENSION(:), INTENT(IN) :: a,b + REAL(WP), DIMENSION(size(a),size(b)) :: outersum + outersum = spread(a,dim=2,ncopies=size(b)) + & + spread(b,dim=1,ncopies=size(a)) + END FUNCTION outersum + +end module outersum_m Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/poly.f90 =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/poly.f90 (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/poly.f90 (revision 1795) @@ -0,0 +1,199 @@ +MODULE poly_m + + IMPLICIT NONE + + INTEGER, PARAMETER, private :: NPAR_POLY=8 + + INTERFACE poly + MODULE PROCEDURE poly_rr,poly_rrv,poly_dd,poly_ddv,& + poly_rc,poly_cc,poly_msk_rrv,poly_msk_ddv + END INTERFACE + + private poly_rr,poly_rrv,poly_dd,poly_ddv,& + poly_rc,poly_cc,poly_msk_rrv,poly_msk_ddv + +CONTAINS + + FUNCTION poly_rr(x,coeffs) + REAL, INTENT(IN) :: x + REAL, DIMENSION(:), INTENT(IN) :: coeffs + REAL :: poly_rr + REAL :: pow + REAL, DIMENSION(:), ALLOCATABLE :: vec + INTEGER :: i,n,nn + n=size(coeffs) + if (n <= 0) then + poly_rr=0.0 + else if (n < NPAR_POLY) then + poly_rr=coeffs(n) + do i=n-1,1,-1 + poly_rr=x*poly_rr+coeffs(i) + end do + else + allocate(vec(n+1)) + pow=x + vec(1:n)=coeffs + do + vec(n+1)=0.0 + nn=ishft(n+1,-1) + vec(1:nn)=vec(1:n:2)+pow*vec(2:n+1:2) + if (nn == 1) exit + pow=pow*pow + n=nn + end do + poly_rr=vec(1) + deallocate(vec) + end if + END FUNCTION poly_rr + !BL + FUNCTION poly_dd(x,coeffs) + DOUBLE PRECISION, INTENT(IN) :: x + DOUBLE PRECISION, DIMENSION(:), INTENT(IN) :: coeffs + DOUBLE PRECISION :: poly_dd + DOUBLE PRECISION :: pow + DOUBLE PRECISION, DIMENSION(:), ALLOCATABLE :: vec + INTEGER :: i,n,nn + n=size(coeffs) + if (n <= 0) then + poly_dd=0.0d0 + else if (n < NPAR_POLY) then + poly_dd=coeffs(n) + do i=n-1,1,-1 + poly_dd=x*poly_dd+coeffs(i) + end do + else + allocate(vec(n+1)) + pow=x + vec(1:n)=coeffs + do + vec(n+1)=0.0d0 + nn=ishft(n+1,-1) + vec(1:nn)=vec(1:n:2)+pow*vec(2:n+1:2) + if (nn == 1) exit + pow=pow*pow + n=nn + end do + poly_dd=vec(1) + deallocate(vec) + end if + END FUNCTION poly_dd + !BL + FUNCTION poly_rc(x,coeffs) + COMPLEX, INTENT(IN) :: x + REAL, DIMENSION(:), INTENT(IN) :: coeffs + COMPLEX :: poly_rc + COMPLEX :: pow + COMPLEX, DIMENSION(:), ALLOCATABLE :: vec + INTEGER :: i,n,nn + n=size(coeffs) + if (n <= 0) then + poly_rc=0.0 + else if (n < NPAR_POLY) then + poly_rc=coeffs(n) + do i=n-1,1,-1 + poly_rc=x*poly_rc+coeffs(i) + end do + else + allocate(vec(n+1)) + pow=x + vec(1:n)=coeffs + do + vec(n+1)=0.0 + nn=ishft(n+1,-1) + vec(1:nn)=vec(1:n:2)+pow*vec(2:n+1:2) + if (nn == 1) exit + pow=pow*pow + n=nn + end do + poly_rc=vec(1) + deallocate(vec) + end if + END FUNCTION poly_rc + !BL + FUNCTION poly_cc(x,coeffs) + COMPLEX, INTENT(IN) :: x + COMPLEX, DIMENSION(:), INTENT(IN) :: coeffs + COMPLEX :: poly_cc + COMPLEX :: pow + COMPLEX, DIMENSION(:), ALLOCATABLE :: vec + INTEGER :: i,n,nn + n=size(coeffs) + if (n <= 0) then + poly_cc=0.0 + else if (n < NPAR_POLY) then + poly_cc=coeffs(n) + do i=n-1,1,-1 + poly_cc=x*poly_cc+coeffs(i) + end do + else + allocate(vec(n+1)) + pow=x + vec(1:n)=coeffs + do + vec(n+1)=0.0 + nn=ishft(n+1,-1) + vec(1:nn)=vec(1:n:2)+pow*vec(2:n+1:2) + if (nn == 1) exit + pow=pow*pow + n=nn + end do + poly_cc=vec(1) + deallocate(vec) + end if + END FUNCTION poly_cc + !BL + FUNCTION poly_rrv(x,coeffs) + REAL, DIMENSION(:), INTENT(IN) :: coeffs,x + REAL, DIMENSION(size(x)) :: poly_rrv + INTEGER :: i,n,m + m=size(coeffs) + n=size(x) + if (m <= 0) then + poly_rrv=0.0 + else if (m < n .or. m < NPAR_POLY) then + poly_rrv=coeffs(m) + do i=m-1,1,-1 + poly_rrv=x*poly_rrv+coeffs(i) + end do + else + do i=1,n + poly_rrv(i)=poly_rr(x(i),coeffs) + end do + end if + END FUNCTION poly_rrv + !BL + FUNCTION poly_ddv(x,coeffs) + DOUBLE PRECISION, DIMENSION(:), INTENT(IN) :: coeffs,x + DOUBLE PRECISION, DIMENSION(size(x)) :: poly_ddv + INTEGER :: i,n,m + m=size(coeffs) + n=size(x) + if (m <= 0) then + poly_ddv=0.0d0 + else if (m < n .or. m < NPAR_POLY) then + poly_ddv=coeffs(m) + do i=m-1,1,-1 + poly_ddv=x*poly_ddv+coeffs(i) + end do + else + do i=1,n + poly_ddv(i)=poly_dd(x(i),coeffs) + end do + end if + END FUNCTION poly_ddv + !BL + FUNCTION poly_msk_rrv(x,coeffs,mask) + REAL, DIMENSION(:), INTENT(IN) :: coeffs,x + LOGICAL, DIMENSION(:), INTENT(IN) :: mask + REAL, DIMENSION(size(x)) :: poly_msk_rrv + poly_msk_rrv=unpack(poly_rrv(pack(x,mask),coeffs),mask,0.0) + END FUNCTION poly_msk_rrv + !BL + FUNCTION poly_msk_ddv(x,coeffs,mask) + DOUBLE PRECISION, DIMENSION(:), INTENT(IN) :: coeffs,x + LOGICAL, DIMENSION(:), INTENT(IN) :: mask + DOUBLE PRECISION, DIMENSION(size(x)) :: poly_msk_ddv + poly_msk_ddv=unpack(poly_ddv(pack(x,mask),coeffs),mask,0.0d0) + END FUNCTION poly_msk_ddv + +END MODULE poly_m Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/poly_term.f90 =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/poly_term.f90 (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/poly_term.f90 (revision 1795) @@ -0,0 +1,51 @@ +MODULE poly_term_m + + IMPLICIT NONE + + INTEGER, PARAMETER, private:: NPAR_POLYTERM=8 + + INTERFACE poly_term + MODULE PROCEDURE poly_term_rr, poly_term_cc + END INTERFACE + + private poly_term_rr,poly_term_cc + +CONTAINS + + RECURSIVE FUNCTION poly_term_rr(a,b) RESULT(u) + REAL, DIMENSION(:), INTENT(IN) :: a + REAL, INTENT(IN) :: b + REAL, DIMENSION(size(a)) :: u + INTEGER n,j + n=size(a) + if (n <= 0) RETURN + u(1)=a(1) + if (n < NPAR_POLYTERM) then + do j=2,n + u(j)=a(j)+b*u(j-1) + end do + else + u(2:n:2)=poly_term_rr(a(2:n:2)+a(1:n-1:2)*b,b*b) + u(3:n:2)=a(3:n:2)+b*u(2:n-1:2) + end if + END FUNCTION poly_term_rr + !BL + RECURSIVE FUNCTION poly_term_cc(a,b) RESULT(u) + COMPLEX, DIMENSION(:), INTENT(IN) :: a + COMPLEX, INTENT(IN) :: b + COMPLEX, DIMENSION(size(a)) :: u + INTEGER n,j + n=size(a) + if (n <= 0) RETURN + u(1)=a(1) + if (n < NPAR_POLYTERM) then + do j=2,n + u(j)=a(j)+b*u(j-1) + end do + else + u(2:n:2)=poly_term_cc(a(2:n:2)+a(1:n-1:2)*b,b*b) + u(3:n:2)=a(3:n:2)+b*u(2:n-1:2) + end if + END FUNCTION poly_term_cc + +END MODULE poly_term_m Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/put_diag.f90 =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/put_diag.f90 (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/put_diag.f90 (revision 1795) @@ -0,0 +1,34 @@ +MODULE put_diag_m + + IMPLICIT NONE + + INTERFACE put_diag + MODULE PROCEDURE put_diag_rv, put_diag_r + END INTERFACE + + private put_diag_rv, put_diag_r + +CONTAINS + + SUBROUTINE put_diag_rv(diagv,mat) + use assert_eq_m, only: assert_eq + REAL, DIMENSION(:), INTENT(IN) :: diagv + REAL, DIMENSION(:,:), INTENT(INOUT) :: mat + INTEGER :: j,n + n=assert_eq(size(diagv),min(size(mat,1),size(mat,2)),'put_diag_rv') + do j=1,n + mat(j,j)=diagv(j) + end do + END SUBROUTINE put_diag_rv + !BL + SUBROUTINE put_diag_r(scal,mat) + REAL, INTENT(IN) :: scal + REAL, DIMENSION(:,:), INTENT(INOUT) :: mat + INTEGER :: j,n + n = min(size(mat,1),size(mat,2)) + do j=1,n + mat(j,j)=scal + end do + END SUBROUTINE put_diag_r + +END MODULE put_diag_m Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/reallocate.f90 =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/reallocate.f90 (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/reallocate.f90 (revision 1795) @@ -0,0 +1,87 @@ +MODULE reallocate_m + + IMPLICIT NONE + + INTERFACE reallocate + MODULE PROCEDURE reallocate_rv, reallocate_rm, reallocate_iv, & + reallocate_im, reallocate_hv + END INTERFACE + + private + public reallocate + +CONTAINS + + FUNCTION reallocate_rv(p, n) + REAL, POINTER:: p(:), reallocate_rv(:) + INTEGER, INTENT(IN):: n + INTEGER nold + !----------------------------------------------- + allocate(reallocate_rv(n)) + if (.not. associated(p)) RETURN + nold = size(p) + reallocate_rv(:min(nold, n)) = p(:min(nold, n)) + deallocate(p) + END FUNCTION reallocate_rv + + !*********************************************************** + + FUNCTION reallocate_iv(p, n) + INTEGER, POINTER:: p(:), reallocate_iv(:) + INTEGER, INTENT(IN):: n + INTEGER nold + !----------------------------------------------- + allocate(reallocate_iv(n)) + if (.not. associated(p)) RETURN + nold = size(p) + reallocate_iv(:min(nold, n)) = p(:min(nold, n)) + deallocate(p) + END FUNCTION reallocate_iv + + !*********************************************************** + + FUNCTION reallocate_hv(p, n) + CHARACTER, POINTER:: p(:), reallocate_hv(:) + INTEGER, INTENT(IN):: n + INTEGER nold + !----------------------------------------------- + allocate(reallocate_hv(n)) + if (.not. associated(p)) RETURN + nold = size(p) + reallocate_hv(:min(nold, n)) = p(:min(nold, n)) + deallocate(p) + END FUNCTION reallocate_hv + + !*********************************************************** + + FUNCTION reallocate_rm(p, n, m) + REAL, DIMENSION(:, :), POINTER:: p, reallocate_rm + INTEGER, INTENT(IN):: n, m + INTEGER nold, mold + !----------------------------------------------- + allocate(reallocate_rm(n, m)) + if (.not. associated(p)) RETURN + nold = size(p, 1) + mold = size(p, 2) + reallocate_rm(:min(nold, n), :min(mold, m)) & + = p(:min(nold, n), :min(mold, m)) + deallocate(p) + END FUNCTION reallocate_rm + + !*********************************************************** + + FUNCTION reallocate_im(p, n, m) + INTEGER, DIMENSION(:, :), POINTER:: p, reallocate_im + INTEGER, INTENT(IN):: n, m + INTEGER nold, mold + !----------------------------------------------- + allocate(reallocate_im(n, m)) + if (.not. associated(p)) RETURN + nold = size(p, 1) + mold = size(p, 2) + reallocate_im(:min(nold, n), :min(mold, m)) & + = p(:min(nold, n), :min(mold, m)) + deallocate(p) + END FUNCTION reallocate_im + +END MODULE reallocate_m Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/scatter_add.f90 =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/scatter_add.f90 (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/scatter_add.f90 (revision 1795) @@ -0,0 +1,43 @@ +MODULE scatter_add_m + + ! Bug correction: intent of "dest". + + IMPLICIT NONE + + INTERFACE scatter_add + MODULE PROCEDURE scatter_add_r,scatter_add_d + END INTERFACE + + private scatter_add_r,scatter_add_d + +CONTAINS + + SUBROUTINE scatter_add_r(dest,source,dest_index) + use assert_eq_m, only: assert_eq + REAL, DIMENSION(:), INTENT(inOUT) :: dest + REAL, DIMENSION(:), INTENT(IN) :: source + INTEGER, DIMENSION(:), INTENT(IN) :: dest_index + INTEGER :: m,n,j,i + n=assert_eq(size(source),size(dest_index),'scatter_add_r') + m=size(dest) + do j=1,n + i=dest_index(j) + if (i > 0 .and. i <= m) dest(i)=dest(i)+source(j) + end do + END SUBROUTINE scatter_add_r + + SUBROUTINE scatter_add_d(dest,source,dest_index) + use assert_eq_m, only: assert_eq + DOUBLE PRECISION, DIMENSION(:), INTENT(inOUT) :: dest + DOUBLE PRECISION, DIMENSION(:), INTENT(IN) :: source + INTEGER, DIMENSION(:), INTENT(IN) :: dest_index + INTEGER :: m,n,j,i + n=assert_eq(size(source),size(dest_index),'scatter_add_d') + m=size(dest) + do j=1,n + i=dest_index(j) + if (i > 0 .and. i <= m) dest(i)=dest(i)+source(j) + end do + END SUBROUTINE scatter_add_d + +END MODULE scatter_add_m Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/scatter_max.f90 =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/scatter_max.f90 (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/scatter_max.f90 (revision 1795) @@ -0,0 +1,49 @@ +MODULE scatter_max_m + + ! Bug correction: intent of "dest". + + IMPLICIT NONE + + INTERFACE scatter_max + MODULE PROCEDURE scatter_max_r, scatter_max_d + END INTERFACE + + private scatter_max_r, scatter_max_d + +CONTAINS + + SUBROUTINE scatter_max_r(dest, source, dest_index) + use assert_eq_m, only: assert_eq + REAL, DIMENSION(:), INTENT(inOUT) :: dest + REAL, DIMENSION(:), INTENT(IN) :: source + INTEGER, DIMENSION(:), INTENT(IN) :: dest_index + + INTEGER :: m, n, j, i + !-------------------------------------- + n=assert_eq(size(source), size(dest_index), 'scatter_max_r') + m=size(dest) + do j=1, n + i=dest_index(j) + if (i > 0 .and. i <= m) dest(i)=max(dest(i), source(j)) + end do + END SUBROUTINE scatter_max_r + + !******************************************************* + + SUBROUTINE scatter_max_d(dest, source, dest_index) + use assert_eq_m, only: assert_eq + DOUBLE PRECISION, DIMENSION(:), INTENT(inOUT) :: dest + DOUBLE PRECISION, DIMENSION(:), INTENT(IN) :: source + INTEGER, DIMENSION(:), INTENT(IN) :: dest_index + + INTEGER :: m, n, j, i + !-------------------------------------- + n=assert_eq(size(source), size(dest_index), 'scatter_max_d') + m=size(dest) + do j=1, n + i=dest_index(j) + if (i > 0 .and. i <= m) dest(i)=max(dest(i), source(j)) + end do + END SUBROUTINE scatter_max_d + +END MODULE scatter_max_m Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/swap.f90 =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/swap.f90 (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/swap.f90 (revision 1795) @@ -0,0 +1,143 @@ +MODULE swap_m + + IMPLICIT NONE + + INTERFACE swap + MODULE PROCEDURE swap_i,swap_r,swap_rv,swap_c, & + swap_cv,swap_cm,swap_z,swap_zv,swap_zm, & + masked_swap_rs,masked_swap_rv,masked_swap_rm + END INTERFACE + + private + public swap + +CONTAINS + + SUBROUTINE swap_i(a,b) + INTEGER, INTENT(INOUT) :: a,b + INTEGER :: dum + dum=a + a=b + b=dum + END SUBROUTINE swap_i + + !************************************************ + + SUBROUTINE swap_r(a,b) + REAL, INTENT(INOUT) :: a,b + REAL :: dum + dum=a + a=b + b=dum + END SUBROUTINE swap_r + + !************************************************ + + SUBROUTINE swap_rv(a,b) + REAL, DIMENSION(:), INTENT(INOUT) :: a,b + REAL, DIMENSION(SIZE(a)) :: dum + dum=a + a=b + b=dum + END SUBROUTINE swap_rv + + !************************************************ + + SUBROUTINE swap_c(a,b) + COMPLEX, INTENT(INOUT) :: a,b + COMPLEX :: dum + dum=a + a=b + b=dum + END SUBROUTINE swap_c + + !************************************************ + + SUBROUTINE swap_cv(a,b) + COMPLEX, DIMENSION(:), INTENT(INOUT) :: a,b + COMPLEX, DIMENSION(SIZE(a)) :: dum + dum=a + a=b + b=dum + END SUBROUTINE swap_cv + + !************************************************ + + SUBROUTINE swap_cm(a,b) + COMPLEX, DIMENSION(:,:), INTENT(INOUT) :: a,b + COMPLEX, DIMENSION(size(a,1),size(a,2)) :: dum + dum=a + a=b + b=dum + END SUBROUTINE swap_cm + + !************************************************ + + SUBROUTINE swap_z(a,b) + COMPLEX(KIND(0D0)), INTENT(INOUT) :: a,b + COMPLEX(KIND(0D0)) :: dum + dum=a + a=b + b=dum + END SUBROUTINE swap_z + + !************************************************ + + SUBROUTINE swap_zv(a,b) + COMPLEX(KIND(0D0)), DIMENSION(:), INTENT(INOUT) :: a,b + COMPLEX(KIND(0D0)), DIMENSION(SIZE(a)) :: dum + dum=a + a=b + b=dum + END SUBROUTINE swap_zv + + !************************************************ + + SUBROUTINE swap_zm(a,b) + COMPLEX(KIND(0D0)), DIMENSION(:,:), INTENT(INOUT) :: a,b + COMPLEX(KIND(0D0)), DIMENSION(size(a,1),size(a,2)) :: dum + dum=a + a=b + b=dum + END SUBROUTINE swap_zm + + !************************************************ + + SUBROUTINE masked_swap_rs(a,b,mask) + REAL, INTENT(INOUT) :: a,b + LOGICAL, INTENT(IN) :: mask + REAL :: swp + if (mask) then + swp=a + a=b + b=swp + end if + END SUBROUTINE masked_swap_rs + + !************************************************ + + SUBROUTINE masked_swap_rv(a,b,mask) + REAL, DIMENSION(:), INTENT(INOUT) :: a,b + LOGICAL, DIMENSION(:), INTENT(IN) :: mask + REAL, DIMENSION(size(a)) :: swp + where (mask) + swp=a + a=b + b=swp + end where + END SUBROUTINE masked_swap_rv + + !************************************************ + + SUBROUTINE masked_swap_rm(a,b,mask) + REAL, DIMENSION(:,:), INTENT(INOUT) :: a,b + LOGICAL, DIMENSION(:,:), INTENT(IN) :: mask + REAL, DIMENSION(size(a,1),size(a,2)) :: swp + where (mask) + swp=a + a=b + b=swp + end where + END SUBROUTINE masked_swap_rm + +END MODULE swap_m Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/unit_matrix.f90 =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/unit_matrix.f90 (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/unit_matrix.f90 (revision 1795) @@ -0,0 +1,18 @@ +module unit_matrix_m + + implicit none + +contains + + SUBROUTINE unit_matrix(mat) + USE nrtype, only: wp + REAL(WP), DIMENSION(:,:), INTENT(OUT) :: mat + INTEGER i,n + n=min(size(mat,1),size(mat,2)) + mat(:,:)=0.0_wp + do i=1,n + mat(i,i)=1.0_wp + end do + END SUBROUTINE unit_matrix + +end module unit_matrix_m Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/upper_triangle.f90 =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/upper_triangle.f90 (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/upper_triangle.f90 (revision 1795) @@ -0,0 +1,19 @@ +module upper_triangle_m + + implicit none + +contains + + FUNCTION upper_triangle(j,k,extra) + use arth_m, only: arth + use outerdiff_m, only: outerdiff + INTEGER, INTENT(IN) :: j,k + INTEGER, OPTIONAL, INTENT(IN) :: extra + LOGICAL, DIMENSION(j,k) :: upper_triangle + INTEGER :: n + n=0 + if (present(extra)) n=extra + upper_triangle=(outerdiff(arth(1,1,j),arth(1,1,k)) < n) + END FUNCTION upper_triangle + +end module upper_triangle_m Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/vabs.f90 =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/vabs.f90 (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/vabs.f90 (revision 1795) @@ -0,0 +1,14 @@ +module vabs_m + + implicit none + +contains + + FUNCTION vabs(v) + USE nrtype, only: wp + REAL(WP), DIMENSION(:), INTENT(IN) :: v + REAL(WP) :: vabs + vabs=sqrt(dot_product(v,v)) + END FUNCTION vabs + +end module vabs_m Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/zroots_unity.f90 =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/zroots_unity.f90 (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NR_util/zroots_unity.f90 (revision 1795) @@ -0,0 +1,35 @@ +module zroots_unity_m + + implicit none + +contains + + FUNCTION zroots_unity(n, nn) + + ! Returns nn powers of the nth root of unity. + + USE nrtype, only: wp, twopi + + INTEGER, INTENT(IN) :: n, nn + COMPLEX(WP) zroots_unity(nn) + + ! Local: + INTEGER k + REAL(WP) :: theta + + !------------------------------------------------- + + zroots_unity(1) = 1. + theta = TWOPI / n + k=1 + do + if (k >= nn) exit + zroots_unity(k + 1) = cmplx(cos(k * theta), sin(k * theta), WP) + zroots_unity(k + 2: min(2 * k, nn)) = zroots_unity(k + 1) & + * zroots_unity(2: min(k, nn - k)) + k = 2 * k + end do + + END FUNCTION zroots_unity + +end module zroots_unity_m Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/GNUmakefile =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/GNUmakefile (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/GNUmakefile (revision 1795) @@ -0,0 +1,35 @@ +# This is a makefile for GNU make. +# This makefile builds "NetCDF95". + +# 1. Source files + +sources = nf95_inq_dimid.f90 nf95_inquire_dimension.f90 nf95_inq_varid.f90 nf95_inquire_variable.f90 nf95_create.f90 nf95_def_dim.f90 nf95_redef.f90 nf95_enddef.f90 nf95_close.f90 nf95_copy_att.f90 nf95_inquire_attribute.f90 nf95_inquire.f90 netcdf95.f90 nf95_def_var.f90 nf95_gw_var.f90 nf95_put_var.f90 nf95_put_att.f90 handle_err.f90 nf95_get_att.f90 nf95_get_var.f90 find_coord.f90 nf95_open.f90 + +# 2. Objects and library + +objects := $(sources:.f90=.o) +lib = libnetcdf95.a + +# 3. Compiler-dependent part + +override FFLAGS += $(if ${NETCDF_INC_DIR}, -I${NETCDF_INC_DIR}) + +# 4. Rules + +# Extend known suffixes: +%.o: %.f90 + $(COMPILE.f) $(OUTPUT_OPTION) $< + +.PHONY: all clean depend + +all: ${lib} +${lib}: ${lib}(${objects}) + +depend depend.mk: + makedepf90 -Wmissing -Wconfused -nosrc -u netcdf -u typesizes ${sources} >depend.mk + +clean: + rm -f ${lib} ${objects} + +# Dependencies between object files and include files: +include depend.mk Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/Read_me =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/Read_me (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/Read_me (revision 1795) @@ -0,0 +1,2 @@ +See: +http://www.lmd.jussieu.fr/~lglmd/NetCDF95 Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/depend.mk =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/depend.mk (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/depend.mk (revision 1795) @@ -0,0 +1,21 @@ +find_coord.o : nf95_inquire_variable.o nf95_inquire.o nf95_inquire_dimension.o nf95_inq_varid.o nf95_get_att.o +netcdf95.o : nf95_redef.o nf95_put_var.o nf95_put_att.o nf95_open.o nf95_inquire_variable.o nf95_inquire.o nf95_inquire_dimension.o nf95_inquire_attribute.o nf95_inq_varid.o nf95_inq_dimid.o nf95_gw_var.o nf95_get_var.o nf95_get_att.o nf95_enddef.o nf95_def_var.o nf95_def_dim.o nf95_create.o nf95_copy_att.o nf95_close.o handle_err.o find_coord.o +nf95_close.o : handle_err.o +nf95_copy_att.o : handle_err.o +nf95_create.o : handle_err.o +nf95_def_dim.o : handle_err.o +nf95_def_var.o : handle_err.o +nf95_enddef.o : handle_err.o +nf95_get_att.o : nf95_inquire_attribute.o handle_err.o +nf95_get_var.o : handle_err.o +nf95_gw_var.o : nf95_inquire_dimension.o nf95_inquire_variable.o nf95_get_var.o +nf95_inq_dimid.o : handle_err.o +nf95_inq_varid.o : handle_err.o +nf95_inquire.o : handle_err.o +nf95_inquire_attribute.o : handle_err.o +nf95_inquire_dimension.o : handle_err.o +nf95_inquire_variable.o : handle_err.o +nf95_open.o : handle_err.o +nf95_put_att.o : handle_err.o +nf95_put_var.o : handle_err.o +nf95_redef.o : handle_err.o Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/find_coord.f90 =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/find_coord.f90 (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/find_coord.f90 (revision 1795) @@ -0,0 +1,103 @@ +module find_coord_m + + implicit none + +contains + + subroutine find_coord(ncid, name, dimid, varid, std_name) + + ! This procedure returns the name, dimension id or variable id of + ! the NetCDF coordinate with standard name "std_name", if such a + ! coordinate exists. The standard name is only used to know what + ! to search, it is not used for the search itself. The search + ! itself is done via a string match on the attribute "units". + + use netcdf, only: NF90_MAX_NAME, NF90_NOERR + use nf95_get_att_m, only: nf95_get_att + use nf95_inq_varid_m, only: nf95_inq_varid + use nf95_inquire_dimension_m, only: nf95_inquire_dimension + use nf95_inquire_m, only: nf95_inquire + use nf95_inquire_variable_m, only: nf95_inquire_variable + + integer, intent(in):: ncid + + character(len=*), intent(out), optional:: name ! blanks if not found + ! The actual character argument should normally have the length + ! "NF90_MAX_NAME". + + integer, intent(out), optional:: dimid ! 0 if not found + integer, intent(out), optional:: varid ! 0 if not found + + character(len=*), intent(in):: std_name + ! standard name : "latitude", "longitude" or "time" + + ! Variables local to the procedure: + + character(len=13) units + logical exact ! "units" must be matched exactly + + integer ncerr, nDimensions, dimid_local, varid_local + character(len=NF90_MAX_NAME) name_local + integer, pointer:: dimids(:) + character(len=80) values + logical found + + !---------------------------------------------- + + select case (std_name) + case("longitude") + units="degrees_east" + exact=.true. + case("latitude") + units="degrees_north" + exact=.true. + case("time") + units=" since" + exact=.false. + case default + print *, "find_coord: bad value of std_name" + print *, "std_name = ", std_name + stop 1 + end select + + call nf95_inquire(ncid, nDimensions) + dimid_local = 0 + found = .false. + + ! Loop on dimensions: + do while (.not. found .and. dimid_local < nDimensions) + dimid_local = dimid_local + 1 + call nf95_inquire_dimension(ncid, dimid_local, name_local) + call nf95_inq_varid(ncid, name_local, varid_local, ncerr) + if (ncerr == NF90_NOERR) then + call nf95_inquire_variable(ncid, varid_local, dimids=dimids) + if (size(dimids) == 1) then + if (dimids(1) == dimid_local) then + ! We have found a coordinate + call nf95_get_att(ncid, varid_local, "units", values, ncerr) + if (ncerr == NF90_NOERR)then + if (exact) then + found = values == units + else + found = index(values, trim(units)) /= 0 + end if + end if + end if + end if + deallocate(dimids) ! pointer + end if + end do + + if (found) then + if (present(name)) name = name_local + if (present(dimid)) dimid = dimid_local + if (present(varid)) varid = varid_local + else + if (present(name)) name = "" + if (present(dimid)) dimid = 0 + if (present(varid)) varid = 0 + end if + + end subroutine find_coord + +end module find_coord_m Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/handle_err.f90 =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/handle_err.f90 (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/handle_err.f90 (revision 1795) @@ -0,0 +1,45 @@ +module handle_err_m + + implicit none + +contains + + subroutine handle_err(message, ncerr, ncid, varid) + + use netcdf, only: nf90_strerror, nf90_noerr, nf90_close + + character(len=*), intent(in):: message + ! (should include name of calling procedure) + + integer, intent(in):: ncerr + + integer, intent(in), optional :: ncid + ! (Provide this argument if you want "handle_err" to try to close + ! the file.) + + integer, intent(in), optional :: varid + + ! Variable local to the procedure: + integer ncerr_close + + !------------------- + + if (ncerr /= nf90_noerr) then + print *, message, ":" + if (present(varid)) print *, "varid = ", varid + print *, trim(nf90_strerror(ncerr)) + if (present(ncid)) then + ! Try to close, to leave the file in a consistent state: + ncerr_close = nf90_close(ncid) + ! (do not call "nf95_close", we do not want to recurse) + if (ncerr_close /= nf90_noerr) then + print *, "nf90_close:" + print *, trim(nf90_strerror(ncerr_close)) + end if + end if + stop 1 + end if + + end subroutine handle_err + +end module handle_err_m Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/netcdf95.f90 =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/netcdf95.f90 (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/netcdf95.f90 (revision 1795) @@ -0,0 +1,27 @@ +module netcdf95 + + ! Author: Lionel GUEZ + + use find_coord_m + use handle_err_m + use nf95_close_m + use nf95_copy_att_m + use nf95_create_m + use nf95_def_dim_m + use nf95_def_var_m + use nf95_enddef_m + use nf95_get_att_m + use nf95_get_var_m + use nf95_gw_var_m + use nf95_inq_dimid_m + use nf95_inq_varid_m + use nf95_inquire_attribute_m + use nf95_inquire_dimension_m + use nf95_inquire_m + use nf95_inquire_variable_m + use nf95_open_m + use nf95_put_att_m + use nf95_put_var_m + use nf95_redef_m + +end module netcdf95 Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/nf95_close.f90 =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/nf95_close.f90 (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/nf95_close.f90 (revision 1795) @@ -0,0 +1,29 @@ +module nf95_close_m + + implicit none + +contains + + subroutine nf95_close(ncid, ncerr) + + use handle_err_m, only: handle_err + use netcdf, only: nf90_close + + integer, intent( in) :: ncid + integer, intent(out), optional :: ncerr + + ! Variable local to the procedure: + integer ncerr_not_opt + + !------------------- + + ncerr_not_opt = nf90_close(ncid) + if (present(ncerr)) then + ncerr = ncerr_not_opt + else + call handle_err("nf95_close", ncerr_not_opt) + end if + + end subroutine nf95_close + +end module nf95_close_m Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/nf95_copy_att.f90 =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/nf95_copy_att.f90 (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/nf95_copy_att.f90 (revision 1795) @@ -0,0 +1,32 @@ +module nf95_copy_att_m + + implicit none + +contains + + + subroutine nf95_copy_att(ncid_in, varid_in, name, ncid_out, varid_out, ncerr) + + use handle_err_m, only: handle_err + use netcdf, only: nf90_copy_att + + integer, intent( in):: ncid_in, varid_in + character(len=*), intent( in):: name + integer, intent( in):: ncid_out, varid_out + integer, intent(out), optional:: ncerr + + ! Variable local to the procedure: + integer ncerr_not_opt + + !------------------- + + ncerr_not_opt = nf90_copy_att(ncid_in, varid_in, name, ncid_out, varid_out) + if (present(ncerr)) then + ncerr = ncerr_not_opt + else + call handle_err("nf95_copy_att " // name, ncerr_not_opt, ncid_out) + end if + + end subroutine nf95_copy_att + +end module nf95_copy_att_m Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/nf95_create.f90 =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/nf95_create.f90 (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/nf95_create.f90 (revision 1795) @@ -0,0 +1,33 @@ +module nf95_create_m + + implicit none + +contains + + subroutine nf95_create(path, cmode, ncid, initialsize, chunksize, ncerr) + + use handle_err_m, only: handle_err + use netcdf, only: nf90_create + + character (len = *), intent(in ) :: path + integer, intent(in ) :: cmode + integer, intent( out) :: ncid + integer, optional, intent(in ) :: initialsize + integer, optional, intent(inout) :: chunksize + integer, intent(out), optional :: ncerr + + ! Variable local to the procedure: + integer ncerr_not_opt + + !------------------- + + ncerr_not_opt = nf90_create(path, cmode, ncid, initialsize, chunksize) + if (present(ncerr)) then + ncerr = ncerr_not_opt + else + call handle_err("nf95_create " // path, ncerr_not_opt) + end if + + end subroutine nf95_create + +end module nf95_create_m Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/nf95_def_dim.f90 =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/nf95_def_dim.f90 (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/nf95_def_dim.f90 (revision 1795) @@ -0,0 +1,32 @@ +module nf95_def_dim_m + + implicit none + +contains + + subroutine nf95_def_dim(ncid, name, nclen, dimid, ncerr) + + use handle_err_m, only: handle_err + use netcdf, only: nf90_def_dim + + integer, intent( in) :: ncid + character (len = *), intent( in) :: name + integer, intent( in) :: nclen + integer, intent(out) :: dimid + integer, intent(out), optional :: ncerr + + ! Variable local to the procedure: + integer ncerr_not_opt + + !------------------- + + ncerr_not_opt = nf90_def_dim(ncid, name, nclen, dimid) + if (present(ncerr)) then + ncerr = ncerr_not_opt + else + call handle_err("nf95_def_dim " // name, ncerr_not_opt, ncid) + end if + + end subroutine nf95_def_dim + +end module nf95_def_dim_m Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/nf95_def_var.f90 =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/nf95_def_var.f90 (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/nf95_def_var.f90 (revision 1795) @@ -0,0 +1,101 @@ +module nf95_def_var_m + + ! The generic procedure name "nf90_def_var" applies to + ! "nf90_def_var_Scalar" but we cannot apply the generic procedure name + ! "nf95_def_var" to "nf95_def_var_scalar" because of the additional + ! optional argument. + ! "nf95_def_var_scalar" cannot be distinguished from "nf95_def_var_oneDim". + + implicit none + + interface nf95_def_var + module procedure nf95_def_var_oneDim, nf95_def_var_ManyDims + end interface + + private + public nf95_def_var, nf95_def_var_scalar + +contains + + subroutine nf95_def_var_scalar(ncid, name, xtype, varid, ncerr) + + use netcdf, only: nf90_def_var + use handle_err_m, only: handle_err + + integer, intent( in) :: ncid + character (len = *), intent( in) :: name + integer, intent( in) :: xtype + integer, intent(out) :: varid + integer, intent(out), optional:: ncerr + + ! Variable local to the procedure: + integer ncerr_not_opt + + !------------------- + + ncerr_not_opt = nf90_def_var(ncid, name, xtype, varid) + if (present(ncerr)) then + ncerr = ncerr_not_opt + else + call handle_err("nf95_def_var_scalar " // name, ncerr_not_opt, ncid) + end if + + end subroutine nf95_def_var_scalar + + !*********************** + + subroutine nf95_def_var_oneDim(ncid, name, xtype, dimids, varid, ncerr) + + use netcdf, only: nf90_def_var + use handle_err_m, only: handle_err + + integer, intent( in) :: ncid + character (len = *), intent( in) :: name + integer, intent( in) :: xtype + integer, intent( in) :: dimids + integer, intent(out) :: varid + integer, intent(out), optional:: ncerr + + ! Variable local to the procedure: + integer ncerr_not_opt + + !------------------- + + ncerr_not_opt = nf90_def_var(ncid, name, xtype, dimids, varid) + if (present(ncerr)) then + ncerr = ncerr_not_opt + else + call handle_err("nf95_def_var_oneDim " // name, ncerr_not_opt, ncid) + end if + + end subroutine nf95_def_var_oneDim + + !*********************** + + subroutine nf95_def_var_ManyDims(ncid, name, xtype, dimids, varid, ncerr) + + use netcdf, only: nf90_def_var + use handle_err_m, only: handle_err + + integer, intent( in) :: ncid + character (len = *), intent( in) :: name + integer, intent( in) :: xtype + integer, dimension(:), intent( in) :: dimids + integer, intent(out) :: varid + integer, intent(out), optional:: ncerr + + ! Variable local to the procedure: + integer ncerr_not_opt + + !------------------- + + ncerr_not_opt = nf90_def_var(ncid, name, xtype, dimids, varid) + if (present(ncerr)) then + ncerr = ncerr_not_opt + else + call handle_err("nf95_def_var_ManyDims " // name, ncerr_not_opt, ncid) + end if + + end subroutine nf95_def_var_ManyDims + +end module nf95_def_var_m Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/nf95_enddef.f90 =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/nf95_enddef.f90 (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/nf95_enddef.f90 (revision 1795) @@ -0,0 +1,30 @@ +module nf95_enddef_m + + implicit none + +contains + + subroutine nf95_enddef(ncid, h_minfree, v_align, v_minfree, r_align, ncerr) + + use handle_err_m, only: handle_err + use netcdf, only: nf90_enddef + + integer, intent( in) :: ncid + integer, optional, intent( in) :: h_minfree, v_align, v_minfree, r_align + integer, intent(out), optional :: ncerr + + ! Variable local to the procedure: + integer ncerr_not_opt + + !------------------- + + ncerr_not_opt = nf90_enddef(ncid, h_minfree, v_align, v_minfree, r_align) + if (present(ncerr)) then + ncerr = ncerr_not_opt + else + call handle_err("nf95_enddef", ncerr_not_opt, ncid) + end if + + end subroutine nf95_enddef + +end module nf95_enddef_m Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/nf95_get_att.f90 =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/nf95_get_att.f90 (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/nf95_get_att.f90 (revision 1795) @@ -0,0 +1,106 @@ +module nf95_get_att_m + + use handle_err_m, only: handle_err + use netcdf, only: nf90_get_att, nf90_noerr + use nf95_inquire_attribute_m, only: nf95_inquire_attribute + + implicit none + + interface nf95_get_att + module procedure nf95_get_att_text, nf95_get_att_one_FourByteInt + + ! The difference between the specific procedures is the type of + ! argument "values". + end interface + + private + public nf95_get_att + +contains + + subroutine nf95_get_att_text(ncid, varid, name, values, ncerr) + + integer, intent( in) :: ncid, varid + character(len = *), intent( in) :: name + character(len = *), intent(out) :: values + integer, intent(out), optional:: ncerr + + ! Variables local to the procedure: + integer ncerr_not_opt + integer att_len + + !------------------- + + ! Check that the length of "values" is large enough: + call nf95_inquire_attribute(ncid, varid, name, nclen=att_len, & + ncerr=ncerr_not_opt) + if (ncerr_not_opt == nf90_noerr) then + if (len(values) < att_len) then + print *, "nf95_get_att_text" + print *, "varid = ", varid + print *, "attribute name: ", name + print *, 'length of "values" is not large enough' + print *, "len(values) = ", len(values) + print *, "number of characters in attribute: ", att_len + stop 1 + end if + end if + + values = "" ! useless in NetCDF version 3.6.2 or better + ncerr_not_opt = nf90_get_att(ncid, varid, name, values) + if (present(ncerr)) then + ncerr = ncerr_not_opt + else + call handle_err("nf95_get_att_text " // trim(name), ncerr_not_opt, & + ncid, varid) + end if + + if (att_len >= 1 .and. ncerr_not_opt == nf90_noerr) then + ! Remove null terminator, if any: + if (iachar(values(att_len:att_len)) == 0) values(att_len:att_len) = " " + end if + + end subroutine nf95_get_att_text + + !*********************** + + subroutine nf95_get_att_one_FourByteInt(ncid, varid, name, values, ncerr) + + use typesizes, only: FourByteInt + + integer, intent( in) :: ncid, varid + character(len = *), intent( in) :: name + integer (kind = FourByteInt), intent(out) :: values + integer, intent(out), optional:: ncerr + + ! Variables local to the procedure: + integer ncerr_not_opt + integer att_len + + !------------------- + + ! Check that the attribute contains a single value: + call nf95_inquire_attribute(ncid, varid, name, nclen=att_len, & + ncerr=ncerr_not_opt) + if (ncerr_not_opt == nf90_noerr) then + if (att_len /= 1) then + print *, "nf95_get_att_one_FourByteInt" + print *, "varid = ", varid + print *, "attribute name: ", name + print *, 'the attribute does not contain a single value' + print *, "number of values in attribute: ", att_len + stop 1 + end if + end if + + ncerr_not_opt = nf90_get_att(ncid, varid, name, values) + if (present(ncerr)) then + ncerr = ncerr_not_opt + else + call handle_err("nf95_get_att_one_FourByteInt " // trim(name), & + ncerr_not_opt, ncid, varid) + end if + + end subroutine nf95_get_att_one_FourByteInt + +end module nf95_get_att_m Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/nf95_get_var.f90 =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/nf95_get_var.f90 (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/nf95_get_var.f90 (revision 1795) @@ -0,0 +1,409 @@ +module nf95_get_var_m + + use netcdf, only: nf90_get_var + use handle_err_m, only: handle_err + + implicit none + + interface nf95_get_var + module procedure nf95_get_var_FourByteReal, nf95_get_var_FourByteInt, & + nf95_get_var_1D_FourByteReal, nf95_get_var_1D_FourByteInt, & + nf95_get_var_1D_EightByteReal, nf95_get_var_2D_FourByteReal, & + nf95_get_var_2D_EightByteReal, nf95_get_var_3D_FourByteInt, & + nf95_get_var_3D_FourByteReal, nf95_get_var_3D_EightByteReal, & + nf95_get_var_4D_FourByteReal, nf95_get_var_4D_EightByteReal, & + nf95_get_var_5D_FourByteReal, nf95_get_var_5D_EightByteReal + end interface + + private + public nf95_get_var + +contains + + subroutine nf95_get_var_FourByteReal(ncid, varid, values, start, ncerr) + + use typesizes, only: FourByteReal + + integer, intent(in) :: ncid, varid + real(kind = FourByteReal), intent(out) :: values + integer, dimension(:), optional, intent(in) :: start + integer, intent(out), optional:: ncerr + + ! Variable local to the procedure: + integer ncerr_not_opt + + !------------------- + + ncerr_not_opt = nf90_get_var(ncid, varid, values, start) + if (present(ncerr)) then + ncerr = ncerr_not_opt + else + call handle_err("nf95_get_var_FourByteReal", ncerr_not_opt, ncid, & + varid) + end if + + end subroutine nf95_get_var_FourByteReal + + !*********************** + + subroutine nf95_get_var_FourByteInt(ncid, varid, values, start, ncerr) + + use typesizes, only: FourByteInt + + integer, intent(in) :: ncid, varid + integer(kind = FourByteInt), intent(out) :: values + integer, dimension(:), optional, intent(in) :: start + integer, intent(out), optional:: ncerr + + ! Variable local to the procedure: + integer ncerr_not_opt + + !------------------- + + ncerr_not_opt = nf90_get_var(ncid, varid, values, start) + if (present(ncerr)) then + ncerr = ncerr_not_opt + else + call handle_err("nf95_get_var_FourByteInt", ncerr_not_opt, ncid, & + varid) + end if + + end subroutine nf95_get_var_FourByteInt + + !*********************** + + subroutine nf95_get_var_1D_FourByteReal(ncid, varid, values, start, & + count_nc, stride, map, ncerr) + + use typesizes, only: FourByteReal + + integer, intent(in) :: ncid, varid + real(kind = FourByteReal), intent(out) :: values(:) + integer, dimension(:), optional, intent(in) :: start, count_nc, stride, map + integer, intent(out), optional:: ncerr + + ! Variable local to the procedure: + integer ncerr_not_opt + + !------------------- + + ncerr_not_opt = nf90_get_var(ncid, varid, values, start, count_nc, & + stride, map) + if (present(ncerr)) then + ncerr = ncerr_not_opt + else + call handle_err("nf95_get_var_1D_FourByteReal", ncerr_not_opt, ncid, & + varid) + end if + + end subroutine nf95_get_var_1D_FourByteReal + + !*********************** + + subroutine nf95_get_var_1D_FourByteInt(ncid, varid, values, start, & + count_nc, stride, map, ncerr) + + use typesizes, only: FourByteInt + + integer, intent(in) :: ncid, varid + integer(kind = FourByteInt), intent(out) :: values(:) + integer, dimension(:), optional, intent(in) :: start, count_nc, stride, map + integer, intent(out), optional:: ncerr + + ! Variable local to the procedure: + integer ncerr_not_opt + + !------------------- + + ncerr_not_opt = nf90_get_var(ncid, varid, values, start, count_nc, & + stride, map) + if (present(ncerr)) then + ncerr = ncerr_not_opt + else + call handle_err("nf95_get_var_1D_FourByteInt", ncerr_not_opt, ncid, & + varid) + end if + + end subroutine nf95_get_var_1D_FourByteInt + + !*********************** + + subroutine nf95_get_var_1D_EightByteReal(ncid, varid, values, start, & + count_nc, stride, map, ncerr) + + use typesizes, only: eightByteReal + + integer, intent(in) :: ncid, varid + real (kind = EightByteReal), intent(out) :: values(:) + integer, dimension(:), optional, intent(in):: start, count_nc, stride, map + integer, intent(out), optional:: ncerr + + ! Variable local to the procedure: + integer ncerr_not_opt + + !------------------- + + ncerr_not_opt = nf90_get_var(ncid, varid, values, start, count_nc, & + stride, map) + if (present(ncerr)) then + ncerr = ncerr_not_opt + else + call handle_err("nf95_get_var_1D_eightByteReal", ncerr_not_opt, ncid, & + varid) + end if + + end subroutine nf95_get_var_1D_EightByteReal + + !*********************** + + subroutine nf95_get_var_2D_FourByteReal(ncid, varid, values, start, & + count_nc, stride, map, ncerr) + + use typesizes, only: FourByteReal + + integer, intent(in) :: ncid, varid + real (kind = FourByteReal), intent(out) :: values(:, :) + integer, dimension(:), optional, intent(in) :: start, count_nc, stride, map + integer, intent(out), optional:: ncerr + + ! Variable local to the procedure: + integer ncerr_not_opt + + !------------------- + + ncerr_not_opt = nf90_get_var(ncid, varid, values, start, count_nc, & + stride, map) + if (present(ncerr)) then + ncerr = ncerr_not_opt + else + call handle_err("nf95_get_var_2D_FourByteReal", ncerr_not_opt, ncid, & + varid) + end if + + end subroutine nf95_get_var_2D_FourByteReal + + !*********************** + + subroutine nf95_get_var_2D_EightByteReal(ncid, varid, values, start, & + count_nc, stride, map, ncerr) + + use typesizes, only: EightByteReal + + integer, intent(in) :: ncid, varid + real (kind = EightByteReal), intent(out) :: values(:, :) + integer, dimension(:), optional, intent(in) :: start, count_nc, stride, map + integer, intent(out), optional:: ncerr + + ! Variable local to the procedure: + integer ncerr_not_opt + + !------------------- + + ncerr_not_opt = nf90_get_var(ncid, varid, values, start, count_nc, & + stride, map) + if (present(ncerr)) then + ncerr = ncerr_not_opt + else + call handle_err("nf95_get_var_2D_EightByteReal", ncerr_not_opt, ncid, & + varid) + end if + + end subroutine nf95_get_var_2D_EightByteReal + + !*********************** + + subroutine nf95_get_var_3D_FourByteInt(ncid, varid, values, start, & + count_nc, stride, map, ncerr) + + use typesizes, only: FourByteInt + + integer, intent(in):: ncid, varid + integer(kind = FourByteInt), intent(out):: values(:, :, :) + integer, dimension(:), optional, intent(in):: start, count_nc, stride, map + integer, intent(out), optional:: ncerr + + ! Variable local to the procedure: + integer ncerr_not_opt + + !------------------- + + ncerr_not_opt = nf90_get_var(ncid, varid, values, start, count_nc, & + stride, map) + if (present(ncerr)) then + ncerr = ncerr_not_opt + else + call handle_err("nf95_get_var_3D_FourByteInt", ncerr_not_opt, ncid, & + varid) + end if + + end subroutine nf95_get_var_3D_FourByteInt + + !*********************** + + subroutine nf95_get_var_3D_FourByteReal(ncid, varid, values, start, & + count_nc, stride, map, ncerr) + + use typesizes, only: FourByteReal + + integer, intent(in) :: ncid, varid + real (kind = FourByteReal), intent(out) :: values(:, :, :) + integer, dimension(:), optional, intent(in) :: start, count_nc, stride, map + integer, intent(out), optional:: ncerr + + ! Variable local to the procedure: + integer ncerr_not_opt + + !------------------- + + ncerr_not_opt = nf90_get_var(ncid, varid, values, start, count_nc, & + stride, map) + if (present(ncerr)) then + ncerr = ncerr_not_opt + else + call handle_err("nf95_get_var_3D_FourByteReal", ncerr_not_opt, ncid, & + varid) + end if + + end subroutine nf95_get_var_3D_FourByteReal + + !*********************** + + subroutine nf95_get_var_3D_EightByteReal(ncid, varid, values, start, & + count_nc, stride, map, ncerr) + + use typesizes, only: eightByteReal + + integer, intent(in) :: ncid, varid + real (kind = EightByteReal), intent(out) :: values(:, :, :) + integer, dimension(:), optional, intent(in) :: start, count_nc, stride, map + integer, intent(out), optional:: ncerr + + ! Variable local to the procedure: + integer ncerr_not_opt + + !------------------- + + ncerr_not_opt = nf90_get_var(ncid, varid, values, start, count_nc, & + stride, map) + if (present(ncerr)) then + ncerr = ncerr_not_opt + else + call handle_err("nf95_get_var_3D_eightByteReal", ncerr_not_opt, ncid, & + varid) + end if + + end subroutine nf95_get_var_3D_EightByteReal + + !*********************** + + subroutine nf95_get_var_4D_FourByteReal(ncid, varid, values, start, & + count_nc, stride, map, ncerr) + + use typesizes, only: FourByteReal + + integer, intent(in) :: ncid, varid + real (kind = FourByteReal), intent(out) :: values(:, :, :, :) + integer, dimension(:), optional, intent(in) :: start, count_nc, stride, map + integer, intent(out), optional:: ncerr + + ! Variable local to the procedure: + integer ncerr_not_opt + + !------------------- + + ncerr_not_opt = nf90_get_var(ncid, varid, values, start, count_nc, & + stride, map) + if (present(ncerr)) then + ncerr = ncerr_not_opt + else + call handle_err("nf95_get_var_4D_FourByteReal", ncerr_not_opt, ncid, & + varid) + end if + + end subroutine nf95_get_var_4D_FourByteReal + + !*********************** + + subroutine nf95_get_var_4D_EightByteReal(ncid, varid, values, start, & + count_nc, stride, map, ncerr) + + use typesizes, only: EightByteReal + + integer, intent(in):: ncid, varid + real(kind = EightByteReal), intent(out):: values(:, :, :, :) + integer, dimension(:), optional, intent(in):: start, count_nc, stride, map + integer, intent(out), optional:: ncerr + + ! Variable local to the procedure: + integer ncerr_not_opt + + !------------------- + + ncerr_not_opt = nf90_get_var(ncid, varid, values, start, count_nc, & + stride, map) + if (present(ncerr)) then + ncerr = ncerr_not_opt + else + call handle_err("nf95_get_var_4D_EightByteReal", ncerr_not_opt, ncid, & + varid) + end if + + end subroutine nf95_get_var_4D_EightByteReal + + !*********************** + + subroutine nf95_get_var_5D_FourByteReal(ncid, varid, values, start, & + count_nc, stride, map, ncerr) + + use typesizes, only: FourByteReal + + integer, intent(in):: ncid, varid + real(kind = FourByteReal), intent(out):: values(:, :, :, :, :) + integer, dimension(:), optional, intent(in):: start, count_nc, stride, map + integer, intent(out), optional:: ncerr + + ! Variable local to the procedure: + integer ncerr_not_opt + + !------------------- + + ncerr_not_opt = nf90_get_var(ncid, varid, values, start, count_nc, & + stride, map) + if (present(ncerr)) then + ncerr = ncerr_not_opt + else + call handle_err("nf95_get_var_5D_FourByteReal", ncerr_not_opt, ncid, & + varid) + end if + + end subroutine nf95_get_var_5D_FourByteReal + + !*********************** + + subroutine nf95_get_var_5D_EightByteReal(ncid, varid, values, start, & + count_nc, stride, map, ncerr) + + use typesizes, only: EightByteReal + + integer, intent(in):: ncid, varid + real(kind = EightByteReal), intent(out):: values(:, :, :, :, :) + integer, dimension(:), optional, intent(in):: start, count_nc, stride, map + integer, intent(out), optional:: ncerr + + ! Variable local to the procedure: + integer ncerr_not_opt + + !------------------- + + ncerr_not_opt = nf90_get_var(ncid, varid, values, start, count_nc, & + stride, map) + if (present(ncerr)) then + ncerr = ncerr_not_opt + else + call handle_err("nf95_get_var_5D_EightByteReal", ncerr_not_opt, ncid, & + varid) + end if + + end subroutine nf95_get_var_5D_EightByteReal + +end module nf95_get_var_m Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/nf95_gw_var.f90 =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/nf95_gw_var.f90 (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/nf95_gw_var.f90 (revision 1795) @@ -0,0 +1,403 @@ +module nf95_gw_var_m + + use nf95_get_var_m, only: NF95_GET_VAR + use nf95_inquire_variable_m, only: nf95_inquire_variable + use nf95_inquire_dimension_m, only: nf95_inquire_dimension + + implicit none + + interface nf95_gw_var + ! "nf95_gw_var" stands for "NetCDF 1995 get whole variable". + ! These procedures read a whole NetCDF variable (coordinate or + ! primary) into an array. + ! The difference between the procedures is the rank and type of + ! argument "values". + ! The procedures do not check the type of the NetCDF variable. + + module procedure nf95_gw_var_real_1d, nf95_gw_var_real_2d, & + nf95_gw_var_real_3d, nf95_gw_var_real_4d, nf95_gw_var_real_5d, & + nf95_gw_var_dble_1d, nf95_gw_var_dble_2d, nf95_gw_var_dble_3d, & + nf95_gw_var_dble_4d, nf95_gw_var_int_1d, nf95_gw_var_int_3d + end interface + + private + public nf95_gw_var + +contains + + subroutine nf95_gw_var_real_1d(ncid, varid, values) + + ! Real type, the array has rank 1. + + integer, intent(in):: ncid + integer, intent(in):: varid + real, pointer:: values(:) + + ! Variables local to the procedure: + integer nclen + integer, pointer:: dimids(:) + + !--------------------- + + call nf95_inquire_variable(ncid, varid, dimids=dimids) + + if (size(dimids) /= 1) then + print *, "nf95_gw_var_real_1d:" + print *, "varid = ", varid + print *, "rank of NetCDF variable is ", size(dimids), ", not 1" + stop 1 + end if + + call nf95_inquire_dimension(ncid, dimids(1), nclen=nclen) + deallocate(dimids) ! pointer + + allocate(values(nclen)) + if (nclen /= 0) call NF95_GET_VAR(ncid, varid, values) + + end subroutine nf95_gw_var_real_1d + + !************************************ + + subroutine nf95_gw_var_real_2d(ncid, varid, values) + + ! Real type, the array has rank 2. + + integer, intent(in):: ncid + integer, intent(in):: varid + real, pointer:: values(:, :) + + ! Variables local to the procedure: + integer nclen1, nclen2 + integer, pointer:: dimids(:) + + !--------------------- + + call nf95_inquire_variable(ncid, varid, dimids=dimids) + + if (size(dimids) /= 2) then + print *, "nf95_gw_var_real_2d:" + print *, "varid = ", varid + print *, "rank of NetCDF variable is ", size(dimids), ", not 2" + stop 1 + end if + + call nf95_inquire_dimension(ncid, dimids(1), nclen=nclen1) + call nf95_inquire_dimension(ncid, dimids(2), nclen=nclen2) + deallocate(dimids) ! pointer + + allocate(values(nclen1, nclen2)) + if (nclen1 /= 0 .and. nclen2 /= 0) call NF95_GET_VAR(ncid, varid, values) + + end subroutine nf95_gw_var_real_2d + + !************************************ + + subroutine nf95_gw_var_real_3d(ncid, varid, values) + + ! Real type, the array has rank 3. + + integer, intent(in):: ncid + integer, intent(in):: varid + real, pointer:: values(:, :, :) + + ! Variables local to the procedure: + integer nclen1, nclen2, nclen3 + integer, pointer:: dimids(:) + + !--------------------- + + call nf95_inquire_variable(ncid, varid, dimids=dimids) + + if (size(dimids) /= 3) then + print *, "nf95_gw_var_real_3d:" + print *, "varid = ", varid + print *, "rank of NetCDF variable is ", size(dimids), ", not 3" + stop 1 + end if + + call nf95_inquire_dimension(ncid, dimids(1), nclen=nclen1) + call nf95_inquire_dimension(ncid, dimids(2), nclen=nclen2) + call nf95_inquire_dimension(ncid, dimids(3), nclen=nclen3) + deallocate(dimids) ! pointer + + allocate(values(nclen1, nclen2, nclen3)) + if (nclen1 * nclen2 * nclen3 /= 0) call NF95_GET_VAR(ncid, varid, values) + + end subroutine nf95_gw_var_real_3d + + !************************************ + + subroutine nf95_gw_var_real_4d(ncid, varid, values) + + ! Real type, the array has rank 4. + + integer, intent(in):: ncid + integer, intent(in):: varid + real, pointer:: values(:, :, :, :) + + ! Variables local to the procedure: + integer len_dim(4), i + integer, pointer:: dimids(:) + + !--------------------- + + call nf95_inquire_variable(ncid, varid, dimids=dimids) + + if (size(dimids) /= 4) then + print *, "nf95_gw_var_real_4d:" + print *, "varid = ", varid + print *, "rank of NetCDF variable is ", size(dimids), ", not 4" + stop 1 + end if + + do i = 1, 4 + call nf95_inquire_dimension(ncid, dimids(i), nclen=len_dim(i)) + end do + deallocate(dimids) ! pointer + + allocate(values(len_dim(1), len_dim(2), len_dim(3), len_dim(4))) + if (all(len_dim /= 0)) call NF95_GET_VAR(ncid, varid, values) + + end subroutine nf95_gw_var_real_4d + + !************************************ + + subroutine nf95_gw_var_real_5d(ncid, varid, values) + + ! Real type, the array has rank 5. + + integer, intent(in):: ncid + integer, intent(in):: varid + real, pointer:: values(:, :, :, :, :) + + ! Variables local to the procedure: + integer len_dim(5), i + integer, pointer:: dimids(:) + + !--------------------- + + call nf95_inquire_variable(ncid, varid, dimids=dimids) + + if (size(dimids) /= 5) then + print *, "nf95_gw_var_real_5d:" + print *, "varid = ", varid + print *, "rank of NetCDF variable is ", size(dimids), ", not 5" + stop 1 + end if + + do i = 1, 5 + call nf95_inquire_dimension(ncid, dimids(i), nclen=len_dim(i)) + end do + deallocate(dimids) ! pointer + + allocate(values(len_dim(1), len_dim(2), len_dim(3), len_dim(4), len_dim(5))) + if (all(len_dim /= 0)) call NF95_GET_VAR(ncid, varid, values) + + end subroutine nf95_gw_var_real_5d + + !************************************ + + subroutine nf95_gw_var_dble_1d(ncid, varid, values) + + ! Double precision, the array has rank 1. + + integer, intent(in):: ncid + integer, intent(in):: varid + double precision, pointer:: values(:) + + ! Variables local to the procedure: + integer nclen + integer, pointer:: dimids(:) + + !--------------------- + + call nf95_inquire_variable(ncid, varid, dimids=dimids) + + if (size(dimids) /= 1) then + print *, "nf95_gw_var_dble_1d:" + print *, "varid = ", varid + print *, "rank of NetCDF variable is ", size(dimids), ", not 1" + stop 1 + end if + + call nf95_inquire_dimension(ncid, dimids(1), nclen=nclen) + deallocate(dimids) ! pointer + + allocate(values(nclen)) + if (nclen /= 0) call NF95_GET_VAR(ncid, varid, values) + + end subroutine nf95_gw_var_dble_1d + + !************************************ + + subroutine nf95_gw_var_dble_2d(ncid, varid, values) + + ! Double precision, the array has rank 2. + + integer, intent(in):: ncid + integer, intent(in):: varid + double precision, pointer:: values(:, :) + + ! Variables local to the procedure: + integer nclen1, nclen2 + integer, pointer:: dimids(:) + + !--------------------- + + call nf95_inquire_variable(ncid, varid, dimids=dimids) + + if (size(dimids) /= 2) then + print *, "nf95_gw_var_dble_2d:" + print *, "varid = ", varid + print *, "rank of NetCDF variable is ", size(dimids), ", not 2" + stop 1 + end if + + call nf95_inquire_dimension(ncid, dimids(1), nclen=nclen1) + call nf95_inquire_dimension(ncid, dimids(2), nclen=nclen2) + deallocate(dimids) ! pointer + + allocate(values(nclen1, nclen2)) + if (nclen1 /= 0 .and. nclen2 /= 0) call NF95_GET_VAR(ncid, varid, values) + + end subroutine nf95_gw_var_dble_2d + + !************************************ + + subroutine nf95_gw_var_dble_3d(ncid, varid, values) + + ! Double precision, the array has rank 3. + + integer, intent(in):: ncid + integer, intent(in):: varid + double precision, pointer:: values(:, :, :) + + ! Variables local to the procedure: + integer nclen1, nclen2, nclen3 + integer, pointer:: dimids(:) + + !--------------------- + + call nf95_inquire_variable(ncid, varid, dimids=dimids) + + if (size(dimids) /= 3) then + print *, "nf95_gw_var_dble_3d:" + print *, "varid = ", varid + print *, "rank of NetCDF variable is ", size(dimids), ", not 3" + stop 1 + end if + + call nf95_inquire_dimension(ncid, dimids(1), nclen=nclen1) + call nf95_inquire_dimension(ncid, dimids(2), nclen=nclen2) + call nf95_inquire_dimension(ncid, dimids(3), nclen=nclen3) + deallocate(dimids) ! pointer + + allocate(values(nclen1, nclen2, nclen3)) + if (nclen1 * nclen2 * nclen3 /= 0) call NF95_GET_VAR(ncid, varid, values) + + end subroutine nf95_gw_var_dble_3d + + !************************************ + + subroutine nf95_gw_var_dble_4d(ncid, varid, values) + + ! Double precision, the array has rank 4. + + integer, intent(in):: ncid + integer, intent(in):: varid + double precision, pointer:: values(:, :, :, :) + + ! Variables local to the procedure: + integer len_dim(4), i + integer, pointer:: dimids(:) + + !--------------------- + + call nf95_inquire_variable(ncid, varid, dimids=dimids) + + if (size(dimids) /= 4) then + print *, "nf95_gw_var_dble_4d:" + print *, "varid = ", varid + print *, "rank of NetCDF variable is ", size(dimids), ", not 4" + stop 1 + end if + + do i = 1, 4 + call nf95_inquire_dimension(ncid, dimids(i), nclen=len_dim(i)) + end do + deallocate(dimids) ! pointer + + allocate(values(len_dim(1), len_dim(2), len_dim(3), len_dim(4))) + if (all(len_dim /= 0)) call NF95_GET_VAR(ncid, varid, values) + + end subroutine nf95_gw_var_dble_4d + + !************************************ + + subroutine nf95_gw_var_int_1d(ncid, varid, values) + + ! Integer type, the array has rank 1. + + integer, intent(in):: ncid + integer, intent(in):: varid + integer, pointer:: values(:) + + ! Variables local to the procedure: + integer nclen + integer, pointer:: dimids(:) + + !--------------------- + + call nf95_inquire_variable(ncid, varid, dimids=dimids) + + if (size(dimids) /= 1) then + print *, "nf95_gw_var_int_1d:" + print *, "varid = ", varid + print *, "rank of NetCDF variable is ", size(dimids), ", not 1" + stop 1 + end if + + call nf95_inquire_dimension(ncid, dimids(1), nclen=nclen) + deallocate(dimids) ! pointer + + allocate(values(nclen)) + if (nclen /= 0) call NF95_GET_VAR(ncid, varid, values) + + end subroutine nf95_gw_var_int_1d + + !************************************ + + subroutine nf95_gw_var_int_3d(ncid, varid, values) + + ! Integer type, the array has rank 3. + + integer, intent(in):: ncid + integer, intent(in):: varid + integer, pointer:: values(:, :, :) + + ! Variables local to the procedure: + integer nclen1, nclen2, nclen3 + integer, pointer:: dimids(:) + + !--------------------- + + call nf95_inquire_variable(ncid, varid, dimids=dimids) + + if (size(dimids) /= 3) then + print *, "nf95_gw_var_int_3d:" + print *, "varid = ", varid + print *, "rank of NetCDF variable is ", size(dimids), ", not 3" + stop 1 + end if + + call nf95_inquire_dimension(ncid, dimids(1), nclen=nclen1) + call nf95_inquire_dimension(ncid, dimids(2), nclen=nclen2) + call nf95_inquire_dimension(ncid, dimids(3), nclen=nclen3) + deallocate(dimids) ! pointer + + allocate(values(nclen1, nclen2, nclen3)) + if (nclen1 * nclen2 * nclen3 /= 0) call NF95_GET_VAR(ncid, varid, values) + + end subroutine nf95_gw_var_int_3d + +end module nf95_gw_var_m Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/nf95_inq_dimid.f90 =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/nf95_inq_dimid.f90 (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/nf95_inq_dimid.f90 (revision 1795) @@ -0,0 +1,31 @@ +module nf95_inq_dimid_m + + implicit none + +contains + + subroutine nf95_inq_dimid(ncid, name, dimid, ncerr) + + use handle_err_m, only: handle_err + use netcdf, only: nf90_inq_dimid + + integer, intent(in) :: ncid + character (len = *), intent(in) :: name + integer, intent(out) :: dimid + integer, intent(out), optional:: ncerr + + ! Variable local to the procedure: + integer ncerr_not_opt + + !------------------- + + ncerr_not_opt = nf90_inq_dimid(ncid, name, dimid) + if (present(ncerr)) then + ncerr = ncerr_not_opt + else + call handle_err("nf95_inq_dimid " // name, ncerr_not_opt, ncid) + end if + + end subroutine nf95_inq_dimid + +end module nf95_inq_dimid_m Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/nf95_inq_varid.f90 =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/nf95_inq_varid.f90 (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/nf95_inq_varid.f90 (revision 1795) @@ -0,0 +1,31 @@ +module nf95_inq_varid_m + + implicit none + +contains + + subroutine nf95_inq_varid(ncid, name, varid, ncerr) + + use handle_err_m, only: handle_err + use netcdf, only: nf90_inq_varid + + integer, intent(in) :: ncid + character(len=*), intent(in):: name + integer, intent(out) :: varid + integer, intent(out), optional:: ncerr + + ! Variable local to the procedure: + integer ncerr_not_opt + + !------------------- + + ncerr_not_opt = nf90_inq_varid(ncid, name, varid) + if (present(ncerr)) then + ncerr = ncerr_not_opt + else + call handle_err("nf95_inq_varid, name = " // name, ncerr_not_opt, ncid) + end if + + end subroutine nf95_inq_varid + +end module nf95_inq_varid_m Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/nf95_inquire.f90 =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/nf95_inquire.f90 (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/nf95_inquire.f90 (revision 1795) @@ -0,0 +1,34 @@ +module nf95_inquire_m + + implicit none + +contains + + + subroutine nf95_inquire(ncid, nDimensions, nVariables, nAttributes, & + unlimitedDimId, formatNum, ncerr) + + use handle_err_m, only: handle_err + use netcdf, only: nf90_inquire + + integer, intent( in) :: ncid + integer, optional, intent(out) :: nDimensions, nVariables, nAttributes + integer, optional, intent(out) :: unlimitedDimId, formatNum + integer, intent(out), optional:: ncerr + + ! Variable local to the procedure: + integer ncerr_not_opt + + !------------------- + + ncerr_not_opt = nf90_inquire(ncid, nDimensions, nVariables, nAttributes, & + unlimitedDimId, formatNum) + if (present(ncerr)) then + ncerr = ncerr_not_opt + else + call handle_err("nf95_inquire", ncerr_not_opt, ncid) + end if + + end subroutine nf95_inquire + +end module nf95_inquire_m Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/nf95_inquire_attribute.f90 =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/nf95_inquire_attribute.f90 (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/nf95_inquire_attribute.f90 (revision 1795) @@ -0,0 +1,35 @@ +module nf95_inquire_attribute_m + + implicit none + +contains + + + subroutine nf95_inquire_attribute(ncid, varid, name, xtype, nclen, attnum, & + ncerr) + + use handle_err_m, only: handle_err + use netcdf, only: nf90_inquire_attribute + + integer, intent( in) :: ncid, varid + character (len = *), intent( in) :: name + integer, intent(out), optional :: xtype, nclen, attnum + integer, intent(out), optional:: ncerr + + ! Variable local to the procedure: + integer ncerr_not_opt + + !------------------- + + ncerr_not_opt = nf90_inquire_attribute(ncid, varid, name, xtype, nclen, & + attnum) + if (present(ncerr)) then + ncerr = ncerr_not_opt + else + call handle_err("nf95_inquire_attribute " // name, ncerr_not_opt, & + ncid, varid) + end if + + end subroutine nf95_inquire_attribute + +end module nf95_inquire_attribute_m Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/nf95_inquire_dimension.f90 =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/nf95_inquire_dimension.f90 (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/nf95_inquire_dimension.f90 (revision 1795) @@ -0,0 +1,31 @@ +module nf95_inquire_dimension_m + + implicit none + +contains + + subroutine nf95_inquire_dimension(ncid, dimid, name, nclen, ncerr) + + use handle_err_m, only: handle_err + use netcdf, only: nf90_inquire_dimension + + integer, intent( in) :: ncid, dimid + character (len = *), optional, intent(out) :: name + integer, optional, intent(out) :: nclen + integer, intent(out), optional:: ncerr + + ! Variable local to the procedure: + integer ncerr_not_opt + + !------------------- + + ncerr_not_opt = nf90_inquire_dimension(ncid, dimid, name, nclen) + if (present(ncerr)) then + ncerr = ncerr_not_opt + else + call handle_err("nf95_inquire_dimension", ncerr_not_opt, ncid) + end if + + end subroutine nf95_inquire_dimension + +end module nf95_inquire_dimension_m Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/nf95_inquire_variable.f90 =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/nf95_inquire_variable.f90 (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/nf95_inquire_variable.f90 (revision 1795) @@ -0,0 +1,53 @@ +module nf95_inquire_variable_m + + implicit none + +contains + + subroutine nf95_inquire_variable(ncid, varid, name, xtype, ndims, dimids, & + nAtts, ncerr) + + ! In "nf90_inquire_variable", "dimids" is an assumed-size array. + ! This is not optimal. + ! We are in the classical case of an array the size of which is + ! unknown in the calling procedure, before the call. + ! Here we use a better solution: a pointer argument array. + ! This procedure associates and defines "dimids" if it is present. + + use handle_err_m, only: handle_err + use netcdf, only: nf90_inquire_variable, nf90_max_var_dims + + integer, intent(in):: ncid, varid + character(len = *), optional, intent(out):: name + integer, optional, intent(out) :: xtype, ndims + integer, dimension(:), optional, pointer :: dimids + integer, optional, intent(out) :: nAtts + integer, intent(out), optional :: ncerr + + ! Variable local to the procedure: + integer ncerr_not_opt + integer dimids_local(nf90_max_var_dims) + integer ndims_not_opt + + !------------------- + + if (present(dimids)) then + ncerr_not_opt = nf90_inquire_variable(ncid, varid, name, xtype, & + ndims_not_opt, dimids_local, nAtts) + allocate(dimids(ndims_not_opt)) ! also works if ndims_not_opt == 0 + dimids = dimids_local(:ndims_not_opt) + if (present(ndims)) ndims = ndims_not_opt + else + ncerr_not_opt = nf90_inquire_variable(ncid, varid, name, xtype, ndims, & + nAtts=nAtts) + end if + + if (present(ncerr)) then + ncerr = ncerr_not_opt + else + call handle_err("nf95_inquire_variable", ncerr_not_opt, ncid, varid) + end if + + end subroutine nf95_inquire_variable + +end module nf95_inquire_variable_m Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/nf95_open.f90 =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/nf95_open.f90 (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/nf95_open.f90 (revision 1795) @@ -0,0 +1,32 @@ +module nf95_open_m + + implicit none + +contains + + subroutine nf95_open(path, mode, ncid, chunksize, ncerr) + + use handle_err_m, only: handle_err + use netcdf, only: nf90_open + + character(len=*), intent(in):: path + integer, intent(in):: mode + integer, intent(out):: ncid + integer, intent(inout), optional:: chunksize + integer, intent(out), optional:: ncerr + + ! Variable local to the procedure: + integer ncerr_not_opt + + !------------------- + + ncerr_not_opt = nf90_open(path, mode, ncid, chunksize) + if (present(ncerr)) then + ncerr = ncerr_not_opt + else + call handle_err("nf95_open " // path, ncerr_not_opt) + end if + + end subroutine nf95_open + +end module nf95_open_m Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/nf95_put_att.f90 =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/nf95_put_att.f90 (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/nf95_put_att.f90 (revision 1795) @@ -0,0 +1,96 @@ +module nf95_put_att_m + + implicit none + + interface nf95_put_att + module procedure nf95_put_att_text, nf95_put_att_one_FourByteInt, & + nf95_put_att_one_FourByteReal + end interface + + private + public nf95_put_att + +contains + + subroutine nf95_put_att_text(ncid, varid, name, values, ncerr) + + use netcdf, only: nf90_put_att + use handle_err_m, only: handle_err + + integer, intent(in) :: ncid, varid + character(len = *), intent(in) :: name + character(len = *), intent(in) :: values + integer, intent(out), optional:: ncerr + + ! Variable local to the procedure: + integer ncerr_not_opt + + !------------------- + + ncerr_not_opt = nf90_put_att(ncid, varid, name, values) + if (present(ncerr)) then + ncerr = ncerr_not_opt + else + call handle_err("nf95_put_att_text " // trim(name), ncerr_not_opt, & + ncid, varid) + end if + + end subroutine nf95_put_att_text + + !************************************ + + subroutine nf95_put_att_one_FourByteInt(ncid, varid, name, values, ncerr) + + use netcdf, only: nf90_put_att + use handle_err_m, only: handle_err + use typesizes, only: FourByteInt + + integer, intent(in) :: ncid, varid + character(len = *), intent(in) :: name + integer(kind = FourByteInt), intent(in) :: values + integer, intent(out), optional:: ncerr + + ! Variable local to the procedure: + integer ncerr_not_opt + + !------------------- + + ncerr_not_opt = nf90_put_att(ncid, varid, name, values) + if (present(ncerr)) then + ncerr = ncerr_not_opt + else + call handle_err("nf95_put_att_one_FourByteInt " // trim(name), & + ncerr_not_opt, ncid, varid) + end if + + end subroutine nf95_put_att_one_FourByteInt + + !************************************ + + subroutine nf95_put_att_one_FourByteReal(ncid, varid, name, values, ncerr) + + use netcdf, only: nf90_put_att + use handle_err_m, only: handle_err + use typesizes, only: FourByteReal + + integer, intent(in) :: ncid, varid + character(len = *), intent(in) :: name + real(kind = FourByteReal), intent(in) :: values + integer, intent(out), optional:: ncerr + + ! Variable local to the procedure: + integer ncerr_not_opt + + !------------------- + + ncerr_not_opt = nf90_put_att(ncid, varid, name, values) + if (present(ncerr)) then + ncerr = ncerr_not_opt + else + call handle_err("nf95_put_att_one_FourByteReal " // trim(name), & + ncerr_not_opt, ncid, varid) + end if + + end subroutine nf95_put_att_one_FourByteReal + +end module nf95_put_att_m Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/nf95_put_var.f90 =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/nf95_put_var.f90 (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/nf95_put_var.f90 (revision 1795) @@ -0,0 +1,343 @@ +module nf95_put_var_m + + implicit none + + interface nf95_put_var + module procedure nf95_put_var_FourByteReal, nf95_put_var_FourByteInt, & + nf95_put_var_1D_FourByteReal, nf95_put_var_1D_FourByteInt, & + nf95_put_var_1D_EightByteReal, nf95_put_var_2D_FourByteReal, & + nf95_put_var_2D_EightByteReal, nf95_put_var_3D_FourByteReal, & + nf95_put_var_3D_EightByteReal, nf95_put_var_4D_FourByteReal, & + nf95_put_var_4D_EightByteReal + end interface + + private + public nf95_put_var + +contains + + subroutine nf95_put_var_FourByteReal(ncid, varid, values, start, ncerr) + + use typesizes, only: FourByteReal + use netcdf, only: nf90_put_var + use handle_err_m, only: handle_err + + integer, intent(in) :: ncid, varid + real(kind = FourByteReal), intent(in) :: values + integer, dimension(:), optional, intent(in) :: start + integer, intent(out), optional:: ncerr + + ! Variable local to the procedure: + integer ncerr_not_opt + + !------------------- + + ncerr_not_opt = nf90_put_var(ncid, varid, values, start) + if (present(ncerr)) then + ncerr = ncerr_not_opt + else + call handle_err("nf95_put_var_FourByteReal", ncerr_not_opt, ncid, & + varid) + end if + + end subroutine nf95_put_var_FourByteReal + + !*********************** + + subroutine nf95_put_var_FourByteInt(ncid, varid, values, start, ncerr) + + use typesizes, only: FourByteInt + use netcdf, only: nf90_put_var + use handle_err_m, only: handle_err + + integer, intent(in) :: ncid, varid + integer(kind = FourByteInt), intent(in) :: values + integer, dimension(:), optional, intent(in) :: start + integer, intent(out), optional:: ncerr + + ! Variable local to the procedure: + integer ncerr_not_opt + + !------------------- + + ncerr_not_opt = nf90_put_var(ncid, varid, values, start) + if (present(ncerr)) then + ncerr = ncerr_not_opt + else + call handle_err("nf95_put_var_FourByteInt", ncerr_not_opt, ncid, & + varid) + end if + + end subroutine nf95_put_var_FourByteInt + + !*********************** + + subroutine nf95_put_var_1D_FourByteReal(ncid, varid, values, start, & + count_nc, stride, map, ncerr) + + use typesizes, only: FourByteReal + use netcdf, only: nf90_put_var + use handle_err_m, only: handle_err + + integer, intent(in) :: ncid, varid + real(kind = FourByteReal), intent(in) :: values(:) + integer, dimension(:), optional, intent(in) :: start, count_nc, stride, map + integer, intent(out), optional:: ncerr + + ! Variable local to the procedure: + integer ncerr_not_opt + + !------------------- + + ncerr_not_opt = nf90_put_var(ncid, varid, values, start, count_nc, & + stride, map) + if (present(ncerr)) then + ncerr = ncerr_not_opt + else + call handle_err("nf95_put_var_1D_FourByteReal", ncerr_not_opt, ncid, & + varid) + end if + + end subroutine nf95_put_var_1D_FourByteReal + + !*********************** + + subroutine nf95_put_var_1D_FourByteInt(ncid, varid, values, start, & + count_nc, stride, map, ncerr) + + use typesizes, only: FourByteInt + use netcdf, only: nf90_put_var + use handle_err_m, only: handle_err + + integer, intent(in) :: ncid, varid + integer(kind = FourByteInt), intent(in) :: values(:) + integer, dimension(:), optional, intent(in) :: start, count_nc, stride, map + integer, intent(out), optional:: ncerr + + ! Variable local to the procedure: + integer ncerr_not_opt + + !------------------- + + ncerr_not_opt = nf90_put_var(ncid, varid, values, start, count_nc, & + stride, map) + if (present(ncerr)) then + ncerr = ncerr_not_opt + else + call handle_err("nf95_put_var_1D_FourByteInt", ncerr_not_opt, ncid, & + varid) + end if + + end subroutine nf95_put_var_1D_FourByteInt + + !*********************** + + subroutine nf95_put_var_1D_EightByteReal(ncid, varid, values, start, & + count_nc, stride, map, ncerr) + + use typesizes, only: eightByteReal + use netcdf, only: nf90_put_var + use handle_err_m, only: handle_err + + integer, intent(in) :: ncid, varid + real (kind = EightByteReal), intent(in) :: values(:) + integer, dimension(:), optional, intent(in):: start, count_nc, stride, map + integer, intent(out), optional:: ncerr + + ! Variable local to the procedure: + integer ncerr_not_opt + + !------------------- + + ncerr_not_opt = nf90_put_var(ncid, varid, values, start, count_nc, & + stride, map) + if (present(ncerr)) then + ncerr = ncerr_not_opt + else + call handle_err("nf95_put_var_1D_eightByteReal", ncerr_not_opt, ncid, & + varid) + end if + + end subroutine nf95_put_var_1D_EightByteReal + + !*********************** + + subroutine nf95_put_var_2D_FourByteReal(ncid, varid, values, start, & + count_nc, stride, map, ncerr) + + use typesizes, only: FourByteReal + use netcdf, only: nf90_put_var + use handle_err_m, only: handle_err + + integer, intent(in) :: ncid, varid + real (kind = FourByteReal), intent(in) :: values(:, :) + integer, dimension(:), optional, intent(in) :: start, count_nc, stride, map + integer, intent(out), optional:: ncerr + + ! Variable local to the procedure: + integer ncerr_not_opt + + !------------------- + + ncerr_not_opt = nf90_put_var(ncid, varid, values, start, count_nc, & + stride, map) + if (present(ncerr)) then + ncerr = ncerr_not_opt + else + call handle_err("nf95_put_var_2D_FourByteReal", ncerr_not_opt, ncid, & + varid) + end if + + end subroutine nf95_put_var_2D_FourByteReal + + !*********************** + + subroutine nf95_put_var_2D_EightByteReal(ncid, varid, values, start, & + count_nc, stride, map, ncerr) + + use typesizes, only: EightByteReal + use netcdf, only: nf90_put_var + use handle_err_m, only: handle_err + + integer, intent(in) :: ncid, varid + real (kind = EightByteReal), intent(in) :: values(:, :) + integer, dimension(:), optional, intent(in) :: start, count_nc, stride, map + integer, intent(out), optional:: ncerr + + ! Variable local to the procedure: + integer ncerr_not_opt + + !------------------- + + ncerr_not_opt = nf90_put_var(ncid, varid, values, start, count_nc, & + stride, map) + if (present(ncerr)) then + ncerr = ncerr_not_opt + else + call handle_err("nf95_put_var_2D_EightByteReal", ncerr_not_opt, ncid, & + varid) + end if + + end subroutine nf95_put_var_2D_EightByteReal + + !*********************** + + subroutine nf95_put_var_3D_FourByteReal(ncid, varid, values, start, & + count_nc, stride, map, ncerr) + + use typesizes, only: FourByteReal + use netcdf, only: nf90_put_var + use handle_err_m, only: handle_err + + integer, intent(in) :: ncid, varid + real (kind = FourByteReal), intent(in) :: values(:, :, :) + integer, dimension(:), optional, intent(in) :: start, count_nc, stride, map + integer, intent(out), optional:: ncerr + + ! Variable local to the procedure: + integer ncerr_not_opt + + !------------------- + + ncerr_not_opt = nf90_put_var(ncid, varid, values, start, count_nc, & + stride, map) + if (present(ncerr)) then + ncerr = ncerr_not_opt + else + call handle_err("nf95_put_var_3D_FourByteReal", ncerr_not_opt, ncid, & + varid) + end if + + end subroutine nf95_put_var_3D_FourByteReal + + !*********************** + + subroutine nf95_put_var_3D_EightByteReal(ncid, varid, values, start, & + count_nc, stride, map, ncerr) + + use typesizes, only: eightByteReal + use netcdf, only: nf90_put_var + use handle_err_m, only: handle_err + + integer, intent(in) :: ncid, varid + real (kind = EightByteReal), intent(in) :: values(:, :, :) + integer, dimension(:), optional, intent(in) :: start, count_nc, stride, map + integer, intent(out), optional:: ncerr + + ! Variable local to the procedure: + integer ncerr_not_opt + + !------------------- + + ncerr_not_opt = nf90_put_var(ncid, varid, values, start, count_nc, & + stride, map) + if (present(ncerr)) then + ncerr = ncerr_not_opt + else + call handle_err("nf95_put_var_3D_eightByteReal", ncerr_not_opt, ncid, & + varid) + end if + + end subroutine nf95_put_var_3D_EightByteReal + + !*********************** + + subroutine nf95_put_var_4D_FourByteReal(ncid, varid, values, start, & + count_nc, stride, map, ncerr) + + use typesizes, only: FourByteReal + use netcdf, only: nf90_put_var + use handle_err_m, only: handle_err + + integer, intent(in) :: ncid, varid + real (kind = FourByteReal), intent(in) :: values(:, :, :, :) + integer, dimension(:), optional, intent(in) :: start, count_nc, stride, map + integer, intent(out), optional:: ncerr + + ! Variable local to the procedure: + integer ncerr_not_opt + + !------------------- + + ncerr_not_opt = nf90_put_var(ncid, varid, values, start, count_nc, & + stride, map) + if (present(ncerr)) then + ncerr = ncerr_not_opt + else + call handle_err("nf95_put_var_4D_FourByteReal", ncerr_not_opt, ncid, & + varid) + end if + + end subroutine nf95_put_var_4D_FourByteReal + + !*********************** + + subroutine nf95_put_var_4D_EightByteReal(ncid, varid, values, start, & + count_nc, stride, map, ncerr) + + use typesizes, only: EightByteReal + use netcdf, only: nf90_put_var + use handle_err_m, only: handle_err + + integer, intent(in):: ncid, varid + real(kind = EightByteReal), intent(in):: values(:, :, :, :) + integer, dimension(:), optional, intent(in):: start, count_nc, stride, map + integer, intent(out), optional:: ncerr + + ! Variable local to the procedure: + integer ncerr_not_opt + + !------------------- + + ncerr_not_opt = nf90_put_var(ncid, varid, values, start, count_nc, & + stride, map) + if (present(ncerr)) then + ncerr = ncerr_not_opt + else + call handle_err("nf95_put_var_4D_EightByteReal", ncerr_not_opt, ncid, & + varid) + end if + + end subroutine nf95_put_var_4D_EightByteReal + +end module nf95_put_var_m Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/nf95_redef.f90 =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/nf95_redef.f90 (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/NetCDF95/nf95_redef.f90 (revision 1795) @@ -0,0 +1,29 @@ +module nf95_redef_m + + implicit none + +contains + + subroutine nf95_redef(ncid, ncerr) + + use handle_err_m, only: handle_err + use netcdf, only: nf90_redef + + integer, intent( in) :: ncid + integer, intent(out), optional :: ncerr + + ! Variable local to the procedure: + integer ncerr_not_opt + + !------------------- + + ncerr_not_opt = nf90_redef(ncid) + if (present(ncerr)) then + ncerr = ncerr_not_opt + else + call handle_err("nf95_redef", ncerr_not_opt, ncid) + end if + + end subroutine nf95_redef + +end module nf95_redef_m Index: LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/Read_me =================================================================== --- LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/Read_me (revision 1795) +++ LMDZ5/branches/testing/tools/Max_diff_nc_with_lib/Read_me (revision 1795) @@ -0,0 +1,2 @@ +See: +www.lmd.jussieu.fr/~lglmd/Max_diff_nc Index: LMDZ5/branches/testing/tools/install_1d_src.sh =================================================================== --- LMDZ5/branches/testing/tools/install_1d_src.sh (revision 1750) +++ LMDZ5/branches/testing/tools/install_1d_src.sh (revision 1795) @@ -54,5 +54,5 @@ ! Declarations specifiques pour le 1D. A reprendre \\ INCLUDE \"flux_arp.h\" - }; /^.*REAL *::.*fsens *, *flat/d; /^.*LOGICAL *::.*ok_flux_surf/d; /^.*COMMON.*flux_arp.*/d " pbl_surface_mod.F90 + }; /^.*REAL *::.*fsens *, *flat/d; /^.*LOGICAL *::.*ok_flux_surf/d; /^.*COMMON.*flux_arp.*/d " pbl_surface_mod.F90 ln -s ../dyn3d/mod_const_mpi.F90 . @@ -62,2 +62,3 @@ ln -s ../dyn3d/control_mod.F90 . ln -s ../dyn3d/q_sat.F . + ln -s ../dyn3d/disvert.F90 .